[Dune] Boundary conditions, segments and indexing

[Dedner, Andreas]

Hi Max.
Looks alright to me but I might have missed something.
Could you please add which grid manager you are using since it might be a bug in the grid implementation you are using.
Best
Andreas
________________________________
From: Dune <dune-bounces at lists.dune-project.org> on behalf of Firmbach, Max <max.firmbach at tum.de>
Sent: 21 August 2020 17:58
To: dune at lists.dune-project.org <dune at lists.dune-project.org>
Subject: [Dune] Boundary conditions, segments and indexing


Hello DUNE-developers,


I have a question regarding the boundary segment indexing and

connecting tags corresponding to the segments with boundary conditions.


Right now, I load a gmsh-file and obtain a vector with boundary tags.

Those tags should trigger some kind of boundary condition for example

set some values in a block-vector or sparse-matrix.

After reading the book from Mr. Sander and some how-to's I think I

understood the principle concept of intersections and so, but still

something is not working right for me.

My actual code loops over all elements of the gridview and over

the corresponding intersections. Inside the loop i check if

the intersection is in contact with the boundary (because I

want to set boundary conditions). If that is true, a switch

should select the right b.c. code, depending on the tag that is

stored inside the boundaryEntitity vector I obtain from loading the

gmsh-file. The array-entry of the boundaryEntitiy vector is

selected via the intersection boundary segment index (is this the right

approach ?).


And here starts my actual question ... in two dimensions that is working

fine (even for complex geometries and several tags), but for three dimensions

it's not working properly even for a single element (for example a hexa-element

with different tags on each of the boundary faces).

So is my description the common approach to handle the boundary situation

or is there a more efficient or better way ? Or am I missing something ?


My code looks like this (the 3 lines after the double loop and the

code inside case 1 should set the corresponding entry in a block-

vector, by looping over the vertices of the intersection geometry and

getting the "global" index for that position):


for( const auto& element : elements(gridView)) {
    for( const auto& isect : intersections(gridView, element)) {


      GeometryType elementGeometry = element.type();
      using Factory = ReferenceElements<double, dim>;
      const auto &ref = Factory::general(elementGeometry);


      if(isect.boundary()) {

        switch(boundaryEntity[isect.boundarySegmentIndex()]) {
          case 0:
            ...
            break;
          case 1:

            for(int i=0; i<ref.size(isect.indexInInside(), 1, dim); i++) {
              int row = indexSet.subIndex(element, ref.subEntity(isect.indexInInside(), 1, i, dim), dim);
              loadVector[row] = Load;
            }

           ...
        }
      }
    }
  }

As written before, for 2d meshes it works pretty well, but for 3d it's the complete
opposite and I don't really understand why.

Thanks for your help and best regards,
Max

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.dune-project.org/pipermail/dune/attachments/20200821/e8f43272/attachment.htm>