fsolve with two variables in a loop

2 visualizzazioni (ultimi 30 giorni)
Mei Cheng
Mei Cheng il 30 Gen 2023
Commentato: Mei Cheng il 31 Gen 2023
Hi, I thank you in advance.
I try to solve two equation with two variables (x and y) with a changing parameters T, i try to use a fsolve method. I met a problem in the loop. Below is my code.
clear all
% effective normal stress (Pa)
p.sigma=5e7;p.To=273.15+20;p.R=8.314;p.K=6e10;p.Lt=8e-4;p.lambda=0.0625;
p.phio=0.32;p.dsb=2e-6;p.dbulk=2e-5;p.qsb=0.4;p.phic=0.2;p.phio=0.02;
p.qbulk=0.7;p.H=0.577;p.n=1.7;p.m=3;p.An=0.0759;p.At=4.4*p.An;
p.Ea=213e3;p.omega=3.69e-7;p.muo=0.6;p.ao=0.006;p.xo = 0.1813;
p.V0=1e-6;p.V1=0.1e-6;p.gammao=0.02;p.M = 2;p.N = 1;
% a changing parameter T with two variables x and y
T=linspace(273.15,600,100);
x=linspace(0,0.2,100);
y=linspace(0,0.2,100);
% Equation 1: At reference (prestep) velocity of V0
fx=@(x,T) p.Lt.*p.lambda./p.V0*p.An.*exp(-p.Ea./p.R./T).*...
(p.sigma).^p.n./(p.dsb).^p.m.*((p.phic-x)./(x-p.phio)).^(-p.M)-p.H.*(p.qsb-2.*x).^p.N;
% Equation 2: At a velocity (poststep) of V1
fy=@(Y,T) p.Lt.*p.lambda./p.V1*p.An.*exp(-p.Ea./p.R./T).*...
(p.sigma).^p.n./(p.dsb).^p.m.*((p.phic-y)./(y-p.phio)).^(-p.M)-p.H.*(p.qsb-2.*y).^p.N;
musi=@(T) p.muo+p.ao*T./p.To.*log(p.V1./p.V0);
muso=p.muo;
tanso=@(x) p.H.*(p.qsb-2.*x).^p.N;
tansi=@(y) p.H.*(p.qsb-2.*y).^p.N;
aa=@(x,T) p.ao.*T./p.To.*(1+tanso(x).^2)./(1-muso.*tanso(x))./(1-musi(T).*tanso(x));
bb=@(y,T) -tansi(y).*(1+musi(T).^2)./(1-musi(T).*tansi(y))./...
(1-musi(T).*tansi(y)).*(1-(p.V1./p.V0).^(1/(p.M+p.N)))./log(p.V1./p.V0);
% loop
options=odeset('Refine',1,'RelTol',1e-15,'InitialStep',1e-6,'MaxStep',3e6);
for j=1:length(T)
fun=@(x,y)[fx(x,T(j));fy(y,T(j))];
z(:,j)=fsolve(fun,[0.199;0.199],options); %plotted to get initial guess close.
end
Not enough input arguments.

Error in solution (line 33)
fun=@(x,y)[fx(x,T(j));fy(y,T(j))];

Error in fsolve (line 264)
fuser = feval(funfcn{3},x,varargin{:});

Caused by:
Failure in initial objective function evaluation. FSOLVE cannot continue.
% figure
plot(T-273.15,aa(z,T),'k-',T-273.15,aa(z,T),'r-');
  2 Commenti
Torsten
Torsten il 30 Gen 2023
If the solutions from the symbolic code did not satisfy your constraints, using "fsolve" won't help.
The three solutions you got for x and y from the symbolic approach are the only ones your equations have.
Mei Cheng
Mei Cheng il 30 Gen 2023
@Torsten thanks for your remind. Previously i used the fsolve to get a good result, but it is for one variable. I guess it should work for two variables.

Accedi per commentare.

Risposta accettata

Matt J
Matt J il 30 Gen 2023
fun=@(z)[fx(z(1),T(j));fy(z(2),T(j))];
  5 Commenti
Matt J
Matt J il 30 Gen 2023
Since z(:,i) are column vectors, perhaps T should be a column-vector as well?
plot(T-273.15,aa(z(:,1),T(:)),'k-',T-273.15,bb(z(:,2),T(:)),'r-');
Mei Cheng
Mei Cheng il 31 Gen 2023
@Matt J thank you very much for your suggestions. I have reivsed my code accordingly.

Accedi per commentare.

Più risposte (0)

Categorie

Scopri di più su Stress and Strain in Help Center e File Exchange

Prodotti


Release

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by