scqjdtfx.cpp

这是本人两年前兼职为某个公司做的石油钻进设计软件

// Scqjdtfx.cpp: implementation of the CScqjdtfx class.
//
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
#include "cvenus.h"
#include "Scqjdtfx.h"

#include "yqsxz.h"
#include "math.h"
#include "qjyljs.h"
#include "LZXJD.h"
#include "direct.h"//用于得到系统路径
#include "Common.h"

#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
#define new DEBUG_NEW
#endif


#define MAXLEN	152		// 200 


extern	char	cGraphColor[][16];
extern	int		iGraphColorNum;
extern	int		g_iProgress ;	// 计算进度
extern	bool	bFinishThread;


//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

CScqjdtfx::CScqjdtfx()
{
	g_iProgress = 0;
}

CScqjdtfx::~CScqjdtfx()
{

}

void CScqjdtfx::Dtfx()
{

	double dblCSZ[4][11];
	int iFXCount[4];

	//D地层深度, P地层, T地层, P井口, T井口, P回压, AF, BF, CF, NF
	double dblDDCSD,dblPDC,dblTDC,dblPJK,dblTJK,dblPHY,dblAF,dblBF,dblCF,dblNF;
	// 地层有效厚度, 气井供给边界半径, 井底半径, 视表皮系数, 有效渗透率, Qtest, Pwftest, Ptftest, Qgmax
	double dblDCYXHD,dblQJGGBJBJ,dblJDBJ,dblSBPXS,dblYXSTL,dblQTest,dblPwftest,dblPtftest,dblQgmax;
	// DIC, DIT1, DVT1, DIT2, DVT2, RR, GG, CO2, H2S, N2, LGR, GL
	double dblDIC,dblDIT1,dblDVT1,dblDIT2,dblDVT2,dblRR,dblGG,dblCO2,dblH2S,dblN2,dblLGR,dblGL;
	
	CString strTemp;
	int	iTemp,iTemp01,iTemp02,iTemp03;

	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(1, 5);
	dblDDCSD=atof(strTemp);

	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(2, 5);
	dblPDC=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(3, 5);
	dblTDC=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(4, 5);
	dblPJK=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(5, 5);
	dblTJK=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(6, 5);
	dblPHY=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(7, 5);
	dblAF=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(8, 5);
	dblBF=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(9, 5);
	dblCF=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(10, 5);
	dblNF=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(11, 5);
	dblDCYXHD=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(12, 5);
	dblQJGGBJBJ=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(13, 5);
	dblJDBJ=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(14, 5);
	dblSBPXS=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(15, 5);   
	dblYXSTL=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(16, 5);  
	dblQTest=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(17, 5);   
	dblPwftest=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(18, 5);   
	dblPtftest=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(19, 5);   
	dblQgmax=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(20, 5);   
	dblDIC=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(21, 5);   
	dblDIT1=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(22, 5);   
	dblDVT1=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(23, 5);   
	dblDIT2=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(24, 5);   
	dblDVT2=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(25, 5);   
	dblRR=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(26, 5);   
	dblGG=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(27, 5);   
	dblCO2=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(28, 5);   
	dblH2S=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(29, 5);   
	dblN2=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(30, 5);   
	dblLGR=atof(strTemp);
	strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(31, 5);   
	dblGL=atof(strTemp);
	

	iTemp = 31;
    for(iTemp01 = 1; iTemp01 <= 3;iTemp01++)
	{
		for(iTemp02 = 1; iTemp02 <= 10; iTemp02++)
		{
			if(iTemp01 == 3 && iTemp02 > 4)
				break;

			 iTemp++;
			 strTemp = ((CCVenusApp*)AfxGetApp())->m_pMainDlg->m_grid.GetTextRC(iTemp, 5);
			 dblCSZ[iTemp01][iTemp02]=atof(strTemp);
		}
	}

	
	for(iTemp01 = 1; iTemp01 <= 3; iTemp01++)
	{
		for(iTemp02 = 1; iTemp02 <= 9; iTemp02++)
		{
			if(dblCSZ[iTemp01][iTemp02] == 0)
			{
				for(iTemp03 = iTemp02 + 1; iTemp03 <= 10; iTemp03++)
				{
					if(dblCSZ[iTemp01][iTemp03] >= 0)
					{
						dblCSZ[iTemp01][iTemp02] = dblCSZ[iTemp01][iTemp03];
						dblCSZ[iTemp01][iTemp03] = 0;

						break;
					}
				}
			}
		}
	}
	
	g_iProgress = 10;

	int iNKK;
	double dblAXY = 0;

	for(iTemp01 = 1; iTemp01 <= 3; iTemp01++)
	{
		if(iTemp01 == 3)
			iNKK = 5;
		else
			iNKK = 10;

		for(iTemp02 = 1; iTemp02 <= iNKK - 1; iTemp02++)
		{
			for(iTemp03 = iTemp02 + 1; iTemp03 <= iNKK; iTemp03++)
			{
				if((iTemp01 == 1 && dblCSZ[iTemp01][iTemp02] > dblCSZ[iTemp01][iTemp03] && dblCSZ[iTemp01][iTemp03] > 0)
					|| (iTemp01 != 1 && dblCSZ[iTemp01][iTemp02] < dblCSZ[iTemp01][iTemp03]))
				{
					dblAXY = dblCSZ[iTemp01][iTemp03];
					dblCSZ[iTemp01][iTemp03] = dblCSZ[iTemp01][iTemp02];
					dblCSZ[iTemp01][iTemp02] = dblAXY;
				}
			}
		}
	}

	g_iProgress = 20;

	if(dblCSZ[1][1] == 0)	
		dblCSZ[1][1] = dblPJK;
	if(dblCSZ[2][1] == 0)	
		dblCSZ[2][1] = dblPDC;

	for(iTemp01 = 1; iTemp01 <= 3; iTemp01++)
	{
		iFXCount[iTemp01] = 0;
		if(iTemp01 == 3)
			iNKK = 5;
		else
			iNKK = 10;
		for(iTemp02 = 1; iTemp02 <= iNKK; iTemp02++)
		{
			if(dblCSZ[iTemp01][iTemp02] > 0)
				iFXCount[iTemp01]++;
		}

	}
	iFXCount[3]++;

	g_iProgress = 30;

  // 估算: 分析参数个数(1), 分析参数个数(2) <= 21,取21, 分析参数个数(3) <= 21,取21
	double dblQQL[MAXLEN], dblPLR[21][MAXLEN],dblPLC[21][21][MAXLEN],dblQL[21][21][21];

	double dblPDC00 = 0;
	double dblLC = 0;
	double dblTG = 0;

	// 参数变换及计算有关参数
	dblTJK += 273.15;
	dblTDC += 273.15;
	dblH2S /= 100;
	dblCO2 /= 100;
	dblN2 /= 100;
	dblLGR /= 10000;
	dblPDC00 = dblPDC;
	dblDIT1 /= 1000;
	dblDIT1 /= 1000;
	dblDIC /= 1000;
	dblLC = dblDDCSD - dblDVT1 - dblDVT2;
	dblTG = (dblTDC - dblTJK) / dblDDCSD;

	double dblPC, dblTC;
    double dblDIT10,dblDVT10,dblDIT20,dblDVT20,dblDIT30,dblDVT30,dblQG;
	double dblTR, dblPR,dblUG,dblZ,dblBT;
	double dblGPX = 0;
	double dblGP = 0;
	double dblTun = 0;
	qjyljs qjyljsz;

	dblPC = 0;
	dblTC = 0;
	yqsxz csjs;
	//天然气的临界参数(PC, TC, GG, H2S, CO2, "干气")
	csjs.trqnjcs(dblGG, dblH2S, dblCO2,1);//临界csz
	dblPC = csjs.Pc;//对比压力
	dblTC = csjs.Tc;//对比温度
	
	if(abs(int(dblAF * dblBF)) > 0)
		dblQgmax = (-dblAF + pow((pow(dblAF,2) + 4 * dblBF * (pow(dblPDC,2) - pow(0.0101,2))),0.5)) / 2 / dblBF;
		//dblQgmax = (-dblAF + pow((pow(dblAF,2) + 4* dblBF * (pow(dblPDC,2) - pow(0.0101,2)), 0.5) / 2 / dblBF;
	else if(dblQTest * dblPwftest  >0 || dblQTest * dblPtftest  >0)
	{
		if(dblPwftest <= 0)
		{
			dblDIT10 = dblDIT1;
			dblDVT10 = dblDVT1;
			dblDIT20 = dblDIT2;
			dblDVT20 = dblDVT2;
			dblDIT30 = dblDIC;
			dblDVT30 = dblLC;
			dblQG = dblQTest * 10000;
			
			//Call 气井压力(G压力X, Ptftest, T井口, TG, DVT10, DIT10, 0, RR, QG, GG, LGR, GL, H2S, N2, CO2, 0, 1, 1)
			dblGPX = qjyljsz.qjyl(dblPtftest, dblTJK, dblTG, dblDVT10, dblDIT10, 0, dblRR, dblQG, dblGG, dblLGR, dblGL, dblH2S, dblN2, dblCO2,0, 0, 1, 1);
			if( dblDVT20 == 0)
			{
				dblGP = dblGPX;
			}
			else
			{
				dblTun = dblTJK + dblTG * dblDVT10;
				//Call 气井压力(Pwftest, G压力X, Tun, TG, DVT20, DIT20, 0, RR, QG, GG, LGR, GL, H2S, N2, CO2, 0, 1, 1)
				dblPwftest=qjyljsz.qjyl(dblGPX, dblGP, dblTG, dblDVT20, dblDIT20, 0, dblRR, dblQG, dblGG, dblLGR, dblGL, dblH2S, dblN2, dblCO2,0, 0, 1, 1);
			}

			if(dblLC == 0)
				dblPwftest = dblGP;
			else
			{
				dblTun = dblTJK + dblTG * (dblDVT10 + dblDVT20);
				// Call 气井压力(Pwftest, G压力, Tun, TG, DVT30, DIT30, 0, RR, QG, GG, LGR, GL, H2S, N2, CO2, 0, 1, 1)
				dblPwftest=qjyljsz.qjyl(dblGP, dblTun, dblTG, dblDVT30, dblDIT30, 0, dblRR, dblQG, dblGG, dblLGR, dblGL, dblH2S, dblN2, dblCO2,0, 0, 1, 1);
			}

		}
		
		if(dblQgmax <= 0)
		{
			double dblPD = 0;
			dblPD = (pow(dblPDC,2) - pow(dblPwftest,2)) / pow(dblPDC,2);
			double dblQgmax1,dblQgmax2,dblQgmax3;
			dblQgmax1 = dblQTest / (1.8 * dblPD - 0.8 * pow(dblPD,2));
			dblQgmax2 = 6 * dblQTest / (pow((1 + 48 * dblPD),0.5) - 1);
			dblQgmax3 = dblQTest / 1.0434 / pow(dblPD,0.6594);
			dblQgmax = (dblQgmax1 + dblQgmax2 + dblQgmax3) / 3;
			
		}
		dblBF = ((pow(dblPDC,2) - pow(0.1,2)) / dblQgmax - (pow(dblPDC,2) - pow(dblPwftest,2)) / dblQTest) / (dblQgmax - dblQTest);
		dblAF = (pow(dblPDC,2) - pow(dblPwftest,2) - dblBF * pow(dblQTest,2)) / dblQTest;

	}
	else if(dblCF * dblNF > 0)
	{
		dblQgmax = dblCF * pow((pow(dblPDC,2) - pow(0.1,2)),dblNF);
	}
	else if(dblDCYXHD * dblQJGGBJBJ * dblJDBJ * dblYXSTL > 0)
	{
		
		dblTR = dblTDC / dblTC;
		dblPR = dblPDC / dblPC;
		//Call 天然气的粘度(UG, PR, TR, CO2, H2S, N2, GG, T地层 - 273.15)
        dblUG=csjs.trqnd(dblPR, dblTR, dblCO2, dblH2S, dblN2, dblGG, dblTC - 273.15); 	
		//Call 天然气偏差系数(Z, PR, TR)
        dblZ=csjs.trqpcxs(dblPR, dblTR);
		dblAF = 0.000127 * dblUG * dblZ * dblTDC * (log(dblQJGGBJBJ / dblJDBJ) / log(10) + dblSBPXS) / dblYXSTL / dblDCYXHD;
		dblBT = 76440000000.0 / pow(dblYXSTL,1.5);
		dblBF = 2.3336*pow(10,-21) * dblBT * dblGG * dblTDC * dblZ * (1 / dblJDBJ - 1 / dblQJGGBJBJ) / pow(dblDCYXHD,2);
		dblQgmax = (-dblAF + pow((pow(dblAF,2) + 4 * dblBF * (pow(dblPDC,2) - pow(0.0101,2))),0.5)) / 2 / dblBF;

	}
	else
	{
		MessageBox(NULL, "生产气井动态分析IPR参数值错误，请更正后再计算！","参数值错误",MB_ICONWARNING);
		return;
	}
	
	g_iProgress = 40;

	dblQgmax = int(1000 * dblQgmax) / 1000;
	int iJSDS;
	// hss
	iJSDS = 150;
	double dblQQLZZ;
	dblQQLZZ = int(1000 * dblQgmax / iJSDS) / 1000;
	if(dblQQLZZ < 0.5)
		dblQQLZZ = 0.5;
	iJSDS = 1 + int(dblQgmax / dblQQLZZ);
	
	// 打开保存计算结果的文件
	CStdioFile stdfRes,stdfGraph;
	CFileException fe;
	BOOL bStatus;
	CString strPath;

	if(((CCVenusApp*)AfxGetApp())->m_cBaseFileLeft.IsEmpty())
	{	
		strTemp = ((CCVenusApp*)AfxGetApp())->GetOpenFileName(1);
		if(strTemp.IsEmpty()) return;
	}
	else
	{
		strTemp = ((CCVenusApp*)AfxGetApp())->m_cBaseFileLeft +  _T(".RGS");
	}

	bStatus = stdfRes.Open(strTemp,CFile::modeCreate|CFile::modeWrite);
	if(!bStatus)
	{
	   TCHAR szMsg[100];
	   sprintf(szMsg, "Error opening file for reading. Code:%d",
		   fe.m_cause);
	   AfxMessageBox(szMsg);
	   return;

	}
	//打开保存计算的图形结果的文件
	if(((CCVenusApp*)AfxGetApp())->m_cBaseFileLeft.IsEmpty())
	{	
		strTemp = ((CCVenusApp*)AfxGetApp())->GetOpenFileName(2);
		if(strTemp.IsEmpty()) return;
	}
	else
	{
		strTemp = ((CCVenusApp*)AfxGetApp())->m_cBaseFileLeft +  _T(".GGS");
	}

	bStatus = stdfGraph.Open(strTemp,CFile::modeCreate|CFile::modeWrite);
	if(!bStatus)
	{
	   TCHAR szMsg[100];
	   sprintf(szMsg, "Error opening file for reading. Code:%d",
		   fe.m_cause);
	   AfxMessageBox(szMsg);
	   return;

	}

	stdfRes.WriteString(_T("\n"));
	stdfRes.WriteString(_T("\t\t$$生产气井动态分析计算结果$0\n"));
	stdfRes.WriteString(_T("\n"));
	stdfRes.WriteString(_T("井口压力$7\t 油管直径$7"));
	
	int iPDC;
	for(iPDC = 1; iPDC <= iFXCount[2]; iPDC++)
	{
		stdfRes.WriteString(_T("\t地层压力$7"));
	}
	stdfRes.WriteString(_T("\n"));

	stdfRes.WriteString(_T("  MPa$11\t   mm$11"));
	for(iPDC = 1; iPDC <= iFXCount[2]; iPDC++)
	{
		strTemp.Format(_T("\t%.2f$11"), dblCSZ[2][iPDC]);
		stdfRes.WriteString(strTemp);
	}
	stdfRes.WriteString(_T("\n"));

	
	int iPercent;
	iPercent = 2;
	// 主窗口.SSPanel1.FloodPercent = 百分比
	
	// 流入关系计算
	double dblPrs,dblT,dblZ1,dblUG1,dblZ0,dblUG0,dblAF1,dblBF1,dblPWF;
	for(iPDC = 1; iPDC <= iFXCount[2]; iPDC++)
	{
		dblPrs = dblCSZ[2][iPDC];;
		for(iTemp02 = 0; iTemp02 <= iJSDS; iTemp02++)
		{
			dblQG = int(1000 * dblQQLZZ * iTemp02) / 1000;