[dune-pdelab] Dune PDELab (2.5) Constraints questions

Steffen Müthing steffen.muething at iwr.uni-heidelberg.de
Tue Jun 12 13:18:41 CEST 2018

Hi Michael,

> Am 12.06.2018 um 11:07 schrieb Michael Wenske <m_wens01 at wwu.de>:
> Dear all!
> I am struggling with the Constraints implementation and the assembly of
> such in Dune::PDELab.
> I want to constrain all degrees of freedom within a subvolume(tumor) in
> the inside of my domain(brain).
> My degrees of freedom are an the edges of a quadrilateral Grid and are
> assembled, via the alpha_volume(...)
> methods in the local operator. I've studied files like
> pdelab/constraints/conforming.hh and their interface seems to imply that
> only fluxes over an -interface-  are
> to be dirichlet constrained. Now theoretically there should be no
> difference wether a dof is restricted by a loop
> which is iterating over the interfaces, or if it is iterating the volumes.

The currently implemented constraints assemblers all work on intersections
because that is what you need most of the time. No one has needed a volume
constraints assembler until now, so there isn’t one in PDELab.

> Question 1: can I use the skeleton constraints and volume constraints
> interchangeably as long as i
> hit the right dof's and set their transformation to empty?

But you are right, you can write your own constraints assembler (like you did)
and have that constrain DOFs on a cell instead of DOFs on an intersection.

> I see that there is an interface to pass an object to the constraints
> assembler which might be used
> to let the user define the constraints via a simple lambda. something like:
> using CON = Dune::PDELab::OverlappingConformingDirichletConstraints;
> using CC = typename GFS::template
> ConstraintsContainer<RangeFieldType>::Type;
> auto bc_lambda = [&](const auto& i,const auto& x) {return
> modelLOP_p->b(i,x);};
> auto b = Dune::PDELab::makeBoundaryConditionFromCallable(*gv_p,bc_lambda);
> Dune::PDELab::constraints(b,*gfs_p, *cc_p);
> My problem is that the above function will derive from
> DirichletConstraintsParameters, again
> being restricted to the outer boundary and not the inside. On a wider
> scale I see the ConstraintsParameters
> interfaces which I take are there to define the type of the object and
> signature of the methods contained in b.
> This object is then passed all the way down to the actual method
> deciding isDirichlet or isNeumann etc.
> 2. Is the XXXParameters concept to have the user write these
> interfaces...and then implement the constraints class using it him/herself?
> How would I use the Parameter Interface so that my constraints would
> work with existing implementations?

The existing parameter classes are written for the more common concept of boundary
constraints, so they will not work in your case. So it’s like you guessed: You have to write
your own constraints assembler, and you can design the parameter interface any way
you want.

> I did'nt give up yet. I wrote my own constraints class. It's a little
> class fulfilling (afaik) the interface for the
> VolumeConstraints given in pdelab/constraints/common/constraints.hh (see
> file)
> it will get an elementGeometry from the assembler, decide if this dof is
> to be restricted or not by the position of
> the dof's (isTumor(pos)) and then set the transformation to empty. I can
> debug right into the line where the
> transformation is set to empty.
>             // check if this corner-dof is supposed to be constrained
>             if (medData_p->isTumor(pos)) trafo[lfs.dofIndex(i)] = empty;
> I can not grasp where i went wrong, and I guess there is some conceptual
> thing
> i am missing about the way these interfaces are written. When i run the
> program it tells me that there are -zero-
> constrained dof's in my gridfunctionspace even though, when I debug the
> assembly procedure seems to take
> place correctly. Also the way I now do it, I do not use the
> ConstraintsParameter Object passed into my class at all.

That’s strange, looking at the attached code, that should do the right thing. What exactly
are you doing in your main driver function? Are you sure that the isTumor() call actually returns

> 3. If anybody could give a quick overview about the -idea- behind the
> setup of constraints in pdelab , or what i am currently missing, i would
> be most
> grateful. What would be the cleanest way for a user to implement the
> constriction of dof's at different positions in the domain?

I think you already managed to figure it out pretty much on your own… :-)

> I would also be gratefull for hits to any code examples where people
> have implemented interior boundary conditions before.

As I said: I’m personally not aware of such examples.


> thank you in advance,
> Michael
> <SegmentationConstraints.hh>_______________________________________________
> dune-pdelab mailing list
> dune-pdelab at lists.dune-project.org
> https://lists.dune-project.org/mailman/listinfo/dune-pdelab

-------------- next part --------------
A non-text attachment was scrubbed...
Name: signature.asc
Type: application/pgp-signature
Size: 488 bytes
Desc: Message signed with OpenPGP
URL: <https://lists.dune-project.org/pipermail/dune-pdelab/attachments/20180612/af47da7b/attachment.sig>

More information about the dune-pdelab mailing list