6 |
|
|
7 |
private llm |
private llm |
8 |
|
|
9 |
real ap(llm+1) ! in Pa |
real ap(llm+1), pa ! in Pa |
10 |
real bp(llm+1), dpres(llm) |
real bp(llm+1) |
11 |
real presnivs(llm), pa ! in Pa |
real presnivs(llm) ! pressions approximatives des milieux de couches, en Pa |
12 |
real, parameter:: preff = 101325. ! in Pa |
real, parameter:: preff = 101325. ! in Pa |
13 |
real nivsigs(llm), nivsig(llm+1) |
real nivsigs(llm), nivsig(llm+1) |
14 |
|
|
22 |
! Auteur : P. Le Van |
! Auteur : P. Le Van |
23 |
|
|
24 |
! This procedure sets the vertical grid. |
! This procedure sets the vertical grid. |
25 |
! It defines the host variables "ap", "bp", "dpres", "presnivs", |
! It defines the host variables "ap", "bp", "presnivs", "nivsigs" |
26 |
! "nivsigs" and "nivsig". |
! and "nivsig". |
27 |
! "pa" should be defined before this procedure is called. |
! "pa" should be defined before this procedure is called. |
28 |
|
|
29 |
use comconst, only: pi |
use nr_util, only: pi |
30 |
|
use jumble, only: new_unit |
31 |
|
use unit_nml_m, only: unit_nml |
32 |
|
|
33 |
REAL s(llm+1) |
REAL s(llm+1) |
34 |
! (atmospheric hybrid sigma-pressure coordinate at the interface |
! "s(l)" is the atmospheric hybrid sigma-pressure coordinate at |
35 |
! between layers "l" and "l-1") |
! the interface between layers "l" and "l-1" |
36 |
|
|
37 |
real ds(llm) |
real ds(llm) |
38 |
! "ds(l)" : épaisseur de la couche "l" dans la coordonnée "s" |
! "ds(l)" : épaisseur de la couche "l" dans la coordonnée "s" |
39 |
|
|
40 |
INTEGER l |
INTEGER l, unit |
41 |
REAL alpha, x(llm) |
REAL alpha, x(llm) |
42 |
|
|
43 |
character(len=7):: s_sampling = "LMD5" |
character(len=7):: s_sampling = "tropo" |
44 |
! (other allowed values are "param", "strato1" and "strato2") |
! (other allowed values are "param", "strato1", "strato" and "read") |
45 |
|
|
46 |
real:: h = 7. ! scale height, in km |
real:: h = 7. ! scale height, in km |
47 |
! (used only if "s_sampling" == "param" or "strato1") |
! (used only if "s_sampling" == "param" or "strato1") |
65 |
|
|
66 |
! Compute "s": |
! Compute "s": |
67 |
|
|
|
s(1) = 1. |
|
|
s(llm+1) = 0. |
|
|
|
|
68 |
print *, "Enter namelist 'disvert_nml'." |
print *, "Enter namelist 'disvert_nml'." |
69 |
read(unit=*, nml=disvert_nml) |
read(unit=*, nml=disvert_nml) |
70 |
write(unit=*, nml=disvert_nml) |
write(unit_nml, nml=disvert_nml) |
71 |
|
|
72 |
select case (s_sampling) |
select case (s_sampling) |
73 |
case ("param") |
case ("param") |
74 |
|
s(1) = 1. |
75 |
|
s(llm+1) = 0. |
76 |
alpha = deltaz / tanh(1./k0) * 2. |
alpha = deltaz / tanh(1./k0) * 2. |
77 |
forall (l = 2: llm) s(l) & |
forall (l = 2: llm) s(l) & |
78 |
= cosh((l - 1) / k0) **(- alpha * k0 / h) & |
= cosh((l - 1) / k0) **(- alpha * k0 / h) & |
79 |
* exp(- alpha / h * tanh((llm - k1) / k0) & |
* exp(- alpha / h * tanh((llm - k1) / k0) & |
80 |
* beta **(l - 1 - (llm - k1)) / log(beta)) |
* beta **(l - 1 - (llm - k1)) / log(beta)) |
81 |
case ("LMD5") |
case ("tropo") |
82 |
! Ancienne discrétisation |
s(1) = 1. |
83 |
|
s(llm+1) = 0. |
84 |
forall (l = 1: llm) ds(l) & |
forall (l = 1: llm) ds(l) & |
85 |
= 1. + 7. * SIN(pi * (REAL(l)-0.5) / real(llm+1))**2 |
= 1. + 7. * SIN(pi * (REAL(l)-0.5) / real(llm+1))**2 |
86 |
ds = ds / sum(ds) |
ds = ds / sum(ds) |
90 |
ENDDO |
ENDDO |
91 |
case ("strato1") |
case ("strato1") |
92 |
! F. Lott 70 niveaux et plus |
! F. Lott 70 niveaux et plus |
93 |
|
s(1) = 1. |
94 |
|
s(llm+1) = 0. |
95 |
forall (l = 1: llm) dz(l) = 1.56 + TANH(REAL(l - 12) / 5.) & |
forall (l = 1: llm) dz(l) = 1.56 + TANH(REAL(l - 12) / 5.) & |
96 |
+ TANH(REAL(l - llm) / 10.) / 2. |
+ TANH(REAL(l - llm) / 10.) / 2. |
97 |
|
|
102 |
|
|
103 |
s(2:llm) = (exp(- zz(2:llm) / h) - exp(- zz(llm + 1) / h)) & |
s(2:llm) = (exp(- zz(2:llm) / h) - exp(- zz(llm + 1) / h)) & |
104 |
/ (1. - exp(- zz(llm + 1) / h)) |
/ (1. - exp(- zz(llm + 1) / h)) |
105 |
case ("strato2") |
case ("strato") |
106 |
! Recommended by F. Lott for a domain including the stratosphere |
! Recommended by F. Lott for a domain including the stratosphere |
107 |
|
s(1) = 1. |
108 |
|
s(llm+1) = 0. |
109 |
forall (l = 1: llm) x(l) = pi * (l - 0.5) / (llm + 1) |
forall (l = 1: llm) x(l) = pi * (l - 0.5) / (llm + 1) |
110 |
|
|
111 |
ds = (1. + 7. * SIN(x)**2) * (1. - tanh(2 * x / pi - 1.))**2 / 4. |
ds = (1. + 7. * SIN(x)**2) * (1. - tanh(2 * x / pi - 1.))**2 / 4. |
114 |
DO l = llm, 2, -1 |
DO l = llm, 2, -1 |
115 |
s(l) = s(l+1) + ds(l) |
s(l) = s(l+1) + ds(l) |
116 |
ENDDO |
ENDDO |
117 |
|
case("read") |
118 |
|
call new_unit(unit) |
119 |
|
open(unit, file="hybrid.csv", status="old", action="read", & |
120 |
|
position="rewind") |
121 |
|
read(unit, fmt=*) ! skip title line |
122 |
|
do l = 1, llm + 1 |
123 |
|
read(unit, fmt=*) s(l) |
124 |
|
end do |
125 |
|
close(unit) |
126 |
|
! Quick check: |
127 |
|
if (s(1) /= 1. .or. s(llm + 1) /= 0. .or. s(1) <= s(2)) then |
128 |
|
print *, '"s" should be in descending order, from 1 to 0.' |
129 |
|
stop 1 |
130 |
|
end if |
131 |
case default |
case default |
132 |
print *, 'Wrong value for "s_sampling"' |
print *, 'Wrong value for "s_sampling"' |
133 |
stop 1 |
stop 1 |
138 |
bp(llm + 1) = 0. |
bp(llm + 1) = 0. |
139 |
ap = pa * (s - bp) |
ap = pa * (s - bp) |
140 |
|
|
141 |
forall (l = 1: llm) |
forall (l = 1: llm) presnivs(l) = 0.5 & |
142 |
dpres(l) = bp(l) - bp(l+1) |
* (ap(l) + bp(l) * preff + ap(l+1) + bp(l+1) * preff) |
|
presnivs(l) = 0.5 * (ap(l) + bp(l) * preff + ap(l+1) + bp(l+1) * preff) |
|
|
end forall |
|
143 |
|
|
144 |
END SUBROUTINE disvert |
END SUBROUTINE disvert |
145 |
|
|