You cannot select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
nmWTAI-Platform/3rd/Pebi/include/pch.h

479 lines
16 KiB
C++

This file contains ambiguous Unicode characters!

This file contains ambiguous Unicode characters that may be confused with others in your current locale. If your use case is intentional and legitimate, you can safely ignore this warning. Use the Escape button to highlight these characters.

#pragma once
#ifndef PCH_H
#define PCH_H
#include "framework.h"
#endif //PCH_H
#define HX_API extern "C" _declspec(dllexport)
#include <iostream>
#include <vector>
#include <cmath>
#include <algorithm>
#include <fstream>
#include <ctime>
#include <iomanip>
#include <unordered_set>
#include <Windows.h>
#ifndef M_PI
const double M_PI = acos(-1.0);
#endif
typedef std::vector<std::vector<std::vector<double>>>dVec3; //三维数组:double
typedef std::vector<std::vector<double>>dVec2; //二维数组:double
typedef std::vector<std::vector<int>>iVec2; //二维数组:int
typedef std::vector<double>dVec1; //一维数组:double
typedef std::vector<int>iVec1; //一维数组:int
template<class T> void HX_copy(std::vector<std::vector<std::vector<T>>>& p1, const std::vector<std::vector<std::vector<T>>>& p0)
{
int m = p0.size(); p1.resize(m);
for (int i = 0; i < m; ++i)
{
int n = p0[i].size(); p1[i].resize(n);
for (int j = 0; j < n; ++j)
{
int l = p0[i][j].size(); p1[i][j].resize(l);
for (int k = 0; k < l; ++k)
{
p1[i][j][k] = p0[i][j][k];
}
}
}
}
template<class T> void HX_copy(std::vector<std::vector<T>>& p1, const std::vector<std::vector<T>>& p0)
{
int m = p0.size(); p1.resize(m);
for (int i = 0; i < m; ++i)
{
int n = p0[i].size(); p1[i].resize(n);
for (int j = 0; j < n; ++j)
{
p1[i][j] = p0[i][j];
}
}
}
template<class T> void HX_copy(std::vector<T>& p1, const std::vector<T>& p0)
{
int m = p0.size(); p1.resize(m);
for (int i = 0; i < m; ++i)
{
p1[i] = p0[i];
}
}
template<class T> void HX_copy(std::vector<std::vector<std::vector<T>>>& p1, T*** p0, int m, int* n, int l)
{
p1.resize(m);
for (int i = 0; i < m; ++i)
{
p1[i].resize(n[i]);
for (int j = 0; j < n[i]; ++j)
{
p1[i][j].resize(l);
for (int k = 0; k < l; ++k)
{
p1[i][j][k] = p0[i][j][k];
}
}
}
}
template<class T> void HX_copy(std::vector<std::vector<T>>& p1, T** p0, int m, int n)
{
p1.resize(m);
for (int i = 0; i < m; ++i)
{
p1[i].resize(n);
for (int j = 0; j < n; ++j)
{
p1[i][j] = p0[i][j];
}
}
}
template<class T> void HX_copy(std::vector<T>& p1, T* p0, int m)
{
p1.resize(m);
for (int i = 0; i < m; ++i)
{
p1[i] = p0[i];
}
}
//点结构体
struct point
{
//点结构体
double x; double y; //点坐标
point() { x = 0; y = 0; }
~point() {}
void set(const double& x_ = 0, const double& y_ = 0) { x = x_; y = y_; }
void set(const point& p) { x = p.x; y = p.y; }
};
struct point3
{
double x;
double y;
double z;
point3() { x = 0; y = 0; z = 0; }
~point3() {}
point3(const dVec1&p) { x = p[0]; y = p[1]; z = p[2]; }
point3(const double& x_ = 0, const double& y_ = 0, const double&z_ = 0) { x = x_; y = y_; z = z_;}
void set(const double& x_ = 0, const double& y_ = 0, const double& z_ = 0) { x = x_; y = y_; z = z_; }
void set(const point3&p) { x = p.x; y = p.y; z = p.z; }
};
//网格结构体
struct cell
{
//网格单元结构体
std::vector<point> p;
iVec2 pindex;
iVec1 isplot;
cell() {}
~cell() {}
};
//网格算法输入参数结构体
struct HX_NWTM_GRID_INPUT
{
// 网格划分算法输入参数结构体
dVec2 Boundary; //3D:{x0, y0, z0, x1, y1, z1} //2D:{x0, y0, x1, y1} 边界数据
dVec2 VerticalWell; //3D:{x0, y0, z0, x1, y1, z1, rw} //2D:{x0, y0, x1, y1, rw} 直井数据
dVec2 HorizontalWell; //3D:{x0, y0, z0, x1, y1, z1, rw} 水平井数据
dVec2 FractureVerticalWell; //3D:{x0, y0, z0, x1, y1, z1, wf} //2D:{x0, y0, x1, y1, wf, FC} 压裂直井数据(wf裂缝半宽,m,FC,裂缝导流能力,mD.m(FC为0时为无限导流大于0时为有限导流))
dVec3 MultistageFracturedHorizontalWell; //3D:{x0, y0, z0, x1, y1, z1, wf} //2D:{x0, y0, x1, y1, wf, FC} 多级压裂水平井数据(wf裂缝半宽,m,FC,裂缝导流能力,mD.m(FC为0时为无限导流大于0时为有限导流))
dVec2 InclinedWell; //3D:{x0, y0, z0, x1, y1, z1, rw} 斜井数据
dVec2 Fault; //3D:{x0, y0, z0, x1, y1, z1} //2D:{x0, y0, x1, y1} 断层数据
double GridControl; // 网格大小控制参数
int D; // 维数
//默认初始化
HX_NWTM_GRID_INPUT()
{
dVec1 a(3), b(4), c(5), d(6);
Boundary.resize(4);
b[0] = -1500.0; b[1] = -1500.0; b[2] = -1500.0; b[3] = 1500.0; Boundary[0] = b;
b[0] = -1500.0; b[1] = 1500.0; b[2] = 1500.0; b[3] = 1500.0; Boundary[1] = b;
b[0] = 1500.0; b[1] = 1500.0; b[2] = 1500.0; b[3] = -1500.0; Boundary[2] = b;
b[0] = 1500.0; b[1] = -1500.0; b[2] = -1500.0; b[3] = -1500.0; Boundary[3] = b;
VerticalWell.resize(3);
a[0] = 0; a[1] = 0; a[2] = 0.1; VerticalWell[0] = a;
a[0] = 1000; a[1] = 1000; a[2] = 0.1; VerticalWell[1] = a;
a[0] = -1000; a[1] = -1000; a[2] = 0.1; VerticalWell[2] = a;
HorizontalWell.resize(0);
FractureVerticalWell.resize(1);
d[0] = -200; d[1] = -200; d[2] = 200; d[3] = -200; d[4] = 0.05; d[5] = 0; FractureVerticalWell[0] = d;
MultistageFracturedHorizontalWell.resize(1);
MultistageFracturedHorizontalWell[0].resize(3, dVec1(5));
d[0] = -600; d[1] = 600; d[2] = -400; d[3] = 600; d[4] = 0.1; d[5] = 0; MultistageFracturedHorizontalWell[0][0] = d;
d[0] = -600; d[1] = 400; d[2] = -400; d[3] = 400; d[4] = 0.1; d[5] = 0; MultistageFracturedHorizontalWell[0][1] = d;
d[0] = -600; d[1] = 200; d[2] = -400; d[3] = 200; d[4] = 0.1; d[5] = 0; MultistageFracturedHorizontalWell[0][2] = d;
InclinedWell.resize(0);
Fault.resize(1);
b[0] = -500; b[1] = 1000; b[2] = 500; b[3] = 500; Fault[0] = b;
GridControl = 150.0;
D = 2;
}
~HX_NWTM_GRID_INPUT() {}
};
//网格算法输出参数结构体(绘图用)
struct HX_NWTM_GRID_OUTPUT1
{
cell TRI_cell; //三角形网格
cell PEBI_cell; //PEBI网格
HX_NWTM_GRID_OUTPUT1() {}
~HX_NWTM_GRID_OUTPUT1() {}
};
//网格算法输出参数结构体(模型用)
struct HX_NWTM_GRID_OUTPUT2
{
dVec2 Trinodexy;
dVec1 Area;
dVec2 D;
struct {
int n;
dVec2 XiLinw;
dVec2 lw;
dVec2 dw;
dVec1 rw;
iVec2 inwell;
dVec1 dwell;
} ZhiJingNeiBianJie;
struct {
int n;
dVec2 XiLinf;
dVec2 lf;
dVec2 df;
dVec1 xf;
iVec2 infra;
iVec1 nf;
iVec1 jjf;
iVec1 jjfl;
dVec1 lfcd;
iVec2 infra1;
dVec2 lf1;
dVec2 df1;
} LieFengJingNeiBianJie;
struct {
int n;
dVec2 XiLinh;
dVec2 lh;
dVec2 dh;
dVec2 dsxf;
iVec2 inhor;
iVec1 nhor;
dVec2 areah;
iVec2 inhor1;
iVec2 nh;
iVec2 jjh;
dVec2 hfcd;
iVec2 jjhl;
iVec2 jjh2;
iVec2 inhor2;
dVec2 lh1;
dVec2 dh1;
} DuoJiYaLieShuiPingJingNeiBianJie;
struct {
int n;
dVec2 WaiBianh;
dVec2 WaiBianl;
dVec2 WaiBiand;
} WaiBianJie;
struct {
int n;
dVec2 faultb1;
dVec2 faultb2;
dVec2 faultl1;
dVec2 faultd1;
} NeiBuDuanCeng;
struct {
iVec1 ia;
iVec1 ja;
iVec2 nzeros;
int numk;
} YuChuLiJuZhen;
HX_NWTM_GRID_OUTPUT2() {}
~HX_NWTM_GRID_OUTPUT2() {}
};
//KRINGING插值输入参数结构体
struct HX_KRING_INPUT
{
double nugget; //块金值:表示空间点在零距离处的变异程度,即测量误差和小于采样尺度的随机变异之和
double sill; //基台值:表示变差函数随距离增加而趋于稳定的极限值,反映区域化变量的总变异程度
double range; //变程:表示空间相关性的有效距离。当两点间距离超过 range 时,它们之间不再具有空间相关性
double model; //变差函数模型:指定变差函数的数学形式,描述空间相关性随距离的变化规律{, , }
//高斯模型SPHERICAL(0):相关性随距离增加呈指数衰减,适用于连续性较强的变量
//指数模型EXPONENTIAL(1):相关性快速衰减,适用于局部变异性较大的变量
//球状模型GAUSSIAN(2):在变程内呈抛物线变化,超过变程后相关性为零
dVec2 p; //插值点
dVec2 v; //
~HX_KRING_INPUT() {}
HX_KRING_INPUT(double nugget0, double sill0, double range0, double model0 , const dVec2& p0, const dVec2& v0)
{
nugget = nugget0; sill = sill0; range = range0; model = model0;
p = p0;
v = v0;
}
};
//KRINGING插值输出参数结构体
struct HX_KRING_OUTPUT
{
dVec1 v;
HX_KRING_OUTPUT() {}
~HX_KRING_OUTPUT() {}
};
//数值试井模型求解器输入参数结构体
struct HX_NWTM_MODEL_INPUT
{
int T; //1:油单相常数pvt; 2:油单相变化pvt; 3:水单相常数pvt; 4:水单相变化pvt; 5:气单相变化pvt; 6:气单相拟压力; 7:油气两相; 8:油水两相; 9:气水两相; 10:油气水三相
HX_NWTM_GRID_OUTPUT2 GRID;
struct Rate //流量数据
{
dVec2 t; //时间, h [一口井一组数]
dVec2 qo; //油流量,m^3/d [一口井一组数]
dVec2 qg; //气流量,m^3/d [一口井一组数]
dVec2 qw; //水流量,m^3/d [一口井一组数]
}Rate;
struct Pressure //压力数据
{
dVec2 t; //时间, h [一口井一组数]
dVec2 p; //压力, MPa [一口井一组数]
}Pressure;
struct CS //井储表皮数据
{
dVec1 C; //井储, m^3/MPa [一口井一个数]
dVec1 S; //表皮, [一口井一个数]
}CS;
struct PVT //流体性质数据
{
dVec1 p; //压力, MPa
double pb; //饱和压力, MPa
dVec1 Rso; //溶解气油比, m^3/m^3
dVec1 Bo; //油体积系数, m^3/m^3
dVec1 Co; //油压缩系数, 1/MPa
dVec1 miuo; //油粘度, mPa·s
dVec1 rouo; //油密度, kg/m^3
dVec1 Rv; //凝析油气比, m^3/m^3
dVec1 Bg; //气体积系数, m^3/m^3
dVec1 Cg; //气压缩系数, 1/MPa
dVec1 miug; //气粘度, mPa·s
dVec1 roug; //气密度, kg/m^3
dVec1 Z; //气偏差因子, 1
dVec1 Rsw; //溶解气水比, m^3/m^3
dVec1 Bw; //水体积系数, m^3/m^3
dVec1 Cw; //水压缩系数, 1/MPa
dVec1 miuw; //水粘度, mPa·s
dVec1 rouw; //水密度, kg/m^3
dVec1 V; //吸附气量, m^3/kg
dVec1 k_kinitial; //渗透率比, 1
dVec1 Cf_Cfinitial; //岩石压缩系数比, 1
dVec1 So; //油饱和度
dVec1 Kro; //油相对渗透率
dVec1 Sg; //气饱和度
dVec1 Krg; //气相对渗透率
dVec1 Sw; //水饱和度
dVec1 Krw; //水相对渗透率
}PVT;
struct Base //基础数据
{
double Pi; //初始压力, MPa
double Cti; //综合压缩系数, 1/MPa
double Cf; //岩石压缩系数, 1/MPa
double Soi; //初始含油饱和度
double Sgi; //初始含气饱和度
double Swi; //初始含水饱和度
dVec1 k; //渗透率, D [一个网格单元一个值]
dVec1 phi; //孔隙度, 1 [一个网格单元一个值]
dVec1 h; //储层厚度, m [一个网格单元一个值]
double d; //时间增长指数
double dt_Min; //最小时间间隔, h
double dt_Max; //最大时间间隔, h
}Base;
//初始化
HX_NWTM_MODEL_INPUT() {}
~HX_NWTM_MODEL_INPUT() {}
HX_NWTM_MODEL_INPUT(const HX_NWTM_GRID_OUTPUT2& p0)
{
T =1;
GRID = p0;
Rate.t.resize(5);
Rate.qo.resize(5);
Rate.qw.resize(5);
Rate.qg.resize(5);
Rate.t[0].resize(2); Rate.t[0][0] = 2000; Rate.t[0][1] = 500;
Rate.qo[0].resize(2); Rate.qo[0][0] = 10; Rate.qo[0][1] = 0;
Rate.qw[0].resize(2); Rate.qw[0][0] = 2; Rate.qw[0][1] = 0;
Rate.qg[0].resize(2); Rate.qg[0][0] = 20000; Rate.qg[0][1] = 0;
Rate.t[1].resize(0);
Rate.qo[1].resize(0);
Rate.qw[1].resize(0);
Rate.qg[1].resize(0);
Rate.t[2].resize(3); Rate.t[2][0] = 1000; Rate.t[2][1] = 1000; Rate.t[2][2] = 500;
Rate.qo[2].resize(3); Rate.qo[2][0] = 30; Rate.qo[2][1] = 40; Rate.qo[2][2] = 20;
Rate.qw[2].resize(3); Rate.qw[2][0] = 3; Rate.qw[2][1] = 4; Rate.qw[2][2] = 2;
Rate.qg[2].resize(3); Rate.qg[2][0] = 30000; Rate.qg[2][1] = 40000; Rate.qg[2][2] = 20000;
Rate.t[3].resize(2); Rate.t[3][0] = 1500; Rate.t[3][1] = 1000;
Rate.qo[3].resize(2); Rate.qo[3][0] = 30; Rate.qo[3][1] = 20;
Rate.qw[3].resize(2); Rate.qw[3][0] = 5; Rate.qw[3][1] = 2;
Rate.qg[3].resize(2); Rate.qg[3][0] = 50000; Rate.qg[3][1] = 20000;
Rate.t[4].resize(2); Rate.t[4][0] = 1000; Rate.t[4][1] = 1500;
Rate.qo[4].resize(2); Rate.qo[4][0] = -50; Rate.qo[4][1] = -60;
Rate.qw[4].resize(2); Rate.qw[4][0] = -2; Rate.qw[4][1] = -5;
Rate.qg[4].resize(2); Rate.qg[4][0] = -20000; Rate.qg[4][1] = -50000;
Pressure.t.resize(0);
Pressure.p.resize(0);
CS.C.resize(5);
CS.C[0] = 0.1; CS.C[1] = 0.1; CS.C[2] = 0.1; CS.C[3] = 0.1; CS.C[4] = 0.1;
CS.S.resize(5);
CS.S[0] = 0.1; CS.S[1] = 0.1; CS.S[2] = 0.1; CS.S[3] = 0.1; CS.S[4] = 0.1;
PVT.p = dVec1(200, 0); for (int i = 0; i < 200; ++i) { PVT.p[i] = (i + 1.0); }
PVT.pb = 40.0;
PVT.Rso = dVec1(200, 0);
PVT.Bo = dVec1(200, 1.2);
PVT.Co = dVec1(200, 5e-4);
PVT.miuo = dVec1(200, 0.5);
PVT.rouo = dVec1(200, 800);
PVT.Rv = dVec1(200, 0);
PVT.Bg = dVec1(200, 5e-3);
PVT.Cg = dVec1(200, 2e-2);
PVT.miug = dVec1(200, 2e-2);
PVT.roug = dVec1(200, 200);
PVT.Z = dVec1(200, 1);
PVT.Rsw = dVec1(200, 0);
PVT.Bw = dVec1(200, 1.05);
PVT.Cw = dVec1(200, 1e-4);
PVT.miuw = dVec1(200, 0.8);
PVT.rouw = dVec1(200, 1000);
PVT.V = dVec1(200, 0);
PVT.k_kinitial = dVec1(200, 1);
PVT.Cf_Cfinitial = dVec1(200, 1);
PVT.So = dVec1(81, 0); for (int i = 0; i < 81; ++i) { PVT.So[i] = (0.1+i*0.01); }
PVT.Kro = dVec1(81, 0); for (int i = 0; i < 81; ++i) { PVT.Kro[i] = (i*0.0125); }
PVT.Krw = dVec1(81, 0); for (int i = 0; i < 81; ++i) { PVT.Krw[i] = (1 - i * 0.0125); }
PVT.Sg = dVec1(100, 0);
PVT.Krg = dVec1(100, 0);
PVT.Sw = dVec1(100, 0);
Base.Pi = 40.0;
Base.Cti = 1e-3;
Base.Cf = 1e-4;
Base.Soi = 0.8;
Base.Sgi = 0.0;
Base.Swi = 0.2;
Base.k = dVec1(p0.Trinodexy.size(), 0.001);
Base.phi = dVec1(p0.Trinodexy.size(), 0.1);
Base.h = dVec1(p0.Trinodexy.size(), 10);
Base.d = 1.05;
Base.dt_Min = 0.0025;
Base.dt_Max = 12.5;
}
};
//数值试井模型求解器输出参数结构体
struct HX_NWTM_MODEL_OUTPUT
{
dVec1 t; //时间, h
dVec2 pw; //井底压力, MPa [一口井一组数]
dVec2 p; //压力分布, MPa [一个时间一组数]
dVec2 So; //油饱和度分布 [一个时间一组数]
dVec2 Sg; //气饱和度分布 [一个时间一组数]
dVec2 Sw; //水饱和度分布 [一个时间一组数]
dVec2 k; //渗透率分布,mD [一个时间一组数]
HX_NWTM_MODEL_OUTPUT() {}
~HX_NWTM_MODEL_OUTPUT() {}
};
HX_API void HX_NWTM_GRID(HX_NWTM_GRID_OUTPUT1& p1, HX_NWTM_GRID_OUTPUT2& p2, const HX_NWTM_GRID_INPUT& p0, std::string LIC); //数值试井网格接口
HX_API void HX_NWTM_KRINGING(HX_KRING_OUTPUT& p1, const HX_KRING_INPUT p0, std::string LIC); //数值试井非均质性计算接口
HX_API void HX_NWTM_MODEL(HX_NWTM_MODEL_OUTPUT& p1, const HX_NWTM_MODEL_INPUT& p0, std::string LIC); //数值试井模型求解器接口