0
$\begingroup$

I have a data named Hc2Tdata where its plot is shown in the first figure. The x-axis is $T$ while the y-axis is $H_{c2}(T)$. I want to fit this with the equation named WHHeqn[T,Hc2T,a,b] where it is an implicit equation of $H_{c2}(T)$ and $T$ while $\{a,b\}$ are fitting parameters. Since the equation is implicit, I followed the answer in this post in order to do the fitting.

In that post, I would have to use FindFit where T is the variable. Since Hc2T ($H_{c2}(T)$) is implicit in the equation, I would have to employ FindRoot to extract Hc2T, then FindFit can do the fitting for $\{a,b\}$. Based on the data plot, the $H_{c2}(T)$ data is in the range $[0,8]$, so I chose {Hc2T,4} as the starting point for FindRoot.

After executing FindFit it gave values for $\{a,b\}$ but there are errors related to line search, so I'm not sure if it is correct or not. In addition, I already know that the fitting curve has to follow the data for larger $T$ then deviate away from the data points as it goes to lower $T$, i.e., the data should be a little bit higher than the fitting curve for lower $T$. Using the values returned, the fitting curve looks as shown in the second figure. Any help or guidance on this?

data

data1

Temp = {7, 8, 9, 9.5, 10, 10.5, 11}; Hc2 = {8.1, 6.2, 4.7, 4.1, 2.7, 0.6, 0}; Hc2Tdata = Transpose[{Temp, Hc2}]; eulerconst = 0.5772; Tc = 11; slope = (Hc2[[6]] - Hc2[[7]])/(Temp[[6]] - Temp[[7]]); Hc20 = -Pi^2/(8 Exp[eulerconst]) (Tc) (slope); hs[Hc2T_] = Hc2T/(-slope Tc); hb[Hc2T_] = hs[Hc2T]/(Pi^2/4); gam[Hc2T_, a_, b_] = Sqrt[(a hb[Hc2T])^2 - (0.5 b)^2]; WHHeqn[T_, Hc2T_, a_, b_] = -Log[1/(T/Tc)] + (0.5 + (I b)/(4 gam[Hc2T, a, b])) PolyGamma[0.5 + (hb[Hc2T] + 0.5 b + I gam[Hc2T, a, b])/(2 (T/Tc))] + (0.5 - (I b)/(4 gam[Hc2T, a, b])) PolyGamma[0.5 + (hb[Hc2T] + 0.5 b - I gam[Hc2T, a, b])/(2 (T/Tc))] - PolyGamma[0.5]; fitfunc[T_?NumericQ, a_?NumericQ, b_?NumericQ] := Hc2T /. FindRoot[Re[SetPrecision[WHHeqn[T, Hc2T, a, b], 35]] == 0, {Hc2T, 4}, AccuracyGoal -> 5, PrecisionGoal -> 5, WorkingPrecision -> 30] FindFit[Hc2Tdata, fitfunc[T, a, b], {a, b}, T] FindRoot::lstol: 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 30.` digits of working precision to meet these tolerances. FindRoot::lstol: 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 30.` digits of working precision to meet these tolerances. FindRoot::lstol: 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 30.` digits of working precision to meet these tolerances. General::stop: Further output of FindRoot::lstol will be suppressed during this calculation. FindFit::sszero: The step size in the search has become less than the tolerance prescribed by the PrecisionGoal option, but the gradient is larger than the tolerance specified by the AccuracyGoal option. There is a possibility that the method has stalled at a point that is not a local minimum. {a -> 3.112012071*10^-8, b -> 10.36533923}

Update: The slope that I need to calculate is the slope of the data near $T=11$. Originally, I calculated the slope using points 6 & 7 ($T=10.5$ & $T=11$). Due to a small fluctuation of the data at $T=10.5$, when I calculate the slope there for sure it will greatly affect the final result as compared to say, when I calculate the slope at points 5 & 7. I think maybe I can just use points 5 & 7 so that the slope would give a more consistent fitting.

Improve this question
$\endgroup$
8
  • $\begingroup$ Why not try plotting the fit result in the same plot as the data to see if the fit is good? $\endgroup$ Commented May 6 at 13:30
  • $\begingroup$ @SjoerdSmit I have updated the post including the fitting curve. As you see, the fitting curve only covered two points. I'm not saying this is wrong or correct. It's just that I received errors when calculating the fitting parameters, so I'm unsure whether the result is correct or not. $\endgroup$ Commented May 6 at 13:54
  • $\begingroup$ @Bill Good point about the code depending strongly on the slope! As for the fitting function, that is the required function to be used. On the other hand, I just wanted to know why the errors are coming out, and whether there are fixes needed to be done. Or the errors are just complaints by MMA such that even if the errors are fixed it will give almost the same parameters? $\endgroup$ Commented May 6 at 14:20
  • $\begingroup$ The three-dot button will show the stack trace. This reveals the values of a and b being tested. Are they reasonable, or should you restrict the search range in FindFit? $\endgroup$ Commented May 6 at 14:38
  • 1
    $\begingroup$ I wouldn't be concerned about errors for parameter values for a and b that are large negative numbers, and therefore irrelevant with regard to your reasonable range, not to mention that FindFit[] does not return those values for a and b. FWIW, here are two examples giving the same error, one having a root and the other not: FindRoot[CubeRoot[x^2 - 2] == 0, {x, 1}] and FindRoot[CubeRoot[1*^16 x^2 - 2*^16]^2 + 10 == 0, {x, 1}]. You might examine your equation and see if something similar is going on, if you are interested in parameter values like {a -> -76.9678, b -> -257.561}. $\endgroup$ Commented May 6 at 19:31

1 Answer 1

Reset to default
4
$\begingroup$

Assuming, as you mention in the comments {a,b}>0 and using data 5 and 7 for the slope, we may define an error function and search for a and b that minimize this function.

err[a_, b_] := Abs[Total[ WHHeqn[#[[1]], #[[2]], a, b]^2 & /@ Hc2Tdata]];

Now we can minimize. This gives some messages about encountering 1/0 during searching, but nevertheless we get a result:

sol = {a, b} /. NMinimize[{err[a, b], a > 0, b > 0}, {a, b}][[2]] {1.44783, 5.9091*10^-9}

With this we can plot the intrinsic function together with the data:

gr1 = ContourPlot[ WHHeqn[T, Hc2T, sol[[1]], sol[[2]]] == 0, {T, 5, 11}, {Hc2T, 0, 8.5}]; gr2 = Graphics[Point[Hc2Tdata]]; Show[gr1, gr2]

enter image description here

Improve this answer
$\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.