<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1">
<style type="text/css" style="display:none;"><!-- P {margin-top:0;margin-bottom:0;} --></style>
</head>
<body dir="ltr">
<div id="divtagdefaultwrapper" style="font-size: 12pt; color: rgb(0, 0, 0); font-family: Calibri, Helvetica, sans-serif, "EmojiFont", "Apple Color Emoji", "Segoe UI Emoji", NotoColorEmoji, "Segoe UI Symbol", "Android Emoji", EmojiSymbols;" dir="ltr">
<p><br>
</p>
<meta content="text/html; charset=UTF-8">
<div dir="ltr">
<div id="x_divtagdefaultwrapper" dir="ltr" style="font-size: 12pt; color: rgb(0, 0, 0); font-family: Calibri, Helvetica, sans-serif, "EmojiFont", "Apple Color Emoji", "Segoe UI Emoji", NotoColorEmoji, "Segoe UI Symbol", "Android Emoji", EmojiSymbols;">
<p>Hey Carsten,</p>
<p><br>
</p>
<p>ok, then i will use the insertionIndex() to be on the save side.</p>
<p>Right now, Load is a vector with fixed values so I guess that</p>
<p>should work fine.</p>
<p><br>
</p>
<p>So for a minimal example, let's consider the domain given in the</p>
<p>attached picture. This could be a slender beam or plate or something</p>
<p>similar. On one side it is fixed at a wall (tag 1 for that B.C.), on the other</p>
<p>end a loadvector, which could be a force is applied (tag 3 for that B.C.).</p>
<p>The other boundary faces are tagged with 0, meaning that no B.C. should</p>
<p>be applied there.</p>
<p><br>
</p>
<p>The domain is meshed with one hexaeder element and let's consider linear</p>
<p>shape functions. The mesh is loaded via the gmsh-reader of DUNE taking</p>
<p>the boundary-segments into consideration. Printing out the vector with the</p>
<p>boundary tags they are given as:</p>
<p><br>
</p>
<p></p>
<div>B.C. tag: 0<br>
B.C. tag: 1<br>
B.C. tag: 3<br>
B.C. tag: 0<br>
B.C. tag: 0<br>
B.C. tag: 0</div>
<p></p>
<p><br>
</p>
<p>so i guess they are just stored in the way appearing in the .msh file (</p>
<div>consecutively ordered). Now I want to check how DUNE sees the boundary</div>
<div>segments and which face sees which tag.</div>
<p></p>
<p><br>
</p>
<p></p>
<div>for( const auto& element : elements(gridView)) {<br>
for( const auto& isect : intersections(gridView, element)) {<br>
if(isect.boundary()) {<br>
<br>
auto ref = referenceElement<double, dim>(element.type());<br>
<br>
cout << "B.C. seg. index: " << isect.boundarySegmentIndex() << endl;<br>
<br>
for(int i=0; i<ref.size(isect.indexInInside(), 1, dim); i++) {<br>
cout << "Coordinate of " << i << " is: "<br>
<< <span>isect.inside().geometry().corner(ref.subEntity(isect.indexInInside(), 1, i, dim))</span><br>
<< endl;<br>
}</div>
<div> }<br>
</div>
}
<p></p>
<p><br>
</p>
<p>This generates the following output listing which index</p>
<p>the actual boundary-segment has and lists the coordinates</p>
<p>of the points that span the boundary face (hopefully they are</p>
<p>printed in the right way):</p>
<p><br>
</p>
<p></p>
<div>
<div>
<div>
<div>B.C. seg. index: 2<br>
Coordinate of 0 is: 0 0 0<br>
Coordinate of 1 is: 5 0 0<br>
Coordinate of 2 is: 0 0.1 0<br>
Coordinate of 3 is: 5 0.1 0<br>
B.C. seg. index: 0<br>
Coordinate of 0 is: 0 0 0<br>
Coordinate of 1 is: 5 0 0<br>
Coordinate of 2 is: 0 0 1<br>
Coordinate of 3 is: 5 0 1<br>
B.C. seg. index: 4<br>
Coordinate of 0 is: 5 0 0<br>
Coordinate of 1 is: 5 0.1 0<br>
Coordinate of 2 is: 5 0 1<br>
Coordinate of 3 is: 5 0.1 1<br>
B.C. seg. index: 5<br>
Coordinate of 0 is: 0 0.1 0<br>
Coordinate of 1 is: 5 0.1 0<br>
Coordinate of 2 is: 0 0.1 0<br>
Coordinate of 3 is: 5 0.1 1<br>
B.C. seg. index: 1<br>
Coordinate of 0 is: 0 0 0<br>
Coordinate of 1 is: 0 0.1 0<br>
Coordinate of 2 is: 0 0 1<br>
Coordinate of 3 is: 0 0.1 1<br>
B.C. seg. index: 3<br>
Coordinate of 0 is: 0 0 1<br>
Coordinate of 1 is: 5 0 1<br>
Coordinate of 2 is: 0 0.1 1<br>
Coordinate of 3 is: 5 0.1 1<br>
</div>
<br>
</div>
</div>
</div>
Checking with the domain from above, the insertion index chooses
<p></p>
<p>the tag for the boundary faces. So boundary segment face</p>
<p>with tag 1 corresponds to insertion index 1 that fits with the domain, <br>
</p>
<p>looking at the corresponding coordinates.</p>
<p>Now, tag 3 should correspond to insertion index 2, but the coordinates</p>
<p>behind that boundary face are not right, they don't fit with the picture</p>
<p>of the domain. So the insertion index, the boundary tag and the corresponding</p>
<p>boundary face coordinates don't match in the right way.</p>
<p>This results in wrong B.C. for the overall simulation, due to the mismatch.<br>
</p>
<p><br>
</p>
<p>Coming back to the "not working", this is what i mean/meant with that. <br>
</p>
<p>In 2d, using the same code just with the different dimension, i didn't face such a problem.</p>
<p><br>
</p>
<p>This is quite a lengthy answer, but i hope it clarifies my problem.</p>
<p>I'm not able to attach the msh.file somehow, if the file is needed i will post it.</p>
<p><br>
</p>
<p>Thanks for the help.</p>
<p>Max<br>
</p>
<p><br>
</p>
<p><br>
</p>
</div>
<hr style="display:inline-block; width:98%" tabindex="-1">
<div id="x_divRplyFwdMsg" dir="ltr"><font style="font-size:11pt" face="Calibri, sans-serif" color="#000000"><b>Von:</b> Carsten Gräser <graeser@mi.fu-berlin.de><br>
<b>Gesendet:</b> Montag, 24. August 2020 10:34:50<br>
<b>An:</b> Firmbach, Max; dune@lists.dune-project.org<br>
<b>Betreff:</b> Re: [Dune] Boundary conditions, segments and indexing</font>
<div> </div>
</div>
</div>
<font size="2"><span style="font-size:10pt">
<div class="PlainText">Hi Max,<br>
<br>
Am 24.08.20 um 09:46 schrieb Firmbach, Max:<br>
[...]<br>
> For the simple one element, hexaeder case the tags are read into<br>
> the boundaryEntity vector like [0, 0, 3, 0, 1, 0]. So on four faces<br>
> there is just a 0, which means no B.C. is applied. For front and back<br>
> either B.C. 1 or 3 should be applied. Looping over the boundary<br>
> intersections my insertion index for the segments is something like<br>
> [0, 1, 2, 4, 3, 5], which is identical to the boundarySegmentIndex()<br>
> method (somewhere it's said that they don't need to be identical<br>
> somehow, for this case they are).<br>
in general they are not, because a grid may renumber them<br>
during creation. The GridFactory::insertionIndex() method<br>
exists, to keep track of such renumberings. However, UGGrid<br>
does not do any renumbering, as you observed. So that's not<br>
the problem here (but you may still adjust your code to not<br>
depend on implementation details).<br>
<br>
There's another potential issue in your code. But since<br>
you've only send a code fragment, it's not clear if it<br>
really is a problem:<br>
<br>
> for( const auto& element : elements(gridView)) {<br>
> for( const auto& isect : intersections(gridView, element)) {<br>
> <br>
> GeometryType elementGeometry = element.type();<br>
> using Factory = ReferenceElements<double, dim>;<br>
> const auto &ref = Factory::general(elementGeometry);<br>
[...]<br>
> for(int i=0; i<ref.size(isect.indexInInside(), 1, dim); i++) {<br>
> int row = indexSet.subIndex(element, ref.subEntity(isect.indexInInside(), 1, i, dim), dim);
<br>
> loadVector[row] = Load; <br>
If Load is a fixed value, that's OK. If it's computed depending<br>
on i, then this is probably a bug, because the subentity<br>
orientation in the reference element and grid elements<br>
does in general not coincide. To avoid ever using i<br>
directly better use re.subEntities(isect.indexInInside(), 1, dim)<br>
instead.<br>
<br>
<br>
> As written before, for 2d meshes it works pretty well,<br>
> but for 3d it's the complete opposite and I don't really<br>
> understand why.<br>
Can you be a little more precise and describe what<br>
"it does not work" means? A minimal complete example<br>
would also be helpful.<br>
<br>
Best,<br>
Carsten<br>
</div>
</span></font></div>
</body>
</html>