2
$\begingroup$

How do I find the zero of this function:

Ef=5.53;
e=1.6*10^-19;
qc[w1_?NumericQ] := 
  w1 (1/(2 a[w1]) (1 + Sqrt[1 - 4/(27 a[w1])])^(1/3) + 
     1/(2 a[w1]) (1 - Sqrt[1 - 4/(27 a[w1])])^(1/3));

a[w1_?NumericQ] := 
  Piecewise[{{2/3, 
     w1 < 0.62}, {0.9069 + 0.3577 w1^(-2/3) - 1.0565 w1^(-1) + 
      1.478 w1^(-4/3) - 0.4524 w1^(-5/3), w1 > 0.62}}];
wmax2[w1_?NumericQ, T_?NumericQ] := Min[T - Ef, T/2];
s4[w1_?NumericQ, T_?NumericQ] := (Pi e^4)/
   T (Log[wmax2[w1, T]/qc[w1]] + T/(T - wmax2[w1, T]) - T/(
     T - qc[w1]) + 2 Log[(T - wmax2[w1, T])/(T - qc[w1])]);

If I plot it it crosses the x axis, so there is a zero, but NSolve[s4[100,T]==0,T] gives a complex number as a solution and says to use Refine, and Refine can't find the zero.

| improve this question | |
$\endgroup$
  • 1
    $\begingroup$ FindRoot perhaps? You haven't defined qc for us, so we can't investigate more. $\endgroup$ – John Doty Jul 4 '18 at 14:50
  • 2
    $\begingroup$ Use FindRoot. Can't say more, because you didn't give full definitions of your symbols (e.g. wmax2, qc). You have asked 48 questions on this site so far, you should know how to ask a good question by now. $\endgroup$ – MarcoB Jul 4 '18 at 14:50
  • $\begingroup$ @MarcoB added the rest of the function I am sorry. $\endgroup$ – mattiav27 Jul 4 '18 at 15:01
  • 2
    $\begingroup$ @mattiav27 a is still undefined. The best way to check for these things is to copy back your posted code in a clean MMA session (see this Q) and try to execute it. It has to work under those conditions, because that's what we will be trying on our own computers. $\endgroup$ – MarcoB Jul 4 '18 at 15:26
  • $\begingroup$ @MarcoB Added the definition of a $\endgroup$ – mattiav27 Jul 4 '18 at 16:53
2
$\begingroup$

With the full code, we can explore the plot to find an initial estimate of the root. COnsidering the extremely small values of the function, I used a LogPlot on its absolute value:

LogPlot[Abs@s4[100, T], {T, 100, 250}, PlotRange -> All]

Mathematica graphics

We confirm that with a regular plot:

Plot[s4[100, T], {T, 120, 250}, PlotRange -> All]

Mathematica graphics

A good estimate for the root appears to be around $180$. We can feed that into FindRoot as a starting point for root search:

FindRoot[s4[100, T], {T, 180}]
(* Out: {T -> 183.495} *)
| improve this answer | |
$\endgroup$
0
$\begingroup$

With a large value for PlotPoints we can also use purely graphical approach to get the roots: 

Normal[Plot[s4[100, T], {T, 120, 250}, PlotRange -> All, 
   MeshFunctions -> {#2 &}, Mesh -> {{0}}, 
   MeshStyle -> PointSize[Large], PlotPoints -> 300]] /. 
 Point[x_] :> {Point[x], Text[Framed @ Style[x, 14], Offset[{50,-30}, x]]}

enter image description here

Use 

Text[Framed@Style[root = T/. FindRoot[s4[100, T], {T, x[[1]]}, WorkingPrecision->20 ], 12]

in place of Text[...] above to get

enter image description here

| improve this answer | |
$\endgroup$

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.