dolfin team mailing list archive

[Kent Andre <kent-and@xxxxxxxxx>]

On ti., 2009-04-07 at 18:09 +0200, Jed Brown wrote:
> On Tue 2009-04-07 14:19, Anders Logg wrote:
> > On Tue, Apr 07, 2009 at 01:51:54PM +0200, Jed Brown wrote:
> > > On Mon 2009-04-06 07:42, Matthew Knepley wrote:
> > > > On Mon, Apr 6, 2009 at 1:09 AM, Anders Logg <logg@xxxxxxxxx> wrote:
> > > > 
> > > >     On Sun, Apr 05, 2009 at 05:11:51PM -0500, Matthew Knepley wrote:
> > > >     > On Sun, Apr 5, 2009 at 2:54 PM, Anders Logg <logg@xxxxxxxxx> wrote:
> > > >     >
> > > >     >     Dirichlet BC would need to be added in/after each multiplication and
> > > >     >     it should be possible to build it into the mult() operator and make
> > > >     it
> > > >     >     efficient.
> > > >     >
> > > >     >
> > > >     > I still think the best way to handle this is to eliminate them, as I
> > > >     talked
> > > >     > about
> > > >     > last time at Simula.
> > > >     >
> > > >     >   Matt
> > > > 
> > > >     I think we're in position to do that now with the new FunctionSpace
> > > >     class. One could allow restrictions to be imposed, for example
> > > > 
> > > >      V.restrict(bc);
> > > > 
> > > >     But I haven't seen any good numbers to indicate it's worth the effort.
> > > >     How much better is it to eliminate, condition numbers etc for letting
> > > >     the constrained dofs just sit along?
> > > > 
> > > > 
> > > > Its not about performance, so much as
> > > > 
> > > >   1) Programming ease, ESPECIALLY with nonlinear problems
> > > > 
> > > >   2) Separation of storage organization from traversal organization
> > > > 
> > > >   3) Cache performance
> > > > 
> > > >   4) Interaction with external packages
> > > > 
> > > > I think I must have done a terrible job explaining this. I have always coded
> > > > both ways and always found that elimination is better.
> > > 
> > > I agree with all your points, but with caveats on 3.
> > > 
> > > There are much bigger gains in cache performance by interlacing the dofs
> > > so as to enable the use of inodes and blocked matrix formats (BAIJ).
> > > 
> > > If we want to use blocked formats with slip boundary conditions, we need
> > > to be able to leave the condition in.  A slip boundary condition imposes
> > > a Dirichlet condition on the normal component and stress conditions on
> > > the tangent components.  The global degrees of freedom can be written in
> > > a rotated coordinate frame so that the Dirichlet condition can be
> > > completely removed, but when using BAIJ we have to put something there.
> > > This can be done, complete with "zeroing the column", but it requires
> > > manipulations in function evaluation and when setting values in the
> > > matrix.  In some of my experiments, the performance benefits of BAIJ
> > > over AIJ+inodes justify this strategy.
> > > 
> > > If you have Dirichlet conditions on all components, then removing these
> > > dofs from the system really is better.
> > 
> > How much faster?
> 
> Interlacing and blocking often mean a factor of 2 or more, depending on
> the number of degrees of freedom per node and the preconditioner.  Note
> that relaxation preconditioners (such as a SOR smoother within
> multigrid) can be stronger when fields are interlaced because it relaxes
> the whole block rather than each component individually.
> 
> For example, see Table 1 of
> 
> @inproceedings{gropp2000pmt,
>   author = {Gropp, William D. and Kaushik, Dinesh K. and Keyes, David E. and Smith, Barry},
>   title = {Performance modeling and tuning of an unstructured mesh CFD application},
>   booktitle = {Supercomputing '00: Proceedings of the 2000 ACM/IEEE conference on Supercomputing (CDROM)},
>   year = {2000},
>   isbn = {0-7803-9802-5},
>   pages = {34},
>   location = {Dallas, Texas, United States},
>   publisher = {IEEE Computer Society},
>   address = {Washington, DC, USA},
> }
> 
> 
> Removing the Dirichlet conditions completely won't have much impact on
> the solver (in terms of number of iterations, the cost per iteration and
> memory usage is affected but not by huge factors) provided you are using
> a "normal" preconditioner.  With "physics-based" preconditioners,
> keeping Dirichlet dofs around means you have to think about how their
> presence manifests itself in the preconditioner.
> 
> > >  Also, if you are not using blocked matrix formats, you might as
> > > well remove all Dirichlet dofs for all the reasons Matt states.
> > 
> > It still doesn't make sense to me.
> > 
> > > >   1) Programming ease, ESPECIALLY with nonlinear problems
> > 
> > I think keeping the dofs is easier.
> 
> You still have to store the Dirichlet values somewhere and you have to
> manipulate the current iterate when you start function evaluation and
> Jacobian assembly, unless you are satisfied with just zeroing rows
> (which ruins symmetry and often slows convergence).  In addition, you
> have to manipulate the element matrices prior to insertion (manipulating
> an assembled matrix is a terrible thing to want to do) and set the
> diagonal entries for Dirichlet conditions afterward.  It is more steps
> to leave them in.
> 
> 
> Jed

Removing the bc will give a small performance boost. It is unlikely 
to be large since it is a matter of removing a small diagonal matrix
from the much larger matrix. 

As it is in Dolin now, the solution lives on a function space
that lives on the whole mesh. The class FunctionSpace can be used
to represent both unknowns and the Dirichlet conditions, when 
eliminating the Dirichlet dofs from the linear system. Nobody has coded
up this yet, but the abstractions are in place and therefore this
elimination should not affect the dolfin demos or the end-user.  
Hence, anyone that wants to implement support for this may do it and the
patch should not be difficult to include in Dolfin. 

Kent