0
|
933 Views
|
7 Replies
|
0 Total Likes
View groups...
Share

# Plot the solution of a differential equation

Posted 4 months ago
 Consider: $$\frac{1}{(1-e)^2 w}=v-z$$ where $v=f(e,i,L,W,a,n)$, $z=g(e,i,L,W,a,w,n)$, and $L = h(e,w,a,b)$ are given functions. From the above equation, $e$ can be derived as a function of $w$ along with the other parameters, i.e., $e=e(w;\dots)$. Using the latter, we can consider the following differential equation: $$\frac{\partial e}{\partial w}=\frac{e}{w}$$ I would like find $e$ and $w$ that satisfies the above differential equation. Obviously, both $e$ and $w$ will come out as a function of $a$ along with the other parameters. Finally, I would like to create four separate plots for $e$, $w$, $\frac{e}{w}$, and $\frac{ae}{w}$ each of them against $a$ with the parameter values of $n=1$, $W=60$, $i=0.1$, and $b=0.6$. By closely referring to Michael E2's answer to this post, I tried to come up with the following code: Clear["Global*"]; n = 1; W = 60; L = (w/(b (a e)^b))^(1/(b - 1)); v = -(((-1 + e - i) (-i + (-1 + e + e i) L n) (1/((-1 + e) (-1 + e - i) w) + w/(-1 + e - i) + ((-1 + e) (1 - e L n) (-(1 + i)^2 (-1 + (1 + i)^((1 - e L n)/(L n - e L n)))^2 + a^2 i^2 (-1 + (1 + i)^(1 + 1/(L n - e L n)))^2 W^2))/(a i (1 + i) (1 - e + i) (1 - (1 + i)^((1 - e L n)/(L n - e L n))) (1 - (1 + i)^(1 + 1/(L n - e L n))) (i - (-1 + e + e i) L n) W)))/(i (1 + i + L n + e (-1 + (-2 + e - i) L n)))); z = -((a i (1 + i)^((1 + e)/(1 - e)) (-1 + (1 + i)^((1 - e L n)/(L n - e L n))) (-1 + (1 + i)^(1 + 1/(L n - e L n))) L n W + a (-1 + e) i (1 + i)^((1 + e)/(1 - e)) (-1 + (1 + i)^((1 - e L n)/(L n - e L n))) (-1 + (1 + i)^(1 + 1/(L n - e L n))) L n w^2 W + (-1 + e - i) (1 + i)^(-((2 e)/(-1 + e))) (-1 + e L n) w (1 + 2 i + i^2 - a^2 i^2 W^2 - 2 (1 + i)^(1 + 1/(L n - e L n)) ((1 + i)^(2 + 1/(-1 + e)) - a^2 i^2 W^2) + (1 + i)^(2 + 2/(L n - e L n)) ((1 + i)^((2 e)/(-1 + e)) -a^2 i^2 W^2)))/(a i^2 (1 + i) ((1 + i)^(e/(1 - e)) - (1 + i)^(1/(L n - e L n))) ((1 + i)^(e/(1 - e)) - (1 + i)^((1 + L n)/(L n - e L n))) (1 + i + L n + e (-1 + (-2 + e - i) L n)) w W)); myEQ = 1/((1 - e)^2 w) == v - z; Block[{n = 1, W = 60, i = 1/10, b = 6/10, a = 1/2}, {#, Dt[#, w]} &@myEQ] /. {e -> e[w]} /. {e'[w] -> e[w]/w} /. {e[w] -> e} // Simplify; icEQ = % /. Equal -> Subtract // Simplify; icSOL = NSolve[icEQ == 0 && 0 < e < 1 && 0 < w < 2, {e, w}] NDSolveValue[{ode = Solve[Block[{n = 1, W = 60, i = 1/10, b = 6/10}, D[{myEQ /. Equal -> Subtract, {e, w} . D[myEQ /. Equal -> Subtract, {{e, w}}]} /. {w -> w[a], e -> e[a]}, a] == 0], {e'[a], w'[a]}] /. Rule -> Equal, e[1/2] == (e /. First@icSOL), w[1/2] == (w /. First@icSOL)}, {e, w}, {a, $MachineEpsilon, 1}] ListLinePlot[{e\[FivePointedStar], w\[FivePointedStar]}, PlotLegends -> Block[{e\[FivePointedStar], w\[FivePointedStar]}, HoldForm /@ {e\[FivePointedStar], w\[FivePointedStar]}], PlotRange -> All]  The code runs forever. Also, I failed to come up with a code that creates the four separate plots. 7 Replies Sort By: Posted 4 months ago  From the following plot myTwoEQs = Block[{n = 1, W = 60, i = 1/10, b = 6/10, a = 1/2}, {#, Dt[#, w]} &@myEQ] /. {e -> e[w]} /. {e'[w] -> e[w]/w} /. {e[w] -> e}; ContourPlot[myTwoEQs, {e, 0, 1}, {w, 0, 2}, PlotPoints -> 200, ContourStyle -> {Directive[Blue, Thickness[.01], Opacity[.2]], Directive[Darker@Red, Thin]}] it would seem that the system of the two equations has infinitely many solutions, along a whole line (near the left side of the plot). If this is correct, then the command icEQ = myTwoEQs /. Equal -> Subtract; icSOL = NSolve[icEQ == 0 && 0 < e < 1 && 0 < w < 2, {e, w}] cannot possibly work. I would rather give a fixed value for w, for example w==1.A closer look suggests that the situation is complicated: Plot[Evaluate[ myTwoEQs[[1]] /. Equal -> Subtract /. w -> 1], {e, .0191, .0192}, PlotRange -> 10 {-1, 1}]  Posted 3 months ago  Thanks for the comment. In fact, the code is exactly the same as the code suggested in Michael E2's answer to the linked post except for the equation in myEQ. And that code worked just fine. This is just for your information. In any case, I was trying to following your comments but failed. May I ask if you could provide a full code to work? I would really appreciate that. Thank you so much for your help. Posted 3 months ago  Here is the full code Clear["Global*"]; n = 1; W = 60; i = 1/10; b = 6/10; a = 1/2; L = (w/(b (a e)^b))^(1/(b - 1)); v = -(((-1 + e - i) (-i + (-1 + e + e i) L n) (1/((-1 + e) (-1 + e - i) w) + w/(-1 + e - i) + ((-1 + e) (1 - e L n) (-(1 + i)^2 (-1 + (1 + i)^((1 - e L n)/(L n - e L n)))^2 + a^2 i^2 (-1 + (1 + i)^(1 + 1/(L n - e L n)))^2 W^2))/(a i (1 + i) (1 - e + i) (1 - (1 + i)^((1 - e L n)/(L n - e L n))) (1 - (1 + i)^(1 + 1/(L n - e L n))) (i - (-1 + e + e i) L n) W)))/(i (1 + i + L n + e (-1 + (-2 + e - i) L n)))); z = -((a i (1 + i)^((1 + e)/(1 - e)) (-1 + (1 + i)^((1 - e L n)/(L n - e L n))) (-1 + (1 + i)^(1 + 1/(L n - e L n))) L n W + a (-1 + e) i (1 + i)^((1 + e)/(1 - e)) (-1 + (1 + i)^((1 - e L n)/(L n - e L n))) (-1 + (1 + i)^(1 + 1/(L n - e L n))) L n w^2 W + (-1 + e - i) (1 + i)^(-((2 e)/(-1 + e))) (-1 + e L n) w (1 + 2 i + i^2 - a^2 i^2 W^2 - 2 (1 + i)^(1 + 1/(L n - e L n)) ((1 + i)^(2 + 1/(-1 + e)) - a^2 i^2 W^2) + (1 + i)^(2 + 2/(L n - e L n)) ((1 + i)^((2 e)/(-1 + e)) - a^2 i^2 W^2)))/(a i^2 (1 + i) ((1 + i)^(e/(1 - e)) - (1 + i)^(1/(L n - e L n))) ((1 + i)^(e/(1 - e)) - (1 + i)^((1 + L n)/(L n - e L n))) (1 + i + L n + e (-1 + (-2 + e - i) L n)) w W)); toBeZero = 1/((1 - e)^2 w) - v + z // Simplify; Plot[Evaluate[toBeZero /. w -> 1], {e, 0, 1}] Plot[Evaluate[toBeZero /. w -> 1], {e, .0191, .01915}, PlotRange -> 1] Your function is very wiggly. It seems a very hard problem. What makes you think that there is a solution at all? Posted 3 months ago  Thanks! But another equation of the system is missing if you code which is:$\frac{\partial e}{\partial w}=\frac{e}{w}$. By the way, the$e$function here,$e=e(w)$, is obtained from the first equation,$\frac{1}{(1-e)^2 w}=v-z\$.
Posted 3 months ago
 Your differential equation is easily solved, and it means that e(w)=c*w for some constant c. You are looking for solutions of the very complicated nonlinear equation toBeZero ==0 that are straight lines through the origin. Are you sure that there is such a straight line? I tried for a while, and could not find it.
Posted 4 months ago
 On a quick run, it seems to me that those Simplify take forever, on expressions half a Megabyte long.Then icSOL = NSolve[icEQ == 0 && 0 < e < 1 && 0 < w < 2, {e, w}] is unlikely to work, because the equations icEQ contain a third parameter beta.
Posted 4 months ago
 Thanks for catching it. beta was a typo. I corrected it to b. Can you please help what I can do then? Thank!