L2 = zeros(14);
that should probably be
L2 = zeros(14,1);
like the other variables.
how to perform data fit like excel? and plot
4 次查看（过去 30 天）
显示 更早的评论
- I have observed array of data ( y_obs) and predicted data (y_pred)
- Predicted data is obtained from an equation
- How do I fit the observed data to the predicted data by minimizing the coefficient "d" in the equation? ( This is possible in excel, but I could not find a suitable method in matlab
Below is my code for steps 1 and 2:
% observed data
y_obs = [0.3 0.2 0.28 0.318 0.421 0.492 0.572 0.55 0.63 0.61 0.73 0.8 0.81 0.84 0.93 0.91]'; % If y_obs should equal to predicted, I can have more data. J us fo rthe code, I am providing limited observed data
t1 = [300:300:21600]';
a=0.0011;
gama = 0.01005;
d=0.000000000302;
n=1;
t=300;
L2 = zeros(14,1);
L3 = zeros(14,1);
L4 = zeros(14,1);
At = zeros(14,1);
t = 300;
n =1;
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
format shortE
for t= 300:300:21600
for n=1:1:14
L2(n) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n) = ((L2(n).*L3(n))/L4(n));
end
S(t/300) = sum(L5);
y_pred(t/300) = 1 -L1*S(t/300); % predicted data
end
2 个评论
Anand Ra
2021-6-16
Thanks for the response, I can update it.
Can you please guide me on how to perform the data fitting in the fashion I described in bullet point 3?
采纳的回答
Walter Roberson
2021-6-16
编辑：Walter Roberson
2021-6-16
format shortE
% observed data
y_obs = [0.3 0.2 0.28 0.318 0.421 0.492 0.572 0.55 0.63 0.61 0.73 0.8 0.81 0.84 0.93 0.91]'; % If y_obs should equal to predicted, I can have more data. J us fo rthe code, I am providing limited observed data
t1 = [300:300:21600]';
T1 = t1(1:length(y_obs)).';
a=0.0011;
gama = 0.01005;
d0 = 0.000000000302;
syms d
n=1;
t=300;
L2 = sym(zeros(14,1));
L3 = sym(zeros(14,1));
L4 = sym(zeros(14,1));
At = sym(zeros(14,1));
t = 300;
n =1;
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
y_pred = sym(zeros(length(T1),1));
for t = T1
for n=1:1:14
L2(n) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n) = ((L2(n).*L3(n))/L4(n));
end
S(t/300) = sum(L5);
y_pred(t/300) = 1 -L1*S(t/300); % predicted data
end
sse = expand(sum((y_pred - y_obs(:)).^2));
f = matlabFunction(sse)
ans = 
opt1 = fmincon(f, d0)
opt2 = fminsearch(f, d0)
35 个评论
Walter Roberson
2021-6-16
My test system got really slow suddenly, so I have not tested the fmincon / fminsearch steps yet.
Walter Roberson
2021-6-16
Please note that you are building 72 y_pred but your y_obs only has 16 entries. I coded above to only go as far in predictions as you have y_obs for.
Walter Roberson
2021-6-16
fmincon suggests d about 1e-11 and fminsearch suggests d about 3e-13 . But d = 0 is better than 1e-11, and the value at -3e-13 is a little better than that. -3.3508267261899e-13 checks out as a local minima when I plot.
Anand Ra
2021-6-16
Thanks a lot Walter!!
However, I am encountering the following error when trying to execute your code. Been trying to debug, but no luck. Following is the error:
Unable to perform assignment because value of type 'sym' is not convertible to 'double'.
Error in att1 (line 36)
L5(n) = ((L2(n).*L3(n))/L4(n));
Caused by:
Error using symengine
Unable to convert expression containing symbolic variables into double array. Apply 'subs' function first to substitute values for variables.
I am not sure how to use the subs function to the coefficient "d"
Anand Ra
2021-6-16
Sorry, one last request. Trying to plot the observed vs predicted fit following the minimization.
I am unable to pull the minimized data. How do I obtain "y_pred" data array?
I am looking to plot the observed vsp redicted after minimization.
Anand Ra
2021-6-17
Sorry, one last request. Trying to plot the observed vs predicted fit following the minimization.
I am unable to pull the minimized data. How do I obtain "y_pred" data array?
I am looking to plot the observed vsp redicted after minimization.
Walter Roberson
2021-6-17
%assuming that opt is the optimal d value:
y_pred_n = double(subs(y_pred, d, opt));
plot(T1, y_obs, 'k*', T1, y_pred_n, 'b')
legend({'observed', 'modeled'})
Anand Ra
2021-6-17
Thanks, but when plotted why does y_pred_n ends up as a straight line =1?
Isnt it suppose to fit the y_obs?
Anand Ra
2021-6-17
And more over, the predicted minmum is same as the inital value of d.
opt1 =opt2 = inital assignment d0. Not sure why.
Walter Roberson
2021-6-17
format long g
% observed data
y_obs = [0.3 0.2 0.28 0.318 0.421 0.492 0.572 0.55 0.63 0.61 0.73 0.8 0.81 0.84 0.93 0.91]'; % If y_obs should equal to predicted, I can have more data. J us fo rthe code, I am providing limited observed data
t1 = [300:300:21600]';
T1 = t1(1:length(y_obs)).';
a=0.0011;
gama = 0.01005;
d0 = 0.000000000302;
syms d
n=1;
t=300;
L2 = sym(zeros(14,1));
L3 = sym(zeros(14,1));
L4 = sym(zeros(14,1));
At = sym(zeros(14,1));
t = 300;
n =1;
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
y_pred = sym(zeros(length(T1),1));
for t = T1
for n=1:1:14
L2(n) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n) = ((L2(n).*L3(n))/L4(n));
end
S(t/300) = sum(L5);
y_pred(t/300) = 1 -L1*S(t/300); % predicted data
end
sse = expand(sum((y_pred - y_obs(:)).^2));
f = matlabFunction(sse)
f = function_handle with value:
@(d)exp(d.*pi.^2.*(-1.146818181818182e+12)).*9.289221334321387e-7+exp(d.*pi.^2.*(-9.94090909090909e+11)).*1.236275809772296e-6+exp(d.*pi.^2.*(-8.522727272727272e+11)).*1.681943238408239e-6+exp(d.*pi.^2.*(-7.213636363636363e+11)).*2.347788685302854e-6+exp(d.*pi.^2.*(-6.013636363636363e+11)).*3.378269437248533e-6-exp(d.*pi.^2.*(-5.734090909090909e+11)).*5.204552758007278e-4-exp(d.*pi.^2.*(-4.970454545454545e+11)).*6.004148783379719e-4+exp(d.*pi.^2.*(-4.922727272727272e+11)).*5.041454287362399e-6-exp(d.*pi.^2.*(-4.261363636363636e+11)).*7.003246735049698e-4+exp(d.*pi.^2.*(-3.940909090909091e+11)).*7.866397979720083e-6-exp(d.*pi.^2.*(-3.606818181818182e+11)).*8.274147573220533e-4+exp(d.*pi.^2.*(-3.068181818181818e+11)).*1.297791811534305e-5-exp(d.*pi.^2.*(-3.006818181818182e+11)).*9.925237366943283e-4-exp(d.*pi.^2.*(-2.461363636363636e+11)).*1.212471884290467e-3+exp(d.*pi.^2.*(-2.4e+11)).*4.487443864485274e-5+exp(d.*pi.^2.*(-2.304545454545454e+11)).*2.300373797956269e-5+exp(d.*pi.^2.*(-2.25e+11)).*4.487443864485274e-5+exp(d.*pi.^2.*(-2.1e+11)).*4.487443864485274e-5-exp(d.*pi.^2.*(-1.970454545454545e+11)).*1.514543380325033e-3+exp(d.*pi.^2.*(-1.95e+11)).*4.487443864485274e-5+exp(d.*pi.^2.*(-1.8e+11)).*4.487443864485274e-5+exp(d.*pi.^2.*(-1.65e+11)).*4.487443864485274e-5-exp(d.*pi.^2.*(-1.534090909090909e+11)).*1.945343394476675e-3+exp(d.*pi.^2.*(-1.5e+11)).*4.487443864485274e-5+exp(d.*pi.^2.*(-1.35e+11)).*4.487443864485274e-5-exp(d.*pi.^2.*(-1.2e+11)).*1.160916510107554e-3-exp(d.*pi.^2.*(-1.152272727272727e+11)).*2.589959458146108e-3-exp(d.*pi.^2.*(-1.125e+11)).*9.378374045852051e-4+exp(d.*pi.^2.*(-1.104545454545454e+11)).*1.001387582346529e-4-exp(d.*pi.^2.*(-1.05e+11)).*2.098753914692759e-3-exp(d.*pi.^2.*(-9.75e+10)).*2.545558669588414e-3-exp(d.*pi.^2.*(-9.0e+10)).*2.634661003027162e-3-exp(d.*pi.^2.*(-8.25e+10)).*3.61737284625722e-3-exp(d.*pi.^2.*(-7.5e+10)).*5.180219672615576e-3-exp(d.*pi.^2.*(-6.75e+10)).*4.957140567093227e-3+exp(d.*pi.^2.*(-6.681818181818181e+10)).*2.736383570858638e-4-exp(d.*pi.^2.*(-6.0e+10)).*5.98408030511718e-3-exp(d.*pi.^2.*(-5.522727272727272e+10)).*5.403745173601804e-3-exp(d.*pi.^2.*(-5.25e+10)).*5.734205845178111e-3-exp(d.*pi.^2.*(-4.5e+10)).*6.761145583202065e-3-exp(d.*pi.^2.*(-3.75e+10)).*7.757255103640483e-3+exp(d.*pi.^2.*(-3.409090909090909e+10)).*1.051219258989808e-3-exp(d.*pi.^2.*(-3.340909090909091e+10)).*8.932689680395143e-3-exp(d.*pi.^2.*(-3.0e+10)).*9.092341417456717e-3-exp(d.*pi.^2.*(-2.25e+10)).*9.646327590019253e-3-exp(d.*pi.^2.*(-1.704545454545454e+10)).*1.75081562684775e-2-exp(d.*pi.^2.*(-1.5e+10)).*1.067326732804321e-2+exp(d.*pi.^2.*(-1.227272727272727e+10)).*8.111137929522739e-3-exp(d.*pi.^2.*(-7.5e+9)).*9.378374045852051e-3-exp(d.*pi.^2.*(-6.136363636363636e+9)).*4.863340231002588e-2+exp(d.*pi.^2.*(-1.668099173553719e+12)).*9.289221334321387e-7+exp(d.*pi.^2.*(-1.563842975206611e+12)).*9.289221334321387e-7+exp(d.*pi.^2.*(-1.459586776859504e+12)).*9.289221334321387e-7+exp(d.*pi.^2.*(-1.44595041322314e+12)).*1.236275809772296e-6+exp(d.*pi.^2.*(-1.355578512396694e+12)).*1.236275809772296e-6+exp(d.*pi.^2.*(-1.355330578512397e+12)).*9.289221334321387e-7+exp(d.*pi.^2.*(-1.265206611570248e+12)).*1.236275809772296e-6+exp(d.*pi.^2.*(-1.251074380165289e+12)).*9.289221334321387e-7+exp(d.*pi.^2.*(-1.239669421487603e+12)).*1.681943238408239e-6+exp(d.*pi.^2.*(-1.174834710743802e+12)).*1.236275809772296e-6+exp(d.*pi.^2.*(-1.162190082644628e+12)).*1.681943238408239e-6+exp(d.*pi.^2.*(-1.084710743801653e+12)).*1.681943238408239e-6+exp(d.*pi.^2.*(-1.084462809917355e+12)).*1.236275809772296e-6+exp(d.*pi.^2.*(-1.049256198347107e+12)).*2.347788685302854e-6+exp(d.*pi.^2.*(-1.042561983471074e+12)).*9.289221334321387e-7+exp(d.*pi.^2.*(-1.007231404958678e+12)).*1.681943238408239e-6+exp(d.*pi.^2.*(-9.836776859504132e+11)).*2.347788685302854e-6+exp(d.*pi.^2.*(-9.383057851239669e+11)).*9.289221334321387e-7+exp(d.*pi.^2.*(-9.297520661157024e+11)).*1.681943238408239e-6+exp(d.*pi.^2.*(-9.18099173553719e+11)).*2.347788685302854e-6+exp(d.*pi.^2.*(-9.037190082644628e+11)).*1.236275809772296e-6+exp(d.*pi.^2.*(-8.747107438016528e+11)).*3.378269437248533e-6+exp(d.*pi.^2.*(-8.525206611570247e+11)).*2.347788685302854e-6-exp(d.*pi.^2.*(-8.340495867768595e+11)).*1.725561698001438e-4+exp(d.*pi.^2.*(-8.200413223140495e+11)).*3.378269437248533e-6+exp(d.*pi.^2.*(-8.133471074380165e+11)).*1.236275809772296e-6+exp(d.*pi.^2.*(-7.869421487603305e+11)).*2.347788685302854e-6-exp(d.*pi.^2.*(-7.819214876033057e+11)).*1.349328492816702e-4+exp(d.*pi.^2.*(-7.74793388429752e+11)).*1.681943238408239e-6+exp(d.*pi.^2.*(-7.653719008264462e+11)).*3.378269437248533e-6-exp(d.*pi.^2.*(-7.29793388429752e+11)).*3.07489019081814e-4-exp(d.*pi.^2.*(-7.229752066115702e+11)).*1.989020169695517e-4+exp(d.*pi.^2.*(-7.160330578512396e+11)).*5.041454287362399e-6+exp(d.*pi.^2.*(-7.107024793388429e+11)).*3.378269437248533e-6+exp(d.*pi.^2.*(-6.973140495867768e+11)).*1.681943238408239e-6-exp(d.*pi.^2.*(-6.777892561983471e+11)).*1.556631166061409e-4-exp(d.*pi.^2.*(-6.776652892561983e+11)).*3.662463051931047e-4+exp(d.*pi.^2.*(-6.712809917355372e+11)).*5.041454287362399e-6+exp(d.*pi.^2.*(-6.560330578512396e+11)).*3.378269437248533e-6+exp(d.*pi.^2.*(-6.557851239669421e+11)).*2.347788685302854e-6-exp(d.*pi.^2.*(-6.326033057851239e+11)).*3.545651335756926e-4+exp(d.*pi.^2.*(-6.265289256198347e+11)).*5.041454287362399e-6-exp(d.*pi.^2.*(-6.255371900826446e+11)).*3.845935043856255e-4-exp(d.*pi.^2.*(-6.198347107438016e+11)).*2.317596145965817e-4+exp(d.*pi.^2.*(-5.902066115702479e+11)).*2.347788685302854e-6-exp(d.*pi.^2.*(-5.874173553719008e+11)).*4.225141736452395e-4+exp(d.*pi.^2.*(-5.817768595041322e+11)).*5.041454287362399e-6-exp(d.*pi.^2.*(-5.81095041322314e+11)).*1.81565656093881e-4+exp(d.*pi.^2.*(-5.732231404958677e+11)).*7.866397979720083e-6+exp(d.*pi.^2.*(-5.46694214876033e+11)).*3.378269437248533e-6-exp(d.*pi.^2.*(-5.423553719008264e+11)).*4.133252706904627e-4-exp(d.*pi.^2.*(-5.422314049586777e+11)).*4.435154859220588e-4+exp(d.*pi.^2.*(-5.37396694214876e+11)).*7.866397979720083e-6+exp(d.*pi.^2.*(-5.370247933884297e+11)).*5.041454287362399e-6-exp(d.*pi.^2.*(-5.246280991735537e+11)).*2.734571304220482e-4-exp(d.*pi.^2.*(-5.212809917355372e+11)).*7.508398095787302e-4-exp(d.*pi.^2.*(-5.036157024793388e+11)).*4.928210665405343e-4+exp(d.*pi.^2.*(-5.015702479338843e+11)).*7.866397979720083e-6+exp(d.*pi.^2.*(-4.920247933884297e+11)).*3.378269437248533e-6-exp(d.*pi.^2.*(-4.918388429752066e+11)).*2.145149370834953e-4-exp(d.*pi.^2.*(-4.691528925619834e+11)).*7.132164890602566e-4+exp(d.*pi.^2.*(-4.657438016528925e+11)).*7.866397979720083e-6-exp(d.*pi.^2.*(-4.648760330578512e+11)).*5.170770741726805e-4-exp(d.*pi.^2.*(-4.590495867768595e+11)).*4.879720675055435e-4-exp(d.*pi.^2.*(-4.518595041322314e+11)).*8.660296595672983e-4+exp(d.*pi.^2.*(-4.475206611570248e+11)).*5.041454287362399e-6+exp(d.*pi.^2.*(-4.462809917355372e+11)).*1.297791811534305e-5-exp(d.*pi.^2.*(-4.373553719008264e+11)).*3.274629761275276e-4+exp(d.*pi.^2.*(-4.299173553719008e+11)).*7.866397979720083e-6-exp(d.*pi.^2.*(-4.262603305785124e+11)).*5.822548292266301e-4+exp(d.*pi.^2.*(-4.183884297520661e+11)).*1.297791811534305e-5-exp(d.*pi.^2.*(-4.170247933884297e+11)).*8.664965375344475e-4-exp(d.*pi.^2.*(-4.100206611570248e+11)).*2.573209687726036e-4-exp(d.*pi.^2.*(-4.066735537190082e+11)).*8.227907592038875e-4+exp(d.*pi.^2.*(-4.027685950413223e+11)).*5.041454287362399e-6-exp(d.*pi.^2.*(-3.934710743801653e+11)).*6.105520315532551e-4+exp(d.*pi.^2.*(-3.90495867768595e+11)).*1.297791811534305e-5-exp(d.*pi.^2.*(-3.87396694214876e+11)).*1.009898140713215e-3-exp(d.*pi.^2.*(-3.826859504132231e+11)).*5.847839449001312e-4-exp(d.*pi.^2.*(-3.64896694214876e+11)).*8.250179927507833e-4+exp(d.*pi.^2.*(-3.626033057851239e+11)).*1.297791811534305e-5-exp(d.*pi.^2.*(-3.614876033057851e+11)).*9.994551880868476e-4+exp(d.*pi.^2.*(-3.582644628099173e+11)).*7.866397979720083e-6-exp(d.*pi.^2.*(-3.580165289256198e+11)).*3.991158404761266e-4-exp(d.*pi.^2.*(-3.553512396694215e+11)).*6.984426295256384e-4-exp(d.*pi.^2.*(-3.486570247933884e+11)).*9.597041822105141e-4-exp(d.*pi.^2.*(-3.356404958677686e+11)).*3.143445625938248e-4+exp(d.*pi.^2.*(-3.352066115702479e+11)).*2.300373797956269e-5+exp(d.*pi.^2.*(-3.347107438016529e+11)).*1.297791811534305e-5-exp(d.*pi.^2.*(-3.280165289256198e+11)).*7.318244984844761e-4-exp(d.*pi.^2.*(-3.278925619834711e+11)).*1.192806860779885e-3+exp(d.*pi.^2.*(-3.224380165289256e+11)).*7.866397979720083e-6-exp(d.*pi.^2.*(-3.16301652892562e+11)).*9.517687701061185e-4+exp(d.*pi.^2.*(-3.142561983471074e+11)).*2.300373797956269e-5-exp(d.*pi.^2.*(-3.132644628099173e+11)).*7.134604030699514e-4-exp(d.*pi.^2.*(-3.127685950413223e+11)).*9.782980412249742e-4-exp(d.*pi.^2.*(-3.099173553719008e+11)).*1.165525845936541e-3-exp(d.*pi.^2.*(-2.95103305785124e+11)).*1.133864667441332e-3+exp(d.*pi.^2.*(-2.933057851239669e+11)).*2.300373797956269e-5-exp(d.*pi.^2.*(-2.908884297520661e+11)).*8.532209556118101e-4-exp(d.*pi.^2.*(-2.866115702479339e+11)).*4.969813954619574e-4+exp(d.*pi.^2.*(-2.789256198347107e+11)).*1.297791811534305e-5-exp(d.*pi.^2.*(-2.733471074380165e+11)).*1.430267128010115e-3+exp(d.*pi.^2.*(-2.723553719008264e+11)).*2.300373797956269e-5-exp(d.*pi.^2.*(-2.711776859504132e+11)).*1.11014429725973e-3-exp(d.*pi.^2.*(-2.711157024793388e+11)).*1.128433198989079e-3-exp(d.*pi.^2.*(-2.68698347107438e+11)).*3.926593948990825e-4-exp(d.*pi.^2.*(-2.685123966942149e+11)).*8.930858674092798e-4-exp(d.*pi.^2.*(-2.623140495867768e+11)).*1.376676806851453e-3-exp(d.*pi.^2.*(-2.606404958677686e+11)).*1.116087424772672e-3+exp(d.*pi.^2.*(-2.514049586776859e+11)).*2.300373797956269e-5+exp(d.*pi.^2.*(-2.510330578512397e+11)).*1.297791811534305e-5-exp(d.*pi.^2.*(-2.507851239669421e+11)).*8.896407903610399e-4-exp(d.*pi.^2.*(-2.460123966942149e+11)).*1.360125120655191e-3-exp(d.*pi.^2.*(-2.328719008264463e+11)).*1.065789786154653e-3-exp(d.*pi.^2.*(-2.324380165289256e+11)).*1.315965960985757e-3-exp(d.*pi.^2.*(-2.295247933884297e+11)).*1.311605615310514e-3-exp(d.*pi.^2.*(-2.259297520661157e+11)).*1.287556350213651e-3-exp(d.*pi.^2.*(-2.237603305785124e+11)).*1.74630682302109e-3-exp(d.*pi.^2.*(-2.231404958677686e+11)).*6.354698800435486e-4-exp(d.*pi.^2.*(-2.186776859504132e+11)).*1.650827958386632e-3-exp(d.*pi.^2.*(-2.149586776859504e+11)).*1.11401758744623e-3+exp(d.*pi.^2.*(-2.095041322314049e+11)).*2.300373797956269e-5-exp(d.*pi.^2.*(-2.09194214876033e+11)).*5.043482874569158e-4-exp(d.*pi.^2.*(-2.085123966942149e+11)).*1.313702552296554e-3-exp(d.*pi.^2.*(-2.013842975206611e+11)).*1.661535545138788e-3-exp(d.*pi.^2.*(-1.967355371900826e+11)).*1.554417754720634e-3-exp(d.*pi.^2.*(-1.952479338842975e+11)).*1.139818167500464e-3-exp(d.*pi.^2.*(-1.93698347107438e+11)).*1.501807355405102e-3-exp(d.*pi.^2.*(-1.913429752066116e+11)).*1.573333923352491e-3+exp(d.*pi.^2.*(-1.885537190082645e+11)).*2.300373797956269e-5-exp(d.*pi.^2.*(-1.81301652892562e+11)).*1.368945351668771e-3-exp(d.*pi.^2.*(-1.807438016528926e+11)).*1.515367231695772e-3-exp(d.*pi.^2.*(-1.791322314049587e+11)).*2.179807373600882e-3-exp(d.*pi.^2.*(-1.790082644628099e+11)).*2.015745019530083e-3-exp(d.*pi.^2.*(-1.67603305785124e+11)).*8.403160814024732e-4-exp(d.*pi.^2.*(-1.673553719008264e+11)).*1.428017188904416e-3-exp(d.*pi.^2.*(-1.640082644628099e+11)).*1.864036761083932e-3-exp(d.*pi.^2.*(-1.639462809917355e+11)).*1.774344979590625e-3-exp(d.*pi.^2.*(-1.612190082644628e+11)).*2.075485373038008e-3+exp(d.*pi.^2.*(-1.606611570247934e+11)).*1.001387582346529e-4-exp(d.*pi.^2.*(-1.571280991735537e+11)).*6.714709706304724e-4-exp(d.*pi.^2.*(-1.566322314049587e+11)).*1.921992468430814e-3-exp(d.*pi.^2.*(-1.563842975206611e+11)).*1.387880735468607e-3-exp(d.*pi.^2.*(-1.549586776859504e+11)).*1.767286306133404e-3+exp(d.*pi.^2.*(-1.506198347107438e+11)).*1.001387582346529e-4-exp(d.*pi.^2.*(-1.466528925619835e+11)).*1.511787052032946e-3-exp(d.*pi.^2.*(-1.433057851239669e+11)).*2.516372569228667e-3+exp(d.*pi.^2.*(-1.405785123966942e+11)).*1.001387582346529e-4-exp(d.*pi.^2.*(-1.394628099173554e+11)).*2.796962540573188e-3-exp(d.*pi.^2.*(-1.366735537190083e+11)).*2.128412013133393e-3-exp(d.*pi.^2.*(-1.361776859504132e+11)).*1.822564063139854e-3-exp(d.*pi.^2.*(-1.355578512396694e+11)).*1.601106342234592e-3-exp(d.*pi.^2.*(-1.342561983471074e+11)).*2.276201942822109e-3-exp(d.*pi.^2.*(-1.311570247933884e+11)).*2.08764059832818e-3+exp(d.*pi.^2.*(-1.305371900826446e+11)).*1.001387582346529e-4-exp(d.*pi.^2.*(-1.25702479338843e+11)).*1.895484749536073e-3-exp(d.*pi.^2.*(-1.255165289256198e+11)).*2.665840947986555e-3-exp(d.*pi.^2.*(-1.253925619834711e+11)).*2.400831728811533e-3+exp(d.*pi.^2.*(-1.20495867768595e+11)).*1.001387582346529e-4-exp(d.*pi.^2.*(-1.162190082644628e+11)).*1.867532462679919e-3-exp(d.*pi.^2.*(-1.118801652892562e+11)).*2.600078596311779e-3-exp(d.*pi.^2.*(-1.115702479338843e+11)).*3.229261072679115e-3-exp(d.*pi.^2.*(-1.093388429752066e+11)).*2.503663169175833e-3-exp(d.*pi.^2.*(-1.074793388429752e+11)).*2.841718925002193e-3-exp(d.*pi.^2.*(-1.047520661157025e+11)).*3.718048812675926e-3-exp(d.*pi.^2.*(-1.042561983471074e+11)).*1.541160783942798e-3+exp(d.*pi.^2.*(-1.004132231404959e+11)).*1.001387582346529e-4-exp(d.*pi.^2.*(-9.836776859504132e+10)).*2.206439352858809e-3-exp(d.*pi.^2.*(-9.762396694214876e+10)).*3.083729529022285e-3+exp(d.*pi.^2.*(-9.719008264462809e+10)).*2.736383570858638e-4-exp(d.*pi.^2.*(-9.427685950413223e+10)).*3.549203701903925e-3+exp(d.*pi.^2.*(-9.111570247933884e+10)).*2.736383570858638e-4-exp(d.*pi.^2.*(-9.037190082644628e+10)).*1.677632012882899e-3-exp(d.*pi.^2.*(-8.956611570247933e+10)).*3.247854137808125e-3-exp(d.*pi.^2.*(-8.950413223140495e+10)).*3.057572712698188e-3+exp(d.*pi.^2.*(-8.504132231404958e+10)).*2.736383570858638e-4-exp(d.*pi.^2.*(-8.380165289256198e+10)).*4.293595358930617e-3-exp(d.*pi.^2.*(-8.367768595041322e+10)).*3.647149653714846e-3-exp(d.*pi.^2.*(-8.200413223140495e+10)).*2.646729964518209e-3-exp(d.*pi.^2.*(-8.033057851239669e+10)).*1.701109632965949e-3+exp(d.*pi.^2.*(-7.896694214876033e+10)).*2.736383570858638e-4-exp(d.*pi.^2.*(-7.74793388429752e+10)).*2.073354126405947e-3-exp(d.*pi.^2.*(-7.53099173553719e+10)).*1.400970970933801e-3-exp(d.*pi.^2.*(-7.332644628099173e+10)).*4.10556536328346e-3+exp(d.*pi.^2.*(-7.289256198347107e+10)).*2.736383570858638e-4-exp(d.*pi.^2.*(-7.165289256198347e+10)).*3.817757992322769e-3-exp(d.*pi.^2.*(-7.02892561983471e+10)).*3.10208060389975e-3-exp(d.*pi.^2.*(-6.973140495867768e+10)).*4.171680834822203e-3-exp(d.*pi.^2.*(-6.712809917355372e+10)).*3.233258358107912e-3-exp(d.*pi.^2.*(-6.557851239669421e+10)).*2.449251492268929e-3-exp(d.*pi.^2.*(-6.526859504132231e+10)).*3.802635492534603e-3-exp(d.*pi.^2.*(-6.285123966942148e+10)).*4.849957020310151e-3+exp(d.*pi.^2.*(-6.074380165289256e+10)).*2.736383570858638e-4-exp(d.*pi.^2.*(-6.024793388429752e+10)).*3.90263544443335e-3-exp(d.*pi.^2.*(-5.578512396694215e+10)).*4.900815396822036e-3-exp(d.*pi.^2.*(-5.46694214876033e+10)).*2.663794445163786e-3-exp(d.*pi.^2.*(-5.37396694214876e+10)).*4.03878234753342e-3-exp(d.*pi.^2.*(-5.237603305785124e+10)).*5.554024171357764e-3-exp(d.*pi.^2.*(-5.212809917355372e+10)).*1.349328492816702e-3-exp(d.*pi.^2.*(-5.020661157024793e+10)).*7.705270936967953e-3+exp(d.*pi.^2.*(-4.958677685950413e+10)).*1.051219258989808e-3-exp(d.*pi.^2.*(-4.859504132231405e+10)).*2.703924869712517e-3+exp(d.*pi.^2.*(-4.648760330578512e+10)).*1.051219258989808e-3-exp(d.*pi.^2.*(-4.555785123966942e+10)).*2.315882509732074e-3-exp(d.*pi.^2.*(-4.518595041322314e+10)).*8.961763440997215e-3-exp(d.*pi.^2.*(-4.475206611570248e+10)).*3.587467832499206e-3+exp(d.*pi.^2.*(-4.338842975206611e+10)).*1.051219258989808e-3-exp(d.*pi.^2.*(-4.252066115702479e+10)).*5.019807379444591e-3-exp(d.*pi.^2.*(-4.190082644628099e+10)).*6.519042004448754e-3-exp(d.*pi.^2.*(-4.183884297520661e+10)).*5.187582385271133e-3+exp(d.*pi.^2.*(-4.028925619834711e+10)).*1.051219258989808e-3-exp(d.*pi.^2.*(-4.016528925619835e+10)).*8.906103197768354e-3-exp(d.*pi.^2.*(-3.948347107438016e+10)).*6.285966812129915e-3-exp(d.*pi.^2.*(-3.87396694214876e+10)).*1.81565656093881e-3+exp(d.*pi.^2.*(-3.71900826446281e+10)).*1.051219258989808e-3-exp(d.*pi.^2.*(-3.644628099173554e+10)).*6.343168813577205e-3-exp(d.*pi.^2.*(-3.582644628099173e+10)).*4.47966954372408e-3-exp(d.*pi.^2.*(-3.514462809917355e+10)).*8.565936793709527e-3-exp(d.*pi.^2.*(-3.278925619834711e+10)).*2.145149370834953e-3-exp(d.*pi.^2.*(-3.142561983471074e+10)).*6.906558555056287e-3+exp(d.*pi.^2.*(-3.099173553719008e+10)).*1.051219258989808e-3-exp(d.*pi.^2.*(-3.037190082644628e+10)).*1.262913562570712e-2-exp(d.*pi.^2.*(-3.012396694214876e+10)).*1.006690771654208e-2-exp(d.*pi.^2.*(-2.789256198347107e+10)).*4.699783250973886e-3-exp(d.*pi.^2.*(-2.733471074380165e+10)).*1.481430295345271e-2-exp(d.*pi.^2.*(-2.510330578512397e+10)).*1.158803131672387e-2-exp(d.*pi.^2.*(-2.479338842975206e+10)).*4.784832830502693e-3-exp(d.*pi.^2.*(-2.429752066115702e+10)).*1.461417777690604e-2-exp(d.*pi.^2.*(-2.324380165289256e+10)).*4.539151625160834e-3-exp(d.*pi.^2.*(-2.237603305785124e+10)).*3.143445625938248e-3-exp(d.*pi.^2.*(-2.169421487603306e+10)).*9.323984455663526e-3-exp(d.*pi.^2.*(-2.12603305785124e+10)).*1.415996734521897e-2-exp(d.*pi.^2.*(-2.095041322314049e+10)).*7.650950212082979e-3-exp(d.*pi.^2.*(-2.014462809917355e+10)).*1.232055441115083e-2-exp(d.*pi.^2.*(-2.008264462809917e+10)).*1.354932127286324e-2-exp(d.*pi.^2.*(-1.859504132231405e+10)).*1.191778538432686e-2-exp(d.*pi.^2.*(-1.822314049586777e+10)).*1.653305185639833e-2-exp(d.*pi.^2.*(-1.791322314049587e+10)).*3.926593948990825e-3+exp(d.*pi.^2.*(-1.785123966942149e+10)).*8.111137929522739e-3+exp(d.*pi.^2.*(-1.673553719008264e+10)).*8.111137929522739e-3+exp(d.*pi.^2.*(-1.56198347107438e+10)).*8.111137929522739e-3-exp(d.*pi.^2.*(-1.549586776859504e+10)).*2.423833979547769e-2-exp(d.*pi.^2.*(-1.518595041322314e+10)).*1.915565675906958e-2-exp(d.*pi.^2.*(-1.506198347107438e+10)).*1.440998712960481e-2+exp(d.*pi.^2.*(-1.450413223140496e+10)).*8.111137929522739e-3-exp(d.*pi.^2.*(-1.394628099173554e+10)).*2.903614146470499e-2+exp(d.*pi.^2.*(-1.338842975206612e+10)).*8.111137929522739e-3-exp(d.*pi.^2.*(-1.239669421487603e+10)).*2.81290411884727e-2-exp(d.*pi.^2.*(-1.214876033057851e+10)).*2.22896740948752e-2+exp(d.*pi.^2.*(-1.115702479338843e+10)).*8.111137929522739e-3-exp(d.*pi.^2.*(-1.084710743801653e+10)).*2.775366993669767e-2-exp(d.*pi.^2.*(-1.047520661157025e+10)).*6.714709706304724e-3-exp(d.*pi.^2.*(-1.004132231404959e+10)).*7.799820122914621e-3-exp(d.*pi.^2.*(-9.297520661157024e+9)).*3.189005253503453e-2-exp(d.*pi.^2.*(-9.111570247933884e+9)).*2.382050581438705e-2-exp(d.*pi.^2.*(-8.925619834710743e+9)).*8.099996173819221e-3-exp(d.*pi.^2.*(-8.367768595041322e+9)).*1.260865985815486e-2-exp(d.*pi.^2.*(-7.8099173553719e+9)).*2.070865603197408e-2-exp(d.*pi.^2.*(-7.74793388429752e+9)).*3.754526844240175e-2-exp(d.*pi.^2.*(-7.252066115702479e+9)).*3.422350532927747e-2-exp(d.*pi.^2.*(-6.694214876033058e+9)).*2.791360452234828e-2-exp(d.*pi.^2.*(-6.198347107438016e+9)).*4.317308657472003e-2-exp(d.*pi.^2.*(-6.074380165289256e+9)).*2.619359032556641e-2-exp(d.*pi.^2.*(-5.578512396694215e+9)).*6.213710985162576e-2-exp(d.*pi.^2.*(-5.020661157024793e+9)).*8.065548324529941e-2-exp(d.*pi.^2.*(-4.648760330578512e+9)).*4.668841671594001e-2-exp(d.*pi.^2.*(-4.462809917355372e+9)).*7.294453258718706e-2-exp(d.*pi.^2.*(-3.90495867768595e+9)).*7.709294884700399e-2-exp(d.*pi.^2.*(-3.347107438016529e+9)).*8.339170789822966e-2-exp(d.*pi.^2.*(-3.099173553719008e+9)).*5.082479931427686e-2-exp(d.*pi.^2.*(-3.037190082644628e+9)).*2.315882509732074e-2-exp(d.*pi.^2.*(-2.789256198347107e+9)).*1.042916293981666e-1-exp(d.*pi.^2.*(-2.231404958677686e+9)).*1.147332338313575e-1-exp(d.*pi.^2.*(-1.673553719008264e+9)).*1.296890728267357e-1-exp(d.*pi.^2.*(-1.549586776859504e+9)).*4.539151625160834e-2-exp(d.*pi.^2.*(-1.115702479338843e+9)).*1.359878318779614e-1-exp(d.*pi.^2.*(-5.578512396694215e+8)).*1.260865985815486e-1+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-7.229752066115702e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-6.198347107438016e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-5.246280991735537e+11)).*2.953582823893289e-6+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-4.373553719008264e+11)).*3.542964438408886e-6+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-3.580165289256198e+11)).*4.328102804794397e-6+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-2.866115702479339e+11)).*5.406392954178542e-6+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-2.231404958677686e+11)).*6.944199128254334e-6+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-1.67603305785124e+11)).*9.245254211948807e-6+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-8.033057851239669e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-4.859504132231405e+10)).*3.188659445303408e-5+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-2.479338842975206e+10)).*6.249802674535459e-5+exp(d.*pi.^2.*(-8.340495867768595e+11)).*exp(d.*pi.^2.*(-8.925619834710743e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-6.777892561983471e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-5.81095041322314e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-4.918388429752066e+11)).*2.953582823893289e-6+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-4.100206611570248e+11)).*3.542964438408886e-6+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-3.356404958677686e+11)).*4.328102804794397e-6+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-2.68698347107438e+11)).*5.406392954178542e-6+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-2.09194214876033e+11)).*6.944199128254334e-6+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-1.571280991735537e+11)).*9.245254211948807e-6+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-7.53099173553719e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-4.555785123966942e+10)).*3.188659445303408e-5+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-2.324380165289256e+10)).*6.249802674535459e-5+exp(d.*pi.^2.*(-7.819214876033057e+11)).*exp(d.*pi.^2.*(-8.367768595041322e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-6.326033057851239e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-5.423553719008264e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-4.590495867768595e+11)).*2.953582823893289e-6+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-3.826859504132231e+11)).*3.542964438408886e-6+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-3.132644628099173e+11)).*4.328102804794397e-6+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-2.507851239669421e+11)).*5.406392954178542e-6+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-1.952479338842975e+11)).*6.944199128254334e-6+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-1.466528925619835e+11)).*9.245254211948807e-6+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-7.02892561983471e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-4.252066115702479e+10)).*3.188659445303408e-5+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-2.169421487603306e+10)).*6.249802674535459e-5+exp(d.*pi.^2.*(-7.29793388429752e+11)).*exp(d.*pi.^2.*(-7.8099173553719e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-6.198347107438016e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-5.246280991735537e+11)).*3.407353435202765e-6+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-4.373553719008264e+11)).*4.08728407829132e-6+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-3.580165289256198e+11)).*4.993046357300908e-6+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-2.866115702479339e+11)).*6.236998487211489e-6+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-2.231404958677686e+11)).*8.01106390617483e-6+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-1.67603305785124e+11)).*1.066563918356019e-5+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-8.033057851239669e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-4.859504132231405e+10)).*3.678545807740265e-5+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-2.479338842975206e+10)).*7.209984578779232e-5+exp(d.*pi.^2.*(-7.229752066115702e+11)).*exp(d.*pi.^2.*(-8.925619834710743e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-5.81095041322314e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-4.918388429752066e+11)).*3.407353435202765e-6+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-4.100206611570248e+11)).*4.08728407829132e-6+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-3.356404958677686e+11)).*4.993046357300908e-6+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-2.68698347107438e+11)).*6.236998487211489e-6+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-2.09194214876033e+11)).*8.01106390617483e-6+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-1.571280991735537e+11)).*1.066563918356019e-5+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-7.53099173553719e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-4.555785123966942e+10)).*3.678545807740265e-5+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-2.324380165289256e+10)).*7.209984578779232e-5+exp(d.*pi.^2.*(-6.777892561983471e+11)).*exp(d.*pi.^2.*(-8.367768595041322e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-5.874173553719008e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-5.036157024793388e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-4.262603305785124e+11)).*2.953582823893289e-6+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-3.553512396694215e+11)).*3.542964438408886e-6+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-2.908884297520661e+11)).*4.328102804794397e-6+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-2.328719008264463e+11)).*5.406392954178542e-6+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-1.81301652892562e+11)).*6.944199128254334e-6+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-1.361776859504132e+11)).*9.245254211948807e-6+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-6.526859504132231e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-3.948347107438016e+10)).*3.188659445303408e-5+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-2.014462809917355e+10)).*6.249802674535459e-5+exp(d.*pi.^2.*(-6.776652892561983e+11)).*exp(d.*pi.^2.*(-7.252066115702479e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-5.423553719008264e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-4.590495867768595e+11)).*3.407353435202765e-6+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-3.826859504132231e+11)).*4.08728407829132e-6+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-3.132644628099173e+11)).*4.993046357300908e-6+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-2.507851239669421e+11)).*6.236998487211489e-6+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-1.952479338842975e+11)).*8.01106390617483e-6+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-1.466528925619835e+11)).*1.066563918356019e-5+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-7.02892561983471e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-4.252066115702479e+10)).*3.678545807740265e-5+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-2.169421487603306e+10)).*7.209984578779232e-5+exp(d.*pi.^2.*(-6.326033057851239e+11)).*exp(d.*pi.^2.*(-7.8099173553719e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-5.422314049586777e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-4.648760330578512e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-3.934710743801653e+11)).*2.953582823893289e-6+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-3.280165289256198e+11)).*3.542964438408886e-6+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-2.685123966942149e+11)).*4.328102804794397e-6+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-2.149586776859504e+11)).*5.406392954178542e-6+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-1.673553719008264e+11)).*6.944199128254334e-6+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-1.25702479338843e+11)).*9.245254211948807e-6+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-6.024793388429752e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-3.644628099173554e+10)).*3.188659445303408e-5+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-1.859504132231405e+10)).*6.249802674535459e-5+exp(d.*pi.^2.*(-6.255371900826446e+11)).*exp(d.*pi.^2.*(-6.694214876033058e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-5.246280991735537e+11)).*3.974341356479841e-6+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-4.373553719008264e+11)).*4.767413318562332e-6+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-3.580165289256198e+11)).*5.823895586417536e-6+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-2.866115702479339e+11)).*7.274842924110025e-6+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-2.231404958677686e+11)).*9.344115072646962e-6+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-1.67603305785124e+11)).*1.244041503372621e-5+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-8.033057851239669e+10)).*2.595594015368134e-5+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-4.859504132231405e+10)).*4.290660482814774e-5+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-2.479338842975206e+10)).*8.409735131958502e-5+exp(d.*pi.^2.*(-6.198347107438016e+11)).*exp(d.*pi.^2.*(-8.925619834710743e+9)).*2.336019999602518e-4+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-5.036157024793388e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-4.262603305785124e+11)).*3.407353435202765e-6+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-3.553512396694215e+11)).*4.08728407829132e-6+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-2.908884297520661e+11)).*4.993046357300908e-6+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-2.328719008264463e+11)).*6.236998487211489e-6+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-1.81301652892562e+11)).*8.01106390617483e-6+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-1.361776859504132e+11)).*1.066563918356019e-5+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-6.526859504132231e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-3.948347107438016e+10)).*3.678545807740265e-5+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-2.014462809917355e+10)).*7.209984578779232e-5+exp(d.*pi.^2.*(-5.874173553719008e+11)).*exp(d.*pi.^2.*(-7.252066115702479e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-4.918388429752066e+11)).*3.974341356479841e-6+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-4.100206611570248e+11)).*4.767413318562332e-6+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-3.356404958677686e+11)).*5.823895586417536e-6+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-2.68698347107438e+11)).*7.274842924110025e-6+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-2.09194214876033e+11)).*9.344115072646962e-6+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-1.571280991735537e+11)).*1.244041503372621e-5+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-7.53099173553719e+10)).*2.595594015368134e-5+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-4.555785123966942e+10)).*4.290660482814774e-5+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-2.324380165289256e+10)).*8.409735131958502e-5+exp(d.*pi.^2.*(-5.81095041322314e+11)).*exp(d.*pi.^2.*(-8.367768595041322e+9)).*2.336019999602518e-4+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-4.590495867768595e+11)).*3.974341356479841e-6+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-3.826859504132231e+11)).*4.767413318562332e-6+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-3.132644628099173e+11)).*5.823895586417536e-6+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-2.507851239669421e+11)).*7.274842924110025e-6+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-1.952479338842975e+11)).*9.344115072646962e-6+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-1.466528925619835e+11)).*1.244041503372621e-5+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-7.02892561983471e+10)).*2.595594015368134e-5+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-4.252066115702479e+10)).*4.290660482814774e-5+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-2.169421487603306e+10)).*8.409735131958502e-5+exp(d.*pi.^2.*(-5.423553719008264e+11)).*exp(d.*pi.^2.*(-7.8099173553719e+9)).*2.336019999602518e-4+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-4.648760330578512e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-3.934710743801653e+11)).*3.407353435202765e-6+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-3.280165289256198e+11)).*4.08728407829132e-6+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-2.685123966942149e+11)).*4.993046357300908e-6+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-2.149586776859504e+11)).*6.236998487211489e-6+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-1.673553719008264e+11)).*8.01106390617483e-6+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-1.25702479338843e+11)).*1.066563918356019e-5+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-6.024793388429752e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-3.644628099173554e+10)).*3.678545807740265e-5+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-1.859504132231405e+10)).*7.209984578779232e-5+exp(d.*pi.^2.*(-5.422314049586777e+11)).*exp(d.*pi.^2.*(-6.694214876033058e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-5.246280991735537e+11)).*exp(d.*pi.^2.*(-4.373553719008264e+11)).*5.632570553726442e-6+exp(d.*pi.^2.*(-5.246280991735537e+11)).*exp(d.*pi.^2.*(-3.580165289256198e+11)).*6.88077592524012e-6+exp(d.*pi.^2.*(-5.246280991735537e+11)).*exp(d.*pi.^2.*(-2.866115702479339e+11)).*8.595031162451022e-6+exp(d.*pi.^2.*(-5.246280991735537e+11)).*exp(d.*pi.^2.*(-2.231404958677686e+11)).*1.103982052571315e-5+exp(d.*pi.^2.*(-5.246280991735537e+11)).*exp(d.*pi.^2.*(-1.67603305785124e+11)).*1.469801561410095e-5+exp(d.*pi.^2.*(-5.246280991735537e+11)).*exp(d.*pi.^2.*(-8.033057851239669e+10)).*3.066624486588446e-5+exp(d.*pi.^2.*(-5.246280991735537e+11)).*exp(d.*pi.^2.*(-4.859504132231405e+10)).*5.069299906815745e-5+exp(d.*pi.^2.*(-5.246280991735537e+11)).*exp(d.*pi.^2.*(-2.479338842975206e+10)).*9.935875768202262e-5+exp(d.*pi.^2.*(-5.246280991735537e+11)).*exp(d.*pi.^2.*(-8.925619834710743e+9)).*2.759944771611512e-4+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-4.518595041322314e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-3.87396694214876e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-3.278925619834711e+11)).*2.953582823893289e-6+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-2.733471074380165e+11)).*3.542964438408886e-6+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-2.237603305785124e+11)).*4.328102804794397e-6+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-1.791322314049587e+11)).*5.406392954178542e-6+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-1.394628099173554e+11)).*6.944199128254334e-6+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-1.047520661157025e+11)).*9.245254211948807e-6+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-5.020661157024793e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-3.037190082644628e+10)).*3.188659445303408e-5+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-1.549586776859504e+10)).*6.249802674535459e-5+exp(d.*pi.^2.*(-5.212809917355372e+11)).*exp(d.*pi.^2.*(-5.578512396694215e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-4.262603305785124e+11)).*3.974341356479841e-6+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-3.553512396694215e+11)).*4.767413318562332e-6+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-2.908884297520661e+11)).*5.823895586417536e-6+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-2.328719008264463e+11)).*7.274842924110025e-6+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-1.81301652892562e+11)).*9.344115072646962e-6+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-1.361776859504132e+11)).*1.244041503372621e-5+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-6.526859504132231e+10)).*2.595594015368134e-5+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-3.948347107438016e+10)).*4.290660482814774e-5+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-2.014462809917355e+10)).*8.409735131958502e-5+exp(d.*pi.^2.*(-5.036157024793388e+11)).*exp(d.*pi.^2.*(-7.252066115702479e+9)).*2.336019999602518e-4+exp(d.*pi.^2.*(-4.918388429752066e+11)).*exp(d.*pi.^2.*(-4.100206611570248e+11)).*5.632570553726442e-6+exp(d.*pi.^2.*(-4.918388429752066e+11)).*exp(d.*pi.^2.*(-3.356404958677686e+11)).*6.88077592524012e-6+exp(d.*pi.^2.*(-4.918388429752066e+11)).*exp(d.*pi.^2.*(-2.68698347107438e+11)).*8.595031162451022e-6+exp(d.*pi.^2.*(-4.918388429752066e+11)).*exp(d.*pi.^2.*(-2.09194214876033e+11)).*1.103982052571315e-5+exp(d.*pi.^2.*(-4.918388429752066e+11)).*exp(d.*pi.^2.*(-1.571280991735537e+11)).*1.469801561410095e-5+exp(d.*pi.^2.*(-4.918388429752066e+11)).*exp(d.*pi.^2.*(-7.53099173553719e+10)).*3.066624486588446e-5+exp(d.*pi.^2.*(-4.918388429752066e+11)).*exp(d.*pi.^2.*(-4.555785123966942e+10)).*5.069299906815745e-5+exp(d.*pi.^2.*(-4.918388429752066e+11)).*exp(d.*pi.^2.*(-2.324380165289256e+10)).*9.935875768202262e-5+exp(d.*pi.^2.*(-4.918388429752066e+11)).*exp(d.*pi.^2.*(-8.367768595041322e+9)).*2.759944771611512e-4+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-4.066735537190082e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-3.486570247933884e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-2.95103305785124e+11)).*2.953582823893289e-6+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-2.460123966942149e+11)).*3.542964438408886e-6+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-2.013842975206611e+11)).*4.328102804794397e-6+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-1.612190082644628e+11)).*5.406392954178542e-6+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-1.255165289256198e+11)).*6.944199128254334e-6+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-9.427685950413223e+10)).*9.245254211948807e-6+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-4.518595041322314e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-2.733471074380165e+10)).*3.188659445303408e-5+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*6.249802674535459e-5+exp(d.*pi.^2.*(-4.691528925619834e+11)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-3.934710743801653e+11)).*3.974341356479841e-6+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-3.280165289256198e+11)).*4.767413318562332e-6+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-2.685123966942149e+11)).*5.823895586417536e-6+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-2.149586776859504e+11)).*7.274842924110025e-6+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-1.673553719008264e+11)).*9.344115072646962e-6+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-1.25702479338843e+11)).*1.244041503372621e-5+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-6.024793388429752e+10)).*2.595594015368134e-5+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-3.644628099173554e+10)).*4.290660482814774e-5+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-1.859504132231405e+10)).*8.409735131958502e-5+exp(d.*pi.^2.*(-4.648760330578512e+11)).*exp(d.*pi.^2.*(-6.694214876033058e+9)).*2.336019999602518e-4+exp(d.*pi.^2.*(-4.590495867768595e+11)).*exp(d.*pi.^2.*(-3.826859504132231e+11)).*5.632570553726442e-6+exp(d.*pi.^2.*(-4.590495867768595e+11)).*exp(d.*pi.^2.*(-3.132644628099173e+11)).*6.88077592524012e-6+exp(d.*pi.^2.*(-4.590495867768595e+11)).*exp(d.*pi.^2.*(-2.507851239669421e+11)).*8.595031162451022e-6+exp(d.*pi.^2.*(-4.590495867768595e+11)).*exp(d.*pi.^2.*(-1.952479338842975e+11)).*1.103982052571315e-5+exp(d.*pi.^2.*(-4.590495867768595e+11)).*exp(d.*pi.^2.*(-1.466528925619835e+11)).*1.469801561410095e-5+exp(d.*pi.^2.*(-4.590495867768595e+11)).*exp(d.*pi.^2.*(-7.02892561983471e+10)).*3.066624486588446e-5+exp(d.*pi.^2.*(-4.590495867768595e+11)).*exp(d.*pi.^2.*(-4.252066115702479e+10)).*5.069299906815745e-5+exp(d.*pi.^2.*(-4.590495867768595e+11)).*exp(d.*pi.^2.*(-2.169421487603306e+10)).*9.935875768202262e-5+exp(d.*pi.^2.*(-4.590495867768595e+11)).*exp(d.*pi.^2.*(-7.8099173553719e+9)).*2.759944771611512e-4+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-3.87396694214876e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-3.278925619834711e+11)).*3.407353435202765e-6+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-2.733471074380165e+11)).*4.08728407829132e-6+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-2.237603305785124e+11)).*4.993046357300908e-6+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-1.791322314049587e+11)).*6.236998487211489e-6+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-1.394628099173554e+11)).*8.01106390617483e-6+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-1.047520661157025e+11)).*1.066563918356019e-5+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-5.020661157024793e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-3.037190082644628e+10)).*3.678545807740265e-5+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-1.549586776859504e+10)).*7.209984578779232e-5+exp(d.*pi.^2.*(-4.518595041322314e+11)).*exp(d.*pi.^2.*(-5.578512396694215e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-4.373553719008264e+11)).*exp(d.*pi.^2.*(-3.580165289256198e+11)).*8.253821160718708e-6+exp(d.*pi.^2.*(-4.373553719008264e+11)).*exp(d.*pi.^2.*(-2.866115702479339e+11)).*1.031015264215268e-5+exp(d.*pi.^2.*(-4.373553719008264e+11)).*exp(d.*pi.^2.*(-2.231404958677686e+11)).*1.324279489053236e-5+exp(d.*pi.^2.*(-4.373553719008264e+11)).*exp(d.*pi.^2.*(-1.67603305785124e+11)).*1.763097557809323e-5+exp(d.*pi.^2.*(-4.373553719008264e+11)).*exp(d.*pi.^2.*(-8.033057851239669e+10)).*3.678563341458988e-5+exp(d.*pi.^2.*(-4.373553719008264e+11)).*exp(d.*pi.^2.*(-4.859504132231405e+10)).*6.080868683344755e-5+exp(d.*pi.^2.*(-4.373553719008264e+11)).*exp(d.*pi.^2.*(-2.479338842975206e+10)).*1.191856013869514e-4+exp(d.*pi.^2.*(-4.373553719008264e+11)).*exp(d.*pi.^2.*(-8.925619834710743e+9)).*3.310686295535354e-4+exp(d.*pi.^2.*(-4.262603305785124e+11)).*exp(d.*pi.^2.*(-3.553512396694215e+11)).*5.632570553726442e-6+exp(d.*pi.^2.*(-4.262603305785124e+11)).*exp(d.*pi.^2.*(-2.908884297520661e+11)).*6.88077592524012e-6+exp(d.*pi.^2.*(-4.262603305785124e+11)).*exp(d.*pi.^2.*(-2.328719008264463e+11)).*8.595031162451022e-6+exp(d.*pi.^2.*(-4.262603305785124e+11)).*exp(d.*pi.^2.*(-1.81301652892562e+11)).*1.103982052571315e-5+exp(d.*pi.^2.*(-4.262603305785124e+11)).*exp(d.*pi.^2.*(-1.361776859504132e+11)).*1.469801561410095e-5+exp(d.*pi.^2.*(-4.262603305785124e+11)).*exp(d.*pi.^2.*(-6.526859504132231e+10)).*3.066624486588446e-5+exp(d.*pi.^2.*(-4.262603305785124e+11)).*exp(d.*pi.^2.*(-3.948347107438016e+10)).*5.069299906815745e-5+exp(d.*pi.^2.*(-4.262603305785124e+11)).*exp(d.*pi.^2.*(-2.014462809917355e+10)).*9.935875768202262e-5+exp(d.*pi.^2.*(-4.262603305785124e+11)).*exp(d.*pi.^2.*(-7.252066115702479e+9)).*2.759944771611512e-4+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-3.614876033057851e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-3.099173553719008e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-2.623140495867768e+11)).*2.953582823893289e-6+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-2.186776859504132e+11)).*3.542964438408886e-6+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-1.790082644628099e+11)).*4.328102804794397e-6+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-1.433057851239669e+11)).*5.406392954178542e-6+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-1.115702479338843e+11)).*6.944199128254334e-6+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-8.380165289256198e+10)).*9.245254211948807e-6+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-4.016528925619835e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-2.429752066115702e+10)).*3.188659445303408e-5+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-1.239669421487603e+10)).*6.249802674535459e-5+exp(d.*pi.^2.*(-4.170247933884297e+11)).*exp(d.*pi.^2.*(-4.462809917355372e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-4.100206611570248e+11)).*exp(d.*pi.^2.*(-3.356404958677686e+11)).*8.253821160718708e-6+exp(d.*pi.^2.*(-4.100206611570248e+11)).*exp(d.*pi.^2.*(-2.68698347107438e+11)).*1.031015264215268e-5+exp(d.*pi.^2.*(-4.100206611570248e+11)).*exp(d.*pi.^2.*(-2.09194214876033e+11)).*1.324279489053236e-5+exp(d.*pi.^2.*(-4.100206611570248e+11)).*exp(d.*pi.^2.*(-1.571280991735537e+11)).*1.763097557809323e-5+exp(d.*pi.^2.*(-4.100206611570248e+11)).*exp(d.*pi.^2.*(-7.53099173553719e+10)).*3.678563341458988e-5+exp(d.*pi.^2.*(-4.100206611570248e+11)).*exp(d.*pi.^2.*(-4.555785123966942e+10)).*6.080868683344755e-5+exp(d.*pi.^2.*(-4.100206611570248e+11)).*exp(d.*pi.^2.*(-2.324380165289256e+10)).*1.191856013869514e-4+exp(d.*pi.^2.*(-4.100206611570248e+11)).*exp(d.*pi.^2.*(-8.367768595041322e+9)).*3.310686295535354e-4+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-3.486570247933884e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-2.95103305785124e+11)).*3.407353435202765e-6+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-2.460123966942149e+11)).*4.08728407829132e-6+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-2.013842975206611e+11)).*4.993046357300908e-6+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-1.612190082644628e+11)).*6.236998487211489e-6+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-1.255165289256198e+11)).*8.01106390617483e-6+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-9.427685950413223e+10)).*1.066563918356019e-5+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-4.518595041322314e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-2.733471074380165e+10)).*3.678545807740265e-5+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*7.209984578779232e-5+exp(d.*pi.^2.*(-4.066735537190082e+11)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-3.934710743801653e+11)).*exp(d.*pi.^2.*(-3.280165289256198e+11)).*5.632570553726442e-6+exp(d.*pi.^2.*(-3.934710743801653e+11)).*exp(d.*pi.^2.*(-2.685123966942149e+11)).*6.88077592524012e-6+exp(d.*pi.^2.*(-3.934710743801653e+11)).*exp(d.*pi.^2.*(-2.149586776859504e+11)).*8.595031162451022e-6+exp(d.*pi.^2.*(-3.934710743801653e+11)).*exp(d.*pi.^2.*(-1.673553719008264e+11)).*1.103982052571315e-5+exp(d.*pi.^2.*(-3.934710743801653e+11)).*exp(d.*pi.^2.*(-1.25702479338843e+11)).*1.469801561410095e-5+exp(d.*pi.^2.*(-3.934710743801653e+11)).*exp(d.*pi.^2.*(-6.024793388429752e+10)).*3.066624486588446e-5+exp(d.*pi.^2.*(-3.934710743801653e+11)).*exp(d.*pi.^2.*(-3.644628099173554e+10)).*5.069299906815745e-5+exp(d.*pi.^2.*(-3.934710743801653e+11)).*exp(d.*pi.^2.*(-1.859504132231405e+10)).*9.935875768202262e-5+exp(d.*pi.^2.*(-3.934710743801653e+11)).*exp(d.*pi.^2.*(-6.694214876033058e+9)).*2.759944771611512e-4+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-3.278925619834711e+11)).*3.974341356479841e-6+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-2.733471074380165e+11)).*4.767413318562332e-6+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-2.237603305785124e+11)).*5.823895586417536e-6+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-1.791322314049587e+11)).*7.274842924110025e-6+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-1.394628099173554e+11)).*9.344115072646962e-6+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-1.047520661157025e+11)).*1.244041503372621e-5+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-5.020661157024793e+10)).*2.595594015368134e-5+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-3.037190082644628e+10)).*4.290660482814774e-5+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-1.549586776859504e+10)).*8.409735131958502e-5+exp(d.*pi.^2.*(-3.87396694214876e+11)).*exp(d.*pi.^2.*(-5.578512396694215e+9)).*2.336019999602518e-4+exp(d.*pi.^2.*(-3.826859504132231e+11)).*exp(d.*pi.^2.*(-3.132644628099173e+11)).*8.253821160718708e-6+exp(d.*pi.^2.*(-3.826859504132231e+11)).*exp(d.*pi.^2.*(-2.507851239669421e+11)).*1.031015264215268e-5+exp(d.*pi.^2.*(-3.826859504132231e+11)).*exp(d.*pi.^2.*(-1.952479338842975e+11)).*1.324279489053236e-5+exp(d.*pi.^2.*(-3.826859504132231e+11)).*exp(d.*pi.^2.*(-1.466528925619835e+11)).*1.763097557809323e-5+exp(d.*pi.^2.*(-3.826859504132231e+11)).*exp(d.*pi.^2.*(-7.02892561983471e+10)).*3.678563341458988e-5+exp(d.*pi.^2.*(-3.826859504132231e+11)).*exp(d.*pi.^2.*(-4.252066115702479e+10)).*6.080868683344755e-5+exp(d.*pi.^2.*(-3.826859504132231e+11)).*exp(d.*pi.^2.*(-2.169421487603306e+10)).*1.191856013869514e-4+exp(d.*pi.^2.*(-3.826859504132231e+11)).*exp(d.*pi.^2.*(-7.8099173553719e+9)).*3.310686295535354e-4+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-3.16301652892562e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-2.711776859504132e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-2.295247933884297e+11)).*2.953582823893289e-6+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-1.913429752066116e+11)).*3.542964438408886e-6+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-1.566322314049587e+11)).*4.328102804794397e-6+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-1.253925619834711e+11)).*5.406392954178542e-6+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-9.762396694214876e+10)).*6.944199128254334e-6+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-7.332644628099173e+10)).*9.245254211948807e-6+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-3.514462809917355e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-2.12603305785124e+10)).*3.188659445303408e-5+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-1.084710743801653e+10)).*6.249802674535459e-5+exp(d.*pi.^2.*(-3.64896694214876e+11)).*exp(d.*pi.^2.*(-3.90495867768595e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-3.099173553719008e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-2.623140495867768e+11)).*3.407353435202765e-6+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-2.186776859504132e+11)).*4.08728407829132e-6+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-1.790082644628099e+11)).*4.993046357300908e-6+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-1.433057851239669e+11)).*6.236998487211489e-6+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-1.115702479338843e+11)).*8.01106390617483e-6+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-8.380165289256198e+10)).*1.066563918356019e-5+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-4.016528925619835e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-2.429752066115702e+10)).*3.678545807740265e-5+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-1.239669421487603e+10)).*7.209984578779232e-5+exp(d.*pi.^2.*(-3.614876033057851e+11)).*exp(d.*pi.^2.*(-4.462809917355372e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-3.580165289256198e+11)).*exp(d.*pi.^2.*(-2.866115702479339e+11)).*1.259493323856204e-5+exp(d.*pi.^2.*(-3.580165289256198e+11)).*exp(d.*pi.^2.*(-2.231404958677686e+11)).*1.617746345057028e-5+exp(d.*pi.^2.*(-3.580165289256198e+11)).*exp(d.*pi.^2.*(-1.67603305785124e+11)).*2.153808658747825e-5+exp(d.*pi.^2.*(-3.580165289256198e+11)).*exp(d.*pi.^2.*(-8.033057851239669e+10)).*4.493751092498283e-5+exp(d.*pi.^2.*(-3.580165289256198e+11)).*exp(d.*pi.^2.*(-4.859504132231405e+10)).*7.428419127963449e-5+exp(d.*pi.^2.*(-3.580165289256198e+11)).*exp(d.*pi.^2.*(-2.479338842975206e+10)).*1.455977175671664e-4+exp(d.*pi.^2.*(-3.580165289256198e+11)).*exp(d.*pi.^2.*(-8.925619834710743e+9)).*4.044350681638781e-4+exp(d.*pi.^2.*(-3.553512396694215e+11)).*exp(d.*pi.^2.*(-2.908884297520661e+11)).*8.253821160718708e-6+exp(d.*pi.^2.*(-3.553512396694215e+11)).*exp(d.*pi.^2.*(-2.328719008264463e+11)).*1.031015264215268e-5+exp(d.*pi.^2.*(-3.553512396694215e+11)).*exp(d.*pi.^2.*(-1.81301652892562e+11)).*1.324279489053236e-5+exp(d.*pi.^2.*(-3.553512396694215e+11)).*exp(d.*pi.^2.*(-1.361776859504132e+11)).*1.763097557809323e-5+exp(d.*pi.^2.*(-3.553512396694215e+11)).*exp(d.*pi.^2.*(-6.526859504132231e+10)).*3.678563341458988e-5+exp(d.*pi.^2.*(-3.553512396694215e+11)).*exp(d.*pi.^2.*(-3.948347107438016e+10)).*6.080868683344755e-5+exp(d.*pi.^2.*(-3.553512396694215e+11)).*exp(d.*pi.^2.*(-2.014462809917355e+10)).*1.191856013869514e-4+exp(d.*pi.^2.*(-3.553512396694215e+11)).*exp(d.*pi.^2.*(-7.252066115702479e+9)).*3.310686295535354e-4+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-2.95103305785124e+11)).*3.974341356479841e-6+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-2.460123966942149e+11)).*4.767413318562332e-6+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-2.013842975206611e+11)).*5.823895586417536e-6+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-1.612190082644628e+11)).*7.274842924110025e-6+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-1.255165289256198e+11)).*9.344115072646962e-6+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-9.427685950413223e+10)).*1.244041503372621e-5+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-4.518595041322314e+10)).*2.595594015368134e-5+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-2.733471074380165e+10)).*4.290660482814774e-5+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*8.409735131958502e-5+exp(d.*pi.^2.*(-3.486570247933884e+11)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*2.336019999602518e-4+exp(d.*pi.^2.*(-3.356404958677686e+11)).*exp(d.*pi.^2.*(-2.68698347107438e+11)).*1.259493323856204e-5+exp(d.*pi.^2.*(-3.356404958677686e+11)).*exp(d.*pi.^2.*(-2.09194214876033e+11)).*1.617746345057028e-5+exp(d.*pi.^2.*(-3.356404958677686e+11)).*exp(d.*pi.^2.*(-1.571280991735537e+11)).*2.153808658747825e-5+exp(d.*pi.^2.*(-3.356404958677686e+11)).*exp(d.*pi.^2.*(-7.53099173553719e+10)).*4.493751092498283e-5+exp(d.*pi.^2.*(-3.356404958677686e+11)).*exp(d.*pi.^2.*(-4.555785123966942e+10)).*7.428419127963449e-5+exp(d.*pi.^2.*(-3.356404958677686e+11)).*exp(d.*pi.^2.*(-2.324380165289256e+10)).*1.455977175671664e-4+exp(d.*pi.^2.*(-3.356404958677686e+11)).*exp(d.*pi.^2.*(-8.367768595041322e+9)).*4.044350681638781e-4+exp(d.*pi.^2.*(-3.280165289256198e+11)).*exp(d.*pi.^2.*(-2.685123966942149e+11)).*8.253821160718708e-6+exp(d.*pi.^2.*(-3.280165289256198e+11)).*exp(d.*pi.^2.*(-2.149586776859504e+11)).*1.031015264215268e-5+exp(d.*pi.^2.*(-3.280165289256198e+11)).*exp(d.*pi.^2.*(-1.673553719008264e+11)).*1.324279489053236e-5+exp(d.*pi.^2.*(-3.280165289256198e+11)).*exp(d.*pi.^2.*(-1.25702479338843e+11)).*1.763097557809323e-5+exp(d.*pi.^2.*(-3.280165289256198e+11)).*exp(d.*pi.^2.*(-6.024793388429752e+10)).*3.678563341458988e-5+exp(d.*pi.^2.*(-3.280165289256198e+11)).*exp(d.*pi.^2.*(-3.644628099173554e+10)).*6.080868683344755e-5+exp(d.*pi.^2.*(-3.280165289256198e+11)).*exp(d.*pi.^2.*(-1.859504132231405e+10)).*1.191856013869514e-4+exp(d.*pi.^2.*(-3.280165289256198e+11)).*exp(d.*pi.^2.*(-6.694214876033058e+9)).*3.310686295535354e-4+exp(d.*pi.^2.*(-3.278925619834711e+11)).*exp(d.*pi.^2.*(-2.733471074380165e+11)).*5.632570553726442e-6+exp(d.*pi.^2.*(-3.278925619834711e+11)).*exp(d.*pi.^2.*(-2.237603305785124e+11)).*6.88077592524012e-6+exp(d.*pi.^2.*(-3.278925619834711e+11)).*exp(d.*pi.^2.*(-1.791322314049587e+11)).*8.595031162451022e-6+exp(d.*pi.^2.*(-3.278925619834711e+11)).*exp(d.*pi.^2.*(-1.394628099173554e+11)).*1.103982052571315e-5+exp(d.*pi.^2.*(-3.278925619834711e+11)).*exp(d.*pi.^2.*(-1.047520661157025e+11)).*1.469801561410095e-5+exp(d.*pi.^2.*(-3.278925619834711e+11)).*exp(d.*pi.^2.*(-5.020661157024793e+10)).*3.066624486588446e-5+exp(d.*pi.^2.*(-3.278925619834711e+11)).*exp(d.*pi.^2.*(-3.037190082644628e+10)).*5.069299906815745e-5+exp(d.*pi.^2.*(-3.278925619834711e+11)).*exp(d.*pi.^2.*(-1.549586776859504e+10)).*9.935875768202262e-5+exp(d.*pi.^2.*(-3.278925619834711e+11)).*exp(d.*pi.^2.*(-5.578512396694215e+9)).*2.759944771611512e-4+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-2.711776859504132e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-2.295247933884297e+11)).*3.407353435202765e-6+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-1.913429752066116e+11)).*4.08728407829132e-6+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-1.566322314049587e+11)).*4.993046357300908e-6+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-1.253925619834711e+11)).*6.236998487211489e-6+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-9.762396694214876e+10)).*8.01106390617483e-6+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-7.332644628099173e+10)).*1.066563918356019e-5+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-3.514462809917355e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-2.12603305785124e+10)).*3.678545807740265e-5+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-1.084710743801653e+10)).*7.209984578779232e-5+exp(d.*pi.^2.*(-3.16301652892562e+11)).*exp(d.*pi.^2.*(-3.90495867768595e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-3.132644628099173e+11)).*exp(d.*pi.^2.*(-2.507851239669421e+11)).*1.259493323856204e-5+exp(d.*pi.^2.*(-3.132644628099173e+11)).*exp(d.*pi.^2.*(-1.952479338842975e+11)).*1.617746345057028e-5+exp(d.*pi.^2.*(-3.132644628099173e+11)).*exp(d.*pi.^2.*(-1.466528925619835e+11)).*2.153808658747825e-5+exp(d.*pi.^2.*(-3.132644628099173e+11)).*exp(d.*pi.^2.*(-7.02892561983471e+10)).*4.493751092498283e-5+exp(d.*pi.^2.*(-3.132644628099173e+11)).*exp(d.*pi.^2.*(-4.252066115702479e+10)).*7.428419127963449e-5+exp(d.*pi.^2.*(-3.132644628099173e+11)).*exp(d.*pi.^2.*(-2.169421487603306e+10)).*...
[opt_mc, err_mc] = fmincon(f, d0)
Local minimum possible. Constraints satisfied.
fmincon stopped because the size of the current step is less than
the value of the step size tolerance and constraints are
satisfied to within the value of the constraint tolerance.
opt_mc =
1.82376530468465e-11
err_mc =
2.75102571512948
[opt_ms, err_ms] = fminsearch(f, d0)
opt_ms =
-3.31787109375621e-13
err_ms =
1.68163295678904
y_pred_n = double(subs(y_pred, d, opt_mc));
plot(T1, y_obs, 'k*', T1, y_pred_n, 'b')
legend({'observed', 'modeled'})
title('fmincon predictions')
y_pred_n = double(subs(y_pred, d, opt_ms));
plot(T1, y_obs, 'k*', T1, y_pred_n, 'b')
legend({'observed', 'modeled'})
title('fminsearch predictions')
These are not straight lines, and the predicted minimum is not the same as d0.
The visual fit is not good in either case.
Walter Roberson
2021-6-17
If you continue on with
trypoints = -logspace(-14,-5, 10000);
predictions = f(trypoints);
[bestprediction, bestidx] = min(predictions)
bestd = trypoints(bestidx)
trypointsp = logspace(-14,-5, 10000);
predictions = f(trypointsp);
[bestpredictionp, bestidx] = min(predictions)
bestdp = trypointsp(bestidx)
then the bestd is still around -3.34774581070573e-13
Either your calculation is quite sensitive to numeric values (that is, if you were to proceed step by step numerically rather than creating a formula and substituting into that, then you might get a different result due to round-off happening much more often)... or else the model is not a good fit for the data.
Anand Ra
2021-6-18
Thanks you for the guidance.
Since the model isn’t yielding the correct initial condition with the fitting algorithm. Will adding constraints help?
Also, wanted to let you know that I expanded to y_obs to 72 entries. Still produces the same result as yours.
Walter Roberson
2021-6-18
Constraints will not help much. All a constraint could help you do would be to avoid going below about -4e-13, which is the point at which it becomes NaN.
My tests show that at least for the data that you posted earlier, that the shape goes:
- almost everything negative gives NaN
- between roughly -4e-13 and -3.35e-13, the error calculation exists and decreases like a quadratic
- above roughly -3.35e-13 and 0, the value increases again like a quadratic
- above 0 the value keeps increasing
When I took the derivative of y_pred and solved for it being 0, the solution was the approximately -3.35 value.
If your assignment statements are correct, then the model is simply wrong for the data. Which implies you should check your assignment statements in detail.
Anand Ra
2021-6-18
Thank you Mr.Walter! very much appreciated. I had certain level of confidence with the assignments( y_pred data) because I cross verified with my excel model to ensure consistency. Also, I can try again with different set of y_obs to see if the trend changes.
As far as the solver is concerned, the literature I am referring to suggests that this can be performed using lsqnonlin. But, I have spent a lot of time on it, and I find it very difficlt to condense my complex equation to accomodate to the form that is compatible to execute using lsqnonlin solver. Any thoughts on what solver I can use for my objective? Thanks for all your guidance!!
Walter Roberson
2021-6-18
format long g
% observed data
y_obs = [0.3 0.2 0.28 0.318 0.421 0.492 0.572 0.55 0.63 0.61 0.73 0.8 0.81 0.84 0.93 0.91]'; % If y_obs should equal to predicted, I can have more data. J us fo rthe code, I am providing limited observed data
t1 = [300:300:21600]';
T1 = t1(1:length(y_obs)).';
a=0.0011;
gama = 0.01005;
d0 = 0.000000000302;
syms d
n=1;
t=300;
L2 = sym(zeros(14,1));
L3 = sym(zeros(14,1));
L4 = sym(zeros(14,1));
At = sym(zeros(14,1));
t = 300;
n =1;
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
y_pred = sym(zeros(length(T1),1));
for t = T1
for n=1:1:14
L2(n) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n) = ((L2(n).*L3(n))/L4(n));
end
S(t/300) = sum(L5);
y_pred(t/300) = 1 -L1*S(t/300); % predicted data
end
flsq = matlabFunction(y_pred(:))
flsq = function_handle with value:
@(d)[exp(d.*pi.^2.*(-7.5e+9)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-5.212809917355372e+10)).*9.638060662976441e-4-exp(d.*pi.^2.*(-4.518595041322314e+10)).*1.111879404329578e-3-exp(d.*pi.^2.*(-3.87396694214876e+10)).*1.296897543527722e-3-exp(d.*pi.^2.*(-3.278925619834711e+10)).*1.532249550596395e-3-exp(d.*pi.^2.*(-2.733471074380165e+10)).*1.838006919804312e-3-exp(d.*pi.^2.*(-2.237603305785124e+10)).*2.245318304241606e-3-exp(d.*pi.^2.*(-1.791322314049587e+10)).*2.804709963564875e-3-exp(d.*pi.^2.*(-1.394628099173554e+10)).*3.602487767549398e-3-exp(d.*pi.^2.*(-1.047520661157025e+10)).*4.796221218789088e-3-exp(d.*pi.^2.*(-5.020661157024793e+9)).*1.000693550667001e-2-exp(d.*pi.^2.*(-3.037190082644628e+9)).*1.654201792665767e-2-exp(d.*pi.^2.*(-1.549586776859504e+9)).*3.242251160829167e-2-exp(d.*pi.^2.*(-5.578512396694215e+8)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-1.5e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-1.042561983471074e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-9.037190082644628e+10)).*1.111879404329578e-3-exp(d.*pi.^2.*(-7.74793388429752e+10)).*1.296897543527722e-3-exp(d.*pi.^2.*(-6.557851239669421e+10)).*1.532249550596395e-3-exp(d.*pi.^2.*(-5.46694214876033e+10)).*1.838006919804312e-3-exp(d.*pi.^2.*(-4.475206611570248e+10)).*2.245318304241606e-3-exp(d.*pi.^2.*(-3.582644628099173e+10)).*2.804709963564875e-3-exp(d.*pi.^2.*(-2.789256198347107e+10)).*3.602487767549398e-3-exp(d.*pi.^2.*(-2.095041322314049e+10)).*4.796221218789088e-3-exp(d.*pi.^2.*(-1.004132231404959e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-6.074380165289256e+9)).*1.654201792665767e-2-exp(d.*pi.^2.*(-3.099173553719008e+9)).*3.242251160829167e-2-exp(d.*pi.^2.*(-1.115702479338843e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-2.25e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-1.563842975206611e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-1.355578512396694e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-1.162190082644628e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-9.836776859504132e+10)).*1.532249550596395e-3-exp(d.*pi.^2.*(-8.200413223140495e+10)).*1.838006919804312e-3-exp(d.*pi.^2.*(-6.712809917355372e+10)).*2.245318304241606e-3-exp(d.*pi.^2.*(-5.37396694214876e+10)).*2.804709963564875e-3-exp(d.*pi.^2.*(-4.183884297520661e+10)).*3.602487767549398e-3-exp(d.*pi.^2.*(-3.142561983471074e+10)).*4.796221218789088e-3-exp(d.*pi.^2.*(-1.506198347107438e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-9.111570247933884e+9)).*1.654201792665767e-2-exp(d.*pi.^2.*(-4.648760330578512e+9)).*3.242251160829167e-2-exp(d.*pi.^2.*(-1.673553719008264e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-3.0e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-2.085123966942149e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-1.807438016528926e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-1.549586776859504e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-1.311570247933884e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-1.093388429752066e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-8.950413223140495e+10)).*2.245318304241606e-3-exp(d.*pi.^2.*(-7.165289256198347e+10)).*2.804709963564875e-3-exp(d.*pi.^2.*(-5.578512396694215e+10)).*3.602487767549398e-3-exp(d.*pi.^2.*(-4.190082644628099e+10)).*4.796221218789088e-3-exp(d.*pi.^2.*(-2.008264462809917e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-1.214876033057851e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-6.198347107438016e+9)).*3.242251160829167e-2-exp(d.*pi.^2.*(-2.231404958677686e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-3.75e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-2.606404958677686e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-2.259297520661157e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-1.93698347107438e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-1.639462809917355e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-1.366735537190083e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-1.118801652892562e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-8.956611570247933e+10)).*2.804709963564875e-3-exp(d.*pi.^2.*(-6.973140495867768e+10)).*3.602487767549398e-3-exp(d.*pi.^2.*(-5.237603305785124e+10)).*4.796221218789088e-3-exp(d.*pi.^2.*(-2.510330578512397e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-1.518595041322314e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-7.74793388429752e+9)).*3.242251160829167e-2-exp(d.*pi.^2.*(-2.789256198347107e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-4.5e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-3.127685950413223e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-2.711157024793388e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-2.324380165289256e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-1.967355371900826e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-1.640082644628099e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-1.342561983471074e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-1.074793388429752e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-8.367768595041322e+10)).*3.602487767549398e-3-exp(d.*pi.^2.*(-6.285123966942148e+10)).*4.796221218789088e-3-exp(d.*pi.^2.*(-3.012396694214876e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-1.822314049586777e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-9.297520661157024e+9)).*3.242251160829167e-2-exp(d.*pi.^2.*(-3.347107438016529e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-5.25e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-3.64896694214876e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-3.16301652892562e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-2.711776859504132e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-2.295247933884297e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-1.913429752066116e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-1.566322314049587e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-1.253925619834711e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-9.762396694214876e+10)).*3.602487767549398e-3-exp(d.*pi.^2.*(-7.332644628099173e+10)).*4.796221218789088e-3-exp(d.*pi.^2.*(-3.514462809917355e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-2.12603305785124e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-1.084710743801653e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-3.90495867768595e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-6.0e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-4.170247933884297e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-3.614876033057851e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-3.099173553719008e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-2.623140495867768e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-2.186776859504132e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-1.790082644628099e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-1.433057851239669e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-1.115702479338843e+11)).*3.602487767549398e-3-exp(d.*pi.^2.*(-8.380165289256198e+10)).*4.796221218789088e-3-exp(d.*pi.^2.*(-4.016528925619835e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-2.429752066115702e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-1.239669421487603e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-4.462809917355372e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-6.75e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-4.691528925619834e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-4.066735537190082e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-3.486570247933884e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-2.95103305785124e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-2.460123966942149e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-2.013842975206611e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-1.612190082644628e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-1.255165289256198e+11)).*3.602487767549398e-3-exp(d.*pi.^2.*(-9.427685950413223e+10)).*4.796221218789088e-3-exp(d.*pi.^2.*(-4.518595041322314e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-2.733471074380165e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-1.394628099173554e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-5.020661157024793e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-7.5e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-5.212809917355372e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-4.518595041322314e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-3.87396694214876e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-3.278925619834711e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-2.733471074380165e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-2.237603305785124e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-1.791322314049587e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-1.394628099173554e+11)).*3.602487767549398e-3-exp(d.*pi.^2.*(-1.047520661157025e+11)).*4.796221218789088e-3-exp(d.*pi.^2.*(-5.020661157024793e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-3.037190082644628e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-1.549586776859504e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-5.578512396694215e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-5.734090909090909e+11)).*(-9.638060662976441e-4)-exp(d.*pi.^2.*(-4.970454545454545e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-4.261363636363636e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-3.606818181818182e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-3.006818181818182e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-2.461363636363636e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-1.970454545454545e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-1.534090909090909e+11)).*3.602487767549398e-3-exp(d.*pi.^2.*(-1.152272727272727e+11)).*4.796221218789088e-3-exp(d.*pi.^2.*(-8.25e+10)).*6.698838604180037e-3-exp(d.*pi.^2.*(-5.522727272727272e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-3.340909090909091e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-1.704545454545454e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-6.136363636363636e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-9.0e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-6.255371900826446e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-5.422314049586777e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-4.648760330578512e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-3.934710743801653e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-3.280165289256198e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-2.685123966942149e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-2.149586776859504e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-1.673553719008264e+11)).*3.602487767549398e-3-exp(d.*pi.^2.*(-1.25702479338843e+11)).*4.796221218789088e-3-exp(d.*pi.^2.*(-6.024793388429752e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-3.644628099173554e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-1.859504132231405e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-6.694214876033058e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-9.75e+10)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-6.776652892561983e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-5.874173553719008e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-5.036157024793388e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-4.262603305785124e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-3.553512396694215e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-2.908884297520661e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-2.328719008264463e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-1.81301652892562e+11)).*3.602487767549398e-3-exp(d.*pi.^2.*(-1.361776859504132e+11)).*4.796221218789088e-3-exp(d.*pi.^2.*(-6.526859504132231e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-3.948347107438016e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-2.014462809917355e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-7.252066115702479e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-1.05e+11)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-7.29793388429752e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-6.326033057851239e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-5.423553719008264e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-4.590495867768595e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-3.826859504132231e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-3.132644628099173e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-2.507851239669421e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-1.952479338842975e+11)).*3.602487767549398e-3-exp(d.*pi.^2.*(-1.466528925619835e+11)).*4.796221218789088e-3-exp(d.*pi.^2.*(-7.02892561983471e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-4.252066115702479e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-2.169421487603306e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-7.8099173553719e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-1.125e+11)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-7.819214876033057e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-6.777892561983471e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-5.81095041322314e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-4.918388429752066e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-4.100206611570248e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-3.356404958677686e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-2.68698347107438e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-2.09194214876033e+11)).*3.602487767549398e-3-exp(d.*pi.^2.*(-1.571280991735537e+11)).*4.796221218789088e-3-exp(d.*pi.^2.*(-7.53099173553719e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-4.555785123966942e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-2.324380165289256e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-8.367768595041322e+9)).*9.006185612967756e-2+1.0;exp(d.*pi.^2.*(-1.2e+11)).*(-6.698838604180037e-3)-exp(d.*pi.^2.*(-8.340495867768595e+11)).*9.638060662976441e-4-exp(d.*pi.^2.*(-7.229752066115702e+11)).*1.111879404329578e-3-exp(d.*pi.^2.*(-6.198347107438016e+11)).*1.296897543527722e-3-exp(d.*pi.^2.*(-5.246280991735537e+11)).*1.532249550596395e-3-exp(d.*pi.^2.*(-4.373553719008264e+11)).*1.838006919804312e-3-exp(d.*pi.^2.*(-3.580165289256198e+11)).*2.245318304241606e-3-exp(d.*pi.^2.*(-2.866115702479339e+11)).*2.804709963564875e-3-exp(d.*pi.^2.*(-2.231404958677686e+11)).*3.602487767549398e-3-exp(d.*pi.^2.*(-1.67603305785124e+11)).*4.796221218789088e-3-exp(d.*pi.^2.*(-8.033057851239669e+10)).*1.000693550667001e-2-exp(d.*pi.^2.*(-4.859504132231405e+10)).*1.654201792665767e-2-exp(d.*pi.^2.*(-2.479338842975206e+10)).*3.242251160829167e-2-exp(d.*pi.^2.*(-8.925619834710743e+9)).*9.006185612967756e-2+1.0]
[opt_lsq, err_lsq] = lsqnonlin(flsq, d0)
Local minimum possible.
lsqnonlin stopped because the size of the current step is less than
the value of the step size tolerance.
opt_lsq =
3.02e-10
err_lsq =
15.9575171934164
y_pred_n = double(subs(y_pred, d, opt_lsq));
plot(T1, y_obs, 'k*', T1, y_pred_n, 'b')
legend({'observed', 'modeled'})
title('lsqnonlin predictions')
Walter Roberson
2021-6-18
I recommend that you take an approach like I showed in https://www.mathworks.com/matlabcentral/answers/857320-how-to-perform-data-fit-like-excel-and-plot#comment_1589825 where you construct a vector of trial d values to evaluate f at, and plot the result. You will be able to see visually that there isn't any hidden low-residue area, that you are not dealing with multiple minima and maxima. Not unless mabye you start going for much higher d values... I can't say that I plotted for (say) d = 317.
Anand Ra
2021-6-19
I will re-examine my assignment and observed data as you suggested. Thanks a lot for your help and gudiance!! Very much appreciated.
Walter Roberson
2021-6-19
Your model cannot even get close if you confine to two observations :(
If you confine to one observation, then the fminsearch() version can get quite close for that one point, but the fmincon() version doesn't get close.
Taylor series doesn't help it seems.
format long g
% observed data
y_obs = [0.3 0.2 0.28 0.318 0.421 0.492 0.572 0.55 0.63 0.61 0.73 0.8 0.81 0.84 0.93 0.91]'; % If y_obs should equal to predicted, I can have more data. J us fo rthe code, I am providing limited observed data
t1 = [300:300:21600]';
T1 = t1(1:length(y_obs)).';
y_obsX = y_obs([1 3]);
T1X = T1([1 3]);
a=0.0011;
gama = 0.01005;
d0 = 0.000000000302;
syms d
n=1;
t=300;
L2 = sym(zeros(14,1));
L3 = sym(zeros(14,1));
L4 = sym(zeros(14,1));
At = sym(zeros(14,1));
t = 300;
n =1;
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
y_pred = sym(zeros(length(T1X),1));
for tidx = 1 : length(T1X)
t = T1X(tidx);
for n=1:1:14
L2(n) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n) = ((L2(n).*L3(n))/L4(n));
end
S(tidx) = sum(L5);
y_pred(tidx) = 1 - L1*S(tidx); % predicted data
end
sse = expand(sum((y_pred - y_obsX(:)).^2));
f = matlabFunction(sse)
f = function_handle with value:
@(d)exp(d.*pi.^2.*(-4.5e+10)).*4.487443864485274e-5-exp(d.*pi.^2.*(-2.25e+10)).*9.646327590019253e-3+exp(d.*pi.^2.*(-1.5e+10)).*4.487443864485274e-5-exp(d.*pi.^2.*(-7.5e+9)).*9.378374045852051e-3+exp(d.*pi.^2.*(-3.127685950413223e+11)).*9.289221334321387e-7+exp(d.*pi.^2.*(-2.711157024793388e+11)).*1.236275809772296e-6+exp(d.*pi.^2.*(-2.324380165289256e+11)).*1.681943238408239e-6+exp(d.*pi.^2.*(-1.967355371900826e+11)).*2.347788685302854e-6+exp(d.*pi.^2.*(-1.640082644628099e+11)).*3.378269437248533e-6-exp(d.*pi.^2.*(-1.563842975206611e+11)).*1.387880735468607e-3-exp(d.*pi.^2.*(-1.355578512396694e+11)).*1.601106342234592e-3+exp(d.*pi.^2.*(-1.342561983471074e+11)).*5.041454287362399e-6-exp(d.*pi.^2.*(-1.162190082644628e+11)).*1.867532462679919e-3+exp(d.*pi.^2.*(-1.074793388429752e+11)).*7.866397979720083e-6+exp(d.*pi.^2.*(-1.042561983471074e+11)).*9.289221334321387e-7-exp(d.*pi.^2.*(-9.836776859504132e+10)).*2.206439352858809e-3+exp(d.*pi.^2.*(-9.037190082644628e+10)).*1.236275809772296e-6+exp(d.*pi.^2.*(-8.367768595041322e+10)).*1.297791811534305e-5-exp(d.*pi.^2.*(-8.200413223140495e+10)).*2.646729964518209e-3+exp(d.*pi.^2.*(-7.74793388429752e+10)).*1.681943238408239e-6-exp(d.*pi.^2.*(-6.712809917355372e+10)).*3.233258358107912e-3+exp(d.*pi.^2.*(-6.557851239669421e+10)).*2.347788685302854e-6+exp(d.*pi.^2.*(-6.285123966942148e+10)).*2.300373797956269e-5+exp(d.*pi.^2.*(-5.46694214876033e+10)).*3.378269437248533e-6-exp(d.*pi.^2.*(-5.37396694214876e+10)).*4.03878234753342e-3-exp(d.*pi.^2.*(-5.212809917355372e+10)).*1.349328492816702e-3-exp(d.*pi.^2.*(-4.518595041322314e+10)).*1.556631166061409e-3+exp(d.*pi.^2.*(-4.475206611570248e+10)).*5.041454287362399e-6-exp(d.*pi.^2.*(-4.183884297520661e+10)).*5.187582385271133e-3-exp(d.*pi.^2.*(-3.87396694214876e+10)).*1.81565656093881e-3+exp(d.*pi.^2.*(-3.582644628099173e+10)).*7.866397979720083e-6-exp(d.*pi.^2.*(-3.278925619834711e+10)).*2.145149370834953e-3-exp(d.*pi.^2.*(-3.142561983471074e+10)).*6.906558555056287e-3+exp(d.*pi.^2.*(-3.012396694214876e+10)).*1.001387582346529e-4+exp(d.*pi.^2.*(-2.789256198347107e+10)).*1.297791811534305e-5-exp(d.*pi.^2.*(-2.733471074380165e+10)).*2.573209687726036e-3-exp(d.*pi.^2.*(-2.237603305785124e+10)).*3.143445625938248e-3+exp(d.*pi.^2.*(-2.095041322314049e+10)).*2.300373797956269e-5+exp(d.*pi.^2.*(-1.822314049586777e+10)).*2.736383570858638e-4-exp(d.*pi.^2.*(-1.791322314049587e+10)).*3.926593948990825e-3-exp(d.*pi.^2.*(-1.506198347107438e+10)).*1.440998712960481e-2-exp(d.*pi.^2.*(-1.394628099173554e+10)).*5.043482874569158e-3-exp(d.*pi.^2.*(-1.047520661157025e+10)).*6.714709706304724e-3+exp(d.*pi.^2.*(-1.004132231404959e+10)).*1.001387582346529e-4+exp(d.*pi.^2.*(-9.297520661157024e+9)).*1.051219258989808e-3-exp(d.*pi.^2.*(-9.111570247933884e+9)).*2.382050581438705e-2+exp(d.*pi.^2.*(-6.074380165289256e+9)).*2.736383570858638e-4-exp(d.*pi.^2.*(-5.020661157024793e+9)).*1.400970970933801e-2-exp(d.*pi.^2.*(-4.648760330578512e+9)).*4.668841671594001e-2+exp(d.*pi.^2.*(-3.347107438016529e+9)).*8.111137929522739e-3+exp(d.*pi.^2.*(-3.099173553719008e+9)).*1.051219258989808e-3-exp(d.*pi.^2.*(-3.037190082644628e+9)).*2.315882509732074e-2-exp(d.*pi.^2.*(-1.673553719008264e+9)).*1.296890728267357e-1-exp(d.*pi.^2.*(-1.549586776859504e+9)).*4.539151625160834e-2+exp(d.*pi.^2.*(-1.115702479338843e+9)).*8.111137929522739e-3-exp(d.*pi.^2.*(-5.578512396694215e+8)).*1.260865985815486e-1+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-1.355578512396694e+11)).*2.143272229768516e-6+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-1.162190082644628e+11)).*2.499915439637062e-6+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-9.836776859504132e+10)).*2.953582823893289e-6+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-8.200413223140495e+10)).*3.542964438408886e-6+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-6.712809917355372e+10)).*4.328102804794397e-6+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-5.37396694214876e+10)).*5.406392954178542e-6+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-4.183884297520661e+10)).*6.944199128254334e-6+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-3.142561983471074e+10)).*9.245254211948807e-6+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*1.928949029275568e-5+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*3.188659445303408e-5+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*6.249802674535459e-5+exp(d.*pi.^2.*(-1.563842975206611e+11)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*1.736043265596178e-4+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-1.162190082644628e+11)).*2.883987336348192e-6+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-9.836776859504132e+10)).*3.407353435202765e-6+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-8.200413223140495e+10)).*4.08728407829132e-6+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-6.712809917355372e+10)).*4.993046357300908e-6+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-5.37396694214876e+10)).*6.236998487211489e-6+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-4.183884297520661e+10)).*8.01106390617483e-6+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-3.142561983471074e+10)).*1.066563918356019e-5+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*2.22530109806415e-5+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*3.678545807740265e-5+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*7.209984578779232e-5+exp(d.*pi.^2.*(-1.355578512396694e+11)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*2.00275845892564e-4+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-9.836776859504132e+10)).*3.974341356479841e-6+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-8.200413223140495e+10)).*4.767413318562332e-6+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-6.712809917355372e+10)).*5.823895586417536e-6+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-5.37396694214876e+10)).*7.274842924110025e-6+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-4.183884297520661e+10)).*9.344115072646962e-6+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-3.142561983471074e+10)).*1.244041503372621e-5+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*2.595594015368134e-5+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*4.290660482814774e-5+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*8.409735131958502e-5+exp(d.*pi.^2.*(-1.162190082644628e+11)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*2.336019999602518e-4+exp(d.*pi.^2.*(-9.836776859504132e+10)).*exp(d.*pi.^2.*(-8.200413223140495e+10)).*5.632570553726442e-6+exp(d.*pi.^2.*(-9.836776859504132e+10)).*exp(d.*pi.^2.*(-6.712809917355372e+10)).*6.88077592524012e-6+exp(d.*pi.^2.*(-9.836776859504132e+10)).*exp(d.*pi.^2.*(-5.37396694214876e+10)).*8.595031162451022e-6+exp(d.*pi.^2.*(-9.836776859504132e+10)).*exp(d.*pi.^2.*(-4.183884297520661e+10)).*1.103982052571315e-5+exp(d.*pi.^2.*(-9.836776859504132e+10)).*exp(d.*pi.^2.*(-3.142561983471074e+10)).*1.469801561410095e-5+exp(d.*pi.^2.*(-9.836776859504132e+10)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*3.066624486588446e-5+exp(d.*pi.^2.*(-9.836776859504132e+10)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*5.069299906815745e-5+exp(d.*pi.^2.*(-9.836776859504132e+10)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*9.935875768202262e-5+exp(d.*pi.^2.*(-9.836776859504132e+10)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*2.759944771611512e-4+exp(d.*pi.^2.*(-8.200413223140495e+10)).*exp(d.*pi.^2.*(-6.712809917355372e+10)).*8.253821160718708e-6+exp(d.*pi.^2.*(-8.200413223140495e+10)).*exp(d.*pi.^2.*(-5.37396694214876e+10)).*1.031015264215268e-5+exp(d.*pi.^2.*(-8.200413223140495e+10)).*exp(d.*pi.^2.*(-4.183884297520661e+10)).*1.324279489053236e-5+exp(d.*pi.^2.*(-8.200413223140495e+10)).*exp(d.*pi.^2.*(-3.142561983471074e+10)).*1.763097557809323e-5+exp(d.*pi.^2.*(-8.200413223140495e+10)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*3.678563341458988e-5+exp(d.*pi.^2.*(-8.200413223140495e+10)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*6.080868683344755e-5+exp(d.*pi.^2.*(-8.200413223140495e+10)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*1.191856013869514e-4+exp(d.*pi.^2.*(-8.200413223140495e+10)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*3.310686295535354e-4+exp(d.*pi.^2.*(-6.712809917355372e+10)).*exp(d.*pi.^2.*(-5.37396694214876e+10)).*1.259493323856204e-5+exp(d.*pi.^2.*(-6.712809917355372e+10)).*exp(d.*pi.^2.*(-4.183884297520661e+10)).*1.617746345057028e-5+exp(d.*pi.^2.*(-6.712809917355372e+10)).*exp(d.*pi.^2.*(-3.142561983471074e+10)).*2.153808658747825e-5+exp(d.*pi.^2.*(-6.712809917355372e+10)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*4.493751092498283e-5+exp(d.*pi.^2.*(-6.712809917355372e+10)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*7.428419127963449e-5+exp(d.*pi.^2.*(-6.712809917355372e+10)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*1.455977175671664e-4+exp(d.*pi.^2.*(-6.712809917355372e+10)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*4.044350681638781e-4+exp(d.*pi.^2.*(-5.37396694214876e+10)).*exp(d.*pi.^2.*(-4.183884297520661e+10)).*2.020786667053276e-5+exp(d.*pi.^2.*(-5.37396694214876e+10)).*exp(d.*pi.^2.*(-3.142561983471074e+10)).*2.690401887959805e-5+exp(d.*pi.^2.*(-5.37396694214876e+10)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*5.613310344061699e-5+exp(d.*pi.^2.*(-5.37396694214876e+10)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*9.27911249927311e-5+exp(d.*pi.^2.*(-5.37396694214876e+10)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*1.818714827031469e-4+exp(d.*pi.^2.*(-5.37396694214876e+10)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*5.051947704481059e-4+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-4.518595041322314e+10)).*2.143272229768516e-6+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-3.87396694214876e+10)).*2.499915439637062e-6+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-3.278925619834711e+10)).*2.953582823893289e-6+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-2.733471074380165e+10)).*3.542964438408886e-6+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-2.237603305785124e+10)).*4.328102804794397e-6+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-1.791322314049587e+10)).*5.406392954178542e-6+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*6.944199128254334e-6+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-1.047520661157025e+10)).*9.245254211948807e-6+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*1.928949029275568e-5+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*3.188659445303408e-5+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*6.249802674535459e-5+exp(d.*pi.^2.*(-5.212809917355372e+10)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*1.736043265596178e-4+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-3.87396694214876e+10)).*2.883987336348192e-6+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-3.278925619834711e+10)).*3.407353435202765e-6+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-2.733471074380165e+10)).*4.08728407829132e-6+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-2.237603305785124e+10)).*4.993046357300908e-6+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-1.791322314049587e+10)).*6.236998487211489e-6+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*8.01106390617483e-6+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-1.047520661157025e+10)).*1.066563918356019e-5+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*2.22530109806415e-5+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*3.678545807740265e-5+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*7.209984578779232e-5+exp(d.*pi.^2.*(-4.518595041322314e+10)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*2.00275845892564e-4+exp(d.*pi.^2.*(-4.183884297520661e+10)).*exp(d.*pi.^2.*(-3.142561983471074e+10)).*3.455665654229711e-5+exp(d.*pi.^2.*(-4.183884297520661e+10)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*7.209972550686889e-5+exp(d.*pi.^2.*(-4.183884297520661e+10)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*1.191848344627342e-4+exp(d.*pi.^2.*(-4.183884297520661e+10)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*2.336034029241982e-4+exp(d.*pi.^2.*(-4.183884297520661e+10)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*6.488934700599144e-4+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-3.278925619834711e+10)).*3.974341356479841e-6+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-2.733471074380165e+10)).*4.767413318562332e-6+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-2.237603305785124e+10)).*5.823895586417536e-6+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-1.791322314049587e+10)).*7.274842924110025e-6+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*9.344115072646962e-6+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-1.047520661157025e+10)).*1.244041503372621e-5+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*2.595594015368134e-5+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*4.290660482814774e-5+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*8.409735131958502e-5+exp(d.*pi.^2.*(-3.87396694214876e+10)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*2.336019999602518e-4+exp(d.*pi.^2.*(-3.278925619834711e+10)).*exp(d.*pi.^2.*(-2.733471074380165e+10)).*5.632570553726442e-6+exp(d.*pi.^2.*(-3.278925619834711e+10)).*exp(d.*pi.^2.*(-2.237603305785124e+10)).*6.88077592524012e-6+exp(d.*pi.^2.*(-3.278925619834711e+10)).*exp(d.*pi.^2.*(-1.791322314049587e+10)).*8.595031162451022e-6+exp(d.*pi.^2.*(-3.278925619834711e+10)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*1.103982052571315e-5+exp(d.*pi.^2.*(-3.278925619834711e+10)).*exp(d.*pi.^2.*(-1.047520661157025e+10)).*1.469801561410095e-5+exp(d.*pi.^2.*(-3.278925619834711e+10)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*3.066624486588446e-5+exp(d.*pi.^2.*(-3.278925619834711e+10)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*5.069299906815745e-5+exp(d.*pi.^2.*(-3.278925619834711e+10)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*9.935875768202262e-5+exp(d.*pi.^2.*(-3.278925619834711e+10)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*2.759944771611512e-4+exp(d.*pi.^2.*(-3.142561983471074e+10)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*9.599095282428926e-5+exp(d.*pi.^2.*(-3.142561983471074e+10)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*1.5867835476285e-4+exp(d.*pi.^2.*(-3.142561983471074e+10)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*3.110110762842481e-4+exp(d.*pi.^2.*(-3.142561983471074e+10)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*8.639131707453793e-4+exp(d.*pi.^2.*(-2.733471074380165e+10)).*exp(d.*pi.^2.*(-2.237603305785124e+10)).*8.253821160718708e-6+exp(d.*pi.^2.*(-2.733471074380165e+10)).*exp(d.*pi.^2.*(-1.791322314049587e+10)).*1.031015264215268e-5+exp(d.*pi.^2.*(-2.733471074380165e+10)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*1.324279489053236e-5+exp(d.*pi.^2.*(-2.733471074380165e+10)).*exp(d.*pi.^2.*(-1.047520661157025e+10)).*1.763097557809323e-5+exp(d.*pi.^2.*(-2.733471074380165e+10)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*3.678563341458988e-5+exp(d.*pi.^2.*(-2.733471074380165e+10)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*6.080868683344755e-5+exp(d.*pi.^2.*(-2.733471074380165e+10)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*1.191856013869514e-4+exp(d.*pi.^2.*(-2.733471074380165e+10)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*3.310686295535354e-4+exp(d.*pi.^2.*(-2.237603305785124e+10)).*exp(d.*pi.^2.*(-1.791322314049587e+10)).*1.259493323856204e-5+exp(d.*pi.^2.*(-2.237603305785124e+10)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*1.617746345057028e-5+exp(d.*pi.^2.*(-2.237603305785124e+10)).*exp(d.*pi.^2.*(-1.047520661157025e+10)).*2.153808658747825e-5+exp(d.*pi.^2.*(-2.237603305785124e+10)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*4.493751092498283e-5+exp(d.*pi.^2.*(-2.237603305785124e+10)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*7.428419127963449e-5+exp(d.*pi.^2.*(-2.237603305785124e+10)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*1.455977175671664e-4+exp(d.*pi.^2.*(-2.237603305785124e+10)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*4.044350681638781e-4+exp(d.*pi.^2.*(-1.791322314049587e+10)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*2.020786667053276e-5+exp(d.*pi.^2.*(-1.791322314049587e+10)).*exp(d.*pi.^2.*(-1.047520661157025e+10)).*2.690401887959805e-5+exp(d.*pi.^2.*(-1.791322314049587e+10)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*5.613310344061699e-5+exp(d.*pi.^2.*(-1.791322314049587e+10)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*9.27911249927311e-5+exp(d.*pi.^2.*(-1.791322314049587e+10)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*1.818714827031469e-4+exp(d.*pi.^2.*(-1.791322314049587e+10)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*5.051947704481059e-4+exp(d.*pi.^2.*(-1.506198347107438e+10)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*3.310698130844849e-4+exp(d.*pi.^2.*(-1.506198347107438e+10)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*6.488999652568689e-4+exp(d.*pi.^2.*(-1.506198347107438e+10)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*1.802486371801353e-3+exp(d.*pi.^2.*(-1.394628099173554e+10)).*exp(d.*pi.^2.*(-1.047520661157025e+10)).*3.455665654229711e-5+exp(d.*pi.^2.*(-1.394628099173554e+10)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*7.209972550686889e-5+exp(d.*pi.^2.*(-1.394628099173554e+10)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*1.191848344627342e-4+exp(d.*pi.^2.*(-1.394628099173554e+10)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*2.336034029241982e-4+exp(d.*pi.^2.*(-1.394628099173554e+10)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*6.488934700599144e-4+exp(d.*pi.^2.*(-1.047520661157025e+10)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*9.599095282428926e-5+exp(d.*pi.^2.*(-1.047520661157025e+10)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*1.5867835476285e-4+exp(d.*pi.^2.*(-1.047520661157025e+10)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*3.110110762842481e-4+exp(d.*pi.^2.*(-1.047520661157025e+10)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*8.639131707453793e-4+exp(d.*pi.^2.*(-9.111570247933884e+9)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*1.072667536503255e-3+exp(d.*pi.^2.*(-9.111570247933884e+9)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*2.979609677210381e-3+exp(d.*pi.^2.*(-5.020661157024793e+9)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*3.310698130844849e-4+exp(d.*pi.^2.*(-5.020661157024793e+9)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*6.488999652568689e-4+exp(d.*pi.^2.*(-5.020661157024793e+9)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*1.802486371801353e-3+exp(d.*pi.^2.*(-4.648760330578512e+9)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*5.84006315165753e-3+exp(d.*pi.^2.*(-3.037190082644628e+9)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*1.072667536503255e-3+exp(d.*pi.^2.*(-3.037190082644628e+9)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*2.979609677210381e-3+exp(d.*pi.^2.*(-1.549586776859504e+9)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*5.84006315165753e-3+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-1.563842975206611e+11)).*1.291276256771512e-5+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-1.355578512396694e+11)).*1.489660135383136e-5+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-1.162190082644628e+11)).*1.737541466049952e-5+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-9.836776859504132e+10)).*2.052858488154529e-5+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-8.200413223140495e+10)).*2.462502341827033e-5+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-6.712809917355372e+10)).*3.008204987025145e-5+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-5.37396694214876e+10)).*3.757659875491354e-5+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-4.183884297520661e+10)).*4.826496825669254e-5+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-3.142561983471074e+10)).*6.425822370922354e-5+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-1.506198347107438e+10)).*1.340696917632419e-4+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-9.111570247933884e+9)).*2.216246165562652e-4+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-4.648760330578512e+9)).*4.343863448121992e-4+exp(d.*pi.^2.*(-2.25e+10)).*exp(d.*pi.^2.*(-1.673553719008264e+9)).*1.206619677211185e-3+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-5.212809917355372e+10)).*1.291276256771512e-5+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-4.518595041322314e+10)).*1.489660135383136e-5+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-3.87396694214876e+10)).*1.737541466049952e-5+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-3.278925619834711e+10)).*2.052858488154529e-5+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-2.733471074380165e+10)).*2.462502341827033e-5+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-2.237603305785124e+10)).*3.008204987025145e-5+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-1.791322314049587e+10)).*3.757659875491354e-5+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-1.394628099173554e+10)).*4.826496825669254e-5+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-1.047520661157025e+10)).*6.425822370922354e-5+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-5.020661157024793e+9)).*1.340696917632419e-4+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-3.037190082644628e+9)).*2.216246165562652e-4+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-1.549586776859504e+9)).*4.343863448121992e-4+exp(d.*pi.^2.*(-7.5e+9)).*exp(d.*pi.^2.*(-5.578512396694215e+8)).*1.206619677211185e-3+1.0084
[opt_mc, err_mc] = fmincon(f, d0)
Local minimum possible. Constraints satisfied.
fmincon stopped because the size of the current step is less than
the value of the step size tolerance and constraints are
satisfied to within the value of the constraint tolerance.
opt_mc =
1.82376530468464e-11
err_mc =
0.737074171900867
[opt_ms, err_ms] = fminsearch(f, d0)
opt_ms =
-3.48007812500063e-12
err_ms =
0.237225768747742
y_pred_n = double(subs(y_pred, d, opt_mc));
plot(T1, y_obs, 'k*', T1X, y_pred_n, 'r-v')
legend({'observed', 'modeled'})
title('fmincon predictions')
y_pred_n = double(subs(y_pred, d, opt_ms));
plot(T1, y_obs, 'k*', T1X, y_pred_n, 'r-v')
legend({'observed', 'modeled'})
title('fminsearch predictions')
best_d = vpasolve(diff(sse, d))
best_d =
1.0
f(double(best_d))
ans =
1.0084
FT = taylor(sse, d, 0, 'order', 10)
FT =
trial_sol = vpasolve(FT)
trial_sol =
f(double(trial_sol))
ans =
320.450229336017 + 0i
0.343266356293643 + 0i
0.654561380806186 + 0i
0.562484132868365 + 0.0790092483031968i
0.562484132868365 - 0.0790092483031968i
0.630971889213914 + 0.068321190380106i
0.630971889213914 - 0.068321190380106i
0.648768811395252 + 0.0351822086775853i
0.648768811395252 - 0.0351822086775853i
Anand Ra
2021-6-19
编辑：Anand Ra
2021-6-19
Oh ok.
Btw, wanted to share the full set of data for y_obs to you. The previous "y_obs" data that I gave could probably be the issue. I am sorry about that. Below is the complete set of observed data
format long g
% observed data
y_obs = [ 0
1.6216e-01
2.9813e-01
4.0805e-01
4.9338e-01
5.5928e-01
6.0894e-01
6.5506e-01
6.8876e-01
7.1773e-01
7.4284e-01
7.6433e-01
7.8448e-01
7.9912e-01
8.1484e-01
8.2950e-01
8.3760e-01
8.4707e-01
8.5767e-01
8.6299e-01
8.6862e-01
8.7433e-01
8.7879e-01
8.8736e-01
8.9152e-01
8.9903e-01
9.0343e-01
9.0984e-01
9.1344e-01
9.1615e-01
9.2046e-01
9.2313e-01
9.2640e-01
9.2992e-01
9.3134e-01
9.3242e-01
9.3616e-01
9.3986e-01
9.4201e-01
9.4145e-01
9.4434e-01
9.4252e-01
9.4131e-01
9.4249e-01
9.4283e-01
9.4395e-01
9.4355e-01
9.4690e-01
9.4919e-01
9.5255e-01
9.5626e-01
9.6282e-01
9.6283e-01
9.6672e-01
9.6439e-01
9.6887e-01
9.7099e-01
9.7529e-01
9.8034e-01
9.8302e-01
9.8494e-01
9.8828e-01
9.8769e-01
9.9130e-01
9.9108e-01
9.9422e-01
9.9561e-01
9.9737e-01
1.0002e+00
1.0034e+00
1.0044e+00]
t1 = [300:300:21600]';
T1 = t1(1:length(y_obs)).';
a=0.0011;
gama = 0.01005;
d0 = 0.000000000302;
syms d
n=1;
t=300;
L2 = sym(zeros(14,1));
L3 = sym(zeros(14,1));
L4 = sym(zeros(14,1));
At = sym(zeros(14,1));
t = 300;
n =1;
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
y_pred = sym(zeros(length(T1),1));
for t = T1
for n=1:1:14
L2(n) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n) = ((L2(n).*L3(n))/L4(n));
end
S(t/300) = sum(L5);
y_pred(t/300) = 1 -L1*S(t/300); % predicted data
end
sse = expand(sum((y_pred - y_obs(:)).^2));
f = matlabFunction(sse)
[opt_mc, err_mc] = fmincon(f, d0)
[opt_ms, err_ms] = fminsearch(f, d0)
y_pred_n = double(subs(y_pred, d, opt_mc));
plot(T1, y_obs, 'k*', T1, y_pred_n, 'b')
legend({'observed', 'modeled'})
title('fmincon predictions')
y_pred_n = double(subs(y_pred, d, opt_ms));
plot(T1, y_obs, 'k*', T1, y_pred_n, 'b')
legend({'observed', 'modeled'})
title('fminsearch predictions')
Walter Roberson
2021-6-20
Above d=1e-9, your function pretty much predicts 1.0 for all values
There just isn't anywhere left to go. Your y_pred do not model your data.
Anand Ra
2021-6-20
编辑：Anand Ra
2021-6-21
I went and coded the model assignment again and verified the two sets of data. The problem was that when the assignments were changed to syms, I am unable to verify the array coming out of the equation..
Below is the code:
y_obs = [
0
1.6216e-01
2.9813e-01
4.0805e-01
4.9338e-01
5.5928e-01
6.0894e-01
6.5506e-01
6.8876e-01
7.1773e-01
7.4284e-01
7.6433e-01
7.8448e-01
7.9912e-01
8.1484e-01
8.2950e-01
8.3760e-01
8.4707e-01
8.5767e-01
8.6299e-01
8.6862e-01
8.7433e-01
8.7879e-01
8.8736e-01
8.9152e-01
8.9903e-01
9.0343e-01
9.0984e-01
9.1344e-01
9.1615e-01
9.2046e-01
9.2313e-01
9.2640e-01
9.2992e-01
9.3134e-01
9.3242e-01
9.3616e-01
9.3986e-01
9.4201e-01
9.4145e-01
9.4434e-01
9.4252e-01
9.4131e-01
9.4249e-01
9.4283e-01
9.4395e-01
9.4355e-01
9.4690e-01
9.4919e-01
9.5255e-01
9.5626e-01
9.6282e-01
9.6283e-01
9.6672e-01
9.6439e-01
9.6887e-01
9.7099e-01
9.7529e-01
9.8034e-01
9.8302e-01
9.8494e-01
9.8828e-01
9.8769e-01
9.9130e-01
9.9108e-01
9.9422e-01
9.9561e-01
9.9737e-01
1.0002e+00
1.0034e+00
1.0044e+00
1.00329
1.0];
t1 = [0:300:21600]';
a=0.0011;
gama = 0.01005;
d=0.000000000302;
n=1;
L2 = zeros(14,1);
L3 = zeros(14,1);
L4 = zeros(14,1);
L5 = zeros(14,1);
S= zeros(73,1);
y_pred = zeros(73,1);
At = zeros(73,1);
% t = 0;
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
format longE
j=1;
k =1;
for t= 0:300:21600
for n=0:1:14
L2(n+1) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n+1) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n+1)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n+1) = ((L2(n+1)*L3(n+1))/L4(n+1));
end
S((t/300) +1) = sum(L5);
y_pred((t/300)+1)= 1 -(L1*S((t/300) +1)); % predicted data
end
plot(t1, y_obs, 'k*', t1, y_pred, 'b')
legend({'observed', 'model data from equation'})
title('Plot before data fitting')
Anand Ra
2021-6-20
The y_obs and y_pred are similar. But, I am unable to perform fitting. If I use syms, the model predicted values are not verifiable. Is it possible to do what I am looking for, now that I have two clean data sets.
Walter Roberson
2021-6-21
You can perform fitting with numeric data.
t1 = [0:300:21600]';
y_obs = [
0
1.6216e-01
2.9813e-01
4.0805e-01
4.9338e-01
etc];
d0 = 0.000000000302;
fun = @(d) ypred(d, t1) - y_obs;
best_d = lsqnonlin(fun, d0);
predicted_y = ypred(best_d, t1);
plot(t1, y_obs, '*', t1, predicted_y, '-');
function y_pred = ypred(d, t1)
%HERE
end
at the point marked HERE insert your code that takes a numeric d and numeric t1 and calculates y_pred -- basically just moving some of your code around, such as moving a and gama inside of the new function, and otherwise mostly doing the same thing
Anand Ra
2021-6-21
Thanks, and its executing. However, its not minimizing the d or performing the fitting process. It throws back the same initial assignment. The reason is as you mentioned earlier? the data is very sensitive? I am wondering why its not working because, I looking to replicate the fitting in the literature.
t1 = [0:300:21600]';
y_obs = [
0
1.6216e-01
2.9813e-01
4.0805e-01
4.9338e-01
5.5928e-01
6.0894e-01
6.5506e-01
6.8876e-01
7.1773e-01
7.4284e-01
7.6433e-01
7.8448e-01
7.9912e-01
8.1484e-01
8.2950e-01
8.3760e-01
8.4707e-01
8.5767e-01
8.6299e-01
8.6862e-01
8.7433e-01
8.7879e-01
8.8736e-01
8.9152e-01
8.9903e-01
9.0343e-01
9.0984e-01
9.1344e-01
9.1615e-01
9.2046e-01
9.2313e-01
9.2640e-01
9.2992e-01
9.3134e-01
9.3242e-01
9.3616e-01
9.3986e-01
9.4201e-01
9.4145e-01
9.4434e-01
9.4252e-01
9.4131e-01
9.4249e-01
9.4283e-01
9.4395e-01
9.4355e-01
9.4690e-01
9.4919e-01
9.5255e-01
9.5626e-01
9.6282e-01
9.6283e-01
9.6672e-01
9.6439e-01
9.6887e-01
9.7099e-01
9.7529e-01
9.8034e-01
9.8302e-01
9.8494e-01
9.8828e-01
9.8769e-01
9.9130e-01
9.9108e-01
9.9422e-01
9.9561e-01
9.9737e-01
1.0002e+00
1.0034e+00
1.0044e+00
1.00329
1.0];
d0 = 0.000000000302;
fun = @(d) ypred(d, t1) - y_obs;
best_d = lsqnonlin(fun, d0)
predicted_y = ypred(best_d, t1);
plot(t1, y_obs, '*', t1, predicted_y, '-');
legend({'observed', 'predicted'})
title('lsqnonlin data fitting')
function y_pred = ypred(d, t1)
t1 = [0:300:21600]';
a=0.0011;
gama = 0.01005;
d=0.000000000302;
L2 = zeros(14,1);
L3 = zeros(100,1);
L4 = zeros(100,1);
L5 = zeros(100,1);
S= zeros(73,1);
y_pred = zeros(73,1);
% t = 0;
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
format longE
k =1;
for t= 0:300:21600
for n=0:1:100
L2(n+1) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n+1) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n+1)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n+1) = ((L2(n+1)*L3(n+1))/L4(n+1));
end
S((t/300) +1) = sum(L5);
y_pred((t/300)+1)= 1 -(L1*S((t/300) +1)); % predicted data
end
end
Anand Ra
2021-6-21
I plotted the equation predicted vs solver predicted. Just overlaps. No minimation is perfpormed by the solver
solver predicted
solver predicted
Walter Roberson
2021-6-21
编辑：Walter Roberson
2021-6-21
You overwrote d and t1 in the function and you hardcoded t1 in the function.
t1 = [0:300:21600]';
y_obs = [
0
1.6216e-01
2.9813e-01
4.0805e-01
4.9338e-01
5.5928e-01
6.0894e-01
6.5506e-01
6.8876e-01
7.1773e-01
7.4284e-01
7.6433e-01
7.8448e-01
7.9912e-01
8.1484e-01
8.2950e-01
8.3760e-01
8.4707e-01
8.5767e-01
8.6299e-01
8.6862e-01
8.7433e-01
8.7879e-01
8.8736e-01
8.9152e-01
8.9903e-01
9.0343e-01
9.0984e-01
9.1344e-01
9.1615e-01
9.2046e-01
9.2313e-01
9.2640e-01
9.2992e-01
9.3134e-01
9.3242e-01
9.3616e-01
9.3986e-01
9.4201e-01
9.4145e-01
9.4434e-01
9.4252e-01
9.4131e-01
9.4249e-01
9.4283e-01
9.4395e-01
9.4355e-01
9.4690e-01
9.4919e-01
9.5255e-01
9.5626e-01
9.6282e-01
9.6283e-01
9.6672e-01
9.6439e-01
9.6887e-01
9.7099e-01
9.7529e-01
9.8034e-01
9.8302e-01
9.8494e-01
9.8828e-01
9.8769e-01
9.9130e-01
9.9108e-01
9.9422e-01
9.9561e-01
9.9737e-01
1.0002e+00
1.0034e+00
1.0044e+00
1.00329
1.0];
d0 = 1e-10;
fun = @(d) ypred(d, t1) - y_obs;
best_d = lsqnonlin(fun, d0)
Local minimum possible.
lsqnonlin stopped because the size of the current step is less than
the value of the step size tolerance.
best_d =
1.000000000000000e-10
predicted_y = ypred(best_d, t1);
plot(t1, y_obs, '*', t1, predicted_y, '-');
legend({'observed', 'predicted'})
title('lsqnonlin data fitting')
function y_pred = ypred(d, t1)
a=0.0011;
gama = 0.01005;
L2 = zeros(14,1);
L3 = zeros(100,1);
L4 = zeros(100,1);
L5 = zeros(100,1);
S= zeros(73,1);
y_pred = zeros(73,1);
% t = 0;
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
format longE
k =1;
for t = t1(:).'
for n=0:1:100
L2(n+1) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n+1) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n+1)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n+1) = ((L2(n+1)*L3(n+1))/L4(n+1));
end
S((t/300) +1) = sum(L5);
y_pred((t/300)+1)= 1 -(L1*S((t/300) +1)); % predicted data
end
end
Walter Roberson
2021-6-21
But still notice that the fitting is giving the same d result as input. The fitting can only progress down to about 1e-8 with the default fitting options.
Anand Ra
2021-6-25
Hello,
All of a sudden, the curvesa are fitting. I have no idea what I did.
t1 = [0:300:21600]';
y_obs = [
0
1.6216e-01
2.9813e-01
4.0805e-01
4.9338e-01
5.5928e-01
6.0894e-01
6.5506e-01
6.8876e-01
7.1773e-01
7.4284e-01
7.6433e-01
7.8448e-01
7.9912e-01
8.1484e-01
8.2950e-01
8.3760e-01
8.4707e-01
8.5767e-01
8.6299e-01
8.6862e-01
8.7433e-01
8.7879e-01
8.8736e-01
8.9152e-01
8.9903e-01
9.0343e-01
9.0984e-01
9.1344e-01
9.1615e-01
9.2046e-01
9.2313e-01
9.2640e-01
9.2992e-01
9.3134e-01
9.3242e-01
9.3616e-01
9.3986e-01
9.4201e-01
9.4145e-01
9.4434e-01
9.4252e-01
9.4131e-01
9.4249e-01
9.4283e-01
9.4395e-01
9.4355e-01
9.4690e-01
9.4919e-01
9.5255e-01
9.5626e-01
9.6282e-01
9.6283e-01
9.6672e-01
9.6439e-01
9.6887e-01
9.7099e-01
9.7529e-01
9.8034e-01
9.8302e-01
9.8494e-01
9.8828e-01
9.8769e-01
9.9130e-01
9.9108e-01
9.9422e-01
9.9561e-01
9.9737e-01
1.0002e+00
1.0034e+00
1.0044e+00
1.00329
1.0];
d0 = 3.0e-10;
fun = @(d) ypred(d, t1) - y_obs;
best_d = lsqnonlin(fun, d0)
predicted_y = ypred(best_d, t1);
plot(t1, y_pred, '*', t1, predicted_y, '-');
legend({'observed', 'predicted'})
%title('lsqnonlin data fitting')
function y_pred = ypred(d, t1)
a=0.0011;
gama = 0.01005;
L2 = zeros(14,1);
L3 = zeros(100,1);
L4 = zeros(100,1);
L5 = zeros(100,1);
S= zeros(73,1);
y_pred = zeros(73,1);
L1 = ((8*gama)/((pi*(1-exp(-2*gama*a)))));
format longE
for t= t1(:).'
for n=0:1:100
L2(n+1) = exp((((2*n + 1)^2)*-d*pi*pi*t)/(4*a*a));
L3(n+1) = (((-1)^n)*2*gama)+(((2*n+1)*pi)*exp(-2*gama*a))/(2*a);
L4(n+1)= ((2*n)+1)*((4*gama*gama)+((((2*n)+1)*pi)/(2*a))^2);
L5(n+1) = ((L2(n+1)*L3(n+1))/L4(n+1));
end
S((t/300) +1) = sum(L5);
y_pred((t/300)+1)= 1 -(L1*S((t/300) +1)); % predicted data
end
end
Anand Ra
2021-6-25
However, I had to keep attempting the inital value to get the right number that would produce a fit. The optimized coefficient is same as my initial assumption.
When I tried with different datab set for y_obs, I am unable to find that perfect inital guess that would produce me a good fit.
Not sure what is going wrong.
Anand Ra
2021-6-26
Did I make any mistake like earlier with the code? Is there a way to get a good fit with an arbitrary initial guess?
更多回答（0 个）
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Nonspherical Models 的更多信息
标签
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!发生错误
由于页面发生更改，无法完成操作。请重新加载页面以查看其更新后的状态。
Translated by 
您也可以从以下列表中选择网站：
如何获得最佳网站性能
选择中国网站（中文或英文）以获得最佳网站性能。其他 MathWorks 国家/地区网站并未针对您所在位置的访问进行优化。
美洲
- América Latina (Español)
- Canada (English)
- United States (English)
欧洲
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)
亚太
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)