!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !! !! !! GNU General Public License !! !! !! !! This file is part of the Flexible Modeling System (FMS). !! !! !! !! FMS is free software; you can redistribute it and/or modify !! !! it and are expected to follow the terms of the GNU General Public !! !! License as published by the Free Software Foundation. !! !! !! !! FMS is distributed in the hope that it will be useful, !! !! but WITHOUT ANY WARRANTY; without even the implied warranty of !! !! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the !! !! GNU General Public License for more details. !! !! !! !! You should have received a copy of the GNU General Public License !! !! along with FMS; if not, write to: !! !! Free Software Foundation, Inc. !! !! 59 Temple Place, Suite 330 !! !! Boston, MA 02111-1307 USA !! !! or see: !! !! http://www.gnu.org/licenses/gpl.txt !! !! !! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! #include ! ! ! Richard D. Slater ! ! ! John P. Dunne ! ! ! ! Ocean Carbon Model Intercomparison Study II: Biotic module ! ! ! ! Implementation of routines to solve the OCMIP-2 Biotic ! simulations as outlined in the Biotic-HOWTO documentation, ! revision 1.7, 1999/10/05. ! ! ! ! ! http://www.ipsl.jussieu.fr/OCMIP/phase2/simulations/Biotic/HOWTO-Biotic.html ! ! ! ! Press, W. H., S. A. Teukosky, W. T. Vetterling, B. P. Flannery, 1992. ! Numerical Recipes in FORTRAN, Second Edition, Cambridge University Press. ! ! ! ! Enting, I.G., T. M. L. Wigley, M. Heimann, 1994. Future Emissions ! and concentrations of carbon dioxide: key ocean / atmosphere / ! land analyses, CSIRO Aust. Div. Atmos. Res. Tech. Pap. No. 31, ! 118 pp. ! ! ! ! !------------------------------------------------------------------ ! ! Module ocmip2_biotic_mod ! ! Implementation of routines to solve the OCMIP-2 Biotic ! simulations as outlined in the Biotic-HOWTO documentation, ! revision 1.7, 1999/10/05. ! !------------------------------------------------------------------ ! module ocmip2_biotic_mod !{ ! !------------------------------------------------------------------ ! ! Global definitions ! !------------------------------------------------------------------ ! ! !---------------------------------------------------------------------- ! ! Modules ! !---------------------------------------------------------------------- ! use time_manager_mod, only: time_type use diag_manager_mod, only: send_data use field_manager_mod, only: fm_field_name_len, fm_path_name_len, fm_string_len use field_manager_mod, only: fm_get_length, fm_get_value, fm_new_value use fms_mod, only: field_exist, file_exist use fms_io_mod, only: register_restart_field, save_restart, restore_state use fms_io_mod, only: restart_file_type use mpp_mod, only: stdout, stdlog, mpp_error, mpp_sum, FATAL use time_manager_mod, only: get_date use time_interp_external_mod, only: time_interp_external use mpp_domains_mod, only: domain2d use constants_mod, only: WTMCO2, WTMO2 use ocean_tpm_util_mod, only: otpm_set_tracer_package, otpm_set_prog_tracer use fm_util_mod, only: fm_util_check_for_bad_fields, fm_util_set_value use fm_util_mod, only: fm_util_get_string, fm_util_get_logical, fm_util_get_integer, fm_util_get_real use fm_util_mod, only: fm_util_get_logical_array, fm_util_get_real_array, fm_util_get_string_array use fm_util_mod, only: fm_util_start_namelist, fm_util_end_namelist use coupler_types_mod, only: ind_alpha, ind_csurf, coupler_2d_bc_type use ocean_types_mod, only: ocean_prog_tracer_type use ocmip2_co2calc_mod, only: ocmip2_co2calc ! !---------------------------------------------------------------------- ! ! force all variables to be "typed" ! !---------------------------------------------------------------------- ! implicit none ! !---------------------------------------------------------------------- ! ! Make all routines and variables private by default ! !---------------------------------------------------------------------- ! private ! !---------------------------------------------------------------------- ! ! Public routines ! !---------------------------------------------------------------------- ! public :: ocmip2_biotic_bbc public :: ocmip2_biotic_end public :: ocmip2_biotic_init public :: ocmip2_biotic_flux_init public :: ocmip2_biotic_sbc public :: ocmip2_biotic_source public :: ocmip2_biotic_start public :: ocmip2_biotic_init_sfc public :: ocmip2_biotic_avg_sfc public :: ocmip2_biotic_sum_sfc public :: ocmip2_biotic_zero_sfc public :: ocmip2_biotic_sfc_end public :: ocmip2_biotic_restart ! !---------------------------------------------------------------------- ! ! Private routines ! !---------------------------------------------------------------------- ! private :: allocate_arrays private :: locate private :: set_array !public :: km_c !public :: compensation_depth !public :: ind_biotic_dop !public :: bc_ptr_prev_mo !public :: bc_ptr_next_mo !public :: bc_frac_prev_mo !public :: bc_frac_next_mo !public :: biotic_flux_dop_a !public :: biotic_flux_dop_d !public :: impvd_dop !public :: conv_dop ! !---------------------------------------------------------------------- ! ! Private parameters ! !---------------------------------------------------------------------- ! character(len=fm_field_name_len), parameter :: package_name = 'ocmip2_biotic' character(len=48), parameter :: mod_name = 'ocmip2_biotic_mod' character(len=48), parameter :: diag_name = 'ocean_ocmip2_biotic' character(len=fm_string_len), parameter :: default_restart_file = 'ocmip2_biotic.res.nc' character(len=fm_string_len), parameter :: default_local_restart_file = 'ocmip2_biotic_local.res.nc' character(len=fm_string_len), parameter :: default_ice_restart_file = 'ice_ocmip2_biotic.res.nc' character(len=fm_string_len), parameter :: default_ocean_restart_file = 'ocmip2_biotic_airsea_flux.res.nc' ! ! coefficients for O2 saturation ! real, parameter :: a_0 = 2.00907 real, parameter :: a_1 = 3.22014 real, parameter :: a_2 = 4.05010 real, parameter :: a_3 = 4.94457 real, parameter :: a_4 = -2.56847e-01 real, parameter :: a_5 = 3.88767 real, parameter :: b_0 = -6.24523e-03 real, parameter :: b_1 = -7.37614e-03 real, parameter :: b_2 = -1.03410e-02 real, parameter :: b_3 = -8.17083e-03 real, parameter :: c_0 = -4.88682e-07 ! !---------------------------------------------------------------------- ! ! Private types ! !---------------------------------------------------------------------- ! ! ! sio2_const = SiO2 concentration (mol/m^3) ! add_phosphate : if true, then add sufficient PO4 to keep ! the predicted PO4 the same as if no depletion ! or changed uptake rate were in effect ! type mask_region_type real, dimension(:,:,:), pointer :: mask => NULL() real, dimension(:), pointer :: elon => NULL() real, dimension(:), pointer :: nlat => NULL() real, dimension(:), pointer :: slat => NULL() real, dimension(:), pointer :: wlon => NULL() logical :: coastal_only real :: factor logical, dimension(:), pointer :: t_mask => NULL() end type mask_region_type type biotic_type !{ logical :: soft_tissue_pump logical :: add_phosphate real :: bio_tau real :: c_2_p real :: ca_remin_depth real :: caco3_2_c real :: compensation_depth real, _ALLOCATABLE, dimension(:,:) :: alpha _NULL real, _ALLOCATABLE, dimension(:,:) :: csurf _NULL real, _ALLOCATABLE, dimension(:,:,:) :: fc _NULL real, _ALLOCATABLE, dimension(:,:,:) :: fca _NULL character(len=fm_string_len) :: local_restart_file real, _ALLOCATABLE, dimension(:,:) :: flux_caco3 _NULL real, _ALLOCATABLE, dimension(:,:) :: flux_poc _NULL real, _ALLOCATABLE, dimension(:,:) :: htotal _NULL integer :: id_alpha = -1 integer :: id_csurf = -1 integer :: id_fc = -1 integer :: id_fca = -1 integer :: id_flux_caco3 = -1 integer :: id_flux_poc = -1 integer :: id_htotal = -1 integer :: id_jca = -1 integer :: id_jdop = -1 integer :: id_jo2 = -1 integer :: id_jpo4 = -1 integer :: id_jpo4_add = -1 integer :: id_jprod = -1 integer :: id_pco2surf = -1 integer :: id_sfc_flux_co2 = -1 integer :: id_sfc_flux_o2 = -1 integer :: id_vstf_alk = -1 integer :: id_vstf_dic = -1 integer :: id_vstf_dop = -1 integer :: id_vstf_o2 = -1 integer :: id_vstf_po4 = -1 integer :: ind_alk integer :: ind_dic integer :: ind_dop integer :: ind_o2 integer :: ind_po4 integer :: ind_co2_flux integer :: ind_o2_flux real, _ALLOCATABLE, dimension(:,:,:) :: jca _NULL real, _ALLOCATABLE, dimension(:,:,:) :: jdop _NULL real, _ALLOCATABLE, dimension(:,:,:) :: jo2 _NULL real, _ALLOCATABLE, dimension(:,:,:) :: jpo4 _NULL real, _ALLOCATABLE, dimension(:,:,:) :: jpo4_add _NULL real, _ALLOCATABLE, dimension(:,:,:) :: jprod _NULL real, _ALLOCATABLE, dimension(:,:,:) :: jprod_norm _NULL real :: kappa integer :: km_c real :: martin_coeff real :: martin_coeff_alt real :: n_2_p character(len=fm_field_name_len) :: name real :: o2_min real :: o_2_p real, _ALLOCATABLE, dimension(:,:) :: pco2surf _NULL real :: r_bio_tau type(mask_region_type) :: r_bio_tau_a type(mask_region_type) :: nut_depl type(mask_region_type) :: norm_remin type(mask_region_type) :: no_caco3 integer :: id_sc_co2 = -1 integer :: id_sc_o2 = -1 real, _ALLOCATABLE, dimension(:,:) :: sc_co2 _NULL real :: sc_co2_0 real :: sc_co2_1 real :: sc_co2_2 real :: sc_co2_3 real, _ALLOCATABLE, dimension(:,:) :: sc_o2 _NULL real :: sc_o2_0 real :: sc_o2_1 real :: sc_o2_2 real :: sc_o2_3 real :: sigma real, _ALLOCATABLE, dimension(:,:) :: sio2 _NULL real :: sio2_const real :: stp_alkalinity real :: stp_salinity real :: stp_temperature real, _ALLOCATABLE, dimension(:,:) :: vstf_alk _NULL real, _ALLOCATABLE, dimension(:,:) :: vstf_dic _NULL real, _ALLOCATABLE, dimension(:,:) :: vstf_dop _NULL real, _ALLOCATABLE, dimension(:,:) :: vstf_o2 _NULL real, _ALLOCATABLE, dimension(:,:) :: vstf_po4 _NULL real, _ALLOCATABLE, dimension(:) :: zfca _NULL real, _ALLOCATABLE, dimension(:) :: zforg _NULL real, _ALLOCATABLE, dimension(:) :: zforg_alt _NULL end type biotic_type !} ! !---------------------------------------------------------------------- ! ! Public variables ! !---------------------------------------------------------------------- ! logical, public :: do_ocmip2_biotic ! !---------------------------------------------------------------------- ! ! Private variables ! !---------------------------------------------------------------------- ! !logical :: no_sbc !logical :: no_source integer :: indsal integer :: indtemp integer :: package_index logical :: module_initialized = .false. character(len=128) :: version = '$Id: ocmip2_biotic.F90,v 16.0.4.1.2.1.2.1.84.1 2009/10/10 00:41:46 nnz Exp $' character(len=128) :: tagname = '$Name: mom4p1_pubrel_dec2009_nnz $' ! !---------------------------------------------------------------------- ! ! Input parameters: ! ! htotal_in = default value for htotal for an initial run ! htotal_scale_lo = scaling parameter to chose htotallo ! htotal_scale_hi = scaling parameter to chose htotalhi ! !---------------------------------------------------------------------- ! real :: htotal_in real, allocatable, dimension(:,:) :: htotal_scale_hi real :: htotal_scale_hi_in real, allocatable, dimension(:,:) :: htotal_scale_lo real :: htotal_scale_lo_in ! !---------------------------------------------------------------------- ! ! Calculated parameters (with possible initial input values): ! !---------------------------------------------------------------------- ! integer :: id_o2_sat integer :: km_c_max integer :: po4_star_id real, allocatable, dimension(:,:,:) :: po4_star_t character*128 :: po4_star_file character*32 :: po4_star_name real, allocatable, dimension(:,:) :: sc_no_term type(biotic_type), allocatable, dimension(:) :: biotic integer :: instances real, allocatable, dimension(:,:) :: o2_saturation real, allocatable, dimension(:) :: tk real, allocatable, dimension(:) :: ts real, allocatable, dimension(:) :: ts2 real, allocatable, dimension(:) :: ts3 real, allocatable, dimension(:) :: ts4 real, allocatable, dimension(:) :: ts5 real, allocatable, dimension(:) :: tt ! for restart integer :: num_restart = 0 type(restart_file_type), allocatable :: restart(:) ! !----------------------------------------------------------------------- ! ! Subroutine and function definitions ! !----------------------------------------------------------------------- ! contains !####################################################################### ! ! ! ! Dynamically allocate arrays ! ! subroutine allocate_arrays(isc, iec, jsc, jec, nk, isd, ied, jsd, jed) !{ implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isc integer, intent(in) :: iec integer, intent(in) :: jsc integer, intent(in) :: jec integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed integer, intent(in) :: nk ! ! local variables ! integer :: i integer :: j integer :: k integer :: l integer :: n ! !----------------------------------------------------------------------- ! start executable code !----------------------------------------------------------------------- ! ! ! global variables ! allocate( sc_no_term(isc:iec,jsc:jec) ) allocate( htotal_scale_lo(isc:iec,jsc:jec) ) allocate( htotal_scale_hi(isc:iec,jsc:jec) ) allocate( o2_saturation(isc:iec,jsc:jec) ) allocate( tt(isc:iec) ) allocate( tk(isc:iec) ) allocate( ts(isc:iec) ) allocate( ts2(isc:iec) ) allocate( ts3(isc:iec) ) allocate( ts4(isc:iec) ) allocate( ts5(isc:iec) ) !allocate( po4_star_t(isd:ied,jsd:jed,km_c_max) ) ! this should be dimensioned as above, but the time_interp routine ! requires that the array dimensions match the datasets dimensions allocate( po4_star_t(isd:ied,jsd:jed,nk) ) ! ! initialize some arrays ! sc_no_term(:,:) = 0.0 htotal_scale_lo(:,:) = 0.0 htotal_scale_hi(:,:) = 0.0 o2_saturation(:,:) = 0.0 tt(:) = 0.0 tk(:) = 0.0 ts(:) = 0.0 ts2(:) = 0.0 ts3(:) = 0.0 ts4(:) = 0.0 ts5(:) = 0.0 po4_star_t(:,:,:) = 0.0 ! ! allocate biotic array elements ! do n = 1, instances !{ allocate( biotic(n)%sc_co2(isc:iec,jsc:jec) ) allocate( biotic(n)%sc_o2(isc:iec,jsc:jec) ) allocate( biotic(n)%htotal(isc:iec,jsc:jec) ) allocate( biotic(n)%csurf(isc:iec,jsc:jec) ) allocate( biotic(n)%alpha(isc:iec,jsc:jec) ) allocate( biotic(n)%pco2surf(isc:iec,jsc:jec) ) allocate( biotic(n)%sio2(isd:ied,jsd:jed) ) allocate( biotic(n)%zforg(nk) ) allocate( biotic(n)%zforg_alt(nk) ) allocate( biotic(n)%zfca(nk) ) allocate( biotic(n)%jprod(isc:iec,jsc:jec,nk) ) allocate( biotic(n)%jpo4(isc:iec,jsc:jec,nk) ) allocate( biotic(n)%jdop(isc:iec,jsc:jec,nk) ) allocate( biotic(n)%jo2(isc:iec,jsc:jec,nk) ) allocate( biotic(n)%jca(isc:iec,jsc:jec,nk) ) allocate( biotic(n)%fc(isc:iec,jsc:jec,nk) ) allocate( biotic(n)%fca(isc:iec,jsc:jec,nk) ) allocate( biotic(n)%flux_poc(isc:iec,jsc:jec) ) allocate( biotic(n)%flux_caco3(isc:iec,jsc:jec) ) allocate( biotic(n)%nut_depl%mask(isc:iec,jsc:jec,12) ) allocate( biotic(n)%no_caco3%mask(isc:iec,jsc:jec,12) ) allocate( biotic(n)%norm_remin%mask(isc:iec,jsc:jec,12) ) allocate( biotic(n)%r_bio_tau_a%mask(isc:iec,jsc:jec,12) ) if (biotic(n)%add_phosphate) then !{ allocate( biotic(n)%jprod_norm(isc:iec,jsc:jec,nk) ) allocate( biotic(n)%jpo4_add(isc:iec,jsc:jec,nk) ) endif !} enddo !} ! ! initialize biotic array elements ! do n = 1, instances !{ do j = jsd, jed !{ do i = isd, ied !{ biotic(n)%sio2(i,j) = 0.0 biotic(n)%sc_co2(i,j) = 0.0 biotic(n)%sc_o2(i,j) = 0.0 enddo !} i enddo !} j do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%htotal(i,j) = 0.0 biotic(n)%csurf(i,j) = 0.0 biotic(n)%alpha(i,j) = 0.0 biotic(n)%pco2surf(i,j) = 0.0 enddo !} i enddo !} j do k = 1, nk !{ biotic(n)%zforg(k) = 0.0 biotic(n)%zforg_alt(k) = 0.0 biotic(n)%zfca(k) = 0.0 enddo !} k do j = jsc, jec !{ do i = isc, iec !{ do k = 1, nk !{ biotic(n)%jprod(i,j,k) = 0.0 biotic(n)%jpo4(i,j,k) = 0.0 biotic(n)%jdop(i,j,k) = 0.0 biotic(n)%jo2(i,j,k) = 0.0 biotic(n)%jca(i,j,k) = 0.0 biotic(n)%fc(i,j,k) = 0.0 biotic(n)%fca(i,j,k) = 0.0 enddo !} k enddo !} i enddo !} j do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%flux_poc(i,j) = 0.0 biotic(n)%flux_caco3(i,j) = 0.0 enddo !} i enddo !} j do j = jsc, jec !{ do i = isc, iec !{ do l = 1, 12 !{ biotic(n)%nut_depl%mask(i,j,l) = 0.0 biotic(n)%no_caco3%mask(i,j,l) = 0.0 biotic(n)%norm_remin%mask(i,j,l) = 0.0 biotic(n)%r_bio_tau_a%mask(i,j,l) = 0.0 enddo !} l enddo !} i enddo !} j if (biotic(n)%add_phosphate) then !{ do j = jsc, jec !{ do i = isc, iec !{ do k = 1, nk !{ biotic(n)%jprod_norm(i,j,k) = 0.0 biotic(n)%jpo4_add(i,j,k) = 0.0 enddo !} k enddo !} i enddo !} j endif !} enddo !} n return end subroutine allocate_arrays !} ! NAME="allocate_arrays" !####################################################################### ! ! ! ! After Numerical recipes: ! ! Given an array XX of length N, and a given value of X, returns a ! value of J such that X is between XX(J) and XX(J+1). XX must be ! monotonic, either increasing or decreasing. J=0 or J=N is ! returned to indicate that X is out of range. ! New features: ! ! If "period" is specified, then the array, xx, is considered ! to be periodic with a period of "period". If "x_in" is out ! of range, then add or subtract "period" once to attempt to ! make "x_in" be in range. ! ! If "nearest" is specified, and true, then return "j" such ! that it is the element of "xx" which is nearest to the value ! of "x_in" (where "x_in" may have been modified by the value ! "period", above). With this option, "j" will be in ! the range 1 <= j <= n. ! ! subroutine locate(xx , n, x_in, j, period, nearest) !{ implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: n real, intent(in) :: x_in real, dimension(n), intent(in) :: xx integer, intent(out) :: j real, optional, intent(in) :: period logical, optional, intent(in) :: nearest ! ! local variables ! integer :: jl, ju, jm real :: x, xt logical :: increasing increasing = xx(1) .lt. xx(n) if (present(period)) then !{ if (increasing) then !{ ! increasing array if (x_in .lt. xx(1)) then !{ ! original value less than start, therefore add period xt = x_in + period if (xt .gt. xx(n)) then !{ ! new value greater than end if (abs(x_in - xx(1)) .gt. abs(xt - xx(n))) then !{ ! new value closer to end than original value to start ! use new value x = xt else !}{ ! original value closer to start than new value to end ! use original value x = x_in endif !} else !}{ ! new value in range ! use new value x = xt endif !} elseif (x_in .gt. xx(n)) then !}{ ! original value greater than end, therefore subtract period xt = x_in - period if (xt .lt. xx(1)) then !{ ! new value less than start if (abs(xt - xx(1)) .lt. abs(x_in - xx(n))) then !{ ! new value closer to start than original value to end ! use new value x = xt else !}{ ! original value closer to end than new value to start ! use original value x = x_in endif !} else !}{ ! new value in range ! use new value x = xt endif !} else !}{ ! original value in range ! use original value x = x_in endif !} else !}{ ! decreasing array if (x_in .gt. xx(1)) then !{ ! original value greater than start, therefore subtract period xt = x_in - period if (xt .lt. xx(n)) then !{ ! new value less than end if (abs(x_in - xx(1)) .gt. abs(xt - xx(n))) then !{ ! new value closer to end than original value to start ! use new value x = xt else !}{ ! original value closer to start than new value to end ! use original value x = x_in endif !} else !}{ ! new value in range ! use new value x = xt endif !} elseif (x_in .lt. xx(n)) then !}{ ! original value less than end, therefore add period xt = x_in + period if (xt .gt. xx(1)) then !{ ! new value greater than start if (abs(xt - xx(1)) .lt. abs(x_in - xx(n))) then !{ ! new value closer to start than original value to end ! use new value x = xt else !}{ ! original value closer to end than new value to start ! use original value x = x_in endif !} else !}{ ! new value in range ! use new value x = xt endif !} else !}{ ! original value in range ! use original value x = x_in endif !} endif !} else !}{ ! no period specified ! use original value x = x_in endif !} jl = 0 ju = n+1 10 continue if (ju - jl .gt. 1) then jm = (ju + jl) / 2 if (increasing .eqv. (x .gt. xx(jm))) then jl = jm else ju = jm endif go to 10 endif j = jl if (present(nearest)) then !{ if (nearest) then !{ if (j .eq. 0) then !{ j = 1 elseif (j .lt. n) then !}{ if (abs(x - xx(j)) .gt. abs(x - xx(j+1))) then !{ j = j + 1 endif !} endif !} endif !} endif !} return end subroutine locate !} ! NAME="locate" !####################################################################### ! ! ! ! calculate the surface boundary conditions ! ! subroutine ocmip2_biotic_bbc(isc, iec, jsc, jec, isd, ied, jsd, jed, T_prog, grid_kmt) !{ ! !----------------------------------------------------------------------- ! modules (have to come first) !----------------------------------------------------------------------- ! implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isc integer, intent(in) :: iec integer, intent(in) :: jsc integer, intent(in) :: jec integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed type(ocean_prog_tracer_type), intent(inout), dimension(:) :: T_prog integer, dimension(isd:,jsd:), intent(in) :: grid_kmt ! !----------------------------------------------------------------------- ! local parameters !----------------------------------------------------------------------- ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_bbc' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: i, j, n, kz ! ! ===================================================================== ! begin executable code ! ===================================================================== ! ! ! set the bottom flux of the column for phosphate to reflect a ! regenerative flux from the sediments where the compensation ! depth is greater than the bottom depth ! do n = 1, instances !{ do j = jsc, jec !{ do i = isc, iec !{ kz = grid_kmt(i,j) if (kz .lt. biotic(n)%km_c .and. kz .gt. 0) then !{ t_prog(biotic(n)%ind_po4)%btf(i,j) = & -biotic(n)%flux_poc(i,j) t_prog(biotic(n)%ind_o2)%btf(i,j) = & biotic(n)%o_2_p * biotic(n)%flux_poc(i,j) t_prog(biotic(n)%ind_dic)%btf(i,j) = & -(biotic(n)%c_2_p * (1.0 + biotic(n)%caco3_2_c) * & biotic(n)%flux_poc(i,j)) t_prog(biotic(n)%ind_alk)%btf(i,j) = & (biotic(n)%n_2_p - 2.0 * biotic(n)%c_2_p * & biotic(n)%caco3_2_c) * biotic(n)%flux_poc(i,j) endif !} enddo !} i enddo !} j enddo !} n return end subroutine ocmip2_biotic_bbc !} ! NAME="ocmip2_biotic_bbc" !####################################################################### ! ! ! ! Clean up various BIOTIC quantities for this run. ! ! subroutine ocmip2_biotic_end(isc, iec, jsc, jec, nk, isd, ied, jsd, jed, & T_prog, grid_dat, grid_tmask, mpp_domain2d, rho_dzt, taup1) !{ ! !----------------------------------------------------------------------- ! modules (have to come first) !----------------------------------------------------------------------- ! implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isc integer, intent(in) :: iec integer, intent(in) :: jsc integer, intent(in) :: jec integer, intent(in) :: nk integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed type(ocean_prog_tracer_type), intent(in), dimension(:) :: T_prog integer, intent(in) :: taup1 real, dimension(isd:,jsd:), intent(in) :: grid_dat real, dimension(isd:,jsd:,:), intent(in) :: grid_tmask type(domain2d), intent(in) :: mpp_domain2d real, dimension(isd:,jsd:,:,:), intent(in) :: rho_dzt ! !----------------------------------------------------------------------- ! local parameters !----------------------------------------------------------------------- ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_end' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: i integer :: j integer :: k integer :: n character(len=fm_field_name_len+1) :: suffix real :: total_alkalinity real :: total_dic real :: total_dop real :: total_o2 real :: total_phosphate integer :: stdoutunit stdoutunit=stdout() ! !----------------------------------------------------------------------- ! statement functions !----------------------------------------------------------------------- ! ! ! ===================================================================== ! begin executable code ! ===================================================================== ! ! ! integrate the total concentrations of some tracers ! for the end of the run ! total_phosphate = 0.0 total_dop = 0.0 total_o2 = 0.0 total_dic = 0.0 total_alkalinity = 0.0 ! ! Use taup1 time index for the start of a run, and taup1 time ! index for the end of a run so that we are integrating the ! same time level and should therefore get identical results ! write (stdoutunit,*) trim(note_header), & 'Global integrals at end of run' do n = 1, instances !{ do k = 1,nk !{ do j = jsc, jec !{ do i = isc, iec !{ total_phosphate = total_phosphate + & t_prog(biotic(n)%ind_po4)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) total_dop = total_dop + & t_prog(biotic(n)%ind_dop)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) total_o2 = total_o2 + & t_prog(biotic(n)%ind_o2)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) total_dic = total_dic + & t_prog(biotic(n)%ind_dic)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) total_alkalinity = total_alkalinity + & t_prog(biotic(n)%ind_alk)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) enddo !} i enddo !} j enddo !} k call mpp_sum(total_phosphate) call mpp_sum(total_dop) call mpp_sum(total_o2) call mpp_sum(total_dic) call mpp_sum(total_alkalinity) write (stdoutunit,*) ' Instance ', trim(biotic(n)%name) write (stdoutunit, & '(/'' Total phosphate = '',es19.12,'' Gmol'')') & total_phosphate * 1.0e-09 write (stdoutunit, & '(/'' Total DOP = '',es19.12,'' Gmol'')') & total_DOP * 1.0e-09 write (stdoutunit, & '(/'' Total O2 = '',es19.12,'' Gmol'')') & total_o2 * 1.0e-09 write (stdoutunit, & '(/'' Total DIC = '',es19.12,'' Gmol'')') & total_dic * 1.0e-09 write (stdoutunit, & '(/'' Total alkalinity = '',es19.12,'' Geq'')') & total_alkalinity * 1.0e-09 write (stdoutunit, & '(/'' Total phosphorous = '',es19.12,'' Gmol'')') & (total_phosphate + total_dop) * 1.0e-09 write (stdoutunit, & '(/'' Total real O2 = '',es19.12,'' Gmol'')') & (total_o2 + biotic(n)%o_2_p * total_phosphate) * 1.0e-09 write (stdoutunit, & '(/'' Total Carbon = '',es19.12,'' Gmol'')') & (total_dic + biotic(n)%c_2_p * total_dop) * 1.0e-09 write (stdoutunit, & '(/'' Total real alkalinity = '',es19.12,'' Geq'')') & (total_alkalinity + biotic(n)%n_2_p * total_phosphate) * 1.0e-09 enddo !} n ! !----------------------------------------------------------------------- ! save out additional information for a restart !----------------------------------------------------------------------- ! write(stdoutunit,*) call ocmip2_biotic_restart do n = 1, instances !{ write(stdoutunit,*) trim(note_header), & 'Writing additional restart information for instance ', & trim(biotic(n)%name) write (stdoutunit,*) trim(note_header), & 'Done writing additional restart information for instance ',& trim(biotic(n)%name) enddo !} n return end subroutine ocmip2_biotic_end !} ! NAME="ocmip2_biotic_end" !####################################################################### ! ! ! Write out restart files registered through register_restart_file ! subroutine ocmip2_biotic_restart(time_stamp) character(len=*), intent(in), optional :: time_stamp integer :: n do n=1, num_restart call save_restart(restart(n), time_stamp) end do end subroutine ocmip2_biotic_restart ! NAME="ocmip2_biotic_restart" !####################################################################### ! ! ! ! Calculate the surface boundary conditions ! ! subroutine ocmip2_biotic_sbc(isc, iec, jsc, jec, nk, isd, ied, jsd, jed, & isc_bnd, iec_bnd, jsc_bnd, jec_bnd, & T_prog, taum1, model_time, grid_tmask, ice_ocean_boundary_fluxes) !{ ! !----------------------------------------------------------------------- ! modules (have to come first) !----------------------------------------------------------------------- ! use coupler_types_mod, only : coupler_2d_bc_type, ind_flux use mpp_mod, only : mpp_sum implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isc integer, intent(in) :: iec integer, intent(in) :: jsc integer, intent(in) :: jec integer, intent(in) :: nk integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed integer, intent(in) :: isc_bnd integer, intent(in) :: iec_bnd integer, intent(in) :: jsc_bnd integer, intent(in) :: jec_bnd type(ocean_prog_tracer_type), intent(inout), dimension(:) :: T_prog integer, intent(in) :: taum1 type(time_type), intent(in) :: model_time real, dimension(isd:,jsd:,:), intent(in) :: grid_tmask type(coupler_2d_bc_type), intent(in) :: ice_ocean_boundary_fluxes ! !----------------------------------------------------------------------- ! local parameters !----------------------------------------------------------------------- ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_sbc' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' real, parameter :: po4_molwt = 31.0 real, parameter :: dop_molwt = 31.0 real, parameter :: o2_molwt = 32.0 real, parameter :: dic_molwt = 12.0 real, parameter :: alk_molwt = 1.0 ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: i integer :: j integer :: n logical :: used ! ! ===================================================================== ! begin executable code ! ===================================================================== ! ! !--------------------------------------------------------------------- ! use the surface fluxes from the coupler ! stf is in mol/m^2/s, flux from coupler is positive upwards !--------------------------------------------------------------------- ! do n = 1, instances !{ t_prog(biotic(n)%ind_dic)%stf(isc:iec,jsc:jec) = & -ice_ocean_boundary_fluxes%bc(biotic(n)%ind_co2_flux)%field(ind_flux)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) t_prog(biotic(n)%ind_o2)%stf(isc:iec,jsc:jec) = & -ice_ocean_boundary_fluxes%bc(biotic(n)%ind_o2_flux)%field(ind_flux)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) enddo !} n ! !----------------------------------------------------------------------- ! Save variables for diagnostics !----------------------------------------------------------------------- ! do n = 1, instances !{ if (biotic(n)%id_sfc_flux_co2 .gt. 0) then !{ used = send_data(biotic(n)%id_sfc_flux_co2, & t_prog(biotic(n)%ind_dic)%stf(:,:), & model_time, rmask = grid_tmask(:,:,1), & is_in=isc, js_in=jsc, ie_in=iec, je_in=jec) endif !} if (biotic(n)%id_sfc_flux_o2 .gt. 0) then !{ used = send_data(biotic(n)%id_sfc_flux_o2, & t_prog(biotic(n)%ind_o2)%stf(:,:), & model_time, rmask = grid_tmask(:,:,1), & is_in=isc, js_in=jsc, ie_in=iec, je_in=jec) endif !} enddo !} n return end subroutine ocmip2_biotic_sbc !} ! NAME="ocmip2_biotic_sbc" !####################################################################### ! ! ! ! Set up any extra fields needed by the ocean-atmosphere gas fluxes ! subroutine ocmip2_biotic_flux_init !{ use atmos_ocean_fluxes_mod, only: aof_set_coupler_flux implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! ! ! local parameters ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_flux_init' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: n character(len=fm_field_name_len) :: name character(len=fm_path_name_len) :: path_to_names character(len=fm_field_name_len+1) :: suffix character(len=256) :: caller_str integer :: stdoutunit stdoutunit=stdout() ! ! First, perform some initialization if this module has not been ! initialized because the normal initialization routine will ! not have been called as part of the normal ocean model ! initialization if this is an Atmosphere pe of a coupled ! model running in concurrent mode ! if (.not. module_initialized) then !{ ! ! Initialize the package ! package_index = otpm_set_tracer_package(package_name, & restart_file = default_restart_file, & caller = trim(mod_name) // '(' // trim(sub_name) // ')') ! ! Check whether to use this package ! path_to_names = '/ocean_mod/tracer_packages/' // trim(package_name) // '/names' instances = fm_get_length(path_to_names) if (instances .lt. 0) then !{ call mpp_error(FATAL, trim(error_header) // ' Could not get number of instances') endif !} ! ! Check some things ! write (stdoutunit,*) if (instances .eq. 0) then !{ write (stdoutunit,*) trim(note_header), ' No instances' do_ocmip2_biotic = .false. else !}{ if (instances .eq. 1) then !{ write (stdoutunit,*) trim(note_header), ' ', instances, ' instance' else !}{ write (stdoutunit,*) trim(note_header), ' ', instances, ' instances' endif !} do_ocmip2_biotic = .true. endif !} module_initialized = .true. endif !} ! ! Return if we don't want to use this package ! if (.not. do_ocmip2_biotic) then !{ return endif !} if (.not. allocated(biotic)) then !{ ! ! allocate storage for biotic array ! allocate ( biotic(instances) ) ! ! loop over the names, saving them into the biotic array ! do n = 1, instances !{ if (fm_get_value(path_to_names, name, index = n)) then !{ biotic(n)%name = name else !}{ write (name,*) n call mpp_error(FATAL, trim(error_header) // & 'Bad field name for index ' // trim(name)) endif !} enddo !} endif !} ! ! Set up the ocean-atmosphere gas flux fields ! caller_str = trim(mod_name) // '(' // trim(sub_name) // ')' do n = 1, instances !{ name = biotic(n)%name if (name(1:1) .eq. '_') then !{ suffix = ' ' else !}{ suffix = '_' // name endif !} ! ! Coupler fluxes ! biotic(n)%ind_co2_flux = aof_set_coupler_flux('co2_flux' // suffix, & flux_type = 'air_sea_gas_flux', implementation = 'ocmip2', & mol_wt = WTMCO2, param = (/ 9.36e-07, 9.7561e-06 /), & ice_restart_file = default_ice_restart_file, & ocean_restart_file = default_ocean_restart_file, & caller = caller_str) biotic(n)%ind_o2_flux = aof_set_coupler_flux('o2_flux' // suffix, & flux_type = 'air_sea_gas_flux', implementation = 'ocmip2', & mol_wt = WTMO2, param = (/ 9.36e-07, 9.7561e-06 /), & ice_restart_file = default_ice_restart_file, & ocean_restart_file = default_ocean_restart_file, & caller = caller_str) ! ! Coupler fields ! enddo !} n return end subroutine ocmip2_biotic_flux_init !} ! !####################################################################### ! ! ! ! Set up any extra fields needed by the tracer packages ! ! Save pointers to various "types", such as Grid and Domains. ! subroutine ocmip2_biotic_init !{ implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! ! ! local parameters ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_init' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' real, parameter :: rho_avg = 1024.5 real, parameter :: sperd = 24.0 * 3600.0 real, parameter :: spery = 365.25 * sperd ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: n character(len=fm_field_name_len) :: name character(len=fm_path_name_len) :: path_to_names character(len=fm_field_name_len+1) :: suffix character(len=fm_field_name_len+3) :: long_suffix logical, dimension(12) :: t_mask character(len=256) :: caller_str character(len=fm_string_len), pointer, dimension(:) :: good_list integer :: stdoutunit stdoutunit=stdout() ! ! Initialize the package ! package_index = otpm_set_tracer_package(package_name, & restart_file = default_restart_file, & caller = trim(mod_name) // '(' // trim(sub_name) // ')') ! ! Check whether to use this package ! path_to_names = '/ocean_mod/tracer_packages/' // trim(package_name) // '/names' instances = fm_get_length(path_to_names) if (instances .lt. 0) then !{ call mpp_error(FATAL, trim(error_header) // ' Could not get number of instances') endif !} ! ! Check some things ! write (stdoutunit,*) if (instances .eq. 0) then !{ write (stdoutunit,*) trim(note_header), ' No instances' do_ocmip2_biotic = .false. else !}{ if (instances .eq. 1) then !{ write (stdoutunit,*) trim(note_header), ' ', instances, ' instance' else !}{ write (stdoutunit,*) trim(note_header), ' ', instances, ' instances' endif !} do_ocmip2_biotic = .true. endif !} module_initialized = .true. ! ! Return if we don't want to use this package ! if (.not. do_ocmip2_biotic) then !{ return endif !} ! ! allocate storage for biotic array ! allocate ( biotic(instances) ) ! ! loop over the names, saving them into the biotic array ! do n = 1, instances !{ if (fm_get_value(path_to_names, name, index = n)) then !{ biotic(n)%name = name else !}{ write (name,*) n call mpp_error(FATAL, trim(error_header) // & 'Bad field name for index ' // trim(name)) endif !} enddo !} ! ! Set up the field input ! caller_str = trim(mod_name) // '(' // trim(sub_name) // ')' do n = 1, instances !{ name = biotic(n)%name if (name(1:1) .eq. '_') then !{ suffix = ' ' long_suffix = ' ' else !}{ suffix = '_' // name long_suffix = ' (' // trim(name) // ')' endif !} ! ! PO4 ! biotic(n)%ind_po4 = otpm_set_prog_tracer('po4' // suffix, & package_name, & longname = 'Phosphate' // trim(long_suffix), & units = 'mol/m^3', flux_units = 'mol/m^2/s', & caller = caller_str) ! ! DOP ! biotic(n)%ind_dop = otpm_set_prog_tracer('dop' // suffix, & package_name, & longname = 'DOP' // trim(long_suffix), & units = 'mol/m^3', flux_units = 'mol/m^2/s', & caller = caller_str) ! ! DIC ! biotic(n)%ind_dic = otpm_set_prog_tracer('dic' // suffix, & package_name, & longname = 'DIC' // trim(long_suffix), & units = 'mol/m^3', flux_units = 'mol/m^2/s', & caller = caller_str) ! ! O2 ! biotic(n)%ind_o2 = otpm_set_prog_tracer('o2' // suffix, & package_name, & longname = 'Oxygen' // trim(long_suffix), & units = 'mol/m^3', flux_units = 'mol/m^2/s', & caller = caller_str) ! ! ALK ! biotic(n)%ind_alk = otpm_set_prog_tracer('alk' // suffix, & package_name, & longname = 'Alkalinity' // trim(long_suffix), & units = 'mol/m^3', flux_units = 'mol/m^2/s', & caller = caller_str) enddo !} n ! !----------------------------------------------------------------------- ! Process the namelists !----------------------------------------------------------------------- ! ! ! Add the package name to the list of good namelists, to be used ! later for a consistency check ! if (fm_new_value('/ocean_mod/GOOD/good_namelists', package_name, append = .true.) .le. 0) then !{ call mpp_error(FATAL, trim(error_header) // & ' Could not add ' // trim(package_name) // ' to "good_namelists" list') endif !} ! !----------------------------------------------------------------------- ! Set up the *global* namelist !----------------------------------------------------------------------- ! call fm_util_start_namelist(package_name, '*global*', caller = caller_str, no_overwrite = .true., & check = .true.) call fm_util_set_value('po4_star_file', 'INPUT/po4_star_ocmip2.nc') call fm_util_set_value('po4_star_name', 'po4_star') call fm_util_set_value('htotal_scale_lo_in', 0.01 ) call fm_util_set_value('htotal_scale_hi_in', 100.0) call fm_util_set_value('htotal_in', 1.0e-08) call fm_util_end_namelist(package_name, '*global*', caller = caller_str, check = .true.) ! !----------------------------------------------------------------------- ! Set up the instance namelists !----------------------------------------------------------------------- ! t_mask(:) = .true. do n = 1, instances !{ ! ! create the instance namelist ! call fm_util_start_namelist(package_name, biotic(n)%name, caller = caller_str, no_overwrite = .true., & check = .true.) call fm_util_set_value('compensation_depth', 75.0) ! m call fm_util_set_value('add_phosphate', .false.) call fm_util_set_value('martin_coeff', 0.9) call fm_util_set_value('martin_coeff_alt', 0.0) call fm_util_set_value('ca_remin_depth', 3500.0) call fm_util_set_value('soft_tissue_pump', .false.) call fm_util_set_value('stp_temperature', 10.0) call fm_util_set_value('stp_salinity', 34.7) call fm_util_set_value('stp_alkalinity', 2370.0 * 1024.5 * 1.0e-06) ! alkalinity is in ueq/kg, converted to eq/m^3 call fm_util_set_value('sio2_const', 7.7e-03) ! mol/m^3 call fm_util_set_value('local_restart_file', default_local_restart_file) call fm_util_set_value('n_2_p', 16.0) call fm_util_set_value('c_2_p', 117.0) call fm_util_set_value('o_2_p', 170.0) call fm_util_set_value('o2_min', 4.0 * rho_avg * 1.0e-06) call fm_util_set_value('bio_tau', 30.0 * sperd) call fm_util_set_value('sigma', 0.67) call fm_util_set_value('kappa', 1.0 / 0.5 / spery) call fm_util_set_value('caco3_2_c', 0.07) !Old Wanninkhof numbers ! call fm_util_set_value('sc_co2_0', 2073.1) ! call fm_util_set_value('sc_co2_1', -125.62) ! call fm_util_set_value('sc_co2_2', 3.6276) ! call fm_util_set_value('sc_co2_3', -0.043219) ! call fm_util_set_value('sc_o2_0', 1638.0) ! call fm_util_set_value('sc_o2_1', -81.83) ! call fm_util_set_value('sc_o2_2', 1.483) ! call fm_util_set_value('sc_o2_3', -0.008004) !New Wanninkhof numbers call fm_util_set_value('sc_co2_0', 2068.9) call fm_util_set_value('sc_co2_1', -118.63) call fm_util_set_value('sc_co2_2', 2.9311) call fm_util_set_value('sc_co2_3', -0.027) call fm_util_set_value('sc_o2_0', 1929.7) call fm_util_set_value('sc_o2_1', -117.46) call fm_util_set_value('sc_o2_2', 3.116) call fm_util_set_value('sc_o2_3', -0.0306) call fm_util_end_namelist(package_name, biotic(n)%name, check = .true., caller = caller_str) ! ! create some sub-namelists ! call fm_util_start_namelist(trim(package_name), trim(biotic(n)%name) // '+norm_remin', & caller = caller_str, no_overwrite = .true., & check = .true.) call fm_util_set_value('factor', 0.0) call fm_util_set_value('coastal_only', .false.) call fm_util_set_value('t_mask', t_mask, size(t_mask)) call fm_util_set_value('wlon', 0.0, index = 0) call fm_util_set_value('elon', 0.0, index = 0) call fm_util_set_value('slat', 0.0, index = 0) call fm_util_set_value('nlat', 0.0, index = 0) call fm_util_end_namelist(trim(package_name), trim(biotic(n)%name) // '+norm_remin', caller = caller_str) call fm_util_start_namelist(trim(package_name), trim(biotic(n)%name) // '+no_caco3', & caller = caller_str, no_overwrite = .true., & check = .true.) call fm_util_set_value('factor', 0.0) call fm_util_set_value('coastal_only', .false.) call fm_util_set_value('t_mask', t_mask, size(t_mask)) call fm_util_set_value('wlon', 0.0, index = 0) call fm_util_set_value('elon', 0.0, index = 0) call fm_util_set_value('slat', 0.0, index = 0) call fm_util_set_value('nlat', 0.0, index = 0) call fm_util_end_namelist(trim(package_name), trim(biotic(n)%name) // '+no_caco3', caller = caller_str) call fm_util_start_namelist(trim(package_name), trim(biotic(n)%name) // '+nut_depl', & caller = caller_str, no_overwrite = .true., & check = .true.) call fm_util_set_value('factor', 0.0) call fm_util_set_value('coastal_only', .false.) call fm_util_set_value('t_mask', t_mask, size(t_mask)) call fm_util_set_value('wlon', 0.0, index = 0) call fm_util_set_value('elon', 0.0, index = 0) call fm_util_set_value('slat', 0.0, index = 0) call fm_util_set_value('nlat', 0.0, index = 0) call fm_util_end_namelist(trim(package_name), trim(biotic(n)%name) // '+nut_depl', caller = caller_str) call fm_util_start_namelist(trim(package_name), trim(biotic(n)%name) // '+r_bio_tau_a', & caller = caller_str, no_overwrite = .true., & check = .true.) call fm_util_set_value('factor', 0.0) call fm_util_set_value('coastal_only', .false.) call fm_util_set_value('t_mask', t_mask, size(t_mask)) call fm_util_set_value('wlon', 0.0, index = 0) call fm_util_set_value('elon', 0.0, index = 0) call fm_util_set_value('slat', 0.0, index = 0) call fm_util_set_value('nlat', 0.0, index = 0) call fm_util_end_namelist(trim(package_name), trim(biotic(n)%name) // '+r_bio_tau_a', caller = caller_str) enddo !} n ! ! Check for any errors in the number of fields in the namelists for this package ! good_list => fm_util_get_string_array('/ocean_mod/GOOD/namelists/' // trim(package_name) // '/good_values', & caller = trim(mod_name) // '(' // trim(sub_name) // ')') if (associated(good_list)) then !{ call fm_util_check_for_bad_fields('/ocean_mod/namelists/' // trim(package_name), good_list, & caller = trim(mod_name) // '(' // trim(sub_name) // ')') deallocate(good_list) else !}{ call mpp_error(FATAL,trim(error_header) // ' Empty "' // trim(package_name) // '" list') endif !} return end subroutine ocmip2_biotic_init !} ! NAME="ocmip2_biotic_init" !####################################################################### ! ! ! ! Initialize surface fields for flux calculations ! ! Note: this subroutine should be merged into ocmip2_biotic_start ! subroutine ocmip2_biotic_init_sfc(isc, iec, jsc, jec, nk, isd, ied, jsd, jed, & isc_bnd, iec_bnd, jsc_bnd, jec_bnd, & Ocean_fields, T_prog, rho, taum1, model_time, grid_tmask) !{ ! !----------------------------------------------------------------------- ! modules (have to come first) !----------------------------------------------------------------------- ! implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isc integer, intent(in) :: iec integer, intent(in) :: jsc integer, intent(in) :: jec integer, intent(in) :: nk integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed integer, intent(in) :: isc_bnd integer, intent(in) :: iec_bnd integer, intent(in) :: jsc_bnd integer, intent(in) :: jec_bnd type(coupler_2d_bc_type), intent(inout) :: Ocean_fields type(ocean_prog_tracer_type), dimension(:), intent(in) :: T_prog real, dimension(isd:,jsd:,:,:), intent(in) :: rho integer, intent(in) :: taum1 type(time_type), intent(in) :: model_time real, dimension(isd:,jsd:,:), intent(in) :: grid_tmask ! ! local parameters ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_init_sfc' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: i integer :: j integer :: n integer :: ind real :: epsln=1.0e-30 do n = 1, instances !{ ! ! CO2 flux ! ind = biotic(n)%ind_co2_flux if (.not. field_exist('INPUT/'//trim(Ocean_fields%bc(ind)%ocean_restart_file), & Ocean_fields%bc(ind)%field(ind_alpha)%name)) then !{ call ocmip2_co2calc(isd, ied, jsd, jed, isc, iec, jsc, jec, & grid_tmask(isd:ied,jsd:jed,1), & t_prog(indtemp)%field(isd:ied,jsd:jed,1,taum1), & t_prog(indsal)%field(isd:ied,jsd:jed,1,taum1), & t_prog(biotic(n)%ind_dic)%field(isd:ied,jsd:jed,1,taum1), & t_prog(biotic(n)%ind_alk)%field(isd:ied,jsd:jed,1,taum1), & t_prog(biotic(n)%ind_po4)%field(isd:ied,jsd:jed,1,taum1), & biotic(n)%sio2, & htotal_scale_lo, htotal_scale_hi, biotic(n)%htotal, & co2star = biotic(n)%csurf, alpha = biotic(n)%alpha, & pco2surf = biotic(n)%pco2surf) ! !--------------------------------------------------------------------- ! Compute the Schmidt number of CO2 in seawater using the ! formulation presented by Wanninkhof (1992, J. Geophys. Res., 97, ! 7373-7382). !--------------------------------------------------------------------- ! do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%sc_co2(i,j) = & biotic(n)%sc_co2_0 + t_prog(indtemp)%field(i,j,1,taum1) * & (biotic(n)%sc_co2_1 + t_prog(indtemp)%field(i,j,1,taum1) * & (biotic(n)%sc_co2_2 + t_prog(indtemp)%field(i,j,1,taum1) * & biotic(n)%sc_co2_3)) * grid_tmask(i,j,1) sc_no_term(i,j) = sqrt(660.0 / (biotic(n)%sc_co2(i,j) + epsln)) enddo !} i enddo !} j Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & biotic(n)%alpha(isc:iec,jsc:jec) * rho(isc:iec,jsc:jec,1,taum1) * sc_no_term(isc:iec,jsc:jec) Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & biotic(n)%csurf(isc:iec,jsc:jec) * rho(isc:iec,jsc:jec,1,taum1) * sc_no_term(isc:iec,jsc:jec) endif !} ! ! O2 flux ! ind = biotic(n)%ind_o2_flux if (.not. field_exist('INPUT/'//trim(Ocean_fields%bc(ind)%ocean_restart_file), & Ocean_fields%bc(ind)%field(ind_alpha)%name)) then !{ ! !--------------------------------------------------------------------- ! Compute the oxygen saturation concentration at 1 atm total ! pressure in mol/m^3 given the temperature (t, in deg C) and ! the salinity (s, in permil) ! ! From Garcia and Gordon (1992), Limnology and Oceonography. ! The formula used is from page 1310, eq (8). ! ! *** Note: the "a3*ts^2" term (in the paper) is incorrect. *** ! *** It shouldn't be there. *** ! ! o2_saturation is defined between T(freezing) <= T <= 40 deg C and ! 0 permil <= S <= 42 permil ! ! check value: T = 10 deg C, S = 35 permil, ! o2_saturation = 0.282015 mol/m^3 !--------------------------------------------------------------------- ! do j = jsc, jec !{ do i = isc, iec !{ tt(i) = 298.15 - t_prog(indtemp)%field(i,j,1,taum1) tk(i) = 273.15 + t_prog(indtemp)%field(i,j,1,taum1) ts(i) = log(tt(i) / tk(i)) ts2(i) = ts(i) * ts(i) ts3(i) = ts2(i) * ts(i) ts4(i) = ts3(i) * ts(i) ts5(i) = ts4(i) * ts(i) o2_saturation(i,j) = & exp(a_0 + a_1*ts(i) + a_2*ts2(i) + & a_3*ts3(i) + a_4*ts4(i) + a_5*ts5(i) + & t_prog(indsal)%field(i,j,1,taum1) * & (b_0 + b_1*ts(i) + b_2*ts2(i) + b_3*ts3(i) + & c_0*t_prog(indsal)%field(i,j,1,taum1))) enddo !} i enddo !} j ! ! convert from ml/l to mol/m^3 ! do j = jsc, jec !{ do i = isc, iec !{ o2_saturation(i,j) = o2_saturation(i,j) * (1000.0/22391.6) enddo !} i enddo !} j ! !--------------------------------------------------------------------- ! Compute the Schmidt number of O2 in seawater using the ! formulation proposed by Keeling et al. (1998, Global Biogeochem. ! Cycles, 12, 141-163). !--------------------------------------------------------------------- ! do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%sc_o2(i,j) = & biotic(n)%sc_o2_0 + t_prog(indtemp)%field(i,j,1,taum1) * & (biotic(n)%sc_o2_1 + t_prog(indtemp)%field(i,j,1,taum1) * & (biotic(n)%sc_o2_2 + t_prog(indtemp)%field(i,j,1,taum1) * & biotic(n)%sc_o2_3)) * grid_tmask(i,j,1) sc_no_term(i,j) = sqrt(660.0 / (biotic(n)%sc_o2(i,j) + epsln)) enddo !} i enddo !} j Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & o2_saturation(isc:iec,jsc:jec) * sc_no_term(isc:iec,jsc:jec) Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & t_prog(biotic(n)%ind_o2)%field(isc:iec,jsc:jec,1,taum1) * rho(isc:iec,jsc:jec,1,taum1) * sc_no_term(isc:iec,jsc:jec) endif !} enddo !} n return end subroutine ocmip2_biotic_init_sfc !} ! NAME="ocmip2_biotic_init_sfc" !####################################################################### ! ! ! ! Sum surface fields for flux calculations ! subroutine ocmip2_biotic_sum_sfc(isc, iec, jsc, jec, nk, isd, ied, jsd, jed, & isc_bnd, iec_bnd, jsc_bnd, jec_bnd, & Ocean_fields, T_prog, rho, taum1, model_time, grid_tmask) !{ ! !----------------------------------------------------------------------- ! modules (have to come first) !----------------------------------------------------------------------- ! implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isc integer, intent(in) :: iec integer, intent(in) :: jsc integer, intent(in) :: jec integer, intent(in) :: nk integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed integer, intent(in) :: isc_bnd integer, intent(in) :: iec_bnd integer, intent(in) :: jsc_bnd integer, intent(in) :: jec_bnd type(coupler_2d_bc_type), intent(inout) :: Ocean_fields type(ocean_prog_tracer_type), intent(in), dimension(:) :: T_prog real, dimension(isd:,jsd:,:,:), intent(in) :: rho integer, intent(in) :: taum1 type(time_type), intent(in) :: model_time real, dimension(isd:,jsd:,:), intent(in) :: grid_tmask ! ! local parameters ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_sum_sfc' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: i integer :: j integer :: n integer :: ind real :: epsln=1.0e-30 do n = 1, instances !{ ind = biotic(n)%ind_co2_flux call ocmip2_co2calc(isd, ied, jsd, jed, isc, iec, jsc, jec, & grid_tmask(isd:ied,jsd:jed,1), & t_prog(indtemp)%field(isd:ied,jsd:jed,1,taum1), & t_prog(indsal)%field(isd:ied,jsd:jed,1,taum1), & t_prog(biotic(n)%ind_dic)%field(isd:ied,jsd:jed,1,taum1), & t_prog(biotic(n)%ind_alk)%field(isd:ied,jsd:jed,1,taum1), & t_prog(biotic(n)%ind_po4)%field(isd:ied,jsd:jed,1,taum1), & biotic(n)%sio2, & htotal_scale_lo, htotal_scale_hi, biotic(n)%htotal, & co2star = biotic(n)%csurf, alpha = biotic(n)%alpha, & pco2surf = biotic(n)%pco2surf) ! !--------------------------------------------------------------------- ! Compute the Schmidt number of CO2 in seawater using the ! formulation presented by Wanninkhof (1992, J. Geophys. Res., 97, ! 7373-7382). !--------------------------------------------------------------------- ! do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%sc_co2(i,j) = & biotic(n)%sc_co2_0 + t_prog(indtemp)%field(i,j,1,taum1) * & (biotic(n)%sc_co2_1 + t_prog(indtemp)%field(i,j,1,taum1) * & (biotic(n)%sc_co2_2 + t_prog(indtemp)%field(i,j,1,taum1) * & biotic(n)%sc_co2_3)) * grid_tmask(i,j,1) sc_no_term(i,j) = sqrt(660.0 / (biotic(n)%sc_co2(i,j) + epsln)) enddo !} i enddo !} j Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) + & biotic(n)%alpha(isc:iec,jsc:jec) * rho(isc:iec,jsc:jec,1,taum1) * sc_no_term(isc:iec,jsc:jec) Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) + & biotic(n)%csurf(isc:iec,jsc:jec) * rho(isc:iec,jsc:jec,1,taum1) * sc_no_term(isc:iec,jsc:jec) ind = biotic(n)%ind_o2_flux ! !--------------------------------------------------------------------- ! Compute the oxygen saturation concentration at 1 atm total ! pressure in mol/m^3 given the temperature (t, in deg C) and ! the salinity (s, in permil) ! ! From Garcia and Gordon (1992), Limnology and Oceonography. ! The formula used is from page 1310, eq (8). ! ! *** Note: the "a3*ts^2" term (in the paper) is incorrect. *** ! *** It shouldn't be there. *** ! ! o2_saturation is defined between T(freezing) <= T <= 40 deg C and ! 0 permil <= S <= 42 permil ! ! check value: T = 10 deg C, S = 35 permil, ! o2_saturation = 0.282015 mol/m^3 !--------------------------------------------------------------------- ! do j = jsc, jec !{ do i = isc, iec !{ tt(i) = 298.15 - t_prog(indtemp)%field(i,j,1,taum1) tk(i) = 273.15 + t_prog(indtemp)%field(i,j,1,taum1) ts(i) = log(tt(i) / tk(i)) ts2(i) = ts(i) * ts(i) ts3(i) = ts2(i) * ts(i) ts4(i) = ts3(i) * ts(i) ts5(i) = ts4(i) * ts(i) o2_saturation(i,j) = & exp(a_0 + a_1*ts(i) + a_2*ts2(i) + & a_3*ts3(i) + a_4*ts4(i) + a_5*ts5(i) + & t_prog(indsal)%field(i,j,1,taum1) * & (b_0 + b_1*ts(i) + b_2*ts2(i) + b_3*ts3(i) + & c_0*t_prog(indsal)%field(i,j,1,taum1))) enddo !} i enddo !} j ! ! convert from ml/l to mol/m^3 ! do j = jsc, jec !{ do i = isc, iec !{ o2_saturation(i,j) = o2_saturation(i,j) * (1000.0/22391.6) enddo !} i enddo !} j ! !--------------------------------------------------------------------- ! Compute the Schmidt number of O2 in seawater using the ! formulation proposed by Keeling et al. (1998, Global Biogeochem. ! Cycles, 12, 141-163). !--------------------------------------------------------------------- ! do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%sc_o2(i,j) = & biotic(n)%sc_o2_0 + t_prog(indtemp)%field(i,j,1,taum1) * & (biotic(n)%sc_o2_1 + t_prog(indtemp)%field(i,j,1,taum1) * & (biotic(n)%sc_o2_2 + t_prog(indtemp)%field(i,j,1,taum1) * & biotic(n)%sc_o2_3)) * grid_tmask(i,j,1) sc_no_term(i,j) = sqrt(660.0 / (biotic(n)%sc_o2(i,j) + epsln)) enddo !} i enddo !} j Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) + & o2_saturation(isc:iec,jsc:jec) * sc_no_term(isc:iec,jsc:jec) Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) + & t_prog(biotic(n)%ind_o2)%field(isc:iec,jsc:jec,1,taum1) * rho(isc:iec,jsc:jec,1,taum1) * sc_no_term(isc:iec,jsc:jec) enddo !} n return end subroutine ocmip2_biotic_sum_sfc !} ! NAME="ocmip2_biotic_sum_sfc" !####################################################################### ! ! ! ! Sum surface fields for flux calculations ! subroutine ocmip2_biotic_zero_sfc(Ocean_fields) !{ implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! type(coupler_2d_bc_type), intent(inout) :: Ocean_fields ! ! local parameters ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_zero_sfc' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: n integer :: ind do n = 1, instances !{ ind = biotic(n)%ind_co2_flux Ocean_fields%bc(ind)%field(ind_alpha)%values = 0.0 Ocean_fields%bc(ind)%field(ind_csurf)%values = 0.0 ind = biotic(n)%ind_o2_flux Ocean_fields%bc(ind)%field(ind_alpha)%values = 0.0 Ocean_fields%bc(ind)%field(ind_csurf)%values = 0.0 enddo !} n return end subroutine ocmip2_biotic_zero_sfc !} ! NAME="ocmip2_biotic_zero_sfc" !####################################################################### ! ! ! ! Sum surface fields for flux calculations ! subroutine ocmip2_biotic_avg_sfc(isc, iec, jsc, jec, nk, isd, ied, jsd, jed, & isc_bnd, iec_bnd, jsc_bnd, jec_bnd, Ocean_fields, Ocean_avg_kount, grid_tmask) !{ implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isc integer, intent(in) :: iec integer, intent(in) :: jsc integer, intent(in) :: jec integer, intent(in) :: nk integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed integer, intent(in) :: isc_bnd integer, intent(in) :: iec_bnd integer, intent(in) :: jsc_bnd integer, intent(in) :: jec_bnd type(coupler_2d_bc_type), intent(inout) :: Ocean_fields integer :: Ocean_avg_kount real, dimension(isd:,jsd:,:), intent(in) :: grid_tmask ! ! local parameters ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_avg_sfc' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! arguments !----------------------------------------------------------------------- ! ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: n integer :: ind real :: divid divid = 1./float(Ocean_avg_kount) do n = 1, instances !{ ind = biotic(n)%ind_co2_flux where (grid_tmask(isc:iec,jsc:jec,1) == 1.0) !{ Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) * divid Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) * divid endwhere !} ind = biotic(n)%ind_o2_flux where (grid_tmask(isc:iec,jsc:jec,1) == 1.0) !{ Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & Ocean_fields%bc(ind)%field(ind_alpha)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) * divid Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = & Ocean_fields%bc(ind)%field(ind_csurf)%values(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) * divid endwhere !} enddo !} n return end subroutine ocmip2_biotic_avg_sfc !} ! NAME="ocmip2_biotic_avg_sfc" !####################################################################### ! ! ! ! Initialize surface fields for flux calculations ! subroutine ocmip2_biotic_sfc_end !{ implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! ! ! local parameters ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_sfc_end' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! return end subroutine ocmip2_biotic_sfc_end !} ! NAME="ocmip2_biotic_sfc_end" !####################################################################### ! ! ! ! compute the source terms for the BIOTICs, including boundary ! conditions (not done in setvbc, to minimize number ! of hooks required in MOM base code) ! ! subroutine ocmip2_biotic_source(isc, iec, jsc, jec, nk, isd, ied, jsd, jed, & T_prog, taum1, model_time, grid_zw, grid_ht, grid_tmask, rho_dzt) !{ ! !----------------------------------------------------------------------- ! modules (have to come first) !----------------------------------------------------------------------- ! implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isc integer, intent(in) :: iec integer, intent(in) :: jsc integer, intent(in) :: jec integer, intent(in) :: nk integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed type(ocean_prog_tracer_type), intent(inout), dimension(:) :: T_prog integer, intent(in) :: taum1 type(time_type), intent(in) :: model_time real, dimension(nk), intent(in) :: grid_zw real, dimension(isd:,jsd:), intent(in) :: grid_ht real, dimension(isd:,jsd:,:), intent(in) :: grid_tmask real, dimension(isd:,jsd:,:,:), intent(in) :: rho_dzt ! !----------------------------------------------------------------------- ! local parameters !----------------------------------------------------------------------- ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_source' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! integer :: i integer :: j integer :: k integer :: n integer :: ind_po4 integer :: ind_dop integer :: ind_dic integer :: ind_alk integer :: ind_o2 integer :: km_c logical :: used integer :: day integer :: month integer :: year integer :: hour integer :: minute integer :: second real, dimension(isc:iec,jsc:jec,nk) :: grid_tmask_comp integer :: stdoutunit stdoutunit=stdout() ! ! ===================================================================== ! begin executable code ! ===================================================================== ! !if (no_source) then !{ !write (stdoutunit,*) trim(warn_header), 'Skipping calculations' !return !else !}{ !write (stdoutunit,*) ', trim(note_header), 'Doing calculations' !endif !} ! ! get the model month ! call get_date(model_time, year, month, day, & hour, minute, second) ! !----------------------------------------------------------------------- ! calculate the source terms for BIOTICs !----------------------------------------------------------------------- ! ! !--------------------------------------------------------------------- ! calculate interpolated PO4_star !--------------------------------------------------------------------- ! call time_interp_external(po4_star_id, model_time, po4_star_t) ! ! Loop over multiple instances ! do n = 1, instances !{ ind_po4 = biotic(n)%ind_po4 ind_dic = biotic(n)%ind_dic ind_dop = biotic(n)%ind_dop ind_o2 = biotic(n)%ind_o2 ind_alk = biotic(n)%ind_alk km_c = biotic(n)%km_c ! !----------------------------------------------------------------------- ! Production !----------------------------------------------------------------------- ! ! ! compute PO4 restoring term and correct for partial ! production in the bottom box ! do k = 1, km_c !{ do j = jsc, jec !{ do i = isc, iec !{ if (t_prog(ind_po4)%field(i,j,k,taum1) .gt. & po4_star_t(i,j,k) * & biotic(n)%nut_depl%mask(i,j,month)) then !{ biotic(n)%jprod(i,j,k) = & (t_prog(ind_po4)%field(i,j,k,taum1) - & po4_star_t(i,j,k) * & biotic(n)%nut_depl%mask(i,j,month)) * & biotic(n)%r_bio_tau_a%mask(i,j,month) * grid_tmask(i,j,k) else !}{ biotic(n)%jprod(i,j,k) = 0.0 endif !} enddo !} i enddo !} j enddo !} k do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jprod(i,j,km_c) = biotic(n)%jprod(i,j,km_c) * & (min(grid_ht(i,j),biotic(n)%compensation_depth) - grid_zw(km_c-1)) / & rho_dzt(i,j,km_c,taum1) enddo !} i enddo !} j ! !----------------------------------------------------------------------- ! Normal production (used to maintain PO4 concentrations assuming no ! other changes) !----------------------------------------------------------------------- ! if (biotic(n)%add_phosphate) then !{ ! ! compute PO4 restoring term and correct for partial ! production in the bottom box ! do k = 1, km_c !{ do j = jsc, jec !{ do i = isc, iec !{ if (t_prog(ind_po4)%field(i,j,k,taum1) .gt. & po4_star_t(i,j,k)) then !{ biotic(n)%jprod_norm(i,j,k) = & (t_prog(ind_po4)%field(i,j,k,taum1) - & po4_star_t(i,j,k)) * & biotic(n)%r_bio_tau * grid_tmask(i,j,k) else !}{ biotic(n)%jprod_norm(i,j,k) = 0.0 endif !} enddo !} i enddo !} j enddo !} k do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jprod_norm(i,j,km_c) = & biotic(n)%jprod_norm(i,j,km_c) * & (min(grid_ht(i,j),biotic(n)%compensation_depth) - grid_zw(km_c-1)) / & rho_dzt(i,j,km_c,taum1) enddo !} i enddo !} j endif !} ! !----------------------------------------------------------------------- ! Particle flux !----------------------------------------------------------------------- ! do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%flux_poc(i,j) = (1.0 - biotic(n)%sigma) * & biotic(n)%jprod(i,j,1) * rho_dzt(i,j,1,taum1) enddo !} i enddo !} j do k = 2, km_c !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%flux_poc(i,j) = biotic(n)%flux_poc(i,j) + & (1.0 - biotic(n)%sigma) * biotic(n)%jprod(i,j,k) * & rho_dzt(i,j,k,taum1) enddo !} i enddo !} j enddo !} k ! ! calculate the flux at the base of each layer below the ! compensation depth ! do k = km_c, nk !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%fc(i,j,k) = biotic(n)%flux_poc(i,j) * & grid_tmask(i,j,k) * & (biotic(n)%norm_remin%mask(i,j,month) * & biotic(n)%zforg(k) + & (1.0 - biotic(n)%norm_remin%mask(i,j,month)) * & biotic(n)%zforg_alt(k)) enddo !} i enddo !} j enddo !} k ! !----------------------------------------------------------------------- ! Calcium Carbonate !----------------------------------------------------------------------- ! ! ! calculate the formation of CaCO3 above the compensation ! depth ! do k = 1, km_c !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jca(i,j,k) = & -biotic(n)%caco3_2_c * biotic(n)%c_2_p * & (1.0 - biotic(n)%sigma) * & biotic(n)%jprod(i,j,k) * grid_tmask(i,j,k) * & biotic(n)%no_caco3%mask(i,j,month) enddo !} i enddo !} j enddo !} k ! ! calculate the flux of CaCO3 at compensation depth and at ! bottom of each level below ! do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%flux_caco3(i,j) = & biotic(n)%caco3_2_c * biotic(n)%c_2_p * & biotic(n)%flux_poc(i,j) * biotic(n)%no_caco3%mask(i,j,month) enddo !} i enddo !} j do k=km_c,nk !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%fca(i,j,k) = biotic(n)%flux_caco3(i,j) * & biotic(n)%zfca(k) * grid_tmask(i,j,k) enddo !} i enddo !} j enddo !} k ! ! calculate the dissolution of CaCO3 below the compensation ! depth ! do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jca(i,j,km_c) = biotic(n)%jca(i,j,km_c) + & (biotic(n)%flux_caco3(i,j) * grid_tmask(i,j,km_c) - & biotic(n)%fca(i,j,km_c) * grid_tmask(i,j,km_c+1)) / & rho_dzt(i,j,km_c,taum1) enddo !} i enddo !} j do k = km_c + 1, nk - 1 !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jca(i,j,k) = & (biotic(n)%fca(i,j,k-1) * grid_tmask(i,j,k) - & biotic(n)%fca(i,j,k) * grid_tmask(i,j,k+1)) / rho_dzt(i,j,k,taum1) enddo !} i enddo !} j enddo !} k do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jca(i,j,nk) = biotic(n)%fca(i,j,nk-1) * & grid_tmask(i,j,nk) / rho_dzt(i,j,nk,taum1) enddo !} i enddo !} j ! !----------------------------------------------------------------------- ! PO4 !----------------------------------------------------------------------- ! if (biotic(n)%add_phosphate) then !{ do k = 1, km_c !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jpo4_add(i,j,k) = & biotic(n)%jprod(i,j,k) - biotic(n)%jprod_norm(i,j,k) enddo !} i enddo !} j enddo !} k endif !} do k = 1, km_c !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jpo4(i,j,k) = -biotic(n)%jprod(i,j,k) + & biotic(n)%kappa * t_prog(ind_dop)%field(i,j,k,taum1) enddo !} i enddo !} j enddo !} k do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jpo4(i,j,km_c) = biotic(n)%jpo4(i,j,km_c) + & (biotic(n)%flux_poc(i,j) * grid_tmask(i,j,km_c) - & biotic(n)%fc(i,j,km_c) * grid_tmask(i,j,km_c+1)) / & rho_dzt(i,j,km_c,taum1) enddo !} i enddo !} j do k = km_c + 1, nk - 1 !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jpo4(i,j,k) = & (biotic(n)%fc(i,j,k-1) * grid_tmask(i,j,k) - & biotic(n)%fc(i,j,k) * grid_tmask(i,j,k+1)) / & rho_dzt(i,j,k,taum1) + & biotic(n)%kappa * t_prog(ind_dop)%field(i,j,k,taum1) enddo !} i enddo !} j enddo !} k do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jpo4(i,j,nk) = & biotic(n)%fc(i,j,nk-1) * grid_tmask(i,j,nk) / & rho_dzt(i,j,nk,taum1) + & biotic(n)%kappa * t_prog(ind_dop)%field(i,j,nk,taum1) enddo !} i enddo !} j do k = 1, nk !{ do j = jsc, jec !{ do i = isc, iec !{ t_prog(ind_po4)%th_tendency(i,j,k) = & t_prog(ind_po4)%th_tendency(i,j,k) + & biotic(n)%jpo4(i,j,k) * rho_dzt(i,j,k,taum1) enddo !} i enddo !} j enddo !} k if (biotic(n)%add_phosphate) then !{ do k = 1, km_c !{ do j = jsc, jec !{ do i = isc, iec !{ t_prog(ind_po4)%th_tendency(i,j,k) = & t_prog(ind_po4)%th_tendency(i,j,k) + & biotic(n)%jpo4_add(i,j,k) * rho_dzt(i,j,k,taum1) enddo !} i enddo !} j enddo !} k endif !} ! !----------------------------------------------------------------------- ! DOP !----------------------------------------------------------------------- ! do k = 1, km_c !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jdop(i,j,k) = & biotic(n)%sigma * biotic(n)%jprod(i,j,k) - & biotic(n)%kappa * t_prog(ind_dop)%field(i,j,k,taum1) enddo !} i enddo !} j enddo !} k do k = km_c + 1, nk !{ do j = jsc, jec !{ do i = isc, iec !{ biotic(n)%jdop(i,j,k) = -biotic(n)%kappa * & t_prog(ind_dop)%field(i,j,k,taum1) enddo !} i enddo !} j enddo !} k do k = 1, nk !{ do j = jsc, jec !{ do i = isc, iec !{ t_prog(ind_dop)%th_tendency(i,j,k) = & t_prog(ind_dop)%th_tendency(i,j,k) + & biotic(n)%jdop(i,j,k) * rho_dzt(i,j,k,taum1) enddo !} i enddo !} j enddo !} k ! !----------------------------------------------------------------------- ! O2 !----------------------------------------------------------------------- ! do k = 1, nk !{ do j = jsc, jec !{ do i = isc, iec !{ if (t_prog(ind_o2)%field(i,j,k,taum1) .gt. & biotic(n)%o2_min) then !{ biotic(n)%jo2(i,j,k) = -biotic(n)%o_2_p * & biotic(n)%jpo4(i,j,k) else !}{ biotic(n)%jo2(i,j,k) = 0.0 endif !} enddo !} i enddo !} j enddo !} k do k = 1, nk !{ do j = jsc, jec !{ do i = isc, iec !{ t_prog(ind_o2)%th_tendency(i,j,k) = & t_prog(ind_o2)%th_tendency(i,j,k) + & biotic(n)%jo2(i,j,k) * rho_dzt(i,j,k,taum1) enddo !} i enddo !} j enddo !} k ! !----------------------------------------------------------------------- ! DIC !----------------------------------------------------------------------- ! do k = 1, nk !{ do j = jsc, jec !{ do i = isc, iec !{ t_prog(ind_dic)%th_tendency(i,j,k) = & t_prog(ind_dic)%th_tendency(i,j,k) + & (biotic(n)%c_2_p * biotic(n)%jpo4(i,j,k) + & biotic(n)%jca(i,j,k)) * rho_dzt(i,j,k,taum1) enddo !} i enddo !} j enddo !} k ! !----------------------------------------------------------------------- ! ALK !----------------------------------------------------------------------- ! do k = 1, nk !{ do j = jsc, jec !{ do i = isc, iec !{ t_prog(ind_alk)%th_tendency(i,j,k) = & t_prog(ind_alk)%th_tendency(i,j,k) + & (-biotic(n)%n_2_p * biotic(n)%jpo4(i,j,k) + & 2.0 * biotic(n)%jca(i,j,k)) * rho_dzt(i,j,k,taum1) enddo !} i enddo !} j enddo !} k enddo !} n ! !----------------------------------------------------------------------- ! Save variables for diagnostics !----------------------------------------------------------------------- ! ! ! set up the grid mask on the computational grid so that we ! will not need to implicitly copy arrays in the following ! subroutine calls ! grid_tmask_comp = grid_tmask(isc:iec,jsc:jec,:) if (id_o2_sat .gt. 0) then used = send_data(id_o2_sat, o2_saturation(:,:), & model_time, rmask = grid_tmask_comp(:,:,1)) endif do n = 1, instances !{ if (biotic(n)%id_sc_co2 .gt. 0) then used = send_data(biotic(n)%id_sc_co2, & biotic(n)%sc_co2(:,:), & model_time, rmask = grid_tmask_comp(:,:,1)) endif if (biotic(n)%id_sc_o2 .gt. 0) then used = send_data(biotic(n)%id_sc_o2, & biotic(n)%sc_o2(:,:), & model_time, rmask = grid_tmask_comp(:,:,1)) endif if (biotic(n)%id_alpha .gt. 0) then used = send_data(biotic(n)%id_alpha, & biotic(n)%alpha(:,:), & model_time, rmask = grid_tmask_comp(:,:,1)) endif if (biotic(n)%id_csurf .gt. 0) then used = send_data(biotic(n)%id_csurf, & biotic(n)%csurf(:,:), & model_time, rmask = grid_tmask_comp(:,:,1)) endif if (biotic(n)%id_pco2surf .gt. 0) then used = send_data(biotic(n)%id_pco2surf, & biotic(n)%pco2surf(:,:), & model_time, rmask = grid_tmask_comp(:,:,1)) endif if (biotic(n)%id_flux_poc .gt. 0) then used = send_data(biotic(n)%id_flux_poc, & biotic(n)%flux_poc(:,:), & model_time, rmask = grid_tmask_comp(:,:,1)) endif if (biotic(n)%id_flux_caco3 .gt. 0) then used = send_data(biotic(n)%id_flux_caco3, & biotic(n)%flux_caco3(:,:), & model_time, rmask = grid_tmask_comp(:,:,1)) endif if (biotic(n)%id_htotal .gt. 0) then used = send_data(biotic(n)%id_htotal, & biotic(n)%htotal(:,:), & model_time, rmask = grid_tmask_comp(:,:,1)) endif if (biotic(n)%id_jprod .gt. 0) then used = send_data(biotic(n)%id_jprod, & biotic(n)%jprod(:,:,:), & model_time, rmask = grid_tmask_comp(:,:,:)) endif if (biotic(n)%id_jca .gt. 0) then used = send_data(biotic(n)%id_jca, & biotic(n)%jca(:,:,:), & model_time, rmask = grid_tmask_comp(:,:,:)) endif if (biotic(n)%id_jpo4 .gt. 0) then used = send_data(biotic(n)%id_jpo4, & biotic(n)%jpo4(:,:,:), & model_time, rmask = grid_tmask_comp(:,:,:)) endif if (biotic(n)%id_jdop .gt. 0) then used = send_data(biotic(n)%id_jdop, & biotic(n)%jdop(:,:,:), & model_time, rmask = grid_tmask_comp(:,:,:)) endif if (biotic(n)%id_jo2 .gt. 0) then used = send_data(biotic(n)%id_jo2, & biotic(n)%jo2(:,:,:), & model_time, rmask = grid_tmask_comp(:,:,:)) endif if (biotic(n)%id_fc .gt. 0) then used = send_data(biotic(n)%id_fc, & biotic(n)%fc(:,:,:), & model_time, rmask = grid_tmask_comp(:,:,:)) endif if (biotic(n)%id_fca .gt. 0) then used = send_data(biotic(n)%id_fca, & biotic(n)%fca(:,:,:), & model_time, rmask = grid_tmask_comp(:,:,:)) endif if (biotic(n)%add_phosphate) then !{ if (biotic(n)%id_jpo4_add .gt. 0) then used = send_data(biotic(n)%id_jpo4_add, & biotic(n)%jpo4_add(:,:,1:km_c), & model_time, rmask = grid_tmask_comp(:,:,1:km_c)) endif endif !} enddo !} n return end subroutine ocmip2_biotic_source !} ! NAME="ocmip2_biotic_source" !####################################################################### ! ! ! ! Initialize variables, read in namelists, calculate constants ! for a given run and allocate diagnostic arrays ! ! subroutine ocmip2_biotic_start(isc, iec, jsc, jec, nk, isd, ied, jsd, jed, & T_prog, taup1, model_time, grid_dat, grid_tmask, grid_kmt, & grid_xt, grid_yt, grid_zt, grid_zw, grid_dzt, grid_name, grid_tracer_axes, & mpp_domain2d, rho_dzt) !{ ! !----------------------------------------------------------------------- ! modules (have to come first) !----------------------------------------------------------------------- ! use time_interp_external_mod, only: init_external_field use diag_manager_mod, only : register_diag_field, diag_axis_init use field_manager_mod, only : fm_get_index implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isc integer, intent(in) :: iec integer, intent(in) :: jsc integer, intent(in) :: jec integer, intent(in) :: nk integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed type(ocean_prog_tracer_type), dimension(:), intent(in) :: T_prog integer, intent(in) :: taup1 type(time_type), intent(in) :: model_time real, dimension(isd:,jsd:), intent(in) :: grid_dat real, dimension(isd:,jsd:,:), intent(in) :: grid_tmask integer, dimension(isd:,jsd:), intent(in) :: grid_kmt real, dimension(isd:,jsd:), intent(in) :: grid_xt real, dimension(isd:,jsd:), intent(in) :: grid_yt real, dimension(:), intent(in) :: grid_zt real, dimension(:), intent(in) :: grid_zw real, dimension(:), intent(in) :: grid_dzt character(len=*), intent(in) :: grid_name integer, dimension(3), intent(in) :: grid_tracer_axes type(domain2d), intent(in) :: mpp_domain2d real, dimension(isd:,jsd:,:,:), intent(in) :: rho_dzt ! !----------------------------------------------------------------------- ! local parameters !----------------------------------------------------------------------- ! character(len=64), parameter :: sub_name = 'ocmip2_biotic_start' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! !----------------------------------------------------------------------- ! local variables !----------------------------------------------------------------------- ! ! !---------------------------------------------------------------------- ! ! Global values to apply the following inhibitions ! and depletions ! ! coastal_only : if true, then only apply the changes in ! coastal boxes ! t_mask_len : parameter giving the number of elements in ! the time mask per year (eg., 12 would ! imply monthly) ! t_mask_array : logical array controlling whether to apply ! the following inhibitions and depletions to ! each time-period (true means set the masks, ! false means use the defaults everywhere) ! num_reg : number of regions ! factor : factor by which to scale the field ! : in the selected regions ! wlon : western longitude of region ! elon : eastern longitude of region ! slat : southern latitude of region ! nlat : northern latitude of region ! mask(imt,jmt) : mask array (0.0 - alternate, 1.0 - normal) ! ! Set up a mask array using wlon,elon,nlat,slat ! (any box with its lon,lat inside the box bounded by ! wlon,elon,nlat,slat value in mask set to factor). ! !---------------------------------------------------------------------- ! integer :: done integer :: i integer :: j integer :: k integer :: l integer :: n character(len=fm_field_name_len+1) :: suffix character(len=fm_field_name_len+3) :: long_suffix character(len=256) :: caller_str integer :: len_w integer :: len_e integer :: len_s integer :: len_n real :: total_alkalinity real :: total_dic real :: total_dop real :: total_o2 real :: total_phosphate character(len=fm_string_len), allocatable :: local_restart_file(:) integer :: ind logical :: fld_exist integer :: id_restart integer :: stdoutunit stdoutunit=stdout() ! ! ===================================================================== ! begin of executable code ! ===================================================================== ! ! !----------------------------------------------------------------------- ! give info !----------------------------------------------------------------------- ! write(stdoutunit,*) write(stdoutunit,*) trim(note_header), & 'Starting ', trim(package_name), ' module' ! ! Determine indices for temperature and salinity ! indtemp = fm_get_index('/ocean_mod/prog_tracers/temp') if (indtemp .le. 0) then !{ call mpp_error(FATAL,trim(error_header) // ' Could not get the temperature index') endif !} indsal = fm_get_index('/ocean_mod/prog_tracers/salt') if (indsal .le. 0) then !{ call mpp_error(FATAL,trim(error_header) // ' Could not get the salinity index') endif !} ! !----------------------------------------------------------------------- ! dynamically allocate the global BIOTIC arrays !----------------------------------------------------------------------- ! call allocate_arrays(isc, iec, jsc, jec, nk, isd, ied, jsd, jed) ! !----------------------------------------------------------------------- ! save the *global* namelist values !----------------------------------------------------------------------- ! caller_str = trim(mod_name) // '(' // trim(sub_name) // ')' call fm_util_start_namelist(package_name, '*global*', caller = caller_str) po4_star_file = fm_util_get_string ('po4_star_file', scalar = .true.) po4_star_name = fm_util_get_string ('po4_star_name', scalar = .true.) htotal_scale_lo_in = fm_util_get_real ('htotal_scale_lo_in', scalar = .true.) htotal_scale_hi_in = fm_util_get_real ('htotal_scale_hi_in', scalar = .true.) htotal_in = fm_util_get_real ('htotal_in', scalar = .true.) call fm_util_end_namelist(package_name, '*global*', caller = caller_str) ! !----------------------------------------------------------------------- ! Open up the PO4 file for restoring !----------------------------------------------------------------------- ! po4_star_id = init_external_field(po4_star_file, & po4_star_name, & domain = mpp_domain2d) if (po4_star_id .eq. 0) then !{ call mpp_error(FATAL, trim(error_header) // & 'Could not open po4_star file: ' // & trim(po4_star_file)) endif !} ! ! set default values for htotal_scale bounds ! htotal_scale_lo(:,:) = htotal_scale_lo_in htotal_scale_hi(:,:) = htotal_scale_hi_in ! !----------------------------------------------------------------------- ! read in the namelists for each instance !----------------------------------------------------------------------- ! do n = 1, instances !{ call fm_util_start_namelist(package_name, biotic(n)%name, caller = caller_str) biotic(n)%compensation_depth = fm_util_get_real ('compensation_depth', scalar = .true.) biotic(n)%add_phosphate = fm_util_get_logical('add_phosphate', scalar = .true.) biotic(n)%martin_coeff = fm_util_get_real ('martin_coeff', scalar = .true.) biotic(n)%martin_coeff_alt = fm_util_get_real ('martin_coeff_alt', scalar = .true.) biotic(n)%ca_remin_depth = fm_util_get_real ('ca_remin_depth', scalar = .true.) biotic(n)%soft_tissue_pump = fm_util_get_logical('soft_tissue_pump', scalar = .true.) biotic(n)%stp_temperature = fm_util_get_real ('stp_temperature', scalar = .true.) biotic(n)%stp_salinity = fm_util_get_real ('stp_salinity', scalar = .true.) biotic(n)%stp_alkalinity = fm_util_get_real ('stp_alkalinity', scalar = .true.) biotic(n)%sio2_const = fm_util_get_real ('sio2_const', scalar = .true.) biotic(n)%local_restart_file = fm_util_get_string ('local_restart_file', scalar = .true.) biotic(n)%n_2_p = fm_util_get_real ('n_2_p', scalar = .true.) biotic(n)%c_2_p = fm_util_get_real ('c_2_p', scalar = .true.) biotic(n)%o_2_p = fm_util_get_real ('o_2_p', scalar = .true.) biotic(n)%o2_min = fm_util_get_real ('o2_min', scalar = .true.) biotic(n)%bio_tau = fm_util_get_real ('bio_tau', scalar = .true.) biotic(n)%sigma = fm_util_get_real ('sigma', scalar = .true.) biotic(n)%kappa = fm_util_get_real ('kappa', scalar = .true.) biotic(n)%caco3_2_c = fm_util_get_real ('caco3_2_c', scalar = .true.) biotic(n)%sc_co2_0 = fm_util_get_real ('sc_co2_0', scalar = .true.) biotic(n)%sc_co2_1 = fm_util_get_real ('sc_co2_1', scalar = .true.) biotic(n)%sc_co2_2 = fm_util_get_real ('sc_co2_2', scalar = .true.) biotic(n)%sc_co2_3 = fm_util_get_real ('sc_co2_3', scalar = .true.) biotic(n)%sc_o2_0 = fm_util_get_real ('sc_o2_0', scalar = .true.) biotic(n)%sc_o2_1 = fm_util_get_real ('sc_o2_1', scalar = .true.) biotic(n)%sc_o2_2 = fm_util_get_real ('sc_o2_2', scalar = .true.) biotic(n)%sc_o2_3 = fm_util_get_real ('sc_o2_3', scalar = .true.) call fm_util_end_namelist(package_name, biotic(n)%name, caller = caller_str) biotic(n)%r_bio_tau = 1.0 / biotic(n)%bio_tau call fm_util_start_namelist(trim(package_name), trim(biotic(n)%name) // '+norm_remin', caller = caller_str) biotic(n)%norm_remin%factor = fm_util_get_real ('factor', scalar = .true.) biotic(n)%norm_remin%coastal_only = fm_util_get_logical ('coastal_only', scalar = .true.) biotic(n)%norm_remin%wlon => fm_util_get_real_array ('wlon') biotic(n)%norm_remin%elon => fm_util_get_real_array ('elon') biotic(n)%norm_remin%slat => fm_util_get_real_array ('slat') biotic(n)%norm_remin%nlat => fm_util_get_real_array ('nlat') biotic(n)%norm_remin%t_mask => fm_util_get_logical_array ('t_mask') call fm_util_end_namelist(trim(package_name), trim(biotic(n)%name) // '+norm_remin', caller = caller_str) call fm_util_start_namelist(trim(package_name), trim(biotic(n)%name) // '+no_caco3', caller = caller_str) biotic(n)%no_caco3%factor = fm_util_get_real ('factor', scalar = .true.) biotic(n)%no_caco3%coastal_only = fm_util_get_logical ('coastal_only', scalar = .true.) biotic(n)%no_caco3%wlon => fm_util_get_real_array ('wlon') biotic(n)%no_caco3%elon => fm_util_get_real_array ('elon') biotic(n)%no_caco3%slat => fm_util_get_real_array ('slat') biotic(n)%no_caco3%nlat => fm_util_get_real_array ('nlat') biotic(n)%no_caco3%t_mask => fm_util_get_logical_array ('t_mask') call fm_util_end_namelist(trim(package_name), trim(biotic(n)%name) // '+no_caco3', caller = caller_str) call fm_util_start_namelist(trim(package_name), trim(biotic(n)%name) // '+nut_depl', caller = caller_str) biotic(n)%nut_depl%factor = fm_util_get_real ('factor', scalar = .true.) biotic(n)%nut_depl%coastal_only = fm_util_get_logical ('coastal_only', scalar = .true.) biotic(n)%nut_depl%wlon => fm_util_get_real_array ('wlon') biotic(n)%nut_depl%elon => fm_util_get_real_array ('elon') biotic(n)%nut_depl%slat => fm_util_get_real_array ('slat') biotic(n)%nut_depl%nlat => fm_util_get_real_array ('nlat') biotic(n)%nut_depl%t_mask => fm_util_get_logical_array ('t_mask') call fm_util_end_namelist(trim(package_name), trim(biotic(n)%name) // '+nut_depl', caller = caller_str) call fm_util_start_namelist(trim(package_name), trim(biotic(n)%name) // '+r_bio_tau_a', caller = caller_str) biotic(n)%r_bio_tau_a%factor = fm_util_get_real ('factor', scalar = .true.) biotic(n)%r_bio_tau_a%coastal_only = fm_util_get_logical ('coastal_only', scalar = .true.) biotic(n)%r_bio_tau_a%wlon => fm_util_get_real_array ('wlon') biotic(n)%r_bio_tau_a%elon => fm_util_get_real_array ('elon') biotic(n)%r_bio_tau_a%slat => fm_util_get_real_array ('slat') biotic(n)%r_bio_tau_a%nlat => fm_util_get_real_array ('nlat') biotic(n)%r_bio_tau_a%t_mask => fm_util_get_logical_array ('t_mask') call fm_util_end_namelist(trim(package_name), trim(biotic(n)%name) // '+r_bio_tau_a', caller = caller_str) enddo !} n ! !----------------------------------------------------------------------- ! calculate the index for the box containing the compensation depth !----------------------------------------------------------------------- ! km_c_max = 0 do n = 1, instances !{ call locate(grid_zw, nk, biotic(n)%compensation_depth, & biotic(n)%km_c, nearest = .true.) if (grid_zw(biotic(n)%km_c) .lt. & biotic(n)%compensation_depth) then !{ biotic(n)%km_c = biotic(n)%km_c + 1 endif !} write (stdoutunit,*) trim(note_header), & 'The compensation depth for instance ', & n, ', ', biotic(n)%compensation_depth, & ' m, occurs in box ', biotic(n)%km_c , & ' between depths ', & grid_zw(biotic(n)%km_c-1), & ' m and ', grid_zw(biotic(n)%km_c), ' m' km_c_max = max(km_c_max, biotic(n)%km_c) enddo !} n ! !----------------------------------------------------------------------- ! read in the norm_remin namelist data !----------------------------------------------------------------------- ! do n = 1, instances !{ ! ! Check some things ! if (associated(biotic(n)%norm_remin%wlon)) then !{ len_w = size(biotic(n)%norm_remin%wlon) else !}{ len_w = 0 endif !} if (associated(biotic(n)%norm_remin%elon)) then !{ len_e = size(biotic(n)%norm_remin%elon) else !}{ len_e = 0 endif !} if (associated(biotic(n)%norm_remin%slat)) then !{ len_s = size(biotic(n)%norm_remin%slat) else !}{ len_s = 0 endif !} if (associated(biotic(n)%norm_remin%nlat)) then !{ len_n = size(biotic(n)%norm_remin%nlat) else !}{ len_n = 0 endif !} if (len_e .ne. len_w .or. len_w .ne. len_s .or. len_s .ne. len_n) then !{ call mpp_error(FATAL, trim(error_header) // ' Region sizes are not equal for ' // trim(biotic(n)%name)) endif !} if (size(biotic(n)%norm_remin%t_mask) .ne. 12) then !{ call mpp_error(FATAL, trim(error_header) // ' t_mask size is not 12 for ' // trim(biotic(n)%name)) endif !} ! ! set all of the values to the default ! biotic(n)%norm_remin%mask(:,:,:) = 1.0 if (len_w .gt. 0) then !{ write (stdoutunit,*) write (stdoutunit,*) trim(note_header), 'Process norm_remin array for ', trim(biotic(n)%name) write (stdoutunit,*) ! ! set values for this time-level ! done = 0 do l = 1, 12 !{ if (biotic(n)%norm_remin%t_mask(l)) then !{ if (done .eq. 0) then !{ ! ! set the values via the input values, saving this time index ! afterwards ! write (stdoutunit,*) 'Assigning month ', l call set_array(biotic(n)%norm_remin%mask(:,:,l), & isd, ied, jsd, jed, grid_xt, grid_yt, grid_kmt, & len_w, biotic(n)%norm_remin%wlon, & biotic(n)%norm_remin%elon, & biotic(n)%norm_remin%slat, & biotic(n)%norm_remin%nlat, & biotic(n)%norm_remin%factor, 1.0, & 'Normal remineralization', biotic(n)%norm_remin%coastal_only) done = l else !}{ ! ! Duplicate the values for a previous time-level ! write (stdoutunit,*) 'Duplicating month ', done, ' as ', l biotic(n)%norm_remin%mask(:,:,l) = biotic(n)%norm_remin%mask(:,:,done) endif !} endif !} enddo !} l endif !} enddo !} n ! !----------------------------------------------------------------------- ! read in the no_caco3 namelist data !----------------------------------------------------------------------- ! do n = 1, instances !{ ! ! Check some things ! if (associated(biotic(n)%no_caco3%wlon)) then !{ len_w = size(biotic(n)%no_caco3%wlon) else !}{ len_w = 0 endif !} if (associated(biotic(n)%no_caco3%elon)) then !{ len_e = size(biotic(n)%no_caco3%elon) else !}{ len_e = 0 endif !} if (associated(biotic(n)%no_caco3%slat)) then !{ len_s = size(biotic(n)%no_caco3%slat) else !}{ len_s = 0 endif !} if (associated(biotic(n)%no_caco3%nlat)) then !{ len_n = size(biotic(n)%no_caco3%nlat) else !}{ len_n = 0 endif !} if (len_e .ne. len_w .or. len_w .ne. len_s .or. len_s .ne. len_n) then !{ call mpp_error(FATAL, trim(error_header) // ' Region sizes are not equal for ' // trim(biotic(n)%name)) endif !} if (size(biotic(n)%no_caco3%t_mask) .ne. 12) then !{ call mpp_error(FATAL, trim(error_header) // ' t_mask size is not 12 for ' // trim(biotic(n)%name)) endif !} ! ! set all of the values to the default ! biotic(n)%no_caco3%mask(:,:,:) = 1.0 if (len_w .gt. 0) then !{ write (stdoutunit,*) write (stdoutunit,*) trim(note_header), 'Process no_caco3 array for ', trim(biotic(n)%name) write (stdoutunit,*) ! ! set values for this time-level ! done = 0 do l = 1, 12 !{ if (biotic(n)%no_caco3%t_mask(l)) then !{ if (done .eq. 0) then !{ ! ! set the values via the input values, saving this time index ! afterwards ! write (stdoutunit,*) 'Assigning month ', l call set_array(biotic(n)%no_caco3%mask(:,:,l), & isd, ied, jsd, jed, grid_xt, grid_yt, grid_kmt, & len_w, biotic(n)%no_caco3%wlon, & biotic(n)%no_caco3%elon, & biotic(n)%no_caco3%slat, & biotic(n)%no_caco3%nlat, & biotic(n)%no_caco3%factor, 1.0, & 'Carbonate inhibition', biotic(n)%no_caco3%coastal_only) done = l else !}{ ! ! Duplicate the values for a previous time-level ! write (stdoutunit,*) 'Duplicating month ', done, ' as ', l biotic(n)%no_caco3%mask(:,:,l) = biotic(n)%no_caco3%mask(:,:,done) endif !} endif !} enddo !} l endif !} enddo !} n ! !----------------------------------------------------------------------- ! read in the nut_depl namelist data !----------------------------------------------------------------------- ! do n = 1, instances !{ ! ! Check some things ! if (associated(biotic(n)%nut_depl%wlon)) then !{ len_w = size(biotic(n)%nut_depl%wlon) else !}{ len_w = 0 endif !} if (associated(biotic(n)%nut_depl%elon)) then !{ len_e = size(biotic(n)%nut_depl%elon) else !}{ len_e = 0 endif !} if (associated(biotic(n)%nut_depl%slat)) then !{ len_s = size(biotic(n)%nut_depl%slat) else !}{ len_s = 0 endif !} if (associated(biotic(n)%nut_depl%nlat)) then !{ len_n = size(biotic(n)%nut_depl%nlat) else !}{ len_n = 0 endif !} if (len_e .ne. len_w .or. len_w .ne. len_s .or. len_s .ne. len_n) then !{ call mpp_error(FATAL, trim(error_header) // ' Region sizes are not equal for ' // trim(biotic(n)%name)) endif !} if (size(biotic(n)%nut_depl%t_mask) .ne. 12) then !{ call mpp_error(FATAL, trim(error_header) // ' t_mask size is not 12 for ' // trim(biotic(n)%name)) endif !} ! ! set all of the values to the default ! biotic(n)%nut_depl%mask(:,:,:) = 1.0 if (len_w .gt. 0) then !{ write (stdoutunit,*) write (stdoutunit,*) trim(note_header), 'Process nut_depl array for ', trim(biotic(n)%name) write (stdoutunit,*) ! ! set values for this time-level ! done = 0 do l = 1, 12 !{ if (biotic(n)%nut_depl%t_mask(l)) then !{ if (done .eq. 0) then !{ ! ! set the values via the input values, saving this time index ! afterwards ! write (stdoutunit,*) 'Assigning month ', l call set_array(biotic(n)%nut_depl%mask(:,:,l), & isd, ied, jsd, jed, grid_xt, grid_yt, grid_kmt, & len_w, biotic(n)%nut_depl%wlon, & biotic(n)%nut_depl%elon, & biotic(n)%nut_depl%slat, & biotic(n)%nut_depl%nlat, & biotic(n)%nut_depl%factor, 1.0, & 'Nutrient depletion', biotic(n)%nut_depl%coastal_only) done = l else !}{ ! ! Duplicate the values for a previous time-level ! write (stdoutunit,*) 'Duplicating month ', done, ' as ', l biotic(n)%nut_depl%mask(:,:,l) = biotic(n)%nut_depl%mask(:,:,done) endif !} endif !} enddo !} l endif !} enddo !} n ! !----------------------------------------------------------------------- ! read in the r_bio_tau_a namelist data !----------------------------------------------------------------------- ! do n = 1, instances !{ ! ! Check some things ! if (associated(biotic(n)%r_bio_tau_a%wlon)) then !{ len_w = size(biotic(n)%r_bio_tau_a%wlon) else !}{ len_w = 0 endif !} if (associated(biotic(n)%r_bio_tau_a%elon)) then !{ len_e = size(biotic(n)%r_bio_tau_a%elon) else !}{ len_e = 0 endif !} if (associated(biotic(n)%r_bio_tau_a%slat)) then !{ len_s = size(biotic(n)%r_bio_tau_a%slat) else !}{ len_s = 0 endif !} if (associated(biotic(n)%r_bio_tau_a%nlat)) then !{ len_n = size(biotic(n)%r_bio_tau_a%nlat) else !}{ len_n = 0 endif !} if (len_e .ne. len_w .or. len_w .ne. len_s .or. len_s .ne. len_n) then !{ call mpp_error(FATAL, trim(error_header) // ' Region sizes are not equal for ' // trim(biotic(n)%name)) endif !} if (size(biotic(n)%r_bio_tau_a%t_mask) .ne. 12) then !{ call mpp_error(FATAL, trim(error_header) // ' t_mask size is not 12 for ' // trim(biotic(n)%name)) endif !} ! ! set all of the values to the default ! biotic(n)%r_bio_tau_a%mask(:,:,:) = 1.0 if (len_w .gt. 0) then !{ write (stdoutunit,*) write (stdoutunit,*) trim(note_header), 'Process r_bio_tau_a array for ', trim(biotic(n)%name) write (stdoutunit,*) ! ! set values for this time-level ! done = 0 do l = 1, 12 !{ if (biotic(n)%r_bio_tau_a%t_mask(l)) then !{ if (done .eq. 0) then !{ ! ! set the values via the input values, saving this time index ! afterwards ! write (stdoutunit,*) 'Assigning month ', l call set_array(biotic(n)%r_bio_tau_a%mask(:,:,l), & isd, ied, jsd, jed, grid_xt, grid_yt, grid_kmt, & len_w, biotic(n)%r_bio_tau_a%wlon, & biotic(n)%r_bio_tau_a%elon, & biotic(n)%r_bio_tau_a%slat, & biotic(n)%r_bio_tau_a%nlat, & biotic(n)%r_bio_tau_a%factor, 1.0, & 'Primary production limitation', biotic(n)%r_bio_tau_a%coastal_only) done = l else !}{ ! ! Duplicate the values for a previous time-level ! write (stdoutunit,*) 'Duplicating month ', done, ' as ', l biotic(n)%r_bio_tau_a%mask(:,:,l) = biotic(n)%r_bio_tau_a%mask(:,:,done) endif !} endif !} enddo !} l endif !} enddo !} n ! ! multiply by the restoring factor ! do n = 1, instances !{ biotic(n)%r_bio_tau_a%mask(:,:,:) = & biotic(n)%r_bio_tau * biotic(n)%r_bio_tau_a%mask(:,:,:) enddo !} n ! !----------------------------------------------------------------------- ! initialize special arrays for remineralization !----------------------------------------------------------------------- ! do n = 1, instances !{ do k = 1, nk !{ biotic(n)%zforg(k) = & (grid_zw(k) / & biotic(n)%compensation_depth) ** (-biotic(n)%martin_coeff) biotic(n)%zfca(k) = & exp(-(grid_zw(k) - & biotic(n)%compensation_depth) / biotic(n)%ca_remin_depth) enddo !} k enddo !} n ! !----------------------------------------------------------------------- ! initialize special arrays for alternate remineralization !----------------------------------------------------------------------- ! do n = 1, instances !{ do k = 1, nk !{ biotic(n)%zforg_alt(k) = (grid_zw(k) / & biotic(n)%compensation_depth) ** (-biotic(n)%martin_coeff_alt) enddo !} k enddo !} n ! !----------------------------------------------------------------------- ! Read in additional information for a restart. ! ! We must process all of the instances before restoring any files ! as all fields must be registered before the fields are ! restored, and fields from different instances may be in the ! same file. ! ! Note that the restart file names here must be different from ! those for the tracer values. !----------------------------------------------------------------------- ! allocate(restart(instances)) allocate(local_restart_file(instances)) write(stdoutunit,*) do n = 1, instances !{ ! ! Set the suffix for this instance (if instance name is "_", ! then use a blank suffix). ! if (biotic(n)%name(1:1) .eq. '_') then !{ suffix = ' ' else !}{ suffix = '_' // biotic(n)%name endif !} ! ! Check whether we are already using this restart file, if so, ! we do not want to duplicate it in the list of restart files ! since we only read each restart file once. ! ind = 0 do l = 1, num_restart if (biotic(n)%local_restart_file == local_restart_file(l)) then ind = l exit endif end do if (ind .eq. 0) then num_restart = num_restart + 1 ind = num_restart local_restart_file(ind) = trim(biotic(n)%local_restart_file) end if ! ! Check whether the field already exists in the restart file. ! If not, then set a default value. ! fld_exist = field_exist('INPUT/' // trim(biotic(n)%local_restart_file), 'htotal' // trim(suffix) ) if ( fld_exist ) then !{ write (stdoutunit,*) trim(note_header), & 'Reading additional information for instance ', & ': Initializing instance ', trim(biotic(n)%name) else !}{ write (stdoutunit,*) trim(note_header), & 'Initializing instance ', trim(biotic(n)%name) biotic(n)%htotal(:,:) = htotal_in endif !} ! ! Register the field for restart ! id_restart = register_restart_field(restart(ind), biotic(n)%local_restart_file, & 'htotal' // trim(suffix), biotic(n)%htotal, & domain=mpp_domain2d, mandatory=fld_exist ) enddo !} ! ! Restore the restart fields if the file exists ! do l = 1, num_restart if (file_exist('INPUT/' // trim(local_restart_file(l)))) then call restore_state(restart(l)) end if end do deallocate(local_restart_file) ! !----------------------------------------------------------------------- ! ! initialize some arrays which are held constant for this ! simulation ! !----------------------------------------------------------------------- ! do n = 1, instances !{ biotic(n)%sio2(:,:) = biotic(n)%sio2_const enddo !} n ! !----------------------------------------------------------------------- ! Set up analyses !----------------------------------------------------------------------- ! ! ! register the global fields ! suffix = '_' // package_name long_suffix = ' (' // trim(package_name) // ')' id_o2_sat = register_diag_field(trim(diag_name), & 'o2_saturation' // trim(suffix), grid_tracer_axes(1:2), & model_time, 'O2 saturation' // trim(long_suffix), 'mol/m^3/atm', & missing_value = -1.0e+10) ! ! register the instance fields ! do n = 1, instances !{ if (instances .eq. 1) then !{ suffix = ' ' long_suffix = ' ' else !}{ suffix = '_' // biotic(n)%name long_suffix = ' (' // trim(biotic(n)%name) // ')' endif !} biotic(n)%id_sc_co2 = register_diag_field(trim(diag_name), & 'sc_co2' // trim(suffix), grid_tracer_axes(1:2), & model_time, 'Schmidt number - CO2' // trim(long_suffix), ' ', & missing_value = -1.0e+10) biotic(n)%id_sc_o2 = register_diag_field(trim(diag_name), & 'sc_o2' // trim(suffix), grid_tracer_axes(1:2), & model_time, 'Schmidt number - O2' // trim(long_suffix), ' ', & missing_value = -1.0e+10) biotic(n)%id_alpha = register_diag_field(trim(diag_name), & 'alpha' // trim(suffix), grid_tracer_axes(1:2), & model_time, 'Alpha CO2' // trim(long_suffix), 'mol/kg/atm', & missing_value = -1.0e+10) biotic(n)%id_csurf = register_diag_field(trim(diag_name), & 'csurf' // trim(suffix), grid_tracer_axes(1:2), & model_time, 'CO2* water' // trim(long_suffix), 'mol/kg', & missing_value = -1.0e+10) biotic(n)%id_pco2surf = register_diag_field(trim(diag_name), & 'pco2surf' // trim(suffix), grid_tracer_axes(1:2), & model_time, 'Oceanic pCO2' // trim(long_suffix), 'ppm', & missing_value = -1.0e+10) biotic(n)%id_flux_poc = register_diag_field(trim(diag_name), & 'flux_poc' // trim(suffix), grid_tracer_axes(1:2), & model_time, 'POC flux' // trim(long_suffix), 'mol/m^2/s', & missing_value = -1.0e+10) biotic(n)%id_flux_caco3 = register_diag_field(trim(diag_name), & 'flux_caco3' // trim(suffix), grid_tracer_axes(1:2), & model_time, 'CaCO3 flux' // trim(long_suffix), 'mol/m^2/s', & missing_value = -1.0e+10) biotic(n)%id_htotal = register_diag_field(trim(diag_name), & 'htotal' // trim(suffix), grid_tracer_axes(1:2), & model_time, 'H+ ion concentration' // trim(long_suffix), ' ', & missing_value = -1.0e+10) biotic(n)%id_jprod = register_diag_field(trim(diag_name), & 'jprod' // trim(suffix), grid_tracer_axes(1:3), & model_time, 'Restoring production' // trim(long_suffix), 'mol/kg/s', & missing_value = -1.0e+10) biotic(n)%id_jca = register_diag_field(trim(diag_name), & 'jca' // trim(suffix), grid_tracer_axes(1:3), & model_time, 'CaCO3 change' // trim(long_suffix), 'mol/kg/s', & missing_value = -1.0e+10) biotic(n)%id_jpo4 = register_diag_field(trim(diag_name), & 'jpo4' // trim(suffix), grid_tracer_axes(1:3), & model_time, 'PO4 source' // trim(long_suffix), 'mol/kg/s', & missing_value = -1.0e+10) biotic(n)%id_jdop = register_diag_field(trim(diag_name), & 'jdop' // trim(suffix), grid_tracer_axes(1:3), & model_time, 'DOP source' // trim(long_suffix), 'mol/kg/s', & missing_value = -1.0e+10) biotic(n)%id_jo2 = register_diag_field(trim(diag_name), & 'jo2' // trim(suffix), grid_tracer_axes(1:3), & model_time, 'O2 source' // trim(long_suffix), 'mol/kg/s', & missing_value = -1.0e+10) biotic(n)%id_fc = register_diag_field(trim(diag_name), & 'fc' // trim(suffix), grid_tracer_axes(1:3), & model_time, 'POP change' // trim(long_suffix), 'mol/m^2/s', & missing_value = -1.0e+10) biotic(n)%id_fca = register_diag_field(trim(diag_name), & 'fca' // trim(suffix), grid_tracer_axes(1:3), & model_time, 'CaCO3 change' // trim(long_suffix), 'mol/m^2/s', & missing_value = -1.0e+10) if (biotic(n)%add_phosphate) then !{ biotic(n)%id_jpo4_add = register_diag_field(trim(diag_name), & 'jpo4_add' // trim(suffix), grid_tracer_axes(1:3), & model_time, 'Additional phosphate' // trim(long_suffix), 'mol/kg/s', & missing_value = -1.0e+10) endif !} enddo !} n ! ! integrate the total concentrations of some tracers ! for the start of the run ! total_phosphate = 0.0 total_dop = 0.0 total_o2 = 0.0 total_dic = 0.0 total_alkalinity = 0.0 ! ! Use taup1 time index for the start of a run, and taup1 time ! index for the end of a run so that we are integrating the ! same time level and should therefore get identical results ! write (stdoutunit,*) trim(note_header), & 'Global integrals at start of run' do n = 1, instances !{ do k = 1,nk !{ do j = jsc, jec !{ do i = isc, iec !{ total_phosphate = total_phosphate + & t_prog(biotic(n)%ind_po4)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) total_dop = total_dop + & t_prog(biotic(n)%ind_dop)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) total_o2 = total_o2 + & t_prog(biotic(n)%ind_o2)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) total_dic = total_dic + & t_prog(biotic(n)%ind_dic)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) total_alkalinity = total_alkalinity + & t_prog(biotic(n)%ind_alk)%field(i,j,k,taup1) * & grid_dat(i,j) * grid_tmask(i,j,k) * rho_dzt(i,j,k,taup1) enddo !} i enddo !} j enddo !} k call mpp_sum(total_phosphate) call mpp_sum(total_dop) call mpp_sum(total_o2) call mpp_sum(total_dic) call mpp_sum(total_alkalinity) write (stdoutunit,*) ' Instance ', trim(biotic(n)%name) write (stdoutunit, & '(/'' Total phosphate = '',es19.12,'' Gmol'')') & total_phosphate * 1.0e-09 write (stdoutunit, & '(/'' Total DOP = '',es19.12,'' Gmol'')') & total_DOP * 1.0e-09 write (stdoutunit, & '(/'' Total O2 = '',es19.12,'' Gmol'')') & total_o2 * 1.0e-09 write (stdoutunit, & '(/'' Total DIC = '',es19.12,'' Gmol'')') & total_dic * 1.0e-09 write (stdoutunit, & '(/'' Total alkalinity = '',es19.12,'' Geq'')') & total_alkalinity * 1.0e-09 write (stdoutunit, & '(/'' Total phosphorous = '',es19.12,'' Gmol'')') & (total_phosphate + total_dop) * 1.0e-09 write (stdoutunit, & '(/'' Total real O2 = '',es19.12,'' Gmol'')') & (total_o2 + biotic(n)%o_2_p * total_phosphate) * 1.0e-09 write (stdoutunit, & '(/'' Total Carbon = '',es19.12,'' Gmol'')') & (total_dic + biotic(n)%c_2_p * total_dop) * 1.0e-09 write (stdoutunit, & '(/'' Total real alkalinity = '',es19.12,'' Geq'')') & (total_alkalinity + biotic(n)%n_2_p * total_phosphate) * 1.0e-09 enddo !} n ! !----------------------------------------------------------------------- ! give info !----------------------------------------------------------------------- ! write(stdoutunit,*) write(stdoutunit,*) trim(note_header), 'Tracer runs initialized' write(stdoutunit,*) return end subroutine ocmip2_biotic_start !} ! NAME="ocmip2_biotic_start" !####################################################################### ! ! ! ! Set up an array covering the model domain with a user-specified ! value, in user-specified regions. There are a given number of ! 2-d regions specified by the values slat, nlat, wlon and elon. ! The longitudes are for a cyclic domain, and if wlon and elon ! are on opposite sides of the cut, the correct thing will ! be done. Elon is considered to be east of wlon, so if elon is ! less than wlon, then the region east of elon to the cut will be ! filled, and the region from the cut to wlon will be filled. ! ! After setting up the array in this routine, it may prove useful ! to allow fine-tuning the settings via an array in a namelist. ! ! Arguments: ! Input: ! num_regions = number of user-specified regions which will be ! filled ! ! wlon = 1-d array of western (starting) longitudes for the ! rectangular regions ! ! elon = 1-d array of eastern (ending) longitudes for the ! rectangular regions ! ! slat = 1-d array of southern (starting) latitudes for the ! rectangular regions ! ! nlat = 1-d array of northern (ending) latitudes for the ! rectangular regions ! ! Note: if slat >= nlat, then nothing is done ! for that region ! ! set_value = the value to assign to array in the user-specified ! regions ! ! unset_value = the value to assign to array outside of the ! user-specified regions ! ! name = character variable used in informative messages ! ! coastal_only = true to limit changes only to coastal points ! (i.e., at least one bordering point is land) ! ! Output: ! ! array = 2-d array which will contain the set- and unset- ! values. The array is assumed to have a border ! one unit wide on all edges, ala MOM. A cyclic ! boundary condition will be set if requested. ! ! subroutine set_array(array, isd, ied, jsd, jed, & xt, yt, kmt, & num_regions, wlon_in, elon_in, slat, nlat, & set_value, unset_value, name, & coastal_only) !{ implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! integer, intent(in) :: isd integer, intent(in) :: ied integer, intent(in) :: jsd integer, intent(in) :: jed integer, intent(in) :: num_regions real, dimension(isd:ied,jsd:jed), intent(out) :: array logical, intent(in) :: coastal_only real, dimension(num_regions), intent(in) :: elon_in integer, dimension(isd:,jsd:), intent(in) :: kmt character(len=*), intent(in) :: name real, dimension(num_regions), intent(in) :: nlat real, intent(in) :: set_value real, dimension(num_regions), intent(in) :: slat real, intent(in) :: unset_value real, dimension(num_regions), intent(in) :: wlon_in real, dimension(isd:,jsd:), intent(in) :: xt real, dimension(isd:,jsd:), intent(in) :: yt ! ! local parameters ! character(len=64), parameter :: sub_name = 'set_array' character(len=256), parameter :: error_header = & '==>Error from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: warn_header = & '==>Warning from ' // trim(mod_name) // '(' // trim(sub_name) // '):' character(len=256), parameter :: note_header = & '==>Note from ' // trim(mod_name) // '(' // trim(sub_name) // '):' ! ! local variables ! integer :: i, j, n real, dimension(:), allocatable :: wlon real, dimension(:), allocatable :: elon integer :: stdoutunit stdoutunit=stdout() ! ! save the longitudes in case they need to be modified ! allocate(wlon(num_regions)) allocate(elon(num_regions)) wlon(:) = wlon_in(:) elon(:) = elon_in(:) ! ! loop over the regions, applying changes as necessary ! do n = 1, num_regions !{ if (nlat(n) .ge. slat(n)) then !{ write (stdoutunit,*) write (stdoutunit,*) trim(note_header), & trim(name), ' region: ', n ! ! make sure that all longitudes are in the range [0,360] ! do while (wlon(n) .gt. 360.0) !{ wlon(n) = wlon(n) - 360.0 enddo !} do while (wlon(n) .lt. 0.0) !{ wlon(n) = wlon(n) + 360.0 enddo !} do while (elon(n) .gt. 360.0) !{ elon(n) = elon(n) - 360.0 enddo !} do while (elon(n) .lt. 0.0) !{ elon(n) = elon(n) + 360.0 enddo !} ! ! if the southern and northern latitudes are the same, then ! find the grid box which encompasses them ... ! if (slat(n) .eq. nlat(n)) then !{ call mpp_error(FATAL, trim(error_header) // & 'Equal latitudes not supported') elseif (wlon(n) .eq. elon(n)) then !}{ call mpp_error(FATAL, trim(error_header) // & 'Equal longitudes not supported') else !}{ ! ! ... else find all boxes where the center lies in the ! rectangular region ! do j = jsd, jed !{ do i = isd, ied !{ if (nlat(n) .ge. yt(i,j) .and. & slat(n) .le. yt(i,j) .and. & lon_between(xt(i,j), wlon(n), elon(n))) then !{ array(i,j) = set_value endif !} enddo !} i enddo !} j endif !} endif !} enddo !} n ! ! if desired only apply mask to coastal regions ! if (coastal_only) then !{ do j = jsd, jed !{ do i = isd, ied !{ if (kmt(i,j) .ne. 0 .and. & array(i,j) .eq. set_value) then !{ ! ! if all the surrounding points are ocean, then this is not ! a coastal point, therefore reset the mask ! if (kmt(i-1,j) .ne. 0 .and. & kmt(i+1,j) .ne. 0 .and. & kmt(i,j-1) .ne. 0 .and. & kmt(i,j+1) .ne. 0) then !{ array(i,j) = unset_value endif !} endif !} enddo !} i enddo !} j endif !} ! ! clean up ! deallocate(wlon) deallocate(elon) return contains ! ! Return true if w <= x_in <= e, taking into account the ! periodicity of longitude. ! ! x_in = value to test ! ! w = west longitude of boundary ! ! e = east longitude of boundary ! function lon_between(x_in, w, e) !{ implicit none ! !----------------------------------------------------------------------- ! Arguments !----------------------------------------------------------------------- ! real, intent(in) :: x_in real, intent(in) :: w real, intent(in) :: e ! ! function definition ! logical :: lon_between ! ! local variables ! !real :: w !real :: e real :: x ! ! Save input values so we may modify them safely ! x = x_in ! ! make sure that all longitudes are in the range [0,360] ! do while (x .gt. 360.0) !{ x = x - 360.0 enddo !} do while (x .lt. 0.0) !{ x = x + 360.0 enddo !} if (w .gt. e) then !{ lon_between = w .le. x .or. x .le. e else !}{ lon_between = w .le. x .and. x .le. e endif !} return end function lon_between !} end subroutine set_array !} ! NAME="set_array" end module ocmip2_biotic_mod !}