Local assembly of forms for a posteriori error analysis

I'd like to add to the initial question that often the local domains are patches (e.g. the union of all cells sharing a vertex) rather than single cells. In that case, the DG approach does not work and instead one would need an efficient way to just assemble and solve local problems on a specific set of cells. I had the impression that dolfin currently is not well suited for this kind of calculation which is required for classical implicit error estimators and also for the very efficient recent equilibration error estimators. For the latter, one would determine equilibrated fluxes in an appropriate space (e.g. H(div)) on the patches.

Cheers, Martin

For the intention of the TO, the following questions might prove helpful
https://answers.launchpad.net/dolfin/+question/203757
https://answers.launchpad.net/dolfin/+question/177108

Thanks Martin Sandve Alnæs, that solved my question.

Rerering to Martin Eigels first answer (Post #2) I would also like to see support of patch-wise assembly in dolfin. In some variants of the Dual Weighted Residual error estimator it is also needed to perform patch-wise interpolations of the discrete solution to compute an dual solution which is not an element of the primal solution space.

Thank you Martin for the follow ups. They helped a lot.

Let me just comment that although dolfin does not have a native patch-wise assembler,
the necessary information should be there via the ufl-ffc-ufc toolchain. The patch-wise
assembly could be accomplished by writing your own assembler which iterates
over patches, computes contributes from each cell/facet and adds to a suitable matrix/vector/tensor.

For computing a dual solution of higher-order than the primal solution space,
I would recommend taking a look at the Extrapolation class. This provides functionality for computing an extrapolated
solution based on least-squares curve fitting. This is what we use as the default approach
for computing an improved dual solution for the automated dual-weighted-residual error estimates.

On Wed, Oct 31, 2012 at 08:50:57AM -0000, Marie Rognes wrote:
> Question #212729 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/212729
>
> Marie Rognes posted a new comment:
> Let me just comment that although dolfin does not have a native
> patch-wise assembler, the necessary information should be there via
> the ufl-ffc-ufc toolchain. The patch-wise assembly could be
> accomplished by writing your own assembler which iterates over
> patches, computes contributes from each cell/facet and adds to a
> suitable matrix/vector/tensor.

Good answer. Often our response to a user wanting to do something
slightly more complex than just standard assembly, is it can't be
done, but in fact most things can be done. It's just that we don't
have a function for it. The entire interface of FEniCS (in particular
DOLFIN) has been designed to be easy to work with, which means that
the implementation of the assembler in DOLFIN itself is pretty
short. So take a look at Assembler.cpp, implement your own iteration
over cells, iterate over neighbors etc to solve the local problems.

--
Anders

> For computing a dual solution of higher-order than the primal
> solution space, I would recommend taking a look at the Extrapolation
> class. This provides functionality for computing an extrapolated
> solution based on least-squares curve fitting. This is what we use
> as the default approach for computing an improved dual solution for
> the automated dual-weighted-residual error estimates.
>

I've written a small example how swig can be used for communication between python and cpp (especially with dolfin). Perhaps this is useful to someone. If there are any questions please contact me.

http://scwww.math.uni-augsburg.de/~tommy/python_dolfin_swig_demo.tar.gz

Also have a look at compile_extension_module to do it JIT within a
Python sesion. Try:

  from dolfin import *
  mesh = UnitSquare(4, 4)
  V = FunctionSpace(mesh, "Lagrange", 1)
  u = interpolate(Constant(1.0), V)
  code = """
  void set_to_zero(boost::shared_ptr<dolfin::GenericVector> x) {
   x->zero();
  }
  """
  compiled_module = compile_extension_module(code)
  compiled_module.set_to_zero(u.vector())
  plot(u)
  interactive()

Johan

On 11/05/2012 08:51 AM, Thomas Fraunholz wrote:
> Question #212729 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/212729
>
> Thomas Fraunholz posted a new comment:
> I've written a small example how swig can be used for communication
> between python and cpp (especially with dolfin). Perhaps this is useful
> to someone. If there are any questions please contact me.
>
> http://scwww.math.uni-augsburg.de/~tommy/python_dolfin_swig_demo.tar.gz
>

That's nice! A lot easier...

On 11/05/2012 10:25 PM, Thomas Fraunholz wrote:
> Question #212729 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/212729
>
> Thomas Fraunholz posted a new comment:
> That's nice! A lot easier...

:)

FYI: compiled_expressions and compiled_subdomains all uses this interface.

Johan