3
$\begingroup$

Writing the following code:

Rext = 16;
Rint = 3;
thickness = 1/10;

nholes = 8;
distcenters = 9;
Rholes = 2;

density[x_, y_, z_] := (x^2 + y^2 + z^2) 8 10^-3;

domain = "-thickness/2\[LessEqual]z\[LessEqual]thickness/2" <> 
         "&&(x-0)^2+(y-0)^2\[LessEqual]Rext^2" <> 
         "&&(x-0)^2+(y-0)^2\[GreaterEqual]Rint^2";
Table[domain = domain <> "&&(x-distcenters Cos[2(" <> ToString[i] <> 
         "-1)Pi/nholes])^2+(y-distcenters Sin[2(" <> ToString[i] <>
         "-1)Pi/nholes])^2\[GreaterEqual]Rholes^2", {i, 1, nholes}];

disk = ImplicitRegion[ToExpression[domain], {x, y, z}];
min = NMinValue[{density[x, y, z], {x, y, z} \[Element] disk}, {x, y, z}];
max = NMaxValue[{density[x, y, z], {x, y, z} \[Element] disk}, {x, y, z}];

Magnify[RegionPlot3D[disk,
              ColorFunction -> Function[{x, y, z},
              ColorData["SandyTerrain"][density[x, y, z]]],
              ColorFunctionScaling -> False,
              PlotLegends -> BarLegend[{"SandyTerrain", {min, max}}],
              PlotPoints -> 75], 1.5]

Iz = NIntegrate[(x^2 + y^2) density[x, y, z], {x, y, z} \[Element] disk]

I get:

enter image description here

13474.5

which is wonderfully as desired.

Now, assuming the elastic theory is valid and assuming that this disk is subject to a perpendicular load, I would like to calculate the stress state to which it is subject by plotting as above. For this purpose I tried the following approach:

NDSolveValue[{
 D[w[x, y, z], {x, 4}] + 2 D[w[x, y, z], {x, 2}, {y, 2}] + D[w[x, y, z], {y, 4}] == 1,
 DirichletCondition[w[x, y, z] == 0, x^2 + y^2 == Rint^2],          
 DirichletCondition[D[w[x, y, z], {x, 2}] == 0, x^2 + y^2 == Rint^2],        
 DirichletCondition[D[w[x, y, z], {y, 2}] == 0, x^2 + y^2 == Rint^2]}, 
 w, {x, y, z} \[Element] disk]

but unfortunately I get the following error:

NDSolveValue::femcmsd: The spatial derivative order of the PDE may not exceed two.

Of course before opening this topic for many hours I read other topics and the most similar is this. Could someone help me adapt it to my case? Thank you!

$\endgroup$
  • $\begingroup$ Currently, the FEM can not handle higher than second derivatives. Having said that, you could try to rewrite your equation as a system of two second order equations. $\endgroup$ – user21 Dec 28 '17 at 9:16

Your Answer

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

Browse other questions tagged or ask your own question.