random-simulation-3.cpp

（随机模拟）Stochastic Simulation（已加注）

// Monte Carlo Simulation
// Written by Microsoft Visual C++
// Copyright by UTLab @ Tsinghua University
// http://orsc.edu.cn/UTLab

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include<iostream.h>
#include<fstream.h>
#include<iomanip.h>
#include "UTLab.h"

// 最终得到的是 数组中80%的元素数值 要大于的某一个数 

#define n_unit 32     // 机组总台数
#define p_total 2850  // 峰值负荷

main()
{
    int i,j,n=5000; // 模拟次数
    double sum[5001];
	double   rand_number,effect_cap;  // 单次模拟得到的装机总量
	
	double unit_cap[n_unit+1];   // 机组容量数组 
	double unit_p[n_unit+1];     // 强迫停运率数组 
	double unit_state[n_unit+1]; // 机组状态数组
//  double ic=0.0;               // 装机容量   
    double LOLP=0.0,sys_state;            // 系统停电状态标记 0表示不停电 1表示停电  

	ifstream fin("rts.txt");   // 打开文件
	  
	for(i=1;i<=n_unit;i++)
	{
	     fin>>unit_cap[i]>>unit_p[i]; // 每台机组的容量,故障率
//       ic=ic+unit_cap[i];
	}
//       cout<<"装机容量 ic= "<<ic<<endl;
   
	ofstream fout("chahao.txt");   // 输出保存结果	
	
	for(j=1;j<=n;j++)
{
		effect_cap=0.0; //单次模拟得到的有效容量
  	 
        for(i=1; i<=n_unit; i++)  // 得到 每次抽样后的机组状态和总容量
		{
		   rand_number=(double)rand()/(RAND_MAX);
//    cout<<i<<"机组对应的随机数是"<<rand_number<<"故障率是"<<unit_p[i]<<endl;
	    
		   if(rand_number>unit_p[i])
		       unit_state[i]=1;
	        else
               unit_state[i]=0;	
 	    
			 effect_cap=effect_cap+unit_state[i]*unit_cap[i]; //模拟一次的发电机有效容量
		}
		   if(effect_cap<p_total)		   
			  sys_state=1.0;
		   else
			  sys_state=0.0;
		   
		   LOLP=LOLP+sys_state/double(n);

		   fout<<j<<"次模拟的effect_cap="<<effect_cap<<"系统停运标记"<<sys_state<<endl;	
	       sum[j]=effect_cap;

//        cout<<"j="<<j<<" "<<"sum[j]="<<sum[j]<<endl;

}
           cout<<"LOLP="<<LOLP<<endl;

           printf("L = %3.4f\n",findmaxn(sum,1,5000,4999)); 
 
	// sum是数组的首地址,最小,最大, 整数部分(p*n)
    // findmaxn的目的是找出在sum数组中排在第 2400 位的元素数值是多少?
    
    fin.close();
	fout.close();
	return 1;
}