fsolve with two variables in a loop
5 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
Mei Cheng
el 30 de En. de 2023
Comentada: Mei Cheng
el 31 de En. de 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
% figure
plot(T-273.15,aa(z,T),'k-',T-273.15,aa(z,T),'r-');
2 comentarios
Torsten
el 30 de En. de 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.
Respuesta aceptada
Más respuestas (0)
Ver también
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!