MMA 12: Transient plane stress problem

Question

Based on the numerical example enter link description here which is proposed by @Hugh and @user21, then, I continue to solve transient plane stress problems (u[t, x, y]), howerver, time-dependent loading (i.e. DirichletCondition[v[t, x, y] == ss*t, x == L])does not work in MMA 12.

Description of the numerical problem:

Left side of the plate is fixed;

Right side of the plate is controlled by the time-dependent displacement in y direction, namely Dirichlet BCs.

Code

Needs["NDSolve`FEM`"];
L = 1;
h = 0.125;
(*Shear stress on beam*)
ss = 5;
reg = Rectangle[{0, -h}, {L, h}];
mesh = ToElementMesh[reg];
mesh["Wireframe"]
materialParameters = {Y -> 10^3, \[Nu] -> 33/100};


ps = {Inactive[
      Div][({{-(Y/(1 - \[Nu]^2)), 
         0}, {0, -((Y (1 - \[Nu]))/(2 (1 - \[Nu]^2)))}}.Inactive[
         Grad][u[t, x, y], {x, y}]), {x, y}] + 
    Inactive[
      Div][({{0, -((Y \[Nu])/(1 - \[Nu]^2))}, {-((Y (1 - \[Nu]))/(2 \
(1 - \[Nu]^2))), 0}}.Inactive[Grad][v[t, x, y], {x, y}]), { x, 
      y}], Inactive[
      Div][({{0, -((Y (1 - \[Nu]))/(2 (1 - \[Nu]^2)))}, {-((Y \
\[Nu])/(1 - \[Nu]^2)), 0}}.Inactive[Grad][u[t, x, y], {x, y}]), {
       x, y}] + 
    Inactive[
      Div][({{-((Y (1 - \[Nu]))/(2 (1 - \[Nu]^2))), 
         0}, {0, -(Y/(1 - \[Nu]^2))}}.Inactive[Grad][
        v[t, x, y], {x, y}]), {x, y}]};
 {uu, vv} = 
  NDSolveValue[{ps == {0, 0}, 
     DirichletCondition[v[t, x, y] == ss*t, x == L], 
     DirichletCondition[u[t, x, y] == 0, x == 0], 
     DirichletCondition[v[t, x, y] == 0, x == 0]} /. 
    materialParameters, {u, v}, {t, 0, 1}, {x, y} \[Element] mesh];

user21 · Accepted Answer · 2019-07-30T14:28:51.107

8

You'd need to make the PDE time dependent and give initial conditions:

{uu, vv} = 
  NDSolveValue[{{D[u[t, x, y], t], D[v[t, x, y], t]} + ps == {0, 0}, 
     u[0, x, y] == 0, v[0, x, y] == 0, 
     DirichletCondition[v[t, x, y] == ss*t, x == L], 
     DirichletCondition[u[t, x, y] == 0, x == 0], 
     DirichletCondition[v[t, x, y] == 0, x == 0]} /. 
    materialParameters, {u, v}, {t, 0, 1}, {x, y} \[Element] mesh];

Here is a way to visualize that:

graphics = Function[t,
    dmesh = 
     ElementMeshDeformation[mesh, 
      Transpose[{uu[t, ##], vv[t, ##]} & @@@ mesh["Coordinates"]], 
      "ScalingFactor" -> 0.01];
    Show[{
      mesh["Wireframe"["MeshElement" -> "BoundaryElements"]],
      dmesh[
       "Wireframe"[
        "ElementMeshDirective" -> 
         Directive[EdgeForm[Red], FaceForm[]]]]
      }, PlotRange -> {{0, 1.}, {0.2, -0.2}}]] /@ Range[0, 1, 0.1];
ListAnimate[graphics]

You need to make sure that the material parameters match (I have used a ScaleFactor < 1 to make this specific example work. Use a smaller force or a stronger material)

If you want second order time derivatives, you'd also need to specify derivatives of the initial condition:

Monitor[{uu, vv} = 
  NDSolveValue[{{D[u[t, x, y], {t, 2}], D[v[t, x, y], {t, 2}]} + 
       ps == {0, 0}, u[0, x, y] == 0, v[0, x, y] == 0, 
     Derivative[1, 0, 0][u][0, x, y] == 0, 
     Derivative[1, 0, 0][v][0, x, y] == 0, 
     DirichletCondition[v[t, x, y] == ss*t, x == L], 
     DirichletCondition[u[t, x, y] == 0, x == 0], 
     DirichletCondition[v[t, x, y] == 0, x == 0]} /. 
    materialParameters, {u, v}, {t, 0, 10^-1}, {x, y} \[Element] mesh,
    EvaluationMonitor :> (monitor = 
      Row[{"t = ", CForm[t]}])], monitor]

Also, see the section A Swinging and Dynamically Loaded Beam that talks about Rayleigh damping.

edited Jul 30 '19 at 14:28

answered Jul 30 '19 at 06:46

user21

39,710
8
110
167

thanks very much for your support! If I add D[u[t, x, y], t], D[v[t, x, y], t], then it is a dynamic problem? or? – ABCDEMMM Jul 30 '19 at 06:58
@ABCDEMMM, also add initial conditions. – user21 Jul 30 '19 at 07:01
one more question, how can we significantly speed up this transient simulation test in MMA 12 ? e.g. easy for us to apply GPUs? or MPI? – ABCDEMMM Jul 30 '19 at 07:17
1

@ABCDEMMM, no. And that is not a drawback. GPUs would need a lot of stuff to run on the GPU so putting some part on a GPU is not going to be helpful. And we are not alone (see here) For MPI. no that's not possible. In this case use usage of MPI is for very large problems but to the best of my knowledge it's not faster it's used to be able to solve these problems at all. The CPU time is not bad. You can reduce the time integration tolerance to speed things up a bit. – user21 Jul 30 '19 at 07:28
1

@user21 In the theory of elasticity for beams and plates, equations of the second order in time are used. – Alex Trounev Jul 30 '19 at 14:10
1

@AlexTrounev, see the updated answer. – user21 Jul 30 '19 at 14:29
@user21 can we control many cpus that we can use for this test (pardiso)? – ABCDEMMM Jul 30 '19 at 14:56
@ABCDEMMM, all cpus are used. – user21 Jul 30 '19 at 14:58
@user21 that means we donot need to set which kind of solver and how many cpus (the way like Abaqus....) – ABCDEMMM Jul 30 '19 at 15:00
@ABCDEMMM correct. – user21 Jul 30 '19 at 15:01
@user21 Pardiso is used automatically?? – ABCDEMMM Jul 30 '19 at 15:04
@ABCDEMMM for stationary FEM yes. For time dependent the linear solver is used differently and I am not sure if that uses the Pardiso version. – user21 Jul 30 '19 at 15:12
1

@user21 I understand that this is just a demonstration of the method. But physically, after removing the load, the beam should oscillate. And in your model, the beam just relaxes to its original state. – Alex Trounev Jul 30 '19 at 15:27
@AlexTrounev a physically correct oscillation is shown in the tutorial. – user21 Jul 30 '19 at 15:33
@AlexTrounev, also note that in the question the boundary is forced into it's place via a DirichletCondition - it's not a NeumannValue. – user21 Jul 30 '19 at 15:43

MMA 12: Transient plane stress problem

1 Answers1

Linked