dolfin team mailing list archive

[Jed Brown <jed@xxxxxxxx>]

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

Attachment: pgpg1vgX_6Ise.pgp
Description: PGP signature