6
$\begingroup$

Tell me please to solve the Laplace equation for the ring? Recorded the equation in polar coordinates, set the domain, Dirichlet boundary conditions, but outputs...

sol = NDSolveValue[
  { ρ^2 D[ u[ρ, ϕ], ρ, ρ] + ρ D[ u[ρ, ϕ], ρ] +  D[u[ρ, ϕ], ϕ, ϕ] == 0
  , DirichletCondition[u[ρ, ϕ] == 1000., ρ == .5 && 0 <= ϕ <= 2 π]
  , DirichletCondition[ u[ρ, ϕ] == 0., ρ == 10 && 0 <= ϕ <= 2 π]
  }
, u, {ρ, 0.5, 10}, {ϕ, 0, 2π}  
 ]

DensityPlot[sol[ρ, ϕ], {ρ, 0.5, 10}, {ϕ, 0, 2π }
, Mesh -> All
, ColorFunction -> "Rainbow"
, PlotLegends -> Automatic
]  

Output :

enter image description here

Note that the mesh above is not the mesh used by NDSolveValue. It is a mesh created by DensityPlot. The mesh used by NDSolveValue is :

DensityPlot[sol[\[Rho], \[Phi]], {\[Rho],\[Phi]} \[Element] sol["ElementMesh"]
, Mesh -> All
, ColorFunction -> "Rainbow"
, PlotLegends -> Automatic
]   

enter image description here

The question is : how to transform the coordinates of the mesh(es) so that we retrieve the correct geometry of the domain : a ring ?

$\endgroup$
  • $\begingroup$ Comments are not for extended discussion; this conversation has been moved to chat. $\endgroup$ – Kuba Nov 19 '17 at 16:53
5
$\begingroup$

Here is the solution to get the right geometry :

Mesh from DensityPlot

sol = NDSolveValue[
  { \[Rho]^2 D[ u[\[Rho], \[Phi]], \[Rho], \[Rho]] + \[Rho] D[ u[\[Rho], \[Phi]], \[Rho]] +  D[u[\[Rho], \[Phi]], \[Phi], \[Phi]] == 0
  , DirichletCondition[u[\[Rho], \[Phi]] == 1000., \[Rho] == .5 && 0 <= \[Phi] <= 2 \[Pi]]
  , DirichletCondition[ u[\[Rho], \[Phi]] == 0., \[Rho] == 10 && 0 <= \[Phi] <= 2 \[Pi]]
  }
, u, {\[Rho], 0.5, 10}, {\[Phi], 0, 2\[Pi]}  
 ];

gr00=DensityPlot[sol[\[Rho], \[Phi]], {\[Rho], 0.5, 10}, {\[Phi], 0, 2\[Pi] }
, Mesh -> All
, ColorFunction -> "Rainbow"
, PlotLegends -> Automatic
];  

