1
$\begingroup$

I want to solve the system of non-linear differential equations given below numerically.

$$y''(t)+500y'(t)+100y(t)=-33\cos(500t)-66\cos(1000t)$$ $$300x'(t)=1000y(t)+500y'(t)-35\tanh(50x'(t))$$

Notice that there is a derivative of $x$ inside $\tanh$. I tried the code below, but many errors appeared.

sol1 = 
  NDSolve[{y''[t] + 500 y'[t] + 100 y[t] == -33 Cos[500 t] - 66 Cos[1000 t],
    y[0] == 0, y'[0] == 0}, y, {t, 0, 30}];
sol2 = 
  NDSolve[{300 x'[t] == 
    1000 First[Evaluate[y[t] /. sol1]] + 500 First[Evaluate[y'[t] /. sol1]] - 
    35 Tanh[50 x'[t]], x[0] == 0}, x, {t, 0, 30}];
Plot[{Evaluate[x[t] /. sol1], Evaluate[y[t] /. sol2]}, {t, 0, 30}, PlotRange -> All]

How can I solve my system of differential equations using Mathematica?

Edit

I found the answer! thanks to Nasser. The solution for $x(t)$ is like below. I just need to adjust the axis

enter image description here

Improve this question
$\endgroup$

1 Answer 1

Reset to default
2
$\begingroup$

Is the solution supposed to be something like this? Hopefully I did not break too many things. You should know if the solution looks right or not.

Clear[y, x, t];
eqs = y''[t] + 500 y'[t] + 100 y[t] == -33 Cos[500 t] - 66 Cos[1000 t];
ic = {y[0] == 0, y'[0] == 0};
sol1 = First@NDSolve[{eqs, ic},{y[t], y'[t]}, {t, 0, 30},Method -> "StiffnessSwitching"];
eqs = 300 x'[t] == 1000 Evaluate[y[t] /. sol1] + 500 Evaluate[y'[t] /. sol1] - 
   35 Tanh[50 x'[t]]

Mathematica graphics

sol2 = NDSolve[{eqs, x[0] == 0}, x[t], {t, 0, 30},MaxSteps -> Infinity]

Mathematica graphics

Grid[{{Plot[Evaluate[y[t] /. sol1], {t, 0, 0.5}, PlotRange -> All]},
  {Plot[Evaluate[x[t] /. sol2], {t, 0, 30}, PlotRange -> All]}}]

Mathematica graphics

Improve this answer
$\endgroup$
4
  • $\begingroup$ thank you! i edited my question to include the answer. $\endgroup$
    – Dadan Ari Wibowo
    Nov 11, 2013 at 6:50
  • $\begingroup$ But if it is okay, I would like to know the reason for using StiffnessSwitching and MaxSteps->Infinity... $\endgroup$
    – Dadan Ari Wibowo
    Nov 11, 2013 at 6:56
  • $\begingroup$ "StiffnessSwitching", well, more information on this is here reference.wolfram.com/mathematica/tutorial/…. Physics wise, the natural frequency of your eqs is 10 rad/sec, while the forcing frequencies (RHS) are much higher (you have one at 1000 rad/sec). So it tells one that system will oscillate slowly in the early transient stage, then in steady state, when the load takes over, system will oscillate at much higher frequency (at the driving frequency). Hence you have a system with low and high vibrations in it. Hence this method to try.... $\endgroup$
    – Nasser
    Nov 11, 2013 at 7:10
  • $\begingroup$ ..."MaxSteps->Infinity" this one is easy. Mathematica gave an error saying Maximum number of 10000 steps reached at the point t == 1.5727575203, so you just increase the number of steps and hope for the best. $\endgroup$
    – Nasser
    Nov 11, 2013 at 7:11

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.