dolfin team mailing list archive

["Johan Hoffman" <jhoffman@xxxxxxxxxx>]

> On Mon 2008-05-19 08:41, Johan Hoffman wrote:
>> > 2008/5/19 Johan Hoffman <jhoffman@xxxxxxxxxx>:
>> >>> On Sun 2008-05-18 22:55, Johan Hoffman wrote:
>> >>>> > On Sun 2008-05-18 21:54, Johan Hoffman wrote:
>> >>>> >> > On Sat, May 17, 2008 at 04:40:48PM +0200, Johan Hoffman wrote:
>> >>>> >> >
>> >>>> >> > 1. Solve time may dominate assemble anyway so that's where we
>> >>>> should
>> >>>> >> > optimize.
>> >>>> >>
>> >>>> >> Yes, there may be such cases, in particular for simple forms
>> >>>> (Laplace
>> >>>> >> equation etc.). For more complex forms with more terms and
>> >>>> coefficients,
>> >>>> >> assembly typically dominates, from what I have seen. This is the
>> >>>> case
>> >>>> >> for
>> >>>> >> the flow problems of Murtazo for example.
>> >>>> >
>> >>>> > This probably depends if you use are using a projection method.
>> If
>> >>>> you
>> >>>> > are
>> >>>> > solving the saddle point problem, you can forget about assembly
>> >>>> time.
>> >>>>
>> >>>> Well, this is not what we see. I agree that this is what you would
>> >>>> like,
>> >>>> but this is not the case now. That is why we are now focusing on
>> the
>> >>>> assembly bottleneck.
>> >>>
>> >>> This just occurred to me.  If you have a Newtonian fluid, then the
>> >>> momentum
>> >>> equations are block diagonal, but this is not reflected in the
>> matrix
>> >>> structure.
>> >>> Sure enough, run the stokes demo with -mat_view_draw -draw_pause -1
>> and
>> >>> note
>> >>> that the off-diagonal blocks of the momentum equations are cyan
>> which
>> >>> means they
>> >>> are set, but have value zero.  This almost doubles the number of
>> >>> insertions into
>> >>> the global matrix.
>> >>
>> >> Good. You are right, this piece of information is not used.
>> >>
>> >> I guess the most general thing is to have ffc delete zero matrix
>> entries
>> >> in computing the sparsity pattern. I do not think this is done today?
>> >
>> > We could add it with an appropriate
>> > ufc::dof_map::tabulate_sparsity_pattern(...),
>> > since the form compiler can figure out which entries are always zero.
>> > Currently, this information is hidden from dolfin, and therefore it
>> must
>> > simply use the local blocks.
>>
>> Yes, this would be very useful to add.
>>
>> > But A.add(...) calls a single block-addition function in PETSc
>> > or Trilinos, does anyone know how these will perform if the values
>> > contain zeros that are outside the initialized sparsity pattern?
>>
>> Worst case scenario is that a reallocation is triggered, but maybe it is
>> dealt with in a less drastic way?
>
> It is important to preallocate for all the entries you will be inserting.
> If
> the option MAT_IGNORE_ZERO_ENTRIES is set then zero entries will not be
> inserted.  I don't think you want to rely on this as a general mechanism
> for
> eliminating non-existant coupling since it would also eliminate entries
> that
> just happened to be missing for the current Jacobian.
>

Ok.

> Note that MatSetValues() needs to be called with a ``full'' block so if we
> eliminate the coupling zeros, we must make separate calls for the rows of
> each
> vector component.  However, this should not make much performance
> difference
> since the only overhead from splitting the rows is the function call.  Of
> course, if we use a hashed cache then this is definitely a non-issue.

Ok.

>> >>  Of course, if you really care about speed, you could
>> >>> implement the momentum equations matrix-free with a custom
>> >>> preconditioner
>> >>> that
>> >>> only uses one block (since they are the same).
>> >>>
>> >>> Also, it seems like a strange design decision for general coupled
>> >>> systems
>> >>> for
>> >>> the blocks to be separate like this.  (Actually, ordering the
>> unknowns
>> >>> this way
>> >>> makes some factorizations ``work'' even when the matrix is
>> indefinite.)
>> >>
>> >>  Ok. We should look in to this.
>> >
>> > The reasoning behind this decision was that it makes access to
>> subvectors
>> > for components of general hierarchial mixed elements easy. Vector
>> elements
>> > and mixed elements are treated as basically the same thing.
>> >
>> > Thus, renumbering should preferably happen in dolfin, which will have
>> to
>> > be
>> > done for many occasions anyway.
>>
>> Yes, it seems natural to have the option to renumber freely in dolfin,
>> in
>> particular to switch between "vertex-oriented" and "vector-component
>> oriented" numbering.
>
> I think this would be good.  Note that the non-square coupling matrices
> between
> mixed elements never need to be formed (for most saddle point
> preconditioners)
> so it would be natural to write them in separate form files.
>
>
>> Iterating once is similar to a projection method, but we iterate until
>> convergence in the outer non linear residual. Of course, after some time
>> the number of iterations is typically 1, so the method collapses into a
>> projection type method.
>
> Okay, cool.  Thanks for clarifying.  How fast is the convergence the first
> time
> for the nonlinear problem?  My (limited) experience and impression from
> the
> literature is that you only see quadratic convergence from the Newton
> method if
> you solve the linear saddle point problem instead of replacing it with the
> projection preconditioner.  Of course, if it only takes one iteration
> after some
> time, then $\Delta t \eta$ must be small in which case we expect the
> projection
> method to be good (this is not the case for the slow non-Newtonian fluids
> I work
> with).

For the fixpoint the convergence depends on the time step size (setting
the corresponding contraction mapping). For a Newton iteration I am
unsure, since we haven't used it for the pure flow problems, but only for
fluid-structure interaction.

Our cases are typically high Re flow, so the viscosity is small, and the
dissipative mechanism comes more or less entirely from the numerical
method. After just a few time steps the number of iterations is usually 1.

/Johan

> Jed
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/mailman/listinfo/dolfin-dev
>