selecting the positive result

9 次查看（过去 30 天）
Domenico Fazzari
Domenico Fazzari 2016-12-23
评论： Domenico Fazzari 2017-2-6
Hi everyone! I have a question on fsolve:i did a code for solving with fsolve 2eqs in 2 variables,but one of the variable is squared. What happens is that fsolve gives me the negative root of that variable,instead of the positive one. How can i achieve that? Maybe with some options? thanks.
  2 个评论
Image Analyst
Image Analyst 2016-12-23
So you're basically finding where a line intersects a parabola? Why doesn't it give all solutions? There should be two solutions unless your intersection point is perpendicular or tangent to the vertex of the parabola. Let's see your code - you forgot to include it.
Domenico Fazzari
Domenico Fazzari 2017-2-6
Thanks for your answer and sorry for the late reply! Let me explain better the problem: I've a system of these two eqs:
1) exp(-G/(R*Tout))= -P*(L^(2) / (1-L)*(11+L)) the unks=>Tout,L 2) Hout-Hin=0 , which is a function of the same unk=>Tout,L
As results I have a negative L, which has not a physical meaning, and as a consequence an unexpected Tout.
This is the function code (fsolve), hope you can understand it well (take in account this is one my second code ever, so maybe it is not well written).
global R att Pin nin cp deltahf Tin % Initialization y=[0 0 0]; Tout=a(2); lambda=a(1); Hin=0; Hout=0; hcompostientranti=[0 0 0 0]; %Enthalpies of feeding compounds hcompostiuscenti=[0 0 0]; % Enthaplies of exiiting compounds %Exiting moles nout=[(1-lambda) (lambda) (lambda)]; % total exiting moles ntout=11+lambda; %mol fract in function of lambda for i=1:length(nout) y(i)=nout(i)/ntout; end % Activities , ideal mix of ideal gas, ai=P*yi,species={'EB' 'STY' 'H2'}; att=[Pin*y(1) Pin*y(2) Pin*y(3)];
G=(127619.9-134.02*Tout);
Keq=exp(-G/(R*Tout));
% Stoichiometric vector EB->STY+H2
vettorestechio1=[-1 +1 +1];
% Prod Activities
prodatt=(att(1)^(vettorestechio1(1)))*(att(2)^(vettorestechio1(2)))*(att(3)^(vettorestechio1(3))); % Produttore attività
% Feeding Enthalpies
for i=1:length(nin)
hcompostientranti(i)=nin(i)*(deltahf(i)+(Tin-298)*cp(i));
Hin=Hin+hcompostientranti(i);
end
% Exit Enthalpies
for j=1:length(nout)
hcompostiuscenti(j)=nout(j)*(deltahf(j)+(Tout-298)*cp(j));
Hout=Hout+hcompostiuscenti(j);
end
reattore(2)=Keq-prodatt; reattore(1)=Hin-Hout; end

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

Aylin
Aylin 2016-12-27
编辑：Aylin 2016-12-27
Hello Domenico,
If you are trying to solve a polynomial equation, then you may want to use the roots function instead. However, if the fsolve function works best for you, please refer to the following MATLAB documentation page for some techniques on solving nonlinear systems with constraints:
https://www.mathworks.com/help/optim/ug/nonlinear-systems-with-constraints.html
In particular, if you would like to get only positive zeros from your code, you can run fsolve for a range of starting values, and then use only the starting values that give you positive results.
Rylan
  1 个评论
Domenico Fazzari
Domenico Fazzari 2017-2-6
Thanks for your answer and sorry for the late reply! Let me explain better the problem: I've a system of these two eqs:
1) exp(-G/(R*Tout))= -P*(L^(2) / (1-L)*(11+L)) the unks=>Tout,L 2) Hout-Hin=0 , which is a function of the same unk=>Tout,L
As results I have a negative L, which has not a physical meaning, and as a consequence an unexpected Tout.
This is the function code (fsolve), hope you can understand it well (take in account this is one my second code ever, so maybe it is not well written).
global R att Pin nin cp deltahf Tin % Initialization y=[0 0 0]; Tout=a(2); lambda=a(1); Hin=0; Hout=0; hcompostientranti=[0 0 0 0]; %Enthalpies of feeding compounds hcompostiuscenti=[0 0 0]; % Enthaplies of exiiting compounds %Exiting moles nout=[(1-lambda) (lambda) (lambda)]; % total exiting moles ntout=11+lambda; %mol fract in function of lambda for i=1:length(nout) y(i)=nout(i)/ntout; end % Activities , ideal mix of ideal gas, ai=P*yi,species={'EB' 'STY' 'H2'}; att=[Pin*y(1) Pin*y(2) Pin*y(3)];
G=(127619.9-134.02*Tout);
Keq=exp(-G/(R*Tout));
% Stoichiometric vector EB->STY+H2
vettorestechio1=[-1 +1 +1];
% Prod Activities
prodatt=(att(1)^(vettorestechio1(1)))*(att(2)^(vettorestechio1(2)))*(att(3)^(vettorestechio1(3))); % Produttore attività
% Feeding Enthalpies
for i=1:length(nin)
hcompostientranti(i)=nin(i)*(deltahf(i)+(Tin-298)*cp(i));
Hin=Hin+hcompostientranti(i);
end
% Exit Enthalpies
for j=1:length(nout)
hcompostiuscenti(j)=nout(j)*(deltahf(j)+(Tout-298)*cp(j));
Hout=Hout+hcompostiuscenti(j);
end
reattore(2)=Keq-prodatt; reattore(1)=Hin-Hout; end

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Programming 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by