7
$\begingroup$

I am trying to solve heat equation on a cylinder whose ends are thermally insulated and its circular face is exposed to convection. Therefore I have Neumann boundary condition on all faces of the cylinder. Here is my code:

NDSolveValue[{Laplacian[u[x, y, z], {x, y, z}] == NeumannValue[0, z == 0] + 
 NeumannValue[0, z == 1] + NeumannValue[1 - u[x, y, z],
 x^2 + y^2 == 1]}, u, {x, y, z} \[Element] Cylinder[{{0, 0, 0}, {0, 0, 1}}, 1]]

However I get the error:

NDSolveValue::bcnop: No places were found on the boundary where Coordinate was True, so BoundaryCondition[{Robin,{1,1},{CompiledFunction[{10,11.,5568},{_Real,_Real,_Real},{{3,0,0},{3,0,1},{3,0,2},{3,2,0}},{{{{1.}},{3,2,0}}},{0,0,3,0,1},{{1}},Function[{x,y,z},{{1.}},Listable],Evaluate],CompiledFunction[{10,11.,5568},{_Real,_Real,_Real},{{3,0,0},{3,0,1},{3,0,2},{3,2,0}},{{{{-1.}},{3,2,0}}},{0,0,3,0,1},{{1}},Function[{x,y,z},{{-1.}},Listable],Evaluate]},1},Coordinate,CompiledFunction[{10,11.,5568},{_Real,_Real,_Real},<<5>>,Evaluate],NeumannValue[1-u,x^2+y^2==1]] will effectively be ignored.

My question: what does No places were found on the boundary where Coordinate was True mean? As far I can see I have specified the faces of the cylinder correctly in the NeumannValue boundary conditions. Any help is much appreciated.

$\endgroup$
4
  • 2
    $\begingroup$ Use NeumannValue[1 - u[x, y, z], x^2 + y^2 >= 0.9999]. $\endgroup$ – flinty Oct 10 '20 at 12:40
  • $\begingroup$ @flinty NeumannValue[1 - u[x, y, z], x^2 + y^2 >= 0.99] works for me, but any number greater than 0.99 doesn't work and I get the same error as before. What's going on here? $\endgroup$ – Deep Oct 10 '20 at 12:47
  • 1
    $\begingroup$ I suspect it could be due to the discretization of the cylinder. $\endgroup$ – flinty Oct 10 '20 at 12:49
  • $\begingroup$ @flinty Thanks a lot. How can I refine the discretization? $\endgroup$ – Deep Oct 10 '20 at 12:53
5
$\begingroup$

To avoid having to tweak the boundary condition, load the finite elements package and make an actual mesh:

<< NDSolve`FEM`

mesh = ToElementMesh @ Cylinder[{{0, 0, 0}, {0, 0, 1}}, 1];

NDSolveValue[
  {
    Laplacian[u[x, y, z], {x, y, z}] ==
      + NeumannValue[0, z == 0]
      + NeumannValue[0, z == 1]
      + NeumannValue[1 - u[x, y, z], x^2 + y^2 == 1]
  },
  u,
  Element[{x, y, z}, mesh]
]

By default the mesh will be second-order, and perhaps this is why it is able to handle the curved boundary properly. It seems that ToElementMesh is able to handle curved boundaries much better than the default discretisation method used by NDSolveValue.

$\endgroup$
2
  • $\begingroup$ Thanks a lot. Your solution works. Do you know how I can make the mesh finer? $\endgroup$ – Deep Oct 11 '20 at 4:38
  • 1
    $\begingroup$ @Deep Take a look at the Element Mesh Generation tutorial and the options for ToElementMesh. Most useful to you will be the options MaxCellMeasure and "MaxBoundaryCellMeasure". $\endgroup$ – yawnoc Oct 11 '20 at 12:23
6
$\begingroup$

Update OpenCascadeLink Free Solution

The OP indicated that they are on version 11.0 that does not include the OpenCascadeLink. I do not have version 11, so I do not know if this works, but it does not depend on OpenCascade. Note that the curved surface came out as ElementMarker==3 in this case.

