MoDeNa  1.0
Software framework facilitating sequential multi-scale modelling
pt1.f90
1 !WWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWW
2 !
3 !WWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWW
4 
5 SUBROUTINE f_pt1 ( fres )
6 
7  USE eos_variables, ONLY: pi, t, parame, mseg, rho
8  IMPLICIT NONE
9 
10  !-----------------------------------------------------------------------------
11  REAL, INTENT(OUT) :: fres
12  !-----------------------------------------------------------------------------
13  REAL :: i1_dft, i2_dft, order1
14  !-----------------------------------------------------------------------------
15 
16  order1 = mseg(1)*mseg(1)*parame(1,2)**3 *parame(1,3) / t
17 
18  CALL f_dft (i1_dft,i2_dft)
19 
20  fres = 2.0 * pi * rho * i1_dft * order1
21 
22 END SUBROUTINE f_pt1
23 
24 
25 
26 !WWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWW
27 !
28 !WWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWW
29 
30 SUBROUTINE phipt1 ( my_pt1, p_pt1 )
31 
32  USE eos_variables, ONLY: kbol, t, rho
33  IMPLICIT NONE
34 
35  !-----------------------------------------------------------------------------
36  REAL, INTENT(OUT) :: my_pt1
37  REAL, INTENT(IN OUT) :: p_pt1
38 
39  !-----------------------------------------------------------------------------
40  REAL :: fres, pgesdz, pgesd2, pgesd3
41  REAL :: zres
42  !-----------------------------------------------------------------------------
43 
44 
45  !-----------------------------------------------------------------------------
46  ! density iteration
47  !-----------------------------------------------------------------------------
48  CALL pressure_pt1 ( p_pt1, pgesdz, pgesd2, pgesd3 )
49 
50  !-----------------------------------------------------------------------------
51  ! compressibility factor z = p/(kT*rho)
52  !-----------------------------------------------------------------------------
53  zres = (p_pt1 * 1.e-30) / (kbol*t*rho)
54  ! zres = zges - 1.0
55 
56  !-----------------------------------------------------------------------------
57  ! f res
58  !-----------------------------------------------------------------------------
59  CALL f_pt1 ( fres )
60 
61  my_pt1 = fres + zres
62 
63 END SUBROUTINE phipt1
64 
65 
66 !WWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWW
67 !
68 !WWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWW
69 
70 SUBROUTINE pressure_pt1 ( pges, pgesdz, pgesd2, pgesd3 )
71 
72  USE eos_variables, ONLY: kbol, t, rho, eta, tfr
73  IMPLICIT NONE
74 
75  !-----------------------------------------------------------------------------
76  REAL, INTENT(OUT) :: pges
77  REAL, INTENT(OUT) :: pgesdz
78  REAL, INTENT(OUT) :: pgesd2
79  REAL, INTENT(OUT) :: pgesd3
80 
81  !-----------------------------------------------------------------------------
82  REAL :: dzetdv, dicht, dist, fact, z3t
83  REAL :: fres1, fres2, fres3, fres4, fres5, fres
84  REAL :: df_dr, df_drdr
85  REAL :: tfr_1, tfr_2, tfr_3, tfr_4, tfr_5
86  !-----------------------------------------------------------------------------
87 
88  dzetdv = eta*rho
89  z3t = eta/rho
90 
91  IF (eta > 1d-1) THEN
92  fact = 1.0
93  ELSE IF (eta <= 0.1.AND.eta > 0.01) THEN
94  fact = 10.0
95  ELSE
96  fact = 100.0
97  END IF
98  dist = eta*3.e-3 *fact
99 
100  dicht = eta
101  eta = dicht - 2.0*dist
102  rho = eta/z3t
103  CALL f_pt1 ( fres )
104  fres1 = fres
105  tfr_1 = tfr
106  eta = dicht - dist
107  rho = eta/z3t
108  CALL f_pt1 ( fres )
109  fres2 = fres
110  tfr_2 = tfr
111  eta = dicht + dist
112  rho = eta/z3t
113  CALL f_pt1 ( fres )
114  fres3 = fres
115  tfr_3 = tfr
116  eta = dicht + 2.0*dist
117  rho = eta/z3t
118  CALL f_pt1 ( fres )
119  fres4 = fres
120  tfr_4 = tfr
121  eta = dicht
122  rho = eta/z3t
123  CALL f_pt1 ( fres )
124  fres5 = fres
125  tfr_5 = tfr
126 
127  !-----------------------------------------------------------------------------
128  ! ptfr = (-tfr_4+8.0*tfr_3-8.0*tfr_2+tfr_1)/(12.0*dist)
129  ! & *dzetdv*(KBOL*T)/1.E-30
130  ! ztfr =ptfr /( rho * (KBOL*t) / 1.E-30)
131  ! ptfrdz = (-tfr_4+16.0*tfr_3-3.E1*tfr_5+16.0*tfr_2-tfr_1)
132  ! & /(12.0*(dist**2 ))* dzetdv*(KBOL*T)/1.E-30
133  ! & + (-tfr_4+8.0*tfr_3-8.0*tfr_2+tfr_1)
134  ! & /(12.0*dist) * 2.0 *eta*6.0/PI/D
135  ! & *(KBOL*T)/1.E-30
136  ! ztfrdz=ptfrdz/( rho*(kbol*T)/1.E-30 ) - ztfr/eta
137  ! write (*,*) eta,ztfr,ztfrdz
138 
139  ! dtfr_dr = (-tfr_4+8.0*tfr_3-8.0*tfr_2+tfr_1)/(12.0*dist)
140  ! write (*,*) eta,dtfr_dr
141  ! stop
142  !-----------------------------------------------------------------------------
143 
144  df_dr = (-fres4+8.0*fres3-8.0*fres2+fres1)/(12.0*dist)
145  df_drdr = (-fres4+16.0*fres3-3.e1*fres5+16.0*fres2-fres1) &
146  /(12.0*(dist**2 ))
147 
148  pges = (-fres4+8.0*fres3-8.0*fres2+fres1) &
149  /(12.0*dist) *dzetdv*(kbol*t)/1.e-30
150 
151  pgesdz = (-fres4+16.0*fres3-3.e1*fres5+16.0*fres2-fres1) &
152  /(12.0*(dist**2 ))* dzetdv*(kbol*t)/1.e-30 &
153  + (-fres4+8.0*fres3-8.0*fres2+fres1) /(12.0*dist) * 2.0 *rho &
154  *(kbol*t)/1.e-30
155 
156  pgesd2 = (fres4-2.0*fres3+2.0*fres2-fres1) /(2.0*dist**3 ) &
157  * dzetdv*(kbol*t)/1.e-30 &
158  + (-fres4+16.0*fres3-3.e1*fres5+16.0*fres2-fres1) /(12.0*(dist**2 )) &
159  * 4.0 *rho *(kbol*t)/1.e-30 + (-fres4+8.0*fres3-8.0*fres2+fres1) &
160  /(12.0*dist) * 2.0 /z3t *(kbol*t)/1.e-30
161  pgesd3 = (fres4-4.0*fres3+6.0*fres5-4.0*fres2+fres1) /(dist**4) &
162  * dzetdv*(kbol*t)/1.e-30 + (fres4-2.0*fres3+2.0*fres2-fres1) &
163  /(2.0*dist**3 ) * 6.0 *rho *(kbol*t)/1.e-30 &
164  + (-fres4+16.0*fres3-3.e1*fres5+16.0*fres2-fres1) &
165  /(12.0*dist**2 )* 6.0 /z3t *(kbol*t)/1.e-30
166 
167 
168 END SUBROUTINE pressure_pt1
WWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWWW This module contains paramete...
Definition: modules.f90:120