Optmization: Need to optimize Tg in the following code w.r.t. dtc.
2 views (last 30 days)
Show older comments
Devdatt Thengdi
on 8 Feb 2018
Answered: Devdatt Thengdi
on 13 Feb 2018
Ti = 251.9 + 273.15; %Fluid in degK
T0 = 386.1 + 273.15; %Fluid/gas out degK
cpl = 0.32; %Cp fluid KJ/kgdegK
LHV = 41.05; %Fuel value KJ/kg
flowL = 500; %flowrate fluid kg/s
%Air
ea = 0.15; %excess air
cpair = 0.171; %Kcal/kgdegK
Tairin = 114;
Tref = 41; %Ambient
%Fuel
THD = 129674.5; %Heat released by fuel Kw
Ef = 0.9320; %Fuel efficiency
atfr = 1.344; %in notebook
%Dimensions
ctc = 0.3048; %ctc in m
Lh = 19.236; %Exposed tube length- horizontal m
Lv = 18.431; %Exposed tube length- vertical m
L = Lh + Lv;
d0 = 0.168; %Tube diameter m
pi = 22/7;
Tw = Tref + (Ti+T0)/L; %Tube wall temperature
Tstack = 1212; %Stack exit temperature
qfuel = THD/Ef;
mfuel = qfuel/LHV;
mair = mfuel*atfr*(1+ea);
qair = mair*cpair*(Tairin - Tref);
qloss = 0.05*qfuel;
mstack = mair + mfuel;
cpsteam = 4.1855;
cpco2 = 0.918;
cpn2 = 1.044;
xco2 = 0.054;
xsteam = 0.034; %Total should be 8.8
xn2 = 0.66;
cpgas = xsteam*cpsteam + xco2*cpco2 + xn2*cpn2;
qout = mstack*cpgas*(Tg - Tstack);
qr = qfuel + qair - (qout + qloss);
%Cold plane area
Ar = pi*L*dtc;
%absorptivity
hp = 6.626*10^-34; %m2kg/s
c = 3*10^8; %velocity of light
alpha = 4*pi*qr/hp*c;
%Exchange factor
x = pi*d0/2; %ctc
F = (sqrt(x^2 - 4) - x + 2*sinh(2/x))/2*pi;
%Heat transfer coefficient
k = 0.25*10^-3; %Thermal conductivity
rof = 1320; %density
Gm = 920880; %Dtc. Do I need to opmass flowrate
mu = 130; %Viscosity Pa-s
h0 = (do/k)*(0.023*((Gm*d0/mu*Ai)^0.8)*((cpl*mu/k))^0.3);
sig = 5.67*10^-8; %Stefan-Boltzmann constant
qr = sig*aplha*Ar*F*(Tg^4 - Tw^4) + Ar*h0*(Tg - Tw);
qfuel + qair - (mstack*cpgas*(Tg - Tstack) + qloss) = (pi*L*dtc*h0*(Tg - Tw))/(1-(sig*4*(pi^2)*dtc*L*F*((Tg^4)-(Tw^4)))/hp*c);
I need to optimize Tg w.r.t. dtc or vice versa. How do I do the same without rearranging? (Since, we cannot do that in this case)
1 Comment
Accepted Answer
Walter Roberson
on 8 Feb 2018
Is
qfuel + qair - (mstack*cpgas*(Tg - Tstack) + qloss) = (pi*L*dtc*h0*(Tg - Tw))/(1-(sig*4*(pi^2)*dtc*L*F*((Tg^4)-(Tw^4)))/hp*c);
an equation to be solved for Tg ? If so, then
Tg = roots([4*F*L*Pi^2*c*cpgas*dtc*mstack*sig, -4*F*L*Pi^2*c*dtc*sig*(Tstack*cpgas*mstack+qair+qfuel-qloss), 0, 0, -(4*(F*L*Pi^2*Tw^4*c*dtc*sig+(1/4)*hp))*cpgas*mstack-L*Pi*dtc*h0*hp, (4*(F*c*Tw^3*sig*(Tstack*cpgas*mstack+qair+qfuel-qloss)*Pi+(1/4)*h0*hp))*Pi*Tw*L*dtc+hp*(Tstack*cpgas*mstack+qair+qfuel-qloss)])
It is a degree 5 polynomial, so there will be 5 outputs, and it is likely that some of them will be complex.
6 Comments
Walter Roberson
on 8 Feb 2018
Unfortunately roots() does not accept an array of coefficients, so you would have to loop over all of the dtc values and roots() for each one, storing the result. Store the 5 results as 5 different columns of an array, with the array having the same number of rows as dtc has entries. Next, post-process that to drop out imaginary results,
vals(abs(imag(vals)) > 1e-10) = nan; %arbitrarily say that small enough imaginary components are just round-off error
vals = real(vals);
plot(dtc, vals)
legend({'root 1', 'root 2', 'root 3', 'root 4', 'root 5'})
This would produce 5 different lines, some of which might turn out to have no visible content (because the root had imaginary components).
It is entirely possible that at some point you might appear to get a discontinuity where two of the lines "switch places". That is because the ordering of the roots is not well defined, other than that the first root is going to be real if there are any real roots (and even that is not promised anywhere.)
Torsten
on 8 Feb 2018
For each value of dtc, you will have to figure out which of the 5 values supplied for Tg is the physical correct one.
Best wishes
Torsten.
More Answers (1)
See Also
Categories
Find more on Elements in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!