← Back to team overview

dolfin team mailing list archive

  1. dolfin team
  2. Mailing list archive
  3. Message #19995

Re: Interface for variational problems with goal-oriented error control/adaptivity:

  Thread Previous • Date Previous • • Thread Next 
Johan Hake skrev 2010-10-21 00.28:

On Wednesday October 20 2010 14:27:55 Marie Rognes wrote:

On 20. okt. 2010 23:14, Johan Hake wrote:

On Wednesday October 20 2010 13:24:17 Marie Rognes wrote:

On 20. okt. 2010 21:51, Johan Hake wrote:

Progress has been made towards the "real" implementation (aka in the
c++ core) of automated goal-oriented error control and adaptivity (cf
email thread "Fwd: [Branch ~dolfin-core/dolfin/main] Rev 4734:
Introduce functionality for automated error control and adaptivity"
for background info.)

Summary:
----------------

Marie suggests the following interface:

Same as for standard VariationalProblems, but with

      pde.solve(u, tol, M)

where "tol" should be a (positive) real number and M a Form,
representing the quantity of interest/goal.

Main question:
-----------------------

Could the rest of you live happily with the above suggestion?

Longer explanation:
-------------------------------

For solving a variational problem adaptively such that the error in a
certain goal is controlled, the following algorithm is more or less

standard:
   # Start with initial mesh T_0

   for i in (0, max_iterations):
     (1)  find solution to variational problem on mesh T_i
     (2)  estimate global error of solution

     (2b) stop if error is below tolerance

     (3)  estimate local errors (indicators)
     (4)  refine mesh T_i ->  T_i+1 based on indicators

Steps (2) and (3) depend on the variational problem, and requires the
solving of additional variational problems (which can be derived from
the original variational problem), and the evaluation of specific
additional forms. In other words, these steps require a certain amount
of expert knowledge. However, we can automate this by generating the
required code.

Here is an example of how this could work:

Ex1.ufl:
   V = FiniteElement("CG", "triangle", 1)
   f = Coefficient(V)
   u = TrialFunction(V)
   v = TestFunction(V)
   a = dot(grad(u), grad(v))*dx
   L = f*v*dx
   M = u*dx

Compile with:
   ffc -e -l dolfin Ex1.ufl

(where -e == --with-error-control))

FFC will then generate "DOLFIN code" for the forms and finite elements
involved in a, L and M, additional forms/finite element spaces, and a
(sub-)class

     Ex1::ErrorControl : dolfin::ErrorControl

where dolfin::ErrorControl supplies the methods

     estimate_error(u)
     compute_error_indicators(u)

The functionality of this class will take care of steps 2 and 3.

The main file can then look as follows:

main1.cpp:
  #include "Ex1.h"

   ...

   Ex1::FunctionSpace V(mesh)
   Ex1::BilinearForm a(V, V)
   Ex1::LinearForm L(V);

   VariationalProblem(a, L, ...);

   Function u(V);
   Ex1::GoalFunctional M(V);
   pde.solve(u, tol, M);

Marie thinks this looks pretty clean and knows that it is doable
(there is a (more or less) working prototype in the launchpad branches
dolfin-error-control/ffc-error-control). (The error control
class/object can attached to the GoalFunctional and hence available
for the VariationalProblem.)

Remark: If we want to increase explicitness/flexibility we could do

something like instead:
   Ex1::ErrorControl ec(...);
   pde.solve(u, tol, M, ec);


What kind of guy would ErrorControl be?


The good guy?  ;)


Obviously!


ErrorControl (or maybe ErrorController or ErrorEstimator) should take
care of estimating the error of an approximate solution u (to be used as
a stopping criterion for mesh refinement) and computing error indicators
(to be used for actual mesh refinement).

Take a look at

     dolfin/adaptivity/ErrorControl.h

from

     bzr branch lp:~dolfin-core/dolfin/error-control

With the current plan, each generated sub-class will define the
required, problem-specific forms.


Which a user need to do if she want to implement her own scheme?


I think yes is the answer here (but I'm not sure if I fully understand
what you ask).

If a user want to use the brilliant setup devised by us, FFC/DOLFIN will
take care of everything. If not, the user is on his own, more or less.


:)


