No error but I am getting Complex double!
12 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
clear ALL
close ALL
clc
Area=0.1963;
L=4;
porosity=0.18;
K=61.5;
injrate=0.01;
oilvisc=4.2;
watervisc=0.85;
Boi=1.21;
Swr=0.28;
Sorw=0.35;
KrwMax=0.31;
KroMax=0.65;
Nw=2.25;
No=2.45;
CO=0.000001;
CW=0.000004;
CR=0.000003;
POi=1500;
PWF=1460;
dt=0.2;
Vr=Area*0.4;
for i=1:10
Swi(1,i)=Swr;
PO(1,i)=POi;
PCOW(1,i)=5;
Krw(1,i)=KrwMax*((Swi(1,i)-Swr)/(1-Sorw-Swr))^Nw;
Kro(1,i)=KroMax*((1-Swi(1,i)-Sorw)/(1-Sorw-Swr))^No;
Omob(1,i)= Kro(1,i)/oilvisc;
Wmob(1,i)=Krw(1,i)/watervisc;
end
for n=1:15
for i=1:10
Ctot(n,:)=CR+Swi(n,:)*CW+(1-Swi(n,:))*CO;
if i==1
Ttot=(2.6366e-4)*K*(Wmob(n,i)+Omob(n,i))*(Area/0.4);
Tw=(2.6366e-4)*K*(Wmob(n,i))*(Area/0.4);
A(i,i)=-(Ttot+(Vr*porosity*Ctot(n,i))/(dt)); %E
A(i,i+1)=Ttot; %F
R(i,1)= (-(Vr*porosity*Ctot(n,i)*PO(n,i))/(dt))+Tw*(PCOW(n,i+1)-PCOW(n,i))-injrate; %R
elseif i==10
WI(n,i)=(2.6366e-4)*K*(Wmob(n,i)+Omob(n,i))*(Area/(0.4/2));
TtotL=Ttot;
TwL=Tw;
Ttot=0;
Tw=0;
%Ttot=2.6366*10^-4*K*(Wmob(n,i)+Omob(n,i))*(A/0.4);
%Tw=2.6366*10^-4*K*(Wmob(n,i))*(A/0.4);
A(i,i-1)=TtotL; %D
A(i,i)=-(TtotL+((Vr*porosity*Ctot(n,i))/(dt))+WI(n,i)); %E
R(i,1)=(-(Vr*porosity*Ctot(n,i)*PO(n,i))/(dt))-TwL*(PCOW(n,i)-PCOW(n,i-1))-WI(n,i)*PWF; %R
else
TtotL=Ttot;
TwL=Tw;
Ttot=(2.6366e-4)*K*(Wmob(n,i)+Omob(n,i))*(Area/0.4);
Tw=(2.6366e-4)*K*(Wmob(n,i))*(Area/0.4);
A(i,i)=-(TtotL+Ttot+(Vr*porosity*Ctot(n,i))/dt); %E
A(i,i+1)=Ttot; %F
A(i,i-1)=TtotL; %D
R(i,1)= (-(Vr*porosity*Ctot(n,i)*PO(n,i))/(dt))+Tw*(PCOW(n,i+1)-PCOW(n,i))-TwL*(PCOW(n,i)-PCOW(n,i-1)); %R
end
end
PO(n+1,:)=A\R;
for i=1:10
if i==1
q_t(n,i)=injrate;
Tw=(2.6366e-4)*K*(Wmob(n,i))*(Area/0.4);
Swi(n+1,i)= (1+(CR+CW)*(PO(n+1,i)-PO(n,i)))*Swi(n,i)+...
(dt/(Vr*porosity))*(Tw*(PO(n+1,i+1)-PO(n+1,i))-...
Tw*(PCOW(n,i+1)-PCOW(n,i))+q_t(n,i));
elseif i==10
TwL=Tw;
Tw=0;
WI(n,i)=(2.6366e-4)*K*(Wmob(n,i)+Omob(n,i))*(Area/(0.4/2));
qt(n,i)= -WI(n,i)*(PO(n+1,i)-PWF);
Swi(n+1,i)=(1+(CR+CW)*(PO(n+1,i)-PO(n,i)))*Swi(n,i)+dt/(Vr*porosity)*...
(-TwL*(PO(n+1,i)-PO(n+1,i-1))+TwL*(PCOW(n,i)-PCOW(n,i-1))+qt(n,i));
else
TwL=Tw;
Tw=2.6366e-4*K*(Wmob(n,i))*(Area/0.4);
Swi(n+1,i)=(1+(CR+CW)*(PO(n+1,i)-PO(n,i)))*Swi(n,i)+ ...
(0.2/(Vr*porosity))*Tw*(PO(n+1,i+1)-PO(n+1,i))-...
TwL*(PO(n+1,i)-PO(n+1,i-1))-Tw*(PCOW(n,i+1)-PCOW(n,i))+TwL*(PCOW(n,i)-PCOW(n,i-1));
end
if Swi(n+1,i)<Swr
Swi(n+1,i)=Swr;
end
end
for i=1:10
if Swi(n+1,i)<=Swr
PCOW(n+1,i)=5;
elseif Swi(n+1,i)>=(1-Sorw)
PCOW(n+1,i)=-5;
elseif Swi(n+1,i)>Swr && Swi(n+1,i)<0.5
PCOW(n+1,i)= min(-0.8*log((Swi(n+1,i)-Swr)/(0.5-Swr)),5);
elseif Swi(n+1,i)>0.5 && Swi(n+1)<(1-Sorw)
PCOW(n+1,i)=max(0.8*log((1-Swi(n+1,i)-Sorw)/(1-0.5-Sorw)),-5);
end
end
for i=1:10
Krw(n+1,i)=KrwMax*((Swi(n+1,i)-Swr)/(1-Sorw-Swr))^Nw;
Kro(n+1,i)=KroMax*((1-Swi(n+1,i)-Sorw)/(1-Sorw-Swr))^No;
Omob(n+1,i)= Krw(n+1,i)/watervisc;
Wmob(n+1,i)=Kro(n+1,i)/oilvisc;
end
end
0 commentaires
Réponses (2)
Mark Sherstan
le 1 Mar 2019
You should not use i or j as the names of loop variables as these are both names for the inbuilt imaginary unit. Try changing that and see if the problem persisits.
Walter Roberson
le 1 Mar 2019
When n = 4, i = 1, then 1-Swi(n+1),i)-Sorw becomes negative, so (1-Swi(n+1,i)-Sorw)/(1-Sorw-Swr) becomes negative. You are raising a negative value to a non-integer. A^B with scalar A, B is defined as exp(B*LOG(A)) so when A<0 the log is complex and unless B is an integer you get a complex result from the ^ operation. Everything else gets tainted from there.
0 commentaires
Voir également
Catégories
En savoir plus sur Logical dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!