Needs["NDSolve`FEM`"]
c1 = Cylinder[{{0, 0, 0}, {0, 0, 1}}, 1];
bmesh = ToBoundaryMesh[c1];
groups = bmesh["BoundaryElementMarkerUnion"];
temp = Most[Range[0, 1, 1/(Length[groups])]];
colors = ColorData["BrightBands"][#] & /@ temp
bmesh["Wireframe"["MeshElementStyle" -> FaceForm /@ colors]]
mesh = ToElementMesh[bmesh];
mesh["Wireframe"]
ufun = NDSolveValue[{Laplacian[u[x, y, z], {x, y, z}] == 
     NeumannValue[1 - u[x, y, z], ElementMarker == 3]}, 
   u, {x, y, z} ∈ mesh];
SliceContourPlot3D[
 ufun[x, y, z], "CenterPlanes", {x, y, z} ∈ 
  Cylinder[{{0, 0, 0}, {0, 0, 1}}, 1], PlotLegends -> Automatic]

OpenCascade Free Solution

You can use OpenCascadeLink define the geometry and it will create ElementMarkers to the faces that you may refer to in your boundary condition specification. This will avoid guessing at what discretization is required when the object or scale changes.

Here is an example. Note that the $\color{Red}{Red\ Surface}$ corresponds to the curved face and is ElementMarker==1.

Needs["OpenCascadeLink`"]
Needs["NDSolve`FEM`"]
cyl = OpenCascadeShape[c1 = Cylinder[{{0, 0, 0}, {0, 0, 1}}, 1]];
bmesh = OpenCascadeShapeSurfaceMeshToBoundaryMesh[cyl];
groups = bmesh["BoundaryElementMarkerUnion"];
temp = Most[Range[0, 1, 1/(Length[groups])]];
colors = ColorData["BrightBands"][#] & /@ temp
bmesh["Wireframe"["MeshElementStyle" -> FaceForm /@ colors]]
mesh = ToElementMesh[bmesh];
mesh["Wireframe"]
ufun = NDSolveValue[{Laplacian[u[x, y, z], {x, y, z}] == 
     NeumannValue[1 - u[x, y, z], ElementMarker == 1]}, 
   u, {x, y, z} ∈ mesh];
SliceContourPlot3D[
 ufun[x, y, z], "CenterPlanes", {x, y, z} ∈ 
  Cylinder[{{0, 0, 0}, {0, 0, 1}}, 1], PlotLegends -> Automatic]

Workflow

$\endgroup$
3
  • $\begingroup$ Unfortunately I am unable to load OpenCascadeLink. I am using Mathematica 11.0. $\endgroup$ – Deep Oct 11 '20 at 4:17
  • $\begingroup$ @Deep Look at the update and see if it works. I do not have v11, so I cannot test. $\endgroup$ – Tim Laska Oct 11 '20 at 15:21
  • $\begingroup$ Thank you very much. $\endgroup$ – Deep Oct 12 '20 at 3:21
3
$\begingroup$

I think it has to do with discretizing the region. Consider:

NDSolveValue[{Laplacian[u[x, y, z], {x, y, z}] == 
   NeumannValue[0, z == 0] + NeumannValue[0, z == 1] + 
    NeumannValue[1 - u[x, y, z], 0.999 <= x^2 + y^2 <= 1.001]}
 , u, {x, y, z} \[Element] Cylinder[{{0, 0, 0}, {0, 0, 1}}, 1]]

This produces your error. However, if we soften up the condition x^2+y^2==1 a bit, then it works:

NDSolveValue[{Laplacian[u[x, y, z], {x, y, z}] == 
   NeumannValue[0, z == 0] + NeumannValue[0, z == 1] + 
    NeumannValue[1 - u[x, y, z], 0.99 <= x^2 + y^2 <= 1.01]}
 , u, {x, y, z} \[Element] Cylinder[{{0, 0, 0}, {0, 0, 1}}, 1]]
(*InterpolatingFunction[{{-1., 1.}, {-1., 1.}, {0., 1.}}, <>]*)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.