MATLAB Answers

Will Meg
0

Fit theoretical x-y data set to experimental data set

Will Meg
さんによって質問されました 2019 年 7 月 5 日
最新アクティビティ Star Strider
さんによって 編集されました 2019 年 7 月 8 日
I have a function to model the force-displacement relation between an interface and a probe. The function outputs two data sets (the probe movement and corresponding force) based on three inputs. I have experimental data and essentially want to find what value for the inputs results in the closest fit according to my model. The problem looks like this,
function [f,d] = myFunc(x)
where, f and d are arrays of data points that I want to fit as closely as possible to the experimental data i have, [f_exp, d_exp].
The crux of the problem to me is I am not sure how to constrain my outputs in order to data mine the inputs. I have so far been trying to use the genetic algorithm function in the optimisation tool box but am not sure if this is the correct approach.
Thanks for any help.
Cheers,
Will

  4 件のコメント

Will Meg
2019 年 7 月 6 日
As a starting point I would just like unconstrained optimisation to get a good fit to my data. My parameter 'x' consists of three parameters x(1), x(2) and x(3), which are not constrained. But want the output data sets to be as close a match as possible to the experimental data i have obtained.
Image Analyst
2019 年 7 月 6 日
Will, where are the screenshots? Where did you attach the data? Why not make it easy for people to answer you, not hard?
If you want it to fit your experimental data perfectly, then use spline interpolation. If you want a regression, then there are several ways but you might have to pick a model equation (which we have no guesses at since you didn't attach a screenshot).
Will Meg
2019 年 7 月 8 日
So my function solves the augmented Young-Laplace equation using a 4th order runga kutta numerical method, and based on this solution I can get force displacement data, as shown in my script below;
function [FX0] = RK4AL_func(myin)
%% Augmented Young Laplace
counter = 0;
for d0 = 160:1.5:1558
hold off
counter = counter+1;
A=4e-20;
e=80; % Relative permitivity water
e0=8.85e-12; % Permitivity of free space
k=100; % debeye length nm
dRho = 150; % delta density kg/m^3 (same as mg/cm^3)
zp1=myin(1); % zeta potential 1
zp2=myin(2); % zeta potential 2
%d0=195; % minimum interfacial gap nm
a=10e-6; % bead radius nm
g=myin(3); % surface tension N/m
h=0.1; % step size
x = 0.1:h:25; % Calculates upto y(1)
y = zeros(1,length(x));
z = zeros(1,length(x));
y(1) = d0; % initial condition
z(1) = -0.000000000000001; % initial condition
% F_xy = @(t,r) 3.*exp(-t)-0.4*r; % change the function as you desire
F_xyz = @(x,y,z) z; % change the function as you desire
%G_xyz = @(x,y,z) (2-((20*(((80*8.85e-12*(2*-0.02*-0.02*(exp(-1*y/96.2)+exp(-3*y/96.2))+(((-0.02)^2+(-0.02)^2)*exp(-2*y/96.2))))/((96.2^2)*((1-exp(-2*y/96.2))^2)))*1e18)-(4e-20*1e27/(6*pi*(y^3))))/50))-(z/x);
G_xyz = @(x,y,z) d0*(2-((a*(((e*e0*(2*zp1*zp2*(exp(-1*y/k)+exp(-3*y/k))+(((zp1^2)+(zp2^2))*exp(-2*y/k))))/((k^2)*((1-exp(-2*y/k))^2)))*1e18)-(A*1e27/(6*pi*(y^3))))/g))-(z/x);
%DP(i) = (((e*e0*(2*zp1*zp2*(exp(-1*y(i)/k)+exp(-3*y(i)/k))+((zp1^2+zp2^2)*exp(-2*y(i)/k))))/((k^2)*((1-exp(-2*y(i)/k))^2)))*1e18)-(A*1e27/(6*pi*(y(i)^3)));
for i=1:(length(x)-1) % calculation loop
k_1 = F_xyz(x(i),y(i),z(i));
L_1 = G_xyz(x(i),y(i),z(i));
k_2 = F_xyz(x(i)+0.5*h,y(i)+0.5*h*k_1,z(i)+0.5*h*L_1);
L_2 = G_xyz(x(i)+0.5*h,y(i)+0.5*h*k_1,z(i)+0.5*h*L_1);
k_3 = F_xyz((x(i)+0.5*h),(y(i)+0.5*h*k_2),(z(i)+0.5*h*L_2));
L_3 = G_xyz((x(i)+0.5*h),(y(i)+0.5*h*k_2),(z(i)+0.5*h*L_2));
k_4 = F_xyz((x(i)+h),(y(i)+k_3*h),(z(i)+L_3*h)); % Corrected
L_4 = G_xyz((x(i)+h),(y(i)+k_3*h),(z(i)+L_3*h));
y(i+1) = y(i) + (1/6)*(k_1+2*k_2+2*k_3+k_4)*h; % main equation
z(i+1) = z(i) + (1/6)*(L_1+2*L_2+2*L_3+L_4)*h; % main equation
end
%% Drop Profile
hold on
DrPr = (y-d0)/d0; %DrPr is the drop profile (D(t)-D0)/D0
plot(x,DrPr)
%% Finding G, H, F and X
for j = 1:length(x)
DP(j) = (((e*e0*(2*zp1*zp2*(exp(-1*y(j)/k)+exp(-3*y(j)/k))+((zp1^2+zp2^2)*exp(-2*y(j)/k))))/((k^2)*((1-exp(-2*y(j)/k))^2)))*1e18)-(A*1e27/(6*pi*(y(j)^3)));
dG(j) = a*d0*x(j)*DP(j)/g;
dH(j) = a*d0*x(j)*log(x(j))*DP(j)/g;
end
for jj = 2:length(x);
G(jj) = (dG(jj)+dG(jj-1))*h/2;
H(jj) = (dH(jj)+dH(jj-1))*h/2;
end
Gsum(counter) = sum(G);
Hsum(counter) = sum(H);
B = 0.57721566 + log(sqrt(a)/2*sqrt(g/(dRho*9.81)));
FX0(counter,1) = (d0 - Hsum(counter) - Gsum(counter)*((1/2*log(d0))+B))-160; %Piezo Movement, X0
FX0(counter,2) = (2*pi*g*Gsum(counter)/1000)*1e6; %Force observed by AFM, F
end
end
I then have experimental force displacement data that I would like to get as close a fit as possible to by varying three inputs to my function. The force - displacement data is essentially just x-y data and looks like this.
Capture.JPG
So I want to am trying to use the optimization function so that the exported data from my funciton matches this data as closely as possible.

サインイン to comment.

1 件の回答

Star Strider
回答者: Star Strider
2019 年 7 月 6 日
編集済み: Star Strider
2019 年 7 月 8 日

I have no idea what your actual function is, or what you are doing.
I would do something like this:
function fd = myFunc(x)
fd(:,1) = ...; % Approximate ‘f_exp’
fd(:,2) = ...; % Approximate ‘d_exp’
end
fitnessfcn = @(x) norm([f_exp(:), d_exp(:)] - myFunc(x));
Ideally, if you are fitting data, ‘myFunc’ should be a function of your independent variable and the parameter vector. I have no idea how ‘f_exp’ relates to ‘d_exp’, or to your parameters, so this is illustrative only.
EDIT —
Since you’re fitting differential equations to data, see:

  0 件のコメント

サインイン to comment.



Translated by