水利水エネルギー計画—水力発電所ダムの特徴パラメータ選択——C++プログラム


//           .cpp 
#include "stdafx.h"
#include <fstream>
#include <iostream>
#include <iomanip>
#include <cmath>
const int  NumSeries = 20;//           
double ZSeries[NumSeries], VSeries[NumSeries], qSeries[NumSeries];//       

double ZVqChaZhi(double x[NumSeries], double y[NumSeries], double x0)
{//Z、V、q      (   :  q       )
	double y0 = -1;//      ,        
	for(int i = 0; i < NumSeries - 1; i++)
	{
		if(((x[i]>=x0)&&(x[i+1]<=x0))||((x[i]<=x0)&&(x[i+1]>=x0))) 
			y0 = y[i] + (y[i+1] - y[i])/(x[i+1] - x[i])*(x0 -x[i]);
	}
	return y0;
}
#include "       .h"
#include "       .h"
#include "    .h"
#include "    .h"
#include "          .h"

void main()
{
	using namespace std;
	ifstream infile;
	infile.open("infile_ZVqSeries.txt");
	for(int i = 0; i  < NumSeries; i++)//      
		infile>>ZSeries[i]>>VSeries[i]>>qSeries[i];
	infile.close();
//	YuNiV();//      
	WProcess();
	TiaoHong();
//	ChuLiJS();
//	NandERelation();
}
//       .h
//               ,     
const int StartYear = 1951,//
Y = 50;//            
double rho = 0.237,//       (kg/m3)
m_down = 0.9,//     
gamma = 1650,//     
p_empty = 0.3,//   
alpha = 0.15,//            
MonthQ[Y][12],//            
AverageYearQ[Y] = {0},//             
ZS,//   
YuNiVolume = 0;//            ,  10^6(m3)
void YuNiV()
{
	using namespace std;
	ofstream outfile;
	outfile.open("outfile_AverageYearQ.txt");
	ifstream infile;
	infile.open("infile_MonthQ.txt");
	for(int i = 0; i < Y; i++)
	{
		outfile<<setw(10)<<i+StartYear;
		for(int j = 0; j < 12; j++)
		{
			infile>>MonthQ[i][j];
			YuNiVolume += MonthQ[i][j];
			AverageYearQ[i] += MonthQ[i][j];
		}
		AverageYearQ[i] /= 12;
		outfile<<setw(10)<<AverageYearQ[i]<<endl;
	}
	YuNiVolume *= ((1+alpha)*rho*2.63*m_down/(1-p_empty)/gamma);
	cout<<"       :"<<YuNiVolume<<endl;
}
//       .h
//       3   ,      (m3/s*3h)
//W—  ,T—3  ,J—  ,C—Classical,SJ—  
const int SumT = 52,//       
	J = 3,//         
	NumJ[J] = {1, 8, 24},//           ,
	WProT = 24;//          , NumJ[J - 1]
//NumJ[J]      、    、    3     
//      ,         、        、      
double SJWPro[3][SumT],//                  
//SJWPro[3][NumJ[J - 1]]  
CTW[SumT],//      
MaxSJW[3][J],//        ,3      
MaxCW[J],//        
RatioW[J],//          
temp_W;//   
int date[J];//        ,         

