MATLAB Answers

Why does a for-loop give different outputs when used in different scripts?

Dursman Mchabe 님이 질문을 제출함. 25 Jan 2019
Hi everyone,
I am trying to use a for-loop to get time that is spaced by 600 seconds for the first 4 steps, then spaced by 7200 seconds for the remaining 26 steps to get:
Results =
600
1200
1800
2400
3000
10200
17400
24600
31800
39000
46200
53400
60600
67800
75000
82200
89400
96600
103800
111000
118200
125400
132600
139800
147000
154200
161400
168600
175800
183000
When I use the loop below, I get success:
%% Time
%t0 = 0.0;
t0 = [600;1200;1800;2400;3000;10200;17400;24600;31800;39000;46200;53400;60600;67800;75000;82200;89400;96600;103800;111000;118200;125400;132600;139800;147000;154200;161400;168600;175800;183000];
Delt = 1;
t = t0;
tmax = 6100;
nloop = tmax/Delt;
nstep = 30;
%% Results
% Results = [ t y0' pH pHInter];
Results = [ t ];
%% Loop1
for k = 1:nstep
%y = SO2_OdeDriver(y0,S);
%pH = 3 - log10 (CH) ;
tf = t;
%Results = [Results; tf y' pH pHInter];
Results = [tf];
%y0 = y;
end
Results
However, I need to use the loop inside a script, e.g. Line 91 through Line 115 in the attached script. And, when I do that I get an error message:
Error using horzcat
Dimensions of arrays being concatenated are not consistent.
Error in modeltimedifferent (line 101)
Results = [ t y0' pH pHInter];
In the attached code, I need Results to be spaced at 600 seconds in the first 4 steps, then 7200 in the remaining 26 steps.
How can I do that?
I also paste the code, for those who prefer it that way.
function Results = modeltimedifferent(y,S)
global S_total C_total Ca_total CSO2_bulkslurry CHSO3 CSO3 CCO2_bulkslurry CHCO3 CCO3 CCaCO3 CCaSO3 Trial CH ...
nloop Delt t CH_trial CSO2_inter CHSO3_inter CSO3_inter CCO2_inter CHCO3_inter CCO3_inter pHtrial pi pHInter
format shorte
%% Parameters
S.R = 8.314;
S.HSO2 = 149;
S.HCO2 = 5.15e+3;
S.DCa2 = 1.39e-9;
S.DSO2 = 2.89e-9;
S.DCO2 = 3.53e-9;
S.KCO2 = 1.7e-3;
S.KHCO3 = 6.55e-8;
S.KSO2 = 6.24;
S.KHSO3 = 5.68e-5;
S.Kw = 5.3E-8;
S.KSPCaSO3 = 3.1e-6;
S.BETCaSO3 = 10;
S.ktot = 8.825e-6;
S.BETCaCO3 = 12.54;
S.MWCaCO3 = 100.0869;
S.Kad = 1.03493;
S.rhoCaCO3 = 2703;
S.MWCaSO3 = 258.30;
S.rhoCaSO3 = 2540;
S.kCaS03 = 2e-05;
S.P = 1.01E+05;
S.T = 323.15;
S.V_Headspace = 1.95e-3;
S.F = 1.7e-5;
S.y_so2 = 2e-3;
S.CSO2_in = 7.28e-2;
S.y_co2 = 0.0;
S.CCO2_in = 0.0;
S.kga = 4.23e-6;
S.kLa = 5.78e-3;
S.kLa_CO2 = 9.598e-4;
%% Initial conditions and Relationships
pHtrial = 7;
Ca_total = 4.879E-02;
C_total = 4.879E-02;
S_total = 0.0;
CCaCO3 = 48.624- Ca_total;
CCaSO3 = 0;
fr = 1;
SO2_g = fr * S.CSO2_in;
CO2_g = S.CCO2_in;
CO2_g = 0.2;
CCO2_inter = CO2_g ;
pg = 200;
CSO2_inter = pg/S.HSO2;
CH_trial = 1.e-7;
pH1 = 7.0;
pH = fzero(@(pH)HionpH(pH,S),pH1);
CH = 10^(3-pH);
CSO2_bulkslurry = (S_total*CH^2)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3);
CHSO3 = (S_total*S.KSO2*CH)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3);
CSO3 = (S_total*S.KSO2*S.KHSO3)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3);
CCO2_bulkslurry = (C_total*CH^2)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3);
CHCO3 = (C_total*S.KCO2*CH)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3);
CCO3 = (C_total*S.KCO2*S.KHCO3)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3);
pH2 = 6.0;
pHInter = fzero (@(pH)HionStar(pH,S),pH2);
CH = 10^(3-pH) ;
y0 = [ SO2_g CO2_g S_total C_total Ca_total CCaCO3 CCaSO3 ];
y0 = y0';
Trial(1) = (SO2_g * S.R *S.T ) / S.HSO2;
Trial(2) = 1;
Trial(3) = CH ;
%% Time
t0 = [600;1200;1800;2400;3000;10200;17400;24600;31800;39000;46200;53400;60600;67800;75000;82200;89400;96600;103800;111000;118200;125400;132600;139800;147000;154200;161400;168600;175800;183000];
%t0 = 0.0;
Delt = 1.0;
t = t0;
tmax = 6100;
nloop = tmax/Delt ;
nstep = 30 ;
%% Results
Results = [ t y0' pH pHInter];
%% Loop1
for k = 1:nstep
y = SO2_OdeDriver(y0,S);
pH = 3 - log10 (CH) ;
modeltime = t;
Results = [Results; modeltime y' pH pHInter];
y0 = y;
end
Results
modeltime
end
function ph = HionpH (pH,S)
global S_total C_total Ca_total
CH = 10^(3-pH);
ph = CH + 2 * Ca_total - ((S_total*S.KSO2*CH)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3))...
- 2*((S_total*S.KSO2*S.KHSO3)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3))...
-((C_total*S.KCO2*CH)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3))...
- 2*((C_total*S.KCO2*S.KHCO3)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3))- S.Kw /CH ;
end
function y = SO2_OdeDriver (y0,S)
global S_total C_total Ca_total CSO2_bulkslurry CHSO3 CSO3 CCO2_bulkslurry CHCO3 CCO3 CCaCO3 CCaSO3 CH nloop Delt t CH_trial Tag pHInter
for k = 1:nloop
pH_trial = 8.0;
pH = fzero(@(pH)HionpH(pH,S),pH_trial);
CH = 10^(3-pH);
pH = 3 -log10 (CH) ;
CSO2_bulkslurry = (S_total*CH^2)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3);
CHSO3 = (S_total*S.KSO2*CH)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3);
CSO3 = (S_total*S.KSO2*S.KHSO3)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3);
CCO2_bulkslurry = (C_total*CH^2)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3);
CHCO3 = (C_total*S.KCO2*CH)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3);
CCO3 = (C_total*S.KCO2*S.KHCO3)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3);
dydt = odes (t, y0,S);
y = y0 + dydt' * Delt;
SO2_g = y(1);
CO2_g = y(2);
S_total = y(3) ;
C_total = y(4) ;
Ca_total = y(5);
CCaCO3 = y(6);
CCaSO3 = y(7);
y0 = y ;
t = t +Delt;
CH_trial = CH;
end
end
function dcdt = odes (t,c,S)
global kCaSO3 CSO2_bulkslurry CHSO3 CSO3 CCO2_bulkslurry CHCO3 CCO3 CCaCO3 CCaSO3 CH S_total C_total Ca_total Trial CSO2_inter ...
CHSO3_inter CSO3_inter CCO2_inter CHCO3_inter CCO3_inter pHtrial pi pHInter pg
S_total = c(3) ;
C_total = c(4) ;
Ca_total = c(5) ;
pH1 = 7;
pH = fzero (@(pH)HionpH(pH,S),pH1);
CH = 10^(3-pH) ;
CSO2_bulkslurry = (S_total*CH^2)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3);
CHSO3 = (S_total*S.KSO2*CH)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3);
CSO3 = (S_total*S.KSO2*S.KHSO3)/(CH^2 + S.KSO2*CH + S.KSO2*S.KHSO3);
CCO2_bulkslurry = (C_total*CH^2)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3);
CHCO3 = (C_total*S.KCO2*CH)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3);
CCO3 = (C_total*S.KCO2*S.KHCO3)/(CH^2 + S.KCO2*CH + S.KCO2*S.KHCO3);
CH;
pH;
p_exit = c(1)*S.R*S.T ;
pg = p_exit;
Trial ;
Tnew = fsolve (@(Tr)NSub(Tr,S),Trial);
pHInter = 3-log10(Tnew(3) );
pstar = Tnew(1) * S.HSO2;
NSO2 = S.kga *(pg-pstar) ;
%Trial = Tnew;
dcdt(1) = (1/S.V_Headspace) * (S.F * S.CSO2_in - S.F * c(1) ) - NSO2 ;
E_CO2 = 1;
NCO2 = S.kLa_CO2*E_CO2*( c(2)*S.R*S.T/S.HCO2 - CCO2_bulkslurry);
dcdt(2) = (1/S.V_Headspace) * (S.F * S.CCO2_in - S.F * c(2)) - NCO2 ;
n = 3;
sat = c(5)*CSO3/ S.KSPCaSO3 ;
if ( sat < 1.0 )
RCaSO3 = 0.0;
else
RCaSO3 = 1.0e-08 * ( sat -1)^n ;
end
RCaCO3 = S.ktot*S.BETCaCO3*S.MWCaCO3*CCaCO3*CH *(1 - (S.Kad*CH)/(1+S.Kad*CH));
dcdt(3) = NSO2 - RCaSO3;
dcdt(4) = NCO2 + RCaCO3;
dcdt(5) = RCaCO3 - RCaSO3;
if ( c(6) < 1.0e-3 )
c(6) = 0.0;
dcdt(6) = 0.0;
Time_c = t ;
Tag = 1;
else
dcdt(6) = - RCaCO3 ;
end
dcdt(7) = RCaSO3;
end
function discrep = NSub (Tr,S)
global S_total C_total Ca_total CSO2_bulkslurry CHSO3 CSO3 CSO2_inter CHSO3_inter CSO3_inter CCO2_inter CHCO3_inter CCO3_inter ...
CCO2_bulkslurry CHCO3 CCO3 pg pHInter
Tr;
CSO2_inter = Tr(1) ;
pstar = CSO2_inter * S.HSO2 ;
S1Total = Tr(2) ;
CH = Tr(3) ;
S_total = CHSO3 + CSO3 + CSO2_bulkslurry ;
CHSO3_inter = S.KSO2*CSO2_inter / CH ;
CSO3_inter = S.KHSO3 *CHSO3_inter /CH;
Ctotal = CHCO3 + CCO3 ;
CHCO3_inter = Ctotal * ( 1+ S.KHCO3 / CH) ;
CCO3_inter = S.KHCO3 *CHCO3_inter /CH ;
S_total_inter = CSO2_inter + CHSO3_inter + CSO3_inter;
df = S_total_inter- S_total ;
pstar;
pg;
Rate_gas = S.kga * (pg-pstar) ;
Rate_liq = S.kLa * df;
discrep(1) = S1Total - (CHSO3_inter + CSO3_inter) ;
discrep(2) = Rate_gas - Rate_liq ;
mCH_cal = 2*Ca_total -(CHSO3_inter + 2*CSO3_inter + CHCO3_inter + 2* CCO3_inter + S.Kw /CH ) ;
discrep(3) = CH + mCH_cal ;
end
function ph1 = HionStar (pH,S)
global S_total C_total Ca_total CSO2_bulkslurry CHSO3 CSO3 CSO2_inter CHSO3_inter CSO3_inter CCO2_inter CHCO3_inter ...
CCO3_inter CCO2_bulkslurry CHCO3 CCO3 pg
CH = 10^(3-pH);
CSO2_inter ;
CHSO3_inter = S.KSO2*CSO2_inter / CH ;
CSO3_inter = S.KHSO3 *CHSO3_inter /CH;
Stotal = CHSO3 + CSO3;
Ctotal = CHCO3 + CCO3 ;
CHSO3_inter = Stotal * ( 1+ S.KHSO3 / CH) ;
CSO3_inter = S.KHSO3 * CHSO3_inter /CH;
CHCO3_inter = Ctotal * ( 1+ S.KHCO3 / CH);
CCO3_inter = S.KHCO3 *CHCO3_inter /CH ;
Ca = Ca_total;
ph1 = CH + 2*Ca -(CHSO3_inter + 2*CSO3_inter + CHCO3_inter + 2* CCO3_inter + S.Kw /CH );
end

  댓글 수: 2

Your script loop that you say works overwrites the value of Results every time round the loop, leaving you with just the final value.
The error you point to though is just a simple concatenation error that the things you are concatenating are not of a size that allows them to be horizontally concatenated - e.g. you cannot concatenate a row vector with a column vector.
Just use the debgugger to check the size of each of the components [ t y0' pH pHInter] and you will easily see which is the wrong size. It may just need transposing or it may be completely incompatible.
Thanks alot for the comment. When I use one time step e.g. 6100 seconds, the code works well. I thought it was possible to use two time steps.

로그인 to comment.

태그

답변 수: 0


Translated by