source: trunk/WRF.COMMON/WRFV2/dyn_em/module_init_utilities.F @ 3553

Last change on this file since 3553 was 2088, checked in by mlefevre, 6 years ago

MESOSCALE. Argurment of lmd_driver was called twice and a interpolation function was missing
File size: 8.0 KB
Line 
1MODULE module_init_utilities
2
3CONTAINS
4
5 real function interp_0( v_in,  &
6                         z_in, z_out, nz_in  )
7 implicit none
8 integer nz_in, nz_out
9 real    v_in(nz_in), z_in(nz_in)
10 real    z_out
11
12 integer kp, k, im, ip
13 logical interp, increasing_z
14 real    height, w1, w2
15 logical debug
16 parameter ( debug = .false. )
17
18! does vertical coordinate increase or decrease with increasing k?
19! set offset appropriately
20
21 height = z_out
22
23 if(debug) write(6,*) ' height in interp_0 ',height
24
25 if (z_in(nz_in) .gt. z_in(1)) then
26
27    if(debug) write(6,*) ' monotonic increase in z in interp_0 '
28    IF (height > z_in(nz_in)) then
29      if(debug) write(6,*) ' point 1 in interp_0 '
30      w2 = (z_in(nz_in)-height)/(z_in(nz_in)-z_in(nz_in-1))
31      w1 = 1.-w2
32      interp_0 = w1*v_in(nz_in) + w2*v_in(nz_in-1)
33    ELSE IF (height < z_in(1)) then
34      if(debug) write(6,*) ' point 2 in interp_0 '
35      w2 = (z_in(2)-height)/(z_in(2)-z_in(1))
36      w1 = 1.-w2
37      interp_0 = w1*v_in(2) + w2*v_in(1)
38    ELSE
39      if(debug) write(6,*) ' point 3 in interp_0 '
40      interp = .false.
41      kp = nz_in
42      DO WHILE ( (interp .eqv. .false.) .and. (kp .ge. 2) )
43        IF(   ((z_in(kp)   .ge. height) .and.     &
44               (z_in(kp-1) .le. height))        )   THEN
45          w2 = (height-z_in(kp))/(z_in(kp-1)-z_in(kp))
46          w1 = 1.-w2
47          interp_0 = w1*v_in(kp) + w2*v_in(kp-1)
48          if(debug) write(6,*) ' interp data, kp, w1, w2 ',kp, w1, w2
49          if(debug) write(6,*) ' interp data, v_in(kp), v_in(kp-1), interp_0 ', &
50                     v_in(kp), v_in(kp-1), interp_0
51          interp = .true.
52        END IF
53        kp = kp-1
54      ENDDO
55    ENDIF
56
57 else
58
59    if(debug) write(6,*) ' monotonic decrease in z in interp_0 '
60
61    IF (height < z_in(nz_in)) then
62      if(debug) write(6,*) ' point 1 in interp_0 '
63      w2 = (z_in(nz_in)-height)/(z_in(nz_in)-z_in(nz_in-1))
64      w1 = 1.-w2
65      interp_0 = w1*v_in(nz_in) + w2*v_in(nz_in-1)
66    ELSE IF (height > z_in(1)) then
67      if(debug) write(6,*) ' point 2 in interp_0 '
68      w2 = (z_in(2)-height)/(z_in(2)-z_in(1))
69      w1 = 1.-w2
70      interp_0 = w1*v_in(2) + w2*v_in(1)
71    ELSE
72      if(debug) write(6,*) ' point 3 in interp_0 '
73      interp = .false.
74      kp = nz_in
75      height = z_out
76      DO WHILE ( (interp .eqv. .false.) .and. (kp .ge. 2) )
77        IF(   ((z_in(kp)   .le. height) .and.     &
78               (z_in(kp-1) .ge. height))             )   THEN
79          w2 = (height-z_in(kp))/(z_in(kp-1)-z_in(kp))
80          w1 = 1.-w2
81          interp_0 = w1*v_in(kp) + w2*v_in(kp-1)
82          interp = .true.
83        END IF
84        kp = kp-1
85      ENDDO
86    ENDIF
87
88 end if
89
90 return
91 END FUNCTION interp_0
92
93 real function interp_0_log( v_in,  &
94                         p_in, p_out, nz_in  )
95 implicit none
96 integer nz_in, nz_out
97 real    v_in(nz_in), p_in(nz_in)
98 real    p_out
99
100 integer kp, k, im, ip
101 logical interp, increasing_z
102 real    height, w1, w2
103 logical debug
104 parameter ( debug = .false. )
105
106! does vertical coordinate increase or decrease with increasing k?
107! set offset appropriately
108
109 height = p_out
110
111 if(debug) write(6,*) ' height in interp_0 ',height
112
113 if (p_in(nz_in) .gt. p_in(1)) then
114
115    if(debug) write(6,*) ' monotonic increase in z in interp_0 '
116    IF (height > p_in(nz_in)) then
117      if(debug) write(6,*) ' point 1 in interp_0 '
118      w2 = log(p_in(nz_in)/height)/log(p_in(nz_in)/p_in(nz_in-1))
119      w1 = 1.-w2
120      interp_0_log = w1*v_in(nz_in) + w2*v_in(nz_in-1)
121    ELSE IF (height < p_in(1)) then
122      if(debug) write(6,*) ' point 2 in interp_0 '
123      w2 = log(p_in(2)/height)/log(p_in(2)/p_in(1))
124      w1 = 1.-w2
125      interp_0_log = w1*v_in(2) + w2*v_in(1)
126    ELSE
127      if(debug) write(6,*) ' point 3 in interp_0 '
128      interp = .false.
129      kp = nz_in
130      DO WHILE ( (interp .eqv. .false.) .and. (kp .ge. 2) )
131        IF(   ((p_in(kp)   .ge. height) .and.     &
132               (p_in(kp-1) .le. height))        )   THEN
133          w2 = log(height/p_in(kp))/log(p_in(kp-1)/p_in(kp))
134          w1 = 1.-w2
135          interp_0_log = w1*v_in(kp) + w2*v_in(kp-1)
136          if(debug) write(6,*) ' interp data, kp, w1, w2 ',kp, w1, w2
137          if(debug) write(6,*) ' interp data, v_in(kp), v_in(kp-1), interp_0_p ', &
138                     v_in(kp), v_in(kp-1), interp_0_log
139          interp = .true.
140        END IF
141        kp = kp-1
142      ENDDO
143    ENDIF
144
145 else
146
147    if(debug) write(6,*) ' monotonic decrease in z in interp_0 '
148
149    IF (height < p_in(nz_in)) then
150      if(debug) write(6,*) ' point 1 in interp_0 '
151      w2 = log(p_in(nz_in)/height)/log(p_in(nz_in)/p_in(nz_in-1))
152      w1 = 1.-w2
153      interp_0_log = w1*v_in(nz_in) + w2*v_in(nz_in-1)
154    ELSE IF (height > p_in(1)) then
155      if(debug) write(6,*) ' point 2 in interp_0 '
156      w2 = log(p_in(2)/height)/log(p_in(2)/p_in(1))
157      w1 = 1.-w2
158      interp_0_log = w1*v_in(2) + w2*v_in(1)
159    ELSE
160      if(debug) write(6,*) ' point 3 in interp_0 '
161      interp = .false.
162      kp = nz_in
163      height = p_out
164      DO WHILE ( (interp .eqv. .false.) .and. (kp .ge. 2) )
165        IF(   ((p_in(kp)   .le. height) .and.     &
166               (p_in(kp-1) .ge. height))             )   THEN
167          w2 = log(height/p_in(kp))/log(p_in(kp-1)/p_in(kp))
168          w1 = 1.-w2
169          interp_0_log = w1*v_in(kp) + w2*v_in(kp-1)
170          interp = .true.
171        END IF
172        kp = kp-1
173      ENDDO
174    ENDIF
175
176 end if
177
178 return
179 END FUNCTION interp_0_log
180
181 real function interp_0_log2( v_in,  &
182                         p_in, p_out, nz_in  )
183 implicit none
184 integer nz_in, nz_out
185 real*8  v_in(nz_in)
186 real*8  p_in(nz_in)
187 real*8  p_out
188
189 integer kp, k, im, ip
190 logical interp, increasing_z
191 real    height, w1, w2
192 logical debug
193 parameter ( debug = .false. )
194
195! does vertical coordinate increase or decrease with increasing k?
196! set offset appropriately
197
198 height = p_out
199
200 if(debug) write(6,*) ' height in interp_0 ',height
201
202 if (p_in(nz_in) .gt. p_in(1)) then
203
204    if(debug) write(6,*) ' monotonic increase in z in interp_0 '
205    IF (height > p_in(nz_in)) then
206      if(debug) write(6,*) ' point 1 in interp_0 '
207      w2 = log(p_in(nz_in)/height)/log(p_in(nz_in)/p_in(nz_in-1))
208      w1 = 1.-w2
209      interp_0_log2 = w1*v_in(nz_in) + w2*v_in(nz_in-1)
210    ELSE IF (height < p_in(1)) then
211      if(debug) write(6,*) ' point 2 in interp_0 '
212      w2 = log(p_in(2)/height)/log(p_in(2)/p_in(1))
213      w1 = 1.-w2
214      interp_0_log2 = w1*v_in(2) + w2*v_in(1)
215    ELSE
216      if(debug) write(6,*) ' point 3 in interp_0 '
217      interp = .false.
218      kp = nz_in
219      DO WHILE ( (interp .eqv. .false.) .and. (kp .ge. 2) )
220        IF(   ((p_in(kp)   .ge. height) .and.     &
221               (p_in(kp-1) .le. height))        )   THEN
222          w2 = log(height/p_in(kp))/log(p_in(kp-1)/p_in(kp))
223          w1 = 1.-w2
224          interp_0_log2 = w1*v_in(kp) + w2*v_in(kp-1)
225          if(debug) write(6,*) ' interp data, kp, w1, w2 ',kp, w1, w2
226          if(debug) write(6,*) ' interp data, v_in(kp), v_in(kp-1), interp_0_p', &
227                     v_in(kp), v_in(kp-1), interp_0_log2
228          interp = .true.
229        END IF
230        kp = kp-1
231      ENDDO
232    ENDIF
233
234 else
235
236    if(debug) write(6,*) ' monotonic decrease in z in interp_0 '
237
238    IF (height < p_in(nz_in)) then
239      if(debug) write(6,*) ' point 1 in interp_0 '
240      w2 = log(p_in(nz_in)/height)/log(p_in(nz_in)/p_in(nz_in-1))
241      w1 = 1.-w2
242      interp_0_log2 = w1*v_in(nz_in) + w2*v_in(nz_in-1)
243    ELSE IF (height > p_in(1)) then
244      if(debug) write(6,*) ' point 2 in interp_0 '
245      w2 = log(p_in(2)/height)/log(p_in(2)/p_in(1))
246      w1 = 1.-w2
247      interp_0_log2 = w1*v_in(2) + w2*v_in(1)
248    ELSE
249      if(debug) write(6,*) ' point 3 in interp_0 '
250      interp = .false.
251      kp = nz_in
252      height = p_out
253      DO WHILE ( (interp .eqv. .false.) .and. (kp .ge. 2) )
254        IF(   ((p_in(kp)   .le. height) .and.     &
255               (p_in(kp-1) .ge. height))             )   THEN
256          w2 = log(height/p_in(kp))/log(p_in(kp-1)/p_in(kp))
257          w1 = 1.-w2
258          interp_0_log2 = w1*v_in(kp) + w2*v_in(kp-1)
259          interp = .true.
260        END IF
261        kp = kp-1
262      ENDDO
263    ENDIF
264
265 end if
266
267 return
268 END FUNCTION interp_0_log2
269
270END MODULE module_init_utilities
271
272
Note: See TracBrowser for help on using the repository browser.