Use new MOM_coupler_types interfaces with tracers

Modified the various MOM6 tracer packages to use the MOM_coupler_types module along with new and simpler MOM_coupler_types interfaces. Also modified the module use statements in some of the core and diagnostics modules to access the coupler_types routines via MOM_coupler_types module. All answers are bitwise identical.
NOAA-GFDL · Jan 24, 2021 · 9ad5dcb · 9ad5dcb
1 parent 74f2290
commit 9ad5dcb
Show file tree

Hide file tree

Showing 15 changed files with 79 additions and 116 deletions.
diff --git a/src/core/MOM.F90 b/src/core/MOM.F90
@@ -47,7 +47,6 @@ module MOM
 use MOM_time_manager,         only : operator(-), operator(>), operator(*), operator(/)
 use MOM_time_manager,         only : operator(>=), operator(==), increment_date
 use MOM_unit_tests,           only : unit_tests
-use coupler_types_mod,        only : coupler_type_send_data, coupler_1d_bc_type, coupler_type_spawn
 
 ! MOM core modules
 use MOM_ALE,                   only : ALE_init, ALE_end, ALE_main, ALE_CS, adjustGridForIntegrity

diff --git a/src/core/MOM_forcing_type.F90 b/src/core/MOM_forcing_type.F90
@@ -4,13 +4,15 @@ module MOM_forcing_type
 ! This file is part of MOM6. See LICENSE.md for the license.
 
 use MOM_array_transform, only : rotate_array, rotate_vector, rotate_array_pair
-use MOM_debugging,     only : hchksum, uvchksum
+use MOM_coupler_types, only : coupler_2d_bc_type, coupler_type_destructor
+use MOM_coupler_types, only : coupler_type_increment_data, coupler_type_initialized
 use MOM_cpu_clock,     only : cpu_clock_id, cpu_clock_begin, cpu_clock_end, CLOCK_ROUTINE
+use MOM_debugging,     only : hchksum, uvchksum
 use MOM_diag_mediator, only : post_data, register_diag_field, register_scalar_field
 use MOM_diag_mediator, only : time_type, diag_ctrl, safe_alloc_alloc, query_averaging_enabled
 use MOM_diag_mediator, only : enable_averages, enable_averaging, disable_averaging
-use MOM_error_handler, only : MOM_error, FATAL, WARNING
 use MOM_EOS,           only : calculate_density_derivs, EOS_domain
+use MOM_error_handler, only : MOM_error, FATAL, WARNING
 use MOM_file_parser,   only : get_param, log_param, log_version, param_file_type
 use MOM_grid,          only : ocean_grid_type
 use MOM_opacity,       only : sumSWoverBands, optics_type, extract_optics_slice, optics_nbands
@@ -19,10 +21,6 @@ module MOM_forcing_type
 use MOM_variables,     only : surface, thermo_var_ptrs
 use MOM_verticalGrid,  only : verticalGrid_type
 
-use coupler_types_mod, only : coupler_2d_bc_type, coupler_type_spawn
-use coupler_types_mod, only : coupler_type_increment_data, coupler_type_initialized
-use coupler_types_mod, only : coupler_type_copy_data, coupler_type_destructor
-
 implicit none ; private
 
 #include <MOM_memory.h>

diff --git a/src/core/MOM_variables.F90 b/src/core/MOM_variables.F90
@@ -4,16 +4,14 @@ module MOM_variables
 ! This file is part of MOM6. See LICENSE.md for the license.
 
 use MOM_array_transform, only : rotate_array, rotate_vector
-use MOM_domains, only : MOM_domain_type, get_domain_extent, group_pass_type
-use MOM_debugging, only : hchksum
+use MOM_coupler_types, only : coupler_1d_bc_type, coupler_2d_bc_type
+use MOM_coupler_types, only : coupler_type_spawn, coupler_type_destructor, coupler_type_initialized
+use MOM_debugging,     only : hchksum
+use MOM_domains,       only : MOM_domain_type, get_domain_extent, group_pass_type
+use MOM_EOS,           only : EOS_type
 use MOM_error_handler, only : MOM_error, FATAL
-use MOM_grid, only : ocean_grid_type
-use MOM_verticalGrid, only : verticalGrid_type
-use MOM_EOS, only : EOS_type
-
-use coupler_types_mod, only : coupler_1d_bc_type, coupler_2d_bc_type
-use coupler_types_mod, only : coupler_type_spawn, coupler_type_destructor
-use coupler_types_mod, only : coupler_type_initialized
+use MOM_grid,          only : ocean_grid_type
+use MOM_verticalGrid,  only : verticalGrid_type
 
 implicit none ; private
 
@@ -471,8 +469,7 @@ subroutine rotate_surface_state(sfc_state_in, G_in, sfc_state, G, turns)
 
   ! TODO: tracer field rotation
   if (coupler_type_initialized(sfc_state_in%tr_fields)) &
-    call MOM_error(FATAL, "Rotation of surface state tracers is not yet " &
-        // "implemented.")
+    call MOM_error(FATAL, "Rotation of surface state tracers is not yet implemented.")
 end subroutine rotate_surface_state
 
 !> Allocates the arrays contained within a BT_cont_type and initializes them to 0.

diff --git a/src/diagnostics/MOM_diagnostics.F90 b/src/diagnostics/MOM_diagnostics.F90
@@ -6,6 +6,7 @@ module MOM_diagnostics
 ! This file is part of MOM6. See LICENSE.md for the license.
 
 use MOM_coms,              only : reproducing_sum
+use MOM_coupler_types,     only : coupler_type_send_data
 use MOM_density_integrals, only : int_density_dz
 use MOM_diag_mediator,     only : post_data, get_diag_time_end
 use MOM_diag_mediator,     only : register_diag_field, register_scalar_field
@@ -30,7 +31,6 @@ module MOM_diagnostics
 use MOM_variables,         only : accel_diag_ptrs, cont_diag_ptrs, surface
 use MOM_verticalGrid,      only : verticalGrid_type, get_thickness_units
 use MOM_wave_speed,        only : wave_speed, wave_speed_CS, wave_speed_init
-use coupler_types_mod,     only : coupler_type_send_data
 
 implicit none ; private
 

diff --git a/src/tracer/DOME_tracer.F90 b/src/tracer/DOME_tracer.F90
@@ -3,6 +3,7 @@ module DOME_tracer
 
 ! This file is part of MOM6. See LICENSE.md for the license.
 
+use MOM_coupler_types,   only : set_coupler_type_data, atmos_ocn_coupler_flux
 use MOM_diag_mediator,   only : diag_ctrl
 use MOM_error_handler,   only : MOM_error, FATAL, WARNING
 use MOM_file_parser,     only : get_param, log_param, log_version, param_file_type
@@ -21,9 +22,6 @@ module DOME_tracer
 use MOM_variables,       only : surface
 use MOM_verticalGrid,    only : verticalGrid_type
 
-use coupler_types_mod,      only : coupler_type_set_data, ind_csurf
-use atmos_ocean_fluxes_mod, only : aof_set_coupler_flux
-
 implicit none ; private
 
 #include <MOM_memory.h>
@@ -48,7 +46,7 @@ module DOME_tracer
   real :: land_val(NTR) = -1.0 !< The value of tr used where land is masked out.
   logical :: use_sponge    !< If true, sponges may be applied somewhere in the domain.
 
-  integer, dimension(NTR) :: ind_tr !< Indices returned by aof_set_coupler_flux if it is used and the
+  integer, dimension(NTR) :: ind_tr !< Indices returned by atmos_ocn_coupler_flux if it is used and the
                                     !! surface tracer concentrations are to be provided to the coupler.
 
   type(diag_ctrl), pointer :: diag => NULL() !< A structure that is used to
@@ -130,7 +128,7 @@ function register_DOME_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
     ! values to the coupler (if any).  This is meta-code and its arguments will
     ! currently (deliberately) give fatal errors if it is used.
     if (CS%coupled_tracers) &
-      CS%ind_tr(m) = aof_set_coupler_flux(trim(name)//'_flux', &
+      CS%ind_tr(m) = atmos_ocn_coupler_flux(trim(name)//'_flux', &
           flux_type=' ', implementation=' ', caller="register_DOME_tracer")
   enddo
 
@@ -359,9 +357,8 @@ subroutine DOME_tracer_surface_state(sfc_state, h, G, GV, CS)
     do m=1,NTR
       !   This call loads the surface values into the appropriate array in the
       ! coupler-type structure.
-      call coupler_type_set_data(CS%tr(:,:,1,m), CS%ind_tr(m), ind_csurf, &
-                   sfc_state%tr_fields, idim=(/isd, is, ie, ied/), &
-                   jdim=(/jsd, js, je, jed/) )
+      call set_coupler_type_data(CS%tr(:,:,1,m), CS%ind_tr(m), sfc_state%tr_fields, &
+                   idim=(/isd, is, ie, ied/), jdim=(/jsd, js, je, jed/) )
     enddo
   endif
 

diff --git a/src/tracer/ISOMIP_tracer.F90 b/src/tracer/ISOMIP_tracer.F90
@@ -11,6 +11,7 @@ module ISOMIP_tracer
 ! Adapted to the ISOMIP test case by Gustavo Marques, May 2016
 
 use MOM_coms, only : max_across_PEs
+use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux
 use MOM_diag_mediator, only : diag_ctrl
 use MOM_error_handler, only : MOM_error, FATAL, WARNING
 use MOM_file_parser, only : get_param, log_param, log_version, param_file_type
@@ -28,9 +29,6 @@ module ISOMIP_tracer
 use MOM_variables, only : surface
 use MOM_verticalGrid, only : verticalGrid_type
 
-use coupler_types_mod, only : coupler_type_set_data, ind_csurf
-use atmos_ocean_fluxes_mod, only : aof_set_coupler_flux
-
 implicit none ; private
 
 #include <MOM_memory.h>
@@ -51,7 +49,7 @@ module ISOMIP_tracer
   real :: land_val(NTR) = -1.0 !< The value of tr used where land is masked out.
   logical :: use_sponge    !< If true, sponges may be applied somewhere in the domain.
 
-  integer, dimension(NTR) :: ind_tr !< Indices returned by aof_set_coupler_flux
+  integer, dimension(NTR) :: ind_tr !< Indices returned by atmos_ocn_coupler_flux
              !< if it is used and the surface tracer concentrations are to be
              !< provided to the coupler.
 
@@ -135,7 +133,7 @@ function register_ISOMIP_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
     ! values to the coupler (if any).  This is meta-code and its arguments will
     ! currently (deliberately) give fatal errors if it is used.
     if (CS%coupled_tracers) &
-      CS%ind_tr(m) = aof_set_coupler_flux(trim(name)//'_flux', &
+      CS%ind_tr(m) = atmos_ocn_coupler_flux(trim(name)//'_flux', &
           flux_type=' ', implementation=' ', caller="register_ISOMIP_tracer")
   enddo
 
@@ -345,9 +343,8 @@ subroutine ISOMIP_tracer_surface_state(sfc_state, h, G, GV, CS)
     do m=1,ntr
       !   This call loads the surface values into the appropriate array in the
       ! coupler-type structure.
-      call coupler_type_set_data(CS%tr(:,:,1,m), CS%ind_tr(m), ind_csurf, &
-                   sfc_state%tr_fields, idim=(/isd, is, ie, ied/), &
-                   jdim=(/jsd, js, je, jed/) )
+      call set_coupler_type_data(CS%tr(:,:,1,m), CS%ind_tr(m), sfc_state%tr_fields, &
+                   idim=(/isd, is, ie, ied/), jdim=(/jsd, js, je, jed/) )
     enddo
   endif
 

diff --git a/src/tracer/MOM_OCMIP2_CFC.F90 b/src/tracer/MOM_OCMIP2_CFC.F90
@@ -3,6 +3,8 @@ module MOM_OCMIP2_CFC
 
 ! This file is part of MOM6. See LICENSE.md for the license.
 
+use MOM_coupler_types, only : extract_coupler_type_data, set_coupler_type_data
+use MOM_coupler_types, only : atmos_ocn_coupler_flux
 use MOM_diag_mediator, only : diag_ctrl
 use MOM_error_handler, only : MOM_error, FATAL, WARNING
 use MOM_file_parser, only : get_param, log_param, log_version, param_file_type
@@ -21,10 +23,6 @@ module MOM_OCMIP2_CFC
 use MOM_variables, only : surface
 use MOM_verticalGrid, only : verticalGrid_type
 
-use coupler_types_mod, only : ind_flux, ind_alpha, ind_csurf
-use coupler_types_mod, only : coupler_type_extract_data, coupler_type_set_data
-use atmos_ocean_fluxes_mod, only : aof_set_coupler_flux
-
 implicit none ; private
 
 #include <MOM_memory.h>
@@ -71,9 +69,9 @@ module MOM_OCMIP2_CFC
   character(len=16) :: CFC11_name !< CFC11 variable name
   character(len=16) :: CFC12_name !< CFC12 variable name
 
-  integer :: ind_cfc_11_flux  !< Index returned by aof_set_coupler_flux that is used to
+  integer :: ind_cfc_11_flux  !< Index returned by atmos_ocn_coupler_flux that is used to
                               !! pack and unpack surface boundary condition arrays.
-  integer :: ind_cfc_12_flux  !< Index returned by aof_set_coupler_flux that is used to
+  integer :: ind_cfc_12_flux  !< Index returned by atmos_ocn_coupler_flux that is used to
                               !! pack and unpack surface boundary condition arrays.
 
   type(diag_ctrl), pointer :: diag => NULL() !< A structure that is used to regulate
@@ -127,7 +125,7 @@ function register_OCMIP2_CFC(HI, GV, param_file, CS, tr_Reg, restart_CS)
   ! indicies for the CFC11 and CFC12 flux coupling.
   call flux_init_OCMIP2_CFC(CS, verbosity=3)
   if ((CS%ind_cfc_11_flux < 0) .or. (CS%ind_cfc_12_flux < 0)) then
-    ! This is most likely to happen with the dummy version of aof_set_coupler_flux
+    ! This is most likely to happen with the dummy version of atmos_ocn_coupler_flux
     ! used in ocean-only runs.
     call MOM_ERROR(WARNING, "CFCs are currently only set up to be run in " // &
                    " coupled model configurations, and will be disabled.")
@@ -291,18 +289,18 @@ subroutine flux_init_OCMIP2_CFC(CS, verbosity)
 
   ! These calls obtain the indices for the CFC11 and CFC12 flux coupling.  They
   ! can safely be called multiple times.
-  ind_flux(1) = aof_set_coupler_flux('cfc_11_flux', &
-       flux_type = 'air_sea_gas_flux', implementation = 'ocmip2', &
-       param = (/ 9.36e-07, 9.7561e-06 /), &
+  ind_flux(1) = atmos_ocn_coupler_flux('cfc_11_flux', &
+       flux_type = 'air_sea_gas_flux', implementation='ocmip2', &
+       param=(/ 9.36e-07, 9.7561e-06 /), &
        ice_restart_file = default_ice_restart_file, &
        ocean_restart_file = default_ocean_restart_file, &
        caller = "register_OCMIP2_CFC", verbosity=verbosity)
-  ind_flux(2) = aof_set_coupler_flux('cfc_12_flux', &
-       flux_type = 'air_sea_gas_flux', implementation = 'ocmip2', &
+  ind_flux(2) = atmos_ocn_coupler_flux('cfc_12_flux', &
+       flux_type='air_sea_gas_flux', implementation='ocmip2', &
        param = (/ 9.36e-07, 9.7561e-06 /), &
-       ice_restart_file = default_ice_restart_file, &
-       ocean_restart_file = default_ocean_restart_file, &
-       caller = "register_OCMIP2_CFC", verbosity=verbosity)
+       ice_restart_file=default_ice_restart_file, &
+       ocean_restart_file=default_ocean_restart_file, &
+       caller="register_OCMIP2_CFC", verbosity=verbosity)
 
   if (present(CS)) then ; if (associated(CS)) then
     CS%ind_cfc_11_flux = ind_flux(1)
@@ -459,9 +457,9 @@ subroutine OCMIP2_CFC_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US
   ! These two calls unpack the fluxes from the input arrays.
   !   The -GV%Rho0 changes the sign convention of the flux and changes the units
   ! of the flux from [Conc. m s-1] to [Conc. kg m-2 T-1].
-  call coupler_type_extract_data(fluxes%tr_fluxes, CS%ind_cfc_11_flux, ind_flux, CFC11_flux, &
+  call extract_coupler_type_data(fluxes%tr_fluxes, CS%ind_cfc_11_flux, CFC11_flux, &
                                  scale_factor=-G%US%R_to_kg_m3*GV%Rho0*US%T_to_s, idim=idim, jdim=jdim)
-  call coupler_type_extract_data(fluxes%tr_fluxes, CS%ind_cfc_12_flux, ind_flux, CFC12_flux, &
+  call extract_coupler_type_data(fluxes%tr_fluxes, CS%ind_cfc_12_flux, CFC12_flux, &
                                  scale_factor=-G%US%R_to_kg_m3*GV%Rho0*US%T_to_s, idim=idim, jdim=jdim)
 
   ! Use a tridiagonal solver to determine the concentrations after the
@@ -602,14 +600,14 @@ subroutine OCMIP2_CFC_surface_state(sfc_state, h, G, GV, CS)
 
   !   These calls load these values into the appropriate arrays in the
   ! coupler-type structure.
-  call coupler_type_set_data(CFC11_alpha, CS%ind_cfc_11_flux, ind_alpha, &
-                             sfc_state%tr_fields, idim=idim, jdim=jdim)
-  call coupler_type_set_data(CFC11_Csurf, CS%ind_cfc_11_flux, ind_csurf, &
-                             sfc_state%tr_fields, idim=idim, jdim=jdim)
-  call coupler_type_set_data(CFC12_alpha, CS%ind_cfc_12_flux, ind_alpha, &
-                             sfc_state%tr_fields, idim=idim, jdim=jdim)
-  call coupler_type_set_data(CFC12_Csurf, CS%ind_cfc_12_flux, ind_csurf, &
-                             sfc_state%tr_fields, idim=idim, jdim=jdim)
+  call set_coupler_type_data(CFC11_alpha, CS%ind_cfc_11_flux, sfc_state%tr_fields, &
+                             solubility=.true., idim=idim, jdim=jdim)
+  call set_coupler_type_data(CFC11_Csurf, CS%ind_cfc_11_flux, sfc_state%tr_fields, &
+                             idim=idim, jdim=jdim)
+  call set_coupler_type_data(CFC12_alpha, CS%ind_cfc_12_flux, sfc_state%tr_fields, &
+                             solubility=.true., idim=idim, jdim=jdim)
+  call set_coupler_type_data(CFC12_Csurf, CS%ind_cfc_12_flux, sfc_state%tr_fields, &
+                             idim=idim, jdim=jdim)
 
 end subroutine OCMIP2_CFC_surface_state
 

diff --git a/src/tracer/advection_test_tracer.F90 b/src/tracer/advection_test_tracer.F90
@@ -3,6 +3,7 @@ module advection_test_tracer
 
 ! This file is part of MOM6. See LICENSE.md for the license.
 
+use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux
 use MOM_diag_mediator, only : diag_ctrl
 use MOM_error_handler, only : MOM_error, FATAL, WARNING
 use MOM_file_parser, only : get_param, log_param, log_version, param_file_type
@@ -20,9 +21,6 @@ module advection_test_tracer
 use MOM_variables, only : surface
 use MOM_verticalGrid, only : verticalGrid_type
 
-use coupler_types_mod, only : coupler_type_set_data, ind_csurf
-use atmos_ocean_fluxes_mod, only : aof_set_coupler_flux
-
 implicit none ; private
 
 #include <MOM_memory.h>
@@ -51,7 +49,7 @@ module advection_test_tracer
   real :: y_origin !< Parameters describing the test functions
   real :: y_width  !< Parameters describing the test functions
 
-  integer, dimension(NTR) :: ind_tr !< Indices returned by aof_set_coupler_flux if it is used and
+  integer, dimension(NTR) :: ind_tr !< Indices returned by atmos_ocn_coupler_flux if it is used and
                    !! the surface tracer concentrations are to be provided to the coupler.
 
   type(diag_ctrl), pointer :: diag => NULL() !< A structure that is used to
@@ -153,7 +151,7 @@ function register_advection_test_tracer(HI, GV, param_file, CS, tr_Reg, restart_
     ! values to the coupler (if any).  This is meta-code and its arguments will
     ! currently (deliberately) give fatal errors if it is used.
     if (CS%coupled_tracers) &
-      CS%ind_tr(m) = aof_set_coupler_flux(trim(name)//'_flux', &
+      CS%ind_tr(m) = atmos_ocn_coupler_flux(trim(name)//'_flux', &
           flux_type=' ', implementation=' ', caller="register_advection_test_tracer")
   enddo
 
@@ -337,9 +335,8 @@ subroutine advection_test_tracer_surface_state(sfc_state, h, G, GV, CS)
     do m=1,CS%ntr
       !   This call loads the surface values into the appropriate array in the
       ! coupler-type structure.
-      call coupler_type_set_data(CS%tr(:,:,1,m), CS%ind_tr(m), ind_csurf, &
-                   sfc_state%tr_fields, idim=(/isd, is, ie, ied/), &
-                   jdim=(/jsd, js, je, jed/) )
+      call set_coupler_type_data(CS%tr(:,:,1,m), CS%ind_tr(m), sfc_state%tr_fields, &
+                   idim=(/isd, is, ie, ied/), jdim=(/jsd, js, je, jed/) )
     enddo
   endif
 

diff --git a/src/tracer/boundary_impulse_tracer.F90 b/src/tracer/boundary_impulse_tracer.F90
@@ -3,6 +3,7 @@ module boundary_impulse_tracer
 
 ! This file is part of MOM6. See LICENSE.md for the license.
 
+use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux
 use MOM_diag_mediator, only : diag_ctrl
 use MOM_error_handler, only : MOM_error, FATAL, WARNING
 use MOM_file_parser, only : get_param, log_param, log_version, param_file_type
@@ -21,9 +22,6 @@ module boundary_impulse_tracer
 use MOM_variables, only : surface, thermo_var_ptrs
 use MOM_verticalGrid, only : verticalGrid_type
 
-use coupler_types_mod, only : coupler_type_set_data, ind_csurf
-use atmos_ocean_fluxes_mod, only : aof_set_coupler_flux
-
 implicit none ; private
 
 #include <MOM_memory.h>
@@ -43,7 +41,7 @@ module boundary_impulse_tracer
   type(tracer_registry_type), pointer :: tr_Reg => NULL() !< A pointer to the tracer registry
   real, pointer :: tr(:,:,:,:) => NULL() !< The array of tracers used in this subroutine, in g m-3?
   logical :: tracers_may_reinit  !< If true, boundary_impulse can be initialized if not found in restart file
-  integer, dimension(NTR_MAX) :: ind_tr  !< Indices returned by aof_set_coupler_flux if it is used and the
+  integer, dimension(NTR_MAX) :: ind_tr  !< Indices returned by atmos_ocn_coupler_flux if it is used and the
                                          !! surface tracer concentrations are to be provided to the coupler.
 
   integer :: nkml !< Number of layers in mixed layer
@@ -131,7 +129,7 @@ function register_boundary_impulse_tracer(HI, GV, param_file, CS, tr_Reg, restar
     ! values to the coupler (if any).  This is meta-code and its arguments will
     ! currently (deliberately) give fatal errors if it is used.
     if (CS%coupled_tracers) &
-      CS%ind_tr(m) = aof_set_coupler_flux(trim(var_name)//'_flux', &
+      CS%ind_tr(m) = atmos_ocn_coupler_flux(trim(var_name)//'_flux', &
           flux_type=' ', implementation=' ', caller="register_boundary_impulse_tracer")
   enddo
   ! Register remaining source time as a restart field
@@ -356,9 +354,8 @@ subroutine boundary_impulse_tracer_surface_state(sfc_state, h, G, GV, CS)
     do m=1,CS%ntr
       !   This call loads the surface values into the appropriate array in the
       ! coupler-type structure.
-      call coupler_type_set_data(CS%tr(:,:,1,m), CS%ind_tr(m), ind_csurf, &
-                   sfc_state%tr_fields, idim=(/isd, is, ie, ied/), &
-                   jdim=(/jsd, js, je, jed/) )
+      call set_coupler_type_data(CS%tr(:,:,1,m), CS%ind_tr(m), sfc_state%tr_fields, &
+                   idim=(/isd, is, ie, ied/), jdim=(/jsd, js, je, jed/) )
     enddo
   endif