Chemistry 版 (精华区)
发信人: zjliu (秋天的萝卜), 信区: Chemistry
标 题: boliwa计算器的未经测试的C语言版本.
发信站: BBS 哈工大紫丁香站 (Thu Jul 8 16:37:11 2004)
发信站: 日月光华
//-------------------------------------------------------------
// 配分函数计算器V1.0 原产自：boliwa@chemistry@日月光华@20040706
//
// 以下为翻译的C语言版本
// 按照boliwa的Basic程序依样画葫芦翻译了一下
// 其实是不知道怎么算的，寒:)
// Borland C编译器4.12 windows版编译通过
// 按照标准的c来写，应该可以在linux下被编译
//
// ！！未经调试和测试！！没希望这个能真正工作，仅供大家参考一下
// 我不是搞这个的，调试不来，呵呵
//
// ㊣20040707
//-------------------------------------------------------------
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>
//-------------------------------------------------------------
//-------------------------------------------------------------
const double Boltzmann = 1.380662E-23;
const double Planck = 6.62616E-34;
//-------------------------------------------------------------
struct TCoordinate
{
double x;
double y;
double z;
};
//-------------------------------------------------------------
// 将字符串Str的第Start个位置开始的Count个字符转换为浮点数
double VAL(const char* Str,int Start,int Count)
{
char Buffer[256];
Str = Str + Start -1;
strncpy(Buffer,Str,Count);
return atof(Buffer);
}
//-------------------------------------------------------------
// 以FileName(文件名）和Temp(温度)为输入参数
// 其他的为输出参数
// 如果函数成功执行则返回true，否则返回false
// 如果函数成功执行则返回true，否则返回false
bool printq(const char* FileName, double Temp, double* p_qt,
double* p_qv, double* p_qr, double* p_Epoint,
double* p_qv2, double* p_TotQ, double* p_TotQ2)
{
const double Pi = 3.141592653589793;
const double C = 299792458;
const Pressure = 101325;
const u = 1.66053873E-34;
const MAX_LINE_LENGTH = 1024;
double VibFreq[256];
double MolMass;
TCoordinate RotConst;
int RotSymm, VibFree;
char GetStr[MAX_LINE_LENGTH];
int I=0;
FILE *File=fopen(FileName,"r");
if( File==NULL )
return false; //error open file
while( !feof(File) )
{
{
k=0;
do
{
fscanf(File,"%c",&Char);
if(Char=='\n')
{
GetStr[k++] = 0;
break;
}
GetStr[k++] = Char;
}while(!feof(File));
if( strncmp(GetStr," Molecular mass:",16 )==0 ) //equal
{
MolMass = VAL(GetStr,18,11);
(*p_qt) = (Boltzmann*Temp/Pressure)*
pow(sqrt(2*Pi*MolMass*u*Boltzmann*Temp)/Planck,3);
}else if( strncmp(GetStr," Sum of electronic and \
zero-point Energies=", 43)==0 )
{
(*p_Epoint) = VAL(GetStr,50,20);
}else if( strncmpi(GetStr," ROTATIONAL SYMMETRY NUMBER",
}else if( strncmpi(GetStr," ROTATIONAL SYMMETRY NUMBER",
27)==0 )
{
RotSymm = (int)(VAL(GetStr,29,2));
}else if( strncmpi(GetStr," Rotational constant",20)==0 )
{
RotConst.x = VAL(GetStr,36,11);
RotConst.y = VAL(GetStr,48,11);
RotConst.z = VAL(GetStr,60,11);
if( RotConst.x!=0 )
RotConst.x = Planck /(8 * Pi * Pi * RotConst.x * 1E9);
else
RotConst.x = -1;
if( RotConst.y!=0 )
RotConst.y = Planck /(8 * Pi * Pi * RotConst.y * 1E9);
else
RotConst.y = -1;
if( RotConst.z!=0 )
RotConst.z = Planck /(8 * Pi * Pi * RotConst.z * 1E9);
else
RotConst.z = -1;
printf("\nRotConst:%g,%g,%g\n",
RotConst.x,RotConst.y,RotConst.z);
RotConst.x,RotConst.y,RotConst.z);
}else if( strncmp(GetStr," Frequencies --",16)==0 )
{
I++;
VibFreq[0] = VAL(GetStr,17,10);
if( strlen(GetStr)>26 )
{
I++;
VibFreq[1] = VAL(GetStr,40,10);
if( strlen(GetStr)>49 )
{
I++;
VibFreq[2] = VAL(GetStr,63,10);
}
}
} // end of if
}//end of while( !feof(File) )
fclose(File);
if( (RotConst.x>0) && (RotConst.y>0) )
(*p_qr) = 8*Pi*Pi*pow(sqrt(2*Pi*Boltzmann*Temp),3)*
sqrt(RotConst.x*RotConst.y*RotConst.z)/
RotSymm/pow(Planck,3);
RotSymm/pow(Planck,3);
else
(*p_qr) = 8*Pi*Pi*RotConst.x*Boltzmann*Temp/
RotSymm / pow(Planck,2);
(*p_qv) = 1;
(*p_qv2) = 1;
VibFree = I;
for(I=1;I<=VibFree;I++)
if( VibFreq[I]>0 )
(*p_qv2) = (*p_qv2)*1/
(1 - exp(-Planck * VibFreq[I]*100*C/Boltzmann/Temp));
(*p_TotQ) = (*p_qt)*(*p_qv)*(*p_qr);
(*p_TotQ2) = (*p_qt)*(*p_qv2)*(*p_qr);
return true;
}
//-------------------------------------------------------------
// 计算并输出各个数据
// 四个参数分别为 温度，过渡态文件名，产物一， 产物二
// 成功执行函数返回true，否则返回false
bool CalcPS(double Temperature,const char* TsFileName,
bool CalcPS(double Temperature,const char* TsFileName,
const char* P1FileName, const char* P2FileName)
{
double qt[256], qv[256], qr[256], Epoint[256], qv2[256];
double TotQ[256], TotQ2[256], k, E0;
int i;
//过渡态
if( !printq(TsFileName, Temperature, &qt[0], &qv[0],
&qr[0], &Epoint[0], &qv2[0], &TotQ[0], &TotQ2[0]) )
return false;
//产物一
if( !printq(P1FileName, Temperature, &qt[1], &qv[1],
&qr[1], &Epoint[1], &qv2[1], &TotQ[1], &TotQ2[1]) )
return false;
//产物二
if( !printq(P2FileName, Temperature, &qt[2], &qv[2],
&qr[2], &Epoint[2], &qv2[2], &TotQ[2], &TotQ2[2]) )
return false;
for(i=0;i<3;i++)
{
printf("qt=%g\n",qt[i]);
printf("qt=%g\n",qt[i]);
printf("qv=%g\n",qv[i]);
printf("qr=%g\n",qr[i]);
printf("eq=%g\n",Epoint[i]);
printf("qv2=%g\n",qv2[i]);
printf("totq=%g\n",TotQ[i]);
printf("totq2=%g\n",TotQ2[i]);
}
E0 = (Epoint[0] - (Epoint[1] + Epoint[2])) * 27.2 * 96485;
printf("E0%g\n",E0);
k = Boltzmann * Temperature / Planck * TotQ2[0] / TotQ[1]
/ TotQ[2] * exp( -E0 / 8.31 / Temperature);
printf("k%g\n",k);
return true;
}
//-------------------------------------------------------------
// 程序入口
// 输入三个文件名为参数
int main(int argc, char* argv[])
{
CalcPS(298.15,argv[1],argv[2],argv[3]);
return 0;
}
printf("totq2=%g\n",TotQ2[i]);
//-------------------------------------------------------------
// 程序入口
// 输入三个文件名为参数
int main(int argc, char* argv[])
{
CalcPS(298.15,argv[1],argv[2],argv[3]);
return 0;
}
//-------------------------------------------------------------
// end of source file
//-------------------------------------------------------------
--
╔═══════════════════╗
║★★★★★友谊第一 比赛第二★★★★★║
╚═══════════════════╝
※ 来源:·哈工大紫丁香 http://bbs.hit.edu.cn·[FROM: 202.118.229.*]
Powered by KBS BBS 2.0 (http://dev.kcn.cn)
页面执行时间：2.927毫秒