2006/08/18 | 源岩降解动力学过程在C++的初步实现(1)
类别(程序语言) | 评论(0) | 阅读(48) | 发表于 22:22
 
//源岩降解动力学过程在C++的初步实现
//---------------------------------------------------------
// !!CORPYRIGHT DECLARATION!!
//This program was wirrten by Zeng H.S using std C++ language.
//You can use the program under GNU licence.
// 2006/08/11(YY/MM/DD)
//------------------------------------------------------------
#include <iostream>
#include <math.h>
#include<conio.h>
#include <stdlib.h>
#define R 8.31441
#define T2 2
#define T_initial 25

using namespace std;

int main()
{
//---------------------------------------------------------
  //begin varialbe declaration
  long double delta_t,time1,time2,A,Ea;
//delta_t是累积时间,time是时间
   //A是指前因子,Ea是活化能
  long double T[200],TOC_initial,C_initial,Tx,alpha_x=2,x[200],k;
  //end varialble declaration
  //---------------------------------------------------------
  //begin data initiation
  T[0]=T_initial;   //initial temperature
  Tx=alpha_x+15*log10(alpha_x/2); //heating rate
  Ea=259.1*1E3; //A type kerogen Ea
  A=41.23E17; //frequency
  //----------------------------------------------------------
  //user input
L1: cout<<"Please input Heating Rate(C/Ma):"<<endl
   <<"Heating Rate(C/Ma)=";
  cin>>alpha_x;
  cout<<endl;
  cout<<"Please input Ea(KJ/Mol):"<<endl
   <<"Ea(KJ/Mol)=";
  cin>>Ea;
  Ea=Ea*1000;
  cout<<endl<<"Please input A(frequency factor):"<<endl
   <<"A=";
  cin>>A;
  cout<<endl;  
    
  
  //----------------------------------------------------------
  //end data initiation
  //----------------------------------------------------------
  //begin calculation
  int i=0;
  time1=0;
  x[i]=1;
  while(time1<=100*(1E6)*365*24*3600)
  {
    i++;
    time2=time1+(1E6)*365*24*3600;
    T[i]=T[i-1]+Tx;
    k=A*exp(-Ea/R/(T[i]+273));
    x[i]=x[i-1]*exp(-k*(time2-time1));
    time1=time2;
  }
  
  
  //end calculation
  //----------------------------------------------------------
  //begin output the results
cout<<"Temperature(C) "<<"X(%)"<<endl;
  for(i=0;i<=250;i++)
{
   cout<<" "<<T[i]<<" "<<x[i]<<endl;
if(x[i]<=1E-10) break;
}
  cout<<"Calculate another source rock? Enter Y to continue or N to end the program"<<endl;
  char loop;
  cin>>loop;
  if(loop=='y'||loop=='Y')
  goto L1;
  else
  {
    cout<<"Press any key to quit out";
    getch();
  }
  
  //end output the results
              
}
0

评论Comments

日志分类
首页[227]
+随笔[68]
学习札记[67]
科研工具[44]
文献摘译[2]
程序语言[20]
网络技巧[26]
--[0]