Message Boards

WOLFRAM COMMUNITY

3944 Views

4 Replies

3 Total Likes

View groups...

Follow this post

Share this post:

GROUPS:

Mathematics Calculus Equation Solving Graphics and Visualization Wolfram Language Numerical Computation

Use results of NDSolve for other computations?

Daniel Gotsmann

Posted 5 years ago

I have been working on a model for a pendulum(using a Lagrangian) and I am using NDSolve. The solver gives me outputs for alpha an l yet I can´t us these to recalculate x and y. In[2]:= ClearAll["Global`"] Definition of parameters In[3]:= R = 0 Subscript[l, 0] = 0.5 m = 1 M = 2.394 g = 9.81 Subscript[\[Alpha], 0] = 1/2Pi Expression of position(x,y) as a function of [Alpha] In[9]:= x[\[Alpha]_, l_, t_] := -Cos[\[Alpha][t]]R + Sin[\[Alpha][t]](l[t] - \[Alpha][t]R) y[\[Alpha]_, l_, t_] := Sin[\[Alpha][t]]R + Cos[\[Alpha][t]](l[t] - \[Alpha][t]R) Calculation of the components and the norm of the velocity Computation of kinetic energy, potential energy and Lagrangian In[11]:= dx[\[Alpha]_, l_, t_] := \!\( \SubscriptBox[\(\[PartialD]\), \(t\)]\((x[\[Alpha], l, t])\)\) dy[\[Alpha]_, l_, t_] := \!\( \SubscriptBox[\(\[PartialD]\), \(t\)]\((y[\[Alpha], l, t])\)\) v[\[Alpha]_, l_, t_] := Sqrt[dx[\[Alpha], l, t]^2 + dy[\[Alpha], l, t]^2] // FullSimplify; Ekin[\[Alpha]_, l_, t_] := 1/2mv[\[Alpha], l, t]^2 + 1/2M(l'[t])^2 Epot[\[Alpha]_, l_, t_] := mgy[\[Alpha], l, t] + M*g (l[t] - Subscript[l, 0]) L1[\[Alpha]_, l_, t_] := Ekin[\[Alpha], l, t] - Epot[\[Alpha], l, t] Euler-Lagrange equation In[17]:= eq1[\[Alpha]_, l_, t_] = D[D[L1[\[Alpha], l, t], \[Alpha]'[t]], t] - D[L1[\[Alpha], l, t], \[Alpha][t]] // FullSimplify; eq2[\[Alpha]_, l_, t_] = D[D[L1[\[Alpha], l, t], l'[t]], t] - D[L1[\[Alpha], l, t], l[t]] // FullSimplify; eq = {eq1[\[Alpha], l, t] == 0, eq2[\[Alpha], l, t] == 0} Out[19]= {l[t] (-9.81 Sin[\[Alpha][t]] + 2 Derivative[1][l][t] Derivative[1][\[Alpha]][t] + l[t] (\[Alpha]^\[Prime]\[Prime])[t]) == 0, 23.4851 + 9.81 Cos[\[Alpha][t]] - l[t] Derivative[1][\[Alpha]][t]^2 + 3.394 (l^\[Prime]\[Prime])[t] == 0} In[20]:= Solution of equation of motion [AliasDelimiter]with given initial conditions In[22]:= sol = First@ NDSolve[{eq, \[Alpha][0] == Subscript[\[Alpha], 0], \[Alpha]'[0] == 0, l[0] == Subscript[l, 0], l'[0] == 0}, {\[Alpha][t], l[t]}, {t, 0, 20}] Plot In[23]:= Plot[\[Alpha][t] /. sol, {t, 0, 20}] In[24]:= Plot[l[t] /. sol, {t, 0, 20}] In[29]:= ParametricPlot[{x[\[Alpha], t]/.sol, y[\[Alpha], t]/.sol}, {t, 0, 20}, PlotRange -> All]

I have been working on a model for a pendulum(using a Lagrangian) and I am using NDSolve. The solver gives me outputs for alpha an l yet I can´t us these to recalculate x and y.

