0
$\begingroup$

My function is

(* condição para o valor de b0 *) σ = 0.6 ; minroot[g_?NumericQ, b_?NumericQ] := Module[{rts, r}, rts = r /. Solve[1 - (b/r)^2 - g^-2*(2/15*σ^9 (1/(r - 1)^9 - 1/(r + 1)^9 - 9/(8 r) (1/(r - 1)^8 - 1/(r + 1)^8)) - σ^3 (1/(r - 1)^3 - 1/(r + 1)^3 - 3/(2 r) (1/(r - 1)^2 - 1/(r + 1)^2))) == 0, r]; rts = Select[rts, With[{nval = N[#, 100]}, Im[nval] == 0 && nval > 0] &]; Max[rts]]; rootmin[g_?NumberQ] := Module[{Rrmts, b}, Rrmts = b /. FindRoot[Re[aA[g, b, 5]] == 0, {b, 1, 2}, Method -> "Brent"]; Rrmts]; (* angulo de espalhamento *) aA[g_?NumberQ, b_?NumberQ, i_] := Pi - 2 b NIntegrate[ 1/(r^2*Sqrt[ 1 - (b/r)^2 - g^-2*(2/15*σ^9 (1/(r - 1)^9 - 1/(r + 1)^9 - 9/(8 r) (1/(r - 1)^8 - 1/(r + 1)^8)) - σ^3 (1/(r - 1)^3 - 1/(r + 1)^3 - 3/(2 r) (1/(r - 1)^2 - 1/(r + 1)^2)))]), {r, minroot[g, b], Infinity}, Exclusions -> {r^2* Sqrt[1 - (b/r)^2 - g^-2*(2/ 15*σ^9 (1/(r - 1)^9 - 1/(r + 1)^9 - 9/(8 r) (1/(r - 1)^8 - 1/(r + 1)^8)) - σ^3 (1/(r - 1)^3 - 1/(r + 1)^3 - 3/(2 r) (1/(r - 1)^2 - 1/(r + 1)^2)))] == 0}, MaxRecursion -> i, Method -> {Automatic, "SymbolicProcessing" -> 0}]; 

When I evaluate

rootmin[0.3] 

I get

Solve::ratnz: Solve was unable to solve the system with inexact coefficients. The answer was obtained by solving a corresponding exact system and numericizing the result.

but when I evaluate

rootmin[0.4] 
1.95804 

I still get some errors, but I get a result.

Why? I want to evaluate rootmin[0.1], for example, and I can't.

$\endgroup$
14
  • 3
    $\begingroup$ I get a bunch of different errors as well: NIntegrate::inumr, FindRoot::nlnum etc. Also, in the definition aA, you use minroot[g,b]; Is it the same as rootmin[g,b]? $\endgroup$ Commented May 12, 2015 at 21:12
  • $\begingroup$ @Mahdi, minroot[g,b] is a value of r , and rootmin[g,b] is a value of b , so they are differents ... do you have an ideia of all this errrors ? $\endgroup$ Commented May 13, 2015 at 12:16
  • $\begingroup$ Please edit the question and give a definition for your function minroot. $\endgroup$ Commented Aug 26, 2015 at 8:20
  • $\begingroup$ @m_goldberg , sorry ... now its ok ... any ideas ? $\endgroup$ Commented Aug 26, 2015 at 16:02
  • $\begingroup$ I hope this question will be reopened. Now that you have supplied a definition for minroot, it is pretty clear what your problem is, and I have an answer I would like to post. $\endgroup$ Commented Aug 26, 2015 at 18:54

1 Answer 1

3
$\begingroup$

The first thing I did was to rationalizing all calculations, starting with the defintion of σ and minroot. This stops the Solve::ratnz messages. I also made some other improvements to minroot.

σ = 6/10; minroot[gg_?NumericQ, bb_?NumericQ] := Module[{b, g, rts, r}, b = Rationalize[bb, 0]; g = Rationalize[gg, 0]; rts = r /. Solve[ 1 - (b/r)^2 - g^-2*(2/15*σ^9 (1/(r - 1)^9 - 1/(r + 1)^9 - 9/(8 r) (1/(r - 1)^8 - 1/(r + 1)^8)) - σ^3 (1/(r - 1)^3 - 1/(r + 1)^3 - 3/(2 r) (1/(r - 1)^2 - 1/(r + 1)^2))) == 0, r]; Max[Select[rts, Im[#] == 0 && # > 0 &]]] 

The definition of aA is more or less untouched.

aA[g_?NumberQ, b_?NumberQ, i_] := (Pi - 2 b NIntegrate[ 1/(r^2* Sqrt[ 1 - (b/r)^2 - g^(-2)* (2/15*σ^9 (1/(r - 1)^9 - 1/(r + 1)^9 - 9/(8 r) * (1/(r - 1)^8 - 1/(r + 1)^8)) - σ^3 (1/(r - 1)^3 -1/(r + 1)^3 -3/(2 r) * (1/(r - 1)^2 - 1/(r + 1)^2)))]), {r, minroot[g, b], ∞}, Exclusions -> {r^2*Sqrt[1 - (b/r)^2 - g^(-2)* (2/15*σ^9 (1/(r - 1)^9 - 1/(r + 1)^9 - 9/(8 r) * (1/(r - 1)^8 - 1/(r + 1)^8)) - σ^3 (1/(r - 1)^3 - 1/(r + 1)^3 - 3/(2 r) * (1/(r - 1)^2 - 1/(r + 1)^2)))] == 0}, MaxRecursion -> i, Method -> {Automatic, "SymbolicProcessing" -> False}]) 

A plot of aA over the interval [1, 2] for several values of g shows the problem with small values of g.

Plot[Evaluate[Re[aA[#, b, 5]]& /@ #], {b, 1, 2}, PlotPoints -> 5, PlotRange -> {-.25, 2.25}, PlotLegends -> Evaluate[( Style[Row[{"g = ", #}], 12]&) /@ #], AspectRatio -> 1, ImageSize -> Medium]&[{.1, .3, .36, .4}] 

plots

We now know why rootmin fails for g < .36116; however, using the version of minroot shown above, rootmin produces no messages for values of g >= .36116. I do recommend simplifying it as shown below.

rootmin[g_?NumberQ] := Module[{b}, b /. FindRoot[Re[aA[g, b, 5]] == 0, {b, 1, 2}, Method -> "Brent"]] 

Now we look at

rootmin[.36116] 
2. 

It is the critical point because it is where the zero of the real part of aA falls right at the end of the interval [1, 2].

$\endgroup$

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.