基于DP法的水电站中长期优化调度

//毕设之优化调度程序---动态规划
#include<stdio.h>
#include<stdlib.h>
#include<math.h>

#define z_ss  328.5    //死水位
#define z_qfh  346.30   //防洪限制水位-前汛期
#define z_hfh  346.80   //防洪限制水位-后汛期
#define z_zc  347      //正常蓄水位
#define Y   48        //流量资料年数
#define M   12        //一年月数
#define T   576       //入流流量数据个数
#define NUM1   24        //水位库容关系点数
#define NUM2   12        //下泄流量与水位关系点数
#define NUM3   3         //水头限制出力关系点数
#define  dz    0.05      //水位离散步长
#define NUM4   400  //一个阶段最大离散点数，>=（347-328.5）/0.05
#define fa_xishu 10     //罚系数

float Q[T+1];//来水流量
float vcz[NUM1],zsy_cz[NUM1];//库容水位插值
float qcz[NUM2],zxy_cz[NUM2];//下泄流量水位插值
float z_mao[NUM3],n_yx[NUM3];//水头预想出力插值
float Qfd[T+1];//发电流量
float Qfd1,Qfd2;//发电流量及过渡值
float qs[T+1];//弃水流量
float z[T+1][NUM4];//各时段离散水位值
float zsx[T+1];//各时段水位上限值
float zzy[T+2];//最有水位记录
float N[T+1][NUM4];//各时段各方案出力
float Nzy[T+1];//最优出力记录
float vzy[T+2];
float N_guodu;//时段出力过渡值
float E[T+1][NUM4];//各时段各方案发电量
float Et[T+1][NUM4];//各阶段电量之和
float Ezy[T+1];//最优方案最优出力记录
float qs[T+1];//弃水流量
float Emax;//各时段最大发电量
int numzy[T+1][NUM4];//最优方案点数记录
int num_zy[T+1];//最优点数记录
int n_yt[T+1];//各月天数
float wky[M+1][Y+1];//阶段可用水量
float wkyave[M+1];//可用水量平均值
float Nzyave[M+1];//出力平均值
float wfangcha[M+1]={0};//可用水量方差
float Nfangcha[M+1]={0};//出力方差
float wNlicha[M+1]={0};//可用水量与出力的离差
float wjunfangcha[M+1];//可用水量均方差
float Njunfangcha[M+1];//出力均方差
float wNhe[M+1];//可用水量与出力的方差之积
float r[M+1],b0[M+1],b1[M+1],b2[M+1];//回归方程特征值
float wkyt[M+1]={0};//可用水量之和
float Nzyt[M+1]={0};//出力之和

//插值函数部分
float chazhi(float x,float y1[],float y2[],int n)
{
    int i;
    for(i=0;i<n-1;i++)  //正常值情况处理方式
    {
        if(x<=y1[i+1]&&x>=y1[i])
            return(y2[i]+(y2[i+1]-y2[i])/(y1[i+1]-y1[i])*(x-y1[i]));
    }
    if(x<y1[0])      //上溢出处理方式
            return y2[0];
    if(x>y1[n-1])    //下溢出处理方式
            return (y2[n-2]+(y2[n-1]-y2[n-2])/(y1[n-1]-y1[n-2])*(x-y1[n-2]));
}

//输入函数部分
void data_in()
{
    int i,j;
    int m=1;
 FILE *fQsw,*fv_z,*fq_z,*fn_z,*fmyt;
 if((fQsw=fopen("水文年排列的流量.txt","r"))==NULL||(fv_z=fopen("水位库容.txt","r"))==NULL||(fq_z=fopen("水位流量.txt","r"))==NULL||(fn_z=fopen("水头出力.txt","r"))==NULL||(fmyt=fopen("每月天数.txt","r"))==NULL)
 {
  printf("Can not open data file !\n");
  exit(0);
 }
 for(i=0;i<Y;i++)//输入来流数据
 {
  for(j=0;j<M;j++)
  {
   fscanf(fQsw,"%f",&Q[m++]);
  }
 }
 for(i=0;i<NUM1;i++