void WProcess()
{
	using namespace std;
	ifstream infile;
	infile.open("infile_ClassicalWProcess.txt");
	for(int i = 0; i < SumT; i++)
		infile>>CTW[i];
	infile.close();
	infile.open("infile_MaxSJW.txt");
	for(int k = 0; k<J; k++)
		for(int BZ = 0; BZ < 3; BZ++)
			infile>>MaxSJW[BZ][k];
	infile.close();
	
	//              
	cout<<"      :"<<endl;
	for(int k = 0; k<J; k++)
	{//          
		 MaxCW[k] = 0;//   
		 for(int j = 0; j < SumT - NumJ[k] + 1; j++)
		 {
			 temp_W = 0;//   
			 for(int n = 0; n < NumJ[k]; n++)
				 temp_W += CTW[j + n];
			 if(temp_W >MaxCW[k])
			 {
				MaxCW[k] = temp_W;//(m3/s* )
				date[k] = j;
			 }
		 }
		 if(k==0) cout<<"    :"<<3*date[k];
		 else cout<<"  "<<NumJ[k]<<"          :"
					<<3*date[k]<<" "<<3*(date[k]+NumJ[k] -1);
		 cout<<setw(10)<<"  :"<<MaxCW[k]<<endl;
	}
	cout<<"                     !"
		<<"     1,     !"<<endl;
	cin>>temp_W;
	if(temp_W == 1)
	{
		for(int BZ = 0; BZ < 3; BZ++)
		{//    
			for(int k = 0; k < J; k++)
			{//       
				if( k == 0||k==1) //       
					RatioW[k] = MaxSJW[BZ][k]/MaxCW[k];
				else RatioW[k] = (MaxSJW[BZ][k] - MaxSJW[BZ][k-1])/(MaxCW[k] - MaxCW[k-1]);
				for(int i = date[k]; i < date[k]+NumJ[k]; i++)
				{
					if( k == 0) SJWPro[BZ][i] = CTW[i]*RatioW[k];
					else if(i < date[k-1]|| i>=date[k-1]+NumJ[k-1])
							 SJWPro[BZ][i] = CTW[i]*RatioW[k];
				}
			}
		}
		ofstream outfile;//        
		outfile.open("outfile_WProcess.txt");
		outfile<<"      :"<<endl
			<<setw(10)<<"  "
			<<setw(10)<<"    "
			<<setw(10)<<"    "
			<<setw(10)<<"    "
			<<setw(10)<<"    "<<endl;
		for(int i = date[J - 1]; i < date[J - 1] + WProT; i++)
		{
			outfile<<setw(10)<<3*i<<setw(10)<<CTW[i];
			for(int BZ = 0; BZ < 3; BZ++)
				outfile<<setw(10)<<SJWPro[BZ][i];
			outfile<<endl;
		}
		outfile.close();

		cout<<"        ,  1  !"<<endl;
		cin>>temp_W;//     
	}
}
//    .h
//    ,    :m3/s
//         10^6m3
//         q[i] ,              
//     V[i]               
const double  QToW = 3.0*60.0*60.0/pow(10.0, 6),//         
ZX = 199.3,//     
	SafeQ = 10043,//          
