📄 tesfunc.m
字号:
for I = 4:8
VPR=exp(AVP(I,1)+BVP(I,1)/(TCR+CVP(I,1)));
PPR(I,1)=VPR*XLR(I,1);
PTR=PTR+PPR(I,1);
VPR=exp(AVP(I,1)+BVP(I,1)/(TCS+CVP(I,1)));
PPS(I,1)=VPR*XLS(I,1);
PTS=PTS+PPS(I,1);
end
PTV=UTVV*RG*TKV/VTV;
for I=1:8
XVR(I,1)=PPR(I,1)/PTR;
XVS(I,1)=PPS(I,1)/PTS;
end
UTVR=PTR*VVR/RG/TKR;
UTVS=PTS*VVS/RG/TKS;
for I=4:8
UCVR(I,1)=UTVR*XVR(I,1);
UCVS(I,1)=UTVS*XVS(I,1);
end
% This is where the Matlab RR's show different values than the FORTRAN 77
% RR's. The code below has been modified with the chop command to
% elliminate the stray numbers that are created by the matlab code.
% Matlab: FORTRAN 77:
% R1F
% 1 1.0000000000000
% R2F
% 1 1.0000000000000
% TKR
% 3.935500000017000e+02 393.55000000169
% RR(1,1)
% 3.181229946057000e-09 3.1812215872809D-09
% RR(2,1)
% 1.569591352140000e-10 1.5695892900646D-10
% RR(3,1)
% 7.437635682791000e-11 7.4376432806282D-11
% RR(4,1)
% 5.708298619084000e-11 5.7083044503356D-11
%
% NOTE: I DID A COMPARISON OF THE EXP FUNCTIONS OF MATLAB
% AND FORTRAN 77 AND THERE IS A CONSISTENT DISCREPANCY AT
% 10^-16 BETWEEN THE OUTPUT OF THE TWO LANGUAGES
%
R1F = chop(R1F,13);
R2F = chop(R2F,13);
TKR = chop(TKR,13);
RR(1,1)=chop(exp(31.5859536-40000.0/1.987/TKR)*R1F,13);
RR(2,1)=chop(exp(3.00094014-20000.0/1.987/TKR)*R2F,13);
RR(3,1)=chop(exp(53.4060443-60000.0/1.987/TKR),13);
RR(4,1)=chop(RR(3,1)*0.767488334,13);
% This next section has the following output:
% Matlab: FORTRAN 77:
% PPR(1,1)
% 5.717769731191162e+03 5717.7697311912 * good
% PPR(3,1)
% 4.159943934675613e+03 4159.9439346756 * good
% R1F
% 2.174415112221978e+04 21744.151122220 * good
% R2F
% 2.247652526817436e+01 22.476525268174 * good
% RR(1,1)
% 3.517133031511259e-01 0.35171237901398
% RR(2,1)
% 2.861611337089543e-01 0.28616075776025
if PPR(1,1) > 0.0 & PPR(3,1) > 0.0
R1F=PPR(1,1)^1.1544;
R2F=PPR(3,1)^0.3735;
RR(1,1)=RR(1,1)*R1F*R2F*PPR(4,1);
RR(2,1)=RR(2,1)*R1F*R2F*PPR(5,1);
else
RR(1,1)=0.0;
RR(2,1)=0.0;
end
%
% This next section has the following output:
% Matlab: FORTRAN 77:
% RR(1,1)
% 2.516069042865577e+02 251.60624318195
% RR(2,1)
% 2.047125210635105e+02 204.71225211921
% RR(3,1)
% 1.134874248863408e+00 1.1348754081821
% RR(4,1)
% 5.281888235868223e-02 5.2818936315249D-02
%
RR(3,1)=RR(3,1)*PPR(1,1)*PPR(5,1);
RR(4,1)=RR(4,1)*PPR(1,1)*PPR(4,1);
for I=1:4
RR(I,1)=RR(I,1)*VVR;
end
%
% This next section has the following output:
% Matlab: FORTRAN 77:
% CRXR(1,1)
% -4.574542995989316e+02 -457.45337070934
% CRXR(3,1)
% -4.563194253500682e+02 -456.31849530116
% CRXR(4,1)
% -2.516861326100957e+02 -251.68547158642
% CRXR(5,1)
% -2.058473953123739e+02 -205.84712752739
% CRXR(6,1)
% 1.187693131222090e+00 1.1876943444973
% CRXR(7,1)
% 2.516069042865577e+02 251.60624318195
% CRXR(8,1)
% 2.047125210635105e+02 204.71225211921
% RH
% 2.759494513807820e+01 27.594886078737
CRXR(1,1)=-RR(1,1)-RR(2,1)-RR(3,1);
CRXR(3,1)=-RR(1,1)-RR(2,1);
CRXR(4,1)=-RR(1,1)-1.5*RR(4,1);
CRXR(5,1)=-RR(2,1)-RR(3,1);
CRXR(6,1)=RR(3,1)+RR(4,1);
CRXR(7,1)=RR(1,1);
CRXR(8,1)=RR(2,1);
RH=RR(1,1)*HTR(1,1)+RR(2,1)*HTR(2,1);
XMWS(1,1)=0.0;
XMWS(2,1)=0.0;
XMWS(6,1)=0.0;
XMWS(8,1)=0.0;
XMWS(9,1)=0.0;
XMWS(10,1)=0.0;
%
% XST & XMWS seem to be calculated correcly in the next
% section with little error.
%
for I=1:8
XST(I,6)=XVV(I,1);
XST(I,8)=XVR(I,1);
XST(I,9)=XVS(I,1);
XST(I,10)=XVS(I,1);
XST(I,11)=XLS(I,1);
XST(I,13)=XLC(I,1);
XMWS(1,1)=XMWS(1,1)+XST(I,1)*XMW(I,1);
XMWS(2,1)=XMWS(2,1)+XST(I,2)*XMW(I,1);
XMWS(6,1)=XMWS(6,1)+XST(I,6)*XMW(I,1);
XMWS(8,1)=XMWS(8,1)+XST(I,8)*XMW(I,1);
XMWS(9,1)=XMWS(9,1)+XST(I,9)*XMW(I,1);
XMWS(10,1)=XMWS(10,1)+XST(I,10)*XMW(I,1);
end
TST(6,1)=TCV;
TST(8,1)=TCR;
TST(9,1)=TCS;
TST(10,1)=TCS;
TST(11,1)=TCS;
TST(13,1)=TCC;
% The HST's from the following substitutions for TESUB1 are
% computed correctly.
%
% TESUB1(XST(1,1),TST(1,1),HST(1,1),1);
% TESUB1(Z,T,H,ITY)
% Substitution of TESUB1 for call statement. MWB
Z(1:8,1)=XST(1:8,1);
T=TST(1,1);
H=HST(1,1);
ITY=1;
if ITY == 0
H = 0;
for I = 1:8
HI = T * (AH(I,1) + BH(I,1) * T / 2 + CH(I,1) * T^2 / 3);
HI = 1.8 * HI;
H = H + Z(I,1) * XMW(I,1) * HI;
end
else
H = 0;
for I = 1:8
HI = T * (AG(I,1) + BG(I,1) * T / 2 + CG(I,1) * T^2 / 3);
HI = 1.8 * HI;
HI = HI + AV(I,1);
H = H + Z(I,1) * XMW(I,1) * HI;
end
end
if ITY == 2
R1 = 3.57696 / 1000000;
H = H - R1 * (T + 273.15);
end
HST(1,1) = H;
% end of TESUB1
% TESUB1(XST(1,2),TST(2,1),HST(2,1),1);
% TESUB1(Z,T,H,ITY)
% Substitution of TESUB1 for call statement. MWB
Z(1:8,1)=XST(1:8,2);
T=TST(2,1);
H=HST(2,1);
ITY=1;
if ITY == 0
H = 0;
for I = 1:8
HI = T * (AH(I,1) + BH(I,1) * T / 2 + CH(I,1) * T^2 / 3);
HI = 1.8 * HI;
H = H + Z(I,1) * XMW(I,1) * HI;
end
else
H = 0;
for I = 1:8
HI = T * (AG(I,1) + BG(I,1) * T / 2 + CG(I,1) * T^2 / 3);
HI = 1.8 * HI;
HI = HI + AV(I,1);
H = H + Z(I,1) * XMW(I,1) * HI;
end
end
if ITY == 2
R1 = 3.57696 / 1000000;
H = H - R1 * (T + 273.15);
end
HST(2,1) = H;
% end of TESUB1
% TESUB1(XST(1,3),TST(3,1),HST(3,1),1);
% TESUB1(Z,T,H,ITY)
% Substitution of TESUB1 for call statement. MWB
Z(1:8,1)=XST(1:8,3);
T=TST(3,1);
H=HST(3,1);
ITY=1;
if ITY == 0
H = 0;
for I = 1:8
HI = T * (AH(I,1) + BH(I,1) * T / 2 + CH(I,1) * T^2 / 3);
HI = 1.8 * HI;
H = H + Z(I,1) * XMW(I,1) * HI;
end
else
H = 0;
for I = 1:8
HI = T * (AG(I,1) + BG(I,1) * T / 2 + CG(I,1) * T^2 / 3);
HI = 1.8 * HI;
HI = HI + AV(I,1);
H = H + Z(I,1) * XMW(I,1) * HI;
end
end
if ITY == 2
R1 = 3.57696 / 1000000;
H = H - R1 * (T + 273.15);
end
HST(3,1) = H;
% end of TESUB1
% TESUB1(XST(1,4),TST(4,1),HST(4,1),1);
% TESUB1(Z,T,H,ITY)
% Substitution of TESUB1 for call statement. MWB
Z(1:8,1)=XST(1:8,4);
T=TST(4,1);
H=HST(4,1);
ITY=1;
if ITY == 0
H = 0;
for I = 1:8
HI = T * (AH(I,1) + BH(I,1) * T / 2 + CH(I,1) * T^2 / 3);
HI = 1.8 * HI;
H = H + Z(I,1) * XMW(I,1) * HI;
end
else
H = 0;
for I = 1:8
HI = T * (AG(I,1) + BG(I,1) * T / 2 + CG(I,1) * T^2 / 3);
HI = 1.8 * HI;
HI = HI + AV(I,1);
H = H + Z(I,1) * XMW(I,1) * HI;
end
end
if ITY == 2
R1 = 3.57696 / 1000000;
H = H - R1 * (T + 273.15);
end
HST(4,1) = H;
% end of TESUB1
% TESUB1(XST(1,6),TST(6,1),HST(6,1),1);
% TESUB1(Z,T,H,ITY)
% Substitution of TESUB1 for call statement. MWB
Z(1:8,1)=XST(1:8,6);
T=TST(6,1);
H=HST(6,1);
ITY=1;
if ITY == 0
H = 0;
for I = 1:8
HI = T * (AH(I,1) + BH(I,1) * T / 2 + CH(I,1) * T^2 / 3);
HI = 1.8 * HI;
H = H + Z(I,1) * XMW(I,1) * HI;
end
else
H = 0;
for I = 1:8
HI = T * (AG(I,1) + BG(I,1) * T / 2 + CG(I,1) * T^2 / 3);
HI = 1.8 * HI;
HI = HI + AV(I,1);
H = H + Z(I,1) * XMW(I,1) * HI;
end
end
if ITY == 2
R1 = 3.57696 / 1000000;
H = H - R1 * (T + 273.15);
end
HST(6,1) = H;
% end of TESUB1
% TESUB1(XST(1,8),TST(8,1),HST(8,1),1);
% TESUB1(Z,T,H,ITY)
% Substitution of TESUB1 for call statement. MWB
Z(1:8,1)=XST(1:8,8);
T=TST(8,1);
H=HST(8,1);
ITY=1;
if ITY == 0
H = 0;
for I = 1:8
HI = T * (AH(I,1) + BH(I,1) * T / 2 + CH(I,1) * T^2 / 3);
HI = 1.8 * HI;
H = H + Z(I,1) * XMW(I,1) * HI;
end
else
H = 0;
for I = 1:8
HI = T * (AG(I,1) + BG(I,1) * T / 2 + CG(I,1) * T^2 / 3);
HI = 1.8 * HI;
HI = HI + AV(I,1);
H = H + Z(I,1) * XMW(I,1) * HI;
end
end
if ITY == 2
R1 = 3.57696 / 1000000;
H = H - R1 * (T + 273.15);
end
HST(8,1) = H;
% end of TESUB1
% TESUB1(XST(1,9),TST(9,1),HST(9,1),1);
% TESUB1(Z,T,H,ITY)
% Substitution of TESUB1 for call statement. MWB
Z(1:8,1)=XST(1:8,9);
T=TST(9,1);
H=HST(9,1);
ITY=1;
if ITY == 0
H = 0;
for I = 1:8
HI = T * (AH(I,1) + BH(I,1) * T / 2 + CH(I,1) * T^2 / 3);
HI = 1.8 * HI;
H = H + Z(I,1) * XMW(I,1) * HI;
end
else
H = 0;
for I = 1:8
HI = T * (AG(I,1) + BG(I,1) * T / 2 + CG(I,1) * T^2 / 3);
HI = 1.8 * HI;
HI = HI + AV(I,1);
H = H + Z(I,1) * XMW(I,1) * HI;
end
end
if ITY == 2
R1 = 3.57696 / 1000000;
H = H - R1 * (T + 273.15);
end
HST(9,1) = H;
% end of TESUB1
HST(10,1)=HST(9,1);
% TESUB1(XST(1,11),TST(11,1),HST(11,1),0);
% TESUB1(Z,T,H,ITY)
% Substitution of TESUB1 for call statement. MWB
Z(1:8,1)=XST(1:8,11);
T=TST(11,1);
H=HST(11,1);
ITY=0;
if ITY == 0
H = 0;
for I = 1:8
HI = T * (AH(I,1) + BH(I,1) * T / 2 + CH(I,1) * T^2 / 3);
HI = 1.8 * HI;
H = H + Z(I,1) * XMW(I,1) * HI;
end
else
H = 0;
for I = 1:8
HI = T * (AG(I,1) + BG(I,1) * T / 2 + CG(I,1) * T^2 / 3);
HI = 1.8 * HI;
HI = HI + AV(I,1);
H = H + Z(I,1) * XMW(I,1) * HI;
end
end
if ITY == 2
R1 = 3.57696 / 1000000;
H = H - R1 * (T + 273.15);
end
HST(11,1) = H;
% end of TESUB1
% TESUB1(XST(1,13),TST(13,1),HST(13,1),0);
% TESUB1(Z,T,H,ITY)
% Substitution of TESUB1 for call statement. MWB
Z(1:8,1)=XST(1:8,13);
T=TST(13,1);
H=HST(13,1);
ITY=0;
if ITY == 0
H = 0;
for I = 1:8
HI = T * (AH(I,1) + BH(I,1) * T / 2 + CH(I,1) * T^2 / 3);
HI = 1.8 * HI;
H = H + Z(I,1) * XMW(I,1) * HI;
end
else
H = 0;
for I = 1:8
HI = T * (AG(I,1) + BG(I,1) * T / 2 + CG(I,1) * T^2 / 3);
HI = 1.8 * HI;
HI = HI + AV(I,1);
H = H + Z(I,1) * XMW(I,1) * HI;
end
end
if ITY == 2
R1 = 3.57696 / 1000000;
H = H - R1 * (T + 273.15);
end
HST(13,1) = H;
% end of TESUB1
% This next section is good.
FTM(1,1)=VPOS(1,1)*VRNG(1,1)/100.0;
FTM(2,1)=VPOS(2,1)*VRNG(2,1)/100.0;
FTM(3,1)=VPOS(3,1)*(1-IDV(6,1))*VRNG(3,1)/100.0;
FTM(4,1)=VPOS(4,1)*(1-IDV(7,1)*0.2)*VRNG(4,1)/100.0+1.E-10;
FTM(11,1)=VPOS(7,1)*VRNG(7,1)/100.0;
FTM(13,1)=VPOS(8,1)*VRNG(8,1)/100.0;
UAC=VPOS(9,1)*VRNG(9,1)*(1+TESUB8(9,TIME))/100.0;
FWR=VPOS(10,1)*VRNG(10,1)/100.0;
FWS=VPOS(11,1)*VRNG(11,1)/100.0;
AGSP=(VPOS(12,1)+150.0)/100.0;
DLP=PTV-PTR;
if DLP<0
DLP=0.0;
end
FLMS=1937.6*sqrt(DLP);
FTM(6,1)=FLMS/XMWS(6,1);
DLP=PTR-PTS;
if DLP<0
DLP=0.0;
end
FLMS=4574.21*sqrt(DLP)*(1-0.25*TESUB8(12,TIME));
FTM(8,1)=FLMS/XMWS(8,1);
DLP=PTS-760.0;
if DLP < 0
DLP=0.0;
end
FLMS=VPOS(6,1)*0.151169*sqrt(DLP);
FTM(10,1)=FLMS/XMWS(10,1);
PR=PTV/PTS;
if PR < 1
PR=1.0;
end
if PR > CPPRMX
PR=CPPRMX;
end
FLCOEF=CPFLMX/1.197;
FLMS=CPFLMX+FLCOEF*(1.0-PR^3);
CPDH=FLMS*(TCS+273.15)*1.8E-6*1.9872*(PTV-PTS)/(XMWS(9,1)*PTS);
DLP=PTV-PTS;
if DLP < 0.0
DLP=0.0;
end
FLMS=FLMS-VPOS(5,1)*53.349*sqrt(DLP);
if FLMS < 1.E-3
FLMS=1.E-3;
end
FTM(9,1)=FLMS/XMWS(9,1);
HST(9,1)=HST(9,1)+CPDH/FTM(9,1);
for I=1:8
FCM(I,1)=XST(I,1)*FTM(1,1);
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -