0
$\begingroup$

I'm trying to solve some Lagrange differential equations of motion, where $L = T - U$.

How can I use the solution for each x1, x2, x3, ..., xn?

Clear["Global`*"]
Needs["VariationalMethods`"]
n = 10; 
t0 = 10;
 Subscript[x, 0][t_] := 0;
 Subscript[x, n + 1][t_] := 0;
 m = 1;
 v0 = 0.1; 
  Table[Subscript[k, j] = 1, {j, 0, n + 1}]; 
ic = Table[{Subscript[x, j][0] == 0, Derivative[1][Subscript[x, j]][0] ==RandomReal[{-v0, v0}]}, {j, 1, n}]; 

ue[x_, t_, k_, n_] := 
       (1/2)*Sum[Subscript[k, j]*(Subscript[x, j][t] - Subscript[x, j - 1][t])^2, 
         {j, 1, n + 1}]; 
    te[x_, t_, n_] := (1/2)*m*Sum[Derivative[1][Subscript[x, j]][t]^2, {j, 1, n}];


lg[x_, t_, k_, n_] := te[x, t, n] - ue[x, t, k, n];

eqm[j_] := D[lg[x, t, k, n], {Subscript[x, j][t], 1}] - 
    D[D[lg[x, t, k, n], {Derivative[1][Subscript[x, j]][t], 1}], {t, 1}];
$\endgroup$
2
  • $\begingroup$ 1) There was a question about frequencies, not exact solutions. 2) There are two types of boundary conditions with $k_0=0; 5$. $\endgroup$ Mar 10, 2019 at 18:41
  • 1
    $\begingroup$ Please look at the instructions here to copy readable code from MMA: How to copy code from MMA so it looks good on this site $\endgroup$
    – MarcoB
    Mar 10, 2019 at 20:58

1 Answer 1

5
$\begingroup$

Numerical solution

Clear["Global`*"]
Needs["VariationalMethods`"]
n = 10; t0 = 10; Subscript[x, 0][t_] := 0; Subscript[x, n + 1][t_] := 0; m = 1; v0 = 0.1; 
  Table[Subscript[k, j] = 1, {j, 0, n + 1}]; 
ic = Table[{Subscript[x, j][0] == 0, Derivative[1][Subscript[x, j]][0] == 
      RandomReal[{-v0, v0}]}, {j, 1, n}]; 
ue[x_, t_, k_, n_] := 
   (1/2)*Sum[Subscript[k, j]*(Subscript[x, j][t] - Subscript[x, j - 1][t])^2, 
     {j, 1, n + 1}]; 
te[x_, t_, n_] := (1/2)*m*Sum[Derivative[1][Subscript[x, j]][t]^2, {j, 1, n}]; 
lg[x_, t_, k_, n_] := te[x, t, n] - ue[x, t, k, n]; 
eqm[j_] := D[lg[x, t, k, n], {Subscript[x, j][t], 1}] - 
    D[D[lg[x, t, k, n], {Derivative[1][Subscript[x, j]][t], 1}], {t, 1}]; 
sol = NDSolveValue[Flatten[{Table[eqm[i] == 0, {i, 1, n}], ic}], 
   Table[Subscript[x, i][t], {i, 1, n}], {t, 0, t0}]
Plot[sol, {t, 0, t0}, PlotLegends -> Automatic]

FIG1

Calculation of frequencies for a system of 100 particles with springs of different stiffness

    ks1[i_] := 3 + 2*(-1)^i
ks2[i_] := 3 - 2*(-1)^i
m = 1;
omega[k_, n_] := 
 Block[{k0 = k, P = n}, 
  eq1 = -m*(x[1]'')[t] + k0*(x[0][t] - x[1][t]) - 
    5 (x[1][t] - x[2][t]);
  eqP = -m*(x[P]'')[t] + k0*(x[P + 1] - x[P][t]) - 
    5 (-x[P - 1][t] + x[P][t]);
  eqn = Table[-m*x[i]''[t] - ks1[i]*(x[i][t] - x[i - 1][t]) - 
     ks2[i]*(x[i][t] - x[i + 1][t]), {i, 2, P - 1}];
  Eq1 = Join[{eq1, eqn, eqP}] /. {x[0][t] -> 0, x[P + 1][t] -> 0};
  op = Table[x[i_][t_] := A[i]*Exp[I*om*t], {i, 1, P}];
  Eq2 = Flatten[Eq1 /. t -> 0];
  Eq = Table[Eq2[[i]] == 0, {i, 1, P}];
  matr = CoefficientArrays[Eq, Table[A[i], {i, 1, P}]];
  M = Last[Normal[matr]];
  f = Det[M];
  s = Solve[f == 0, om];
  lst = Table[N[om /. s[[i]], 10], {i, 1, Length[s]}];
  Drop[Sort[lst], P]]

An example of using the function omega[]

lst0 = omega[0, 100];
lst5 = omega[5, 100];
{ListPlot[Abs[lst0 - lst5], PlotRange -> All],
 ListPlot[{lst0, lst5}, 
  PlotLegends -> {Row[{"k0 = ", 0}], Row[{"k0 = ", 5}]}], 
 ListPlot[lst0, PlotLabel -> Row[{"k0 = ", 0}]], 
 ListPlot[lst5, PlotLabel -> Row[{"k0 = ", 5}]]}

fig2

$\endgroup$
1
  • $\begingroup$ @Alrubaie I fixed the code. $\endgroup$ Mar 12, 2019 at 21:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.