Im trying to run fsolve inside a for loop but it seem that my second iteration inside the for loop is not able to read the initial trial point

1 ビュー (過去 30 日間)
John Harvie delos Reyes
John Harvie delos Reyes 2021 年 11 月 25 日
コメント済み: Matt J 2021 年 11 月 25 日
Code below:
clear;
clc;
%%
Disp = [1 10 4 2 3;
2 8 3 2 2;
3 5 4 1.5 3];
options = optimset('TolFun',1e-12,'MaxFunEvals',1e9,'Maxiter',1e9);
options.StepTolerance =1e-12;
%%
N = size(Disp,1);
O = [0,0];
tau = 0.3;
Bdispy = 0;
Cdispy = 0;
GG = zeros(2,N);
exitflags = zeros(1,N);
for i = 1:N
Peakdispx(i) = Disp(i,2);
Peakdispy(i) = Disp(i,3);
Adispy(i) = Peakdispy(i);
mslope(i) = (Disp(i,5)-O(2))/(Disp(i,4)-O(1));
Cdispx(i) = Peakdispx(i)-((Peakdispy(i)-Cdispy)/(mslope(i)));
Bdispx(i) = Cdispx(i)*-1;
Adispx(i) = (Adispy(i)-Bdispy)/(mslope(i))+Bdispx(i);
slopeACdisp(i) = ((Adispy(i)-Cdispy)/(Adispx(i)-Cdispx(i)));
Afull(i) = (0.5*abs(Adispx(i)*Bdispy+Bdispx(i)*Cdispy+Cdispx(i)*Peakdispy(i)+Peakdispx(i)*Adispy(i)-Adispy(i)*Bdispx(i)-Bdispy*Cdispx(i)-Cdispy*Peakdispx(i)-Peakdispy(i)*Adispx(i)));
F = @(GG)[slopeACdisp(i)-((GG(2)-Cdispy)/(GG(1)-Cdispx(i)));
-tau+(0.5*abs((GG(1)*Bdispy+Bdispx(i)*Cdispy+Cdispx(i)*Peakdispy(i)+Peakdispx(i)*GG(2)-GG(2)*Bdispx(i)-Bdispy*Cdispx(i)-Cdispy*Peakdispx(i)-Peakdispy(i)*GG(1))))/Afull]
gg0 = [tau*Disp(i,3);(1-tau)*Disp(i,2)];
GG = fsolve(F,gg0);
%% pinch in x and y
disppinchy(i) = GG(2)/Disp(i,3);
U(i) = ((GG(2)-Cdispy)/mslope(i))+Cdispx(i);
L(i) = ((GG(2)-Bdispy)/mslope(i))+Bdispx(i);
disppinchx(i) = (GG(2)-L(i))/(U(i)-L(i));
end
for i = 1:N
Gdispx(i) = GG(1);
Gdispy(i) = GG(2);
Apin(i) =0.5*abs(Gdispx(i)*Bdispy+Bdispx(i)*Cdispy+Cdispx(i)*Peakdispy(i)+Peakdispx(i)*Gdispy(i)-Gdispy(i)*Bdispx(i)-Bdispy*Cdispx(i)-Cdispy*Peakdispx(i)-Peakdispy(i)*Gdispx(i));
ratio(i) = Apin(i)/Afull(i)
end

回答 (1 件)

Matt J
Matt J 2021 年 11 月 25 日
編集済み: Matt J 2021 年 11 月 25 日
Should the division be by Afull(i)?
F = @(GG)[slopeACdisp(i)-((GG(2)-Cdispy)/(GG(1)-Cdispx(i)));
-tau+(0.5*abs((GG(1)*Bdispy+Bdispx(i)*Cdispy+Cdispx(i)*Peakdispy(i)+Peakdispx(i)*GG(2)-GG(2)*Bdispx(i)-Bdispy*Cdispx(i)-Cdispy*Peakdispx(i)-Peakdispy(i)*GG(1))))/Afull(i)]
  2 件のコメント
Matt J
Matt J 2021 年 11 月 25 日
Since it's a 2D function, you can do a grid search for a more accurate initial guess.

サインインしてコメントする。

製品


リリース

R2018a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by