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 176472 : subroutine setup_grid( nzmax, ngrdcol, sfc_elevation, l_implemented, &
256 176472 : grid_type, deltaz, zm_init, zm_top, &
257 176472 : 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 176472 : 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 176472 : gr%nz = nzmax
351 :
352 : ! Default bounds
353 176472 : begin_height = 1
354 176472 : end_height = gr%nz
355 :
356 : !---------------------------------------------------
357 176472 : 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 3529440 : stat=ierr )
477 :
478 176472 : 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 176472 : gr ) ! intent(inout)
493 :
494 2946672 : do i = 1, ngrdcol
495 2946672 : 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 176472 : subroutine setup_grid_heights( nz, ngrdcol, &
550 : l_implemented, grid_type, &
551 176472 : deltaz, zm_init, momentum_heights, &
552 176472 : 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 176472 : 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 15176592 : do k = 1, nz, 1
738 250643592 : do i = 1, ngrdcol
739 250467120 : 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 15176592 : do k = 1, nz, 1
746 250643592 : do i = 1, ngrdcol
747 250467120 : 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 15000120 : do k=1,nz-1
758 247696920 : do i = 1, ngrdcol
759 247520448 : gr%dzm(i,k) = gr%zt(i,k+1) - gr%zt(i,k)
760 : end do
761 : enddo
762 :
763 2946672 : do i = 1, ngrdcol
764 2946672 : 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 15000120 : do k=2,nz
770 247696920 : do i = 1, ngrdcol
771 247520448 : gr%dzt(i,k) = gr%zm(i,k) - gr%zm(i,k-1)
772 : end do
773 : enddo
774 :
775 2946672 : do i = 1, ngrdcol
776 2946672 : 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 15000120 : do k=1,nz-1
782 247696920 : do i = 1, ngrdcol
783 247520448 : gr%invrs_dzm(i,k) = 1._core_rknd / ( gr%zt(i,k+1) - gr%zt(i,k) )
784 : end do
785 : enddo
786 :
787 2946672 : do i = 1, ngrdcol
788 2946672 : 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 15000120 : do k=2,nz
795 247696920 : do i = 1, ngrdcol
796 247520448 : gr%invrs_dzt(i,k) = 1._core_rknd / ( gr%zm(i,k) - gr%zm(i,k-1) )
797 : end do
798 : enddo
799 :
800 2946672 : do i = 1, ngrdcol
801 2946672 : 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 176472 : 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 176472 : gr ) ! InOut
834 :
835 176472 : 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 )
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 : !---------------------------- Local Variables ----------------------------
1181 : real( kind = core_rknd ), dimension(gr%nz) :: &
1182 0 : redirect_interpolated_azm_col ! Variable when interp. to momentum levels
1183 :
1184 : !---------------------------- Begin Code ----------------------------
1185 :
1186 0 : redirect_interpolated_azm_col = redirect_interpolated_azm_1D(gr, azt)
1187 :
1188 0 : redirect_interpolated_azm_k = redirect_interpolated_azm_col(k)
1189 :
1190 : return
1191 : end function redirect_interpolated_azm_k
1192 :
1193 : !=============================================================================
1194 : ! Wrapped in interface zt2zm
1195 0 : function redirect_interpolated_azm_1D( gr, azt )
1196 :
1197 : ! Description:
1198 : ! Calls the appropriate corresponding function based on l_cubic_temp
1199 : !-----------------------------------------------------------------------
1200 :
1201 : use clubb_precision, only: &
1202 : core_rknd ! Variable(s)
1203 :
1204 : use model_flags, only: &
1205 : l_cubic_interp, & ! Variable(s)
1206 : l_quintic_poly_interp
1207 :
1208 : use constants_clubb, only: &
1209 : fstdout ! Variable
1210 :
1211 : implicit none
1212 :
1213 : !---------------------------- Input Variables ----------------------------
1214 : type (grid), target, intent(in) :: gr
1215 :
1216 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
1217 : azt ! Variable on thermodynamic grid levels [units vary]
1218 :
1219 : !---------------------------- Return Variable ----------------------------
1220 : real( kind = core_rknd ), dimension(gr%nz) :: &
1221 : redirect_interpolated_azm_1D ! Variable when interp. to momentum levels
1222 :
1223 : !---------------------------- Local Variables ----------------------------
1224 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1225 0 : azt_col ! Variable on thermodynamic grid levels [units vary]
1226 :
1227 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1228 0 : redirect_interpolated_azm_col ! Variable when interp. to momentum levels
1229 :
1230 : !---------------------------- Begin Code ----------------------------
1231 :
1232 0 : azt_col(1,:) = azt
1233 :
1234 0 : redirect_interpolated_azm_col = redirect_interpolated_azm_2D(gr%nz, 1, gr, azt_col)
1235 :
1236 0 : redirect_interpolated_azm_1D = redirect_interpolated_azm_col(1,:)
1237 :
1238 0 : return
1239 : end function redirect_interpolated_azm_1D
1240 :
1241 : !=============================================================================
1242 : ! Wrapped in interface zt2zm
1243 15529536 : function redirect_interpolated_azm_2D( nz, ngrdcol, gr, azt )
1244 :
1245 : ! Description:
1246 : ! Calls the appropriate corresponding function based on l_cubic_temp
1247 : !-----------------------------------------------------------------------
1248 :
1249 : use clubb_precision, only: &
1250 : core_rknd ! Variable(s)
1251 :
1252 : use model_flags, only: &
1253 : l_cubic_interp, & ! Variable(s)
1254 : l_quintic_poly_interp
1255 :
1256 : use constants_clubb, only: &
1257 : fstdout ! Variable
1258 :
1259 : implicit none
1260 :
1261 : !---------------------------- Input Variables ----------------------------
1262 : integer, intent(in) :: &
1263 : nz, &
1264 : ngrdcol
1265 :
1266 : type (grid), target, intent(in) :: gr
1267 :
1268 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1269 : azt ! Variable on thermodynamic grid levels [units vary]
1270 :
1271 : !---------------------------- Return Variable ----------------------------
1272 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1273 : redirect_interpolated_azm_2D ! Variable when interp. to momentum levels
1274 :
1275 : !---------------------------- Begin Code ----------------------------
1276 :
1277 : ! Sanity Check
1278 15529536 : if (l_quintic_poly_interp) then
1279 : if (.not. l_cubic_interp) then
1280 : write (fstdout, *) "Error: Model flag l_quintic_poly_interp should not be true if "&
1281 0 : //"l_cubic_interp is false."
1282 0 : error stop
1283 : end if
1284 : end if
1285 :
1286 : ! Redirect
1287 : if ( l_cubic_interp .and. nz >= 3 ) then
1288 : redirect_interpolated_azm_2D = cubic_interpolated_azm_2D( nz, ngrdcol, gr, azt )
1289 : else
1290 : call linear_interpolated_azm_2D( nz, ngrdcol, gr, azt, &
1291 15529536 : redirect_interpolated_azm_2D )
1292 : end if
1293 :
1294 15529536 : return
1295 : end function redirect_interpolated_azm_2D
1296 :
1297 : !=============================================================================
1298 : ! Wrapped in interface zm2zt
1299 0 : function redirect_interpolated_azt_k( gr, azm, k )
1300 :
1301 : ! Description:
1302 : ! Calls the appropriate corresponding function based on l_cubic_temp
1303 : !-----------------------------------------------------------------------
1304 :
1305 : use clubb_precision, only: &
1306 : core_rknd ! Variable(s)
1307 :
1308 : use model_flags, only: &
1309 : l_cubic_interp, & ! Variable(s)
1310 : l_quintic_poly_interp
1311 :
1312 : use constants_clubb, only: &
1313 : fstdout ! Variable
1314 :
1315 : implicit none
1316 :
1317 : type (grid), target, intent(in) :: gr
1318 :
1319 : !---------------------------- Input Variables ----------------------------
1320 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
1321 : azm ! Variable on momentum grid levels [units vary]
1322 :
1323 : integer, intent(in) :: &
1324 : k ! Vertical level
1325 :
1326 : !---------------------------- Return Variable ----------------------------
1327 : real( kind = core_rknd ) :: &
1328 : redirect_interpolated_azt_k ! Variable when interp. to momentum levels
1329 :
1330 : !---------------------------- Local Variables ----------------------------
1331 :
1332 : real( kind = core_rknd ), dimension(gr%nz) :: &
1333 0 : redirect_interpolated_azt_col ! Variable when interp. to momentum levels
1334 :
1335 : !---------------------------- Begin Code ----------------------------
1336 :
1337 0 : redirect_interpolated_azt_col = redirect_interpolated_azt_1D(gr, azm)
1338 :
1339 0 : redirect_interpolated_azt_k = redirect_interpolated_azt_col(k)
1340 :
1341 : return
1342 : end function redirect_interpolated_azt_k
1343 :
1344 :
1345 : !=============================================================================
1346 : ! Wrapped in interface zm2zt
1347 0 : function redirect_interpolated_azt_1D( gr, azm )
1348 :
1349 : ! Description:
1350 : ! Calls the appropriate corresponding function based on l_cubic_temp
1351 : !-----------------------------------------------------------------------
1352 :
1353 : use clubb_precision, only: &
1354 : core_rknd ! Variable(s)
1355 :
1356 : use model_flags, only: &
1357 : l_cubic_interp, & ! Variable(s)
1358 : l_quintic_poly_interp
1359 :
1360 : use constants_clubb, only: &
1361 : fstdout ! Variable
1362 :
1363 : implicit none
1364 :
1365 : type (grid), target, intent(in) :: gr
1366 :
1367 : !---------------------------- Input Variables ----------------------------
1368 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
1369 : azm ! Variable on momentum grid levels [units vary]
1370 :
1371 : !---------------------------- Return Variable ----------------------------
1372 : real( kind = core_rknd ), dimension(gr%nz) :: &
1373 : redirect_interpolated_azt_1D ! Variable when interp. to momentum levels
1374 :
1375 : !---------------------------- Local Variables ----------------------------
1376 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1377 0 : azm_col ! Variable on thermodynamic grid levels [units vary]
1378 :
1379 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
1380 0 : redirect_interpolated_azt_col ! Variable when interp. to momentum levels
1381 :
1382 : !---------------------------- Begin Code ----------------------------
1383 :
1384 0 : azm_col(1,:) = azm
1385 :
1386 0 : redirect_interpolated_azt_col = redirect_interpolated_azt_2D(gr%nz, 1, gr, azm_col)
1387 :
1388 0 : redirect_interpolated_azt_1D = redirect_interpolated_azt_col(1,:)
1389 :
1390 0 : return
1391 : end function redirect_interpolated_azt_1D
1392 :
1393 : !=============================================================================
1394 : ! Wrapped in interface zm2zt
1395 16764840 : function redirect_interpolated_azt_2D( nz, ngrdcol, gr, azm )
1396 :
1397 : ! Description:
1398 : ! Calls the appropriate corresponding function based on l_cubic_temp
1399 : !-----------------------------------------------------------------------
1400 :
1401 : use clubb_precision, only: &
1402 : core_rknd ! Variable(s)
1403 :
1404 : use model_flags, only: &
1405 : l_cubic_interp, & ! Variable(s)
1406 : l_quintic_poly_interp
1407 :
1408 : use constants_clubb, only: &
1409 : fstdout ! Variable
1410 :
1411 : implicit none
1412 :
1413 : !---------------------------- Input Variables ----------------------------
1414 : integer, intent(in) :: &
1415 : nz, &
1416 : ngrdcol
1417 :
1418 : type (grid), target, intent(in) :: gr
1419 :
1420 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1421 : azm ! Variable on momentum grid levels [units vary]
1422 :
1423 : !---------------------------- Return Variable----------------------------
1424 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1425 : redirect_interpolated_azt_2D ! Variable when interp. to momentum levels
1426 :
1427 : !---------------------------- Begin Code ----------------------------
1428 :
1429 : ! Sanity Check
1430 16764840 : if (l_quintic_poly_interp) then
1431 : if (.not. l_cubic_interp) then
1432 : write (fstdout, *) "Error: Model flag l_quintic_poly_interp should not be true if "&
1433 0 : //"l_cubic_interp is false."
1434 0 : error stop
1435 : end if
1436 : end if
1437 :
1438 : ! Redirect
1439 : if (l_cubic_interp .and. nz >= 3 ) then
1440 : redirect_interpolated_azt_2D = cubic_interpolated_azt_2D( nz, ngrdcol, gr, azm )
1441 : else
1442 : call linear_interpolated_azt_2D( nz, ngrdcol, gr, azm, &
1443 16764840 : redirect_interpolated_azt_2D )
1444 : end if
1445 :
1446 16764840 : return
1447 : end function redirect_interpolated_azt_2D
1448 :
1449 : !=============================================================================
1450 15529536 : subroutine linear_interpolated_azm_2D( nz, ngrdcol, gr, azt, &
1451 15529536 : linear_interpolated_azm )
1452 :
1453 : ! Description:
1454 : ! Function to interpolate a variable located on the thermodynamic grid
1455 : ! levels (azt) to the momentum grid levels (azm). This function inputs the
1456 : ! entire azt array and outputs the results as an azm array. The
1457 : ! formulation used is compatible with a stretched (unevenly-spaced) grid.
1458 : !-----------------------------------------------------------------------
1459 :
1460 : use clubb_precision, only: &
1461 : core_rknd ! Variable(s)
1462 :
1463 : use interpolation, only: &
1464 : linear_interp_factor ! Procedure(s)
1465 :
1466 : implicit none
1467 :
1468 : integer, intent(in) :: &
1469 : nz, &
1470 : ngrdcol
1471 :
1472 : type (grid), target, intent(in) :: gr
1473 :
1474 : ! ------------------------------ Input Variable ------------------------------
1475 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1476 : azt ! Variable on thermodynamic grid levels [units vary]
1477 :
1478 : ! ------------------------------ Return Variable ------------------------------
1479 : real( kind = core_rknd ), intent(out), dimension(ngrdcol,nz) :: &
1480 : linear_interpolated_azm ! Variable when interp. to momentum levels
1481 :
1482 : ! ------------------------------ Local Variable ------------------------------
1483 : integer :: i, k ! Grid level loop index
1484 :
1485 : ! ------------------------------ Begin Code ------------------------------
1486 :
1487 : !$acc data copyin( azt, gr, gr%weights_zt2zm, gr%zt, gr%zm ) &
1488 : !$acc copyout( linear_interpolated_azm )
1489 :
1490 : ! Interpolate the value of a thermodynamic-level variable to the central
1491 : ! momentum level, k, between two successive thermodynamic levels using
1492 : ! linear interpolation.
1493 : !$acc parallel loop gang vector collapse(2) default(present)
1494 1320010560 : do k = 1, nz-1
1495 21797328960 : do i = 1, ngrdcol
1496 40954636800 : linear_interpolated_azm(i,k) = gr%weights_zt2zm(i,k,1) &
1497 62736436224 : * ( azt(i,k+1) - azt(i,k) ) + azt(i,k)
1498 : end do
1499 : end do
1500 : !$acc end parallel loop
1501 :
1502 : ! Set the value of the thermodynamic-level variable, azt, at the uppermost
1503 : ! level of the model, which is a momentum level. The name of the variable
1504 : ! when interpolated/extended to momentum levels is azm.
1505 : ! Use a linear extension based on the values of azt at levels gr%nz and
1506 : ! gr%nz-1 to find the value of azm at level gr%nz (the uppermost level
1507 : ! in the model).
1508 : !$acc parallel loop gang vector default(present)
1509 259307136 : do i = 1, ngrdcol
1510 487555200 : linear_interpolated_azm(i,nz) &
1511 487555200 : = ( ( azt(i,nz) - azt(i,nz-1) ) / ( gr%zt(i,nz) - gr%zt(i,nz-1) ) ) &
1512 990639936 : * ( gr%zm(i,nz) - gr%zt(i,nz) ) + azt(i,nz)
1513 : end do
1514 : !$acc end parallel loop
1515 :
1516 : !$acc end data
1517 :
1518 15529536 : return
1519 :
1520 : end subroutine linear_interpolated_azm_2D
1521 :
1522 : !=============================================================================
1523 0 : function zt2zm2zt( nz, ngrdcol, gr, azt )
1524 :
1525 : ! Description:
1526 : ! Function to interpolate a variable located on the thermodynamic grid
1527 : ! levels (azt) to the momentum grid levels (azm), then interpolate back
1528 : ! to thermodynamic grid levels (azt).
1529 : !
1530 : ! Note:
1531 : ! This is intended for smoothing variables.
1532 : !-----------------------------------------------------------------------------
1533 :
1534 : use clubb_precision, only: &
1535 : core_rknd ! Variable(s)
1536 :
1537 : implicit none
1538 :
1539 : ! ------------------------------ Input Variable ------------------------------
1540 : integer, intent(in) :: &
1541 : nz, &
1542 : ngrdcol
1543 :
1544 : type (grid), target, intent(in) :: gr
1545 :
1546 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1547 : azt ! Variable on thermodynamic grid levels [units vary]
1548 :
1549 : ! ------------------------------ Return Variable ------------------------------
1550 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1551 : zt2zm2zt ! Variable when interp. to momentum levels
1552 :
1553 : ! ------------------------------ Local Variable ------------------------------
1554 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1555 0 : azt_zm
1556 :
1557 : ! ------------------------------ Begin Code ------------------------------
1558 :
1559 : !$acc data copyin( azt ) &
1560 : !$acc copyout( zt2zm2zt ) &
1561 : !$acc create( azt_zm )
1562 :
1563 : ! Interpolate azt to momentum levels
1564 0 : azt_zm = zt2zm( nz, ngrdcol, gr, azt )
1565 :
1566 : ! Interpolate back to termodynamic levels
1567 0 : zt2zm2zt = zm2zt( nz, ngrdcol, gr, azt_zm )
1568 :
1569 : !$acc end data
1570 :
1571 0 : return
1572 :
1573 : end function zt2zm2zt
1574 :
1575 : !=============================================================================
1576 1411776 : function zm2zt2zm( nz, ngrdcol, gr, azm )
1577 :
1578 : ! Description:
1579 : ! Function to interpolate a variable located on the momentum grid
1580 : ! levels(azm) to thermodynamic grid levels (azt), then interpolate
1581 : ! back to momentum grid levels (azm).
1582 : !
1583 : ! Note:
1584 : ! This is intended for smoothing variables.
1585 : !-----------------------------------------------------------------------------
1586 :
1587 : use clubb_precision, only: &
1588 : core_rknd ! Variable(s)
1589 :
1590 : implicit none
1591 :
1592 : ! ------------------------------ Input Variable ------------------------------
1593 : integer, intent(in) :: &
1594 : nz, &
1595 : ngrdcol
1596 :
1597 : type (grid), target, intent(in) :: gr
1598 :
1599 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1600 : azm ! Variable on momentum grid levels [units vary]
1601 :
1602 : ! ------------------------------ Return Variable ------------------------------
1603 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1604 : zm2zt2zm ! Variable when interp. to momentum levels
1605 :
1606 : ! ------------------------------ Local Variable ------------------------------
1607 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1608 1411776 : azm_zt
1609 :
1610 : ! ------------------------------ Begin Code ------------------------------
1611 :
1612 : !$acc data copyin( azm ) &
1613 : !$acc copyout( zm2zt2zm ) &
1614 : !$acc create( azm_zt )
1615 :
1616 : ! Interpolate azt to termodynamic levels
1617 1411776 : azm_zt = zm2zt( nz, ngrdcol, gr, azm )
1618 :
1619 : ! Interpolate back to momentum levels
1620 1411776 : zm2zt2zm = zt2zm( nz, ngrdcol, gr, azm_zt )
1621 :
1622 : !$acc end data
1623 :
1624 1411776 : return
1625 :
1626 : end function zm2zt2zm
1627 :
1628 : !=============================================================================
1629 : function cubic_interpolated_azm_2D( nz, ngrdcol, gr, azt )
1630 :
1631 : ! Description:
1632 : ! Function to interpolate a variable located on the momentum grid
1633 : ! levels (azt) to the thermodynamic grid levels (azm). This function outputs the
1634 : ! value of azt at a all grid levels using Steffen's monotonic cubic
1635 : ! interpolation implemented by Tak Yamaguchi.
1636 : !
1637 : ! References:
1638 : ! None
1639 : !-----------------------------------------------------------------------
1640 :
1641 : use interpolation, only: &
1642 : mono_cubic_interp ! Procedure(s)
1643 :
1644 : use clubb_precision, only: &
1645 : core_rknd ! Variable(s)
1646 :
1647 : implicit none
1648 :
1649 : ! Input Variables
1650 : integer, intent(in) :: &
1651 : nz, &
1652 : ngrdcol
1653 :
1654 : type (grid), target, intent(in) :: gr
1655 :
1656 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1657 : azt
1658 :
1659 : ! Return Variable
1660 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
1661 : cubic_interpolated_azm_2D
1662 :
1663 : ! Local Variable(s)
1664 : integer :: &
1665 : k, i
1666 :
1667 : integer :: km1, k00, kp1, kp2
1668 :
1669 : ! ---- Begin Code ----
1670 :
1671 : do k = 1, nz
1672 : do i = 1, ngrdcol
1673 :
1674 : ! k levels are based on Tak's find_indices subroutine -dschanen 24 Oct 2011
1675 : if ( k == nz-1 ) then
1676 : km1 = nz-2
1677 : kp1 = nz
1678 : kp2 = nz
1679 : k00 = nz-1
1680 : else if ( k == nz ) then ! Extrapolation
1681 : km1 = nz
1682 : kp1 = nz
1683 : kp2 = nz
1684 : k00 = nz-1
1685 : else if ( k == 1 ) then
1686 : km1 = 1
1687 : kp1 = 2
1688 : kp2 = 3
1689 : k00 = 1
1690 : else
1691 : km1 = k-1
1692 : kp1 = k+1
1693 : kp2 = k+2
1694 : k00 = k
1695 : end if
1696 :
1697 : ! Do the actual interpolation.
1698 : ! Use a cubic monotonic spline interpolation.
1699 : cubic_interpolated_azm_2D(i,k) = &
1700 : mono_cubic_interp( gr%zm(i,k), km1, k00, kp1, kp2, &
1701 : gr%zt(i,km1), gr%zt(i,k00), gr%zt(i,kp1), gr%zt(i,kp2), &
1702 : azt(i,km1), azt(i,k00), azt(i,kp1), azt(i,kp2) )
1703 : end do
1704 : end do
1705 :
1706 : return
1707 :
1708 : end function cubic_interpolated_azm_2D
1709 :
1710 : !=============================================================================
1711 176472 : subroutine calc_zt2zm_weights( nz, ngrdcol, &
1712 : gr )
1713 :
1714 : ! Description:
1715 : ! Function used to help in an interpolation of a variable (var_zt) located
1716 : ! on the thermodynamic grid levels (azt) to the momentum grid levels (azm).
1717 : ! This function computes a weighting factor for both the upper thermodynamic
1718 : ! level (k+1) and the lower thermodynamic level (k) applied to the central
1719 : ! momentum level (k). For the uppermost momentum grid level (k=gr%nz), a
1720 : ! weighting factor for both the thermodynamic level at gr%nz and the
1721 : ! thermodynamic level at gr%nz-1 are calculated based on the use of a
1722 : ! linear extension. This function outputs the weighting factors at a single
1723 : ! momentum grid level (k). The formulation used is compatible with a
1724 : ! stretched (unevenly-spaced) grid. The weights are defined as follows:
1725 : !
1726 : ! ---var_zt(k+1)------------------------------------------- t(k+1)
1727 : ! azt_weight(t_above) = factor
1728 : ! ===========var_zt(interp)================================ m(k)
1729 : ! azt_weight(t_below) = 1 - factor
1730 : ! ---var_zt(k)--------------------------------------------- t(k)
1731 : !
1732 : ! The vertical indices t(k+1), m(k), and t(k) correspond with altitudes
1733 : ! zt(k+1), zm(k), and zt(k), respectively. The letter "t" is used for
1734 : ! thermodynamic levels and the letter "m" is used for momentum levels.
1735 : !
1736 : ! For all levels k < gr%nz:
1737 : !
1738 : ! The formula for a linear interpolation is given by:
1739 : !
1740 : ! var_zt( interp to zm(k) )
1741 : ! = [ ( var_zt(k+1) - var_zt(k) ) / ( zt(k+1) - zt(k) ) ]
1742 : ! * ( zm(k) - zt(k) ) + var_zt(k);
1743 : !
1744 : ! which can be rewritten as:
1745 : !
1746 : ! var_zt( interp to zm(k) )
1747 : ! = [ ( zm(k) - zt(k) ) / ( zt(k+1) - zt(k) ) ]
1748 : ! * ( var_zt(k+1) - var_zt(k) ) + var_zt(k).
1749 : !
1750 : ! Furthermore, the formula can be rewritten as:
1751 : !
1752 : ! var_zt( interp to zm(k) )
1753 : ! = factor * var_zt(k+1) + ( 1 - factor ) * var_zt(k);
1754 : !
1755 : ! where:
1756 : !
1757 : ! factor = ( zm(k) - zt(k) ) / ( zt(k+1) - zt(k) ).
1758 : !
1759 : ! One of the important uses of this function is in situations where the
1760 : ! variable to be interpolated is being treated IMPLICITLY in an equation.
1761 : ! Usually, the variable to be interpolated is involved in a derivative (such
1762 : ! as d(var_zt)/dz in the diagram below). For the term of the equation
1763 : ! containing the derivative, grid weights are needed for two interpolations,
1764 : ! rather than just one interpolation. Thus, four grid weights (labeled
1765 : ! A(k), B(k), C(k), and D(k) in the diagram below) are needed.
1766 : !
1767 : ! ---var_zt(k+1)------------------------------------------- t(k+1)
1768 : ! A(k)
1769 : ! ===========var_zt(interp)================================ m(k)
1770 : ! B(k) = 1 - A(k)
1771 : ! ---var_zt(k)-----------d(var_zt)/dz---------------------- t(k)
1772 : ! C(k)
1773 : ! ===========var_zt(interp)================================ m(k-1)
1774 : ! D(k) = 1 - C(k)
1775 : ! ---var_zt(k-1)------------------------------------------- t(k-1)
1776 : !
1777 : ! The vertical indices t(k+1), m(k), t(k), m(k-1), and t(k-1) correspond
1778 : ! with altitudes zt(k+1), zm(k), zt(k), zm(k-1), and zt(k-1), respectively.
1779 : ! The letter "t" is used for thermodynamic levels and the letter "m" is used
1780 : ! for momentum levels.
1781 : !
1782 : ! The grid weights, indexed around the central thermodynamic level (k), are
1783 : ! defined as follows:
1784 : !
1785 : ! A(k) = ( zm(k) - zt(k) ) / ( zt(k+1) - zt(k) );
1786 : !
1787 : ! which is the same as "factor" for the interpolation to momentum
1788 : ! level (k). In the code, this interpolation is referenced as
1789 : ! gr%weights_zt2zm(t_above,mk), which can be read as "grid weight in a zt2zm
1790 : ! interpolation of the thermodynamic level above momentum level (k) (applied
1791 : ! to momentum level (k))".
1792 : !
1793 : ! B(k) = 1 - [ ( zm(k) - zt(k) ) / ( zt(k+1) - zt(k) ) ]
1794 : ! = 1 - A(k);
1795 : !
1796 : ! which is the same as "1 - factor" for the interpolation to momentum
1797 : ! level (k). In the code, this interpolation is referenced as
1798 : ! gr%weights_zt2zm(t_below,mk), which can be read as "grid weight in a zt2zm
1799 : ! interpolation of the thermodynamic level below momentum level (k) (applied
1800 : ! to momentum level (k))".
1801 : !
1802 : ! C(k) = ( zm(k-1) - zt(k-1) ) / ( zt(k) - zt(k-1) );
1803 : !
1804 : ! which is the same as "factor" for the interpolation to momentum
1805 : ! level (k-1). In the code, this interpolation is referenced as
1806 : ! gr%weights_zt2zm(t_above,mkm1), which can be read as "grid weight in a
1807 : ! zt2zm interpolation of the thermodynamic level above momentum level (k-1)
1808 : ! (applied to momentum level (k-1))".
1809 : !
1810 : ! D(k) = 1 - [ ( zm(k-1) - zt(k-1) ) / ( zt(k) - zt(k-1) ) ]
1811 : ! = 1 - C(k);
1812 : !
1813 : ! which is the same as "1 - factor" for the interpolation to momentum
1814 : ! level (k-1). In the code, this interpolation is referenced as
1815 : ! gr%weights_zt2zm(t_below,mkm1), which can be read as "grid weight in a
1816 : ! zt2zm interpolation of the thermodynamic level below momentum level (k-1)
1817 : ! (applied to momentum level (k-1))".
1818 : !
1819 : ! Additionally, as long as the central thermodynamic level (k) in the above
1820 : ! scenario is not the uppermost thermodynamic level or the lowermost
1821 : ! thermodynamic level (k /= gr%nz and k /= 1), the four weighting factors
1822 : ! have the following relationships: A(k) = C(k+1) and B(k) = D(k+1).
1823 : !
1824 : !
1825 : ! Special condition for uppermost grid level, k = gr%nz:
1826 : !
1827 : ! The uppermost momentum grid level is above the uppermost thermodynamic
1828 : ! grid level. Thus, a linear extension is used at this level.
1829 : !
1830 : ! For level k = gr%nz:
1831 : !
1832 : ! The formula for a linear extension is given by:
1833 : !
1834 : ! var_zt( extend to zm(k) )
1835 : ! = [ ( var_zt(k) - var_zt(k-1) ) / ( zt(k) - zt(k-1) ) ]
1836 : ! * ( zm(k) - zt(k-1) ) + var_zt(k-1);
1837 : !
1838 : ! which can be rewritten as:
1839 : !
1840 : ! var_zt( extend to zm(k) )
1841 : ! = [ ( zm(k) - zt(k-1) ) / ( zt(k) - zt(k-1) ) ]
1842 : ! * ( var_zt(k) - var_zt(k-1) ) + var_zt(k-1).
1843 : !
1844 : ! Furthermore, the formula can be rewritten as:
1845 : !
1846 : ! var_zt( extend to zm(k) )
1847 : ! = factor * var_zt(k) + ( 1 - factor ) * var_zt(k-1);
1848 : !
1849 : ! where:
1850 : !
1851 : ! factor = ( zm(k) - zt(k-1) ) / ( zt(k) - zt(k-1) ).
1852 : !
1853 : ! Due to the fact that a linear extension is being used, the value of factor
1854 : ! will be greater than 1. The weight of thermodynamic level k = gr%nz on
1855 : ! momentum level k = gr%nz equals the value of factor. The weight of
1856 : ! thermodynamic level k = gr%nz-1 on momentum level k = gr%nz equals
1857 : ! 1 - factor, which is less than 0. However, the sum of the two weights
1858 : ! equals 1.
1859 : !
1860 : !
1861 : ! Brian Griffin; September 12, 2008.
1862 : !
1863 : !-----------------------------------------------------------------------
1864 :
1865 : use constants_clubb, only: &
1866 : one ! Constant(s)
1867 :
1868 : use clubb_precision, only: &
1869 : core_rknd ! Variable(s)
1870 :
1871 : implicit none
1872 :
1873 : ! Constant parameters
1874 : integer, parameter :: &
1875 : t_above = 1, & ! Upper thermodynamic level.
1876 : t_below = 2 ! Lower thermodynamic level.
1877 :
1878 : ! Input Variable
1879 : integer, intent(in) :: &
1880 : nz, &
1881 : ngrdcol
1882 :
1883 : ! Input/Output Variable
1884 : type (grid), target, intent(inout) :: gr
1885 :
1886 : ! Local Variables
1887 : real( kind = core_rknd ) :: factor
1888 :
1889 : integer :: i, k
1890 :
1891 15000120 : do k = 1, nz-1
1892 247696920 : do i = 1, ngrdcol
1893 :
1894 : ! The top model level (gr%nz) is formulated differently because the top
1895 : ! momentum level is above the top thermodynamic level. A linear
1896 : ! extension is required, rather than linear interpolation.
1897 : ! Note: Variable "factor" will be greater than 1 in this situation.
1898 0 : gr%weights_zt2zm(i,k,t_above) = ( gr%zm(i,k) - gr%zt(i,k) ) &
1899 232696800 : / ( gr%zt(i,k+1) - gr%zt(i,k) )
1900 :
1901 : ! Weight of lower thermodynamic level on momentum level.
1902 247520448 : gr%weights_zt2zm(i,k,t_below) = one - gr%weights_zt2zm(i,k,t_above)
1903 :
1904 : end do
1905 : end do
1906 :
1907 : ! At most levels, the momentum level is found in-between two
1908 : ! thermodynamic levels. Linear interpolation is used.
1909 2946672 : do i = 1, ngrdcol
1910 0 : gr%weights_zt2zm(i,nz,t_above) = ( gr%zm(i,nz) - gr%zt(i,nz-1) ) &
1911 2770200 : / ( gr%zt(i,nz) - gr%zt(i,nz-1) )
1912 :
1913 2946672 : gr%weights_zt2zm(i,nz,t_below) = one - gr%weights_zt2zm(i,nz,t_above)
1914 : end do
1915 :
1916 176472 : return
1917 :
1918 : end subroutine calc_zt2zm_weights
1919 :
1920 : !=============================================================================
1921 16764840 : subroutine linear_interpolated_azt_2D( nz, ngrdcol, gr, azm, &
1922 16764840 : linear_interpolated_azt )
1923 :
1924 : ! Description:
1925 : ! Function to interpolate a variable located on the momentum grid levels
1926 : ! (azm) to the thermodynamic grid levels (azt). This function inputs the
1927 : ! entire azm array and outputs the results as an azt array. The formulation
1928 : ! used is compatible with a stretched (unevenly-spaced) grid.
1929 : !-----------------------------------------------------------------------
1930 :
1931 : use clubb_precision, only: &
1932 : core_rknd ! Variable(s)
1933 :
1934 : use interpolation, only: &
1935 : linear_interp_factor ! Procedure(s)
1936 :
1937 : implicit none
1938 :
1939 : integer, intent(in) :: &
1940 : nz, &
1941 : ngrdcol
1942 :
1943 : type (grid), target, intent(in) :: gr
1944 :
1945 : ! ------------------------------ Input Variable ------------------------------
1946 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
1947 : azm ! Variable on momentum grid levels [units vary]
1948 :
1949 : ! ------------------------------ Output Variable ------------------------------
1950 : real( kind = core_rknd ), intent(out), dimension(ngrdcol,nz) :: &
1951 : linear_interpolated_azt ! Variable when interp. to thermodynamic levels
1952 :
1953 : ! ------------------------------ Local Variable ------------------------------
1954 : integer :: i, k ! Grid level loop index
1955 :
1956 : ! ------------------------------ Begin Code ------------------------------
1957 :
1958 : !$acc data copyin( azm, gr, gr%weights_zm2zt, gr%zt, gr%zm ) &
1959 : !$acc copyout( linear_interpolated_azt )
1960 :
1961 : ! Set the value of the momentum-level variable, azm, at the lowermost level
1962 : ! of the model (below the model lower boundary), which is a thermodynamic
1963 : ! level. The name of the variable when interpolated/extended to
1964 : ! thermodynamic levels is azt.
1965 : ! Use a linear extension based on the values of azm at levels 1 and 2 to
1966 : ! find the value of azt at level 1 (the lowermost level in the model).
1967 : !$acc parallel loop gang vector default(present)
1968 279933840 : do i = 1, ngrdcol
1969 526338000 : linear_interpolated_azt(i,1) &
1970 526338000 : = ( ( azm(i,2) - azm(i,1) ) / ( gr%zm(i,2) - gr%zm(i,1) ) ) &
1971 806271840 : * ( gr%zt(i,1) - gr%zm(i,1) ) + azm(i,1)
1972 : end do
1973 : !$acc end parallel loop
1974 :
1975 : ! Interpolate the value of a momentum-level variable to the central
1976 : ! thermodynamic level, k, between two successive momentum levels using
1977 : ! linear interpolation.
1978 : !$acc parallel loop gang vector collapse(2) default(present)
1979 1425011400 : do k = 2, nz
1980 23531207400 : do i = 1, ngrdcol
1981 44212392000 : linear_interpolated_azt(i,k) = gr%weights_zm2zt(i,k,1) &
1982 67726834560 : * ( azm(i,k) - azm(i,k-1) ) + azm(i,k-1)
1983 : end do
1984 : end do
1985 : !$acc end parallel loop
1986 :
1987 : !$acc end data
1988 :
1989 16764840 : return
1990 :
1991 : end subroutine linear_interpolated_azt_2D
1992 :
1993 : !=============================================================================
1994 : function cubic_interpolated_azt_2D( nz, ngrdcol, gr, azm )
1995 :
1996 : ! Description:
1997 : ! Function to interpolate a variable located on the momentum grid
1998 : ! levels (azm) to the thermodynamic grid levels (azt). This function outputs the
1999 : ! value of azt at a all grid levels using Steffen's monotonic cubic
2000 : ! interpolation implemented by Tak Yamaguchi.
2001 : !
2002 : ! References:
2003 : ! None
2004 : !-----------------------------------------------------------------------
2005 :
2006 : use interpolation, only: &
2007 : mono_cubic_interp ! Procedure(s)
2008 :
2009 : use clubb_precision, only: &
2010 : core_rknd ! Variable(s)
2011 :
2012 : implicit none
2013 :
2014 : ! Input Variables
2015 : integer, intent(in) :: &
2016 : nz, &
2017 : ngrdcol
2018 :
2019 : type (grid), target, intent(in) :: gr
2020 :
2021 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
2022 : azm
2023 :
2024 : ! Return Variable
2025 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
2026 : cubic_interpolated_azt_2D
2027 :
2028 : ! Local Variable(s)
2029 : integer :: &
2030 : k, i
2031 :
2032 : integer :: km1, k00, kp1, kp2
2033 :
2034 : ! ---- Begin Code ----
2035 :
2036 : do k = 1, nz
2037 : do i = 1, ngrdcol
2038 :
2039 : ! k levels are based on Tak's find_indices subroutine -dschanen 24 Oct 2011
2040 : if ( k == nz ) then
2041 : km1 = nz-2
2042 : kp1 = nz
2043 : kp2 = nz
2044 : k00 = nz-1
2045 : else if ( k == 2 ) then
2046 : km1 = 1
2047 : kp1 = 2
2048 : kp2 = 3
2049 : k00 = 1
2050 : else if ( k == 1 ) then ! Extrapolation for the ghost point
2051 : km1 = nz
2052 : k00 = 1
2053 : kp1 = 2
2054 : kp2 = 3
2055 : else
2056 : km1 = k-2
2057 : kp1 = k
2058 : kp2 = k+1
2059 : k00 = k-1
2060 : end if
2061 :
2062 : ! Do the actual interpolation.
2063 : ! Use a cubic monotonic spline interpolation.
2064 : cubic_interpolated_azt_2D(i,k) = &
2065 : mono_cubic_interp( gr%zt(i,k), km1, k00, kp1, kp2, &
2066 : gr%zm(i,km1), gr%zm(i,k00), gr%zm(i,kp1), gr%zm(i,kp2), &
2067 : azm(i,km1), azm(i,k00), azm(i,kp1), azm(i,kp2) )
2068 : end do
2069 : end do
2070 :
2071 : return
2072 :
2073 : end function cubic_interpolated_azt_2D
2074 :
2075 : !=============================================================================
2076 176472 : subroutine calc_zm2zt_weights( nz, ngrdcol, &
2077 : gr )
2078 :
2079 : ! Description:
2080 : ! Function used to help in an interpolation of a variable (var_zm) located
2081 : ! on the momentum grid levels (azm) to the thermodynamic grid levels (azt).
2082 : ! This function computes a weighting factor for both the upper momentum
2083 : ! level (k) and the lower momentum level (k-1) applied to the central
2084 : ! thermodynamic level (k). For the lowermost thermodynamic grid
2085 : ! level (k=1), a weighting factor for both the momentum level at 1 and the
2086 : ! momentum level at 2 are calculated based on the use of a linear extension.
2087 : ! This function outputs the weighting factors at a single thermodynamic grid
2088 : ! level (k). The formulation used is compatible with a stretched
2089 : ! (unevenly-spaced) grid. The weights are defined as follows:
2090 : !
2091 : ! ===var_zm(k)============================================= m(k)
2092 : ! azm_weight(m_above) = factor
2093 : ! -----------var_zm(interp)-------------------------------- t(k)
2094 : ! azm_weight(m_below) = 1 - factor
2095 : ! ===var_zm(k-1)=========================================== m(k-1)
2096 : !
2097 : ! The vertical indices m(k), t(k), and m(k-1) correspond with altitudes
2098 : ! zm(k), zt(k), and zm(k-1), respectively. The letter "t" is used for
2099 : ! thermodynamic levels and the letter "m" is used for momentum levels.
2100 : !
2101 : ! For all levels k > 1:
2102 : !
2103 : ! The formula for a linear interpolation is given by:
2104 : !
2105 : ! var_zm( interp to zt(k) )
2106 : ! = [ ( var_zm(k) - var_zm(k-1) ) / ( zm(k) - zm(k-1) ) ]
2107 : ! * ( zt(k) - zm(k-1) ) + var_zm(k-1);
2108 : !
2109 : ! which can be rewritten as:
2110 : !
2111 : ! var_zm( interp to zt(k) )
2112 : ! = [ ( zt(k) - zm(k-1) ) / ( zm(k) - zm(k-1) ) ]
2113 : ! * ( var_zm(k) - var_zm(k-1) ) + var_zm(k-1).
2114 : !
2115 : ! Furthermore, the formula can be rewritten as:
2116 : !
2117 : ! var_zm( interp to zt(k) )
2118 : ! = factor * var_zm(k) + ( 1 - factor ) * var_zm(k-1);
2119 : !
2120 : ! where:
2121 : !
2122 : ! factor = ( zt(k) - zm(k-1) ) / ( zm(k) - zm(k-1) ).
2123 : !
2124 : ! One of the important uses of this function is in situations where the
2125 : ! variable to be interpolated is being treated IMPLICITLY in an equation.
2126 : ! Usually, the variable to be interpolated is involved in a derivative (such
2127 : ! as d(var_zm)/dz in the diagram below). For the term of the equation
2128 : ! containing the derivative, grid weights are needed for two interpolations,
2129 : ! rather than just one interpolation. Thus, four grid weights (labeled
2130 : ! A(k), B(k), C(k), and D(k) in the diagram below) are needed.
2131 : !
2132 : ! ===var_zm(k+1)=========================================== m(k+1)
2133 : ! A(k)
2134 : ! -----------var_zm(interp)-------------------------------- t(k+1)
2135 : ! B(k) = 1 - A(k)
2136 : ! ===var_zm(k)===========d(var_zm)/dz====================== m(k)
2137 : ! C(k)
2138 : ! -----------var_zm(interp)-------------------------------- t(k)
2139 : ! D(k) = 1 - C(k)
2140 : ! ===var_zm(k-1)=========================================== m(k-1)
2141 : !
2142 : ! The vertical indices m(k+1), t(k+1), m(k), t(k), and m(k-1) correspond
2143 : ! with altitudes zm(k+1), zt(k+1), zm(k), zt(k), and zm(k-1), respectively.
2144 : ! The letter "t" is used for thermodynamic levels and the letter "m" is used
2145 : ! for momentum levels.
2146 : !
2147 : ! The grid weights, indexed around the central momentum level (k), are
2148 : ! defined as follows:
2149 : !
2150 : ! A(k) = ( zt(k+1) - zm(k) ) / ( zm(k+1) - zm(k) );
2151 : !
2152 : ! which is the same as "factor" for the interpolation to thermodynamic
2153 : ! level (k+1). In the code, this interpolation is referenced as
2154 : ! gr%weights_zm2zt(m_above,tkp1), which can be read as "grid weight in a
2155 : ! zm2zt interpolation of the momentum level above thermodynamic
2156 : ! level (k+1) (applied to thermodynamic level (k+1))".
2157 : !
2158 : ! B(k) = 1 - [ ( zt(k+1) - zm(k) ) / ( zm(k+1) - zm(k) ) ]
2159 : ! = 1 - A(k);
2160 : !
2161 : ! which is the same as "1 - factor" for the interpolation to thermodynamic
2162 : ! level (k+1). In the code, this interpolation is referenced as
2163 : ! gr%weights_zm2zt(m_below,tkp1), which can be read as "grid weight in a
2164 : ! zm2zt interpolation of the momentum level below thermodynamic
2165 : ! level (k+1) (applied to thermodynamic level (k+1))".
2166 : !
2167 : ! C(k) = ( zt(k) - zm(k-1) ) / ( zm(k) - zm(k-1) );
2168 : !
2169 : ! which is the same as "factor" for the interpolation to thermodynamic
2170 : ! level (k). In the code, this interpolation is referenced as
2171 : ! gr%weights_zm2zt(m_above,tk), which can be read as "grid weight in a zm2zt
2172 : ! interpolation of the momentum level above thermodynamic level (k) (applied
2173 : ! to thermodynamic level (k))".
2174 : !
2175 : ! D(k) = 1 - [ ( zt(k) - zm(k-1) ) / ( zm(k) - zm(k-1) ) ]
2176 : ! = 1 - C(k);
2177 : !
2178 : ! which is the same as "1 - factor" for the interpolation to thermodynamic
2179 : ! level (k). In the code, this interpolation is referenced as
2180 : ! gr%weights_zm2zt(m_below,tk), which can be read as "grid weight in a zm2zt
2181 : ! interpolation of the momentum level below thermodynamic level (k) (applied
2182 : ! to thermodynamic level (k))".
2183 : !
2184 : ! Additionally, as long as the central momentum level (k) in the above
2185 : ! scenario is not the lowermost momentum level or the uppermost momentum
2186 : ! level (k /= 1 and k /= gr%nz), the four weighting factors have the
2187 : ! following relationships: A(k) = C(k+1) and B(k) = D(k+1).
2188 : !
2189 : !
2190 : ! Special condition for lowermost grid level, k = 1:
2191 : !
2192 : ! The lowermost thermodynamic grid level is below the lowermost momentum
2193 : ! grid level. Thus, a linear extension is used at this level. It should
2194 : ! be noted that the thermodynamic level k = 1 is considered to be below the
2195 : ! model lower boundary, which is defined to be at momentum level k = 1.
2196 : ! Thus, the values of most variables at thermodynamic level k = 1 are not
2197 : ! often needed or referenced.
2198 : !
2199 : ! For level k = 1:
2200 : !
2201 : ! The formula for a linear extension is given by:
2202 : !
2203 : ! var_zm( extend to zt(k) )
2204 : ! = [ ( var_zm(k+1) - var_zm(k) ) / ( zm(k+1) - zm(k) ) ]
2205 : ! * ( zt(k) - zm(k) ) + var_zm(k);
2206 : !
2207 : ! which can be rewritten as:
2208 : !
2209 : ! var_zm( extend to zt(k) )
2210 : ! = [ ( zt(k) - zm(k) ) / ( zm(k+1) - zm(k) ) ]
2211 : ! * ( var_zm(k+1) - var_zm(k) ) + var_zm(k).
2212 : !
2213 : ! Furthermore, the formula can be rewritten as:
2214 : !
2215 : ! var_zm( extend to zt(k) )
2216 : ! = factor * var_zm(k+1) + ( 1 - factor ) * var_zm(k);
2217 : !
2218 : ! where:
2219 : !
2220 : ! factor = ( zt(k) - zm(k) ) / ( zm(k+1) - zm(k) ).
2221 : !
2222 : ! Due to the fact that a linear extension is being used, the value of factor
2223 : ! will be less than 0. The weight of the upper momentum level, which is
2224 : ! momentum level k = 2, on thermodynamic level k = 1 equals the value of
2225 : ! factor. The weight of the lower momentum level, which is momentum level
2226 : ! k = 1, on thermodynamic level k = 1 equals 1 - factor, which is greater
2227 : ! than 1. However, the sum of the weights equals 1.
2228 : !
2229 : !
2230 : ! Brian Griffin; September 12, 2008.
2231 : !
2232 : !-----------------------------------------------------------------------
2233 :
2234 : use constants_clubb, only: &
2235 : one ! Constant(s)
2236 :
2237 : use clubb_precision, only: &
2238 : core_rknd ! Variable(s)
2239 :
2240 : implicit none
2241 :
2242 : ! Constant parameters
2243 : integer, parameter :: &
2244 : m_above = 1, & ! Upper momentum level.
2245 : m_below = 2 ! Lower momentum level.
2246 :
2247 : ! Input Variable
2248 : integer, intent(in) :: &
2249 : nz, &
2250 : ngrdcol
2251 :
2252 : ! Input/Output Variable
2253 : type (grid), target, intent(inout) :: gr
2254 :
2255 : integer :: i, k
2256 :
2257 :
2258 : ! At most levels, the momentum level is found in-between two
2259 : ! thermodynamic levels. Linear interpolation is used.
2260 2946672 : do i = 1, ngrdcol
2261 0 : gr%weights_zm2zt(i,1,m_above) = ( gr%zt(i,1) - gr%zm(i,1) ) &
2262 2770200 : / ( gr%zm(i,2) - gr%zm(i,1) )
2263 :
2264 2946672 : gr%weights_zm2zt(i,1,m_below) = one - gr%weights_zm2zt(i,1,m_above)
2265 : end do
2266 :
2267 15000120 : do k = 2, nz
2268 247696920 : do i = 1, ngrdcol
2269 :
2270 : ! The top model level (gr%nz) is formulated differently because the top
2271 : ! momentum level is above the top thermodynamic level. A linear
2272 : ! extension is required, rather than linear interpolation.
2273 : ! Note: Variable "factor" will be greater than 1 in this situation.
2274 0 : gr%weights_zm2zt(i,k,m_above) = ( gr%zt(i,k) - gr%zm(i,k-1) ) &
2275 232696800 : / ( gr%zm(i,k) - gr%zm(i,k-1) )
2276 :
2277 : ! Weight of lower thermodynamic level on momentum level.
2278 247520448 : gr%weights_zm2zt(i,k,m_below) = one - gr%weights_zm2zt(i,k,m_above)
2279 :
2280 : end do
2281 : end do
2282 :
2283 176472 : return
2284 :
2285 : end subroutine calc_zm2zt_weights
2286 :
2287 : !=============================================================================
2288 : ! Wrapped in interface ddzm
2289 0 : function gradzm_2D( nz, ngrdcol, gr, azm )
2290 :
2291 : ! Description:
2292 : ! 2D version of gradzm
2293 : !-----------------------------------------------------------------------
2294 :
2295 : use clubb_precision, only: &
2296 : core_rknd ! Variable(s)
2297 :
2298 : implicit none
2299 :
2300 : integer, intent(in) :: &
2301 : nz, &
2302 : ngrdcol
2303 :
2304 : type (grid), target, intent(in) :: gr
2305 :
2306 : ! Input Variable
2307 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
2308 : azm ! Variable on momentum grid levels [units vary]
2309 :
2310 : ! Return Variable
2311 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
2312 : gradzm_2D ! Vertical derivative of azm [units vary / m]
2313 :
2314 : ! Local Variable
2315 : integer :: i, k ! Grid level loop index
2316 :
2317 : !$acc data copyin( gr, gr%invrs_dzt, azm ) &
2318 : !$acc copyout( gradzm_2D )
2319 :
2320 : !$acc parallel loop gang vector default(present)
2321 0 : do i = 1, ngrdcol
2322 0 : gradzm_2D(i,1) = ( azm(i,2) - azm(i,1) ) * gr%invrs_dzt(i,2)
2323 : end do
2324 : !$acc end parallel loop
2325 :
2326 : !$acc parallel loop gang vector collapse(2) default(present)
2327 0 : do k = 2, nz
2328 0 : do i = 1, ngrdcol
2329 0 : gradzm_2D(i,k) = ( azm(i,k) - azm(i,k-1) ) * gr%invrs_dzt(i,k)
2330 : end do
2331 : end do
2332 : !$acc end parallel loop
2333 :
2334 : !$acc end data
2335 :
2336 0 : return
2337 :
2338 : end function gradzm_2D
2339 :
2340 : !=============================================================================
2341 : ! Wrapped in interface ddzm
2342 0 : function gradzm_1D( gr, azm )
2343 :
2344 : ! Description:
2345 : ! 2D version of gradzm
2346 : !-----------------------------------------------------------------------
2347 :
2348 : use clubb_precision, only: &
2349 : core_rknd ! Variable(s)
2350 :
2351 : implicit none
2352 :
2353 : type (grid), target, intent(in) :: gr
2354 :
2355 : ! Input Variable
2356 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
2357 : azm ! Variable on momentum grid levels [units vary]
2358 :
2359 : ! Return Variable
2360 : real( kind = core_rknd ), dimension(gr%nz) :: &
2361 : gradzm_1D ! Vertical derivative of azm [units vary / m]
2362 :
2363 : ! Local Variables
2364 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
2365 0 : azm_col ! Variable on momentum grid levels [units vary]
2366 :
2367 : ! Return Variable
2368 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
2369 0 : gradzm_1D_col ! Vertical derivative of azm [units vary / m]
2370 :
2371 0 : azm_col(1,:) = azm
2372 :
2373 0 : gradzm_1D_col = gradzm_2D(gr%nz, 1, gr, azm_col)
2374 :
2375 0 : gradzm_1D = gradzm_1D_col(1,:)
2376 :
2377 0 : return
2378 :
2379 : end function gradzm_1D
2380 :
2381 : !=============================================================================
2382 : ! Wrapped in interface ddzt
2383 9882432 : function gradzt_2D( nz, ngrdcol, gr, azt )
2384 :
2385 : ! Description:
2386 : ! 2D version of gradzt
2387 : !-----------------------------------------------------------------------
2388 :
2389 : use clubb_precision, only: &
2390 : core_rknd ! Variable(s)
2391 :
2392 : implicit none
2393 :
2394 : integer, intent(in) :: &
2395 : nz, &
2396 : ngrdcol
2397 :
2398 : type (grid), target, intent(in) :: gr
2399 :
2400 : ! Input Variable
2401 : real( kind = core_rknd ), intent(in), dimension(ngrdcol,nz) :: &
2402 : azt ! Variable on thermodynamic grid levels [units vary]
2403 :
2404 : ! Output Variable
2405 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
2406 : gradzt_2D ! Vertical derivative of azt [units vary / m]
2407 :
2408 : ! Local Variable
2409 : integer :: i, k ! Grid level loop index
2410 :
2411 : !$acc data copyin( gr, gr%invrs_dzm, azt ) &
2412 : !$acc copyout( gradzt_2D )
2413 :
2414 : !$acc parallel loop gang vector default(present)
2415 165013632 : do i = 1, ngrdcol
2416 165013632 : gradzt_2D(i,nz) = ( azt(i,nz) - azt(i,nz-1) ) * gr%invrs_dzm(i,nz-1)
2417 : end do
2418 : !$acc end parallel loop
2419 :
2420 : !$acc parallel loop gang vector collapse(2) default(present)
2421 840006720 : do k = 1, nz-1
2422 13871027520 : do i = 1, ngrdcol
2423 13861145088 : gradzt_2D(i,k) = ( azt(i,k+1) - azt(i,k) ) * gr%invrs_dzm(i,k)
2424 : end do
2425 : end do
2426 : !$acc end parallel loop
2427 :
2428 : !$acc end data
2429 :
2430 9882432 : return
2431 :
2432 : end function gradzt_2D
2433 :
2434 : !=============================================================================
2435 : ! Wrapped in interface ddzt
2436 0 : function gradzt_1D( gr, azt )
2437 :
2438 : ! Description:
2439 : ! 2D version of gradzt
2440 : !-----------------------------------------------------------------------
2441 :
2442 : use clubb_precision, only: &
2443 : core_rknd ! Variable(s)
2444 :
2445 : implicit none
2446 :
2447 : type (grid), target, intent(in) :: gr
2448 :
2449 : ! Input Variable
2450 : real( kind = core_rknd ), intent(in), dimension(gr%nz) :: &
2451 : azt ! Variable on thermodynamic grid levels [units vary]
2452 :
2453 : ! Output Variable
2454 : real( kind = core_rknd ), dimension(gr%nz) :: &
2455 : gradzt_1D ! Vertical derivative of azt [units vary / m]
2456 :
2457 : ! Local Variables
2458 :
2459 : ! Input Variable
2460 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
2461 0 : azt_col ! Variable on thermodynamic grid levels [units vary]
2462 :
2463 : real( kind = core_rknd ), dimension(1,gr%nz) :: &
2464 0 : gradzt_1D_col ! Vertical derivative of azt [units vary / m]
2465 :
2466 0 : azt_col(1,:) = azt
2467 :
2468 0 : gradzt_1D_col = gradzt_2D(gr%nz, 1, gr, azt_col)
2469 :
2470 0 : gradzt_1D = gradzt_1D_col(1,:)
2471 :
2472 0 : return
2473 :
2474 : end function gradzt_1D
2475 :
2476 : !=============================================================================
2477 0 : function flip( x, xdim )
2478 :
2479 : ! Description:
2480 : ! Flips a single dimension array (i.e. a vector), so the first element
2481 : ! becomes the last and vice versa for the whole column. This is a
2482 : ! necessary part of the code because BUGSrad and CLUBB store altitudes in
2483 : ! reverse order.
2484 : !
2485 : ! References:
2486 : ! None
2487 : !-------------------------------------------------------------------------
2488 :
2489 : use clubb_precision, only: &
2490 : dp ! double precision
2491 :
2492 : implicit none
2493 :
2494 : ! Input Variables
2495 : integer, intent(in) :: xdim
2496 :
2497 : real(kind = dp), dimension(xdim), intent(in) :: x
2498 :
2499 : ! Output Variables
2500 : real(kind = dp), dimension(xdim) :: flip
2501 :
2502 : ! Local Variables
2503 0 : real(kind = dp), dimension(xdim) :: tmp
2504 :
2505 : integer :: indx
2506 :
2507 :
2508 : ! Get rid of an annoying compiler warning.
2509 : indx = 1
2510 : indx = indx
2511 :
2512 0 : forall ( indx = 1 : xdim )
2513 0 : tmp(indx) = x((xdim+1) - (indx))
2514 : end forall
2515 :
2516 0 : flip = tmp
2517 :
2518 :
2519 0 : return
2520 :
2521 : end function flip
2522 :
2523 : !===============================================================================
2524 :
2525 0 : end module grid_class
|