0
$\begingroup$

I want to solve the following differential equation. It is about modeling the vibration of a linearized elastic rod.

$$\frac{\partial^2} {\partial x^2} ( E I \frac{\partial^2 w} {\partial x^2}) + \rho S \frac{\partial^2 w} {\partial t^2} =0$$

In the above formula, $E I = 1, \rho S = 1$, $w(x,t)$ is a binary function of $x$ and $t$.

The boundary and initial conditions are as follows:

$$w(x,t) \Big| _{t=0}=\frac{x^2} {6} (3 - x)$$ $$\frac{\partial w} {\partial t}\Big| _{x=0}=0$$

$$\frac{\partial^2 w} {\partial t^2}\Big| _{x=1}=0$$ $$\frac{\partial^3 w} {\partial t^3}\Big| _{x=1}=0$$

I wrote the following code according to the above conditions:

ClearAll["Global`*"]
tau = 10;
L = 1;
Elastic = 1;
Imoment = 1;
ρ = 1;
S = 1;
sol = NDSolveValue[{D[Elastic*Imoment*D[w[x, t], {x, 2}], {x, 2}] + 
     S*ρ*D[w[x, t], {t, 2}] == 0, w[x, 0] == x^2/6 (3 - x), 
   D[w[0, t], {t, 1}] == 0, 
   D[w[L, t], {t, 2}] == 0 D[w[L, t], {t, 3}] == 0}, 
  w[x, t], {x, 0, L}, {t, 0, tau}, 
  Method -> {"MethodOfLines", 
    "DifferentiateBoundaryConditions" -> {True, "ScaleFactor" -> 100},
     "SpatialDiscretization" -> {"TensorProductGrid", 
      "MaxPoints" -> 100, "MinPoints" -> 100, 
      "DifferenceOrder" -> 2}}, MaxSteps -> 10^6]

But I can't get the numerical solution of $w(x,t)$, so I can't draw the vibration image of the first 10 seconds.

What can I do to solve this partial differential equation?

$\endgroup$
7
  • 4
    $\begingroup$ As I said already in a comment your other post, the boundary condition on the third derivative in time does not make much sense. Actually all three boundary conditions in time do not make sense to me. Also, you need to impose initial conditions for the time derivatives. $\endgroup$ – Henrik Schumacher Apr 6 '20 at 8:15
  • 2
    $\begingroup$ Also, for people new to this, you should give at least a bit of context (it is about modelling the vibration of a linearized elastic rod) and link to you previous post. $\endgroup$ – Henrik Schumacher Apr 6 '20 at 8:16
  • 4
    $\begingroup$ hi. It is simple. You have PDE with order 4 in space and order 2 in time. Hence you need 4 boundary conditions and 2 initial conditions. You list above only 4 "conditions". You need a total of "6 conditions". two for time and 4 for space. You can not have a "condition" with equal or higher derivative than the order of derivative in the PDE. So for time, IC can not have higher than first order derivative, since the time in the PDE has order 2 derivative. Similar for space. $\endgroup$ – Nasser Apr 6 '20 at 8:43
  • 2
    $\begingroup$ What I said is: The partial differential equation that you wrote does not make sense. As we do not know which problem you actually want to solve, we cannot tell how to modify the boundary conditions. As for the initial conditions: If you mean the cantilever to be at rest in the beginning, then use D[w[x, t], {t, 1}] == 0 /. t -> 0 is what you should. $\endgroup$ – Henrik Schumacher Apr 6 '20 at 8:46
  • 1
    $\begingroup$ I guess from your other post that you want the left and of the cantilever to be clamped. That would be {w[0, t] == 0, D[w[x, t], {x, 1}] == 0 /. x -> 0}. I also guess that the right end should be free. I am not 100% sure, but I think the correct boundary conditions for that would be {D[w[x, t], {x, 2}] == 0 /. x -> L, D[w[x, t], {x, 3}] == 0 /. x -> L}. Other users will know that better than me. $\endgroup$ – Henrik Schumacher Apr 6 '20 at 8:47
2
$\begingroup$

To extend an answer by Henrik Schumacher and my answer here we conseder numericaly finished code

