!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !! !! !! 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 !! !! !! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! module atmosphere_mod !----------------------------------------------------------------------- ! ! interface for B-grid dynamical core with Held-Suarez forcing ! !----------------------------------------------------------------------- use bgrid_core_driver_mod, only: bgrid_dynam_type, & bgrid_core_driver_init, & bgrid_core_driver, & bgrid_core_time_diff, & bgrid_core_driver_end, & get_bottom_data, & put_bottom_data, & atmosphere_domain use bgrid_prog_var_mod, only: prog_var_type, var_init use bgrid_horiz_mod, only: get_horiz_grid_size, & get_horiz_grid_bound, TGRID use time_manager_mod, only: time_type, get_time, operator(+) use fms_mod, only: file_exist, open_namelist_file, & error_mesg, FATAL, & check_nml_error, stdlog, & write_version_number, & mpp_pe, mpp_root_pe, & close_file, set_domain ! routines used by subroutine bgrid_physics use bgrid_change_grid_mod, only: change_grid, TEMP_GRID, WIND_GRID use bgrid_horiz_mod , only: horiz_grid_type use bgrid_vert_mod , only: vert_grid_type, & compute_pres_full, compute_pres_half use bgrid_halo_mod , only: update_halo, UWND, VWND, TEMP, & NORTH, EAST, WEST, SOUTH use shallow_physics_mod , only: shallow_physics_init, shallow_physics, & shallow_physics_end !----------------------------------------------------------------------- implicit none private public atmosphere, & atmosphere_init, & atmosphere_end, & atmosphere_resolution, & atmosphere_boundary, & get_atmosphere_axes, & atmosphere_domain !----------------------------------------------------------------------- character(len=128) :: version = '$Id: atmosphere.F90,v 14.0 2007/03/15 21:56:36 fms Exp $' character(len=128) :: tag = '$Name: mom4p1_pubrel_dec2009_nnz $' !----------------------------------------------------------------------- !---- namelist (saved in file input.nml) ---- ! ! physics_window The number of "i" by "j" rows processed each time ! the modular physics is called. To process the entire ! domain use physics_window = (/0,0/). ! [integer, default: physics_window = 0,0] integer, dimension(2) :: physics_window = (/0,0/) namelist /atmosphere_nml/ physics_window !----------------------------------------------------------------------- !---- private data ---- type (bgrid_dynam_type), save :: Dynam type (prog_var_type), save :: Var, Var_dt type (time_type) :: Time_step_atmos real :: dt_atmos real, dimension(:,:,:), pointer :: omega =>NULL() integer, dimension(4) :: atmos_axes !----------------------------------------------------------------------- contains !####################################################################### subroutine atmosphere (Time) type (time_type), intent(in) :: Time type(time_type) :: Time_prev, Time_next !----------------------------------------------------------------------- Time_prev = Time ! two time-level scheme Time_next = Time + Time_step_atmos !---- dynamics ----- call bgrid_core_driver ( Time_next, Var, Var_dt, Dynam, omega ) !---- call physics ----- call bgrid_physics ( physics_window, dt_atmos, Time_next, & Dynam%Hgrid, Dynam%Vgrid, Dynam, & Var, Var_dt ) !---- time differencing and diagnostics ----- call bgrid_core_time_diff ( omega, Time_next, Dynam, Var, Var_dt ) !----------------------------------------------------------------------- end subroutine atmosphere !####################################################################### subroutine atmosphere_init ( Time_init, Time, Time_step ) type (time_type), intent(in) :: Time_init, Time, Time_step integer :: unit, sec, ierr, io !----------------------------------------------------------------------- !----- read namelist ----- if (file_exist('input.nml')) then unit = open_namelist_file ( ) ierr=1; do while (ierr /= 0) read (unit, nml=atmosphere_nml, iostat=io, end=10) ierr = check_nml_error (io, 'atmosphere_nml') enddo 10 call close_file (unit) endif !----- write version and namelist to log file ----- call write_version_number ( version, tag ) if ( mpp_pe() == mpp_root_pe() ) write (stdlog(), nml=atmosphere_nml) !---- compute physics/atmos time step in seconds ---- Time_step_atmos = Time_step call get_time (Time_step_atmos, sec) dt_atmos = real(sec) !----- initialize dynamical core ----- call bgrid_core_driver_init ( Time_init, Time, Time_step, & Var, Var_dt, Dynam, atmos_axes ) !----- make sure code is initialized for shallow water version ----- if (Dynam%Vgrid%nlev /= 1) call error_mesg ('atmosphere_init', & 'number of model levels not correct '// & 'for shallow water version', FATAL) !----- initialize storage needed for vert motion ---- omega => var_init (Dynam%Hgrid, Dynam%Vgrid%nlev) !----- initialize physics interface ----- ! (might want to use alm,aph instead of tlm,tph) call shallow_physics_init ( atmos_axes, Time, & Dynam%Hgrid%Tmp%tlm(Dynam%Hgrid%Tmp%is:Dynam%Hgrid%Tmp%ie), & Dynam%Hgrid%Tmp%tph(Dynam%Hgrid%Tmp%js:Dynam%Hgrid%Tmp%je) ) ! ----- use entire grid as window ? ----- if (physics_window(1) <= 0) physics_window(1) = Dynam%Hgrid%Tmp%ie-Dynam%Hgrid%Tmp%is+1 if (physics_window(2) <= 0) physics_window(2) = Dynam%Hgrid%Tmp%je-Dynam%Hgrid%Tmp%js+1 !----------------------------------------------------------------------- end subroutine atmosphere_init !####################################################################### subroutine atmosphere_end integer :: unit call bgrid_core_driver_end ( Var, Dynam ) end subroutine atmosphere_end !####################################################################### ! returns the number of longitude and latitude grid points ! for either the local PEs grid (default) or the global grid subroutine atmosphere_resolution (nlon, nlat, global) integer, intent(out) :: nlon, nlat logical, intent(in), optional :: global !---- return the size of the grid used for physics computations ---- call get_horiz_grid_size (Dynam % Hgrid, TGRID, nlon, nlat, global) end subroutine atmosphere_resolution !####################################################################### ! returns the longitude and latitude grid box edges ! for either the local PEs grid (default) or the global grid subroutine atmosphere_boundary (blon, blat, global) real, intent(out) :: blon(:), blat(:) logical, intent(in), optional :: global !----- return the longitudinal and latitudinal grid box edges ---------- call get_horiz_grid_bound (Dynam % Hgrid, TGRID, blon, blat, global) end subroutine atmosphere_boundary !####################################################################### ! returns the axis indices associated with the coupling grid subroutine get_atmosphere_axes ( axes ) integer, intent(out) :: axes (:) !----- returns the axis indices for the atmospheric (mass) grid ----- if ( size(axes) < 0 .or. size(axes) > 4 ) call error_mesg ( & 'get_atmosphere_axes in atmosphere_mod', & 'size of argument is incorrect', FATAL ) axes (1:size(axes)) = atmos_axes (1:size(axes)) end subroutine get_atmosphere_axes !####################################################################### subroutine bgrid_physics ( window, dt_phys, Time, Hgrid, Vgrid, & Dynam, Var, Var_dt ) !----------------------------------------------------------------------- ! ! Time = current time (time_type, see time manager) ! !----------------------------------------------------------------------- integer, intent(in) :: window(2) real, intent(in) :: dt_phys type(time_type),intent(in) :: Time type (horiz_grid_type),intent(inout) :: Hgrid type (vert_grid_type),intent(in) :: Vgrid type(bgrid_dynam_type),intent(in) :: Dynam type (prog_var_type),intent(in) :: Var type (prog_var_type),intent(inout) :: Var_dt !----------------------------------------------------------------------- integer :: j, k, n, is, ie, js, je, i1, i2, j1, j2, idim, jdim integer :: timelev = +1 !----------------------------------------------------------------------- real, dimension(Hgrid%ilb:Hgrid%iub, & Hgrid%jlb:Hgrid%jub,1) :: uh, vh, uh_dt, vh_dt, & u_dt, v_dt !----------------------------------------------------------------------- !---------------------------- do physics ------------------------------- idim = window(1) jdim = window(2) ! --- momentum and momentum tendency on mass grid --- call update_halo (Hgrid, UWND, Var_dt%u, halos=SOUTH+WEST) call update_halo (Hgrid, VWND, Var_dt%v, halos=SOUTH+WEST) call change_grid (Hgrid, WIND_GRID, TEMP_GRID, Var%u, Var%v, uh , vh ) call change_grid (Hgrid, WIND_GRID, TEMP_GRID, Var_dt%u, Var_dt%v, uh_dt, vh_dt) u_dt = uh_dt ! save copy of tendency before physics v_dt = vh_dt ! --- loop through physics windows --- js = Hgrid%Tmp%js do while ( js <= Hgrid%Tmp%je ) je = min ( js+jdim-1, Hgrid%Tmp%je ) is = Hgrid%Tmp%is do while ( is <= Hgrid%Tmp%ie ) ie = min ( is+idim-1, Hgrid%Tmp%ie ) ! ---- j-axis indices in the global physics grid ---- j1 = js-Hgrid%Tmp%js+1; j2 = j1+(je-js) i1 = is-Hgrid%Tmp%is+1; i2 = i1+(ie-is) !----------------------------------------------------------------------- !-------------------------- call physics ------------------------------- !------------ (need to add leap-frog option for uh,vh) ----------------- !----------------------------------------------------------------------- call shallow_physics ( i1, i2, j1, j2, timelev, dt_phys, Time ,& uh (is:ie,js:je,1) ,& vh (is:ie,js:je,1) ,& Var%ps (is:ie,js:je) ,& uh (is:ie,js:je,1) ,& vh (is:ie,js:je,1) ,& Var%ps (is:ie,js:je) ,& uh_dt (is:ie,js:je,1) ,& vh_dt (is:ie,js:je,1) ,& Var_dt%ps (is:ie,js:je) ) is = is + idim enddo js = js + jdim enddo !----------------------------------------------------------------------- ! compute momentum tendencies on mass grid for physics only ! udt(phys) = udt(current) - udt(before physics) uh_dt = uh_dt - u_dt vh_dt = vh_dt - v_dt ! update halos of momentum tendencies on mass grid ! then move momentum tendencies to momentum grid call update_halo (Hgrid, TEMP, uh_dt, halos=NORTH+EAST) call update_halo (Hgrid, TEMP, vh_dt, halos=NORTH+EAST) call change_grid (Hgrid, TEMP_GRID, WIND_GRID, uh_dt, vh_dt, u_dt, v_dt) u_dt(:,Hgrid%jub,:) = 0.0 ! zero out unused polar halo row v_dt(:,Hgrid%jub,:) = 0.0 ! no harm done when not polar row !---- update momentum tendencies ---- Var_dt%u = Var_dt%u + u_dt Var_dt%v = Var_dt%v + v_dt !---- update all halo rows ---- call update_halo (Hgrid, TEMP, Var_dt%ps) call update_halo (Hgrid, UWND, Var_dt%u) call update_halo (Hgrid, VWND, Var_dt%v) call update_halo (Hgrid, TEMP, Var_dt%t) call update_halo (Hgrid, TEMP, Var_dt%r) Var_dt % pssl = Var_dt % ps !----------------------------------------------------------------------- end subroutine bgrid_physics !####################################################################### end module atmosphere_mod