Message Boards Message Boards

Solving a specific initial value problem

Felix Toperczer, TU Wien
Posted 1 year ago

Hi, I am completely new to Mathematica and I really need some advice on solving a specific initial value problem numerically.

The differential equation is: θ''[x] - p^2 sin[2θ[x]] == 0 , x = {0,1}

for p = π/(2 Sqrt[2])

The initial values are θ'[0] == 0; θ[0] == {0, 2π}.

So I need to solve the equation for various initial angles θ[0] from 0 to 2π

Ι tried it with a For-Loop, but it is not working.
I would like to get a list of angles θ[1] after solving the equation for different initial angles θ[0]
I need to plot the angle θ[0] over the resulting angle θ[1].

Thanks for your help!

Attachments:
Question.nb
POSTED BY: Felix Toperczer
6 Replies
Sort By:
Attachment

Attachments:
Bildschirm­foto ...png
POSTED BY: Felix Toperczer

The following works fine for me, although with a different answer:

Clear[x];
p = Pi/(2 Sqrt[2]);
dgl = \[Theta]''[x] - p^2 Sin[2 \[Theta][x]] == 0;
sol = NDSolve[{dgl, \[Theta][0] == Pi/2, \[Theta]'[0] == 
    0}, \[Theta], {x, 0, 1}]
entry = \[Theta][1] /. sol[[1]]

I cannot reproduce the error you get, but it is similar to what we happens when the independent variable of a differential equation is given a value before NDSolve is called:

x = 0;
NDSolve[{f'[x] == 1, f[0] == 0}, f, {x, 0, 1}]

The variables f, x are not automatically protected from outside interference. However, in your case there must be something else going on.
POSTED BY: Gianluca Gorni

This is the mathematical pendulum equation for the oscillatory case.

Mathematica solves it directly, but the form of the solution is a bit weird, because the Jacobi amplitude is primarily defined for rotational modes with veloctiy :>0 at [Theta] = Pi,. 

    DSolveValue[\[Theta]''[t] - \[Pi]^2/8 Sin[2 \[Theta][t]] ==   0  , \[Theta][t], t] 

 JacobiAmplitude[(I Sqrt[\[Pi]^2-8 C[1] (t+C[2]))/(2 Sqrt[2]),(2 \[Pi]^2)/(\[Pi]^2-8 C[1])]

Its nearly impossible to include the start condtions. As the equation is time independent, one constant C[1] is simply the start time and has be set so to meet the start conditions together with .the second constant,determines the velocity at t=-C[2] , that is in the lower dead point.

      D[JacobiAmplitude[(I Sqrt[\[Pi]^2 - 8 C[1]] (t + C[2]))/(2 Sqrt[2]), ( 2 \[Pi]^2)/(\[Pi]^2 - 8 C[1])], t] /. {t + C[2] -> 0}

      \Omega  =  (I Sqrt[\[Pi]^2 - 8 C[1]])/(2 Sqrt[2])  ->   (Sqrt[8 C[1]]-\[Pi]^2 )/(2 Sqrt[2])

We know that the total energy is conserved

   [Theta]'[t]   [Theta]''[t] - \[Pi]^2/8     [Theta]'[t]  Sin[2 \[Theta][t]] ==   0   -->
    1/2 [Theta]'[t] ^2  +\[Pi]^2/16  (  1- Cos[2 \Theta]  = W = const

From this equation it is possible to solve for C[1] in terms of Theta

Because of the poor state of implementation of the elliptic function complex in Mathematica, most work has to be done by hand. The rather trival solution can be found in Pendel exakte Lösung

POSTED BY: Roland Franzius

In this case it is no pendulum equation, but an equilibrium condition for a rod in a crooked channel with an asymmetric cross-section. There is no oscillation going on.

POSTED BY: Felix Toperczer

Sorry, yes, there is a minus sign. But except for an x -> i t the algebra is all the same. One is in the uncharted area of the elliptic amplitude unforunately.
POSTED BY: Roland Franzius

You cannot use \[Theta]_entry as a variable name, it is a pattern. I would do it this way:

sol[t0_] := NDSolveValue[
   {\[Theta]''[x] - (\[Pi]/(2 Sqrt[2]))^2 Sin[2 \[Theta][x]] == 0,
    \[Theta][0] == t0, \[Theta]'[0] == 0}, \[Theta], {x, 0, 1}];
Table[{t0, sol[t0][1]},
 {t0, 0, 2 Pi, 0.125 Pi}]
POSTED BY: Gianluca Gorni
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Tag limit exceeded
Note: Only the first five people you tag will receive an email notification; the other tagged names will appear as links to their profiles.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract