No error but I am getting Complex double!
12 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
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 Kommentare
Antworten (2)
Mark Sherstan
am 1 Mär. 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
am 1 Mär. 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 Kommentare
Siehe auch
Kategorien
Mehr zu Logical finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!