Message Boards Message Boards

Intersecting boundary problems when solving PDEs

Hi everyone,

the Wolfram Language both in the cloud and in Mathematica 10 offers the possibility to solve PDEs on quite general surfaces. In the documentation it is shown that the Laplace equation can be solved on a 3D model of the space shuttle like this:

mr = DiscretizeGraphics[ExampleData[{"Geometry3D", "SpaceShuttle"}]];

uif = NDSolveValue[{Inactive[Laplacian][u[x, y, z], {x, y, z}] == 0, 
    DirichletCondition[u[x, y, z] == 1, z <= -1.3]}, 
   u, {x, y, z} \[Element] mr];

Needs["NDSolve`FEM`"]
ElementMeshSurfacePlot3D[uif, Boxed -> False, ViewPoint -> {0, -4, 2}]

This is really nice and gives

enter image description here

The problem shows up when we take more irregular objects like 

mr = DiscretizeGraphics[ExampleData[{"Geometry3D", "Cow"}]]

enter image description here

or 

mr = DiscretizeGraphics[ExampleData[{"Geometry3D", "Horse"}]]

When you execute the integration you get:

NDSolveValue::fememib: The input has or generated an intersecting boundary and cannot be processed. >>

Similarly, 

BoundaryDiscretizeGraphics[ExampleData[{"Geometry3D", "Cow"}]]

gives

BoundaryMeshRegion::binsect: "The boundary curves self-intersect or cross each other in BoundaryMeshRegion[{{0.10799399763345718`,0.08466669917106628`,-0.21893799304962158`},{0.1157900020480156`,0.08869930356740952`,-0.22847199440002441`},{0.10719799995422363`,0.09556479752063751`,-0.23278899490833282`},<<46>>,{-0.30769699811935425`,0.093573197722435`,-0.22068199515342712`},<<2853>>},{{},{},{Polygon[{<<1>>}]}}]. \!\(\*ButtonBox[\">>\",
Appearance->{Automatic, None},
BaseStyle->\"Link\",
ButtonData:>\"paclet:ref/BoundaryMeshRegion\",
ButtonNote->\"BoundaryMeshRegion::binsect\"]\)"

So the problem is that there are self-intersections of the surface. Is there any elegant approach to fixing that?

Cheers, Marco
POSTED BY: Marco Thiel
5 Replies
Sort By:

Hi everyone, I have a similar problem that Mathematica fails to solve. Would any options help? I want to solve 3D Laplace equation where one of the boundaries of the region is dynamic, i.e. I can update it and make a next cycle. I have experimented with different shapes of the surface. For rather smooth ones it works. For the ones with many edges it fails by generating self-intersections.
The code below works.

xmin = -1; xmax = +1; ymin = -1; ymax = 1;
shape[x_, y_] := 0.05*Sin[15 x + 0.8]*Sin[15 y + 0.8];
realshape[x_, y_] = Piecewise[{{shape[x, y], shape[x, y] > 0}}, 0];
bottomregion = 
  ImplicitRegion[-0.1 < z < 
    shape[x, y], {{x, xmin, xmax}, {y, ymin, ymax}, {z, -0.1, 
     xcut}}];
ToElementMesh[meshbottom, "MaxCellMeasure" -> 0.003]

but if I replace shape with realshape, then it fails.

POSTED BY: Sergey Slizovskiy

I am working on a similar problem and decided to try and fix the "Cow". Your suggestion of eliminating the self intersecting boundary curves doesn't completely solve the problem with DiscretizeGraphics or BoundaryDiscretizeGraphics on the "Cow" example model.

After using the following routines to iteratively select out the 83 intersecting facets: 

pts=ExampleData[{"Geometry3D","Cow"}, "VertexData"];
facets=Partition[ExampleData[{"Geometry3D","Cow"}, "PolygonData"],1];
getIntersect:=intersect={#}&/@TetGenDetectIntersectingFacets[pts, facets][[2]];
getFacet:=facets=Select[facets,!MemberQ[intersect,#]&];
cow = Graphics3D[GraphicsComplex[pts, Polygon[Flatten[facets, 1]]]]

Attempting to TetGenTetrahedralize cow, we get: TetGenTetrahedralize::reterr: Tetrahedralize returned an error, 3. >>

Some of this problem is simply due to the deletion of facets leaving holes in the mesh. So after an export to .STL and filling the holes (using MeshLab), the DiscretizeGraphics on the import then throws a different error 2.

This seem to have gotten fixed in version 10.2, which now successfully imports the fixed .STL with the BoundaryMeshRegion option: 

mr=Import["cowtest.stl", "BoundaryMeshRegion"]
uif@mr

BTW - there are 2 interesting facets that are not part of the set of 83 intersecting ones. They are in the tail of the cow that throw the following errors (until they are cleaned up): BoundaryMeshRegion::bcsm: There is a closed curve or surface that does not include all cells specified for a boundary. A boundary should consist of a single closed curve or surface. >>

Unfortunately, my 10.2 version doesn't load the ExampleData[{"Geometry3D", "Cow"}] (this may be a problem with my install or a bug, not sure). Nevertheless, I have attached a fixed .STL version of the cow (sans tail).

Attachments:
POSTED BY: J Gregory Moxness

Hi Gregory,

thanks a lot. Indeed, now it works:

mr = Import["~/Desktop/cowfix4.stl", "BoundaryMeshRegion"];
uif = NDSolveValue[{Inactive[Laplacian][
      u[x, y, z], {x, y, z}] == -0.01, 
    DirichletCondition[u[x, y, z] == 1, z <= 0.08]}, 
   u, {x, y, z} \[Element] mr];

Needs["NDSolve`FEM`"]
ElementMeshSurfacePlot3D[uif, Boxed -> False, ViewPoint -> {0, -4, 2}]

enter image description here

I noticed that you had to trim the cow's tail.

Thanks,

Marco
POSTED BY: Marco Thiel

Thank you very much for your help! That is very useful. I will check TetGenDetectIntersectingFacets out and see whether I can make it work. Thanks a lot!

I just wonder one thing. Some time ago I built a little model of cells in a "microfluidic chamber" - the details don't really matter. The idea is that you have "holes" in a rectangular shape. These holes have their own boundary conditions and you want to solve a diffusion equation. The simulation seemed to work fine in the MMA 10 prerelease version (see results here!) and it stopped working in the full version. If you execute the file that you can download using the link, and you execute it in MMA10 you will get an error message that there are intersecting boundaries and the execution fails.

Marco

PS: I am sorry that the notebook is not commented properly, but it should give an idea of the observation.

POSTED BY: Marco Thiel

That is generally a very hard problem and WL is unlikely in the near future to be able to solve 'healing' of intersections in 3D. Supposedly, there are "safe" and very expensive algorithms for doing this.

TetGenDetectIntersectingFacets returns the intersecting facets, then one can try to repair them manually.

POSTED BY: Bruce Miller
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Tag limit exceeded
Note: Only the first five people you tag will receive an email notification; the other tagged names will appear as links to their profiles.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract