Line data Source code
1 : module edyn_init
2 : !
3 : ! Initialize edynamo
4 : !
5 : use cam_logfile, only: iulog
6 : use cam_abortutils, only: endrun
7 : use spmd_utils, only: masterproc
8 : use edynamo, only: alloc_edyn, debug_hist
9 :
10 : implicit none
11 :
12 : private
13 : public :: edynamo_init
14 :
15 : contains
16 : !-----------------------------------------------------------------------
17 1536 : subroutine edynamo_init(mpicomm,ionos_debug_hist)
18 :
19 : !
20 : ! One-time initialization, called from ionosphere_init
21 : ! before dyn_init and phys_init
22 : !
23 : use edyn_maggrid, only: set_maggrid, gmlat, nmlonp1, nmlat, nmlath, nmlev
24 : use edyn_mpi, only: mp_exchange_tasks
25 : use edyn_mpi, only: mp_distribute_mag
26 : use edyn_phys_grid, only: edyn_phys_grid_init
27 : use edyn_solve, only: edyn_solve_init
28 :
29 : !
30 : ! Args:
31 : integer, intent(in) :: mpicomm
32 : logical, intent(in) :: ionos_debug_hist
33 :
34 768 : debug_hist = ionos_debug_hist
35 :
36 768 : if (masterproc) then
37 2 : write(iulog,"('Enter edynamo_init:')")
38 : endif
39 :
40 768 : call set_maggrid () ! set parameter-based global magnetic grid
41 :
42 768 : call edyn_solve_init
43 :
44 768 : call mp_distribute_mag(nmlonp1, nmlat, nmlath, nmlev)
45 :
46 768 : call register_grids()
47 768 : call mp_exchange_tasks(mpicomm, 0, gmlat) ! single arg is iprint
48 :
49 768 : call alloc_edyn() ! allocate dynamo arrays
50 :
51 768 : call edyn_phys_grid_init()
52 :
53 768 : call add_fields() ! add fields to WACCM history master list
54 :
55 768 : end subroutine edynamo_init
56 :
57 : !-----------------------------------------------------------------------
58 768 : subroutine add_fields
59 : use cam_history, only: addfld, horiz_only, add_default
60 : use phys_control, only: phys_getopts !Method used to get flag for waccmx ionosphere output variables
61 :
62 : logical :: history_waccmx
63 :
64 1536 : call addfld ('PED_MAG' ,(/ 'lev' /), 'I', 'S/m ','Pedersen Conductivity' ,gridname='gmag_grid')
65 1536 : call addfld ('HAL_MAG' ,(/ 'lev' /), 'I', 'S/m ','Hall Conductivity' ,gridname='gmag_grid')
66 768 : call addfld ('PHIHM' , horiz_only, 'I', 'VOLTS','High Latitude Electric Potential' ,gridname='gmag_grid')
67 768 : call addfld ('PHIM2D' , horiz_only, 'I', 'VOLTS','PHIM2D: Electric Potential' ,gridname='gmag_grid')
68 768 : call addfld ('ED1' , horiz_only, 'I', 'V/m ','ED1: Eastward Electric Field' ,gridname='gmag_grid')
69 768 : call addfld ('ED2' , horiz_only, 'I', 'V/m ','ED2: Equatorward Electric Field' ,gridname='gmag_grid')
70 1536 : call addfld ('PHIM3D' ,(/ 'lev' /), 'I', 'VOLTS','PHIM3D: 3d Electric Potential' ,gridname='gmag_grid')
71 1536 : call addfld ('ED13D' ,(/ 'lev' /), 'I', 'V/m ','ED13D: Eastward Electric Field' ,gridname='gmag_grid')
72 1536 : call addfld ('ED23D' ,(/ 'lev' /), 'I', 'V/m ','ED23D: Equatorward Electric Field',gridname='gmag_grid')
73 1536 : call addfld ('ZPOT_MAG' ,(/ 'lev' /), 'I', 'cm ','Geopotential on mag grid (h0 min)',gridname='gmag_grid')
74 :
75 768 : call addfld ('HALL_CONDUCTANCE',horiz_only, 'I', 'S','Hall Conductance', gridname='gmag_grid')
76 768 : call addfld ('PED_CONDUCTANCE', horiz_only, 'I', 'S','Pedersen Conductance',gridname='gmag_grid')
77 :
78 1536 : call addfld ('POTEN',(/ 'lev' /), 'I', 'Volts','POTEN: Electric Potential', gridname='geo_grid')
79 :
80 1536 : call addfld ('EX',(/'lev'/),'I','V/m','Geographic zonal component of Electric Field', gridname='geo_grid')
81 1536 : call addfld ('EY',(/'lev'/),'I','V/m','Geographic meridional component of Electric Field', gridname='geo_grid')
82 1536 : call addfld ('EZ',(/'lev'/),'I','V/m','Vertical component of Electric Field', gridname='geo_grid')
83 :
84 768 : call addfld ('BMOD', horiz_only, 'I', 'gauss','magnitude of magnetic field',gridname='geo_grid')
85 768 : call addfld ('XB', horiz_only, 'I', 'gauss','northward component of magnetic field',gridname='geo_grid')
86 768 : call addfld ('YB', horiz_only, 'I', 'gauss','eastward component of magnetic field',gridname='geo_grid')
87 768 : call addfld ('ZB', horiz_only, 'I', 'gauss','downward component of magnetic field',gridname='geo_grid')
88 :
89 768 : if (debug_hist) then
90 0 : call addfld ('EPHI3D' ,(/ 'lev' /), 'I', ' ','EPHI3D' ,gridname='gmag_grid')
91 0 : call addfld ('ELAM3D' ,(/ 'lev' /), 'I', ' ','ELAM3D' ,gridname='gmag_grid')
92 0 : call addfld ('EMZ3D' ,(/ 'lev' /), 'I', ' ','EMZ3D' ,gridname='gmag_grid')
93 :
94 0 : call addfld ('ADOTV1_MAG',(/ 'lev' /), 'I', ' ','ADOTV1 on mag grid' ,gridname='gmag_grid')
95 0 : call addfld ('ADOTV2_MAG',(/ 'lev' /), 'I', ' ','ADOTV2 on mag grid' ,gridname='gmag_grid')
96 0 : call addfld ('adota1_mag_a', horiz_only, 'I', ' ','adota1_mag_a',gridname='gmag_grid')
97 0 : call addfld ('ZIGM11_a', horiz_only, 'I', ' ','ZIGM11_a',gridname='gmag_grid')
98 0 : call addfld ('EDYN_ZIGM11_0', horiz_only, 'I', ' ','EDYN_ZIGM11_0',gridname='gmag_grid')
99 0 : call addfld ('EDYN_ZIGM11', horiz_only, 'I', ' ','EDYN_ZIGM11',gridname='gmag_grid')
100 0 : call addfld ('EDYN_ZIGM22', horiz_only, 'I', ' ','EDYN_ZIGM22',gridname='gmag_grid')
101 0 : call addfld ('EDYN_ZIGMC' , horiz_only, 'I', ' ','EDYN_ZIGMC' ,gridname='gmag_grid')
102 0 : call addfld ('EDYN_ZIGM2' , horiz_only, 'I', ' ','EDYN_ZIGM2' ,gridname='gmag_grid')
103 0 : call addfld ('EDYN_RIM1' , horiz_only, 'I', ' ','EDYN_RIM1' ,gridname='gmag_grid')
104 0 : call addfld ('EDYN_RIM2' , horiz_only, 'I', ' ','EDYN_RIM2' ,gridname='gmag_grid')
105 :
106 : ! rjac: scaled derivatives of geomagnetic coords wrt geographic coordinates.
107 0 : call addfld ('RJAC11',(/'lev'/), 'I', '1','cos(thetas)/cos(theta)*d(lamdas)/d(lamda)' ,gridname='geo_grid')
108 0 : call addfld ('RJAC12',(/'lev'/), 'I', '1','cos(thetas)*d(lamdas)/d(theta)' ,gridname='geo_grid')
109 0 : call addfld ('RJAC21',(/'lev'/), 'I', '1','1./cos(theta)*d(thetas)/d(lamda)' ,gridname='geo_grid')
110 0 : call addfld ('RJAC22',(/'lev'/), 'I', '1','d(thetas)/d(theta)' ,gridname='geo_grid')
111 : endif
112 : !-------------------------------------------------------------------------------
113 : ! Set default values for ionosphere history variables
114 : !-------------------------------------------------------------------------------
115 768 : call phys_getopts(history_waccmx_out=history_waccmx)
116 :
117 768 : if (history_waccmx) then
118 768 : call add_default ('PED_CONDUCTANCE', 1, ' ')
119 768 : call add_default ('HALL_CONDUCTANCE' , 1, ' ')
120 : end if
121 :
122 768 : end subroutine add_fields
123 : !-----------------------------------------------------------------------
124 :
125 768 : subroutine register_grids()
126 :
127 768 : use cam_grid_support, only: horiz_coord_t, horiz_coord_create, iMap
128 : use cam_grid_support, only: cam_grid_register
129 : use edyn_mpi, only: mlat0, mlat1, mlon0, omlon1, ntask, mytid
130 : use edyn_mpi, only: lat0, lat1, lon0, lon1
131 : use edyn_maggrid, only: gmlat, gmlon, nmlat, nmlon
132 : use edyn_geogrid, only: glat, glon, nlat, nlon
133 :
134 : integer, parameter :: mag_decomp = 111 ! Must be unique within CAM
135 : integer, parameter :: geo_decomp = 112 ! Must be unique within CAM
136 :
137 : type(horiz_coord_t), pointer :: lat_coord => null()
138 : type(horiz_coord_t), pointer :: lon_coord => null()
139 : integer(iMap), pointer :: grid_map(:,:) => null()
140 : integer(iMap), pointer :: coord_map(:) => null()
141 : integer :: i, j, ind
142 :
143 768 : if (mytid>=ntask) then
144 :
145 0 : if (mlon0/=1) then
146 0 : call endrun('register_grids: mlon0 needs to be 1 on inactive PEs')
147 : end if
148 0 : if (omlon1/=0) then
149 0 : call endrun('register_grids: omlon1 needs to be 0 on inactive PEs')
150 : end if
151 0 : if (mlat0/=1) then
152 0 : call endrun('register_grids: mlat0 needs to be 1 on inactive PEs')
153 : end if
154 0 : if (mlat1/=0) then
155 0 : call endrun('register_grids: mlat1 needs to be 0 on inactive PEs')
156 : end if
157 :
158 0 : if (lon0/=1) then
159 0 : call endrun('register_grids: lon0 needs to be 1 on inactive PEs')
160 : end if
161 0 : if (lon1/=0) then
162 0 : call endrun('register_grids: lon1 needs to be 0 on inactive PEs')
163 : end if
164 0 : if (lat0/=1) then
165 0 : call endrun('register_grids: lat0 needs to be 1 on inactive PEs')
166 : end if
167 0 : if (lat1/=0) then
168 0 : call endrun('register_grids: lat1 needs to be 0 on inactive PEs')
169 : end if
170 :
171 : endif
172 :
173 2304 : allocate(grid_map(4, ((omlon1 - mlon0 + 1) * (mlat1 - mlat0 + 1))))
174 768 : ind = 0
175 3096 : do i = mlat0, mlat1
176 18616 : do j = mlon0, omlon1
177 15520 : ind = ind + 1
178 15520 : grid_map(1, ind) = j
179 15520 : grid_map(2, ind) = i
180 15520 : grid_map(3, ind) = j
181 17848 : grid_map(4, ind) = i
182 : end do
183 : end do
184 :
185 2304 : allocate(coord_map(mlat1 - mlat0 + 1))
186 768 : if (mlon0==1) then
187 452 : coord_map = (/ (i, i = mlat0, mlat1) /)
188 : else
189 2838 : coord_map = 0
190 : end if
191 : lat_coord => horiz_coord_create('mlat', '', nmlat, 'latitude', &
192 0 : 'degrees_north', mlat0, mlat1, gmlat(mlat0:mlat1), &
193 768 : map=coord_map)
194 768 : nullify(coord_map)
195 :
196 2304 : allocate(coord_map(omlon1 - mlon0 + 1))
197 768 : if (mlat0==1) then
198 344 : coord_map = (/ (i, i = mlon0, omlon1) /)
199 : else
200 5704 : coord_map = 0
201 : end if
202 : lon_coord => horiz_coord_create('mlon', '', nmlon, 'longitude', &
203 0 : 'degrees_east', mlon0, omlon1, gmlon(mlon0:omlon1), &
204 768 : map=coord_map)
205 768 : deallocate(coord_map)
206 : nullify(coord_map)
207 :
208 : call cam_grid_register('gmag_grid', mag_decomp, lat_coord, lon_coord, &
209 768 : grid_map, unstruct=.false.)
210 768 : nullify(grid_map)
211 :
212 :
213 : ! for the Oplus geo grid
214 2304 : allocate(grid_map(4, ((lon1 - lon0 + 1) * (lat1 - lat0 + 1))))
215 768 : ind = 0
216 3072 : do i = lat0, lat1
217 30720 : do j = lon0, lon1
218 27648 : ind = ind + 1
219 27648 : grid_map(1, ind) = j
220 27648 : grid_map(2, ind) = i
221 27648 : grid_map(3, ind) = j
222 29952 : grid_map(4, ind) = i
223 : end do
224 : end do
225 :
226 2304 : allocate(coord_map(lat1 - lat0 + 1))
227 768 : if (lon0==1) then
228 448 : coord_map = (/ (i, i = lat0, lat1) /)
229 : else
230 2816 : coord_map = 0
231 : end if
232 : lat_coord => horiz_coord_create('glat', '', nlat, 'latitude', &
233 0 : 'degrees_north', lat0, lat1, glat(lat0:lat1), &
234 768 : map=coord_map)
235 768 : nullify(coord_map)
236 :
237 2304 : allocate(coord_map(lon1 - lon0 + 1))
238 768 : if (lat0==1) then
239 600 : coord_map = (/ (i, i = lon0, lon1) /)
240 : else
241 9672 : coord_map = 0
242 : end if
243 :
244 : lon_coord => horiz_coord_create('glon', '', nlon, 'longitude', &
245 0 : 'degrees_east', lon0, lon1, glon(lon0:lon1), &
246 768 : map=coord_map)
247 768 : deallocate(coord_map)
248 : nullify(coord_map)
249 :
250 : call cam_grid_register('geo_grid', geo_decomp, lat_coord, lon_coord, &
251 768 : grid_map, unstruct=.false.)
252 768 : nullify(grid_map)
253 :
254 768 : end subroutine register_grids
255 :
256 : !-----------------------------------------------------------------------
257 : end module edyn_init
|