Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Set of 2 nonlinear ODE in mathamtica 9

  1. Aug 4, 2014 #1
    i just signed up here so i hope this is the right place.

    i need to solve a set of 2 non-linear ordinary differencial equations.

    i tryed using NDSolve but it doesnt really work so im not sure whats wrong with my code.

    here is my code (copy paste):

    c = 0.1;
    B] = {2*m[t]*(m'[t])*
    Sin [\[Lambda][t]] + ((m[t])^2 )*(\[Lambda]'[t])*
    Cos[\[Lambda][t]] == -(c/(4*m[t]))*Tan[\[Lambda][t]/2],
    m'[t] == -(c/ (16*(m[t]^2)))*( Tan[\[Lambda][t]/2]^2),
    Sin[\[Lambda][0]] == 0.95, m[0] == 1};
    Subscript[sol, B] = NDSolve[Subscript[sys, B], {m}, t]

    Plot[m[t], {t, 0, 50}]

    as u can see from the code the equations are:

    (1) 2*m*m'*sinλ+m^2*λ'*cosλ = -(c/4m)*tan(λ/2)
    (2) m' = -(c/16m^2)*(tan(λ/2))^2

    where c=0.1 and also its known that:


    i wanna find m(t).
  2. jcsd
  3. Aug 4, 2014 #2
    Small changes to make NDSolve happy and get a plot.

    Code (Text):
    In[1]:= c = 0.1;
       sysB = {2*m[t]*m'[t]*Sin[\[Lambda][t]] + m[t]^2*\[Lambda]'[t]*Cos[\[Lambda][t]] == -c/(4*m[t])*Tan[\[Lambda][t]/2],
       m'[t] == -c/(16*m[t]^2)*Tan[\[Lambda][t]/2]^2, Sin[\[Lambda][0]] == 0.95, m[0] == 1};
       solB = m[t] /. NDSolve[sysB, m[t], {t, 0, 50}]

    Solve::ifun: Inverse functions are being used by Solve, so some solutions may not be found; use Reduce for complete solution information. >>

    Out[3]= {InterpolatingFunction[{{0. 18.8496}}, <>][t]}
    and the result is

    Code (Text):
    In[4]:= Plot[solB, {t, 0, 50}]

    Out[4]= ...PlotSnipped...
    Now you need to study this and make sure it is really correct and the result makes sense. Check the changes character by character until you can see exactly what all the changes were. If you want to you can then try reversing those changes, one at a time, and see if it continues to work or if it breaks. If it breaks then back off, try to understand what went wrong and see if you can find a way to have your original and have it still work.
    Last edited: Aug 5, 2014
  4. Aug 6, 2014 #3
    first of all, thank you for ur help.

    can u explain to me this line the code u wrote:
    solB = m[t] /. NDSolve[sysB, m[t], {t, 0, 50}]

    what exactly is being signed to solB? if i want to plot for m'[t] , how do i do it?
    and what exactly the /. does?

    when i wrote my code i simply used examples i saw for NDSolve and tried to write something similiar but it turns out i did it wrong cuz it didnt work. im trying to understand what was missing and what i did wrong, so ill know the next time i try writing something similiar.
  5. Aug 6, 2014 #4
    NDSolve, and DSolve and several other "solvers" in Mathematica, return "rules", sometimes multiple sets of rules, sometimes single rules. If you change that line to just NDSolve[sysB, m[t], {t, 0, 50}] then you can see the rule nested inside layers of { }. You can compare that with the value that is assigned to solB and see how this is changed.

    Rules are not a number, not like solving x+4==0 for x and getting -4. Understanding how rules work is an important part of learning to use Mathematica.

    /. is Mathematica shorthand for the ReplaceAll function in Mathematica and you can find position the cursor to the right of either of those, tap F1 and read all the documentation. In a sentence, what it is doing here is turning a rule into a definition of the function m[t].

    Correctly interpreting examples, particularly for differential equations and interpolating functions, can be difficult. When I've figured out how to correctly do this for a particular kind of problem that I will probably need to do again someday I save that along with carefully thought out notes to myself in a little Mathematica notebook that I keep on my desktop for future reference.
  6. Aug 7, 2014 #5
    im in abit of a problem now. im trying to write something more general for my problem, using functions instead of manually deriving the 2 equations and just solving them, can u help me get this working?

    basicly i have 2 first order ODE for 2 functions: j[tau] , m[tau] with 2 starting conditions, but i express the ODE with functions of j[tau] and m[tau] that i defined in the code.
    also i have 2 different models for one of those functions of j[tau] and m[tau] so theres 2 different sets of ODE.

    but somehow when i try running it i get alot of messages and i cant get any results. i dont think its because of the subscripts or greek letters because ive tested it on the code u gave earlier and they didnt do harm (as it should be).

    here is the code:

    Code (Text):
    c = 0.1;
    Subscript[\[Lambda], 0] = ArcSin[0.95];
    Subscript[m, 0] =
      2/(Sin[Subscript[\[Lambda], 0]]*Tan[Subscript[\[Lambda], 0]/2]);
    Subscript[j, 0] = 1;
    Subscript[J, 0] = 1; (*arbitrary*)
    Subscript[M, 0] = (Subscript[J, 0]/Sin[Subscript[\[Lambda], 0]])^(1/2);
    Subscript[\[CapitalOmega], H, 0] =
      Tan[Subscript[\[Lambda], 0]/2]/(2*Subscript[M, 0]);

    M[m_] := m*Subscript[J, 0]*Subscript[\[CapitalOmega], H, 0];
    J[j_] := j*Subscript[J, 0];
    \[Lambda][j_, m_] := ArcSin[J[j]/(M[m]^2)];

    Subscript[z, 1][j_, m_] :=
      1 + ((Cos[\[Lambda][j, m]])^2)*( (1 + Sin[\[Lambda][j, m]])^(1/
             3) + (1 - Sin[\[Lambda][j, m]])^(1/3) );
    Subscript[z, 2][j_,
       m_] := (3*(Sin[\[Lambda][j, m]])^2 +
         Subscript[z, 1][\[Lambda][j, m]]^2)^(1/2);
    SubPlus[z][j_, m_] :=
      3 + Subscript[z,
         m]] - ( (3 - Subscript[z, 1][\[Lambda][j, m]])*(3 +
            Subscript[z, 1][\[Lambda][j, m]] +
            2*Subscript[z, 2][\[Lambda][j, m]]) )^(1/2);

    Subscript[\[CapitalOmega], H][j_, m_] :=
      Tan[\[Lambda][j, m]/2]/(2*M[m]);
    Subscript[\[CapitalOmega], ISCO][j_, m_] :=
      1/(M[m]*(SubPlus[z][\[Lambda][j, m]]^(3/2) + Sin[\[Lambda][j, m]]));
    Subscript[\[CapitalOmega], T, A][j_, m_] :=
      Subscript[\[CapitalOmega], ISCO][j, m];
    Subscript[\[CapitalOmega], T, B][j_, m_] := (1/2)*
       Subscript[\[CapitalOmega], H][j, m];

    Subscript[\[Omega], H][j_, m_] :=
      Subscript[\[CapitalOmega], H][j, m]/Subscript[\[CapitalOmega], H, 0];
    Subscript[\[Omega], T, A][j_, m_] :=
      Subscript[\[CapitalOmega], T, A][j, m]/Subscript[\[CapitalOmega], H,
    Subscript[\[Omega], T, B][j_, m_] :=
      Subscript[\[CapitalOmega], T, B][j, m]/Subscript[\[CapitalOmega], H,

    (* the normilized equations of motion using \[Tau]=Subscript[\
    \[CapitalOmega], H,0]*t *)

      A] = {j'[\[Tau]] == -c*(Subscript[\[Omega], H][j[\[Tau]],
            m[\[Tau]]] - Subscript[\[Omega], T, A][j[\[Tau]], m[\[Tau]]]),
        m'[\[Tau]] ==
        Subscript[\[Omega], T, A][j[\[Tau]], m[\[Tau]]]*j'[\[Tau]],
       j[0] == Subscript[j, 0], m[0] == Subscript[m, 0]};
    Subscript[sol, A] =
     m[\[Tau]] /. NDSolve[Subscript[sys, A], m[\[Tau]], {\[Tau], 0, 50}]
    Plot[Subscript[sol, A], {\[Tau], 0, 80}]

      B] = {j'[\[Tau]] == -c*(Subscript[\[Omega], H][j[\[Tau]],
            m[\[Tau]]] - Subscript[\[Omega], T, B][j[\[Tau]], m[\[Tau]]]),
        m'[\[Tau]] ==
        Subscript[\[Omega], T, B][j[\[Tau]], m[\[Tau]]]*j'[\[Tau]],
       j[0] == Subscript[j, 0], m[0] == Subscript[m, 0]};
    Subscript[sol, B] =
     m[\[Tau]] /. NDSolve[Subscript[sys, B], m[\[Tau]], {\[Tau], 0, 50}]
    Plot[Subscript[sol, B], {\[Tau], 0, 80}]

    i know its alot to ask but im kinda hopeless here i dont know what to do...

    and btw, is there a way for me to write here the code so it will show the greek letters and the subscripts directly and not write the long tags for them instead, so it will be like i see it in mathmatica in my computer?
    Last edited: Aug 7, 2014
  7. Aug 7, 2014 #6
    What I would suggest is that you take the time to become an expert at understanding exactly how subscripts really work at all levels inside Mathematica. Hint: They are not variables, really, you might think they are variables, they might look like variables, in very simple cases they can pretend to act like variables, but they are not variables, really. I am told that they are kind of functions, but not in a way that I have ever really been able to understand. And in all the years I have never found tools that really show how they actually work at the deepest levels.

    Different issues:

    You define SubPlus[z] to accept two arguments, j and m, and then you use it SubPlus[z],[Lambda[j,m]]. I'm guessing one or the other of those is an error. Likewise for z2. Likewise for z1. If I guess what those might be, delete all your subscripts, try to correct all the mistakes I made trying to get rid of those and make a few more changes then I can at least get plots, but have no idea whether the result has anything to do with the real world. And the result is nothing like what you are going to want. <delete>

    Good luck. I hope you succeed. If you do then perhaps you can write a brilliantly clear paper that even I could understand and which explains exactly how subscripts really really work deep down inside Mathematica. Someone needs to do it. You might include the more general answer of how people can desktop publish their Mathematica input and still get it to work. If you can do that right then you could probably sell that.
  8. Aug 8, 2014 #7
    im not sure i have the time to become and expert in this.
    if instead of subscript i use something like Ω_T' its gonna be a proper variable? or i cant use greek letters either and i must write something like 'omega_T' ?

    i have corrected the mistake in the functions of z+ z1 z2 and still gotten many messages. some how it drew a plot for the model B set of equations, im not sure if its right, but im sure that theres alot of wrong things in the code somehow due to all those messages i get.

    i will try converting all those subscripts to just simple names, but i need to know - using greek letters is a problem?
    also - for the purpose of defining variables, mathmatica is case sensitive right? so i can use both OMEGA_T and omega_T for 2 different variables right?
  9. Aug 8, 2014 #8
    Underscore _ has a meaning buried deep in Mathematica and you are probably going to have nothing but trouble if you try that.

    That is probably a good guess.

    Anything that already has a defined meaning for Mathematica is going to give you grief. Naming something I or E or N or + or Sin or Solve or ... is going to cause problems. If you think you have a good name you might try looking that up in the help system. If it isn't there then you are less likely to have a problem, but Omega_T isn't going to be in the help system and still isn't going to work. Yes it is case sensitive.

    You might get yourself a good introduction to Mathematica book and learn some of the general knowledge about how to avoid problems. That will save you a lot of time in the future.

    If it is possible then try doing things one step or one function at a time. Write that. Test that. Try to figure out what the errors mean. After you get that correct then add the next layer of complexity. Banging in a whole program, pushing Enter and hoping it somehow will work is usually going to take more time than getting it right one simple step at a time.
  10. Aug 8, 2014 #9
    one thing im still unsure of is whether or not using greek letters without any subscripts is fine, and wont cause the problems subscripts cause?

    ive tried looking for a way to somehow still be able to work with those subscripts. according to here:

    i can use the line :


    if i understand correctly, it invokes the Notation package, and then for example i can use:

    Symbolize[Subscript[m, 0]]

    and after that perphaps Subscript[m, 0] will work properly as i need. are you familiar with this?
    maybe i can try doing it for all the subscripts in my code and see if it works..
  11. Aug 8, 2014 #10
    Try it and see what happens. Unfortunately it isn't like the stone age when we had documents which described exactly in detail what was and was not allowed in a language.

    I know. That is exactly what I said happens to some people several messages ago.

    Try reading SubscriptedVariables101.nb written by one of the Wolfram people after the umpteenth "why is my code broke?" post. That document doesn't explain what any of the problems with subscript really are. All that does is describe how to use one method to get around one kind of problem with subscript. But perhaps that will be enough for you.

    Good luck on your quest.
    Last edited: Aug 8, 2014
Share this great discussion with others via Reddit, Google+, Twitter, or Facebook