tau = 3;
L = 1;
Elastic = 1;
Imoment = 1;
\[Rho] = 1;
S = 1; A = 1/10; Dynamic["time: " <> ToString[CForm[currentTime]]]
AbsoluteTiming[
 sol = NDSolveValue[{S \[Rho] D[w[x, t], {t, 2}] + 
       D[Elastic Imoment D[w[x, t], {x, 2}], {x, 2}] == 0, 
     w[x, t] == A x^2/6 (3 - x) /. t -> 0, 
     D[w[x, t], {t, 1}] == 0 /. t -> 0, w[x, t] == 0 /. x -> 0, 
     D[w[x, t], {x, 1}] == 0 /. x -> 0, 
     D[w[x, t], {x, 2}] == 0 /. x -> L, 
     D[w[x, t], {x, 3}] == If[t <= 10^-4, -A, -A Exp[-10 t]] /. 
      x -> L}, w, {x, 0, L}, {t, 0, tau}, 
    Method -> {"MethodOfLines", 
      "DifferentiateBoundaryConditions" -> False, 
      "SpatialDiscretization" -> {"TensorProductGrid", 
        "MaxPoints" -> 80, "MinPoints" -> 80, 
        "DifferenceOrder" -> 4}}, 
    EvaluationMonitor :> (currentTime = t;)];]

Plot3D[sol[x, t], {x, 0, L}, {t, 0, tau}, ColorFunction -> "Rainbow", 
 AxesLabel -> Automatic, Mesh -> None, PlotPoints -> 50]

Figure 1

$\endgroup$
5
  • $\begingroup$ But I ran this code in Mathematica version 12.1, but I get a lot of error information, and I can't get the image shown in your answer. $\endgroup$ – A little mouse on the pampas Apr 7 '20 at 9:56
  • 1
    $\begingroup$ This code is debugged for version "12.0.0 for Microsoft Windows (64-bit) (April 6, 2019)". What messages did you receive in version 12.1? $\endgroup$ – Alex Trounev Apr 7 '20 at 10:38
  • $\begingroup$ I received the following information in v12.1: NDSolveValue::deqn: Equation or list of equations expected instead of If[t<=1/10000,-A,-A exp(-10 t)] in the first argument {w^(0,2)(x,t)+w^(4,0)(x,t)==0,w(x,0)==1/60 (3-x) x^2,w^(0,1)(x,0)==0,w(0,t)==0,w^(1,0)(0,t)==0,w^(2,0)(1,t)==0,If[t<=1/10000,-A,-A exp(-10 t)]}. NDSolveValue::dsvar: 0.0000204285714 cannot be used as a variable. General::stop: Further output of NDSolveValue::dsvar will be suppressed during this calculation. $\endgroup$ – A little mouse on the pampas Apr 7 '20 at 22:54
  • 1
    $\begingroup$ @PleaseCorrectGrammarMistakes Ah, sorry, there is a missing part of code. Now should be good. $\endgroup$ – Alex Trounev Apr 7 '20 at 23:14
  • $\begingroup$ Thank you very much for your great answer. $\endgroup$ – A little mouse on the pampas Apr 7 '20 at 23:41
1
$\begingroup$

This might work.

sol = NDSolveValue[{
   S \[Rho] D[w[x, t], {t, 2}] + D[Elastic Imoment D[w[x, t], {x, 2}], {x, 2}] == 0,
   w[x, t] == x^2/6 (3 - x) /. t -> 0,
   D[w[x, t], {t, 1}] == 0 /. t -> 0,
   w[x, t] == 0 /. x -> 0,
   D[w[x, t], {x, 1}] == 0 /. x -> 0,
   D[w[x, t], {x, 2}] == 0 /. x -> L,
   D[w[x, t], {x, 3}] == 0 /. x -> L
   },
  w,
  {x, 0, L}, {t, 0, 10},
  Method -> {"MethodOfLines", 
    "DifferentiateBoundaryConditions" -> {True, "ScaleFactor" -> 100},
     "SpatialDiscretization" -> {"TensorProductGrid", 
      "MaxPoints" -> 100, "MinPoints" -> 100, "DifferenceOrder" -> 2}
    },
  MaxSteps -> 10^6
  ]
$\endgroup$
1
  • $\begingroup$ Thank you very much for your answer, but the image of w at t = 0 is obviously not a parabola w=x^2/6 (3 - x). When t = 10, the amplitude of the cantilever beam is obviously too large, which is far greater than the initial amplitude L^2/6 (3 - L) at the first right end. $\endgroup$ – A little mouse on the pampas Apr 6 '20 at 9:06

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.