# How to solve the warning problem and obtain real roots without imaginary part?

I am trying to solve a equation with Newton's method via FindRoot, and the codes are:

Define the functions:

τs[Ωs_,t_]:=2*Ωs*Coth[Ωs/(2*t)]
τb[Ωb_,t_]:=2*Ωb*Coth[Ωb/(2*t)]
Xs[Js_,Jb_,Ωs_,t_,δ_]:=(τs[Ωs,t]-4*Js-4*Js/(Sqrt[2])^δ)/(Jb+4*Js/(Sqrt[2])^δ+4*Js/(Sqrt[3])^δ)
Xb[Jb_,Ωb_,t_,δ_]:=(τb[Ωb,t]-4*Jb-4*Jb/(Sqrt[2])^δ)/(Jb+4*Jb/(Sqrt[2])^δ+4*Jb/(Sqrt[3])^δ)
ϕ[Jb_,Ωb_,t_,δ_]:=ArcCosh[Xb[Jb,Ωb,t,δ]/2]
BN[Jb_,Ωb_,t_,δ_,n_]:=Sinh[(n+1)*ϕ[Jb,Ωb,t,δ]]/Sinh[ϕ[Jb,Ωb,t,δ]]
fumfa[Js_,Jb_,Ωs_,Ωb_,t_,δ_,n_]:=(Xs[Js,Jb,Ωs,t,δ])^2*BN[Jb,Ωb,t,δ,n-2]-2*Xs[Js,Jb,Ωs,t,δ]*BN[Jb,Ωb,t,δ,n-3]+BN[Jb,Ωb,t,δ,n-4]


Solve a series of solutions:

Table[FindRoot[fumfa[2.0,1.0,Ωs,2.0,3.5,t,4]==0.0,{Ωs,4.0}],{t,0.1,4.0,0.2}]


However, I got the following warning message:

The line search decreased the step size to within tolerance specified by AccuracyGoal and PrecisionGoal but was unable to find a sufficient decrease in the merit function. You may need more than MachinePrecision digits of working precision to meet these tolerances.

And some of the solutions have imaginary parts though the coefficients are actually zeros:

{{Ωs->11.6453 +0. I},{Ωs->10.5823 +0. I},{Ωs->9.58167 +0. I},{Ωs->8.63181 +0. I},{Ωs->7.72096 +0. I},{Ωs->6.83656 +0. I},{Ωs->5.96386 +0. I},{Ωs->5.08324 +0. I},{Ωs->4.16381 +0. I},{Ωs->3.14345 +0. I},{Ωs->1.82232},{Ωs->5.70687*10^-6},{Ωs->4.74748},{Ωs->3.29661},{Ωs->1.28711},{Ωs->-1.26412*10^-6},{Ωs->-3.05724*10^-6},{Ωs->-0.0000392875},{Ωs->0.0000453593},{Ωs->0.0000543535}}


My questions are:

1. How to remove the warning message?
2. why there are imaginary parts in the solutions and how to assure Reals only computation?
-
Do you want a "Reals only computation" or is removing the imaginary numbers afterwards acceptable? – Öskå Aug 1 '14 at 12:03
Re[]/.solutions can remove imaginary parts easily; I am wondering how to compute via reals only method. – LCFactorization Aug 1 '14 at 12:05
You can also use Chop to remove very small imaginaries caused by numerical methods applied to inexact numbers, which is your situation. – m_goldberg Aug 1 '14 at 12:11
BTW for the first question, just use Quiet@FindRoot[]. – Öskå Aug 1 '14 at 12:15
Often the warning message indicates that FindRoot got stuck at a non-root. You should check your solutions to see if they are accurate enough. – Michael E2 Aug 1 '14 at 14:57

## 1 Answer

First, you have t in the slot for δ -- that may be a mistake. Be that as it may, the the question about the warning FindRoot::lstol has an explanation.

Second, you're getting complex solutions because the function evaluates to complex numbers:

fumfa[2.0, 1.0, Ωs, 2.0, 3.5, t, 4] /. {t -> 0.1, Ωs -> 4.0}
(* 1.02171 + 0. I *)


Finally, the explanation. The warning FindRoot::lstol means that Mathematica was having trouble making the function(s) close to zero. This can happen at a local extremum, for instance.

We can use Check to catch the values of t that produce warnings. The value of t and the corresponding "solution" sol0 will be stored in err. All solutions together with their corresponding values of t will be stored in sol.

Module[{sol0},
err = {};
sols = Table[
Quiet@Check[
sol0 = Flatten @
{HoldPattern[t] -> t, FindRoot[fumfa[2.0, 1.0, Ωs, 2.0, 3.5, t, 4] == 0.0, {Ωs, 4.0}]},
AppendTo[err, sol0]; sol0,
{FindRoot::lstol}],
{t, 0.1, 4.0, 0.2}]
];


Here are the values of t that produce warnings and the corresponding values of the function:

t /. err
fumfa[2.0, 1.0, Ωs, 2.0, 3.5, t, 4] /. err
(*
{2.3, 3.1, 3.3, 3.5, 3.7, 3.9}
{-0.0643143, 0.298996, 0.960447, 2.04808, 3.70264, 6.09244}
*)


You can have such large errors when the derivative of the objective function is very large, since $dy = f'(x)\;dx$. That is not the case here:

D[fumfa[2.0, 1.0, Ωs, 2.0, 3.5, t, 4], Ωs] /. err
(* {-2.75827*10^-7, -1.67904*10^-7, -7.33596*10^-7, -0.0000152322, 0.0000266041, 0.0000460977} *)


In fact, we see that FindMinimum has gotten stuck at relative extrema:

GraphicsGrid @ Partition[
Hold[Plot[
fumfa[2.0, 1.0, x, 2.0, 3.5, t, 4], {x, Ωs - 0.1, Ωs + 0.1},
PlotLabel -> HoldForm["t" == t],
AxesLabel -> {"Ωs", "fumfa"}]] /. err //
ReleaseHold,
3]


-