Asked by sophp
on 4 Feb 2018

how to find the intersection points of dH_rem and dH_gen within the limits specified below

To=500:1:850; %outlet temp

Ti=570; %inlet temp

y_A=0.003; %proportion of benzene in feed

V=8.5;

P=3e5; %Pa

R=8.3145; %kJ/mol.K

Cp=1.09 %kJ/kgK

Mr_air=29e-3; %kg/mol

dH_1=-1850; %kJ/mol

dH_2=-1423;

dH_3=-3273;

F_ao=0.1; %molar flow rate of benzene

k1 = 1e7.*exp(-12700./To);

k2 = 5e4.*exp(-10800./To);

k3 = 7e7.*exp(-15000./To);

t=(V.*y_A.*P)./(F_ao.*R.*To);

X=(t.*(k1+k3))./(1+t.*(k1+k3));

S=k1./((1+(k2.*t)).*(k1+k3));

Y_B = (k1.*t)./(((k2.*t)+1).*(1+(t.*(k1+k3))));

dH_gen=-((V.*y_A.*P)./(R.*To)).*((((k1.*dH_1)+(k3.*dH_3)).*(1-X))+(k2.*X.*S.*dH_2));

plot(To,dH_gen)

hold on

dH_rem=(Mr_air./y_A).*F_ao.*Cp.*(To-Ti);

plot(To,dH_rem)

hold off

Answer by John D'Errico
on 4 Feb 2018

Standard question: how to find the intersection(s) of two functions.

1. Subtract them. Where the difference is zero, there lies an intersection.

2. Use a root finder. That could be anything from fzero, solve, vpasolve, fsolve, etc.

Note that all standard optimization based root finders will find ONE root, and only one root. It will depend on your starting values. Solve might be able to find the three points of intersection I saw on the plot.

The lazy solution to finding an approximate root is to use a tool like Doug Schwarz's intersections code on the FEX. Evaluate each function at a few hundred points, then call intersections. It will use linear interpolation to find the crossing points. But a virtue of that solution is it will report all crossings between the curves.

Walter Roberson
on 4 Feb 2018

When you set To to numeric at the beginning of the code, you created a numeric vector for dH_gen and dH_rem: the numeric values of To are copied as needed to form the expressions. Then you assign a symbol to To and suddenly expect that you can solve() the numeric vector for a symbol.

You need to do the syms To in place of the To=500:1:850;

If you also syms Ti in place of assigning 370 to it, and then go through the solve() you will get

dH_rem =

(3161*To)/3000 - (3161*Ti)/3000

You can see at a glance that the equation will be satisfied when To == Ti .

Notice that you compute Y_B but never use it.

sophp
on 4 Feb 2018

Thank you!! Okay, so I replaced

To=500:1:850

with

syms To

Matlab then tells me that the equation dH_rem-dH_gen is unsolvable symbolically and gives me a numerical approximation instead.

The solution is correct but how do I alter the code to give the lowest intersection value.

Walter Roberson
on 4 Feb 2018

Sign in to comment.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.