I'm currently trying to create a function which depends on the compressibility factor of species in a mixture. This function uses subprogram of van der waals mixing rules.
The appears to an conistent error there which i can't fix.
Error Message is : The variable appears to change size on every loop iteration
The code used
Subprogram for Mixing rules
function [am, bm, aa] = MIXRULES_VDW(a, b, kappa, eta, x)%The funtion MIXRULES_VDW calculates the am and bm of a mixture using
%VDW MIXING RULES.
%All input/output data are expressed in SI or are dimensionless.
%giving a firt value for c
c=length(x);%initialization the values of am and bm
am=0;bm=0;%calculation of am, bm and aa
for i=1:c for j=1:c aa(i,j)=sqrt(a(i)*a(j))*(1-kappa(i,j)); bb(i,j)=(b(i)+b(j))/2*(1-eta(i,j)); am=am+x(i)*x(j)*aa(i,j); bm=bm+x(i)*x(j)*bb(i,j); endend
Actual program for Calculating compressibilty:
function [Z,V,ZL,ZV,VL,VV] = ZMIX_PRVDW(Pc, Tc, w, kappa, eta, P, T, x)%The funtion ZMIX_PRVDW calculates the compressibility factor z and the specific
%molar volume V for a mixture at given pressure P, temperature T and concentration
%using PR-CEoS and VDW MIXING RULES. This function requires REALROOTS_CARDANO and
% MIXRULES_VDW funtions to run.
%All input/output data are expressed in SI.
%P[Pa], T[K], w[dimensionless], V[m^3/mol], Z[dimensionless]
%universal gas constant
R=8.314;%PR attractive constant and covolume at critical point
b=0.07780*(R*Tc)./Pc;k=0.37464+1.54226*w-0.26992*w.^2;alpha=(1+k.*(1-sqrt(T./Tc))).^2;a=0.45724*alpha.*(R^2.*(Tc.^2))./Pc;%calling MIXRULES_VDW function
[am, bm, aa] = MIXRULES_VDW(a, b, kappa, eta, x);%calculation of Am and Bm
Am=am*P/(R*T)^2;Bm=bm*P/(R*T);%calculation of A and B
a2=-1+Bm;a1=Am-3*Bm^2-2*Bm;a0=-Am*Bm+Bm^2+Bm^3;%calculation of Z and V using REALROOTS_CARDANO function
Z=REALROOTS_CARDANO(a2,a1,a0);V=Z*R*T/P;%selection of Z and V values for liquid and vapor
if min(V)<b ZL=max(Z); ZV=max(Z); VL=max(V); VV=max(V);else ZL=min(Z); ZV=max(Z); VL=min(V); VV=max(V); endend
Both indicate same problem with the use of 'aa'. Same error is shown.
Any help would be appreciated.
Best Answer