【物理应用】基于Matlab模拟井筒多相流附matlab代码

简介: 【物理应用】基于Matlab模拟井筒多相流附matlab代码

 ✅作者简介:热爱科研的Matlab仿真开发者,修心和技术同步精进,matlab项目合作可私信。

🍎个人主页:Matlab科研工作室

🍊个人信条:格物致知。

更多Matlab仿真内容点击👇

智能优化算法神经网络预测雷达通信 无线传感器

信号处理图像处理路径规划元胞自动机无人机

1 内容介绍

采油工程专业的培养目标,是能从事油气田开发方案设计,采油方法优选,油气井增产方法和提高采收率方法研究与应用等方面工作的高级工程技术人才.一名合格的工程技术人员,不仅要掌握扎实的基础理论知识,更要具备一定的生产技术技能和综合研究分析解决实际生产问题的能力

2 部分代码

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%                  Analytic model of SAGD                 %

%% 初始化Matlab,工作空间清理

clc;

clear;

tic;

%% 数据载入

SIP=9.869;                 % 读入参数   % SIP:大气压atm与MPa之间的换算

[Info]=Info(1);            % 水平井水平段信息                            

[WellSlotScr]=WellSlot(1);                

%% 参数读取

[DX,DY,DZ]=deal(Info(1,1),Info(1,2),Info(1,3));    % (1)网格数据读入

[P,T,SO,SW,KX,KY,KZ,VP]=deal(Info(2,1),Info(2,2),Info(2,3),Info(2,4),Info(2,5),Info(2,6),Info(2,7),Info(2,8));      %(2)油藏初始数据读入

[PWFC,TWFC,STQ,S,L,Vsx,PWFCS,TWFCS,STQS]=deal(Info(3,1),Info(3,2),Info(3,3),Info(3,4),Info(3,5),Info(3,6),Info(3,7),Info(3,8),Info(3,9));   %(3)井控数据读入

[D,rlo,rh,hf,ha,lamdace,lamdae,arfa,ls,ws,lu,ts,reh,dert,ms,ns,ngf]=deal(WellSlotScr(1,1),WellSlotScr(1,2),...      %(4)筛管数据读入

   WellSlotScr(1,3),WellSlotScr(1,4),WellSlotScr(1,5),WellSlotScr(1,6),WellSlotScr(1,7) ,WellSlotScr(1,8),...

   WellSlotScr(1,9), WellSlotScr(1,10),WellSlotScr(1,11),WellSlotScr(1,12),WellSlotScr(1,13),WellSlotScr(1,14),...

   WellSlotScr(1,15),WellSlotScr(1,16),WellSlotScr(1,17));

%% 统一压力单位

P=P*SIP;           %※压力相关的要单位转为 atm

PWFC=PWFC*SIP;     %※压力相关的要单位转为 atm

PWFCS=PWFCS*SIP;   %※压力相关的要单位转为 atm

SG=1-SO-SW;

%%

NUM=L/DX;   % 水平段微元段数目

DELT=1;     % 时间步长,day

FT=0;       % 已模拟总时间,day

FTMAX=1.5;   % 模拟总时间,day

%%

TOC=zeros(1,1);         % 总产油量

TWC=zeros(1,1);         % 总注/产水量

Matinj1=[ ];

Matinj2=[ ];

Matinj3=[ ];

%% 大循环

NNMAX=8000000;

Iteration=0;              % 迭代次数记录

for N=1:NNMAX

   FT=FT+DELT;

   if FT<=FTMAX

       PWF(1:NUM)=PWFC;

       TWF(1:NUM)=TWFC;

       SQWF(1:NUM)=STQ;

   elseif  FT>FTMAX

       break

   end

   

   [DDP,DDSQ]=deal(1,0);

   while DDP>=0.01 || DDSQ>0.01    % while判断循环:井筒压力、干度

       Iteration=Iteration+1;      % 迭代次数记录

       

       [QO,QW,QG,QHG]=QRATEsteam(NUM,P,T,SW,SG,S,PWF,TWF,SQWF,DX,DY,DZ,KY,KZ,rlo);          %※  注汽解析模型

%%  注汽解析模型部分

       

    [PWFN,TWFN,SQWFN,TubPWFN,TubTWFN,TubSQWFN,TubLOSS,AunLOSS,QLOSS]=Sagd4LiHeel(NUM,DX,Vsx,QG,T,PWFC,TWFC,STQ,DELT,D,ls,ws,lu,dert,ngf);     % 可运算

       

%% 注汽判断

           DDP=max(abs((PWF(1:NUM)-PWFN(1:NUM))),[],2);    % 求出每行的最大值

           DDT=max(abs((TWF(1:NUM)-TWFN(1:NUM))),[],2);    % PWFN和TWFN为wellbore中计算最新值;应该是计算值与假定值的对比

           DDSQ=max(abs((SQWF(1:NUM)-SQWFN(1:NUM))),[],2);

           

           PWF(1:NUM)=PWFN(1:NUM);

           TWF(1:NUM)=TWFN(1:NUM);

           SQWF(1:NUM)=SQWFN(1:NUM);

           

           if DDP<10 && DDSQ<10          %  该处判断其实可以去掉。误差值0.0001是否太小,调大点可以减少循环次数

               break;

           end

%% 迭代次数达到一定值跳出循环,注意迭代次数越大越接近精确值

       Iteration

       if Iteration>=1000

           break

       end

%%

   end  % 对应while DDP>=0.0001 || DDT>=0.0001 || DDSQ>0.0001

%% 井筒压力单位处理

  TubPWFN(:)=TubPWFN(:)/SIP;       % 输出单位为MPa

  PWFN(:)=PWFN(:)/SIP;             % 输出单位为MPa

 

%% 注汽解析模型数据输出

           TWR=0;

           for M=1:NUM

               PP=P;

               TT=T;

               TIN=TWFN(M);

               TIT=(TT+TIN)/2.0;   % ???是否加权平均

               

               [RG,~,HG]=STEAM(TWFN(M),SQWFN(M));

               [TubRG,~,TubHG]=STEAM(TubTWFN(M),TubSQWFN(M));

               

               HG=HG/1000;            % 焓值单位 10^6 J/kg

               TubHG=TubHG/1000;      % 焓值单位 10^6 J/kg

           

               QGG=QG(M)/1000;                   % QG单位为Kg/D,该处除以水密度则为 m3/D,按地面条件下水的体积来算

               QtonD=QG(M)/1000;                 % 注汽量输出单位为 吨/天

               QKGS=QG(M)/86400;                 % 注汽量输出单位为 kg/s

               

               Entha=QHG(M)/DX/86400;            % QHG单位是 kJ/m.s  

               Tubloss=TubLOSS(M)/DX/86400;      % 单位是 kJ/m.s

               Aunloss=AunLOSS(M)/DX/86400;      % 单位是 kJ/m.s

               Qloss=QLOSS(M)/DX/86400;          % 单位是 kJ/m.s  

               

               HeatTA=Entha+Tubloss+Aunloss;     % 单位是 kJ/m.s ;注入热量加长油管/短油管到地层以及环空到地层的热量;根据井筒中传热方式计算选择地层到底吸收了多少热量

               HeatA=Entha+Aunloss;              % 单位是 kJ/m.s ;注入热量加环空到地层的热量;根据井筒中传热方式计算选择地层到底吸收了多少热量

               

               Matinj1=[Matinj1;FT,M,M*DX,TubPWFN(M),PWFN(M),TubTWFN(M),TWFN(M),TubSQWFN(M),SQWFN(M),QGG,QtonD,QKGS];    % Matinj1指输出的第1个注入井数据矩阵

               Matinj2=[Matinj2;FT,M,M*DX,TubHG,HG,Tubloss,Aunloss,Qloss,Entha,HeatTA,HeatA];    % 之前输出顺序为[TubHG,HG,Tubloss,Aunloss,HeatTA,HeatA,Entha,Qloss]

               TWR=TWR+QtonD;           % 单位为 吨/天

               TWC=TWC+QtonD*DELT;      % 单位为 吨

           end

           Matinj3=[Matinj3;FT,DELT,TWR,TWC];  % Matpro2指输出的第2个注入井数据矩阵

%%

   N

   toc

end   % 对应 for N=1:NNMAX

%% 出图部分

% hold on

% figure;plot(Matinj1(:,3),Matinj1(:,[4,5]),'r+')  figure;plot(Matinj1(:,3),Matinj1(:,4),'r+')

% figure;plot(Matinj1(:,3),Matinj1(:,4),'LineStyle', 'none', 'Marker', '^')

figure;plot(Matinj1(:,3),Matinj1(:,4),'r^')

hold on

plot(Matinj1(:,3),Matinj1(:,5),'bs')

title('井筒蒸汽压力分布图','Fontsize',18)

xlabel('距水平井跟端距离/m','Fontsize',17)

ylabel('井筒蒸汽压力/MPa','Fontsize',17)

h=legend('长管','环空');

set(h,'FontSize',13);     %

set(gca,'FontSize',11);   % 坐标字体大小设置

%

figure;plot(Matinj1(:,3),Matinj1(:,6),'r^')

hold on

plot(Matinj1(:,3),Matinj1(:,7),'bs')

title('井筒蒸汽温度分布图','Fontsize',18)

xlabel('距水平井跟端距离/m','Fontsize',17)

ylabel('井筒蒸汽温度/°C','Fontsize',17)

h=legend('长管','环空');

set(h,'FontSize',13);     %

set(gca,'FontSize',11);   % 坐标字体大小设置

%

figure;plot(Matinj1(:,3),Matinj1(:,8),'r^')

hold on

plot(Matinj1(:,3),Matinj1(:,9),'bs')

title('井筒蒸汽干度分布图','Fontsize',18)

xlabel('距水平井跟端距离/m','Fontsize',17)

ylabel('井筒蒸汽干度/°C','Fontsize',17)

