!*********************************************************************** ! 熱源モジュール Thermal Pool ! 多数の分子間でランダムにエネルギーを交換し、仮想的な熱源を実現する ! Compile: ! g95 -c thpool_N4.f90 ! g95 thpool_N4.o mtfort90.o -o thpool_N4 ! Date: 2007/06/04 MODULE thpool_N4 USE mtmod ! MT乱数モジュール使用 !INTEGER, PARAMETER :: N_MOL = 6 ! 分子の数 INTEGER, PARAMETER :: N_MOL = 4 ! 分子の数 REAL*8, PARAMETER :: MAX_EXCHANGE = 2.5d0 ! 1回に交換する最大値 REAL*8 :: Mol(1:N_MOL) CONTAINS !********************************************************************** ! 初期化 ! 初期平均エネルギー ! 初期の全エネルギーを返す REAL*8 FUNCTION pool_init(e0) IMPLICIT NONE REAL*8, INTENT(IN) :: e0 CALL sgrnd(4357) ! ランダムシード Mol = e0 ! 配列の全要素を初期化 CALL pool_run( 1000 ) ! 十分に安定するまで回す pool_init = N_MOL * e0 ! 最初に与えた全エネルギーを返す return END FUNCTION pool_init !********************************************************************** ! 熱浴内でのエネルギー交換を進める ! 1stepで単純に N_MOLの交換を行う SUBROUTINE pool_run(step) IMPLICIT NONE INTEGER, INTENT(IN) :: step INTEGER :: st, i, j REAL*8 :: chg DO st = 0, step * N_MOL ! 異なる2つの分子 i, j を選び出す i = int(grnd() * N_MOL) + 1 DO j = int(grnd() * N_MOL) + 1 IF (i /= j) THEN EXIT END IF END DO ! i -> j にエネルギーを渡す ! 交換する大きさを決める chg = MAX_EXCHANGE * grnd() IF (chg > Mol(i)) THEN ! その大きさを交換できなければ Mol(j) = Mol(j) + Mol(i) ! iを全てjに渡す Mol(i) = 0.0 ELSE Mol(j) = Mol(j) + chg Mol(i) = Mol(i) - chg END IF END DO END SUBROUTINE pool_run !********************************************************************** ! エネルギーの受け渡しを行う REAL*8 FUNCTION pool_get(i) IMPLICIT NONE INTEGER, INTENT(IN) :: i ! 対象となる分子 pool_get = Mol(i) END FUNCTION pool_get SUBROUTINE pool_set(i,e0) IMPLICIT NONE INTEGER, INTENT(IN) :: i ! 対象となる分子 REAL*8, INTENT(IN) :: e0 Mol(i) = e0 END SUBROUTINE pool_set !********************************************************************** ! 熱浴の全エネルギーを返す REAL*8 FUNCTION pool_total() IMPLICIT NONE INTEGER :: i REAL*8 :: sum sum = 0.0d0 DO i=0, N_MOL sum = sum + Mol(i) END DO pool_total = sum return END FUNCTION pool_total !********************************************************************** ! 熱浴内の様子を出力する SUBROUTINE pool_list IMPLICIT NONE INTEGER, PARAMETER :: N_HIST = 20 INTEGER :: hist(1:N_HIST) INTEGER :: i, pos REAL*8 :: sum REAL*8, PARAMETER :: DELTA = 0.5d0 ! ヒストグラムの刻み幅 hist = 0 ! 配列の0クリアー sum = 0.0 DO i=0, N_MOL pos = int( Mol(i) / DELTA ) + 1 IF (pos > N_HIST) THEN pos = N_HIST END IF hist(pos) = hist(pos) + 1 sum = sum + Mol(i) END DO PRINT *, hist(1:N_HIST), ':', sum END SUBROUTINE pool_list END MODULE thpool_N4 !********************************************************************** ! テスト用メイン !!PROGRAM main !!USE thpool_N4 !!IMPLICIT NONE !! INTEGER :: i !! !! CALL pool_init(10.0) !! CALL pool_list !! !! DO ! i=0,100 !! CALL pool_run(100) !! CALL pool_list !! END DO !! !!END PROGRAM main !*** END of FILE *******************************************************