Is this guy
generated for you by UFL and FFC, or is the one in

   dolfin/adaptivity/ErrorControl .h

always used in your prefered example.


FFC will generate (and the ffc/error-control branch does) a sub-class of
this guy adapted to a given VariationalProblem F (or a/L) and a given
goal M. (Same design as for Forms and FunctionSpaces.)


Ok, but in your prefered example you just instantiate the GoalFunctional.
Where is the subclassed ErrorControl guy, or is the GoalFunctional him in
disguise?


On the other side, if someone wants to supply their own error control,
they probably want to control the entire adaptive algorithms as
well.


Can you give a schematical example of how a user would do that for both
of these suggestions?


With the Marie-preferred suggestion, the code is as in main1.cpp.

With the other (non-preferred suggestion), I say that if the user wants
to do adaptivity in their own way, I suggest the user does adaptivity in
their own way, namely, implements their own adaptive loop, their own
error estimates, marking strategy etc. I would claim that the error
estimation stage is the more difficult, so if the user knows how to do
that, the remainder should be pretty easy. In essence, these are the
steps in

      dolfin/adaptivity/AdaptiveSolver::solve(...)

(in the above mentioned branch)


Agree! It just took some time to grasp the concepts, I am probably not
there yet ;)


So, I suggest not. Other well-founded opinions are of course
welcome ;)

For nonlinear problems, however, a certain amount of differentiations
of forms will be required. We therefore need to know what the unknown
is, so that we can take the derivative of the forms with respect to
it. Anders and Marie have discussed various alternatives quite a
bit. In order not to "change everything", we suggest the following,
namely, that the user is required to state which coefficient the
unknown is:

Ex2.ufl:
   V = FiniteElement("CG", "triangle", 1)
   f = Coefficient(V)
   u = Coefficient(V)
   v = TestFunction(V)
   F = u*dot(grad(u), grad(v))*dx - f*v*dx
   M = u*dx

   # Need to know what the unknown is
   unknown = u (!)


It is probably what we need. Is M the prefered name for an error
functional? It makes me think of a Matrix with a particular name...


Good point. We can easily adjust that.


I also first thought of making unknown a property of F, like:
   F.unknown = u

and thought this might be usefull in the python interface too.


That is another possibility. If so, u would also need to be attached to
M. I would prefer to keep the labelling as terse as possible.


Sure. For now just one error functional can be defined in a form file.
The one with name M, right? Would there be a need to be able to define
several functionals in the same file? If so should there be some sort of
form attribute for this purpose?


Anders and I have been through (but rejected) a couple of other ideas,
including

(1) Making TrialFunction a special sub-class of Function instead of of
Argument

     (Very sensible, consistentifies linear and nonlinear, but way
     radical)


What would a functional of u then be?

   u = TrialFunction(V)
   f = assemble(u*dx) ?


Maybe.


Which means?


(2) Introducing an additional sub-class Unknown of Function, thus
identifying the unknown

     (Looks a bit weird, does not increase consistency with linear

problems.)


I think this could work. The point with a linear problem is that the
problem can be described without the Unknown (=>  Linear). Which means
that a linear problem wont have an Unknown.


Yes, but linear problems have a TrialFunction which plays the same role.


But isn't a TrialFunction only meaningfull in a bilinear form?


But this is what you use the 'u' argument in pde.solve(u) to indicate?


Yes and no. Without ever having looked at the nonlinear solver in
DOLFIN, I assume that the  u argument in pde.solve(u) is (primarily)
needed for the updating in the Newton loop.  The same need is there for
the adaptive solve.


Sure, but you also use it to figure our the Unknown in the pure python
code? How could you otherwise do it?


In the current dolfin Python adaptivity module? We give "u" as an
argument to AdaptiveVariationalProblem and do all the code generation on
the fly there.


Ok.
Just as a side note, we first consider changing the meaning of TrialFunction to be a special subclass of Expression which (a Function in DOLFIN) instead of an Argument. So one would always define a nonlinear problem and it would always be linearized. That would solve some of the issues regarding the need for redefining u when a linear problem is solved: 

u = TrialFunction(V)  # now it's a trial function
...
u = problem.solve()   # now it's something different

But we decided this would require too many changes.

--
Anders

Follow ups

References

  Thread Previous • Date Previous • • Thread Next