# Is Wolfram's example of solving PDEs in 3D inconsistent?

I am slightly puzzled over the example "Solve PDEs over 3D Regions" provided by Wolfram. The example shows how to solve Laplace's equation over a space shuttle mesh. However, if I run the code in 11.0.1.0 the solution produced appears to be completely homogeneous:

mr = BoundaryDiscretizeGraphics[
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} ∈ mr];


Then

vals = uif["ValuesOnGrid"];
{Max[vals], Min[vals]}


produces the output {1.,1.}, suggesting that the temperature field is homogeneous. Somehow, ElementMeshSurfacePlot3D[uif, Boxed -> False, ViewPoint -> {0, -4, 2}] still produces the beautiful coloured plot of the documentation. This confuses me.

Here is a partially successful modified example with slightly altered Dirichlet boundary conditions:

uif = NDSolveValue[{Inactive[Laplacian][u[x, y, z], {x, y, z}] == 0,
DirichletCondition[u[x, y, z] == 1, -1 <= x <= 1],
DirichletCondition[u[x, y, z] == -1, -5 <= x <= 4]},
u, {x, y, z} ∈ mr];


The ElementMeshSurfacePlot3D suggests that a plausible solution was computed:

With[{vals = uif["ValuesOnGrid"]},
Legended[ElementMeshSurfacePlot3D[uif, ViewPoint -> {0, -4, 2},
Axes -> True, AxesLabel -> {X, Y, Z}],
BarLegend[{ColorFunction /. Options[ElementMeshSurfacePlot3D] //
First, MinMax@vals}]]] DensityPlot[uif[X, Y, 0.4], {X, -7, 7}, {Y, -4, 4}, PlotPoints -> 70,
PlotLegends -> Automatic, FrameLabel -> {X, Y, Z},
PlotLabel -> "Slice at Z=0.4"] // Quiet However, a more careful look suggests that my Dirichlet condition is not exactly met at X=0:

Plot[uif[X, 0, 0.4], {X, -7, 7}, AxesLabel -> {X, u},
PlotLabel -> "Slice at Y=0, Z=0.4"] // Quiet I would really appreciate your help fully understanding this example, as I also intend to solve FEM problems on BoundaryMeshRegions imported from STLs.

• I'm also confused by the example from Wolfram. First problem : how is it possible that the Dirichlet condition DirichletCondition[u[x, y, z] == 1, z <= -1.3] is outside mr ??? Aug 20, 2017 at 16:49
• The first part is most likely a typo. The equation should have been a Poisson equation; since the visualization kind of looks right, it was not noticed. The fact that it kind of looks right is because of a numerical curiosity. A perfect example that garbage in does not always mean garbage out.... The second issue you raise, I can not reproduce in 11.1 and I no longer have a 11.0 around to check. I get a flat line at 1. Aug 21, 2017 at 17:32
• @andre, what do you mean with outside? Aug 21, 2017 at 17:40
• @user21 I mean : -1.3 is under the lowest point of the shuttle Aug 21, 2017 at 17:43
• The Example page is updated now. Thanks for the feedback. Sep 11, 2017 at 19:58

You question has two parts. The example on the page you reference was not quite correct (EDIT- It's fixed now.) There are many ways to fix this. One idea is to have two distinct DirichletConditions like this:

mr = BoundaryDiscretizeGraphics[
ExampleData[{"Geometry3D", "SpaceShuttle"}]];
uif = NDSolveValue[{Inactive[Laplacian][u[x, y, z], {x, y, z}] == 0,
DirichletCondition[u[x, y, z] == 1, z <= -1.3],
DirichletCondition[u[x, y, z] == 0, z >= 1]},
u, {x, y, z} \[Element] mr];
Needs["NDSolveFEM"]
ElementMeshSurfacePlot3D[uif, Boxed -> False, ViewPoint -> {0, -4, 2}] Another way to fix this is to use a different PDE, e.g. use Poisson's equations in stead.

The fact that ElementMeshSurfacePlot3D did show a visualization that seemed valid is because the values computed do have slight variations. You will see that if you copy and paste the values in the notebook. What is surprising to me too is that this noise does seem to have a 'direction' that looks like a correct flow... Perfect example where garbage in does not result in garbage out. I'll try to get the web page updated.

For the second point: Maybe I am missing the point but the coordinate {0,0,0.4} is inside the region and not on the boundary.

emsf = ElementMeshSurfacePlot3D[uif, Boxed -> False,
ViewPoint -> {0, -4, 2}];
Show[emsf, Graphics3D[{Red, Ball[{0, 0, 0.4}, 0.5]}],
PlotRange -> {Automatic, {0, 5}, {-2, 5}}]


That, I think, explains why the point does not match the value on the boundary. • This is perfect, thank you, explains everything. Aug 22, 2017 at 15:20