In[2]:= ClearAll["Global`*"]

Definition of parameters

In[3]:= R = 0
Subscript[l, 0] = 0.5
m = 1
M = 2.394
g = 9.81
Subscript[\[Alpha], 0] = 1/2*Pi

Expression of position(x,y) as a function of [Alpha]

In[9]:= x[\[Alpha]_, l_, t_] := -Cos[\[Alpha][t]]*R + 
  Sin[\[Alpha][t]]*(l[t] - \[Alpha][t]*R)
y[\[Alpha]_, l_, t_] := 
 Sin[\[Alpha][t]]*R + Cos[\[Alpha][t]]*(l[t] - \[Alpha][t]*R)

Calculation of the components and the norm of the velocity Computation of kinetic energy, potential energy and Lagrangian

In[11]:= dx[\[Alpha]_, l_, t_] := \!\(
\*SubscriptBox[\(\[PartialD]\), \(t\)]\((x[\[Alpha], l, t])\)\)
dy[\[Alpha]_, l_, t_] := \!\(
\*SubscriptBox[\(\[PartialD]\), \(t\)]\((y[\[Alpha], l, t])\)\)
v[\[Alpha]_, l_, t_] := 
  Sqrt[dx[\[Alpha], l, t]^2 + dy[\[Alpha], l, t]^2] // FullSimplify;
Ekin[\[Alpha]_, l_, t_] := 1/2*m*v[\[Alpha], l, t]^2 + 1/2*M*(l'[t])^2
Epot[\[Alpha]_, l_, t_] := m*g*y[\[Alpha], l, t] + M*g (l[t] - Subscript[l, 0])
L1[\[Alpha]_, l_, t_] := Ekin[\[Alpha], l, t] - Epot[\[Alpha], l, t]

Euler-Lagrange equation

In[17]:= eq1[\[Alpha]_, l_, t_] = 
  D[D[L1[\[Alpha], l, t], \[Alpha]'[t]], t] - 
    D[L1[\[Alpha], l, t], \[Alpha][t]] // FullSimplify;
eq2[\[Alpha]_, l_, t_] = 
  D[D[L1[\[Alpha], l, t], l'[t]], t] - D[L1[\[Alpha], l, t], l[t]] // 
   FullSimplify;
eq = {eq1[\[Alpha], l, t] == 0, eq2[\[Alpha], l, t] == 0}

Out[19]= {l[t] (-9.81 Sin[\[Alpha][t]] + 
     2 Derivative[1][l][t] Derivative[1][\[Alpha]][t] + 
     l[t] (\[Alpha]^\[Prime]\[Prime])[t]) == 0, 
 23.4851 + 9.81 Cos[\[Alpha][t]] - l[t] Derivative[1][\[Alpha]][t]^2 + 
   3.394 (l^\[Prime]\[Prime])[t] == 0}

In[20]:=

Solution of equation of motion [AliasDelimiter]with given initial conditions

In[22]:= sol = First@
  NDSolve[{eq, \[Alpha][0] == Subscript[\[Alpha], 0], \[Alpha]'[0] == 0, 
    l[0] == Subscript[l, 0], l'[0] == 0}, {\[Alpha][t], l[t]}, {t, 0, 20}]

Plot

In[23]:= Plot[\[Alpha][t] /. sol, {t, 0, 20}]

In[24]:= Plot[l[t] /. sol, {t, 0, 20}]

In[29]:= 
ParametricPlot[{x[\[Alpha], t]/.sol, y[\[Alpha], t]/.sol}, {t, 0, 20}, PlotRange -> All]

POSTED BY: Daniel Gotsmann

4 Replies

Sort By:

Moderation Team

Moderation Team, WOLFRAM

Posted 5 years ago

Welcome to Wolfram Community! Please make sure you know the rules: https://wolfr.am/READ-1ST The rules explain how to format your code properly. If you do not format code, it may become corrupted and useless to other members. Please EDIT your post and make sure code blocks start on a new paragraph and look framed and colored like this. int = Integrate[1/(x^3 - 1), x]; Map[Framed, int, Infinity]