!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !! !! !! 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 FV dynamical core with Held-Suarez forcing ! !----------------------------------------------------------------------- use time_manager_mod, only: time_type, get_time, set_time, operator(+) use fms_mod, only: file_exist, open_namelist_file, & error_mesg, FATAL, & check_nml_error, stdlog, & write_version_number, & close_file, set_domain use hs_forcing_mod, only: hs_forcing_init #if defined(MARS_GCM) || defined(VENUS_GCM) || defined(STRAT_GCM) use hs_forcing_mod, only: hs_forcing_end use fms_mod, only: nullify_domain #endif use constants_mod, only: omega, cp_air, rdgas, kappa, radius, grav, rvgas use mpp_domains_mod, only: domain2d !------------------ ! FV specific codes: !------------------ use fv_pack, only: nlon, mlat, nlev, beglat, endlat, beglon, endlon, & rlonb, rlatb, & cold_start, ncnst, pnats, consv_te, ptop, & fv_init, fv_domain, fv_end, change_time, map_dt, & adiabatic, restart_format #if defined(MARS_GCM) || defined(VENUS_GCM) || defined(STRAT_GCM) use fv_pack, only: rlon, rlat, ak, bk, p_ref #endif use fv_diagnostics, only: fv_diag_init, fv_diag, fv_time use timingModule, only: timing_on, timing_off use fv_restart_mod, only: fv_restart, write_fv_rst use fv_dynamics_mod, only: fv_dynamics use fv_phys_mod, only: fv_phys !----------------------------------------------------------------------- implicit none private public atmosphere_init, atmosphere, atmosphere_end, atmosphere_domain !----------------------------------------------------------------------- character(len=128) :: version = '$Id: atmosphere.F90,v 17.0 2009/07/21 02:53:01 fms Exp $' character(len=128) :: tag = '$Name: mom4p1_pubrel_dec2009_nnz $' !----------------------------------------------------------------------- !---- private data ---- type (time_type) :: Time_step_atmos real :: dt_atmos integer :: sec integer days, seconds !----------------------------------------------------------------------- contains !####################################################################### subroutine atmosphere_init ( Time_init, Time, Time_step ) type (time_type), intent(in) :: Time_step type (time_type), intent(in) :: Time_init type (time_type), intent(in) :: Time ! local: integer axes(4) integer ss, ds #if defined(MARS_GCM) || defined(VENUS_GCM) || defined(STRAT_GCM) integer ii, jj, k real, allocatable :: p_std(:) real, allocatable :: rlonb2d(:,:), rlatb2d(:,:) #endif !----- write version and namelist to log file ----- call write_version_number ( version, tag ) !---- compute physics/atmos time step in seconds ---- Time_step_atmos = Time_step call get_time (Time_step_atmos, sec) dt_atmos = real(sec) !----- initialize FV dynamical core ----- call fv_init( sec ) call fv_restart( days, seconds ) if ( .not. cold_start ) then ! ! Check consistency in Time ! fv_time = set_time (seconds, days) call get_time (Time, ss, ds) if( seconds /= ss .or. days /= ds ) then ! write(6,*) 'FMS:', ds, ss ! write(6,*) 'FV:', days, seconds call error_mesg('FV_init:', & 'Time inconsistent between fv_rst and INPUT/atmos_model.res', & FATAL) endif else fv_time = time endif call fv_diag_init( axes, Time ) #if defined(MARS_GCM) || defined(VENUS_GCM) || defined(STRAT_GCM) if( nlev > 1 ) then allocate( p_std(nlev+1) ) DO k= 1, nlev+1 p_std(k)= ak(k) + bk(k)*p_ref ENDDO allocate ( rlatb2d(nlon+1,beglat:endlat+1) ) allocate ( rlonb2d(nlon+1,beglat:endlat+1) ) ! Convert bounding lon and lat arrays to 2-D DO jj= beglat, endlat+1 rlatb2d(beglon:endlon+1,jj)= rlatb(jj) ENDDO DO ii= beglon, endlon+1 rlonb2d(ii,beglat:endlat+1)= rlonb(ii) ENDDO call set_domain ( fv_domain ) call hs_forcing_init ( nlon, mlat, nlev, & rlonb2d(beglon:endlon+1,beglat:endlat+1), & rlatb2d(beglon:endlon+1,beglat:endlat+1), & rlon (beglon:endlon, beglat:endlat), & rlat (beglon:endlon, beglat:endlat), & p_std, axes, Time ) deallocate( p_std ) call nullify_domain( ) endif #else if( nlev > 1 ) call hs_forcing_init ( axes, Time ) #endif !----------------------------------------------------------------------- end subroutine atmosphere_init !####################################################################### subroutine atmosphere (Time) #include type(time_type), intent(in) :: Time ! local: real zvir logical p_map ! Perform only partial remapping #include fv_time = Time + Time_step_atmos call get_time (fv_time, seconds, days) !---- call fv dynamics ----- zvir = 0. ! no virtual effect if not full physics if ( mod(seconds, map_dt) == 0 ) then p_map = .false. else p_map = .true. endif call timing_on('fv_dynamics') #ifndef USE_LIMA call fv_dynamics (nlon, mlat, nlev, beglat, endlat, & ncnst, pnats, p_map, consv_te, & u, v, delp, pt, q, & ps, pe, pk, pkz, phis, & omga, peln, ptop, omega, sec, & zvir, cp_air, rdgas, kappa, radius, ua, va, fv_time ) #else call fv_dynamics (nlon, mlat, nlev, beglat, endlat, & ncnst, pnats, p_map, consv_te, & u, v, delp, pt, q, & ps, pe, pk, pkz, phis, & omga, peln, ptop, omega, sec, & zvir, cp_air, rdgas, kappa, radius, ua, va ) #endif call timing_off('fv_dynamics') if( nlev /=1 .and. .not. adiabatic ) then call timing_on('FV_PHYS') call fv_phys ( fv_time, sec ) call timing_off('FV_PHYS') endif !---- diagnostics for FV dynamics ----- call timing_on('FV_DIAG') !!!!rjw hsphys= .true. ---> uses simple calender: days & seconds call fv_diag(fv_time, nlon, mlat, nlev, beglat, endlat, ncnst, zvir, & dt_atmos, .true.) call timing_off('FV_DIAG') end subroutine atmosphere subroutine atmosphere_end !----- initialize domains for writing global physics data ----- call set_domain ( fv_domain ) call get_time (fv_time, seconds, days) #if defined(MARS_GCM) || defined(VENUS_GCM) || defined(STRAT_GCM) call hs_forcing_end( days ) call nullify_domain( ) #endif call write_fv_rst( 'RESTART/fv_rst.res', days, seconds, grav, & restart_format ) call fv_end(days, seconds) end subroutine atmosphere_end subroutine atmosphere_domain (Domain) type(domain2d), intent(inout) :: Domain Domain = fv_domain end subroutine atmosphere_domain end module atmosphere_mod