1
$\begingroup$

I am trying to solve numerically a difficult implicit equation and plot the solution. The thing that I want to solve numerically is the following:

$$ Q(h)=\frac{0.13}{1-10^{-4}\log{\frac{Y(h)}{100}}} $$

$$ P(h)=2*10^{4}*\sqrt{\frac{1}{1-10^{-4}\log{\frac{Y(h)}{100}}}}-1/6$$

$$ Y(h)=\sqrt{3*Q(h)*h^2-(P(h)+1/6)*10^{12}} $$

So the point is that I want to solve (numerically is the only way) the implicit equation of Y(h) in terms of h, but note that functions Q(h) and P(h) also depend on Y (so implicitly on h).

I want to get numerically the solution of Y(h) in order to use this to plot Q(h) and P(h) in terms of h.

Someone knows how to do this??? Thanks!

Expressions in format code:

Q[h_] := 0.13/(1 - 10^(-4)*Log[Y[h]/100]);

P[h_] := 2*(10^4)*Sqrt[1/(1 - 10^(-4)*Log[Y[h]/100])] - 1/6;

Y[h_] := Sqrt[3*Q[h]*h^2 - (P[h] + 1/6)*10^12];
| improve this question | |
$\endgroup$
  • 3
    $\begingroup$ What is the range of h and Yh? $\endgroup$ – Ulrich Neumann Apr 11 '18 at 11:29
3
$\begingroup$

You can solve your implicit equation as follows:

First I renamed Y[h]->y

Q[h_] := 0.13/(1 - 10^(-4)*Log[y/100])
P[h_] := 2*(10^4)*Sqrt[1/(1 - 10^(-4)*Log[y/100])] - 1/6

which gives the implicit equation 

zero = y - Sqrt[3*Q[h]*h^2 - (P[h] + 1/6)*10^12] /. y -> Exp[logy] //PowerExpand

Thereby I substituted y -> Exp[logy] because of poor scaling. The solution of zero==0 can be evaluated explicitely for h[logy] (not so easy for logy[h]) 

sol = Solve[zero == 0, h ];
ParametricPlot[ {h, Exp[logy]} /. sol, {logy, 0, 25},AxesLabel -> {h, y}, PlotLabel -> zero == 0,PlotRange -> {{-10^9, 10^9}, {0, 10^9 }}]

enter image description here

Q and P can be plotted 

GraphicsRow[{
ParametricPlot[({h, Q[h]} /. y -> Exp[logy] //PowerExpand) /. sol[[2]], {logy, 0, 25}, AspectRatio -> 1,AxesLabel -> {h, "Q[h]"}], 
ParametricPlot[({h, P[h]} /. y -> Exp[logy] // PowerExpand) /.sol[[2]], {logy, 0, 25}, AspectRatio -> 1,AxesLabel -> {h, "P[h]"}]}]

enter image description here

| improve this answer | |
$\endgroup$
  • $\begingroup$ Thank you it's okay but I was searching to get a plot of Q(h) and P(h) in terms of h, once one has solved the equation for y(h) as you have done. Do you know how do that? $\endgroup$ – Joe Apr 11 '18 at 13:40
  • $\begingroup$ Where is the problem? You know the solution pairs {h,y}!!! ParametricPlot[({h, Q[h]} /. y -> Exp[logy] // PowerExpand ) /. sol[[2]], {logy, 0, 25}, AspectRatio -> 1] $\endgroup$ – Ulrich Neumann Apr 11 '18 at 13:53

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.