サンプル C 言語プログラム(H2)

トップページ
水素分子のページに戻る。

下のソースプログラムをそのままテキストエディタ(メモ帳など)にコピー and ペースト すれば、簡単にコンパイルと実行できる。
(この プログラムは 単純なC言語なので、このテキストエディタを "filename.c" とセーブしてコンパイルしてほしい。)
ここでは 1 MM = 1 × 10-14 meter という新しい単位を使用している。

このサンプルプログラムでは、実行すると、最初に H2 の核間距離 ( MM ) を入力する。
その後、電子1の y座標 "r" ( MM )、x 座標 "a" ( MM ) を 入力する。
これら入力値から、このプログラムは H2 の結合エネルギー (eV)、原子核 1 に作用する トータルの力 (= nucforce )、電子1に作用する力の x, y 成分 (= elefx, elefy ) を出力する。
これらの力は 水素原子の基底状態の電子と原子核間の力との比で表されている。

電子に作用する力と その速度から 1軌道に含まれるドブロイ波を計算する。
( ここでは ビリアル定理 E = -T = 1/2 V を使用している。)
"another x" は 電子1が 力の方向へ移動した先の もう1つの x 座標である。
最初の x 座標は +100 まで自動的に計算ごとに増加していく。


#include <stdio.h>
#include <math.h>

int main(void) 
 {
   int i;
   double nnuc,rr,aa,ke,r,nuc,a, rra,rrb,rrc;
   double poten,ppot, kinetic,velo,binding,elefx,elefy,eleforce;
   double radius, debroglie, wave, nucforce, anotherx;
   double me=9.1093826e-31;
   double pai=3.141592653589793; 
   double epsi=8.85418781787346e-12;
   double h=6.62606896e-34; 
   double ele=1.60217653e-19;  
   double boh=5292.0*5292.0*1.0e-28;  /* boh = ( Bohr radius )^2 */ 
     
                                      /* input  nuc, x, y coordinate */
   
    printf(" Internuclear distance (MM) of H2 molecule ? ");  
    scanf("%lf",&nnuc);
 
    printf(" r (MM) ? = y coordinate of electron 1 ");  
    scanf("%lf",&rr);

    printf(" a (MM) ? = x coordinate of electron 1 ");  
    scanf("%lf", &aa);
    printf("                        \n");

    r = rr * 1.0e-14; nuc = nnuc * 1.0e-14;  /*  change MM to meter */ 
    ke = (ele*ele)/(4.0*pai*epsi);
 
   for (i=1; i < 10 ;i++) {      /* repeat until a=initial a+100 */

    a=aa*1.0e-14;
    rra=sqrt(a*a+r*r);       /* rra =distance between electron 1 and n1 */
    rrb=sqrt(r*r+(nuc-a)*(nuc-a));  /* rrb=  between electron 1 and n2 */
    rrc=sqrt(4*r*r+(nuc-2*a)*(nuc-2*a));  /* between two electrons */
 
    poten=ke*(-2.0/rra-2.0/rrb+1.0/nuc+1.0/rrc);   /* poten=potential energy (J) */
    ppot=poten*6.241509*1.0e18;       /* ppot = potential energy (eV) */ 
                                  
    kinetic=-0.5*poten;      /* kinetic = total kinetic energy of two electrons (J) */
    velo=sqrt(kinetic/me);  /* velo = electron's velocity (m/s) */ 
    binding=-ppot*0.5-13.606*2;       /* binding energy (eV) of H2 */

                            /* foeces acting on electron 1 ( x, -y directions ) */
   elefx= -a/(rra*rra*rra) + (nuc-a)/(rrb*rrb*rrb) - (nuc-2*a)/(rrc*rrc*rrc);
   elefy= r/(rra*rra*rra) + r/(rrb*rrb*rrb) - (2*r)/(rrc*rrc*rrc);
   eleforce = ke*sqrt(elefx*elefx + elefy*elefy);  /* total force */

   radius=(me*velo*velo)/eleforce;  /* rotation radius from centrifugal force */
   debroglie = h/(me*velo);       /* de Broglie wavelength of electron */
   wave=(2*pai*radius)/debroglie;  /* wave's number in one orbit */

   nucforce= a/(rra*rra*rra)-1.0/(nuc*nuc)+(nuc-a)/(rrb*rrb*rrb);
   nucforce=nucforce * boh;            /* total force acting on nucleus 1 */

   anotherx=a+(2.0*r*elefx)/elefy; 
   anotherx = anotherx * 1.0e14;      /* another x coordinate of electron 1 */

   elefx=elefx * boh; elefy = elefy * boh;
 
   printf("a:%.1f ", aa);
   printf(" binding-energy: %.4f ", binding); 
   printf(" elefx: %.3f ", elefx);
   printf(" elefy: %.3f \n", elefy);
 
   printf("nucforce: %.3f ", nucforce);
   printf(" another x : %.2f ", anotherx);
   printf(" de Broglie wave: %.3f \n", wave);
   printf("                        \n");

   aa=aa+10;}                             
   
     return 0;
   }