h=legend('长管','环空');

set(h,'FontSize',13);     %

set(gca,'FontSize',11);   % 坐标字体大小设置

%% 数据输出

% Nameinj1={'FT','M','TubPWFN(M)','PWF(M)','TubTWFN(M)','TWF(M)','TubSQWFN(M)','SQWF(M)','QGG','QKGS','RO','RW'};    % 注入井数据输出

% xlswrite('Parameters',Nameinj1,'sheet1','A1');

% xlswrite('Parameters',Matinj1,'sheet1','A2');

%

% Nameinj2={'FT','DELT','日注水TWR','总注水TWC(J)'};

% xlswrite('Parameters',Nameinj2,'sheet1','R1');

% xlswrite('Parameters',Matinj2,'sheet1','R2');

toc

3 运行结果

image.gif编辑

image.gif编辑

4 参考文献

[1]代锋, 孙凯, 厉爽,等. 欠平衡钻井井筒多相流技术研究[J]. 石油矿场机械, 2009, 38(1):4.

博主简介:擅长智能优化算法神经网络预测信号处理元胞自动机图像处理路径规划无人机雷达通信无线传感器等多种领域的Matlab仿真,相关matlab代码问题可私信交流。

部分理论引用网络文献,若有侵权联系博主删除。

相关文章
|
25天前
|
机器学习/深度学习 传感器 算法
【无人车路径跟踪】基于神经网络的数据驱动迭代学习控制(ILC)算法,用于具有未知模型和重复任务的非线性单输入单输出(SISO)离散时间系统的无人车的路径跟踪(Matlab代码实现)
【无人车路径跟踪】基于神经网络的数据驱动迭代学习控制(ILC)算法,用于具有未知模型和重复任务的非线性单输入单输出(SISO)离散时间系统的无人车的路径跟踪(Matlab代码实现)
|
25天前
|
机器学习/深度学习 算法 安全
【图像处理】使用四树分割和直方图移动的可逆图像数据隐藏(Matlab代码实现)
【图像处理】使用四树分割和直方图移动的可逆图像数据隐藏(Matlab代码实现)
108 2
|
25天前
|
canal 算法 vr&ar
【图像处理】基于电磁学优化算法的多阈值分割算法研究(Matlab代码实现)
【图像处理】基于电磁学优化算法的多阈值分割算法研究(Matlab代码实现)
|
25天前
|
机器学习/深度学习 存储 算法
【微电网调度】考虑需求响应的基于改进多目标灰狼算法的微电网优化调度研究(Matlab代码实现)
【微电网调度】考虑需求响应的基于改进多目标灰狼算法的微电网优化调度研究(Matlab代码实现)
|
25天前
|
机器学习/深度学习 存储 算法
【水下机器人建模】基于QLearning自适应强化学习PID控制器在AUV中的应用研究(Matlab代码实现)
【水下机器人建模】基于QLearning自适应强化学习PID控制器在AUV中的应用研究(Matlab代码实现)
206 0
|
25天前
|
传感器 资源调度 算法
【数据融合】【状态估计】基于KF、UKF、EKF、PF、FKF、DKF卡尔曼滤波KF、无迹卡尔曼滤波UKF、拓展卡尔曼滤波数据融合研究(Matlab代码实现)
【数据融合】【状态估计】基于KF、UKF、EKF、PF、FKF、DKF卡尔曼滤波KF、无迹卡尔曼滤波UKF、拓展卡尔曼滤波数据融合研究(Matlab代码实现)
208 0
|
25天前
|
机器学习/深度学习 边缘计算 算法
【无人机】无人机群在三维环境中的碰撞和静态避障仿真(Matlab代码实现)
【无人机】无人机群在三维环境中的碰撞和静态避障仿真(Matlab代码实现)
136 0
|
26天前
|
机器学习/深度学习 分布式计算 算法
【风场景生成与削减】【m-ISODATA、kmean、HAC】无监督聚类算法,用于捕获电力系统中风场景生成与削减研究(Matlab代码实现)
【风场景生成与削减】【m-ISODATA、kmean、HAC】无监督聚类算法,用于捕获电力系统中风场景生成与削减研究(Matlab代码实现)
114 0
|
1月前
|
存储 编解码 算法
【多光谱滤波器阵列设计的最优球体填充】使用MSFA设计方法进行各种重建算法时,图像质量可以提高至多2 dB,并在光谱相似性方面实现了显著提升(Matlab代码实现)
【多光谱滤波器阵列设计的最优球体填充】使用MSFA设计方法进行各种重建算法时,图像质量可以提高至多2 dB,并在光谱相似性方面实现了显著提升(Matlab代码实现)
|
1月前
|
机器学习/深度学习 传感器 算法
【高创新】基于优化的自适应差分导纳算法的改进最大功率点跟踪研究(Matlab代码实现)
【高创新】基于优化的自适应差分导纳算法的改进最大功率点跟踪研究(Matlab代码实现)
150 14

热门文章

最新文章