epsilon_q = 0.1,//    
QFaDian = 1000,//      
QBottomHole = 1280;//        ,            
const int  NumHole = 13;//     
double q[SumT],//      
V[SumT],//         
Z[SumT],//    
///q_ZX,//             
q_ZX,//             
ZFH,//     
ZSJ,//     
ZJH,//     
t_qMax, Q_qMax, q_qMax, V_qMax, Z_qMax,//         
q1, q2;//   、       
int DestroyZFH = 0,//          
ShouWei = 0;//   
double qFunction(int BZ, int i,double QStart, double QEnd, double VStart, double qStart, int t)
{//         
//i——   i    
//t          ,              ,      0
	int time = 0;
	using namespace std;
	double VEnd, qEnd = 0;//        
	q1 = qStart;//   
	while(qEnd == 0)
	{
		VEnd =  VStart + pow(0.5, t)*QToW*(QStart + QEnd - qStart - q1)/2;
		q2 = ZVqChaZhi(VSeries, qSeries, VEnd);//      V[i]
		if(q2 < 0)// V[i]       
		{
			cout<<"    "<<VEnd<<endl//1455.94      
				<<"  !!! BZ = "<<BZ<<" q["<<3*i<<"]        !";
			cin>>q2;
		}
		q2 = QFaDian + NumHole*q2;//          
		if(Z[i - 1]>ZSJ) q2 +=1280;//            
		if(fabs(q2 - q1) < epsilon_q) 
			qEnd = q1;
		else	{q1 = (q1 + q2)/2;time++;if(time == 2000) cout<<"Time" ;}
	}
	return qEnd;
}
void qMax(int BZ, int i, double QStart, double QEnd, double VStart, double qStart)
{//      
	using namespace std;
	double qEnd, deltaT = 0;
	int t = 1;//t          ,              ,      0
    qEnd = qFunction(BZ, i, QStart, (QStart+QEnd)/2, VStart, qStart, t);
     if(DestroyZFH ==0&&qEnd > SafeQ) qEnd = SafeQ;
	while(fabs(qEnd - (QStart+QEnd)/2) >= epsilon_q)
	{
		if(t>100)
		{//   100     ,          !
			cout<<"   BZ = "<<BZ<<"     !"<<endl
				<<"     "<<SJWPro[BZ][i - 1]<<endl
				<<"     "<<q[i - 1]<<endl
				<<"     "<<Z[i - 1]<<endl
				<<"     "<<SJWPro[BZ][i]<<endl
				<<"     "<<q[i]<<endl
				<<"     "<<Z[i]<<endl
				<<"            :"<<qEnd<<endl
				<<"        :"<<(QStart+QEnd)/2<<endl; 
			cin>>t;
		}
		if(qEnd> (QStart+QEnd)/2) 
			QEnd = (QStart+QEnd)/2;
		else
		{
			deltaT += pow(0.5, t);
			VStart += pow(0.5, t)*QToW*(QStart + (QStart+QEnd)/2 - qEnd - qStart);
			 qStart = qEnd;
			 QStart =  (QStart+QEnd)/2;
		}
		t++;
		qEnd = qFunction(BZ, i, QStart, (QStart+QEnd)/2, VStart, qStart, t);
		if(DestroyZFH ==0&&qEnd > SafeQ) qEnd = SafeQ;
	}
	deltaT += pow(0.5, t);
	VStart += pow(0.5, t)*QToW*(QStart + (QStart+QEnd)/2 - qEnd - qStart);
	t_qMax = 3*(i - 1+deltaT);
	Q_qMax = (QStart+QEnd)/2;
	q_qMax = qEnd;
	V_qMax = VStart;
	Z_qMax = ZVqChaZhi(VSeries, ZSeries, V_qMax);
	if(BZ==0) {ZFH = Z_qMax; cout<<"     ZFH ="<<ZFH<<endl;}
	if(BZ==1) {ZSJ = Z_qMax; cout<<"     ZSJ ="<<ZSJ<<endl;}
	if(BZ==2) {ZJH = Z_qMax;cout<<"     ZJH ="<<ZJH<<endl;}
}
void TiaoHong()
{//     
	using namespace std;
	system("cls");
	ifstream infile;
	infile.open("infile_XiuYunWProcess.txt");
	for(int i = date[J - 1]; i < date[J - 1] + WProT; i++)
	{//            
		for(int BZ = 0; BZ < 3; BZ++)
			infile>>SJWPro[BZ][i];
	}
	q_ZX =QFaDian + NumHole*ZVqChaZhi(ZSeries, qSeries, ZX); 
	ZFH = ZSJ = 100000;//               ,      
	infile.close();
	ofstream outfile ;
	outfile.open("outfile_TiaoHong.txt");
	outfile<<"    :"<<ZX<<endl//199.3
		<<"           :"<<q_ZX<<endl//6502.9
		<<"     :"<<ZVqChaZhi(ZSeries, VSeries, ZX)<<endl;
	for(int BZ =0; BZ < 3; BZ++)
	{//                  ,     
		outfile<<endl<<endl<<endl;//(3)
		if(BZ==0) outfile<<"       :";
		if(BZ==1) outfile<<"       :";
		if(BZ==2) outfile<<"      :";
		outfile<<endl<<setw(10)<<"  "<<setw(10)<<"    "
		 <<setw(10)<<"    "<<setw(10)<<"    "<<setw(10)<<"    "<<endl;
		SJWPro[BZ][date[J - 1] - 1] = q[date[J - 1] - 1] = q_ZX;
		V[date[J - 1] - 1] = ZVqChaZhi(ZSeries, VSeries, ZX);
		int qMaxGet = 0;//           
		DestroyZFH = ShouWei = 0;//          
		for(int i = date[J - 1]; i < date[J - 1] + WProT; i++)
		{
			if((SJWPro[BZ][i]<q_ZX&&qMaxGet == 0)||ShouWei == 1)//i = 8~13       SJWPro[0][13] = 5628.24
			{//        
				q[i] = SJWPro[BZ][i];//			
				V[i] = ZVqChaZhi(ZSeries, VSeries, ZX);
				Z[i] = ZX;
			}
			else
			{//      
				q[i] = qFunction(BZ, i, SJWPro[BZ][i-1],SJWPro[BZ][i], V[i-1], q[i-1], 0);
				V[i] =  V[i - 1] + QToW*(SJWPro[BZ][i - 1] + SJWPro[BZ][i] - q[i - 1] - q[i])/2;
				Z[i] = ZVqChaZhi(VSeries, ZSeries, V[i]);
				if(Z[i]>ZFH) DestroyZFH =1;
				if(DestroyZFH ==0&&q[i] > SafeQ) 
				{//        ,        
					q[i] = SafeQ;
					V[i] =  V[i - 1] + QToW*(SJWPro[BZ][i - 1] + SJWPro[BZ][i] - q[i - 1] - q[i])/2;
					Z[i] = ZVqChaZhi(VSeries, ZSeries, V[i]);
				}
				if(qMaxGet == 1&&q[i]<q_ZX) 
				{
//					q[i] = SJWPro[BZ][i];//      ,           
//					V[i] = V[i - 1];
//					Z[i] = Z[i - 1];
					Z[i] = ZX;
					V[i] = ZVqChaZhi(ZSeries, VSeries, ZX);
					q[i] = SJWPro[BZ][i - 1] + SJWPro[BZ][i]  - q[i - 1]- (V[i] - V[i - 1])*2.0/QToW;
					ShouWei = 1;
//					V[i] =  V[i - 1] + QToW*(SJWPro[BZ][i - 1] + SJWPro[BZ][i] - q[i - 1] - q[i])/2;
				}
				if(q[i]>=SJWPro[BZ][i]&&qMaxGet == 0)//        
				{
					qMax(BZ, i, SJWPro[BZ][i -1], SJWPro[BZ][i], V[i-1], q[i -1]);
					qMaxGet = 1;//            
				}
			}
			outfile<<setw(10)<<3*i<<setw(10)<<SJWPro[BZ][i]
			   <<setw(10)<<q[i]<<setw(10)<<V[i]<<setw(10)<<Z[i]<<endl;
		}
		outfile<<"         :"<<endl
			<<setw(10)<<t_qMax<<setw(10)<<Q_qMax
			<<setw(10)<< q_qMax<<setw(10)<< V_qMax<<setw(10)<< Z_qMax;
	}
	outfile.close();
	cout<<"    ,   !"<<endl;
	cin>>temp_W;//     
}
//    .h
//  :  10 ^6m3,  m3/s,    kw
//           
const double N_XZ = 68.4;//    (        ,       )100000;//
const int HLowSeries = 14,//            
HUpSeries = 60;//            
const double kFactor = 8.2,//    
CostH = 0.5,//    
VDead = 185.68;//   ,    
double deltaH,//   QLoss_XZ = N_XZ
ComeQ[12], LossQ[12],//  、    
ThrowQ[12],//  
VAbs[12],//      
 ChuLi[12],//  
