← Back to team overview

dolfin team mailing list archive

Re: Stokes in DOLFIN

 

This looks very nice Anders! It seems that we now in DOLFIN have nice
working modules (Stokes, Navier-Stokes, Elasticity,...) tackling many of
the usual problems when implementing a PDE solver: high order, mixed
element, non-linear problems (as a linearized problem),...

I think the approach with a PDE class holding the pair (a,L) sounds like 
good idea.

For the linear algebra we have Petsc + Hypre giving access to a range of
parallel Krylov based solvers plus algebraic multigrid.

For efficiency the next big challenge is to have all components in DOLFIN
parallelized; that means the mesh and the assembly routine. It seems that 
Matt's Sieve will be the solution.

As for abstraction of the implementation, hiding the linear algebra as you
say is one thing. Another is to automate space-time FEM in FFC. That is,
to give the forms for a space-time domain, instead of just for a space
domain, as is the case today. The functionality of space-time FEM is
already there, but the time discretization today has to be given
explicitly to FFC. This would take the form-file code closer to the PDE
formulation for time dependent problems.

Another thing to think about is to also automate the discretization of non
linear problems. That is to give FFC the non-linear form to either
discretize for the non-linear solvers of Petsc, or to linearize and
dicretize for a linear solver.

/Johan


> Looks like it's working now. The following simple code is the complete
> DOLFIN implementation of the Stokes solver (which looks almost
> identical to the Poisson solver):
>
>     Stokes::BilinearForm a;
>     Stokes::LinearForm L(f);
>
>     Matrix A;
>     Vector x, b;
>     FEM::assemble(a, L, A, b, mesh, bc);
>
>     GMRES solver;
>     solver.solve(A, x, b);
>
>     Function u(x, mesh, a.trial());
>     File file("stokes.m");
>     file << u;
>
> and the following is the FFC implementation of the bilinear form:
>
>     a = (dot(grad(u), grad(v)) - p*div(v) + div(u)*q)*dx
>
> It can't get much simpler. We will probably add an additional set of
> tools to abstract away the linear algebra whenever that is possible,
> so that one can work entirely with Function and never have to look at
> vectors and matrices. The current tools will remain, however, in
> particular the assembly functions in FEM and of course Vector and
> Matrix, but there will be tools that build on these.
>
> One suggestion is to create a class PDE which holds a pair (a, L) of
> bilinear and linear forms. FFC would automatically generate a class
> Stokes that inherits from PDE and just initializes the two forms. Thus,
> instead of accessing the forms directly like Stokes::BilinearForm
> (which will still be possible), one can just do
>
>     Function u;
>     Stokes pde(f, mesh, bc);
>     pde.solve(u);
>     File file("stokes.m");
>     file << u;
>
> Choice of linear solvers etc can be controlled with parameters.
>
> We also need to add some functionality to Function for representing
> slices of vector-valued functions (like velocity and pressure in the
> example above). This is easy (just need to keep an extra integer
> component_offset) but we need to decide on a proper syntax for
> creating slices.
>
> Attached is a plot of the Stokes solution with Taylor-Hood elements.
>
> /Anders
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
>


-- 
Johan Hoffman, PhD
Assistant Professor
School of Computer Science and Communication
Royal Institute of Technology KTH
SE-100 44 Stockholm
Sweden

Email: jhoffman@xxxxxxxxxxx
URL: www.nada.kth.se/~jhoffman




Follow ups

References