Line data Source code
1 : !------------------------------------------------------------------------
2 : ! $Id$
3 : !===============================================================================
4 : module grid_class
5 :
6 : ! Description:
7 : !
8 : ! Definition of a grid class and associated functions
9 : !
10 : ! The grid specification is as follows:
11 : !
12 : ! + ================== zm(nz) =========GP=========
13 : ! |
14 : ! |
15 : ! 1/dzt(nz) + ------------------ zt(nz) ---------GP---------
16 : ! | |
17 : ! | |
18 : ! + 1/dzm(nz-1) ================== zm(nz-1) ==================
19 : ! |
20 : ! |
21 : ! + ------------------ zt(nz-1) ------------------
22 : !
23 : ! .
24 : ! .
25 : ! .
26 : ! .
27 : !
28 : ! ================== zm(k+1) ===================
29 : !
30 : !
31 : ! + ------------------ zt(k+1) -------------------
32 : ! |
33 : ! |
34 : ! + 1/dzm(k) ================== zm(k) =====================
35 : ! | |
36 : ! | |
37 : ! 1/dzt(k) + ------------------ zt(k) ---------------------
38 : ! |
39 : ! |
40 : ! + ================== zm(k-1) ===================
41 : !
42 : !
43 : ! ------------------ zt(k-1) -------------------
44 : !
45 : ! .
46 : ! .
47 : ! .
48 : ! .
49 : !
50 : ! + ================== zm(2) =====================
51 : ! |
52 : ! |
53 : ! 1/dzt(2) + ------------------ zt(2) ---------------------
54 : ! | |
55 : ! | |
56 : ! + 1/dzm(1) ================== zm(1) ============GP======= zm_init
57 : ! | ////////////////////////////////////////////// surface
58 : ! |
59 : ! + ------------------ zt(1) ------------GP-------
60 : !
61 : !
62 : ! The variable zm(k) stands for the momentum level altitude at momentum
63 : ! level k; the variable zt(k) stands for the thermodynamic level altitude at
64 : ! thermodynamic level k; the variable invrs_dzt(k) is the inverse distance
65 : ! between momentum levels (over a central thermodynamic level k); and the
66 : ! variable invrs_dzm(k) is the inverse distance between thermodynamic levels
67 : ! (over a central momentum level k). Please note that in the above diagram,
68 : ! "invrs_dzt" is denoted "dzt", and "invrs_dzm" is denoted "dzm", such that
69 : ! 1/dzt is the distance between successive momentum levels k-1 and k (over a
70 : ! central thermodynamic level k), and 1/dzm is the distance between successive
71 : ! thermodynamic levels k and k+1 (over a central momentum level k).
72 : !
73 : ! The grid setup is compatible with a stretched (unevely-spaced) grid. Thus,
74 : ! the distance between successive grid levels may not always be constant.
75 : !
76 : ! The following diagram is an example of a stretched grid that is defined on
77 : ! momentum levels. The thermodynamic levels are placed exactly halfway
78 : ! between the momentum levels. However, the momentum levels do not fall
79 : ! halfway between the thermodynamic levels.
80 : !
81 : ! =============== zm(k+1) ===============
82 : !
83 : !
84 : !
85 : ! --------------- zt(k+1) ---------------
86 : !
87 : !
88 : !
89 : ! =============== zm(k) ===============
90 : !
91 : ! --------------- zt(k) ---------------
92 : !
93 : ! =============== zm(k-1) ===============
94 : !
95 : ! The following diagram is an example of a stretched grid that is defined on
96 : ! thermodynamic levels. The momentum levels are placed exactly halfway
97 : ! between the thermodynamic levels. However, the thermodynamic levels do not
98 : ! fall halfway between the momentum levels.
99 : !
100 : ! --------------- zt(k+1) ---------------
101 : !
102 : !
103 : !
104 : ! =============== zm(k) ===============
105 : !
106 : !
107 : !
108 : ! --------------- zt(k) ---------------
109 : !
110 : ! =============== zm(k-1) ===============
111 : !
112 : ! --------------- zt(k-1) ---------------
113 : !
114 : ! NOTE: Any future code written for use in the CLUBB parameterization should
115 : ! use interpolation formulas consistent with a stretched grid. The
116 : ! simplest way to do so is to call the appropriate interpolation
117 : ! function from this module. Interpolations should *not* be handled in
118 : ! the form of: ( var_zm(k) + var_zm(k-1) ) / 2; *nor* in the form of:
119 : ! 0.5_core_rknd*( var_zt(k+1) + var_zt(k) ). Rather, all explicit interpolations
120 : ! should call zt2zm or zm2zt; while interpolations for a variable being
121 : ! solved for implicitly in the code should use gr%weights_zt2zm (which
122 : ! refers to interp_weights_zt2zm_imp), or gr%weights_zm2zt (which
123 : ! refers to interp_weights_zm2zt_imp).
124 : !
125 : ! Momentum level 1 is placed at altitude zm_init, which is usually at the
126 : ! surface. However, in general, zm_init can be at any altitude defined by the
127 : ! user.
128 : !
129 : ! GP indicates ghost points. Variables located at those levels are not
130 : ! prognosed, but only used for boundary conditions.
131 : !
132 : ! Chris Golaz, 7/17/99
133 : ! modified 9/10/99
134 : ! schemena, modified 6/11/2014 - Restructered code to add cubic/linear flag
135 :
136 : ! References:
137 : !
138 : ! https://arxiv.org/pdf/1711.03675v1.pdf#nameddest=url:clubb_grid
139 : !
140 : ! Section 3c, p. 3548 /Numerical discretization/ of:
141 : ! ``A PDF-Based Model for Boundary Layer Clouds. Part I:
142 : ! Method and Model Description'' Golaz, et al. (2002)
143 : ! JAS, Vol. 59, pp. 3540--3551.
144 : !-----------------------------------------------------------------------
145 :
146 : use clubb_precision, only: &
147 : core_rknd ! Variable(s)
148 :
149 : implicit none
150 :
151 : public :: grid, zt2zm, zm2zt, zt2zm2zt, zm2zt2zm, &
152 : ddzm, ddzt, &
153 : setup_grid, cleanup_grid, setup_grid_heights, &
154 : read_grid_heights, flip
155 :
156 : private :: t_above, t_below, m_above, m_below
157 :
158 : private ! Default Scoping
159 :
160 : ! Constant parameters
161 : integer, parameter :: &
162 : t_above = 1, & ! Upper thermodynamic level index (gr%weights_zt2zm).
163 : t_below = 2, & ! Lower thermodynamic level index (gr%weights_zt2zm).
164 : m_above = 1, & ! Upper momentum level index (gr%weights_zm2zt).
165 : m_below = 2 ! Lower momentum level index (gr%weights_zm2zt).
166 :
167 : type grid
168 :
169 : integer :: nz ! Number of points in the grid
170 : ! Note: Fortran 90/95 prevents an allocatable array from appearing
171 : ! within a derived type. However, Fortran 2003 does not!!!!!!!!
172 : real( kind = core_rknd ), allocatable, dimension(:,:) :: &
173 : zm, & ! Momentum grid
174 : zt ! Thermo grid
175 :
176 : real( kind = core_rknd ), allocatable, dimension(:,:) :: &
177 : invrs_dzm, & ! The inverse spacing between thermodynamic grid
178 : ! levels; centered over momentum grid levels.
179 : invrs_dzt ! The inverse spacing between momentum grid levels;
180 : ! centered over thermodynamic grid levels.
181 :
182 : real( kind = core_rknd ), allocatable, dimension(:,:) :: &
183 : dzm, & ! Spacing between thermodynamic grid levels; centered over
184 : ! momentum grid levels
185 : dzt ! Spcaing between momentum grid levels; centered over
186 : ! thermodynamic grid levels
187 :
188 : ! These weights are normally used in situations
189 : ! where a momentum level variable is being
190 : ! solved for implicitly in an equation and
191 : ! needs to be interpolated to the thermodynamic grid levels.
192 : real( kind = core_rknd ), allocatable, dimension(:,:,:) :: weights_zm2zt, &
193 : ! These weights are normally used in situations where a
194 : ! thermodynamic level variable is being solved for implicitly in an equation
195 : ! and needs to be interpolated to the momentum grid levels.
196 : weights_zt2zm
197 :
198 : end type grid
199 :
200 : ! The grid is defined here so that it is common throughout the module.
201 : ! The implication is that only one grid can be defined !
202 :
203 :
204 : ! Modification for using CLUBB in a host model (i.e. one grid per column)
205 :
206 :
207 : ! Interfaces provided for function overloading
208 :
209 : interface zt2zm
210 : ! For l_cubic_interp = .true.
211 : ! This version uses cublic spline interpolation of Stefen (1990).
212 : !
213 : ! For l_cubic_interp = .false.
214 : ! This performs a linear extension at the highest grid level and therefore
215 : ! does not guarantee, for positive definite quantities (e.g. wp2), that the
216 : ! extended point is indeed positive definite. Positive definiteness can be
217 : ! ensured with a max statement.
218 : ! In the future, we could add a flag (lposdef) and, when needed, apply the
219 : ! max statement directly within interpolated_azm and interpolated_azmk.
220 : module procedure redirect_interpolated_azm_k ! Works over a single vertical level
221 : module procedure redirect_interpolated_azm_1D ! Works over all vertical levels
222 : module procedure redirect_interpolated_azm_2D ! Works over all vertical levels and columns
223 : end interface
224 :
225 : interface zm2zt
226 : ! For l_cubic_interp = .true.
227 : ! This version uses cublic spline interpolation of Stefen (1990).
228 : !
229 : ! For l_cubic_interp = .false.
230 : ! This performs a linear extension at the lowest grid level and therefore
231 : ! does not guarantee, for positive definite quantities (e.g. wp2), that the
232 : ! extended point is indeed positive definite. Positive definiteness can be
233 : ! ensured with a max statement.
234 : ! In the future, we could add a flag (lposdef) and, when needed, apply the
235 : ! max statement directly within interpolated_azt and interpolated_aztk.
236 : module procedure redirect_interpolated_azt_k ! Works over a single vertical level
237 : module procedure redirect_interpolated_azt_1D ! Works over all vertical levels
238 : module procedure redirect_interpolated_azt_2D ! Works over all vertical levels and columns
239 : end interface
240 :
241 : ! Vertical derivative functions
242 : interface ddzm
243 : module procedure gradzm_1D ! Works over all vertical levels
244 : module procedure gradzm_2D ! Works over all vertical levels and columns
245 : end interface
246 :
247 : interface ddzt
248 : module procedure gradzt_1D ! Works over all vertical levels
249 : module procedure gradzt_2D ! Works over all vertical levels and columns
250 : end interface
251 :
252 : contains
253 :
254 : !=============================================================================
255 4467528 : subroutine setup_grid( nzmax, ngrdcol, sfc_elevation, l_implemented, &
256 4467528 : grid_type, deltaz, zm_init, zm_top, &
257 4467528 : momentum_heights, thermodynamic_heights, &
258 : gr )
259 :
260 : ! Description:
261 : ! Grid Constructor
262 : !
263 : ! This subroutine sets up the CLUBB vertical grid.
264 : !
265 : ! References:
266 : ! ``Equations for CLUBB'', Sec. 8, Grid Configuration.
267 : !-----------------------------------------------------------------------
268 :
269 : use constants_clubb, only: &
270 : fstderr ! Variable(s)
271 :
272 : use error_code, only: &
273 : clubb_at_least_debug_level, & ! Procedure
274 : err_code, & ! Error indicator
275 : clubb_fatal_error, & ! Constant
276 : err_header ! String
277 :
278 : use clubb_precision, only: &
279 : core_rknd ! Variable(s)
280 :
281 : implicit none
282 :
283 : ! Constant parameters
284 : integer, parameter :: &
285 : NWARNING = 250 ! Issue a warning if nzmax exceeds this number.
286 :
287 : ! Input Variables
288 : integer, intent(in) :: &
289 : nzmax, & ! Number of vertical levels in grid [#]
290 : ngrdcol
291 :
292 : type(grid), target, intent(inout) :: gr
293 :
294 : real( kind = core_rknd ), dimension(ngrdcol), intent(in) :: &
295 : sfc_elevation ! Elevation of ground level [m AMSL]
296 :
297 : ! Flag to see if CLUBB is running on it's own,
298 : ! or if it's implemented as part of a host model.
299 : logical, intent(in) :: l_implemented
300 :
301 : ! If CLUBB is running on it's own, this option determines if it is using:
302 : ! 1) an evenly-spaced grid;
303 : ! 2) a stretched (unevenly-spaced) grid entered on the thermodynamic grid
304 : ! levels (with momentum levels set halfway between thermodynamic levels);
305 : ! or
306 : ! 3) a stretched (unevenly-spaced) grid entered on the momentum grid levels
307 : ! (with thermodynamic levels set halfway between momentum levels).
308 : integer, intent(in) :: grid_type
309 :
310 : ! If the CLUBB model is running by itself, and is using an evenly-spaced
311 : ! grid (grid_type = 1), it needs the vertical grid spacing and
312 : ! momentum-level starting altitude as input.
313 : real( kind = core_rknd ), dimension(ngrdcol), intent(in) :: &
314 : deltaz, & ! Vertical grid spacing [m]
315 : zm_init, & ! Initial grid altitude (momentum level) [m]
316 : zm_top ! Maximum grid altitude (momentum level) [m]
317 :
318 : ! If the CLUBB parameterization is implemented in a host model, it needs to
319 : ! use the host model's momentum level altitudes and thermodynamic level
320 : ! altitudes.
321 : ! If the CLUBB model is running by itself, but is using a stretched grid
322 : ! entered on thermodynamic levels (grid_type = 2), it needs to use the
323 : ! thermodynamic level altitudes as input.
324 : ! If the CLUBB model is running by itself, but is using a stretched grid
325 : ! entered on momentum levels (grid_type = 3), it needs to use the momentum
326 : ! level altitudes as input.
327 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nzmax) :: &
328 : momentum_heights, & ! Momentum level altitudes (input) [m]
329 : thermodynamic_heights ! Thermodynamic level altitudes (input) [m]
330 :
331 : ! Local Variables
332 : integer :: &
333 : begin_height, & ! Lower bound for *_heights arrays [-]
334 : end_height ! Upper bound for *_heights arrays [-]
335 :
336 : integer :: ierr, & ! Allocation stat
337 : i, k ! Loop index
338 :
339 :
340 : ! ---- Begin Code ----
341 :
342 : ! Define the grid size
343 :
344 4467528 : if ( nzmax > NWARNING .and. clubb_at_least_debug_level( 1 ) ) then
345 : write(fstderr,*) "Warning: running with vertical grid "// &
346 0 : "which is larger than", NWARNING, "levels."
347 0 : write(fstderr,*) "This may take a lot of CPU time and memory."
348 : end if
349 :
350 4467528 : gr%nz = nzmax
351 :
352 : ! Default bounds
353 4467528 : begin_height = 1
354 4467528 : end_height = gr%nz
355 :
356 : !---------------------------------------------------
357 4467528 : if ( .not. l_implemented ) then
358 :
359 : ! Since l_implemented=.false., we have to calculate the number of vertical levels
360 : ! depending on the grid type. Currently, clubb can only have multiple columns
361 : ! if implemented in a host model, but if in the future we have a standalone
362 : ! version with multiple columns, we may end up with columns that have a
363 : ! different number of vertical levels. CLUBB cannot handle this however, so
364 : ! we use only the first column to determine the number of vertical levels.
365 :
366 0 : if ( grid_type == 1 ) then
367 :
368 : ! Determine the number of grid points given the spacing
369 : ! to fit within the bounds without going over.
370 0 : gr%nz = floor( ( zm_top(1) - zm_init(1) + deltaz(1) ) / deltaz(1) )
371 :
372 0 : else if( grid_type == 2 ) then! Thermo
373 :
374 : ! Find begin_height (lower bound)
375 : k = gr%nz
376 :
377 0 : do while( thermodynamic_heights(1,k) >= zm_init(1) .and. k > 1 )
378 :
379 0 : k = k - 1
380 :
381 : end do
382 :
383 0 : if( thermodynamic_heights(1,k) >= zm_init(1) ) then
384 :
385 0 : write(fstderr,*) err_header, "Stretched zt grid cannot fulfill zm_init requirement"
386 0 : err_code = clubb_fatal_error
387 0 : return
388 :
389 : else
390 :
391 0 : begin_height = k
392 :
393 : end if
394 :
395 : ! Find end_height (upper bound)
396 : k = gr%nz
397 :
398 0 : do while( thermodynamic_heights(1,k) > zm_top(1) .and. k > 1 )
399 :
400 0 : k = k - 1
401 :
402 : end do
403 :
404 0 : if( zm_top(1) < thermodynamic_heights(1,k) ) then
405 :
406 0 : write(fstderr,*) err_header, "Stretched zt grid cannot fulfill zm_top requirement"
407 0 : err_code = clubb_fatal_error
408 0 : return
409 :
410 : else
411 :
412 0 : end_height = k
413 :
414 0 : gr%nz = size( thermodynamic_heights(1,begin_height:end_height ) )
415 :
416 : end if
417 :
418 0 : else if( grid_type == 3 ) then ! Momentum
419 :
420 : ! Find begin_height (lower bound)
421 : k = 1
422 :
423 0 : do while( momentum_heights(1,k) < zm_init(1) .and. k < gr%nz )
424 :
425 0 : k = k + 1
426 :
427 : end do
428 :
429 0 : if( momentum_heights(1,k) < zm_init(1) ) then
430 :
431 0 : write(fstderr,*) err_header, "Stretched zm grid cannot fulfill zm_init requirement"
432 0 : err_code = clubb_fatal_error
433 0 : return
434 :
435 : else
436 :
437 0 : begin_height = k
438 :
439 : end if
440 :
441 : ! Find end_height (lower bound)
442 : k = gr%nz
443 :
444 0 : do while( momentum_heights(1,k) > zm_top(1) .and. k > 1 )
445 :
446 0 : k = k - 1
447 :
448 : end do
449 :
450 0 : if( momentum_heights(1,k) > zm_top(1) ) then
451 :
452 0 : write(fstderr,*) err_header, "Stretched zm grid cannot fulfill zm_top requirement"
453 0 : err_code = clubb_fatal_error
454 0 : return
455 :
456 : else
457 :
458 0 : end_height = k
459 :
460 0 : gr%nz = size( momentum_heights(1,begin_height:end_height ) )
461 :
462 : end if
463 :
464 : endif ! grid_type
465 :
466 : endif ! .not. l_implemented
467 :
468 : !---------------------------------------------------
469 :
470 : ! Allocate memory for the grid levels
471 : allocate( gr%zm(ngrdcol,gr%nz), gr%zt(ngrdcol,gr%nz), &
472 : gr%dzm(ngrdcol,gr%nz), gr%dzt(ngrdcol,gr%nz), &
473 : gr%invrs_dzm(ngrdcol,gr%nz), gr%invrs_dzt(ngrdcol,gr%nz), &
474 : gr%weights_zm2zt(ngrdcol,gr%nz,m_above:m_below), &
475 : gr%weights_zt2zm(ngrdcol,gr%nz,t_above:t_below), &
476 89350560 : stat=ierr )
477 :
478 4467528 : if ( ierr /= 0 ) then
479 0 : write(fstderr,*) err_header, "In setup_grid: allocation of grid variables failed."
480 0 : err_code = clubb_fatal_error
481 0 : return
482 : end if
483 :
484 : ! Set the values for the derived types used for heights, derivatives, and
485 : ! interpolation from the momentum/thermodynamic grid
486 : call setup_grid_heights( &
487 : gr%nz, ngrdcol, & ! intent(in)
488 : l_implemented, grid_type, & ! intent(in)
489 : deltaz, zm_init, & ! intent(in)
490 0 : momentum_heights(:,begin_height:end_height), & ! intent(in)
491 : thermodynamic_heights(:,begin_height:end_height), & ! intent(in)
492 4467528 : gr ) ! intent(inout)
493 :
494 74597328 : do i = 1, ngrdcol
495 74597328 : if ( sfc_elevation(i) > gr%zm(i,1) ) then
496 : write(fstderr,*) "The altitude of the lowest momentum level, " &
497 : // "gr%zm(1,1), must be at or above the altitude of " &
498 : // "the surface, sfc_elevation. The lowest model " &
499 0 : // "momentum level cannot be below the surface."
500 0 : write(fstderr,*) "Altitude of lowest momentum level =", gr%zm(i,1)
501 0 : write(fstderr,*) "Altitude of the surface =", sfc_elevation(i)
502 0 : err_code = clubb_fatal_error
503 0 : return
504 : endif
505 : end do
506 :
507 : return
508 :
509 : end subroutine setup_grid
510 :
511 : !=============================================================================
512 0 : subroutine cleanup_grid( gr )
513 :
514 : ! Description:
515 : ! De-allocates the memory for the grid
516 : !
517 : ! References:
518 : ! None
519 : !------------------------------------------------------------------------------
520 : use constants_clubb, only: &
521 : fstderr ! Constant(s)
522 :
523 : implicit none
524 :
525 : type(grid), target, intent(inout) :: gr
526 :
527 : ! Local Variable(s)
528 : integer :: ierr
529 :
530 : ! ----- Begin Code -----
531 :
532 : ! Allocate memory for grid levels
533 : deallocate( gr%zm, gr%zt, &
534 : gr%dzm, gr%dzt, &
535 : gr%invrs_dzm, gr%invrs_dzt, &
536 : gr%weights_zm2zt, gr%weights_zt2zm, &
537 0 : stat=ierr )
538 :
539 0 : if ( ierr /= 0 ) then
540 0 : write(fstderr,*) "Grid deallocation failed."
541 : end if
542 :
543 :
544 0 : return
545 :
546 : end subroutine cleanup_grid
547 :
548 : !=============================================================================
549 4467528 : subroutine setup_grid_heights( nz, ngrdcol, &
550 : l_implemented, grid_type, &
551 4467528 : deltaz, zm_init, momentum_heights, &
552 4467528 : thermodynamic_heights, &
553 : gr )
554 :
555 : ! Description:
556 : ! Sets the heights and interpolation weights of the column.
557 : ! This is seperated from setup_grid for those host models that have heights
558 : ! that vary with time.
559 : ! References:
560 : ! None
561 : !------------------------------------------------------------------------------
562 :
563 : use constants_clubb, only: &
564 : fstderr ! Constant(s)
565 :
566 : use clubb_precision, only: &
567 : core_rknd ! Variable(s)
568 :
569 : use error_code, only: &
570 : err_code, & ! Error indicator
571 : clubb_fatal_error, & ! Constant
572 : err_header ! String
573 :
574 : implicit none
575 :
576 : ! Input Variables
577 : integer, intent(in) :: &
578 : nz, &
579 : ngrdcol
580 :
581 : type (grid), target, intent(inout) :: gr
582 :
583 : ! Flag to see if CLUBB is running on it's own,
584 : ! or if it's implemented as part of a host model.
585 : logical, intent(in) :: l_implemented
586 :
587 : ! If CLUBB is running on it's own, this option determines if it is using:
588 : ! 1) an evenly-spaced grid;
589 : ! 2) a stretched (unevenly-spaced) grid entered on the thermodynamic grid
590 : ! levels (with momentum levels set halfway between thermodynamic levels);
591 : ! or
592 : ! 3) a stretched (unevenly-spaced) grid entered on the momentum grid levels
593 : ! (with thermodynamic levels set halfway between momentum levels).
594 : integer, intent(in) :: grid_type
595 :
596 : ! If the CLUBB model is running by itself, and is using an evenly-spaced
597 : ! grid (grid_type = 1), it needs the vertical grid spacing and
598 : ! momentum-level starting altitude as input.
599 : real( kind = core_rknd ), dimension(ngrdcol), intent(in) :: &
600 : deltaz, & ! Vertical grid spacing [m]
601 : zm_init ! Initial grid altitude (momentum level) [m]
602 :
603 :
604 : ! If the CLUBB parameterization is implemented in a host model, it needs to
605 : ! use the host model's momentum level altitudes and thermodynamic level
606 : ! altitudes.
607 : ! If the CLUBB model is running by itself, but is using a stretched grid
608 : ! entered on thermodynamic levels (grid_type = 2), it needs to use the
609 : ! thermodynamic level altitudes as input.
610 : ! If the CLUBB model is running by itself, but is using a stretched grid
611 : ! entered on momentum levels (grid_type = 3), it needs to use the momentum
612 : ! level altitudes as input.
613 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
614 : momentum_heights, & ! Momentum level altitudes (input) [m]
615 : thermodynamic_heights ! Thermodynamic level altitudes (input) [m]
616 :
617 : integer :: k, i
618 :
619 : ! ---- Begin Code ----
620 :
621 4467528 : if ( .not. l_implemented ) then
622 :
623 :
624 0 : if ( grid_type == 1 ) then
625 :
626 : ! Evenly-spaced grid.
627 : ! Momentum level altitudes are defined based on the grid starting
628 : ! altitude, zm_init, the constant grid-spacing, deltaz, and the number
629 : ! of grid levels, gr%nz.
630 :
631 : ! Define momentum level altitudes. The first momentum level is at
632 : ! altitude zm_init.
633 0 : do k = 1, nz, 1
634 0 : do i = 1, ngrdcol
635 0 : gr%zm(i,k) = zm_init(i) + real( k-1, kind = core_rknd ) * deltaz(i)
636 : end do
637 : end do
638 :
639 : ! Define thermodynamic level altitudes. Thermodynamic level altitudes
640 : ! are located at the central altitude levels, exactly halfway between
641 : ! momentum level altitudes. The lowermost thermodynamic level is
642 : ! found by taking 1/2 the altitude difference between the bottom two
643 : ! momentum levels and subtracting that value from the bottom momentum
644 : ! level. The first thermodynamic level is below zm_init.
645 0 : do i = 1, ngrdcol
646 0 : gr%zt(i,1) = zm_init(i) - ( 0.5_core_rknd * deltaz(i) )
647 : end do
648 :
649 0 : do k = 2, nz, 1
650 0 : do i = 1, ngrdcol
651 0 : gr%zt(i,k) = 0.5_core_rknd * ( gr%zm(i,k) + gr%zm(i,k-1) )
652 : end do
653 : enddo
654 :
655 :
656 0 : elseif ( grid_type == 2 ) then
657 :
658 : ! Stretched (unevenly-spaced) grid: stretched thermodynamic levels.
659 : ! Thermodynamic levels are defined according to a stretched grid that
660 : ! is entered through the use of an input file. This is similar to a
661 : ! SAM-style stretched grid.
662 :
663 : ! Define thermodynamic level altitudes.
664 0 : do k = 1, nz, 1
665 0 : do i = 1, ngrdcol
666 0 : gr%zt(i,k) = thermodynamic_heights(i,k)
667 : end do
668 : enddo
669 :
670 : ! Define momentum level altitudes. Momentum level altitudes are
671 : ! located at the central altitude levels, exactly halfway between
672 : ! thermodynamic level altitudes. The uppermost momentum level
673 : ! altitude is found by taking 1/2 the altitude difference between the
674 : ! top two thermodynamic levels and adding that value to the top
675 : ! thermodynamic level.
676 0 : do k = 1, nz-1, 1
677 0 : do i = 1, ngrdcol
678 0 : gr%zm(i,k) = 0.5_core_rknd * ( gr%zt(i,k+1) + gr%zt(i,k) )
679 : end do
680 : enddo
681 :
682 0 : do i = 1, ngrdcol
683 0 : gr%zm(i,nz) = gr%zt(i,nz) + &
684 0 : 0.5_core_rknd * ( gr%zt(i,nz) - gr%zt(i,nz-1) )
685 : end do
686 :
687 0 : elseif ( grid_type == 3 ) then
688 :
689 : ! Stretched (unevenly-spaced) grid: stretched momentum levels.
690 : ! Momentum levels are defined according to a stretched grid that is
691 : ! entered through the use of an input file. This is similar to a
692 : ! WRF-style stretched grid.
693 :
694 : ! Define momentum level altitudes.
695 0 : do k = 1, nz, 1
696 0 : do i = 1, ngrdcol
697 0 : gr%zm(i,k) = momentum_heights(i,k)
698 : end do
699 : enddo
700 :
701 : ! Define thermodynamic level altitudes. Thermodynamic level altitudes
702 : ! are located at the central altitude levels, exactly halfway between
703 : ! momentum level altitudes. The lowermost thermodynamic level
704 : ! altitude is found by taking 1/2 the altitude difference between the
705 : ! bottom two momentum levels and subtracting that value from the
706 : ! bottom momentum level.
707 0 : do i = 1, ngrdcol
708 0 : gr%zt(i,1) = gr%zm(i,1) - 0.5_core_rknd * ( gr%zm(i,2) - gr%zm(i,1) )
709 : end do
710 :
711 0 : do k = 2, nz, 1
712 0 : do i = 1, ngrdcol
713 0 : gr%zt(i,k) = 0.5_core_rknd * ( gr%zm(i,k) + gr%zm(i,k-1) )
714 : end do
715 : enddo
716 :
717 :
718 : else
719 :
720 : ! Invalid grid type.
721 0 : write(fstderr,*) err_header, "Invalid grid type: ", grid_type, &
722 0 : ". Valid options are 1, 2, or 3."
723 0 : err_code = clubb_fatal_error
724 0 : return
725 :
726 :
727 : endif
728 :
729 :
730 : else
731 :
732 : ! The CLUBB parameterization is implemented in a host model.
733 : ! Use the host model's momentum level altitudes and thermodynamic level
734 : ! altitudes to set up the CLUBB grid.
735 :
736 : ! Momentum level altitudes from host model.
737 384207408 : do k = 1, nz, 1
738 6345240408 : do i = 1, ngrdcol
739 6340772880 : gr%zm(i,k) = momentum_heights(i,k)
740 : end do
741 : enddo
742 :
743 : ! Thermodynamic level altitudes from host model after possible grid-index
744 : ! adjustment for CLUBB interface.
745 384207408 : do k = 1, nz, 1
746 6345240408 : do i = 1, ngrdcol
747 6340772880 : gr%zt(i,k) = thermodynamic_heights(i,k)
748 : end do
749 : enddo
750 :
751 :
752 : endif ! not l_implemented
753 :
754 :
755 : ! Define dzm, the spacing between thermodynamic grid levels; centered over
756 : ! momentum grid levels
757 379739880 : do k=1,nz-1
758 6270643080 : do i = 1, ngrdcol
759 6266175552 : gr%dzm(i,k) = gr%zt(i,k+1) - gr%zt(i,k)
760 : end do
761 : enddo
762 :
763 74597328 : do i = 1, ngrdcol
764 74597328 : gr%dzm(i,nz) = gr%dzm(i,nz-1)
765 : end do
766 :
767 : ! Define dzt, the spacing between momentum grid levels; centered over
768 : ! thermodynamic grid levels
769 379739880 : do k=2,nz
770 6270643080 : do i = 1, ngrdcol
771 6266175552 : gr%dzt(i,k) = gr%zm(i,k) - gr%zm(i,k-1)
772 : end do
773 : enddo
774 :
775 74597328 : do i = 1, ngrdcol
776 74597328 : gr%dzt(i,1) = gr%dzt(i,2)
777 : end do
778 :
779 : ! Define invrs_dzm, which is the inverse spacing between thermodynamic grid
780 : ! levels; centered over momentum grid levels.
781 379739880 : do k=1,nz-1
782 6270643080 : do i = 1, ngrdcol
783 6266175552 : gr%invrs_dzm(i,k) = 1._core_rknd / ( gr%zt(i,k+1) - gr%zt(i,k) )
784 : end do
785 : enddo
786 :
787 74597328 : do i = 1, ngrdcol
788 74597328 : gr%invrs_dzm(i,nz) = gr%invrs_dzm(i,nz-1)
789 : end do
790 :
791 :
792 : ! Define invrs_dzt, which is the inverse spacing between momentum grid
793 : ! levels; centered over thermodynamic grid levels.
794 379739880 : do k=2,nz
795 6270643080 : do i = 1, ngrdcol
796 6266175552 : gr%invrs_dzt(i,k) = 1._core_rknd / ( gr%zm(i,k) - gr%zm(i,k-1) )
797 : end do
798 : enddo
799 :
800 74597328 : do i = 1, ngrdcol
801 74597328 : gr%invrs_dzt(i,1) = gr%invrs_dzt(i,2)
802 : end do
803 :
804 :
805 : ! Interpolation Weights: zm grid to zt grid.
806 : ! The grid index (k) is the index of the level on the thermodynamic (zt)
807 : ! grid. The result is the weights of the upper and lower momentum levels
808 : ! (that sandwich the thermodynamic level) applied to that thermodynamic
809 : ! level. These weights are normally used in situations where a momentum
810 : ! level variable is being solved for implicitly in an equation, and the
811 : ! aforementioned variable needs to be interpolated from three successive
812 : ! momentum levels (the central momentum level, as well as one momentum level
813 : ! above and below the central momentum level) to the intermediate
814 : ! thermodynamic grid levels that sandwich the central momentum level.
815 : ! For more information, see the comments in function interpolated_aztk_imp.
816 : call calc_zm2zt_weights( nz, ngrdcol, & ! In
817 4467528 : gr ) ! InOut
818 :
819 :
820 : ! Interpolation Weights: zt grid to zm grid.
821 : ! The grid index (k) is the index of the level on the momentum (zm) grid.
822 : ! The result is the weights of the upper and lower thermodynamic levels
823 : ! (that sandwich the momentum level) applied to that momentum level. These
824 : ! weights are normally used in situations where a thermodynamic level
825 : ! variable is being solved for implicitly in an equation, and the
826 : ! aforementioned variable needs to be interpolated from three successive
827 : ! thermodynamic levels (the central thermodynamic level, as well as one
828 : ! thermodynamic level above and below the central thermodynamic level) to
829 : ! the intermediate momentum grid levels that sandwich the central
830 : ! thermodynamic level.
831 : ! For more information, see the comments in function interpolated_azmk_imp.
832 : call calc_zt2zm_weights( nz, ngrdcol, & ! In
833 4467528 : gr ) ! InOut
834 :
835 4467528 : return
836 : end subroutine setup_grid_heights
837 :
838 : !=============================================================================
839 0 : subroutine read_grid_heights( nzmax, grid_type, &
840 : zm_grid_fname, zt_grid_fname, &
841 : file_unit, &
842 0 : momentum_heights, &
843 0 : thermodynamic_heights )
844 :
845 : ! Description:
846 : ! This subroutine is used foremost in cases where the grid_type corresponds
847 : ! with the stretched (unevenly-spaced) grid options (either grid_type = 2 or
848 : ! grid_type = 3). This subroutine reads in the values of the stretched grid
849 : ! altitude levels for either the thermodynamic level grid or the momentum
850 : ! level grid. This subroutine also handles basic error checking for all
851 : ! three grid types.
852 : !------------------------------------------------------------------------
853 :
854 : use constants_clubb, only: &
855 : fstderr ! Variable(s)
856 :
857 : use file_functions, only: &
858 : file_read_1d ! Procedure(s)
859 :
860 : use clubb_precision, only: &
861 : core_rknd ! Variable(s)
862 :
863 : use error_code, only: &
864 : err_code, & ! Error indicator
865 : clubb_fatal_error, & ! Constant
866 : err_header ! String
867 :
868 : implicit none
869 :
870 : ! Input Variables.
871 :
872 : ! Declared number of vertical levels.
873 : integer, intent(in) :: &
874 : nzmax
875 :
876 : ! If CLUBB is running on it's own, this option determines if it is using:
877 : ! 1) an evenly-spaced grid;
878 : ! 2) a stretched (unevenly-spaced) grid entered on the thermodynamic grid
879 : ! levels (with momentum levels set halfway between thermodynamic levels);
880 : ! or
881 : ! 3) a stretched (unevenly-spaced) grid entered on the momentum grid levels
882 : ! (with thermodynamic levels set halfway between momentum levels).
883 : integer, intent(in) :: &
884 : grid_type
885 :
886 : character(len=*), intent(in) :: &
887 : zm_grid_fname, & ! Path and filename of file for momentum level altitudes
888 : zt_grid_fname ! Path and filename of file for thermodynamic level altitudes
889 :
890 : integer, intent(in) :: &
891 : file_unit ! Unit number for zt_grid_fname & zm_grid_fname (based on the OpenMP thread)
892 :
893 : ! Output Variables.
894 :
895 : ! If the CLUBB model is running by itself, but is using a stretched grid
896 : ! entered on thermodynamic levels (grid_type = 2), it needs to use the
897 : ! thermodynamic level altitudes as input.
898 : ! If the CLUBB model is running by itself, but is using a stretched grid
899 : ! entered on momentum levels (grid_type = 3), it needs to use the momentum
900 : ! level altitudes as input.
901 : real( kind = core_rknd ), dimension(nzmax), intent(out) :: &
902 : momentum_heights, & ! Momentum level altitudes (file input) [m]
903 : thermodynamic_heights ! Thermodynamic level altitudes (file input) [m]
904 :
905 : ! Local Variables.
906 :
907 : integer :: &
908 : zt_level_count, & ! Number of altitudes found in zt_grid_fname
909 : zm_level_count ! Number of altitudes found in zm_grid_fname
910 :
911 : integer :: input_status ! Status of file being read:
912 : ! > 0 ==> error reading file.
913 : ! = 0 ==> no error and more file to be read.
914 : ! < 0 ==> end of file indicator.
915 :
916 : ! Generic variable for storing file data while counting the number
917 : ! of file entries.
918 : real( kind = core_rknd ) :: generic_input_item
919 :
920 : integer :: k ! Loop index
921 :
922 : ! ---- Begin Code ----
923 :
924 : ! Declare the momentum level altitude array and the thermodynamic level
925 : ! altitude array to be 0 until overwritten.
926 0 : momentum_heights(1:nzmax) = 0.0_core_rknd
927 0 : thermodynamic_heights(1:nzmax) = 0.0_core_rknd
928 :
929 : ! Avoid uninitialized memory
930 0 : generic_input_item = 0.0_core_rknd
931 :
932 :
933 0 : if ( grid_type == 1 ) then
934 :
935 : ! Evenly-spaced grid.
936 : ! Grid level altitudes are based on a constant distance between them and
937 : ! a starting point for the bottom of the grid.
938 :
939 : ! As a way of error checking, make sure that there isn't any file entry
940 : ! for either momentum level altitudes or thermodynamic level altitudes.
941 0 : if ( zm_grid_fname /= '' ) then
942 0 : write(fstderr,*) err_header, &
943 : "An evenly-spaced grid has been selected. " &
944 0 : // " Please reset zm_grid_fname to ''."
945 0 : err_code = clubb_fatal_error
946 0 : return
947 : endif
948 0 : if ( zt_grid_fname /= '' ) then
949 0 : write(fstderr,*) err_header, &
950 : "An evenly-spaced grid has been selected. " &
951 0 : // " Please reset zt_grid_fname to ''."
952 0 : err_code = clubb_fatal_error
953 0 : return
954 : endif
955 :
956 :
957 0 : elseif ( grid_type == 2 ) then
958 :
959 : ! Stretched (unevenly-spaced) grid: stretched thermodynamic levels.
960 : ! Thermodynamic levels are defined according to a stretched grid that is
961 : ! entered through the use of an input file. Momentum levels are set
962 : ! halfway between thermodynamic levels. This is similar to a SAM-style
963 : ! stretched grid.
964 :
965 : ! As a way of error checking, make sure that there isn't any file entry
966 : ! for momentum level altitudes.
967 0 : if ( zm_grid_fname /= '' ) then
968 0 : write(fstderr,*) err_header, &
969 : "Thermodynamic level altitudes have been selected " &
970 : // "for use in a stretched (unevenly-spaced) grid. " &
971 0 : // " Please reset zm_grid_fname to ''."
972 0 : err_code = clubb_fatal_error
973 0 : return
974 : endif
975 :
976 : !$omp critical
977 : ! Open the file zt_grid_fname.
978 : open( unit=file_unit, file=zt_grid_fname, &
979 0 : status='old', action='read' )
980 :
981 : ! Find the number of thermodynamic level altitudes listed
982 : ! in file zt_grid_fname.
983 0 : zt_level_count = 0
984 0 : do
985 : read( unit=file_unit, fmt=*, iostat=input_status ) &
986 0 : generic_input_item
987 0 : if ( input_status < 0 ) exit ! end of file indicator
988 0 : if ( input_status > 0 ) then
989 0 : write(fstderr,*) err_header, & ! error reading input
990 0 : "Error reading thermodynamic level input file."
991 0 : err_code = clubb_fatal_error
992 0 : exit
993 : end if
994 0 : zt_level_count = zt_level_count + 1
995 : enddo
996 :
997 : ! Close the file zt_grid_fname.
998 0 : close( unit=file_unit )
999 : !$omp end critical
1000 :
1001 0 : if ( err_code == clubb_fatal_error ) return
1002 :
1003 : ! Check that the number of thermodynamic grid altitudes in the input file
1004 : ! matches the declared number of CLUBB grid levels (nzmax).
1005 0 : if ( zt_level_count /= nzmax ) then
1006 : write(fstderr,*) &
1007 : "The number of thermodynamic grid altitudes " &
1008 : // "listed in file " // trim(zt_grid_fname) &
1009 : // " does not match the number of CLUBB grid " &
1010 0 : // "levels specified in the model.in file."
1011 : write(fstderr,*) &
1012 0 : "Number of thermodynamic grid altitudes listed: ", &
1013 0 : zt_level_count
1014 : write(fstderr,*) &
1015 0 : "Number of CLUBB grid levels specified: ", nzmax
1016 0 : err_code = clubb_fatal_error
1017 0 : return
1018 : endif
1019 :
1020 : ! Read the thermodynamic level altitudes from zt_grid_fname.
1021 : call file_read_1d( file_unit, zt_grid_fname, nzmax, 1, & ! intent(in)
1022 0 : thermodynamic_heights ) ! intent(out)
1023 :
1024 : ! Check that each thermodynamic level altitude increases
1025 : ! in height as the thermodynamic level grid index increases.
1026 0 : do k = 2, nzmax, 1
1027 0 : if ( thermodynamic_heights(k) &
1028 0 : <= thermodynamic_heights(k-1) ) then
1029 : write(fstderr,*) &
1030 : "The declared thermodynamic level grid " &
1031 : // "altitudes are not increasing in height " &
1032 0 : // "as grid level index increases."
1033 : write(fstderr,*) &
1034 0 : "Grid index: ", k-1, ";", &
1035 0 : " Thermodynamic level altitude: ", &
1036 0 : thermodynamic_heights(k-1)
1037 : write(fstderr,*) &
1038 0 : "Grid index: ", k, ";", &
1039 0 : " Thermodynamic level altitude: ", &
1040 0 : thermodynamic_heights(k)
1041 0 : err_code = clubb_fatal_error
1042 0 : return
1043 : endif
1044 : enddo
1045 :
1046 :
1047 0 : elseif ( grid_type == 3 ) then
1048 :
1049 : ! Stretched (unevenly-spaced) grid: stretched momentum levels.
1050 : ! Momentum levels are defined according to a stretched grid that is
1051 : ! entered through the use of an input file. Thermodynamic levels are set
1052 : ! halfway between momentum levels. This is similar to a WRF-style
1053 : ! stretched grid.
1054 :
1055 : ! As a way of error checking, make sure that there isn't any file entry
1056 : ! for thermodynamic level altitudes.
1057 0 : if ( zt_grid_fname /= '' ) then
1058 : write(fstderr,*) &
1059 : "Momentum level altitudes have been selected " &
1060 : // "for use in a stretched (unevenly-spaced) grid. " &
1061 0 : // " Please reset zt_grid_fname to ''."
1062 0 : err_code = clubb_fatal_error
1063 0 : return
1064 : endif
1065 :
1066 : ! Open the file zm_grid_fname.
1067 : open( unit=file_unit, file=zm_grid_fname, &
1068 0 : status='old', action='read' )
1069 :
1070 : ! Find the number of momentum level altitudes
1071 : ! listed in file zm_grid_fname.
1072 0 : zm_level_count = 0
1073 0 : do
1074 : read( unit=file_unit, fmt=*, iostat=input_status ) &
1075 0 : generic_input_item
1076 0 : if ( input_status < 0 ) exit ! end of file indicator
1077 0 : if ( input_status > 0 ) then
1078 :
1079 0 : write(fstderr,*) err_header, &! error reading input
1080 0 : "Error reading momentum level input file."
1081 0 : err_code = clubb_fatal_error
1082 0 : return
1083 : end if
1084 0 : zm_level_count = zm_level_count + 1
1085 : enddo
1086 :
1087 : ! Close the file zm_grid_fname.
1088 0 : close( unit=file_unit )
1089 :
1090 : ! Check that the number of momentum grid altitudes in the input file
1091 : ! matches the declared number of CLUBB grid levels (nzmax).
1092 0 : if ( zm_level_count /= nzmax ) then
1093 : write(fstderr,*) &
1094 : "The number of momentum grid altitudes " &
1095 : // "listed in file " // trim(zm_grid_fname) &
1096 : // " does not match the number of CLUBB grid " &
1097 0 : // "levels specified in the model.in file."
1098 : write(fstderr,*) &
1099 0 : "Number of momentum grid altitudes listed: ", &
1100 0 : zm_level_count
1101 : write(fstderr,*) &
1102 0 : "Number of CLUBB grid levels specified: ", nzmax
1103 0 : err_code = clubb_fatal_error
1104 0 : return
1105 : endif
1106 :
1107 : ! Read the momentum level altitudes from zm_grid_fname.
1108 : call file_read_1d( file_unit, zm_grid_fname, nzmax, 1, & ! intent(in)
1109 0 : momentum_heights ) ! intent(out)
1110 :
1111 : ! Check that each momentum level altitude increases in height as the
1112 : ! momentum level grid index increases.
1113 0 : do k = 2, nzmax, 1
1114 0 : if ( momentum_heights(k) &
1115 0 : <= momentum_heights(k-1) ) then
1116 : write(fstderr,*) &
1117 : "The declared momentum level grid " &
1118 : // "altitudes are not increasing in height " &
1119 0 : // "as grid level index increases."
1120 : write(fstderr,*) &
1121 0 : "Grid index: ", k-1, ";", &
1122 0 : " Momentum level altitude: ", &
1123 0 : momentum_heights(k-1)
1124 : write(fstderr,*) &
1125 0 : "Grid index: ", k, ";", &
1126 0 : " Momentum level altitude: ", &
1127 0 : momentum_heights(k)
1128 0 : err_code = clubb_fatal_error
1129 0 : return
1130 : endif
1131 : enddo
1132 :
1133 :
1134 : endif
1135 :
1136 :
1137 : ! The purpose of this if statement is to avoid a compiler warning.
1138 : if ( generic_input_item > 0.0_core_rknd ) then
1139 : ! Do nothing
1140 : endif
1141 : ! Joshua Fasching June 2008
1142 :
1143 : return
1144 :
1145 : end subroutine read_grid_heights
1146 :
1147 : !=============================================================================
1148 : ! Wrapped in interface zt2zm
1149 0 : function redirect_interpolated_azm_k( gr, azt, k, zm_min )
1150 :
1151 : ! Description:
1152 : ! Calls the appropriate corresponding function based on l_cubic_temp
1153 : !-----------------------------------------------------------------------
1154 :
1155 : use clubb_precision, only: &
1156 : core_rknd ! Variable(s)
1157 :
1158 : use model_flags, only: &
1159 : l_cubic_interp, & ! Variable(s)
1160 : l_quintic_poly_interp
1161 :
1162 : use constants_clubb, only: &
1163 : fstdout ! Variable
1164 :
1165 : implicit none
1166 :
1167 : !---------------------------- Input Variables ----------------------------
1168 : type (grid), target, intent(in) :: gr
1169 :
1170 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
1171 : azt ! Variable on thermodynamic grid levels [units vary]
1172 :
1173 : integer, intent(in) :: &
1174 : k
1175 :
1176 : !---------------------------- Return Variable ----------------------------
1177 : real( kind = core_rknd ) :: &
1178 : redirect_interpolated_azm_k ! Variable when interp. to momentum levels
1179 :
1180 : !---------------------------- Optional Input Variable ----------------------------
1181 : real( kind = core_rknd ), optional, intent(in) :: &
1182 : zm_min
1183 :
1184 : !---------------------------- Local Variables ----------------------------
1185 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1186 0 : redirect_interpolated_azm_col ! Variable when interp. to momentum levels
1187 :
1188 : !---------------------------- Begin Code ----------------------------
1189 :
1190 0 : redirect_interpolated_azm_col = redirect_interpolated_azm_2D( gr%nz, 1, gr, azt, zm_min )
1191 :
1192 0 : redirect_interpolated_azm_k = redirect_interpolated_azm_col(1,k)
1193 :
1194 : return
1195 : end function redirect_interpolated_azm_k
1196 :
1197 : !=============================================================================
1198 : ! Wrapped in interface zt2zm
1199 0 : function redirect_interpolated_azm_1D( gr, azt, zm_min )
1200 :
1201 : ! Description:
1202 : ! Calls the appropriate corresponding function based on l_cubic_temp
1203 : !-----------------------------------------------------------------------
1204 :
1205 : use clubb_precision, only: &
1206 : core_rknd ! Variable(s)
1207 :
1208 : use model_flags, only: &
1209 : l_cubic_interp, & ! Variable(s)
1210 : l_quintic_poly_interp
1211 :
1212 : use constants_clubb, only: &
1213 : fstdout ! Variable
1214 :
1215 : implicit none
1216 :
1217 : !---------------------------- Input Variables ----------------------------
1218 : type (grid), target, intent(in) :: gr
1219 :
1220 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
1221 : azt ! Variable on thermodynamic grid levels [units vary]
1222 :
1223 : !---------------------------- Return Variable ----------------------------
1224 : real( kind = core_rknd ), dimension(gr%nz) :: &
1225 : redirect_interpolated_azm_1D ! Variable when interp. to momentum levels
1226 :
1227 : !---------------------------- Optional Input Variable ----------------------------
1228 : real( kind = core_rknd ), optional, intent(in) :: &
1229 : zm_min
1230 :
1231 : !---------------------------- Local Variables ----------------------------
1232 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1233 0 : azt_col ! Variable on thermodynamic grid levels [units vary]
1234 :
1235 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1236 0 : redirect_interpolated_azm_col ! Variable when interp. to momentum levels
1237 :
1238 : !---------------------------- Begin Code ----------------------------
1239 :
1240 0 : azt_col(1,:) = azt
1241 :
1242 0 : redirect_interpolated_azm_col = redirect_interpolated_azm_2D( gr%nz, 1, gr, azt_col, zm_min )
1243 :
1244 0 : redirect_interpolated_azm_1D = redirect_interpolated_azm_col(1,:)
1245 :
1246 0 : return
1247 : end function redirect_interpolated_azm_1D
1248 :
1249 : !=============================================================================
1250 : ! Wrapped in interface zt2zm
1251 393142464 : function redirect_interpolated_azm_2D( nz, ngrdcol, gr, azt, zm_min )
1252 :
1253 : ! Description:
1254 : ! Calls the appropriate corresponding function based on l_cubic_temp
1255 : !-----------------------------------------------------------------------
1256 :
1257 : use clubb_precision, only: &
1258 : core_rknd ! Variable(s)
1259 :
1260 : use model_flags, only: &
1261 : l_cubic_interp, & ! Variable(s)
1262 : l_quintic_poly_interp
1263 :
1264 : use constants_clubb, only: &
1265 : fstdout ! Variable
1266 :
1267 : implicit none
1268 :
1269 : !---------------------------- Input Variables ----------------------------
1270 : integer, intent(in) :: &
1271 : nz, &
1272 : ngrdcol
1273 :
1274 : type (grid), target, intent(in) :: gr
1275 :
1276 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1277 : azt ! Variable on thermodynamic grid levels [units vary]
1278 :
1279 : !---------------------------- Return Variable ----------------------------
1280 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1281 : redirect_interpolated_azm_2D ! Variable when interp. to momentum levels
1282 :
1283 : !---------------------------- Optional Input Variable ----------------------------
1284 : real( kind = core_rknd ), optional, intent(in) :: &
1285 : zm_min
1286 :
1287 : !---------------------------- Begin Code ----------------------------
1288 :
1289 : ! Sanity Check
1290 : if (l_quintic_poly_interp) then
1291 : if (.not. l_cubic_interp) then
1292 : write (fstdout, *) "Error: Model flag l_quintic_poly_interp should not be true if "&
1293 : //"l_cubic_interp is false."
1294 : error stop
1295 : end if
1296 : end if
1297 :
1298 : ! Redirect
1299 : if ( l_cubic_interp .and. nz >= 3 ) then
1300 : redirect_interpolated_azm_2D = cubic_interpolated_azm_2D( nz, ngrdcol, gr, azt )
1301 : else
1302 : call linear_interpolated_azm_2D( nz, ngrdcol, gr, azt, &
1303 : redirect_interpolated_azm_2D, &
1304 393142464 : zm_min )
1305 : end if
1306 :
1307 393142464 : return
1308 : end function redirect_interpolated_azm_2D
1309 :
1310 : !=============================================================================
1311 : ! Wrapped in interface zm2zt
1312 0 : function redirect_interpolated_azt_k( gr, azm, k, zt_min )
1313 :
1314 : ! Description:
1315 : ! Calls the appropriate corresponding function based on l_cubic_temp
1316 : !-----------------------------------------------------------------------
1317 :
1318 : use clubb_precision, only: &
1319 : core_rknd ! Variable(s)
1320 :
1321 : use model_flags, only: &
1322 : l_cubic_interp, & ! Variable(s)
1323 : l_quintic_poly_interp
1324 :
1325 : use constants_clubb, only: &
1326 : fstdout ! Variable
1327 :
1328 : implicit none
1329 :
1330 : type (grid), target, intent(in) :: gr
1331 :
1332 : !---------------------------- Input Variables ----------------------------
1333 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
1334 : azm ! Variable on momentum grid levels [units vary]
1335 :
1336 : integer, intent(in) :: &
1337 : k ! Vertical level
1338 :
1339 : !---------------------------- Return Variable ----------------------------
1340 : real( kind = core_rknd ) :: &
1341 : redirect_interpolated_azt_k ! Variable when interp. to momentum levels
1342 :
1343 : !---------------------------- Optional Input Variable ----------------------------
1344 : real( kind = core_rknd ), optional, intent(in) :: &
1345 : zt_min
1346 :
1347 : !---------------------------- Local Variables ----------------------------
1348 :
1349 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1350 0 : redirect_interpolated_azt_col ! Variable when interp. to momentum levels
1351 :
1352 : !---------------------------- Begin Code ----------------------------
1353 :
1354 0 : redirect_interpolated_azt_col = redirect_interpolated_azt_2D( gr%nz, 1, gr, azm, zt_min )
1355 :
1356 0 : redirect_interpolated_azt_k = redirect_interpolated_azt_col(1,k)
1357 :
1358 : return
1359 : end function redirect_interpolated_azt_k
1360 :
1361 :
1362 : !=============================================================================
1363 : ! Wrapped in interface zm2zt
1364 0 : function redirect_interpolated_azt_1D( gr, azm, zt_min )
1365 :
1366 : ! Description:
1367 : ! Calls the appropriate corresponding function based on l_cubic_temp
1368 : !-----------------------------------------------------------------------
1369 :
1370 : use clubb_precision, only: &
1371 : core_rknd ! Variable(s)
1372 :
1373 : use model_flags, only: &
1374 : l_cubic_interp, & ! Variable(s)
1375 : l_quintic_poly_interp
1376 :
1377 : use constants_clubb, only: &
1378 : fstdout ! Variable
1379 :
1380 : implicit none
1381 :
1382 : type (grid), target, intent(in) :: gr
1383 :
1384 : !---------------------------- Input Variables ----------------------------
1385 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
1386 : azm ! Variable on momentum grid levels [units vary]
1387 :
1388 : !---------------------------- Return Variable ----------------------------
1389 : real( kind = core_rknd ), dimension(gr%nz) :: &
1390 : redirect_interpolated_azt_1D ! Variable when interp. to momentum levels
1391 :
1392 : !---------------------------- Optional Input Variable ----------------------------
1393 : real( kind = core_rknd ), optional, intent(in) :: &
1394 : zt_min
1395 :
1396 : !---------------------------- Local Variables ----------------------------
1397 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1398 0 : azm_col ! Variable on thermodynamic grid levels [units vary]
1399 :
1400 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1401 0 : redirect_interpolated_azt_col ! Variable when interp. to momentum levels
1402 :
1403 : !---------------------------- Begin Code ----------------------------
1404 :
1405 0 : azm_col(1,:) = azm
1406 :
1407 0 : redirect_interpolated_azt_col = redirect_interpolated_azt_2D( gr%nz, 1, gr, azm_col, zt_min )
1408 :
1409 0 : redirect_interpolated_azt_1D = redirect_interpolated_azt_col(1,:)
1410 :
1411 0 : return
1412 : end function redirect_interpolated_azt_1D
1413 :
1414 : !=============================================================================
1415 : ! Wrapped in interface zm2zt
1416 424415160 : function redirect_interpolated_azt_2D( nz, ngrdcol, gr, azm, zt_min )
1417 :
1418 : ! Description:
1419 : ! Calls the appropriate corresponding function based on l_cubic_temp
1420 : !-----------------------------------------------------------------------
1421 :
1422 : use clubb_precision, only: &
1423 : core_rknd ! Variable(s)
1424 :
1425 : use model_flags, only: &
1426 : l_cubic_interp, & ! Variable(s)
1427 : l_quintic_poly_interp
1428 :
1429 : use constants_clubb, only: &
1430 : fstdout ! Variable
1431 :
1432 : implicit none
1433 :
1434 : !---------------------------- Input Variables ----------------------------
1435 : integer, intent(in) :: &
1436 : nz, &
1437 : ngrdcol
1438 :
1439 : type (grid), target, intent(in) :: gr
1440 :
1441 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1442 : azm ! Variable on momentum grid levels [units vary]
1443 :
1444 : !---------------------------- Return Variable----------------------------
1445 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1446 : redirect_interpolated_azt_2D ! Variable when interp. to momentum levels
1447 :
1448 : !---------------------------- Optional Input Variable ----------------------------
1449 : real( kind = core_rknd ), optional, intent(in) :: &
1450 : zt_min
1451 :
1452 : !---------------------------- Begin Code ----------------------------
1453 :
1454 : ! Sanity Check
1455 : if (l_quintic_poly_interp) then
1456 : if (.not. l_cubic_interp) then
1457 : write (fstdout, *) "Error: Model flag l_quintic_poly_interp should not be true if "&
1458 : //"l_cubic_interp is false."
1459 : error stop
1460 : end if
1461 : end if
1462 :
1463 : ! Redirect
1464 : if (l_cubic_interp .and. nz >= 3 ) then
1465 : redirect_interpolated_azt_2D = cubic_interpolated_azt_2D( nz, ngrdcol, gr, azm )
1466 : else
1467 : call linear_interpolated_azt_2D( nz, ngrdcol, gr, azm, &
1468 : redirect_interpolated_azt_2D, &
1469 424415160 : zt_min )
1470 : end if
1471 :
1472 424415160 : return
1473 : end function redirect_interpolated_azt_2D
1474 :
1475 : !=============================================================================
1476 393142464 : subroutine linear_interpolated_azm_2D( nz, ngrdcol, gr, azt, &
1477 393142464 : linear_interpolated_azm, &
1478 : zm_min )
1479 :
1480 : ! Description:
1481 : ! Function to interpolate a variable located on the thermodynamic grid
1482 : ! levels (azt) to the momentum grid levels (azm). This function inputs the
1483 : ! entire azt array and outputs the results as an azm array. The
1484 : ! formulation used is compatible with a stretched (unevenly-spaced) grid.
1485 : !-----------------------------------------------------------------------
1486 :
1487 : use clubb_precision, only: &
1488 : core_rknd ! Variable(s)
1489 :
1490 : use interpolation, only: &
1491 : linear_interp_factor ! Procedure(s)
1492 :
1493 : implicit none
1494 :
1495 : integer, intent(in) :: &
1496 : nz, &
1497 : ngrdcol
1498 :
1499 : type (grid), target, intent(in) :: gr
1500 :
1501 : ! ------------------------------ Input Variable ------------------------------
1502 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1503 : azt ! Variable on thermodynamic grid levels [units vary]
1504 :
1505 : ! ------------------------------ Return Variable ------------------------------
1506 : real( kind = core_rknd ), intent(out), dimension(ngrdcol,nz) :: &
1507 : linear_interpolated_azm ! Variable when interp. to momentum levels
1508 :
1509 : !---------------------------- Optional Input Variable ----------------------------
1510 : real( kind = core_rknd ), optional, intent(in) :: &
1511 : zm_min
1512 :
1513 : ! ------------------------------ Local Variable ------------------------------
1514 : integer :: i, k ! Grid level loop index
1515 :
1516 : ! ------------------------------ Begin Code ------------------------------
1517 :
1518 : !$acc data copyin( azt, gr, gr%weights_zt2zm, gr%zt, gr%zm ) &
1519 : !$acc copyout( linear_interpolated_azm )
1520 :
1521 : !$acc parallel loop gang vector default(present)
1522 6564564864 : do i = 1, ngrdcol
1523 6564564864 : linear_interpolated_azm(i,1) = azt(i,2)
1524 : end do
1525 : !$acc end parallel loop
1526 :
1527 : ! Interpolate the value of a thermodynamic-level variable to the central
1528 : ! momentum level, k, between two successive thermodynamic levels using
1529 : ! linear interpolation.
1530 : !$acc parallel loop gang vector collapse(2) default(present)
1531 33023966976 : do k = 2, nz-1
1532 >54525*10^7 : do i = 1, ngrdcol
1533 >10244*10^8 : linear_interpolated_azm(i,k) = gr%weights_zt2zm(i,k,1) &
1534 >15693*10^8 : * ( azt(i,k+1) - azt(i,k) ) + azt(i,k)
1535 : end do
1536 : end do
1537 : !$acc end parallel loop
1538 :
1539 : ! Set the value of the thermodynamic-level variable, azt, at the uppermost
1540 : ! level of the model, which is a momentum level. The name of the variable
1541 : ! when interpolated/extended to momentum levels is azm.
1542 : ! Use a linear extension based on the values of azt at levels gr%nz and
1543 : ! gr%nz-1 to find the value of azm at level gr%nz (the uppermost level
1544 : ! in the model).
1545 : !$acc parallel loop gang vector default(present)
1546 6564564864 : do i = 1, ngrdcol
1547 12342844800 : linear_interpolated_azm(i,nz) &
1548 12342844800 : = ( ( azt(i,nz) - azt(i,nz-1) ) / ( gr%zt(i,nz) - gr%zt(i,nz-1) ) ) &
1549 25078832064 : * ( gr%zm(i,nz) - gr%zt(i,nz) ) + azt(i,nz)
1550 : end do
1551 : !$acc end parallel loop
1552 :
1553 393142464 : if ( present(zm_min) ) then
1554 : !$acc parallel loop gang vector collapse(2) default(present) copyin(zm_min)
1555 16905125952 : do k = 1, nz
1556 >27919*10^7 : do i = 1, ngrdcol
1557 >27899*10^7 : linear_interpolated_azm(i,k) = max( linear_interpolated_azm(i,k), zm_min )
1558 : end do
1559 : end do
1560 : !$acc end parallel loop
1561 : end if
1562 :
1563 : !$acc end data
1564 :
1565 393142464 : return
1566 :
1567 : end subroutine linear_interpolated_azm_2D
1568 :
1569 : !=============================================================================
1570 0 : function zt2zm2zt( nz, ngrdcol, gr, azt, zt_min )
1571 :
1572 : ! Description:
1573 : ! Function to interpolate a variable located on the thermodynamic grid
1574 : ! levels (azt) to the momentum grid levels (azm), then interpolate back
1575 : ! to thermodynamic grid levels (azt).
1576 : !
1577 : ! Note:
1578 : ! This is intended for smoothing variables.
1579 : !-----------------------------------------------------------------------------
1580 :
1581 : use clubb_precision, only: &
1582 : core_rknd ! Variable(s)
1583 :
1584 : implicit none
1585 :
1586 : ! ------------------------------ Input Variable ------------------------------
1587 : integer, intent(in) :: &
1588 : nz, &
1589 : ngrdcol
1590 :
1591 : type (grid), target, intent(in) :: gr
1592 :
1593 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1594 : azt ! Variable on thermodynamic grid levels [units vary]
1595 :
1596 : ! ------------------------------ Return Variable ------------------------------
1597 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1598 : zt2zm2zt ! Variable when interp. to momentum levels
1599 :
1600 : !---------------------------- Optional Input Variable ----------------------------
1601 : real( kind = core_rknd ), optional, intent(in) :: &
1602 : zt_min
1603 :
1604 : ! ------------------------------ Local Variable ------------------------------
1605 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1606 0 : azt_zm
1607 :
1608 : ! ------------------------------ Begin Code ------------------------------
1609 :
1610 : !$acc data copyin( azt ) &
1611 : !$acc copyout( zt2zm2zt ) &
1612 : !$acc create( azt_zm )
1613 :
1614 : ! Interpolate azt to momentum levels
1615 0 : azt_zm = zt2zm( nz, ngrdcol, gr, azt )
1616 :
1617 : ! Interpolate back to termodynamic levels
1618 0 : zt2zm2zt = zm2zt( nz, ngrdcol, gr, azt_zm, zt_min )
1619 :
1620 : !$acc end data
1621 :
1622 0 : return
1623 :
1624 : end function zt2zm2zt
1625 :
1626 : !=============================================================================
1627 26805168 : function zm2zt2zm( nz, ngrdcol, gr, azm, zm_min )
1628 :
1629 : ! Description:
1630 : ! Function to interpolate a variable located on the momentum grid
1631 : ! levels(azm) to thermodynamic grid levels (azt), then interpolate
1632 : ! back to momentum grid levels (azm).
1633 : !
1634 : ! Note:
1635 : ! This is intended for smoothing variables.
1636 : !-----------------------------------------------------------------------------
1637 :
1638 : use clubb_precision, only: &
1639 : core_rknd ! Variable(s)
1640 :
1641 : implicit none
1642 :
1643 : ! ------------------------------ Input Variable ------------------------------
1644 : integer, intent(in) :: &
1645 : nz, &
1646 : ngrdcol
1647 :
1648 : type (grid), target, intent(in) :: gr
1649 :
1650 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1651 : azm ! Variable on momentum grid levels [units vary]
1652 :
1653 : ! ------------------------------ Return Variable ------------------------------
1654 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1655 : zm2zt2zm ! Variable when interp. to momentum levels
1656 :
1657 : !---------------------------- Optional Input Variable ----------------------------
1658 : real( kind = core_rknd ), optional, intent(in) :: &
1659 : zm_min
1660 :
1661 : ! ------------------------------ Local Variable ------------------------------
1662 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1663 26805168 : azm_zt
1664 :
1665 : ! ------------------------------ Begin Code ------------------------------
1666 :
1667 : !$acc data copyin( azm ) &
1668 : !$acc copyout( zm2zt2zm ) &
1669 : !$acc create( azm_zt )
1670 :
1671 : ! Interpolate azt to termodynamic levels
1672 26805168 : azm_zt = zm2zt( nz, ngrdcol, gr, azm )
1673 :
1674 : ! Interpolate back to momentum levels
1675 26805168 : zm2zt2zm = zt2zm( nz, ngrdcol, gr, azm_zt, zm_min )
1676 :
1677 : !$acc end data
1678 :
1679 26805168 : return
1680 :
1681 : end function zm2zt2zm
1682 :
1683 : !=============================================================================
1684 : function cubic_interpolated_azm_2D( nz, ngrdcol, gr, azt )
1685 :
1686 : ! Description:
1687 : ! Function to interpolate a variable located on the momentum grid
1688 : ! levels (azt) to the thermodynamic grid levels (azm). This function outputs the
1689 : ! value of azt at a all grid levels using Steffen's monotonic cubic
1690 : ! interpolation implemented by Tak Yamaguchi.
1691 : !
1692 : ! References:
1693 : ! None
1694 : !-----------------------------------------------------------------------
1695 :
1696 : use interpolation, only: &
1697 : mono_cubic_interp ! Procedure(s)
1698 :
1699 : use clubb_precision, only: &
1700 : core_rknd ! Variable(s)
1701 :
1702 : implicit none
1703 :
1704 : ! Input Variables
1705 : integer, intent(in) :: &
1706 : nz, &
1707 : ngrdcol
1708 :
1709 : type (grid), target, intent(in) :: gr
1710 :
1711 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1712 : azt
1713 :
1714 : ! Return Variable
1715 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1716 : cubic_interpolated_azm_2D
1717 :
1718 : ! Local Variable(s)
1719 : integer :: &
1720 : k, i
1721 :
1722 : integer :: km1, k00, kp1, kp2
1723 :
1724 : ! ---- Begin Code ----
1725 :
1726 : do k = 1, nz
1727 : do i = 1, ngrdcol
1728 :
1729 : ! k levels are based on Tak's find_indices subroutine -dschanen 24 Oct 2011
1730 : if ( k == nz-1 ) then
1731 : km1 = nz-2
1732 : kp1 = nz
1733 : kp2 = nz
1734 : k00 = nz-1
1735 : else if ( k == nz ) then ! Extrapolation
1736 : km1 = nz
1737 : kp1 = nz
1738 : kp2 = nz
1739 : k00 = nz-1
1740 : else if ( k == 1 ) then
1741 : km1 = 1
1742 : kp1 = 2
1743 : kp2 = 3
1744 : k00 = 1
1745 : else
1746 : km1 = k-1
1747 : kp1 = k+1
1748 : kp2 = k+2
1749 : k00 = k
1750 : end if
1751 :
1752 : ! Do the actual interpolation.
1753 : ! Use a cubic monotonic spline interpolation.
1754 : cubic_interpolated_azm_2D(i,k) = &
1755 : mono_cubic_interp( gr%zm(i,k), km1, k00, kp1, kp2, &
1756 : gr%zt(i,km1), gr%zt(i,k00), gr%zt(i,kp1), gr%zt(i,kp2), &
1757 : azt(i,km1), azt(i,k00), azt(i,kp1), azt(i,kp2) )
1758 : end do
1759 : end do
1760 :
1761 : return
1762 :
1763 : end function cubic_interpolated_azm_2D
1764 :
1765 : !=============================================================================
1766 4467528 : subroutine calc_zt2zm_weights( nz, ngrdcol, &
1767 : gr )
1768 :
1769 : ! Description:
1770 : ! Function used to help in an interpolation of a variable (var_zt) located
1771 : ! on the thermodynamic grid levels (azt) to the momentum grid levels (azm).
1772 : ! This function computes a weighting factor for both the upper thermodynamic
1773 : ! level (k+1) and the lower thermodynamic level (k) applied to the central
1774 : ! momentum level (k). For the uppermost momentum grid level (k=gr%nz), a
1775 : ! weighting factor for both the thermodynamic level at gr%nz and the
1776 : ! thermodynamic level at gr%nz-1 are calculated based on the use of a
1777 : ! linear extension. This function outputs the weighting factors at a single
1778 : ! momentum grid level (k). The formulation used is compatible with a
1779 : ! stretched (unevenly-spaced) grid. The weights are defined as follows:
1780 : !
1781 : ! ---var_zt(k+1)------------------------------------------- t(k+1)
1782 : ! azt_weight(t_above) = factor
1783 : ! ===========var_zt(interp)================================ m(k)
1784 : ! azt_weight(t_below) = 1 - factor
1785 : ! ---var_zt(k)--------------------------------------------- t(k)
1786 : !
1787 : ! The vertical indices t(k+1), m(k), and t(k) correspond with altitudes
1788 : ! zt(k+1), zm(k), and zt(k), respectively. The letter "t" is used for
1789 : ! thermodynamic levels and the letter "m" is used for momentum levels.
1790 : !
1791 : ! For all levels k < gr%nz:
1792 : !
1793 : ! The formula for a linear interpolation is given by:
1794 : !
1795 : ! var_zt( interp to zm(k) )
1796 : ! = [ ( var_zt(k+1) - var_zt(k) ) / ( zt(k+1) - zt(k) ) ]
1797 : ! * ( zm(k) - zt(k) ) + var_zt(k);
1798 : !
1799 : ! which can be rewritten as:
1800 : !
1801 : ! var_zt( interp to zm(k) )
1802 : ! = [ ( zm(k) - zt(k) ) / ( zt(k+1) - zt(k) ) ]
1803 : ! * ( var_zt(k+1) - var_zt(k) ) + var_zt(k).
1804 : !
1805 : ! Furthermore, the formula can be rewritten as:
1806 : !
1807 : ! var_zt( interp to zm(k) )
1808 : ! = factor * var_zt(k+1) + ( 1 - factor ) * var_zt(k);
1809 : !
1810 : ! where:
1811 : !
1812 : ! factor = ( zm(k) - zt(k) ) / ( zt(k+1) - zt(k) ).
1813 : !
1814 : ! One of the important uses of this function is in situations where the
1815 : ! variable to be interpolated is being treated IMPLICITLY in an equation.
1816 : ! Usually, the variable to be interpolated is involved in a derivative (such
1817 : ! as d(var_zt)/dz in the diagram below). For the term of the equation
1818 : ! containing the derivative, grid weights are needed for two interpolations,
1819 : ! rather than just one interpolation. Thus, four grid weights (labeled
1820 : ! A(k), B(k), C(k), and D(k) in the diagram below) are needed.
1821 : !
1822 : ! ---var_zt(k+1)------------------------------------------- t(k+1)
1823 : ! A(k)
1824 : ! ===========var_zt(interp)================================ m(k)
1825 : ! B(k) = 1 - A(k)
1826 : ! ---var_zt(k)-----------d(var_zt)/dz---------------------- t(k)
1827 : ! C(k)
1828 : ! ===========var_zt(interp)================================ m(k-1)
1829 : ! D(k) = 1 - C(k)
1830 : ! ---var_zt(k-1)------------------------------------------- t(k-1)
1831 : !
1832 : ! The vertical indices t(k+1), m(k), t(k), m(k-1), and t(k-1) correspond
1833 : ! with altitudes zt(k+1), zm(k), zt(k), zm(k-1), and zt(k-1), respectively.
1834 : ! The letter "t" is used for thermodynamic levels and the letter "m" is used
1835 : ! for momentum levels.
1836 : !
1837 : ! The grid weights, indexed around the central thermodynamic level (k), are
1838 : ! defined as follows:
1839 : !
1840 : ! A(k) = ( zm(k) - zt(k) ) / ( zt(k+1) - zt(k) );
1841 : !
1842 : ! which is the same as "factor" for the interpolation to momentum
1843 : ! level (k). In the code, this interpolation is referenced as
1844 : ! gr%weights_zt2zm(t_above,mk), which can be read as "grid weight in a zt2zm
1845 : ! interpolation of the thermodynamic level above momentum level (k) (applied
1846 : ! to momentum level (k))".
1847 : !
1848 : ! B(k) = 1 - [ ( zm(k) - zt(k) ) / ( zt(k+1) - zt(k) ) ]
1849 : ! = 1 - A(k);
1850 : !
1851 : ! which is the same as "1 - factor" for the interpolation to momentum
1852 : ! level (k). In the code, this interpolation is referenced as
1853 : ! gr%weights_zt2zm(t_below,mk), which can be read as "grid weight in a zt2zm
1854 : ! interpolation of the thermodynamic level below momentum level (k) (applied
1855 : ! to momentum level (k))".
1856 : !
1857 : ! C(k) = ( zm(k-1) - zt(k-1) ) / ( zt(k) - zt(k-1) );
1858 : !
1859 : ! which is the same as "factor" for the interpolation to momentum
1860 : ! level (k-1). In the code, this interpolation is referenced as
1861 : ! gr%weights_zt2zm(t_above,mkm1), which can be read as "grid weight in a
1862 : ! zt2zm interpolation of the thermodynamic level above momentum level (k-1)
1863 : ! (applied to momentum level (k-1))".
1864 : !
1865 : ! D(k) = 1 - [ ( zm(k-1) - zt(k-1) ) / ( zt(k) - zt(k-1) ) ]
1866 : ! = 1 - C(k);
1867 : !
1868 : ! which is the same as "1 - factor" for the interpolation to momentum
1869 : ! level (k-1). In the code, this interpolation is referenced as
1870 : ! gr%weights_zt2zm(t_below,mkm1), which can be read as "grid weight in a
1871 : ! zt2zm interpolation of the thermodynamic level below momentum level (k-1)
1872 : ! (applied to momentum level (k-1))".
1873 : !
1874 : ! Additionally, as long as the central thermodynamic level (k) in the above
1875 : ! scenario is not the uppermost thermodynamic level or the lowermost
1876 : ! thermodynamic level (k /= gr%nz and k /= 1), the four weighting factors
1877 : ! have the following relationships: A(k) = C(k+1) and B(k) = D(k+1).
1878 : !
1879 : !
1880 : ! Special condition for uppermost grid level, k = gr%nz:
1881 : !
1882 : ! The uppermost momentum grid level is above the uppermost thermodynamic
1883 : ! grid level. Thus, a linear extension is used at this level.
1884 : !
1885 : ! For level k = gr%nz:
1886 : !
1887 : ! The formula for a linear extension is given by:
1888 : !
1889 : ! var_zt( extend to zm(k) )
1890 : ! = [ ( var_zt(k) - var_zt(k-1) ) / ( zt(k) - zt(k-1) ) ]
1891 : ! * ( zm(k) - zt(k-1) ) + var_zt(k-1);
1892 : !
1893 : ! which can be rewritten as:
1894 : !
1895 : ! var_zt( extend to zm(k) )
1896 : ! = [ ( zm(k) - zt(k-1) ) / ( zt(k) - zt(k-1) ) ]
1897 : ! * ( var_zt(k) - var_zt(k-1) ) + var_zt(k-1).
1898 : !
1899 : ! Furthermore, the formula can be rewritten as:
1900 : !
1901 : ! var_zt( extend to zm(k) )
1902 : ! = factor * var_zt(k) + ( 1 - factor ) * var_zt(k-1);
1903 : !
1904 : ! where:
1905 : !
1906 : ! factor = ( zm(k) - zt(k-1) ) / ( zt(k) - zt(k-1) ).
1907 : !
1908 : ! Due to the fact that a linear extension is being used, the value of factor
1909 : ! will be greater than 1. The weight of thermodynamic level k = gr%nz on
1910 : ! momentum level k = gr%nz equals the value of factor. The weight of
1911 : ! thermodynamic level k = gr%nz-1 on momentum level k = gr%nz equals
1912 : ! 1 - factor, which is less than 0. However, the sum of the two weights
1913 : ! equals 1.
1914 : !
1915 : !
1916 : ! Brian Griffin; September 12, 2008.
1917 : !
1918 : !-----------------------------------------------------------------------
1919 :
1920 : use constants_clubb, only: &
1921 : one ! Constant(s)
1922 :
1923 : use clubb_precision, only: &
1924 : core_rknd ! Variable(s)
1925 :
1926 : implicit none
1927 :
1928 : ! Constant parameters
1929 : integer, parameter :: &
1930 : t_above = 1, & ! Upper thermodynamic level.
1931 : t_below = 2 ! Lower thermodynamic level.
1932 :
1933 : ! Input Variable
1934 : integer, intent(in) :: &
1935 : nz, &
1936 : ngrdcol
1937 :
1938 : ! Input/Output Variable
1939 : type (grid), target, intent(inout) :: gr
1940 :
1941 : ! Local Variables
1942 : real( kind = core_rknd ) :: factor
1943 :
1944 : integer :: i, k
1945 :
1946 379739880 : do k = 1, nz-1
1947 6270643080 : do i = 1, ngrdcol
1948 :
1949 : ! The top model level (gr%nz) is formulated differently because the top
1950 : ! momentum level is above the top thermodynamic level. A linear
1951 : ! extension is required, rather than linear interpolation.
1952 : ! Note: Variable "factor" will be greater than 1 in this situation.
1953 0 : gr%weights_zt2zm(i,k,t_above) = ( gr%zm(i,k) - gr%zt(i,k) ) &
1954 5890903200 : / ( gr%zt(i,k+1) - gr%zt(i,k) )
1955 :
1956 : ! Weight of lower thermodynamic level on momentum level.
1957 6266175552 : gr%weights_zt2zm(i,k,t_below) = one - gr%weights_zt2zm(i,k,t_above)
1958 :
1959 : end do
1960 : end do
1961 :
1962 : ! At most levels, the momentum level is found in-between two
1963 : ! thermodynamic levels. Linear interpolation is used.
1964 74597328 : do i = 1, ngrdcol
1965 0 : gr%weights_zt2zm(i,nz,t_above) = ( gr%zm(i,nz) - gr%zt(i,nz-1) ) &
1966 70129800 : / ( gr%zt(i,nz) - gr%zt(i,nz-1) )
1967 :
1968 74597328 : gr%weights_zt2zm(i,nz,t_below) = one - gr%weights_zt2zm(i,nz,t_above)
1969 : end do
1970 :
1971 4467528 : return
1972 :
1973 : end subroutine calc_zt2zm_weights
1974 :
1975 : !=============================================================================
1976 424415160 : subroutine linear_interpolated_azt_2D( nz, ngrdcol, gr, azm, &
1977 424415160 : linear_interpolated_azt, &
1978 : zt_min )
1979 :
1980 : ! Description:
1981 : ! Function to interpolate a variable located on the momentum grid levels
1982 : ! (azm) to the thermodynamic grid levels (azt). This function inputs the
1983 : ! entire azm array and outputs the results as an azt array. The formulation
1984 : ! used is compatible with a stretched (unevenly-spaced) grid.
1985 : !-----------------------------------------------------------------------
1986 :
1987 : use clubb_precision, only: &
1988 : core_rknd ! Variable(s)
1989 :
1990 : use interpolation, only: &
1991 : linear_interp_factor ! Procedure(s)
1992 :
1993 : implicit none
1994 :
1995 : integer, intent(in) :: &
1996 : nz, &
1997 : ngrdcol
1998 :
1999 : type (grid), target, intent(in) :: gr
2000 :
2001 : ! ------------------------------ Input Variable ------------------------------
2002 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
2003 : azm ! Variable on momentum grid levels [units vary]
2004 :
2005 : ! ------------------------------ Output Variable ------------------------------
2006 : real( kind = core_rknd ), intent(out), dimension(ngrdcol,nz) :: &
2007 : linear_interpolated_azt ! Variable when interp. to thermodynamic levels
2008 :
2009 : !---------------------------- Optional Input Variable ----------------------------
2010 : real( kind = core_rknd ), optional, intent(in) :: &
2011 : zt_min
2012 :
2013 : ! ------------------------------ Local Variable ------------------------------
2014 : integer :: i, k ! Grid level loop index
2015 :
2016 : ! ------------------------------ Begin Code ------------------------------
2017 :
2018 : !$acc data copyin( azm, gr, gr%weights_zm2zt, gr%zt, gr%zm ) &
2019 : !$acc copyout( linear_interpolated_azt )
2020 :
2021 : ! Set the value of the momentum-level variable, azm, at the lowermost level
2022 : ! of the model (below the model lower boundary), which is a thermodynamic
2023 : ! level. The name of the variable when interpolated/extended to
2024 : ! thermodynamic levels is azt.
2025 : ! Use a linear extension based on the values of azm at levels 1 and 2 to
2026 : ! find the value of azt at level 1 (the lowermost level in the model).
2027 : !$acc parallel loop gang vector default(present)
2028 7086746160 : do i = 1, ngrdcol
2029 13324662000 : linear_interpolated_azt(i,1) &
2030 13324662000 : = ( ( azm(i,2) - azm(i,1) ) / ( gr%zm(i,2) - gr%zm(i,1) ) ) &
2031 20411408160 : * ( gr%zt(i,1) - gr%zm(i,1) ) + azm(i,1)
2032 : end do
2033 : !$acc end parallel loop
2034 :
2035 : ! Interpolate the value of a momentum-level variable to the central
2036 : ! thermodynamic level, k, between two successive momentum levels using
2037 : ! linear interpolation.
2038 : !$acc parallel loop gang vector collapse(2) default(present)
2039 36075288600 : do k = 2, nz
2040 >59571*10^7 : do i = 1, ngrdcol
2041 >11192*10^8 : linear_interpolated_azt(i,k) = gr%weights_zm2zt(i,k,1) &
2042 >17145*10^8 : * ( azm(i,k) - azm(i,k-1) ) + azm(i,k-1)
2043 : end do
2044 : end do
2045 : !$acc end parallel loop
2046 :
2047 424415160 : if ( present(zt_min) ) then
2048 : !$acc parallel loop gang vector collapse(2) default(present) copyin(zt_min)
2049 19210370400 : do k = 1, nz
2050 >31726*10^7 : do i = 1, ngrdcol
2051 >31703*10^7 : linear_interpolated_azt(i,k) = max( linear_interpolated_azt(i,k), zt_min )
2052 : end do
2053 : end do
2054 : !$acc end parallel loop
2055 : end if
2056 :
2057 : !$acc end data
2058 :
2059 424415160 : return
2060 :
2061 : end subroutine linear_interpolated_azt_2D
2062 :
2063 : !=============================================================================
2064 : function cubic_interpolated_azt_2D( nz, ngrdcol, gr, azm )
2065 :
2066 : ! Description:
2067 : ! Function to interpolate a variable located on the momentum grid
2068 : ! levels (azm) to the thermodynamic grid levels (azt). This function outputs the
2069 : ! value of azt at a all grid levels using Steffen's monotonic cubic
2070 : ! interpolation implemented by Tak Yamaguchi.
2071 : !
2072 : ! References:
2073 : ! None
2074 : !-----------------------------------------------------------------------
2075 :
2076 : use interpolation, only: &
2077 : mono_cubic_interp ! Procedure(s)
2078 :
2079 : use clubb_precision, only: &
2080 : core_rknd ! Variable(s)
2081 :
2082 : implicit none
2083 :
2084 : ! Input Variables
2085 : integer, intent(in) :: &
2086 : nz, &
2087 : ngrdcol
2088 :
2089 : type (grid), target, intent(in) :: gr
2090 :
2091 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
2092 : azm
2093 :
2094 : ! Return Variable
2095 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
2096 : cubic_interpolated_azt_2D
2097 :
2098 : ! Local Variable(s)
2099 : integer :: &
2100 : k, i
2101 :
2102 : integer :: km1, k00, kp1, kp2
2103 :
2104 : ! ---- Begin Code ----
2105 :
2106 : do k = 1, nz
2107 : do i = 1, ngrdcol
2108 :
2109 : ! k levels are based on Tak's find_indices subroutine -dschanen 24 Oct 2011
2110 : if ( k == nz ) then
2111 : km1 = nz-2
2112 : kp1 = nz
2113 : kp2 = nz
2114 : k00 = nz-1
2115 : else if ( k == 2 ) then
2116 : km1 = 1
2117 : kp1 = 2
2118 : kp2 = 3
2119 : k00 = 1
2120 : else if ( k == 1 ) then ! Extrapolation for the ghost point
2121 : km1 = nz
2122 : k00 = 1
2123 : kp1 = 2
2124 : kp2 = 3
2125 : else
2126 : km1 = k-2
2127 : kp1 = k
2128 : kp2 = k+1
2129 : k00 = k-1
2130 : end if
2131 :
2132 : ! Do the actual interpolation.
2133 : ! Use a cubic monotonic spline interpolation.
2134 : cubic_interpolated_azt_2D(i,k) = &
2135 : mono_cubic_interp( gr%zt(i,k), km1, k00, kp1, kp2, &
2136 : gr%zm(i,km1), gr%zm(i,k00), gr%zm(i,kp1), gr%zm(i,kp2), &
2137 : azm(i,km1), azm(i,k00), azm(i,kp1), azm(i,kp2) )
2138 : end do
2139 : end do
2140 :
2141 : return
2142 :
2143 : end function cubic_interpolated_azt_2D
2144 :
2145 : !=============================================================================
2146 4467528 : subroutine calc_zm2zt_weights( nz, ngrdcol, &
2147 : gr )
2148 :
2149 : ! Description:
2150 : ! Function used to help in an interpolation of a variable (var_zm) located
2151 : ! on the momentum grid levels (azm) to the thermodynamic grid levels (azt).
2152 : ! This function computes a weighting factor for both the upper momentum
2153 : ! level (k) and the lower momentum level (k-1) applied to the central
2154 : ! thermodynamic level (k). For the lowermost thermodynamic grid
2155 : ! level (k=1), a weighting factor for both the momentum level at 1 and the
2156 : ! momentum level at 2 are calculated based on the use of a linear extension.
2157 : ! This function outputs the weighting factors at a single thermodynamic grid
2158 : ! level (k). The formulation used is compatible with a stretched
2159 : ! (unevenly-spaced) grid. The weights are defined as follows:
2160 : !
2161 : ! ===var_zm(k)============================================= m(k)
2162 : ! azm_weight(m_above) = factor
2163 : ! -----------var_zm(interp)-------------------------------- t(k)
2164 : ! azm_weight(m_below) = 1 - factor
2165 : ! ===var_zm(k-1)=========================================== m(k-1)
2166 : !
2167 : ! The vertical indices m(k), t(k), and m(k-1) correspond with altitudes
2168 : ! zm(k), zt(k), and zm(k-1), respectively. The letter "t" is used for
2169 : ! thermodynamic levels and the letter "m" is used for momentum levels.
2170 : !
2171 : ! For all levels k > 1:
2172 : !
2173 : ! The formula for a linear interpolation is given by:
2174 : !
2175 : ! var_zm( interp to zt(k) )
2176 : ! = [ ( var_zm(k) - var_zm(k-1) ) / ( zm(k) - zm(k-1) ) ]
2177 : ! * ( zt(k) - zm(k-1) ) + var_zm(k-1);
2178 : !
2179 : ! which can be rewritten as:
2180 : !
2181 : ! var_zm( interp to zt(k) )
2182 : ! = [ ( zt(k) - zm(k-1) ) / ( zm(k) - zm(k-1) ) ]
2183 : ! * ( var_zm(k) - var_zm(k-1) ) + var_zm(k-1).
2184 : !
2185 : ! Furthermore, the formula can be rewritten as:
2186 : !
2187 : ! var_zm( interp to zt(k) )
2188 : ! = factor * var_zm(k) + ( 1 - factor ) * var_zm(k-1);
2189 : !
2190 : ! where:
2191 : !
2192 : ! factor = ( zt(k) - zm(k-1) ) / ( zm(k) - zm(k-1) ).
2193 : !
2194 : ! One of the important uses of this function is in situations where the
2195 : ! variable to be interpolated is being treated IMPLICITLY in an equation.
2196 : ! Usually, the variable to be interpolated is involved in a derivative (such
2197 : ! as d(var_zm)/dz in the diagram below). For the term of the equation
2198 : ! containing the derivative, grid weights are needed for two interpolations,
2199 : ! rather than just one interpolation. Thus, four grid weights (labeled
2200 : ! A(k), B(k), C(k), and D(k) in the diagram below) are needed.
2201 : !
2202 : ! ===var_zm(k+1)=========================================== m(k+1)
2203 : ! A(k)
2204 : ! -----------var_zm(interp)-------------------------------- t(k+1)
2205 : ! B(k) = 1 - A(k)
2206 : ! ===var_zm(k)===========d(var_zm)/dz====================== m(k)
2207 : ! C(k)
2208 : ! -----------var_zm(interp)-------------------------------- t(k)
2209 : ! D(k) = 1 - C(k)
2210 : ! ===var_zm(k-1)=========================================== m(k-1)
2211 : !
2212 : ! The vertical indices m(k+1), t(k+1), m(k), t(k), and m(k-1) correspond
2213 : ! with altitudes zm(k+1), zt(k+1), zm(k), zt(k), and zm(k-1), respectively.
2214 : ! The letter "t" is used for thermodynamic levels and the letter "m" is used
2215 : ! for momentum levels.
2216 : !
2217 : ! The grid weights, indexed around the central momentum level (k), are
2218 : ! defined as follows:
2219 : !
2220 : ! A(k) = ( zt(k+1) - zm(k) ) / ( zm(k+1) - zm(k) );
2221 : !
2222 : ! which is the same as "factor" for the interpolation to thermodynamic
2223 : ! level (k+1). In the code, this interpolation is referenced as
2224 : ! gr%weights_zm2zt(m_above,tkp1), which can be read as "grid weight in a
2225 : ! zm2zt interpolation of the momentum level above thermodynamic
2226 : ! level (k+1) (applied to thermodynamic level (k+1))".
2227 : !
2228 : ! B(k) = 1 - [ ( zt(k+1) - zm(k) ) / ( zm(k+1) - zm(k) ) ]
2229 : ! = 1 - A(k);
2230 : !
2231 : ! which is the same as "1 - factor" for the interpolation to thermodynamic
2232 : ! level (k+1). In the code, this interpolation is referenced as
2233 : ! gr%weights_zm2zt(m_below,tkp1), which can be read as "grid weight in a
2234 : ! zm2zt interpolation of the momentum level below thermodynamic
2235 : ! level (k+1) (applied to thermodynamic level (k+1))".
2236 : !
2237 : ! C(k) = ( zt(k) - zm(k-1) ) / ( zm(k) - zm(k-1) );
2238 : !
2239 : ! which is the same as "factor" for the interpolation to thermodynamic
2240 : ! level (k). In the code, this interpolation is referenced as
2241 : ! gr%weights_zm2zt(m_above,tk), which can be read as "grid weight in a zm2zt
2242 : ! interpolation of the momentum level above thermodynamic level (k) (applied
2243 : ! to thermodynamic level (k))".
2244 : !
2245 : ! D(k) = 1 - [ ( zt(k) - zm(k-1) ) / ( zm(k) - zm(k-1) ) ]
2246 : ! = 1 - C(k);
2247 : !
2248 : ! which is the same as "1 - factor" for the interpolation to thermodynamic
2249 : ! level (k). In the code, this interpolation is referenced as
2250 : ! gr%weights_zm2zt(m_below,tk), which can be read as "grid weight in a zm2zt
2251 : ! interpolation of the momentum level below thermodynamic level (k) (applied
2252 : ! to thermodynamic level (k))".
2253 : !
2254 : ! Additionally, as long as the central momentum level (k) in the above
2255 : ! scenario is not the lowermost momentum level or the uppermost momentum
2256 : ! level (k /= 1 and k /= gr%nz), the four weighting factors have the
2257 : ! following relationships: A(k) = C(k+1) and B(k) = D(k+1).
2258 : !
2259 : !
2260 : ! Special condition for lowermost grid level, k = 1:
2261 : !
2262 : ! The lowermost thermodynamic grid level is below the lowermost momentum
2263 : ! grid level. Thus, a linear extension is used at this level. It should
2264 : ! be noted that the thermodynamic level k = 1 is considered to be below the
2265 : ! model lower boundary, which is defined to be at momentum level k = 1.
2266 : ! Thus, the values of most variables at thermodynamic level k = 1 are not
2267 : ! often needed or referenced.
2268 : !
2269 : ! For level k = 1:
2270 : !
2271 : ! The formula for a linear extension is given by:
2272 : !
2273 : ! var_zm( extend to zt(k) )
2274 : ! = [ ( var_zm(k+1) - var_zm(k) ) / ( zm(k+1) - zm(k) ) ]
2275 : ! * ( zt(k) - zm(k) ) + var_zm(k);
2276 : !
2277 : ! which can be rewritten as:
2278 : !
2279 : ! var_zm( extend to zt(k) )
2280 : ! = [ ( zt(k) - zm(k) ) / ( zm(k+1) - zm(k) ) ]
2281 : ! * ( var_zm(k+1) - var_zm(k) ) + var_zm(k).
2282 : !
2283 : ! Furthermore, the formula can be rewritten as:
2284 : !
2285 : ! var_zm( extend to zt(k) )
2286 : ! = factor * var_zm(k+1) + ( 1 - factor ) * var_zm(k);
2287 : !
2288 : ! where:
2289 : !
2290 : ! factor = ( zt(k) - zm(k) ) / ( zm(k+1) - zm(k) ).
2291 : !
2292 : ! Due to the fact that a linear extension is being used, the value of factor
2293 : ! will be less than 0. The weight of the upper momentum level, which is
2294 : ! momentum level k = 2, on thermodynamic level k = 1 equals the value of
2295 : ! factor. The weight of the lower momentum level, which is momentum level
2296 : ! k = 1, on thermodynamic level k = 1 equals 1 - factor, which is greater
2297 : ! than 1. However, the sum of the weights equals 1.
2298 : !
2299 : !
2300 : ! Brian Griffin; September 12, 2008.
2301 : !
2302 : !-----------------------------------------------------------------------
2303 :
2304 : use constants_clubb, only: &
2305 : one ! Constant(s)
2306 :
2307 : use clubb_precision, only: &
2308 : core_rknd ! Variable(s)
2309 :
2310 : implicit none
2311 :
2312 : ! Constant parameters
2313 : integer, parameter :: &
2314 : m_above = 1, & ! Upper momentum level.
2315 : m_below = 2 ! Lower momentum level.
2316 :
2317 : ! Input Variable
2318 : integer, intent(in) :: &
2319 : nz, &
2320 : ngrdcol
2321 :
2322 : ! Input/Output Variable
2323 : type (grid), target, intent(inout) :: gr
2324 :
2325 : integer :: i, k
2326 :
2327 :
2328 : ! At most levels, the momentum level is found in-between two
2329 : ! thermodynamic levels. Linear interpolation is used.
2330 74597328 : do i = 1, ngrdcol
2331 0 : gr%weights_zm2zt(i,1,m_above) = ( gr%zt(i,1) - gr%zm(i,1) ) &
2332 70129800 : / ( gr%zm(i,2) - gr%zm(i,1) )
2333 :
2334 74597328 : gr%weights_zm2zt(i,1,m_below) = one - gr%weights_zm2zt(i,1,m_above)
2335 : end do
2336 :
2337 379739880 : do k = 2, nz
2338 6270643080 : do i = 1, ngrdcol
2339 :
2340 : ! The top model level (gr%nz) is formulated differently because the top
2341 : ! momentum level is above the top thermodynamic level. A linear
2342 : ! extension is required, rather than linear interpolation.
2343 : ! Note: Variable "factor" will be greater than 1 in this situation.
2344 0 : gr%weights_zm2zt(i,k,m_above) = ( gr%zt(i,k) - gr%zm(i,k-1) ) &
2345 5890903200 : / ( gr%zm(i,k) - gr%zm(i,k-1) )
2346 :
2347 : ! Weight of lower thermodynamic level on momentum level.
2348 6266175552 : gr%weights_zm2zt(i,k,m_below) = one - gr%weights_zm2zt(i,k,m_above)
2349 :
2350 : end do
2351 : end do
2352 :
2353 4467528 : return
2354 :
2355 : end subroutine calc_zm2zt_weights
2356 :
2357 : !=============================================================================
2358 : ! Wrapped in interface ddzm
2359 0 : function gradzm_2D( nz, ngrdcol, gr, azm )
2360 :
2361 : ! Description:
2362 : ! 2D version of gradzm
2363 : !-----------------------------------------------------------------------
2364 :
2365 : use clubb_precision, only: &
2366 : core_rknd ! Variable(s)
2367 :
2368 : implicit none
2369 :
2370 : integer, intent(in) :: &
2371 : nz, &
2372 : ngrdcol
2373 :
2374 : type (grid), target, intent(in) :: gr
2375 :
2376 : ! Input Variable
2377 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
2378 : azm ! Variable on momentum grid levels [units vary]
2379 :
2380 : ! Return Variable
2381 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
2382 : gradzm_2D ! Vertical derivative of azm [units vary / m]
2383 :
2384 : ! Local Variable
2385 : integer :: i, k ! Grid level loop index
2386 :
2387 : !$acc data copyin( gr, gr%invrs_dzt, azm ) &
2388 : !$acc copyout( gradzm_2D )
2389 :
2390 : !$acc parallel loop gang vector default(present)
2391 0 : do i = 1, ngrdcol
2392 0 : gradzm_2D(i,1) = ( azm(i,2) - azm(i,1) ) * gr%invrs_dzt(i,2)
2393 : end do
2394 : !$acc end parallel loop
2395 :
2396 : !$acc parallel loop gang vector collapse(2) default(present)
2397 0 : do k = 2, nz
2398 0 : do i = 1, ngrdcol
2399 0 : gradzm_2D(i,k) = ( azm(i,k) - azm(i,k-1) ) * gr%invrs_dzt(i,k)
2400 : end do
2401 : end do
2402 : !$acc end parallel loop
2403 :
2404 : !$acc end data
2405 :
2406 0 : return
2407 :
2408 : end function gradzm_2D
2409 :
2410 : !=============================================================================
2411 : ! Wrapped in interface ddzm
2412 0 : function gradzm_1D( gr, azm )
2413 :
2414 : ! Description:
2415 : ! 2D version of gradzm
2416 : !-----------------------------------------------------------------------
2417 :
2418 : use clubb_precision, only: &
2419 : core_rknd ! Variable(s)
2420 :
2421 : implicit none
2422 :
2423 : type (grid), target, intent(in) :: gr
2424 :
2425 : ! Input Variable
2426 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
2427 : azm ! Variable on momentum grid levels [units vary]
2428 :
2429 : ! Return Variable
2430 : real( kind = core_rknd ), dimension(gr%nz) :: &
2431 : gradzm_1D ! Vertical derivative of azm [units vary / m]
2432 :
2433 : ! Local Variables
2434 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
2435 0 : azm_col ! Variable on momentum grid levels [units vary]
2436 :
2437 : ! Return Variable
2438 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
2439 0 : gradzm_1D_col ! Vertical derivative of azm [units vary / m]
2440 :
2441 0 : azm_col(1,:) = azm
2442 :
2443 0 : gradzm_1D_col = gradzm_2D(gr%nz, 1, gr, azm_col)
2444 :
2445 0 : gradzm_1D = gradzm_1D_col(1,:)
2446 :
2447 0 : return
2448 :
2449 : end function gradzm_1D
2450 :
2451 : !=============================================================================
2452 : ! Wrapped in interface ddzt
2453 250181568 : function gradzt_2D( nz, ngrdcol, gr, azt )
2454 :
2455 : ! Description:
2456 : ! 2D version of gradzt
2457 : !-----------------------------------------------------------------------
2458 :
2459 : use clubb_precision, only: &
2460 : core_rknd ! Variable(s)
2461 :
2462 : implicit none
2463 :
2464 : integer, intent(in) :: &
2465 : nz, &
2466 : ngrdcol
2467 :
2468 : type (grid), target, intent(in) :: gr
2469 :
2470 : ! Input Variable
2471 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
2472 : azt ! Variable on thermodynamic grid levels [units vary]
2473 :
2474 : ! Output Variable
2475 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
2476 : gradzt_2D ! Vertical derivative of azt [units vary / m]
2477 :
2478 : ! Local Variable
2479 : integer :: i, k ! Grid level loop index
2480 :
2481 : !$acc data copyin( gr, gr%invrs_dzm, azt ) &
2482 : !$acc copyout( gradzt_2D )
2483 :
2484 : !$acc parallel loop gang vector collapse(2) default(present)
2485 21015251712 : do k = 2, nz-1
2486 >34697*10^7 : do i = 1, ngrdcol
2487 >34672*10^7 : gradzt_2D(i,k) = ( azt(i,k+1) - azt(i,k) ) * gr%invrs_dzm(i,k)
2488 : end do
2489 : end do
2490 : !$acc end parallel loop
2491 :
2492 : !$acc parallel loop gang vector default(present)
2493 4177450368 : do i = 1, ngrdcol
2494 3927268800 : gradzt_2D(i,1) = gradzt_2D(i,2)
2495 4177450368 : gradzt_2D(i,nz) = gradzt_2D(i,nz-1)
2496 : end do
2497 : !$acc end parallel loop
2498 :
2499 : !$acc end data
2500 :
2501 250181568 : return
2502 :
2503 : end function gradzt_2D
2504 :
2505 : !=============================================================================
2506 : ! Wrapped in interface ddzt
2507 0 : function gradzt_1D( gr, azt )
2508 :
2509 : ! Description:
2510 : ! 2D version of gradzt
2511 : !-----------------------------------------------------------------------
2512 :
2513 : use clubb_precision, only: &
2514 : core_rknd ! Variable(s)
2515 :
2516 : implicit none
2517 :
2518 : type (grid), target, intent(in) :: gr
2519 :
2520 : ! Input Variable
2521 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
2522 : azt ! Variable on thermodynamic grid levels [units vary]
2523 :
2524 : ! Output Variable
2525 : real( kind = core_rknd ), dimension(gr%nz) :: &
2526 : gradzt_1D ! Vertical derivative of azt [units vary / m]
2527 :
2528 : ! Local Variables
2529 :
2530 : ! Input Variable
2531 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
2532 0 : azt_col ! Variable on thermodynamic grid levels [units vary]
2533 :
2534 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
2535 0 : gradzt_1D_col ! Vertical derivative of azt [units vary / m]
2536 :
2537 0 : azt_col(1,:) = azt
2538 :
2539 0 : gradzt_1D_col = gradzt_2D(gr%nz, 1, gr, azt_col)
2540 :
2541 0 : gradzt_1D = gradzt_1D_col(1,:)
2542 :
2543 0 : return
2544 :
2545 : end function gradzt_1D
2546 :
2547 : !=============================================================================
2548 0 : function flip( x, xdim )
2549 :
2550 : ! Description:
2551 : ! Flips a single dimension array (i.e. a vector), so the first element
2552 : ! becomes the last and vice versa for the whole column. This is a
2553 : ! necessary part of the code because BUGSrad and CLUBB store altitudes in
2554 : ! reverse order.
2555 : !
2556 : ! References:
2557 : ! None
2558 : !-------------------------------------------------------------------------
2559 :
2560 : use clubb_precision, only: &
2561 : dp ! double precision
2562 :
2563 : implicit none
2564 :
2565 : ! Input Variables
2566 : integer, intent(in) :: xdim
2567 :
2568 : real(kind = dp), dimension(xdim), intent(in) :: x
2569 :
2570 : ! Output Variables
2571 : real(kind = dp), dimension(xdim) :: flip
2572 :
2573 : ! Local Variables
2574 0 : real(kind = dp), dimension(xdim) :: tmp
2575 :
2576 : integer :: indx
2577 :
2578 :
2579 : ! Get rid of an annoying compiler warning.
2580 : indx = 1
2581 : indx = indx
2582 :
2583 0 : forall ( indx = 1 : xdim )
2584 0 : tmp(indx) = x((xdim+1) - (indx))
2585 : end forall
2586 :
2587 0 : flip = tmp
2588 :
2589 :
2590 0 : return
2591 :
2592 : end function flip
2593 :
2594 : !===============================================================================
2595 :
2596 0 : end module grid_class
|