# Finite Element Analysis Boundary Condition problem

I'm trying to learn how to use Mathematica to perform FEA. In this case i'm trying to perform a compressive test on a disk. I've started with meshing the component.

Reg1 = Annulus[{0, 0}, {14, 25}];
Mesh1 = ToElementMesh[Reg1, MeshQualityGoal -> 0];


The compressive force is 6.41 Mpa

q = 6.41;


The operator for the plane stress analysis is provided by Mathematica and it's:

op = {Inactive[
Div][({{0, -((Y \[Nu])/(1 - \[Nu]^2))}, {-((Y (1 - \[Nu]))/(
v[x, y], {x, y}]), {x, y}] +
Inactive[
Div][({{-(Y/(1 - \[Nu]^2)),
0}, {0, -((Y (1 - \[Nu]))/(2 (1 - \[Nu]^2)))}}.Inactive[
Grad][u[x, y], {x, y}]), {x, y}],
Inactive[
Div][({{0, -((Y (1 - \[Nu]))/(2 (1 - \[Nu]^2)))}, {-((
u[x, y], {x, y}]), {x, y}] +
Inactive[
Div][({{-((Y (1 - \[Nu]))/(2 (1 - \[Nu]^2))),
v[x, y], {x, y}]), {x, y}]} /. {Y -> 3416, \[Nu] -> 0.33};


And now I have problems with the definition of the boundary conditions:

bc = {0, NeumannValue[-q, x == 0 && y==25],NeumannValue[q, x == 0 && y==-25]};

{uif, vif} =
NDSolveValue[{op == {0, NeumannValue[-q, x == 0]}}, {u, v}, {14^2 < x^2 + y^2 < 25^2}]

Mesh1 = uif["ElementMesh"];
ElementMeshDeformation[Mesh1, {uif, vif}]["Wireframe"]


Obviously this is wrong but I can't understand why and how to set my boundary conditions. **The image above is referring to another type of mesh.

• you declare bc but as far as I can see you are using other boundary conditions when solving the problem Dec 14, 2016 at 0:14

Seems like one needs to have an extended area for the boundary conditions to work properly, so I used z<24.5 instead of z==-25 as you had it. Changed a bit the code but I hope you get the idea:

<< NDSolveFEM
ClearAll[x, y, v, u, mesh1, uif, vif]
Reg1 = Annulus[{0, 0}, {14, 25}];

mesh = ToElementMesh[Reg1, "MaxCellMeasure" -> 0.5]
planeStress = {Inactive[
Div][{{0, -((Y*\[Nu])/(1 - \[Nu]^2))}, {-(Y*(1 - \[Nu]))/(2*(1 \
- \[Nu]^2)), 0}}.Inactive[Grad][v[y, z], {y, z}], {y, z}] +
Inactive[
Div][{{-(Y/(1 - \[Nu]^2)),
0}, {0, -(Y*(1 - \[Nu]))/(2*(1 - \[Nu]^2))}}.Inactive[Grad][
u[y, z], {y, z}], {y, z}],
Inactive[
Div][{{0, -(Y*(1 - \[Nu]))/(2*(1 - \[Nu]^2))}, {-((Y*\[Nu])/(1 \
- \[Nu]^2)), 0}}.Inactive[Grad][u[y, z], {y, z}], {y, z}] +
Inactive[
Div][{{-(Y*(1 - \[Nu]))/(2*(1 - \[Nu]^2)),
v[y, z], {y, z}], {y, z}]} /. {Y -> 10^3, \[Nu] -> 33/100};
{uif, vif} =
NDSolveValue[{planeStress == {0, -NeumannValue[60, z >= 24.5]},
DirichletCondition[u[y, z] == 0, z < -24.5],
DirichletCondition[v[y, z] == 0, z < -24.5]}, {u,
v}, {y, z} \[Element] mesh];
dmesh = ElementMeshDeformation[mesh, {uif, vif}, "ScalingFactor" -> 1];
Show[{mesh["Wireframe"],
dmesh["Wireframe"[
"ElementMeshDirective" -> Directive[EdgeForm[Red], FaceForm[]]]]}]


• Thank you so much. Now I see what the problem was! This is really helpful:)
– user44593
Dec 14, 2016 at 21:28