I have a second order differential equation, in terms of a term that was previously obtained from numerically solving another equation.

Numerical solution of equation before 2nd order equation:

Rgas = Rationalize[8.314]; (* gas constant *)
acoef[tc_, pc_] := (27*Rgas^2*tc^2)/(64*pc);
bcoef[tc_, pc_] := (Rgas*tc)/(8*pc);
Clear[p, t];
latentK = 76.9*10^3*39*1.67*10^-27*6.022*10^23;
tcK = 2223; (* Kelvin *)
pcK = 16*10^6; (*Pa *)

eqn = (p - Rgas*t)*v^3 - bcoef[tcK, pcK]*p*v^2 + acoef[tcK, pcK]*v - acoef[tcK, pcK]*bcoef[tcK, pcK];
asolns = Solve[eqn == 0, v] ;
(eqn /. asolns)[[1]] ;
realroot = v /. asolns[[1]];

And now, my 2nd order differential equation, in terms of realroot, which is a function of p(t) an t:

NDSolve[p''[t] == latentK/(t*realroot), p[t], t]

But this is not working. I guess one of the problems is that realroot is not defined to be a function of p(t) and t...


Two issues need to be addressed. First, boundary conditions need to be specified for NDSolve. Second, p needs to be replaced by p[t] in realroot. In the absence of boundary conditions specified in the Question, I have provide some sample expressions.

s = NDSolve[{p''[t] == latentK/(t*realroot /. p -> p[t]), p[1] == 1, 
    p'[1] == 0}, p[t], {t, 1, 2}, Method -> "StiffnessSwitching"];
Plot[p[t] /. s, {t, 1, 2}]

enter image description here

  • $\begingroup$ @SuperCiocia In response to your comment, now deleted, Plot[Re[p[t]] /. First@s, {t, 273, 293}] will provide a curve, as will Plot[Im[p[t]] /. First@s, {t, 273, 293}], but the results suggests an error in your equations. $\endgroup$ – bbgodfrey Nov 23 '15 at 21:20

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.