Show[gr00 /. GraphicsComplex[array1_, rest___] :>  
                  GraphicsComplex[(#[[1]] {Cos[#[[2]]],Sin[#[[2]]]})& /@ array1, rest],
                  PlotRange -> 10 {{-1,1},{-1,1}}
                  ]   

enter image description here

Mesh from NDSolveValue

EDIT

THERE WAS A ERROR IN MY ANSWER JUST AFTER THIS EDIT

The real mesh used by NDSolve[...] is given by :

Show[sol["ElementMesh"]["Wireframe"] /. GraphicsComplex[array1_, rest___] :>  
                  GraphicsComplex[(#[[1]] {Cos[#[[2]]],Sin[#[[2]]]})& /@ array1, rest],
                  PlotRange -> 10 {{-1,1},{-1,1}}
                  ]    

enter image description here

END EDIT

sol = NDSolveValue[
  { \[Rho]^2 D[ u[\[Rho], \[Phi]], \[Rho], \[Rho]] + \[Rho] D[ u[\[Rho], \[Phi]], \[Rho]] +  D[u[\[Rho], \[Phi]], \[Phi], \[Phi]] == 0
  , DirichletCondition[u[\[Rho], \[Phi]] == 1000., \[Rho] == .5 && 0 <= \[Phi] <= 2 \[Pi]]
  , DirichletCondition[ u[\[Rho], \[Phi]] == 0., \[Rho] == 10 && 0 <= \[Phi] <= 2 \[Pi]]
  }
, u, {\[Rho], 0.5, 10}, {\[Phi], 0, 2\[Pi]}  
 ];

gr01=DensityPlot[sol[\[Rho], \[Phi]], {\[Rho],\[Phi]} \[Element] sol["ElementMesh"]
, Mesh -> All
, ColorFunction -> "Rainbow"
, PlotLegends -> Automatic
];

Show[gr01 /. GraphicsComplex[array1_, rest___] :>  
                  GraphicsComplex[(#[[1]] {Cos[#[[2]]],Sin[#[[2]]]})& /@ array1, rest],
                  PlotRange -> 10 {{-1,1},{-1,1}}
                  ]  

enter image description here

WARNING

The question is restricted to the particular case where the solution doesn't depend of Phi. All derivative related to phi are 0.

Let see what's happening if the diffusion of heat is not radial, for example suppose that the heating is only done in the upper-right part of the inner boundary. Here is the code :

sol = NDSolveValue[
  { \[Rho]^2 D[ u[\[Rho], \[Phi]], \[Rho], \[Rho]] + \[Rho] D[ u[\[Rho], \[Phi]], \[Rho]] +  D[u[\[Rho], \[Phi]], \[Phi], \[Phi]] == 0
  , DirichletCondition[u[\[Rho], \[Phi]] == If[0 <\[Phi]< Pi/2,1000.,200.], \[Rho] == 2 && 0 <= \[Phi] <= 2 \[Pi]]
  , DirichletCondition[ u[\[Rho], \[Phi]] == 0, \[Rho] == 10 && 0 <= \[Phi] <= 2 \[Pi]]
  }
, u, {\[Rho], 2, 10}, {\[Phi], 0, 2\[Pi]}  
 ];

gr02=DensityPlot[sol[\[Rho], \[Phi]], {\[Rho], 2, 10}, {\[Phi], 0, 2\[Pi] }
, Mesh -> All
, ColorFunction -> "Rainbow"
, PlotLegends -> Automatic
,PlotRange -> All
];    

Show[gr02 /. GraphicsComplex[array1_, rest___] :>  
                  GraphicsComplex[(#[[1]] {Cos[#[[2]]],Sin[#[[2]]]})& /@ array1, rest],
                  PlotRange -> 10 {{-1,1},{-1,1}}
                  ]

enter image description here

The heat doesn't diffuse through the line phi=0. This is because the domain has 4 boundaries : R=0.5 , R=10 , but also phi=0 , phi=2 Pi. The boundary conditions are not specified on the boundaries phi=0 , phi=2 Pi. In this case, NDSolve takes the default boudary condition which is Neumann=0 (ie Heat flux=0).The problem can be solved on the most recent versions of mathematica, where periodic Neumann conditions are introduced. (stroken because probably false, not sure)

$\endgroup$
  • $\begingroup$ Maybe it's possible to contruct a mesh with only 2 boundaries (the two circles) . I have not tried. Or maybe "sewing" 2 boundaries of a 4 boundaries domain ? $\endgroup$ – andre314 Nov 19 '17 at 17:58
  • $\begingroup$ maybe this brings clarity into the picture, I want to calculate this kind of configuration: link @ander $\endgroup$ – Alex Nov 19 '17 at 19:33
  • $\begingroup$ Note : The equations of heat are exactly the same than the equations of potential. Just replace °C by Volts. $\endgroup$ – andre314 Nov 19 '17 at 19:58
  • $\begingroup$ seconde case 2) there is "nothing" outside well, I'll look at this virtual charge method $\endgroup$ – Alex Nov 19 '17 at 20:17
2
$\begingroup$

In this configuration, there is no problem with stitching the borders

h = ImplicitRegion[0.25 <= x^2 + y^2 <= 100, {x, y}];


sol = NDSolveValue[{D[u[x, y], x, x] + D[u[x, y], y, y] == 0,    

DirichletCondition[u[x, y] == 10.,x^2 + y^2 == 0.25 &&  N[ArcTan[x, y] <= 0]],   

DirichletCondition[u[x, y] == 0., x^2 + y^2 == 100 && N[ArcTan[x, y] >=  0]]},u, {x, y} ∈ h]

DensityPlot[sol[x, y], {x, -10, 10}, {y, -10, 10}, Mesh -> None, ColorFunction -> "Rainbow", PlotLegends -> Automatic]

Output:

enter image description here

$\endgroup$
  • $\begingroup$ great answer!but why you use N[ArcTan[x, y] <= 0] and N[ArcTan[x, y] >= 0] ,and what about don't use it? $\endgroup$ – dcydhb Aug 9 '18 at 9:38
  • $\begingroup$ And theresult is not symmetrical. $\endgroup$ – dcydhb Aug 9 '18 at 9:56
  • $\begingroup$ I wanted to solve the problem of a cylindrical "capacitor" which has a cut in the wall. Therefore, it is understandable that the result will not be symmetrical @dcydhb $\endgroup$ – Alex Aug 11 '18 at 18:18
  • $\begingroup$ and what about my first question?and why your answer is different from the answer1 ,which is right $\endgroup$ – dcydhb Aug 13 '18 at 3:55

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.