FV3DYCORE  Version 1.1.0
fv_arrays.F90
Go to the documentation of this file.
1 !***********************************************************************
2 !* GNU Lesser General Public License
3 !*
4 !* This file is part of the FV3 dynamical core.
5 !*
6 !* The FV3 dynamical core is free software: you can redistribute it
7 !* and/or modify it under the terms of the
8 !* GNU Lesser General Public License as published by the
9 !* Free Software Foundation, either version 3 of the License, or
10 !* (at your option) any later version.
11 !*
12 !* The FV3 dynamical core is distributed in the hope that it will be
13 !* useful, but WITHOUT ANYWARRANTY; without even the implied warranty
14 !* of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
15 !* See the GNU General Public License for more details.
16 !*
17 !* You should have received a copy of the GNU Lesser General Public
18 !* License along with the FV3 dynamical core.
19 !* If not, see <http://www.gnu.org/licenses/>.
20 !***********************************************************************
23 
25 #include <fms_platform.h>
26  use mpp_domains_mod, only: domain2d
27  use fms_io_mod, only: restart_file_type
28  use time_manager_mod, only: time_type
29  use horiz_interp_type_mod, only: horiz_interp_type
30  use mpp_domains_mod, only: nest_domain_type
31  use mpp_mod, only: mpp_broadcast
32  use platform_mod, only: r8_kind
33  public
34 
35  integer, public, parameter :: r_grid = r8_kind
36 
37  !Several 'auxiliary' structures are introduced here. These are for
38  ! the internal use by certain modules, and although fv_atmos_type
39  ! contains one of each of these structures all memory management
40  ! is performed by the module in question.
41 
42 
43  integer, parameter:: max_step = 1000
44 !--- MAY NEED TO TEST THIS
45 #ifdef OVERLOAD_R4
46  real, parameter:: real_big = 1.e8 ! big enough to cause blowup if used
47  real, parameter:: real_snan=x'FFBFFFFF'
48 #else
49  real, parameter:: real_big = 1.e30 ! big enough to cause blowup if used
50  real, parameter:: real_snan=x'FFF7FFFFFFFFFFFF'
51 #endif
52  real, parameter:: i4_in=-huge(1)
54 
55 
56  integer ::id_ps, id_slp, id_ua, id_va, id_pt, id_omga, id_vort, &
57  id_tm, id_pv, id_zsurf, id_oro, id_sgh, id_divg, id_w, &
58  id_ke, id_te, id_zs, id_ze, id_mq, id_vorts, id_us, id_vs, &
59  id_tq, id_rh, id_c15, id_c25, id_c35, id_c45, &
60  id_f15, id_f25, id_f35, id_f45, id_ctp, &
61  id_ppt, id_ts, id_tb, id_ctt, id_pmask, id_pmaskv2, &
62  id_delp, id_delz, id_zratio, id_ws, id_iw, id_lw, &
63  id_pfhy, id_pfnh, &
64  id_qn, id_qn200, id_qn500, id_qn850, id_qp, id_mdt, &
65  id_qdt, id_aam, id_amdt, &
66  id_acly, id_acl, id_acl2, &
67  id_dbz, id_maxdbz, id_basedbz, id_dbz4km, id_dbztop, id_dbz_m10c, &
68  id_ctz, id_w1km, id_wmaxup, id_wmaxdn, id_cape, id_cin, id_diss
69 
70 ! Selected p-level fields from 3D variables:
71  integer :: id_vort200, id_vort500, id_w500, id_w700
72  integer :: id_vort850, id_w850, id_x850, id_srh25, &
73  id_uh03, id_uh25, id_theta_e, &
74  id_w200, id_s200, id_sl12, id_sl13, id_w5km, id_rain5km, id_w2500m
75  integer :: id_srh1, id_srh3, id_ustm, id_vstm
76 ! NGGPS 31-level diag
77  integer, allocatable :: id_u(:), id_v(:), id_t(:), id_h(:), id_q(:), id_omg(:)
78 
79  integer:: id_u_plev, id_v_plev, id_t_plev, id_h_plev, id_q_plev, id_omg_plev
80 ! IPCC diag
81  integer :: id_rh10, id_rh50, id_rh100, id_rh200, id_rh250, id_rh300, &
82  id_rh500, id_rh700, id_rh850, id_rh925, id_rh1000
83  integer :: id_dp10, id_dp50, id_dp100, id_dp200, id_dp250, id_dp300, &
84  id_dp500, id_dp700, id_dp850, id_dp925, id_dp1000
85 
86  integer :: id_rh1000_cmip, id_rh925_cmip, id_rh850_cmip, id_rh700_cmip, id_rh500_cmip, &
87  id_rh300_cmip, id_rh250_cmip, id_rh100_cmip, id_rh50_cmip, id_rh10_cmip
88 
89  integer :: id_hght
90  integer :: id_u100m, id_v100m, id_w100m
91 
92  ! For initial conditions:
93  integer ic_ps, ic_ua, ic_va, ic_ppt
94  integer ic_sphum
95  integer, allocatable :: id_tracer(:)
96 ! ESM requested diagnostics - dry mass/volume mixing ratios
97  integer, allocatable :: id_tracer_dmmr(:)
98  integer, allocatable :: id_tracer_dvmr(:)
99  real, allocatable :: w_mr(:)
100 
101  real, allocatable :: phalf(:)
102  real, allocatable :: zsurf(:,:)
103  real, allocatable :: zxg(:,:)
104  real, allocatable :: pt1(:)
105 
106 
107  logical :: initialized = .false.
108  real sphum, liq_wat, ice_wat ! GFDL physics
109  real rainwat, snowwat, graupel
110 
111  real :: efx(max_step), efx_sum, efx_nest(max_step), efx_sum_nest, mtq(max_step), mtq_sum
112  integer :: steps
113 
114  end type fv_diag_type
115 
116 
121  real(kind=R_GRID), allocatable, dimension(:,:,:) :: grid_64, agrid_64
122  real(kind=R_GRID), allocatable, dimension(:,:) :: area_64, area_c_64
123  real(kind=R_GRID), allocatable, dimension(:,:) :: sina_64, cosa_64
124  real(kind=R_GRID), allocatable, dimension(:,:) :: dx_64, dy_64
125  real(kind=R_GRID), allocatable, dimension(:,:) :: dxc_64, dyc_64
126  real(kind=R_GRID), allocatable, dimension(:,:) :: dxa_64, dya_64
127 
128  real, allocatable, dimension(:,:,:) :: grid, agrid
129  real, allocatable, dimension(:,:) :: area, area_c
130  real, allocatable, dimension(:,:) :: rarea, rarea_c
131 
132  real, allocatable, dimension(:,:) :: sina, cosa
133  real, allocatable, dimension(:,:,:) :: e1,e2
134  real, allocatable, dimension(:,:) :: dx, dy
135  real, allocatable, dimension(:,:) :: dxc, dyc
136  real, allocatable, dimension(:,:) :: dxa, dya
137  real, allocatable, dimension(:,:) :: rdx, rdy
138  real, allocatable, dimension(:,:) :: rdxc, rdyc
139  real, allocatable, dimension(:,:) :: rdxa, rdya
140 
141  ! Scalars:
142  real(kind=R_GRID), allocatable :: edge_s(:)
143  real(kind=R_GRID), allocatable :: edge_n(:)
144  real(kind=R_GRID), allocatable :: edge_w(:)
145  real(kind=R_GRID), allocatable :: edge_e(:)
146  ! Vector:
147  real(kind=R_GRID), allocatable :: edge_vect_s(:)
148  real(kind=R_GRID), allocatable :: edge_vect_n(:)
149  real(kind=R_GRID), allocatable :: edge_vect_w(:)
150  real(kind=R_GRID), allocatable :: edge_vect_e(:)
151  ! scalar:
152  real(kind=R_GRID), allocatable :: ex_s(:)
153  real(kind=R_GRID), allocatable :: ex_n(:)
154  real(kind=R_GRID), allocatable :: ex_w(:)
155  real(kind=R_GRID), allocatable :: ex_e(:)
156 
157  real, allocatable :: l2c_u(:,:), l2c_v(:,:)
158  ! divergence Damping:
159  real, allocatable :: divg_u(:,:), divg_v(:,:) !
160  ! del6 diffusion:
161  real, allocatable :: del6_u(:,:), del6_v(:,:) !
162  ! Cubed_2_latlon:
163  real, allocatable :: a11(:,:)
164  real, allocatable :: a12(:,:)
165  real, allocatable :: a21(:,:)
166  real, allocatable :: a22(:,:)
167  ! latlon_2_cubed:
168  real, allocatable :: z11(:,:)
169  real, allocatable :: z12(:,:)
170  real, allocatable :: z21(:,:)
171  real, allocatable :: z22(:,:)
172 
173 ! real, allocatable :: w00(:,:)
174 
175  real, allocatable :: cosa_u(:,:)
176  real, allocatable :: cosa_v(:,:)
177  real, allocatable :: cosa_s(:,:)
178  real, allocatable :: sina_u(:,:)
179  real, allocatable :: sina_v(:,:)
180  real, allocatable :: rsin_u(:,:)
181  real, allocatable :: rsin_v(:,:)
182  real, allocatable :: rsina(:,:)
183  real, allocatable :: rsin2(:,:)
184  real(kind=R_GRID), allocatable :: ee1(:,:,:)
185  real(kind=R_GRID), allocatable :: ee2(:,:,:)
186  real(kind=R_GRID), allocatable :: ec1(:,:,:)
187  real(kind=R_GRID), allocatable :: ec2(:,:,:)
188  real(kind=R_GRID), allocatable :: ew(:,:,:,:)
189  real(kind=R_GRID), allocatable :: es(:,:,:,:)
190 
191 
192  !- 3D Super grid to contain all geometrical factors --
193  ! the 3rd dimension is 9
194  real, allocatable :: sin_sg(:,:,:)
195  real, allocatable :: cos_sg(:,:,:)
196  !--------------------------------------------------
197 
198  ! Unit Normal vectors at cell edges:
199  real(kind=R_GRID), allocatable :: en1(:,:,:)
200  real(kind=R_GRID), allocatable :: en2(:,:,:)
201 
202  ! Extended Cubed cross-edge winds
203  real, allocatable :: eww(:,:)
204  real, allocatable :: ess(:,:)
205 
206  ! Unit vectors for lat-lon grid
207  real(kind=R_GRID), allocatable :: vlon(:,:,:), vlat(:,:,:)
208  real, allocatable :: fc(:,:), f0(:,:)
209 
210  integer, dimension(:,:,:), allocatable :: iinta, jinta, iintb, jintb
211 
212  !Scalar data
213 
214  integer :: npx_g, npy_g, ntiles_g ! global domain
215 
216  real(kind=R_GRID) :: global_area
217  logical :: g_sum_initialized = .false.
218  logical:: sw_corner, se_corner, ne_corner, nw_corner
219 
220  real(kind=R_GRID) :: da_min, da_max, da_min_c, da_max_c
221 
222  real :: acapn, acaps
223  real :: globalarea
224 
225  logical :: latlon = .false.
226  logical :: cubed_sphere = .false.
227  logical :: have_south_pole = .false.
228  logical :: have_north_pole = .false.
229  logical :: stretched_grid = .false.
230 
231  logical :: square_domain = .false.
232 
233  integer, pointer :: grid_type
238 
239  logical, pointer :: nested
240 
241  logical, pointer :: regional
242 
243  end type fv_grid_type
244 
246 
247  !! FOR EACH VARIABLE IN FV_FLAGS:
248  !! 1. Must be defined here:
249  !! 2. Must be broadcast in fv_atmos_data
250  !! 3. If a namelist entry, a pointer must
251  !! be defined and associated in fv_control
252  !! 4. Must NOT appear in fv_current_grid_mod.
253  !! (this module will soon be removed)
254  !! 5. Must be referenced through Atm%flagstruct,
255  !! not Atm%, unless a convenience
256  !! pointer is defined
257 
258 !-----------------------------------------------------------------------
259 ! Grid descriptor file setup
260 !-----------------------------------------------------------------------
261  character(len=80) :: grid_name = 'Gnomonic'
262  character(len=120):: grid_file = 'Inline'
263  integer :: grid_type = 0
270 ! -> moved to grid_tools
271 
272 ! Momentum (or KE) options:
273  integer :: hord_mt = 9
289 
290  integer :: kord_mt = 8
294 
295  integer :: kord_wz = 8
298 
300  integer :: hord_vt = 9
302 
303 
304 ! Heat & air mass (delp) transport options:
305  integer :: hord_tm = 9
307 
308  integer :: hord_dp = 9
311 
312  integer :: kord_tm = -8
317  integer :: hord_tr = 12
323 
324  integer :: kord_tr = 8
327 
328  real :: scale_z = 0.
329  real :: w_max = 75.
330  real :: z_min = 0.05
331  real :: lim_fac = 1.0
332 
333  integer :: nord = 1
338 
339  integer :: nord_tr = 0
341 
342  real :: dddmp = 0.0
345 
346  real :: d2_bg = 0.0
349 
350  real :: d4_bg = 0.16
355 
356  real :: vtdm4 = 0.0
362 
363  real :: trdm2 = 0.0
364  real :: d2_bg_k1 = 4.
372 
373  real :: d2_bg_k2 = 2.
376 
377  real :: d2_divg_max_k1 = 0.15
378  real :: d2_divg_max_k2 = 0.08
379  real :: damp_k_k1 = 0.2
380  real :: damp_k_k2 = 0.12
381 
383  integer :: n_zs_filter = 0
390 
391  integer :: nord_zs_filter = 4
398 
399  logical :: full_zs_filter = .false.
408 
409  logical :: rf_fast = .false.
415 
416  logical :: consv_am = .false.
417 
418 
419  logical :: do_sat_adj= .false. !
420  logical :: do_f3d = .false. !
421  logical :: no_dycore = .false.
425 
426  logical :: convert_ke = .false.
431 
432  logical :: do_vort_damp = .false.
443 
444  logical :: use_old_omega = .true.
445 ! PG off centering:
446  real :: beta = 0.0
458 
459 #ifdef SW_DYNAMICS
460  integer :: n_sponge = 0
463 
464  real :: d_ext = 0.
469 
470  integer :: nwat = 0
485 
486  logical :: warm_start = .false.
489 
490  logical :: inline_q = .true.
494 
495  logical :: adiabatic = .true.
498 #else
499  integer :: n_sponge = 1
501 
502  real :: d_ext = 0.02
506 
507  integer :: nwat = 3
522 
523  logical :: warm_start = .true.
526 
527  logical :: inline_q = .false.
531  logical :: adiabatic = .false.
532 #endif
533 !-----------------------------------------------------------
534 ! Grid shifting, rotation, and the Schmidt transformation:
535 !-----------------------------------------------------------
536  real :: shift_fac = 18.
544 
545 ! Defaults for Schmidt transformation:
546  logical :: do_schmidt = .false.
549  real(kind=R_GRID) :: stretch_fac = 1.
557 
558  real(kind=R_GRID) :: target_lat = -90.
564 
565  real(kind=R_GRID) :: target_lon = 0.
567 
568  !-----------------------------------------------------------------------------------------------
569  ! Example #1a: US regional climate simulation, center located over Oklahoma city: (262.4, 35.4)
570  ! stretching factor: 2.5
571  ! Example #1b: US Hurricane model, center over Miami: (279.7, 25.8)
572  ! stretching factor: 3-5
573  ! Example #2a: South-East Asia Regional Climate H*** (SERACH), Central Taiwan: (121.0, 23.5)
574  ! Example #2b: Typhoon Model: (124.0, 22.0)
575  ! stretching factor: 5-10
576  !-----------------------------------------------------------------------------------------------
577 
578  logical :: reset_eta = .false.
579  real :: p_fac = 0.05
588 
589  real :: a_imp = 0.75
599 
600  integer :: n_split = 0
604 
605  real :: fac_n_spl = 1.0
606  real :: fhouri = 0.0
608 
609  integer :: m_split = 0
610  integer :: k_split = 1
612 
613  logical :: use_logp = .false.
618 
619 ! For doubly periodic domain with sim_phys
620 ! 5km 150 20 (7.5 s) 2
621 !
622 ! Estimates for Gnomonic grids:
623  !===================================================
624  ! dx (km) dt (sc) n_split m_split
625  !===================================================
626  ! C1000: ~10 150 16 3
627  ! C2000: ~5 90 18 (5 s) 2
628  !===================================================
629 ! The nonhydrostatic algorithm is described in Lin 2006, QJ, (submitted)
630 ! C2000 should easily scale to at least 6 * 100 * 100 = 60,000 CPUs
631 ! For a 1024 system: try 6 x 13 * 13 = 1014 CPUs
632 
633  integer :: q_split = 0
637 
638  integer :: print_freq = 0
642 
643  logical :: write_3d_diags = .true.
648 !------------------------------------------
649 ! Model Domain parameters
650 !------------------------------------------
651  integer :: npx
654 
655  integer :: npy
658 
659  integer :: npz
662 
663  integer :: npz_rst = 0
667 
668  integer :: ncnst = 0
673 
674  integer :: pnats = 0
678 
679  integer :: dnats = 0
682 
683  integer :: ntiles = 1
687 
688  integer :: ndims = 2 ! Lat-Lon Dims for Grid in Radians
689  integer :: nf_omega = 1
691 
692  integer :: fv_sg_adj = -1
699 
700  integer :: na_init = 0
705 
706 
707  logical :: nudge_dz = .false.
712 
713  real :: p_ref = 1.e5
719 
720  real :: dry_mass = 98290.
724 
725  integer :: nt_prog = 0
726  integer :: nt_phys = 0
727  real :: tau_h2o = 0.
734 
735 
736  real :: delt_max = 1.
742 
743  real :: d_con = 0.
747 
748  real :: ke_bg = 0.
750 
751  real :: consv_te = 0.
761 
762  real :: tau = 0.
770 
771  real :: rf_cutoff = 30.e2
772 
773  logical :: filter_phys = .false.
774  logical :: dwind_2d = .false.
781 
782  logical :: breed_vortex_inline = .false.
785 
786  logical :: range_warn = .false.
791 
792  logical :: fill = .false.
795  logical :: fill_dp = .false.
801 
802 
803  logical :: fill_wz = .false.
804  logical :: check_negative = .false.
805  logical :: non_ortho = .true.
806  logical :: moist_phys = .true.
807  logical :: do_held_suarez = .false.
810  logical :: do_reed_physics = .false.
811  logical :: reed_cond_only = .false.
812  logical :: reproduce_sum = .true.
818 
819  logical :: adjust_dry_mass = .false.
826 
827  logical :: fv_debug = .false.
829  logical :: srf_init = .false.
830  logical :: mountain = .true.
837  logical :: old_divg_damp = .false.
847 
848  logical :: remap_t = .true.
855 
856  logical :: z_tracer = .false.
860 
861  logical :: fv_land = .false.
869 !--------------------------------------------------------------------------------------
870 ! The following options are useful for NWP experiments using datasets on the lat-lon grid
871 !--------------------------------------------------------------------------------------
872  logical :: nudge = .false.
875 
876  logical :: nudge_ic = .false.
879 
880  logical :: ncep_ic = .false.
884 
885  logical :: nggps_ic = .false.
888 
889  logical :: ecmwf_ic = .false.
891 
892  logical :: gfs_phil = .false.
893 
894  logical :: agrid_vel_rst = .false.
898 
899  logical :: use_new_ncep = .false. ! use the NCEP ICs created after 2014/10/22, if want to read CWAT
900  logical :: use_ncep_phy = .false. ! if .T., separate CWAT by weights of liq_wat and liq_ice in FV_IC
901  logical :: fv_diag_ic = .false. ! reconstruct IC from fv_diagnostics on lat-lon grid
902 
903  logical :: external_ic = .false.
910 
911  logical :: external_eta = .false.
915 
916  logical :: read_increment = .false.
917 ! following are namelist parameters for Stochastic Energy Baskscatter dissipation estimate
918  logical :: do_skeb = .false.
919  integer :: skeb_npass = 11
920 ! Default restart files from the "Memphis" latlon FV core:
921  character(len=128) :: res_latlon_dynamics = 'INPUT/fv_rst.res.nc'
923  character(len=128) :: res_latlon_tracers = 'INPUT/atmos_tracers.res.nc'
928 ! The user also needs to copy the "cold start" cubed sphere restart files (fv_core.res.tile1-6)
929 ! to the INPUT dir during runtime
930 !------------------------------------------------
931 ! Parameters related to non-hydrostatic dynamics:
932 !------------------------------------------------
933  logical :: hydrostatic = .true.
935 
936  logical :: phys_hydrostatic = .true.
941 
942  logical :: use_hydro_pressure = .false.
945 
946  logical :: do_uni_zfull = .false.
952 
953  logical :: hybrid_z = .false.
956 
957  logical :: make_nh = .false.
960 
961  logical :: make_hybrid_z = .false.
964 
965  logical :: nudge_qv = .false.
972 
973  real :: add_noise = -1.
976 
977  integer :: a2b_ord = 4
981 
982  integer :: c2l_ord = 4
987 
988 
989  real(kind=R_GRID) :: dx_const = 1000.
992 
993  real(kind=R_GRID) :: dy_const = 1000.
996 
997  real(kind=R_GRID) :: deglat = 15.
1000 
1001  !The following deglat_*, deglon_* options are not used.
1002  real(kind=R_GRID) :: deglon_start = -30., deglon_stop = 30., &
1003  deglat_start = -30., deglat_stop = 30.
1004 
1005  logical :: regional = .false.
1006 
1007  integer :: bc_update_interval = 3
1008 
1010  integer, pointer :: grid_number
1011 
1012  !f1p
1013  logical :: adj_mass_vmr = .false. !TER: This is to reproduce answers for verona patch. This default can be changed
1014  ! to .true. in the next city release if desired
1015 
1016  !integer, pointer :: test_case
1017  !real, pointer :: alpha
1018 
1019  end type fv_flags_type
1020 
1022 
1023  !!! CLEANUP: could we have pointers to np[xyz], nest_domain, and the index/weight arrays?
1024 
1025  real, allocatable, dimension(:,:,:) :: west_t1, east_t1, south_t1, north_t1
1026  real, allocatable, dimension(:,:,:) :: west_t0, east_t0, south_t0, north_t0
1027 
1028  integer :: istag, jstag
1029 
1030  logical :: allocated = .false.
1031  logical :: initialized = .false.
1032 
1033  end type fv_nest_bc_type_3d
1034 
1036 
1037  real, allocatable, dimension(:,:,:,:) :: west_t1, east_t1, south_t1, north_t1
1038  real, allocatable, dimension(:,:,:,:) :: west_t0, east_t0, south_t0, north_t0
1039 
1040  integer :: istag, jstag
1041 
1042  logical :: allocated = .false.
1043  logical :: initialized = .false.
1044 
1045  end type fv_nest_bc_type_4d
1046 
1048 
1049 !nested grid flags:
1050 
1051  integer :: refinement = 3
1057 
1058 
1059  integer :: parent_tile = 1
1065 
1066  logical :: nested = .false.
1067 
1068  integer :: nestbctype = 1
1069  integer :: nsponge = 0
1070  integer :: nestupdate = 0
1072 
1073  logical :: twowaynest = .false.
1076  integer :: ioffset, joffset
1077 
1078  integer :: nest_timestep = 0
1079  integer :: tracer_nest_timestep = 0
1080  real :: s_weight = 1.e-6
1081  logical :: first_step = .true.
1082  integer :: refinement_of_global = 1
1083  integer :: npx_global
1084  integer :: upoff = 1
1085  integer :: isu = -999, ieu = -1000, jsu = -999, jeu = -1000
1086 
1087  type(nest_domain_type) :: nest_domain
1088  type(nest_domain_type), allocatable :: nest_domain_all(:)
1089 
1090  !Interpolation arrays for grid nesting
1091  integer, allocatable, dimension(:,:,:) :: ind_h, ind_u, ind_v, ind_b
1092  real, allocatable, dimension(:,:,:) :: wt_h, wt_u, wt_v, wt_b
1093  integer, allocatable, dimension(:,:,:) :: ind_update_h
1094 
1095  !These arrays are not allocated by allocate_fv_atmos_type; but instead
1096  !allocated for all grids, regardless of whether the grid is
1097  !on a PE of a concurrent run.
1098  logical, allocatable, dimension(:) :: child_grids
1099 
1100  logical :: parent_proc, child_proc
1101  logical :: parent_of_twoway = .false.
1102 
1103  !These are for time-extrapolated BCs
1104  type(fv_nest_bc_type_3d) :: delp_bc, u_bc, v_bc, uc_bc, vc_bc, divg_bc
1105  type(fv_nest_bc_type_3d), allocatable, dimension(:) :: q_bc
1106 #ifndef SW_DYNAMICS
1107  type(fv_nest_bc_type_3d) :: pt_bc, w_bc, delz_bc
1108 #ifdef USE_COND
1109  type(fv_nest_bc_type_3d) :: q_con_bc
1110 #ifdef MOIST_CAPPA
1111  type(fv_nest_bc_type_3d) :: cappa_bc
1112 #endif
1113 #endif
1114 #endif
1115 
1116  !These are for tracer flux BCs
1117  logical :: do_flux_bcs, do_2way_flux_bcs
1118  type(restart_file_type) :: bcfile_ne, bcfile_sw
1119 
1120  end type fv_nest_type
1121 
1127  module procedure allocate_fv_nest_bc_type_3d
1128  module procedure allocate_fv_nest_bc_type_3d_atm
1129  end interface
1130 
1134  module procedure deallocate_fv_nest_bc_type_3d
1135  end interface
1136 
1138 
1139  integer :: is, ie, js, je
1140  integer :: isd, ied, jsd, jed
1141  integer :: isc, iec, jsc, jec
1142 
1143  integer :: ng
1144 
1145  end type fv_grid_bounds_type
1146 
1148 
1149  integer :: is_north ,ie_north ,js_north ,je_north &
1150  ,is_south ,ie_south ,js_south ,je_south &
1151  ,is_east ,ie_east ,js_east ,je_east &
1152  ,is_west ,ie_west ,js_west ,je_west
1153 
1154  integer :: is_north_uvs ,ie_north_uvs ,js_north_uvs ,je_north_uvs &
1155  ,is_south_uvs ,ie_south_uvs ,js_south_uvs ,je_south_uvs &
1156  ,is_east_uvs ,ie_east_uvs ,js_east_uvs ,je_east_uvs &
1157  ,is_west_uvs ,ie_west_uvs ,js_west_uvs ,je_west_uvs
1158 
1159  integer :: is_north_uvw ,ie_north_uvw ,js_north_uvw ,je_north_uvw &
1160  ,is_south_uvw ,ie_south_uvw ,js_south_uvw ,je_south_uvw &
1161  ,is_east_uvw ,ie_east_uvw ,js_east_uvw ,je_east_uvw &
1162  ,is_west_uvw ,ie_west_uvw ,js_west_uvw ,je_west_uvw
1163 
1165 
1167 
1168  logical :: allocated = .false.
1169  logical :: dummy = .false. ! same as grids_on_this_pe(n)
1170  integer :: grid_number = 1
1171 
1172  !Timestep-related variables.
1173 
1174  type(time_type) :: time_init, time, run_length, time_end, time_step_atmos
1175 
1176  logical :: grid_active = .true. !Always active for now
1177 
1178  !This is kept here instead of in neststruct% simply for convenience
1179  type(fv_atmos_type), pointer :: parent_grid => null()
1180 
1181 !-----------------------------------------------------------------------
1182 ! Five prognostic state variables for the f-v dynamics
1183 !-----------------------------------------------------------------------
1184 ! dyn_state:
1185 ! D-grid prognostatic variables: u, v, and delp (and other scalars)
1186 !
1187 ! o--------u(i,j+1)----------o
1188 ! | | |
1189 ! | | |
1190 ! v(i,j)------scalar(i,j)----v(i+1,j)
1191 ! | | |
1192 ! | | |
1193 ! o--------u(i,j)------------o
1194 !
1195 ! The C grid component is "diagnostic" in that it is predicted every time step
1196 ! from the D grid variables.
1197  real, _allocatable :: u(:,:,:) _null
1198  real, _allocatable :: v(:,:,:) _null
1199  real, _allocatable :: pt(:,:,:) _null
1200  real, _allocatable :: delp(:,:,:) _null
1201  real, _allocatable :: q(:,:,:,:) _null
1202  real, _allocatable :: qdiag(:,:,:,:) _null
1203 
1204 !----------------------
1205 ! non-hydrostatic state:
1206 !----------------------------------------------------------------------
1207  real, _allocatable :: w(:,:,:) _null
1208  real, _allocatable :: delz(:,:,:) _null
1209  real, _allocatable :: ze0(:,:,:) _null
1210  real, _allocatable :: q_con(:,:,:) _null
1211 
1212 !-----------------------------------------------------------------------
1213 ! Auxilliary pressure arrays:
1214 ! The 5 vars below can be re-computed from delp and ptop.
1215 !-----------------------------------------------------------------------
1216 ! dyn_aux:
1217  real, _allocatable :: ps (:,:) _null
1218  real, _allocatable :: pe (:,:,: ) _null
1219  real, _allocatable :: pk (:,:,:) _null
1220  real, _allocatable :: peln(:,:,:) _null
1221  real, _allocatable :: pkz (:,:,:) _null
1222 
1223 ! For phys coupling:
1224  real, _allocatable :: u_srf(:,:) _null
1225  real, _allocatable :: v_srf(:,:) _null
1226  real, _allocatable :: sgh(:,:) _null
1227  real, _allocatable :: oro(:,:) _null
1228  real, _allocatable :: ts(:,:) _null
1229 ! For stochastic kinetic energy backscatter (SKEB)
1230  real, _allocatable :: diss_est(:,:,:) _null
1231 
1232 !-----------------------------------------------------------------------
1233 ! Others:
1234 !-----------------------------------------------------------------------
1235  real, _allocatable :: phis(:,:) _null
1236  real, _allocatable :: omga(:,:,:) _null
1237  real, _allocatable :: ua(:,:,:) _null
1238  real, _allocatable :: va(:,:,:) _null
1239  real, _allocatable :: uc(:,:,:) _null
1240  real, _allocatable :: vc(:,:,:) _null
1241 
1242  real, _allocatable :: ak(:) _null
1243  real, _allocatable :: bk(:) _null
1244 
1245  integer :: ks
1246 
1247 ! Accumulated Mass flux arrays
1248  real, _allocatable :: mfx(:,:,:) _null
1249  real, _allocatable :: mfy(:,:,:) _null
1250 ! Accumulated Courant number arrays
1251  real, _allocatable :: cx(:,:,:) _null
1252  real, _allocatable :: cy(:,:,:) _null
1253 
1254  type(fv_flags_type) :: flagstruct
1255 
1256  !! Convenience pointers
1257  integer, pointer :: npx, npy, npz, ncnst, ng
1258 
1259  integer, allocatable, dimension(:) :: pelist
1260 
1262 
1263  type(fv_regional_bc_bounds_type) :: regional_bc_bounds
1264 
1265  type(domain2d) :: domain
1266 #if defined(SPMD)
1267 
1268  type(domain2d) :: domain_for_coupler
1269 
1270  integer :: num_contact, npes_per_tile, tile, npes_this_grid
1271  integer :: layout(2), io_layout(2) = (/ 1,1 /)
1280 
1281 #endif
1282  !These do not actually belong to the grid, but to the process
1283  !integer :: masterproc
1284  !integer :: gid
1285 
1286 !!!!!!!!!!!!!!!!
1287 ! From fv_grid_tools
1288 !!!!!!!!!!!!!!!!
1289 
1290 
1291  real :: ptop
1292 
1293  type(fv_grid_type) :: gridstruct
1294 
1295 
1296 !!!!!!!!!!!!!!!!
1297 !fv_diagnostics!
1298 !!!!!!!!!!!!!!!!
1299 
1300  type(fv_diag_type) :: idiag
1301 
1302 !!!!!!!!!!!!!!
1303 ! From fv_io !
1304 !!!!!!!!!!!!!!
1305  type(restart_file_type) :: fv_restart, sst_restart, fv_tile_restart, &
1306  rsf_restart, mg_restart, lnd_restart, tra_restart
1307 
1308  type(fv_nest_type) :: neststruct
1309 
1310  !Hold on to coarse-grid global grid, so we don't have to waste processor time getting it again when starting to do grid nesting
1311  real(kind=R_GRID), allocatable, dimension(:,:,:,:) :: grid_global
1312 
1313  integer :: atmos_axes(4)
1314 
1315  end type fv_atmos_type
1316 
1317 contains
1318 
1322  subroutine allocate_fv_atmos_type(Atm, isd_in, ied_in, jsd_in, jed_in, is_in, ie_in, js_in, je_in, &
1323  npx_in, npy_in, npz_in, ndims_in, ncnst_in, nq_in, ng_in, dummy, alloc_2d, ngrids_in)
1325  !WARNING: Before calling this routine, be sure to have set up the
1326  ! proper domain parameters from the namelists (as is done in
1327  ! fv_control.F90)
1328 
1329  implicit none
1330  type(fv_atmos_type), intent(INOUT), target :: Atm
1331  integer, intent(IN) :: isd_in, ied_in, jsd_in, jed_in, is_in, ie_in, js_in, je_in
1332  integer, intent(IN) :: npx_in, npy_in, npz_in, ndims_in, ncnst_in, nq_in, ng_in
1333  logical, intent(IN) :: dummy, alloc_2d
1334  integer, intent(IN) :: ngrids_in
1335  integer:: isd, ied, jsd, jed, is, ie, js, je
1336  integer:: npx, npy, npz, ndims, ncnst, nq, ng
1337 
1338  !For 2D utility arrays
1339  integer:: isd_2d, ied_2d, jsd_2d, jed_2d, is_2d, ie_2d, js_2d, je_2d
1340  integer:: npx_2d, npy_2d, npz_2d, ndims_2d, ncnst_2d, nq_2d, ng_2d
1341 
1342  integer :: i,j,k, ns, n
1343 
1344  if (atm%allocated) return
1345 
1346  if (dummy) then
1347  isd = 0
1348  ied= -1
1349  jsd= 0
1350  jed= -1
1351  is= 0
1352  ie= -1
1353  js= 0
1354  je= -1
1355  npx= 1
1356  npy= 1
1357  npz= 1
1358  ndims= 1
1359  ncnst= 1
1360  nq= 1
1361  ng = 1
1362  else
1363  isd = isd_in
1364  ied= ied_in
1365  jsd= jsd_in
1366  jed= jed_in
1367  is= is_in
1368  ie= ie_in
1369  js= js_in
1370  je= je_in
1371  npx= npx_in
1372  npy= npy_in
1373  npz= npz_in
1374  ndims= ndims_in
1375  ncnst= ncnst_in
1376  nq= nq_in
1377  ng = ng_in
1378  endif
1379 
1380  if ((.not. dummy) .or. alloc_2d) then
1381  isd_2d = isd_in
1382  ied_2d= ied_in
1383  jsd_2d= jsd_in
1384  jed_2d= jed_in
1385  is_2d= is_in
1386  ie_2d= ie_in
1387  js_2d= js_in
1388  je_2d= je_in
1389  npx_2d= npx_in
1390  npy_2d= npy_in
1391  npz_2d= npz_in
1392  ndims_2d= ndims_in
1393  ncnst_2d= ncnst_in
1394  nq_2d= nq_in
1395  ng_2d = ng_in
1396  else
1397  isd_2d = 0
1398  ied_2d= -1
1399  jsd_2d= 0
1400  jed_2d= -1
1401  is_2d= 0
1402  ie_2d= -1
1403  js_2d= 0
1404  je_2d= -1
1405  npx_2d= 1
1406  npy_2d= 1
1407  npz_2d= npz_in !for ak, bk, which are 1D arrays and thus OK to allocate
1408  ndims_2d= 1
1409  ncnst_2d= 1
1410  nq_2d= 1
1411  ng_2d = 1
1412  endif
1413 
1414 !This should be set up in fv_mp_mod
1415 !!$ Atm%bd%isd = isd_in
1416 !!$ Atm%bd%ied = ied_in
1417 !!$ Atm%bd%jsd = jsd_in
1418 !!$ Atm%bd%jed = jed_in
1419 !!$
1420 !!$ Atm%bd%is = is_in
1421 !!$ Atm%bd%ie = ie_in
1422 !!$ Atm%bd%js = js_in
1423 !!$ Atm%bd%je = je_in
1424 !!$
1425 !!$ Atm%bd%isc = Atm%bd%is
1426 !!$ Atm%bd%iec = Atm%bd%ie
1427 !!$ Atm%bd%jsc = Atm%bd%js
1428 !!$ Atm%bd%jec = Atm%bd%je
1429 
1430  atm%bd%ng = ng
1431 
1432  !Convenience pointers
1433  atm%npx => atm%flagstruct%npx
1434  atm%npy => atm%flagstruct%npy
1435  atm%npz => atm%flagstruct%npz
1436  atm%ncnst => atm%flagstruct%ncnst
1437 
1438  atm%ng => atm%bd%ng
1439 
1440 !!$ Atm%npx = npx_in
1441 !!$ Atm%npy = npy_in
1442 !!$ Atm%npz = npz_in
1443  atm%flagstruct%ndims = ndims_in
1444 
1445  allocate ( atm%u(isd:ied ,jsd:jed+1,npz) ) ; atm%u=real_snan
1446  allocate ( atm%v(isd:ied+1,jsd:jed ,npz) ) ; atm%v=real_snan
1447 
1448  allocate ( atm%pt(isd:ied ,jsd:jed ,npz) ) ; atm%pt=real_snan
1449  allocate ( atm%delp(isd:ied ,jsd:jed ,npz) ) ; atm%delp=real_snan
1450  allocate ( atm%q(isd:ied ,jsd:jed ,npz, nq) ) ; atm%q=real_snan
1451  allocate (atm%qdiag(isd:ied ,jsd:jed ,npz, nq+1:ncnst) ) ; atm%qdiag=real_snan
1452 
1453  ! Allocate Auxilliary pressure arrays
1454  allocate ( atm%ps(isd:ied ,jsd:jed) ) ; atm%ps=real_snan
1455  allocate ( atm%pe(is-1:ie+1, npz+1,js-1:je+1) ) ; atm%pe=real_snan
1456  allocate ( atm%pk(is:ie ,js:je , npz+1) ) ; atm%pk=real_snan
1457 !!! allocate ( Atm%peln(isd:ied,npz+1,jsd:jed) ) ! Does regional need this or is the following lines okay?
1458  allocate ( atm%peln(is:ie,npz+1,js:je) ) ; atm%peln=real_snan
1459  allocate ( atm%pkz(is:ie,js:je,npz) ) ; atm%pkz=real_snan
1460 
1461  allocate ( atm%u_srf(is:ie,js:je) ) ; atm%u_srf=real_snan
1462  allocate ( atm%v_srf(is:ie,js:je) ) ; atm%v_srf=real_snan
1463 
1464  if ( atm%flagstruct%fv_land ) then
1465  allocate ( atm%sgh(is:ie,js:je) ) ; atm%sgh=real_snan
1466  allocate ( atm%oro(is:ie,js:je) ) ; atm%oro=real_snan
1467  else
1468  allocate ( atm%oro(1,1) ) ; atm%oro=real_snan
1469  endif
1470 
1471  ! Allocate others
1472  allocate ( atm%diss_est(isd:ied ,jsd:jed ,npz) ) ; atm%diss_est=real_snan
1473  allocate ( atm%ts(is:ie,js:je) ) ; atm%ts=real_snan
1474  allocate ( atm%phis(isd:ied ,jsd:jed ) ) ; atm%phis=real_snan
1475  allocate ( atm%omga(isd:ied ,jsd:jed ,npz) ); atm%omga=0.
1476  allocate ( atm%ua(isd:ied ,jsd:jed ,npz) ) ; atm%ua=real_snan
1477  allocate ( atm%va(isd:ied ,jsd:jed ,npz) ) ; atm%va=real_snan
1478  allocate ( atm%uc(isd:ied+1,jsd:jed ,npz) ) ; atm%uc=real_snan
1479  allocate ( atm%vc(isd:ied ,jsd:jed+1,npz) ) ; atm%vc=real_snan
1480  ! For tracer transport:
1481  allocate ( atm%mfx(is:ie+1, js:je, npz) ) ; atm%mfx=real_snan
1482  allocate ( atm%mfy(is:ie , js:je+1,npz) ) ; atm%mfy=real_snan
1483  allocate ( atm%cx(is:ie+1, jsd:jed, npz) ) ; atm%cx=real_snan
1484  allocate ( atm%cy(isd:ied ,js:je+1, npz) ) ; atm%cy=real_snan
1485 
1486  allocate ( atm%ak(npz_2d+1) ) ; atm%ak=real_snan
1487  allocate ( atm%bk(npz_2d+1) ) ; atm%bk=real_snan
1488 
1489  !--------------------------
1490  ! Non-hydrostatic dynamics:
1491  !--------------------------
1492  if ( atm%flagstruct%hydrostatic ) then
1493  !Note length-one initialization if hydrostatic = .true.
1494  allocate ( atm%w(isd:isd, jsd:jsd ,1) ) ; atm%w=real_snan
1495  allocate ( atm%delz(isd:isd, jsd:jsd ,1) ) ; atm%delz=real_snan
1496  allocate ( atm%ze0(is:is, js:js ,1) ) ; atm%ze0=real_snan
1497  else
1498  allocate ( atm%w(isd:ied, jsd:jed ,npz ) ) ; atm%w=real_snan
1499  allocate ( atm%delz(isd:ied, jsd:jed ,npz) ) ; atm%delz=real_snan
1500  if( atm%flagstruct%hybrid_z ) then
1501  allocate ( atm%ze0(is:ie, js:je ,npz+1) ) ; atm%ze0=real_snan
1502  else
1503  allocate ( atm%ze0(is:is, js:js ,1) ) ; atm%ze0=real_snan
1504  endif
1505  ! allocate ( mono(isd:ied, jsd:jed, npz))
1506  endif
1507 
1508 #ifdef USE_COND
1509  allocate ( atm%q_con(isd:ied,jsd:jed,1:npz) ) ; atm%q_con=real_snan; atm%q_con=0.0
1510 #else
1511  allocate ( atm%q_con(isd:isd,jsd:jsd,1) ) ; atm%q_con=real_snan; atm%q_con=0.0
1512 #endif
1513 
1514 #ifndef NO_TOUCH_MEM
1515 ! Notes by SJL
1516 ! Place the memory in the optimal shared mem space
1517 ! This will help the scaling with OpenMP
1518 !$OMP parallel do default(none) shared(isd,ied,jsd,jed,npz,Atm,nq,ncnst)
1519  do k=1, npz
1520  do j=jsd, jed
1521  do i=isd, ied
1522  atm%ua(i,j,k) = real_big
1523  atm%va(i,j,k) = real_big
1524  atm%pt(i,j,k) = real_big
1525  atm%delp(i,j,k) = real_big
1526  enddo
1527  enddo
1528  do j=jsd, jed+1
1529  do i=isd, ied
1530  atm%u(i,j,k) = real_big
1531  atm%vc(i,j,k) = real_big
1532  enddo
1533  enddo
1534  do j=jsd, jed
1535  do i=isd, ied+1
1536  atm%v(i,j,k) = real_big
1537  atm%uc(i,j,k) = real_big
1538  enddo
1539  enddo
1540  if ( .not. atm%flagstruct%hydrostatic ) then
1541  do j=jsd, jed
1542  do i=isd, ied
1543  atm%w(i,j,k) = real_big
1544  atm%delz(i,j,k) = real_big
1545  enddo
1546  enddo
1547  endif
1548  do n=1,nq
1549  do j=jsd, jed
1550  do i=isd, ied
1551  atm%q(i,j,k,n) = real_big
1552  enddo
1553  enddo
1554  enddo
1555  do n=nq+1,ncnst
1556  do j=jsd, jed
1557  do i=isd, ied
1558  atm%qdiag(i,j,k,n) = real_big
1559  enddo
1560  enddo
1561  enddo
1562  enddo
1563 #endif
1564 
1565  allocate ( atm%gridstruct% area(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct% area=real_snan ! Cell Centered
1566  allocate ( atm%gridstruct% area_64(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct% area_64=real_snan ! Cell Centered
1567  allocate ( atm%gridstruct%rarea(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct%rarea=real_snan ! Cell Centered
1568 
1569  allocate ( atm%gridstruct% area_c(isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct% area_c=real_snan ! Cell Corners
1570  allocate ( atm%gridstruct% area_c_64(isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct% area_c_64=real_snan ! Cell Corners
1571  allocate ( atm%gridstruct%rarea_c(isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct%rarea_c=real_snan ! Cell Corners
1572 
1573  allocate ( atm%gridstruct% dx(isd_2d:ied_2d ,jsd_2d:jed_2d+1) ) ; atm%gridstruct% dx=real_snan
1574  allocate ( atm%gridstruct% dx_64(isd_2d:ied_2d ,jsd_2d:jed_2d+1) ) ; atm%gridstruct% dx_64=real_snan
1575  allocate ( atm%gridstruct%rdx(isd_2d:ied_2d ,jsd_2d:jed_2d+1) ) ; atm%gridstruct%rdx=real_snan
1576  allocate ( atm%gridstruct% dy(isd_2d:ied_2d+1,jsd_2d:jed_2d ) ) ; atm%gridstruct% dy=real_snan
1577  allocate ( atm%gridstruct% dy_64(isd_2d:ied_2d+1,jsd_2d:jed_2d ) ) ; atm%gridstruct% dy_64=real_snan
1578  allocate ( atm%gridstruct%rdy(isd_2d:ied_2d+1,jsd_2d:jed_2d ) ) ; atm%gridstruct%rdy=real_snan
1579 
1580  allocate ( atm%gridstruct% dxc(isd_2d:ied_2d+1,jsd_2d:jed_2d ) ) ; atm%gridstruct% dxc=real_snan
1581  allocate ( atm%gridstruct% dxc_64(isd_2d:ied_2d+1,jsd_2d:jed_2d ) ) ; atm%gridstruct% dxc_64=real_snan
1582  allocate ( atm%gridstruct%rdxc(isd_2d:ied_2d+1,jsd_2d:jed_2d ) ) ; atm%gridstruct%rdxc=real_snan
1583  allocate ( atm%gridstruct% dyc(isd_2d:ied_2d ,jsd_2d:jed_2d+1) ) ; atm%gridstruct% dyc=real_snan
1584  allocate ( atm%gridstruct% dyc_64(isd_2d:ied_2d ,jsd_2d:jed_2d+1) ) ; atm%gridstruct% dyc_64=real_snan
1585  allocate ( atm%gridstruct%rdyc(isd_2d:ied_2d ,jsd_2d:jed_2d+1) ) ; atm%gridstruct%rdyc=real_snan
1586 
1587  allocate ( atm%gridstruct% dxa(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct% dxa=real_snan
1588  allocate ( atm%gridstruct% dxa_64(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct% dxa_64=real_snan
1589  allocate ( atm%gridstruct%rdxa(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct%rdxa=real_snan
1590  allocate ( atm%gridstruct% dya(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct% dya=real_snan
1591  allocate ( atm%gridstruct% dya_64(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct% dya_64=real_snan
1592  allocate ( atm%gridstruct%rdya(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct%rdya=real_snan
1593 
1594  allocate ( atm%gridstruct%grid (isd_2d:ied_2d+1,jsd_2d:jed_2d+1,1:ndims_2d) ) ; atm%gridstruct%grid=real_snan
1595  allocate ( atm%gridstruct%grid_64 (isd_2d:ied_2d+1,jsd_2d:jed_2d+1,1:ndims_2d) ) ; atm%gridstruct%grid_64=real_snan
1596  allocate ( atm%gridstruct%agrid(isd_2d:ied_2d ,jsd_2d:jed_2d ,1:ndims_2d) ) ; atm%gridstruct%agrid=real_snan
1597  allocate ( atm%gridstruct%agrid_64(isd_2d:ied_2d ,jsd_2d:jed_2d ,1:ndims_2d) ) ; atm%gridstruct%agrid_64=real_snan
1598  allocate ( atm%gridstruct% sina(isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct% sina=real_snan ! SIN(angle of intersection)
1599  allocate ( atm%gridstruct% sina_64(isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct% sina_64=real_snan ! SIN(angle of intersection)
1600  allocate ( atm%gridstruct%rsina(is_2d:ie_2d+1,js_2d:je_2d+1) ) ; atm%gridstruct%rsina=real_snan ! Why is the size different?
1601  allocate ( atm%gridstruct% cosa(isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct% cosa=real_snan ! COS(angle of intersection)
1602  allocate ( atm%gridstruct% cosa_64(isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct% cosa_64=real_snan ! COS(angle of intersection)
1603 
1604  allocate ( atm%gridstruct% e1(3,isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct%e1=real_snan
1605  allocate ( atm%gridstruct% e2(3,isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct%e2=real_snan
1606 
1607  allocate (atm%gridstruct%iinta(4, isd_2d:ied_2d ,jsd_2d:jed_2d), &
1608  atm%gridstruct%jinta(4, isd_2d:ied_2d ,jsd_2d:jed_2d), &
1609  atm%gridstruct%iintb(4, is_2d:ie_2d+1 ,js_2d:je_2d+1), &
1610  atm%gridstruct%jintb(4, is_2d:ie_2d+1 ,js_2d:je_2d+1) )
1611 
1612  allocate ( atm%gridstruct%edge_s(npx_2d) ) ; atm%gridstruct%edge_s=real_snan
1613  allocate ( atm%gridstruct%edge_n(npx_2d) ) ; atm%gridstruct%edge_n=real_snan
1614  allocate ( atm%gridstruct%edge_w(npy_2d) ) ; atm%gridstruct%edge_w=real_snan
1615  allocate ( atm%gridstruct%edge_e(npy_2d) ) ; atm%gridstruct%edge_e=real_snan
1616 
1617  allocate ( atm%gridstruct%edge_vect_s(isd_2d:ied_2d) ) ; atm%gridstruct%edge_vect_s=real_snan
1618  allocate ( atm%gridstruct%edge_vect_n(isd_2d:ied_2d) ) ; atm%gridstruct%edge_vect_n=real_snan
1619  allocate ( atm%gridstruct%edge_vect_w(jsd_2d:jed_2d) ) ; atm%gridstruct%edge_vect_w=real_snan
1620  allocate ( atm%gridstruct%edge_vect_e(jsd_2d:jed_2d) ) ; atm%gridstruct%edge_vect_e=real_snan
1621 
1622  allocate ( atm%gridstruct%ex_s(npx_2d) ) ; atm%gridstruct%ex_s=real_snan
1623  allocate ( atm%gridstruct%ex_n(npx_2d) ) ; atm%gridstruct%ex_n=real_snan
1624  allocate ( atm%gridstruct%ex_w(npy_2d) ) ; atm%gridstruct%ex_w=real_snan
1625  allocate ( atm%gridstruct%ex_e(npy_2d) ) ; atm%gridstruct%ex_e=real_snan
1626 
1627 
1628  allocate ( atm%gridstruct%l2c_u(is_2d:ie_2d, js_2d:je_2d+1) ) ; atm%gridstruct%l2c_u=real_snan
1629  allocate ( atm%gridstruct%l2c_v(is_2d:ie_2d+1,js_2d:je_2d) ) ; atm%gridstruct%l2c_v=real_snan
1630 
1631  ! For diveregnce damping:
1632  allocate ( atm%gridstruct%divg_u(isd_2d:ied_2d, jsd_2d:jed_2d+1) ) ; atm%gridstruct%divg_u=real_snan
1633  allocate ( atm%gridstruct%divg_v(isd_2d:ied_2d+1,jsd_2d:jed_2d) ) ; atm%gridstruct%divg_v=real_snan
1634  ! For del6 diffusion:
1635  allocate ( atm%gridstruct%del6_u(isd_2d:ied_2d, jsd_2d:jed_2d+1) ) ; atm%gridstruct%del6_u=real_snan
1636  allocate ( atm%gridstruct%del6_v(isd_2d:ied_2d+1,jsd_2d:jed_2d) ) ; atm%gridstruct%del6_v=real_snan
1637 
1638  allocate ( atm%gridstruct%z11(is_2d-1:ie_2d+1,js_2d-1:je_2d+1) ) ; atm%gridstruct%z11=real_snan
1639  allocate ( atm%gridstruct%z12(is_2d-1:ie_2d+1,js_2d-1:je_2d+1) ) ; atm%gridstruct%z12=real_snan
1640  allocate ( atm%gridstruct%z21(is_2d-1:ie_2d+1,js_2d-1:je_2d+1) ) ; atm%gridstruct%z21=real_snan
1641  allocate ( atm%gridstruct%z22(is_2d-1:ie_2d+1,js_2d-1:je_2d+1) ) ; atm%gridstruct%z22=real_snan
1642 
1643 ! if (.not.Atm%flagstruct%hydrostatic) &
1644 ! allocate ( Atm%gridstruct%w00(is_2d-1:ie_2d+1,js_2d-1:je_2d+1) ) ; Atm%gridstruct%w00=real_snan
1645 
1646  allocate ( atm%gridstruct%a11(is_2d-1:ie_2d+1,js_2d-1:je_2d+1) ) ; atm%gridstruct%a11=real_snan
1647  allocate ( atm%gridstruct%a12(is_2d-1:ie_2d+1,js_2d-1:je_2d+1) ) ; atm%gridstruct%a12=real_snan
1648  allocate ( atm%gridstruct%a21(is_2d-1:ie_2d+1,js_2d-1:je_2d+1) ) ; atm%gridstruct%a21=real_snan
1649  allocate ( atm%gridstruct%a22(is_2d-1:ie_2d+1,js_2d-1:je_2d+1) ) ; atm%gridstruct%a22=real_snan
1650  allocate ( atm%gridstruct%vlon(is_2d-2:ie_2d+2,js_2d-2:je_2d+2,3) ) ; atm%gridstruct%vlon=real_snan
1651  allocate ( atm%gridstruct%vlat(is_2d-2:ie_2d+2,js_2d-2:je_2d+2,3) ) ; atm%gridstruct%vlat=real_snan
1652  ! Coriolis parameters:
1653  allocate ( atm%gridstruct%f0(isd_2d:ied_2d ,jsd_2d:jed_2d ) ) ; atm%gridstruct%f0=real_snan
1654  allocate ( atm%gridstruct%fC(isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct%fc=real_snan
1655 
1656  ! Corner unit vectors:
1657  allocate( atm%gridstruct%ee1(3,isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct%ee1=real_snan
1658  allocate( atm%gridstruct%ee2(3,isd_2d:ied_2d+1,jsd_2d:jed_2d+1) ) ; atm%gridstruct%ee2=real_snan
1659 
1660  ! Center unit vectors:
1661  allocate( atm%gridstruct%ec1(3,isd_2d:ied_2d,jsd_2d:jed_2d) ) ; atm%gridstruct%ec1=real_snan
1662  allocate( atm%gridstruct%ec2(3,isd_2d:ied_2d,jsd_2d:jed_2d) ) ; atm%gridstruct%ec2=real_snan
1663 
1664  ! Edge unit vectors:
1665  allocate( atm%gridstruct%ew(3,isd_2d:ied_2d+1,jsd_2d:jed_2d, 2) ) ; atm%gridstruct%ew=real_snan
1666  allocate( atm%gridstruct%es(3,isd_2d:ied_2d ,jsd_2d:jed_2d+1,2) ) ; atm%gridstruct%es=real_snan
1667 
1668  ! Edge unit "Normal" vectors: (for omega computation)
1669  allocate( atm%gridstruct%en1(3,is_2d:ie_2d, js_2d:je_2d+1) ) ; atm%gridstruct%en1=real_snan ! E-W edges
1670  allocate( atm%gridstruct%en2(3,is_2d:ie_2d+1,js_2d:je_2d ) ) ; atm%gridstruct%en2=real_snan ! N-S egdes
1671 
1672  allocate ( atm%gridstruct%cosa_u(isd_2d:ied_2d+1,jsd_2d:jed_2d) ) ; atm%gridstruct%cosa_u=real_snan
1673  allocate ( atm%gridstruct%sina_u(isd_2d:ied_2d+1,jsd_2d:jed_2d) ) ; atm%gridstruct%sina_u=real_snan
1674  allocate ( atm%gridstruct%rsin_u(isd_2d:ied_2d+1,jsd_2d:jed_2d) ) ; atm%gridstruct%rsin_u=real_snan
1675 
1676  allocate ( atm%gridstruct%cosa_v(isd_2d:ied_2d,jsd_2d:jed_2d+1) ) ; atm%gridstruct%cosa_v=real_snan
1677  allocate ( atm%gridstruct%sina_v(isd_2d:ied_2d,jsd_2d:jed_2d+1) ) ; atm%gridstruct%sina_v=real_snan
1678  allocate ( atm%gridstruct%rsin_v(isd_2d:ied_2d,jsd_2d:jed_2d+1) ) ; atm%gridstruct%rsin_v=real_snan
1679 
1680  allocate ( atm%gridstruct%cosa_s(isd_2d:ied_2d,jsd_2d:jed_2d) ) ; atm%gridstruct%rsin_v=real_snan ! cell center
1681 
1682  allocate ( atm%gridstruct%rsin2(isd_2d:ied_2d,jsd_2d:jed_2d) ) ; atm%gridstruct%rsin_v=real_snan ! cell center
1683 
1684 
1685  ! Super (composite) grid:
1686 
1687  ! 9---4---8
1688  ! | |
1689  ! 1 5 3
1690  ! | |
1691  ! 6---2---7
1692 
1693  allocate ( atm%gridstruct%cos_sg(isd_2d:ied_2d,jsd_2d:jed_2d,9) ) ; atm%gridstruct%cos_sg=real_snan
1694  allocate ( atm%gridstruct%sin_sg(isd_2d:ied_2d,jsd_2d:jed_2d,9) ) ; atm%gridstruct%sin_sg=real_snan
1695 
1696  allocate( atm%gridstruct%eww(3,4) ) ; atm%gridstruct%eww=real_snan
1697  allocate( atm%gridstruct%ess(3,4) ) ; atm%gridstruct%ess=real_snan
1698 
1699  if (atm%neststruct%nested) then
1700 
1701  allocate(atm%neststruct%ind_h(isd:ied,jsd:jed,4)) ; atm%neststruct%ind_h=i4_in
1702  allocate(atm%neststruct%ind_u(isd:ied,jsd:jed+1,4)) ; atm%neststruct%ind_u=i4_in
1703  allocate(atm%neststruct%ind_v(isd:ied+1,jsd:jed,4)) ; atm%neststruct%ind_v=i4_in
1704 
1705  allocate(atm%neststruct%wt_h(isd:ied, jsd:jed, 4)) ; atm%neststruct%wt_h=real_snan
1706  allocate(atm%neststruct%wt_u(isd:ied, jsd:jed+1,4)) ; atm%neststruct%wt_u=real_snan
1707  allocate(atm%neststruct%wt_v(isd:ied+1, jsd:jed, 4)) ; atm%neststruct%wt_v=real_snan
1708  allocate(atm%neststruct%ind_b(isd:ied+1,jsd:jed+1,4)) ; atm%neststruct%ind_b=i4_in
1709  allocate(atm%neststruct%wt_b(isd:ied+1, jsd:jed+1,4)) ; atm%neststruct%wt_b=real_snan
1710 
1711  ns = atm%neststruct%nsponge
1712 
1713  call allocate_fv_nest_bc_type(atm%neststruct%delp_BC,atm,ns,0,0,dummy)
1714  call allocate_fv_nest_bc_type(atm%neststruct%u_BC,atm,ns,0,1,dummy)
1715  call allocate_fv_nest_bc_type(atm%neststruct%v_BC,atm,ns,1,0,dummy)
1716  call allocate_fv_nest_bc_type(atm%neststruct%uc_BC,atm,ns,1,0,dummy)
1717  call allocate_fv_nest_bc_type(atm%neststruct%vc_BC,atm,ns,0,1,dummy)
1718  call allocate_fv_nest_bc_type(atm%neststruct%divg_BC,atm,ns,1,1,dummy)
1719 
1720  if (ncnst > 0) then
1721  allocate(atm%neststruct%q_BC(ncnst))
1722  do n=1,ncnst
1723  call allocate_fv_nest_bc_type(atm%neststruct%q_BC(n),atm,ns,0,0,dummy)
1724  enddo
1725  endif
1726 
1727 #ifndef SW_DYNAMICS
1728 
1729  call allocate_fv_nest_bc_type(atm%neststruct%pt_BC,atm,ns,0,0,dummy)
1730 #ifdef USE_COND
1731  call allocate_fv_nest_bc_type(atm%neststruct%q_con_BC,atm,ns,0,0,dummy)
1732 #ifdef MOIST_CAPPA
1733  call allocate_fv_nest_bc_type(atm%neststruct%cappa_BC,atm,ns,0,0,dummy)
1734 #endif
1735 #endif
1736  if (.not.atm%flagstruct%hydrostatic) then
1737  call allocate_fv_nest_bc_type(atm%neststruct%w_BC,atm,ns,0,0,dummy)
1738  call allocate_fv_nest_bc_type(atm%neststruct%delz_BC,atm,ns,0,0,dummy)
1739  endif
1740 
1741 #endif
1742 
1743  if (atm%neststruct%twowaynest) allocate(&
1744  atm%neststruct%ind_update_h( &
1745  atm%parent_grid%bd%isd:atm%parent_grid%bd%ied+1, &
1746  atm%parent_grid%bd%jsd:atm%parent_grid%bd%jed+1,2)); atm%neststruct%ind_update_h=i4_in
1747 
1748  end if
1749 
1750  !--- Do the memory allocation only for nested model
1751  if( ngrids_in > 1 ) then
1752  if (atm%flagstruct%grid_type < 4) then
1753  if (atm%neststruct%nested) then
1754  allocate(atm%grid_global(1-ng_2d:npx_2d +ng_2d,1-ng_2d:npy_2d +ng_2d,2,1)); atm%grid_global=real_snan
1755  else
1756  allocate(atm%grid_global(1-ng_2d:npx_2d +ng_2d,1-ng_2d:npy_2d +ng_2d,2,1:6)); atm%grid_global=real_snan
1757  endif
1758  end if
1759  endif
1760 
1761  atm%ptop = real_snan
1762 
1763  atm%allocated = .true.
1764  if (dummy) atm%dummy = .true.
1765 
1766  end subroutine allocate_fv_atmos_type
1767 
1769  subroutine deallocate_fv_atmos_type(Atm)
1771  implicit none
1772  type(fv_atmos_type), intent(INOUT) :: Atm
1773 
1774  integer :: n
1775 
1776  if (.not.atm%allocated) return
1777 
1778  deallocate ( atm%u )
1779  deallocate ( atm%v )
1780  deallocate ( atm%pt )
1781  deallocate ( atm%delp )
1782  deallocate ( atm%q )
1783  deallocate ( atm%qdiag )
1784  deallocate ( atm%ps )
1785  deallocate ( atm%pe )
1786  deallocate ( atm%pk )
1787  deallocate ( atm%peln )
1788  deallocate ( atm%pkz )
1789  deallocate ( atm%phis )
1790  deallocate ( atm%omga )
1791  deallocate ( atm%ua )
1792  deallocate ( atm%va )
1793  deallocate ( atm%uc )
1794  deallocate ( atm%vc )
1795  deallocate ( atm%mfx )
1796  deallocate ( atm%mfy )
1797  deallocate ( atm%cx )
1798  deallocate ( atm%cy )
1799  deallocate ( atm%ak )
1800  deallocate ( atm%bk )
1801 
1802  deallocate ( atm%u_srf )
1803  deallocate ( atm%v_srf )
1804  if( atm%flagstruct%fv_land ) deallocate ( atm%sgh )
1805  deallocate ( atm%oro )
1806 
1807  deallocate ( atm%w )
1808  deallocate ( atm%delz )
1809  deallocate ( atm%ze0 )
1810  deallocate ( atm%q_con )
1811 
1812  deallocate ( atm%gridstruct% area ) ! Cell Centered
1813  deallocate ( atm%gridstruct%rarea ) ! Cell Centered
1814 
1815  deallocate ( atm%gridstruct% area_c ) ! Cell Corners
1816  deallocate ( atm%gridstruct%rarea_c ) ! Cell Corners
1817 
1818  deallocate ( atm%gridstruct% dx )
1819  deallocate ( atm%gridstruct%rdx )
1820  deallocate ( atm%gridstruct% dy )
1821  deallocate ( atm%gridstruct%rdy )
1822 
1823  deallocate ( atm%gridstruct% dxc )
1824  deallocate ( atm%gridstruct%rdxc )
1825  deallocate ( atm%gridstruct% dyc )
1826  deallocate ( atm%gridstruct%rdyc )
1827 
1828  deallocate ( atm%gridstruct% dxa )
1829  deallocate ( atm%gridstruct%rdxa )
1830  deallocate ( atm%gridstruct% dya )
1831  deallocate ( atm%gridstruct%rdya )
1832 
1833  deallocate ( atm%gridstruct%grid )
1834  deallocate ( atm%gridstruct%agrid )
1835  deallocate ( atm%gridstruct%sina ) ! SIN(angle of intersection)
1836  deallocate ( atm%gridstruct%cosa ) ! COS(angle of intersection)
1837 
1838  deallocate ( atm%gridstruct% e1 )
1839  deallocate ( atm%gridstruct% e2 )
1840 
1841 
1842 
1843 
1844  deallocate (atm%gridstruct%iinta, &
1845  atm%gridstruct%jinta, &
1846  atm%gridstruct%iintb, &
1847  atm%gridstruct%jintb )
1848 
1849  deallocate ( atm%gridstruct%edge_s )
1850  deallocate ( atm%gridstruct%edge_n )
1851  deallocate ( atm%gridstruct%edge_w )
1852  deallocate ( atm%gridstruct%edge_e )
1853 
1854  deallocate ( atm%gridstruct%edge_vect_s )
1855  deallocate ( atm%gridstruct%edge_vect_n )
1856  deallocate ( atm%gridstruct%edge_vect_w )
1857  deallocate ( atm%gridstruct%edge_vect_e )
1858 
1859  deallocate ( atm%gridstruct%ex_s )
1860  deallocate ( atm%gridstruct%ex_n )
1861  deallocate ( atm%gridstruct%ex_w )
1862  deallocate ( atm%gridstruct%ex_e )
1863 
1864 
1865  deallocate ( atm%gridstruct%l2c_u )
1866  deallocate ( atm%gridstruct%l2c_v )
1867  ! For diveregnce damping:
1868  deallocate ( atm%gridstruct%divg_u )
1869  deallocate ( atm%gridstruct%divg_v )
1870  ! For del6 diffusion:
1871 
1872  deallocate ( atm%gridstruct%z11 )
1873  deallocate ( atm%gridstruct%z12 )
1874  deallocate ( atm%gridstruct%z21 )
1875  deallocate ( atm%gridstruct%z22 )
1876 
1877  deallocate ( atm%gridstruct%a11 )
1878  deallocate ( atm%gridstruct%a12 )
1879  deallocate ( atm%gridstruct%a21 )
1880  deallocate ( atm%gridstruct%a22 )
1881  deallocate ( atm%gridstruct%vlon )
1882  deallocate ( atm%gridstruct%vlat )
1883  ! Coriolis parameters:
1884  deallocate ( atm%gridstruct%f0 )
1885  deallocate ( atm%gridstruct%fC )
1886 
1887  ! Corner unit vectors:
1888  deallocate( atm%gridstruct%ee1 )
1889  deallocate( atm%gridstruct%ee2 )
1890 
1891  ! Center unit vectors:
1892  deallocate( atm%gridstruct%ec1 )
1893  deallocate( atm%gridstruct%ec2 )
1894 
1895  ! Edge unit vectors:
1896  deallocate( atm%gridstruct%ew )
1897  deallocate( atm%gridstruct%es )
1898 
1899  ! Edge unit "Normal" vectors: (for omega computation)
1900  deallocate( atm%gridstruct%en1 ) ! E-W edges
1901  deallocate( atm%gridstruct%en2 ) ! N-S egdes
1902 
1903  deallocate ( atm%gridstruct%cosa_u )
1904  deallocate ( atm%gridstruct%sina_u )
1905  deallocate ( atm%gridstruct%rsin_u )
1906 
1907  deallocate ( atm%gridstruct%cosa_v )
1908  deallocate ( atm%gridstruct%sina_v )
1909  deallocate ( atm%gridstruct%rsin_v )
1910 
1911  deallocate ( atm%gridstruct%cosa_s ) ! cell center
1912 
1913  deallocate ( atm%gridstruct%rsin2 ) ! cell center
1914 
1915 
1916  ! Super (composite) grid:
1917 
1918  ! 9---4---8
1919  ! | |
1920  ! 1 5 3
1921  ! | |
1922  ! 6---2---7
1923 
1924  deallocate ( atm%gridstruct%cos_sg )
1925  deallocate ( atm%gridstruct%sin_sg )
1926 
1927  deallocate( atm%gridstruct%eww )
1928  deallocate( atm%gridstruct%ess )
1929 
1930  if (atm%neststruct%nested) then
1931  deallocate(atm%neststruct%ind_h)
1932  deallocate(atm%neststruct%ind_u)
1933  deallocate(atm%neststruct%ind_v)
1934 
1935  deallocate(atm%neststruct%wt_h)
1936  deallocate(atm%neststruct%wt_u)
1937  deallocate(atm%neststruct%wt_v)
1938 
1939  deallocate(atm%neststruct%ind_b)
1940  deallocate(atm%neststruct%wt_b)
1941 
1942  call deallocate_fv_nest_bc_type(atm%neststruct%delp_BC)
1943  call deallocate_fv_nest_bc_type(atm%neststruct%u_BC)
1944  call deallocate_fv_nest_bc_type(atm%neststruct%v_BC)
1945  call deallocate_fv_nest_bc_type(atm%neststruct%uc_BC)
1946  call deallocate_fv_nest_bc_type(atm%neststruct%vc_BC)
1947  call deallocate_fv_nest_bc_type(atm%neststruct%divg_BC)
1948 
1949  if (allocated(atm%neststruct%q_BC)) then
1950  do n=1,size(atm%neststruct%q_BC)
1951  call deallocate_fv_nest_bc_type(atm%neststruct%q_BC(n))
1952  enddo
1953  endif
1954 
1955 #ifndef SW_DYNAMICS
1956  call deallocate_fv_nest_bc_type(atm%neststruct%pt_BC)
1957 #ifdef USE_COND
1958  call deallocate_fv_nest_bc_type(atm%neststruct%q_con_BC)
1959 #ifdef MOIST_CAPPA
1960  call deallocate_fv_nest_bc_type(atm%neststruct%cappa_BC)
1961 #endif
1962 #endif
1963  if (.not.atm%flagstruct%hydrostatic) then
1964  call deallocate_fv_nest_bc_type(atm%neststruct%w_BC)
1965  call deallocate_fv_nest_bc_type(atm%neststruct%delz_BC)
1966  endif
1967 #endif
1968 
1969 
1970  if (atm%neststruct%twowaynest) deallocate(atm%neststruct%ind_update_h)
1971 
1972  end if
1973 
1974  if (atm%flagstruct%grid_type < 4) then
1975  if(allocated(atm%grid_global)) deallocate(atm%grid_global)
1976  end if
1977 
1978  atm%allocated = .false.
1979 
1980  end subroutine deallocate_fv_atmos_type
1981 
1982 
1983 subroutine allocate_fv_nest_bc_type_3d_atm(BC,Atm,ns,istag,jstag,dummy)
1985  type(fv_nest_bc_type_3d), intent(INOUT) :: BC
1986  type(fv_atmos_type), intent(IN) :: Atm
1987  integer, intent(IN) :: ns, istag, jstag
1988  logical, intent(IN) :: dummy
1989 
1990  integer :: is, ie, js, je, isd, ied, jsd, jed, npx, npy, npz, ng
1991 
1992  if (bc%allocated) return
1993 
1994  is = atm%bd%is
1995  ie = atm%bd%ie
1996  js = atm%bd%js
1997  je = atm%bd%je
1998 
1999  isd = atm%bd%isd
2000  ied = atm%bd%ied
2001  jsd = atm%bd%jsd
2002  jed = atm%bd%jed
2003 
2004  npx = atm%npx
2005  npy = atm%npy
2006  npz = atm%npz
2007 
2008  ng = atm%ng
2009 
2010  call allocate_fv_nest_bc_type_3d(bc,is,ie,js,je,isd,ied,jsd,jed,npx,npy,npz,ng,ns,istag,jstag,dummy)
2011 
2012 
2013 end subroutine allocate_fv_nest_bc_type_3d_atm
2014 
2015 subroutine allocate_fv_nest_bc_type_3d(BC,is,ie,js,je,isd,ied,jsd,jed,npx,npy,npz,ng,ns,istag,jstag,dummy)
2017  type(fv_nest_bc_type_3d), intent(INOUT) :: BC
2018  integer, intent(IN) :: ns, istag, jstag
2019  logical, intent(IN) :: dummy
2020 
2021  integer, intent(IN) :: is, ie, js, je, isd, ied, jsd, jed, npx, npy, npz, ng
2022 
2023  if (bc%allocated) return
2024 
2025 
2026  if (ie == npx-1 .and. .not. dummy) then
2027  allocate(bc%east_t1(ie+1-ns+istag:ied+istag,jsd:jed+jstag,npz))
2028  allocate(bc%east_t0(ie+1-ns+istag:ied+istag,jsd:jed+jstag,npz))
2029  do k=1,npz
2030  do j=jsd,jed+jstag
2031  do i=ie+1-ns+istag,ied+istag
2032  bc%east_t1(i,j,k) = 0.
2033  bc%east_t0(i,j,k) = 0.
2034  enddo
2035  enddo
2036  enddo
2037  else
2038  allocate(bc%east_t1(1,1,npz))
2039  allocate(bc%east_t0(1,1,npz))
2040  end if
2041 
2042  if (js == 1 .and. .not. dummy) then
2043  allocate(bc%south_t1(isd:ied+istag,jsd:js-1+ns,npz))
2044  allocate(bc%south_t0(isd:ied+istag,jsd:js-1+ns,npz))
2045  do k=1,npz
2046  do j=jsd,js-1+ns
2047  do i=isd,ied+istag
2048  bc%south_t1(i,j,k) = 0.
2049  bc%south_t0(i,j,k) = 0.
2050  enddo
2051  enddo
2052  enddo
2053  else
2054  allocate(bc%south_t1(1,1,npz))
2055  allocate(bc%south_t0(1,1,npz))
2056  end if
2057 
2058  if (is == 1 .and. .not. dummy) then
2059  allocate(bc%west_t1(isd:is-1+ns,jsd:jed+jstag,npz))
2060  allocate(bc%west_t0(isd:is-1+ns,jsd:jed+jstag,npz))
2061  do k=1,npz
2062  do j=jsd,jed+jstag
2063  do i=isd,is-1+ns
2064  bc%west_t1(i,j,k) = 0.
2065  bc%west_t0(i,j,k) = 0.
2066  enddo
2067  enddo
2068  enddo
2069  else
2070  allocate(bc%west_t1(1,1,npz))
2071  allocate(bc%west_t0(1,1,npz))
2072  end if
2073 
2074  if (je == npy-1 .and. .not. dummy) then
2075  allocate(bc%north_t1(isd:ied+istag,je+1-ns+jstag:jed+jstag,npz))
2076  allocate(bc%north_t0(isd:ied+istag,je+1-ns+jstag:jed+jstag,npz))
2077  do k=1,npz
2078  do j=je+1-ns+jstag,jed+jstag
2079  do i=isd,ied+istag
2080  bc%north_t1(i,j,k) = 0.
2081  bc%north_t0(i,j,k) = 0.
2082  enddo
2083  enddo
2084  enddo
2085  else
2086  allocate(bc%north_t1(1,1,npz))
2087  allocate(bc%north_t0(1,1,npz))
2088  end if
2089 
2090  bc%allocated = .true.
2091 
2092 end subroutine allocate_fv_nest_bc_type_3d
2093 
2094 subroutine deallocate_fv_nest_bc_type_3d(BC)
2096  type(fv_nest_bc_type_3d) :: BC
2097 
2098  if (.not. bc%allocated) return
2099 
2100  deallocate(bc%north_t1)
2101  deallocate(bc%south_t1)
2102  deallocate(bc%west_t1)
2103  deallocate(bc%east_t1)
2104 
2105  if (allocated(bc%north_t0)) then
2106  deallocate(bc%north_t0)
2107  deallocate(bc%south_t0)
2108  deallocate(bc%west_t0)
2109  deallocate(bc%east_t0)
2110  endif
2111 
2112  bc%allocated = .false.
2113 
2114 end subroutine deallocate_fv_nest_bc_type_3d
2115 
2116 
2117 end module fv_arrays_mod
subroutine allocate_fv_nest_bc_type_3d_atm(BC, Atm, ns, istag, jstag, dummy)
Definition: fv_arrays.F90:1984
The type &#39;fv_grid_type&#39; is made up of grid-dependent information from fv_grid_tools and fv_grid_utils...
Definition: fv_arrays.F90:120
real, parameter real_snan
Definition: fv_arrays.F90:50
real, parameter i4_in
Definition: fv_arrays.F90:52
integer, parameter max_step
Definition: fv_arrays.F90:43
integer, parameter, public r_grid
Definition: fv_arrays.F90:35
&#39;allocate_fv_nest_BC_type&#39; is an interface to subroutines that allocate the &#39;fv_nest_BC_type&#39; structu...
Definition: fv_arrays.F90:1126
&#39;deallocate_fv_nest_BC_type&#39; is an interface to a subroutine that deallocates the &#39;fv_nest_BC_type&#39; s...
Definition: fv_arrays.F90:1133
subroutine allocate_fv_nest_bc_type_3d(BC, is, ie, js, je, isd, ied, jsd, jed, npx, npy, npz, ng, ns, istag, jstag, dummy)
Definition: fv_arrays.F90:2016
The module &#39;fv_arrays&#39; contains the &#39;fv_atmos_type&#39; and associated datatypes.
Definition: fv_arrays.F90:24
real, parameter real_big
Definition: fv_arrays.F90:49
subroutine allocate_fv_atmos_type(Atm, isd_in, ied_in, jsd_in, jed_in, is_in, ie_in, js_in, je_in, npx_in, npy_in, npz_in, ndims_in, ncnst_in, nq_in, ng_in, dummy, alloc_2d, ngrids_in)
The subroutine &#39;allocate_fv_atmos_type&#39; allocates the fv_atmos_type.
Definition: fv_arrays.F90:1324
subroutine deallocate_fv_nest_bc_type_3d(BC)
Definition: fv_arrays.F90:2095
subroutine deallocate_fv_atmos_type(Atm)
The subroutine &#39;deallocate_fv_atmos_type&#39; deallocates the fv_atmos_type.
Definition: fv_arrays.F90:1770