clear
clc
% atom_C
filename1 = 'atom_C.mat';
load(filename1);
distances = atom_C(:,1);
energies = atom_C(:,2);
opts = optimoptions('lsqcurvefit','Display','off');
n1 = 4;
%% atom_C, range_cutoff = 10;
range_fit = (1:21).'; range_cutoff = (5:16).';
parameters1_min = cell(length(range_cutoff),1);
parameters1_time = cell(length(range_cutoff),1);
parameters1_rms = cell(length(range_cutoff),1);
%%
tic
for ii = 1:length(range_cutoff)
hfunc = @(x,xdata)(x(1)*(xdata+x(2)).^2+x(3));
% hfunc = @(x,xdata) (x(1)*(1-exp(x(2)*(xdata+x(3)))).^2+x(4));
range_fit1 = range_fit(1:range_cutoff(ii));
p1_min = []; p1_time = []; rmstm1 = inf;
for kk = 1:1e3
p01 = (rand(1,n1)-0.5)*2e-0;
[p11,resnorm1,~,~,~] = lsqcurvefit(hfunc, ...
p01,distance1o(range_fit1),energyo(range_fit1),[],[],opts);
if resnorm1<rmstm1
p1_min = [p1_min; p11];
p1_time = [p1_time; kk];
rmstm1 = resnorm1;
end
end
parameters1_min{ii} = p1_min;
parameters1_time{ii} = p1_time;
end
runtime = toc
%%
resnorms = [];
for ii = 1:length(range_cutoff)
range_fit1 = range_fit(1:range_cutoff(ii));
resnorm_tm1 = mean(((hfunc(parameters1_min{ii}(end,:),distance1o(range_fit1)) ...
-energyo(range_fit1))).^2);
resnorms = [resnorms; resnorm_tm1];
end
p1_ho.hfunc = hfunc;
p1_ho.range_fit = range_fit;
p1_ho.range_cutoff = range_cutoff;
p1_ho.parameters1_min = parameters1_min;
p1_ho.parameters1_time = parameters1_time;
p1_ho.resnorms = resnorms;
% p1_morse.hfunc = hfunc;
% p1_morse.range_fit = range_fit;
% p1_morse.range_cutoff = range_cutoff;
% p1_morse.parameters1_min = parameters1_min;
% p1_morse.parameters1_time = parameters1_time;
% p1_morse.resnorms = resnorms;
%%
% clf
% plot(distance1o,energyo,'o')
% hold on
% ii = 1;
% distance1 = (distance1o(1):1e-2:distance1o(range_cutoff(ii))).';
% distance_diff = abs(distance1-distance1o(range_cutoff(ii)));
% loca = find(distance_diff==min(distance_diff));
% plot(distance1,hfunc(parameters1_min{ii}(end,:),distance1))