HLow_h[HLowSeries],HLow_Q[HLowSeries],
HUp_h[HUpSeries],HUp_V[HUpSeries];
double HLowChaZhi(double x[HLowSeries], double y[HLowSeries], double x0)
{//        
	double y0 = -1;//      ,        
	for(int i = 0; i < HLowSeries - 1; i++)
	{
		if(((x[i]>=x0)&&(x[i+1]<=x0))||((x[i]<=x0)&&(x[i+1]>=x0))) 
			y0 = y[i] + (y[i+1] - y[i])/(x[i+1] - x[i])*(x0 -x[i]);
	}
	return y0;
}
double HUpChaZhi(double x[HUpSeries], double y[HUpSeries], double x0)
{//        ,    10^6m3
	double y0 = -1;//      ,        
	for(int i = 0; i < HUpSeries - 1; i++)
	{
		if(((x[i]>=x0)&&(x[i+1]<=x0))||((x[i]<=x0)&&(x[i+1]>=x0))) 
			y0 = y[i] + (y[i+1] - y[i])/(x[i+1] - x[i])*(x0 -x[i]);
	}
	return y0;
}
void ChuLiJS()
{//    
	using namespace std;
	ifstream infile;
	infile.open("infile_HLowSeries.txt");
	for(int i = 0; i < HLowSeries; i++)
		infile>>HLow_h[i]>>HLow_Q[i];
	infile.close();
	infile.open("infile_HUpSeries.txt");
	for(int i = 0; i < HUpSeries; i++)
		infile>>HUp_h[i]>>HUp_V[i];
	infile.close();
	infile.open("infile_ComeLossQ.txt");
	for(int i = 0; i < 12; i++)
		infile>>ComeQ[i]>> LossQ[i];
	infile.close();
	ofstream outfile;
	outfile.open("outfile_ChuLi.txt");
	outfile<<setw(10)<<"  "<<setw(10)<<"  "<<endl;
	for(int i = 0; i < 12; i++)
	{
		if(i == 0) 
		{
			VAbs[i] = VDead + (ComeQ[i] - LossQ[i])*2.63;
			deltaH = HUpChaZhi(HUp_V, HUp_h, (VAbs[i]+VDead)/2.0);
		} 
		else 
		{
			VAbs[i] = VAbs[i - 1] + (ComeQ[i] - LossQ[i])*2.63;	
			 deltaH= HUpChaZhi(HUp_V, HUp_h, (VAbs[i-1]+VAbs[i])/2.0);
		}
		deltaH -= (HLowChaZhi(HLow_Q, HLow_h, LossQ[i]) + CostH);
		ChuLi[i] =deltaH*kFactor* LossQ[i]/10000.0;//    kw
		ThrowQ[i] = 0;
		 if(ChuLi[i] > N_XZ)
		 {
			 ChuLi[i] = N_XZ;
			ThrowQ[i] = LossQ[i] - 10000.0*N_XZ/kFactor/deltaH;
		 }
		outfile<<setw(10)<<ChuLi[i]<<setw(10)<<ThrowQ[i]<<endl;
	}
}
//          .h
//           
//      : kw      :  kw     :  kw*h
const int NumN = 20;//         
const double MaxWorkLoad[3] = {260, 255, 250};//12、1、2          ,    
double N_g[NumN],//         (         )
E_g[NumN],//                 
DayWorkLoad[24],//                
LowLine;
void NandERelation()
{//               
	using namespace std;
	ifstream infile;
	infile.open("infile_DayWorkLoad.txt");
	for(int i = 0; i < 24; i++)
		infile>>DayWorkLoad[i];
	infile.close();
	ofstream outfile;
	outfile.open("outfile_NandERelation.txt");
	cout<<"           "<<endl;
	cin>>N_g[0];//          
	cout<<"           "<<endl;
	cin>>N_g[NumN  - 1];//          
	for(int i = 0; i < NumN; i++)
	{
		N_g[i] = N_g[0] + (double)(N_g[NumN  - 1] - N_g[0])/(NumN - 1)*i;
		E_g[i] = 0;
		for(int Month = 0; Month < 3; Month++)
		{
			LowLine =MaxWorkLoad[0] - N_g[i];
			for(int Hour = 0; Hour < 24; Hour++)
				if((DayWorkLoad[Hour]/100.0*MaxWorkLoad[Month])>LowLine)
					E_g[i] += ((DayWorkLoad[Hour]/100.0*MaxWorkLoad[Month]) - LowLine);
		}
		E_g[i] *= (30.4/10000.0);
		outfile<<setw(10)<<N_g[i]<<setw(10)<<E_g[i]<<endl;
	}
	outfile.close();
}