4
$\begingroup$

I have defined two functions to fit $\dfrac{a}{\sqrt{x^2-b^2}}$ and $\dfrac{a}{\sqrt{b^2-x^2}}$.

fitSingularity[listFit_] := 
  Module[{res, startVal1, model, a, b, fit},
   startVal1 = listFit[[1]][[1]]; Print[startVal1];
   res = Sqrt[(listFit[[6]][[1]]^2 - startVal1^2)]*listFit[[6]][[2]];
   model = a/Sqrt[(x^2 - b^2)];
   fit = FindFit[listFit[[2 ;; -1]], 
     model, {{a, res}, {b, startVal1}}, x];
   a = a /. fit;
   b = b /. fit;
   {a, b}
   ];

fitSingularity2[listFit_] := Module[{res, startVal1, model, a, b, fit},
   startVal1 = listFit[[-1]][[1]]; Print[N[startVal1]];
   res = Sqrt[(startVal1^2 - listFit[[6]][[1]]^2)]*listFit[[6]][[2]]; 
   Print[N[res]];
   model = a/Sqrt[(b^2 - x^2)];
   fit = FindFit[listFit[[2 ;; -1]], 
     model, {{a, res}, {b, startVal1}}, x];
   a = a /. fit;
   b = b /. fit;
   {a, b}
   ];

To realize the fits I create a list. In the first case with the following code everything is working nicely

g[x_] := a/Sqrt[x^2 - b^2]
a = 1; b = 2;
list1 = Table[{x, g[x]}, {x, 1.999, 4, 0.1}];
Show[Plot[g[x], {x, 2, 4}], ListPlot[list1, PlotStyle -> Green]]
fitSingularity[list1]

But in the second case I obtain an error and I do not know why (it should not be so different than for the first case)

f[x_] := a/Sqrt[b^2 - x^2]
a = 1; b = 2;
list2 = Table[{x, f[x]}, {x, 0, 1.999, 0.01}];
Show[Plot[f[x], {x, 0, 2}], ListPlot[list2, PlotStyle -> Green]]
fitSingularity2[list2]

I would like to mention that I implement this to be sure I can fit the function well, but I need to use this fit in a much more complicated case (inside a numerical function that is supposed to have a singularity of the form $\dfrac{a}{\sqrt{x^2-b^2}}$ or $\dfrac{a}{\sqrt{b^2-x^2}}$)

$\endgroup$
  • $\begingroup$ Welcome to Mathematica.SE! I suggest the following: 0) Browse the common pitfalls question 1) As you receive help, try to give it too, by answering questions in your area of expertise. 2) Read the faq! 3) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! $\endgroup$ – Dr. belisarius Sep 21 '15 at 21:40
5
$\begingroup$

You may work this out by establishing a region where the square roots are positive. I'm using Reduce[ ] here, but any bounding method will work:

sqRootG[x_, b_] := x^2 - b^2
sqRootF[x_, b_] := b^2 - x^2
g[x_] := a/Sqrt[sqRootG[x, b]]
f[x_] := a/Sqrt[sqRootF[x, b]]
listg = Rest@Table[{x, g[x] /. {a -> 1, b -> 2}}, {x, 1.999, 4, 0.01}];
listf = Rest@Table[{x, f[x] /. {a -> 1, b -> 2}}, {x, 0, 1.999, 0.01}];
sg = Reduce[And @@ Thread[sqRootG[#, b] & /@ listg[[All, 1]] > 0], {a, b}, Reals]; 
sf = Reduce[And @@ Thread[sqRootF[#, b] & /@ listf[[All, 1]] > 0], {a, b}, Reals];
NonlinearModelFit[listg, {g[x], sg}, {a, b}, x, 
                  Method -> {"NMinimize", Method -> "NelderMead"}]["BestFitParameters"]
NonlinearModelFit[listf, {f[x], sf}, {a, b}, x, 
                  Method -> {"NMinimize", Method -> "NelderMead"}]["BestFitParameters"]

(*
{a -> 1., b -> 2.}
{a -> 1., b -> 2.}
*)
$\endgroup$
  • $\begingroup$ Thank you for your answer. Is NonlinearModelFit better than FindFit?? I am not sure I completely understand the And @@ Thread[sqRootG[#, b] & /@ listg[[All, 1]] > 0] part, I will try to read the help more carefully. Thanks again $\endgroup$ – lambertmular Sep 21 '15 at 21:36
  • $\begingroup$ As I state below (comment to Jack LaVigne), I try your code with a FindFit[listf, {f[x], sf}, {{a, 1}, {b, 1.99}}, x] and get the same error message as I had initially with basically no value returned. So now I am not sure the problem comes from the constraint but probably on FindFit. What do you think? $\endgroup$ – lambertmular Sep 22 '15 at 9:09
  • $\begingroup$ In fact now I saw that only NonlinearModelFit together with the constraint is working for this simple example. Thank you very much. However, when I implement it in my more complicated numerical function (the one I mention at the end of my question), your suggestion does not work (even with the Reduce etc...). I obtain an error Power::infy: "Infinite expression 1/0.^0.5 encountered." I don't know what to do. $\endgroup$ – lambertmular Sep 22 '15 at 10:07
  • $\begingroup$ My other question related to the more complex case was posted here mathematica.stackexchange.com/questions/95220/… $\endgroup$ – lambertmular Sep 22 '15 at 13:45
3
$\begingroup$

This is not an answer but a comment that is too long to fit in a comment space.

The goal is to elaborate on belisarius's answer.

Copy and execute the code up to NonLinearModelFit[... from belisarius's answer.

sg evaluates to be

-2.001 < b < 2.001

The intent is to use it as a constraint in NonlinearModelFit.

To break down the complete expression

sg = Reduce[And @@ Thread[sqRootG[#, b] & /@ listg[[All, 1]] > 0], {a, b}, Reals]

It is helpful to start with the inner expression and evaluate it step by step. In the example I will only show the first four terms of the result.

listg[[All,1]] takes the first term (i.e., x) from listg pairs.

{2.001, 2.011, 2.021, 2.031, ...}

Next using these input values

sqRootG[#, b] & /@ listg[[All, 1]]

makes a new list of x^2 - b^2 using the x value from above.

{4.004 - b^2, 4.04412 - b^2, 4.08444 - b^2, 4.12496 - b^2, ...}

Thread adds >0 to the list.

Thread[sqRootG[#, b] & /@ listg[[All, 1]] > 0]

and gives

{4.004 - b^2 > 0, 4.04412 - b^2 > 0, 4.08444 - b^2 > 0, 
 4.12496 - b^2 > 0, ...}

And makes it into a single constraint expression.

And @@ Thread[sqRootG[#, b] & /@ listg[[All, 1]] > 0]

4.004 - b^2 > 0 && 4.04412 - b^2 > 0 && 4.08444 - b^2 > 0 && 
 4.12496 - b^2 > 0 ...

Finally then Reduce evaluates this and gives the final constraint expression

-2.001 < b < 2.001

With regard to FindFit the constraint works equally well with that function

FindFit[listg, {g[x], sg}, {a, b}, x]

Used like this you get a warning message and a result of 0.9999 for a and -2 for b.

You can use the first result to refine the search with starting values

FindFit[listg, {g[x], sg}, {{a, 0.99}, {b, 2.}}, x]

Which rapidly gives

{a -> 1.00001, b -> 2.}

Similar remarks apply to the second function.

$\endgroup$
  • $\begingroup$ Waouh thank you!! $\endgroup$ – lambertmular Sep 21 '15 at 23:28
  • $\begingroup$ You are welcome. I think belisarius, along with several other high ranking members on this site, is brilliant. I find that frequently I have to do some work (because of my lack of experience) to break down their answer before I understand the answer in total. The comment was intended to show you a method that I have found helpful. $\endgroup$ – Jack LaVigne Sep 21 '15 at 23:42
  • $\begingroup$ Yes your method is very good. I will apply it next time! Thanks again $\endgroup$ – lambertmular Sep 22 '15 at 8:33
  • $\begingroup$ Sorry but your statement does not work with the f function. I pass FindFit[listf, {f[x], sf}, {a, b}, x] and get a wrong values. After I pass FindFit[listf, {f[x], sf}, {{a, 1}, {b, 2}}, x] and get the good values but with error message $\endgroup$ – lambertmular Sep 22 '15 at 8:47
  • $\begingroup$ Even worse, If I use FindFit[listf, {f[x], sf}, {{a, 1}, {b, 1.99}}, x] (1.99 as a start value instead of 2), I obtain the same error that I had in my initial question with basically no value returned. So finally, I think the problem comes from FindFit not from the constraint. I am still investigating $\endgroup$ – lambertmular Sep 22 '15 at 9:06

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.