7

This problem is axially symmetric. I want to solve using the Finite Element method and taking advantage of the symmetry.

(Note: eleField3D needs to be called with Cartesian coordinates)

Needs["NDSolve`FEM`"]
Clear["Global`*"]

cylinOffset = 0.0125;(*meter*)
cylinRadius = 0.0075;(*meter*)

shellRadius = 0.0233;(*meter*)
scaleLen = 2.10;

cylinVoltage = 100;

shell = ImplicitRegion[0 <= r <= shellRadius && ζ^2 <= (scaleLen shellRadius)^2, {r,ζ}];
halfTorus = ImplicitRegion[(r - cylinOffset)^2 + ζ^2 <= cylinRadius^2 && r <= cylinOffset, {r, ζ}];

reg = ToElementMesh[RegionDifference[shell, halfTorus], AccuracyGoal -> 9, PrecisionGoal -> 9, MaxCellMeasure -> 1*^-8];
RegionPlot[reg, AspectRatio -> scaleLen]

Please check that my 2D (r,z) Cylindrical Laplacian is correct.

solution =  NDSolveValue[{(1/r) D[r D[u[r, ζ], r], r] + D[u[r, ζ], {ζ, 2}] == 0,
   DirichletCondition[u[r, ζ] == 0, r == shellRadius],
   DirichletCondition[u[r, ζ] == 0, ζ^2 == (scaleLen shellRadius)^2],
   DirichletCondition[u[r, ζ] == cylinVoltage, r == cylinOffset && ζ^2 <= cylinRadius^2],
   DirichletCondition[u[r, ζ] == cylinVoltage, (r - cylinOffset)^2 + ζ^2 == cylinRadius^2 && r <= cylinOffset]},
 u, {r, ζ} ∈ reg, Method -> {"FiniteElement"}]

How do I revolve the solution to get the 3D field?

Proposed Answer:

eleFieldCylin[r_, ζ_] = -Grad[solution[r, ζ], {r, ζ}];
eleField3D[x_, y_, z_] = RotationTransform[ArcTan[x, y], {0, 0, 1}][
  Insert[eleFieldCylin[Sqrt[x^2 + y^2], z], 0, 2]];

Show[
 VectorPlot3D[eleField3D[x, y, z],
  {x, -shellRadius, shellRadius}, {y, -shellRadius, shellRadius},
  {z, -(scaleLen shellRadius), (scaleLen shellRadius)}],
 RegionPlot3D[
  DiscretizeRegion[
   ImplicitRegion[(cylinOffset - Sqrt[x^2 + y^2])^2 + z^2 <= 
      cylinRadius^2 && x^2 + y^2 <= cylinOffset^2, {x, y, z}]]]
 ]
xzczd
  • 65,995
  • 9
  • 163
  • 468
Young
  • 7,495
  • 1
  • 20
  • 45
  • 1
    I think you would more likely get an answer if you explain your equations more and in particular phrase a very specific question. For instance, the only question I can see is called "answer". – Felix Jan 31 '17 at 00:10
  • Your Laplacian looks correct to me. The only thing that I found is the double comma after Laplacian == 0. This seems to be a typo and could cause problems. – Felix Feb 02 '17 at 01:54
  • Haven't tested your code yet. (You know Wolfram Cloud is slow… ) But I think something is missing before [Insert……? Also, you added 2 ; in the ToElementMesh[…… line. BTW, you can use this function to copy all the code piece in a post to test it. – xzczd Feb 02 '17 at 03:11
  • If the shell is a rectangle, I'd use Rectangle to represent it, something like Rectangle[{0, -scaleLen shellRadius}, {shellRadius, scaleLen shellRadius}]. You could use NIntegrate[1, {x, y} \[Element] reg] to measure how well the region area is approximated. – user21 Feb 07 '17 at 06:58

1 Answers1

3

I think the proposed answer is correct now, here's just some suggestions.

First, you can deduce Laplace equation with the help of Laplacian to avoid possible mistakes:

Laplacian[u[r, ζ], {r, th, ζ}, "Cylindrical"] == 0
(* Derivative[0, 2][u][r, ζ] + Derivative[1, 0][u][r, ζ]/r + 
     Derivative[2, 0][u][r, ζ] == 0 *)

Second, eleField3D can be defined in a clearer and simpler way with TransformField:

eleFieldCylin[r_, ζ_] = -Grad[solution[r, ζ], {r, th, ζ}, "Cylindrical"];
eleField3D[x_, y_, z_] = 
  TransformedField["Cylindrical" -> "Cartesian", 
   eleFieldCylin[r, ζ], {r, th, ζ} -> {x, y, z}]

Third, it may be better to set a region function for VectorPlot3D, for example:

regionfunc = (Function[{x, y, z}, #1] & )[
       0 <= r <= shellRadius && ζ^2 <= (scaleLen*shellRadius)^2 && 
            ! ((r - cylinOffset)^2 + ζ^2 <= cylinRadius^2 && 
                r <= cylinOffset) /. {r -> Sqrt[x^2 + y^2], ζ -> z}];

VectorPlot3D[
  (Piecewise[{{#1, regionfunc[x, y, z]}}, 
          {0, 0, 0}] & )@eleField3D[x, y, z]//Evaluate, {x, -shellRadius, shellRadius}, 
  {y, -shellRadius, shellRadius}, {z, -(scaleLen shellRadius), scaleLen shellRadius}]

4th, I think PrecisionGoal and AccuracyGoal has no obvious effect here and can be taken away. (Have you read this post?: Is manual adjustment of AccuracyGoal and PrecisionGoal useless?)

Community
  • 1
xzczd
  • 65,995
  • 9
  • 163
  • 468
  • Thank you very much for your answer! I do think AccuracyGoal and PrecisionGoal in the context of ToElementMesh improve the boundary mesh quality based on some limited experimentation... @user21 would need to confirm – Young Feb 04 '17 at 15:24
  • @young Hope user21 will go by here. BTW (maybe you already know that) "@user21" won't remind him because he hasn't appeared in the comment of this answer for so long: http://meta.stackexchange.com/a/43020/284701 – xzczd Feb 05 '17 at 02:52