How to iteratively solve when equations are dependent on each other.

1 次查看（过去 30 天）
John
John 2025-7-27
回答： Sam Chak 2025-7-28
I have two equations, I know the varibles T_oQlast, K_dot_nom, and Sigma_ysrtM which are attached. The top equation is typically used to solve for T_oQlast, but in my case I need to solve for T_o_calc so I rearranged it as you can see. GAMMA is calculated based on T_o_calc.
What funtion can I used to solve for T_o_calc? I have thought about using an interative solution, select a starting GAMMA1 of 50, then calculate T_o_calc for that GAMMA1. Recacluate GAMMA1 for this new T_o_calc. But I can't figure out a way to loop the code until it converges on a solution.
T_o_calc = (((T_oQlast+273.15)*GAMMA1)/(GAMMA1-log(K_dot_nom/0.91)))-273.15;
GAMMA1 = (9.9 * exp((((T_o_calc+273.15)/190)^1.66)+(((Sigma_ysrtM)/722)^1.09)));
Heres the Vpasolve code I tied out. equLeft is the T_o_calc equation solved for GAMMA1, equRight is the equation for GAMMA1. This didn't work for some reason, it gave back no results.
syms T_o_calc
eqnLeft = -(273.15*log(K_dot_nom/0.91)-T_oQlast*log(K_dot_nom/0.91))/(T_o_calc-T_oQlast);
eqnRight = (9.9 * exp((((T_o_calc+273.15)/190)^1.66)+(((Sigma_ysrtM)/722)^1.09)));
T_o = vpasolve(eqnLeft == eqnRight, T_o_calc);
  4 个评论
显示 2更早的评论
Sam Chak
Sam Chak 2025-7-27
Ran in:
Hi @John, Based on your corrected equations, there is an intersection around T_o_calc = -170. You can try out the approaches by @Star Strider and @Torsten.
load('Sigma_ysrtM.mat')
load('T_oQlast.mat')
load('K_dot_nom.mat')
syms T_o_calc
eqnLeft = -(273.15*log(K_dot_nom) - T_oQ*log(K_dot_nom))/(T_o_calc - T_oQ);
eqnRight = (9.9*exp((((T_o_calc + 273.15)/190)^1.66) + (((Sigma_ysrtM)/722)^1.09)));
hold on
fplot(eqnLeft, [-220, -140])
fplot(eqnRight, [-220, -140])
hold off
ylim([0, 100])
legend('eqnL', 'eqnR', 'location', 'northwest')
grid on
John
John 2025-7-28
Third times the charm here, I broke things down to make sure each indivdual term was correct and found another error.
syms T_o_calc
Sigma_ysrtM = sigma_ysrt * 6.895;
K_dot_nom_Metric = table2array(inputs(1,10)); %attached is this value
lnK = log(K_dot_nom_Metric);
Term1 = 273.15 * lnK;
Term2 = T_oQ * lnK;
Term3 = ((Sigma_ysrtM)/722)^1.09;
eqnLeft = -(Term1+Term2)/(T_o_calc-T_oQ);
eqnRight = 9.9*exp((((T_o_calc+273.15)/190)^1.66) + (Term3));
T_o = vpasolve(eqnLeft == eqnRight, T_o_calc);
Answer is -145.3 if you plot it, which matches my excel solver, but vpasolve doesn't work still.

请先登录，再进行评论。

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

采纳的回答

Sam Chak
Sam Chak 2025-7-28
Ran in:
Hi @John
You can specify the range so that 'vpasolve' searches the solution in that region.
load('Sigma_ysrtM.mat')
load('T_oQlast.mat')
load('K_dot_nom_Metric')
syms T_o_calc
lnK = log(K_dot_nom_Metric);
Term1 = 273.15 * lnK;
Term2 = T_oQ * lnK;
Term3 = ((Sigma_ysrtM)/722)^1.09;
eqnLeft = -(Term1+Term2)/(T_o_calc-T_oQ);
eqnRight = 9.9*exp((((T_o_calc+273.15)/190)^1.66) + (Term3));
hold on
fplot(eqnLeft, [-200, -130])
fplot(eqnRight, [-200, -130])
hold off
ylim([0, 100])
legend('eqnL', 'eqnR', 'location', 'northwest')
grid on
●
sol = vpasolve(eqnLeft == eqnRight, T_o_calc, [-150 140])
sol = 

更多回答（2 个）

Torsten
Torsten 2025-7-27
编辑：Torsten 2025-7-27
Ran in:
Sigma_ysrtM = load("Sigma_ysrtM.mat").Sigma_ysrtM;
T_oQlast = load("T_oQlast.mat").T_oQ;
K_dot_nom = load("K_dot_nom.mat").K_dot_nom;
T_o_calc = fsolve(@(T_o_calc)fun(T_o_calc,Sigma_ysrtM,T_oQlast,K_dot_nom),300)
Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient.
T_o_calc = -123.4789
function res = fun(T_o_calc,Sigma_ysrtM,T_oQlast,K_dot_nom)
GAMMA1 = (9.9 * exp((((T_o_calc+273.15)/190)^1.66)+(((Sigma_ysrtM)/722)^1.09)));
res = T_o_calc - ((((T_oQlast+273.15)*GAMMA1)/(GAMMA1-log(K_dot_nom/0.91)))-273.15);
end
  1 个评论
John
John 2025-7-28
T_o_calc = fsolve(@(T_o_calc)fun(T_o_calc,Sigma_ysrtM,T_oQlast,K_dot_nom_Metric),300)
function res = fun(T_o_calc,Sigma_ysrtM,T_oQlast,K_dot_nom_Metric)
GAMMA1 = (9.9 * exp((((T_o_calc+273.15)/190)^1.66)+(((Sigma_ysrtM)/722)^1.09)));
res = T_o_calc - ((((T_oQlast+273.15)*GAMMA1)/(GAMMA1-log(K_dot_nom_Metric)))-273.15);
end
I get -123.4790, but it should be aroun -145.3

请先登录，再进行评论。

John
John 2025-7-28
编辑：John 2025-7-28
syms T_o_calc
Sigma_ysrtM = sigma_ysrt * 6.895;
K_dot_nom_Metric = table2array(inputs(1,10)); %Load from attached
lnK = log(K_dot_nom_Metric);
Term1 = 273.15 * lnK;
Term2 = T_oQ * lnK;
Term3 = ((Sigma_ysrtM)/722)^1.09;
eqnLeft = -(Term1+Term2)/(T_o_calc-T_oQ);
eqnRight = 9.9*exp((((T_o_calc+273.15)/190)^1.66) + (Term3));
T_o = vpasolve(eqnLeft == eqnRight, T_o_calc,[-500,500]); %ADDED bounds
Returns correct value and NOT an imaginary number which is the issue I had with the origianl code.

类别

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

标签

产品


版本

R2023b

Community Treasure Hunt

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

Start Hunting!

Translated by