! Stefan Petri ! ! ! ! !----------------------------------------------------------------------- ! ! Driver for the atmospheric model Aeolus2, contains routines to advance the ! atmospheric model state by one time step. ! ! This version of atmos_model_mod is going to implement ! a coupling interface to the atmosphere model Aeolus2, which is ! basedon a pseudo-spectral moist-convective Thermal Rotating Shallow ! Water model (mcTRSW). ! It requires two routines to advance ! the atmospheric model one time step into the future, up and and down. ! The boundary variables needed by other component models for coupling ! are contained in a derived data type. A variable of this derived type ! is returned when initializing the atmospheric model. It is used by other ! routines in this module and by coupling routines. The contents of ! this derived type should only be modified by the atmospheric model. ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !! !! !! 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 module atmosphere_mod use time_manager_mod, only: time_type, get_time, set_time, get_date, set_date, operator(+), & operator(/), & days_in_year, get_calendar_type, THIRTY_DAY_MONTHS, NOLEAP use fms_mod, only: & file_exist, open_namelist_file, & error_mesg, FATAL, & check_nml_error, stdlog, stdout, & write_version_number, & mpp_pe, mpp_root_pe, & close_file, set_domain, & mpp_clock_id, mpp_clock_begin, & mpp_clock_end, CLOCK_SUBCOMPONENT, & clock_flag_default, & mpp_npes !use fms_io_mod, only: & ! return_domain use mpp_domains_mod, only: domain1D, domain2D, & mpp_domains_init, mpp_define_domains, & mpp_get_compute_domain, mpp_get_compute_domains, & mpp_get_domain_components, mpp_get_pelist use field_manager_mod, only: MODEL_ATMOS use constants_mod, only: omega, cp_air, rdgas, grav, rvgas, kappa, radius, pstd_mks, PI, stefan use stock_constants_mod, only : ISTOCK_WATER, ISTOCK_HEAT, ISTOCK_SALT use diag_manager_mod, only: diag_axis_init, register_diag_field use diag_manager_mod, only: register_static_field, send_data use diag_manager_mod, only: diag_manager_end use mpp_mod, only: mpp_error, mpp_get_current_pelist, mpp_sum use tracer_manager_mod, only: get_number_tracers, get_tracer_index, NO_TRACER use xgrid_mod, only: grid_box_type use grid_mod, only: get_grid_cell_area, get_grid_cell_centers use topography_mod, only: get_topog_mean, get_topog_stdev ! we now have two alternate ways to get the cosine of the zenith angle use astronomy_mod, only: astronomy_init, daily_mean_solar !, annual_mean_solar, get_period !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler use ifcore use aeolus2fms_mod, only: aeolus2fms_interface_check use aeolus2fms_mod, only: aeolus2fms_init_grid, aeolus2fms_get_axes_coords, aeolus2fms_init use aeolus2fms_mod, only: aeolus2fms_atmosphere_down, aeolus2fms_atmosphere_up use aeolus2fms_mod, only: aeolus2fms_finish, aeolus2fms_restart !use aeolus2fms_mod, only: aeolus2fms_get_boundary use aeolus2fms_mod, only: aeolus2fms_get_bottom_mass, aeolus2fms_get_bottom_wind, aeolus2fms_get_stock_pe !----------------------------------------------------------------------- implicit none private public atmosphere_down, atmosphere_up, & atmosphere_init, atmosphere_end, & atmosphere_resolution, atmosphere_boundary, & atmosphere_cell_area, atmosphere_restart, & get_atmosphere_axes, atmosphere_domain public get_bottom_mass, get_bottom_wind, get_stock_pe public surf_diff_type !----------------------------------------------------------------------- ! type surf_diff_type real, pointer, dimension(:,:) :: dtmass => NULL() ! dt/mass, where dt = atmospheric time step (sec) ! mass = mass per unit area of lowest atmospheric layer (Kg/m2) real, pointer, dimension(:,:) :: dflux_t => NULL() ! derivative of implicit part of downward temperature flux at the top of the lowest ! atmospheric layer with respect to the temperature ! of that layer (J/(m2 K)) real, pointer, dimension(:,:,:) :: dflux_tr => NULL() ! derivative of implicit part of downward flux of specific humidity ! at the top of the lowest atmospheric layer with respect to ! the specific humidity of that layer (kg/(m2 s)) real, pointer, dimension(:,:) :: delta_t => NULL() ! the increment in temperature in the lowest atmospheric ! layer (((i+1)-(i-1) if atmos model is leapfrog) (K) ! (defined in gcm_vert_diff_down as the increment computed up ! to this point in model, including effect of vertical ! diffusive flux at top of lowest model level, presumed ! to be modified to include effects of surface fluxes ! outside of this module, then used to start upward ! tridiagonal sweep, real, pointer, dimension(:,:,:):: delta_tr => NULL() ! similarly for the increment in specific humidity ! (non-dimensional = Kg/Kg) real, pointer, dimension(:,:) :: delta_u => NULL() real, pointer, dimension(:,:) :: delta_v => NULL() real, pointer, dimension(:,:) :: sst_miz => NULL() ! only used if /flux_exchange_nml/do_forecast is .true. end type surf_diff_type ! !----------------------------------------------------------------------- character(len=128) :: version = '$Id: atmosphere.F90 3862 2021-12-06 15:51:04Z petri $' character(len=128) :: tagname = '$Name: mom4p1_pubrel_dec2009_nnz $' !---- namelist (saved in file input.nml) ---- ! ! ! Longitudinal resolution of atmosphere grid (at the equator). ! Must match the data in the grid specification file. ! ! ! Latitudinal resolution of atmosphere grid. ! Must match the data in the grid specification file. ! !! !! Number of vertical layers of atmosphere grid in the Troposphere !! !! !! Number of vertical layers of atmosphere grid in the Stratosphere !! ! ! If true, use flat topography in Aeolus2, else use real topography heights ! ! ! If true, update Aeolus2 topo.land_frac in every _down and _up step. ! Also update the surface type fractions. ! If false, do it only in the very first update step. ! To update in every step might be useful if the land or ocean might decide to change the ! land sea mask. With current implementations, that seems never to be the case. ! ! ! Freeze time for obtaining astronomical data. ! This can be useful to study spin-up behaviour / stability / equilibrium without interference of seasonal (or diurnal) cycle. ! ! ! If fix_astro is true, this specifies the time to use ! ! ! Solar constant ! ! namelist /atmos_aeolus2_nml/ & NoNx, NoNy, & ! Number of grid cells in lon, lat !NoNz, NoNz_Strato, & ! Number of grid cells in height (Tropo+Strato) lon_0_360, & ! whether longitudes run 0 to 360 in the grid spec file without_topo, & ! choose between flat and real topography update_land_frac_always, & ! wether to update land_frac (and with it the surface type fractions) ! in every _down and _up step, or only in the first one fix_astro, & ! freeze astronomical time, i.e. keep time for obtaining astronomical data constant astro_date, & ! if fix_astro, this gives the time to use solar_constant !----------------------------------------------------------------------- !---- private data ---- integer :: NoNx, NoNy !, NoNz, NoNz_Strato ! values from namelist logical :: lon_0_360 = .false. ! .true. was the default for old Aeolus logical :: without_topo = .false. logical :: update_land_frac_always = .false. type (time_type) :: Time_step_atmos real :: dt_atmos ! Time_step_atmos converted to seconds integer, dimension(4) :: atmos_axes ! dimension(4) is required by atmos_coupled_mod logical :: fix_astro = .false. ! freeze astronomical time, i.e. keep time for solar data constant integer, dimension(6) :: astro_date = (/ 0, 0, 0, 0, 0, 0 /) ! year, month, day, hour, minute, second type (time_type) :: Time_fix_astro ! if fix_astro, this gives the time to use integer :: fixdayoftheyear = -1 ! if fix_astro, this the day of the year to use for sinsol real :: solar_constant = 1365.0 ! clock ids for timing integer :: clock_id_aeolus2 integer :: clockid_init, clockid_up, clockid_down, clockid_diag logical :: module_is_initialized =.false. integer :: nr_tracers = 0 ! total number of tracers in surf_diff_type fields integer :: q_ind = NO_TRACER ! index of humidity in tracers array integer :: co2_ind = NO_TRACER ! index of CO2 in tracers array type(domain2D) :: grid_domain integer :: is, ie, js, je ! local domain boundary indices real, dimension(:), allocatable :: lonb, latb, lon, lat ! , level, levelb ! for the _global_ axes ! diag var ids integer :: id_t = -1, id_astro_cosz = -1, id_astro_solar = -1 integer :: atmos_stock_water_id = -1, atmos_stock_heat_id = -1 integer :: id_flux_lw_sf_up = -1 ! to diagnose problems with interpolation to reduced grid ! for the FMS astronomy_mod integer :: DaysPerYear real, dimension(:), allocatable :: rad_lat ! grid-cell midpoint coordinates in meridional direction, in radians !------------------------- axis names ---------------------------------- character(len=8) :: axiset = 'aeolus2' character(len=8) :: mod_name = 'aeolus2' !----------------------------------------------------------------------- !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler integer(4) :: original_fpe_flags integer(4) :: strict_fpe_flags = FPE_M_TRAP_UND + FPE_M_TRAP_OVF + FPE_M_TRAP_DIV0 & + FPE_M_TRAP_INV + FPE_M_ABRUPT_UND + FPE_M_ABRUPT_DMZ !! Some Fortran-operatiosn apparently abuse floating point registers for fast (?) transfer !! of non-floating-point data. This can trigger unwanted exceptions. !! Lets try suppress some of that. integer(4) :: loose_fpe_flags = FPE_M_TRAP_DIV0 + FPE_M_ABRUPT_UND + FPE_M_ABRUPT_DMZ contains ! ================================================================================== ! ! ! ! public routine required for atmospheric components of coupled models ! Read in restart files and initialize arrays ! ! ! ! public routine required for atmospheric components of coupled models ! Read in restart files and initialize arrays ! ! ! ! ! ! The base (or initial) time of the experiment. ! ! ! The current time. ! ! ! The atmospheric model/physics time step. ! ! ! The surface terms for vertical diffusion that are exchanged ! with other component models. On input fields have not been ! allocated, on output all arrays are allocated. ! ! ! Contains information about the grid cells ! In the 18dec2009 release, all other atmos models set the fields ! of this parameter simply to zero, or even just ignore it entirely. ! So do we, for now. ! subroutine atmosphere_init (Time_init, Time, Time_step, Surf_diff, Grid_box) type (time_type), intent(in) :: Time_init, Time, Time_step type(surf_diff_type), intent(inout) :: Surf_diff type(grid_box_type), intent(inout) :: Grid_box integer :: Time_init_year, Time_init_month, Time_init_day, Time_init_hours, Time_init_minutes, Time_init_seconds integer :: Time_year, Time_month, Time_days, Time_hours, Time_seconds integer :: Time_step_days, Time_step_seconds integer :: unit, ierr, io integer :: log_unit integer :: pe, npes, layout(2) integer :: atmos_comm ! the MPI communicator for the atmosphere integer, allocatable, dimension(:) :: atmos_pelist ! needed to obtain the peset integer :: id_lonb, id_latb !, id_levelb integer, dimension(1:2) :: axes_2d logical :: used ! global axis info retrieved from exchange grid, for comparison with ! the axis info retrieved from the python side. real, dimension(:), allocatable :: exglon, exglat, londiff, latdiff ! 2-D grid cell boundary arrays, for calls to topography_mod real, dimension(:,:), allocatable :: lon_boundaries, lat_boundaries logical :: got_topo real, _ALLOCATABLE :: SIGORO(:,:) ! standard deviation of topography on local domain real, _ALLOCATABLE :: HORO(:,:) ! Orography height, averaged over grid cell real, _ALLOCATABLE :: land_frac(:,:) ! provisional land fraction mask. The real one is delivered by the coupler in every down and up step. real, _ALLOCATABLE :: area(:,:) ! grid cell area (obtained from the exchange grid) integer id_horo, id_sigoro, id_area !character(len=64) :: filename integer(4) :: previous_fpe_flags if (module_is_initialized) return !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler call clearstatusfpqq() ! do not inherit signalled FP exceptions from FMS original_fpe_flags = FOR_SET_FPE(strict_fpe_flags) clock_id_aeolus2 = mpp_clock_id('Aeolus2', & flags=clock_flag_default, grain=CLOCK_SUBCOMPONENT ) clockid_init = mpp_clock_id('Aeolus2_init', & flags=clock_flag_default, grain=CLOCK_SUBCOMPONENT ) clockid_up = mpp_clock_id('Aeolus2_update_up', & flags=clock_flag_default, grain=CLOCK_SUBCOMPONENT ) clockid_down = mpp_clock_id('Aeolus2_update_down', & flags=clock_flag_default, grain=CLOCK_SUBCOMPONENT ) clockid_diag = mpp_clock_id('Aeolus2_diagnostic', & flags=clock_flag_default, grain=CLOCK_SUBCOMPONENT ) call mpp_clock_begin(clockid_init) !----- read namelist ----- if (file_exist('input.nml')) then unit = open_namelist_file ( ) ierr=1; do while (ierr /= 0) read (unit, nml=atmos_aeolus2_nml, iostat=io, end=10) ierr = check_nml_error (io, 'atmos_aeolus2_nml') enddo 10 call close_file (unit) endif !----- write version and namelist to log file ----- call write_version_number ( version, tagname ) if ( mpp_pe() == mpp_root_pe() ) then log_unit = stdlog() write (log_unit, nml=atmos_aeolus2_nml) write (*, nml=atmos_aeolus2_nml) endif ! some sanity checks if (0 == NoNx .or. 0 == NoNy) & write(*,*) 'Warning: strange grid size in atmos_aeolus2_nml:',NoNx,NoNy !---- compute physics/atmos time step in seconds ---- Time_step_atmos = Time_step ! remember Time_step call get_time (Time_step, Time_step_seconds) dt_atmos = real(Time_step_seconds) ! initialize domain decomposition ! heavily inspired by atmos_spectral/tools/spec_mpp.F90 call mpp_domains_init() pe = mpp_pe() npes = mpp_npes() !grid domain: atmos_fv and EBM do, by default, 1D decomposition along Y layout = (/1,npes/) !if( PRESENT(grid_layout) ) layout = grid_layout call mpp_define_domains( (/1,NoNx,1,NoNy/), layout, grid_domain ) if(pe == mpp_root_pe()) call print_decomp (npes, layout, grid_domain ) !!requirement of equal domains: can be generalized to retain mirror symmetry between N/S if unequal. !!the equal-domains requirement permits us to eliminate one buffer/unbuffer in the transpose_fourier routines. !if( mod(NoNy,layout(2)).NE.0 ) then ! ! call mpp_error( FATAL, 'atmos_aeolus2fms: currently requires equal grid domains on all PEs.' ) ! write(chtmp1,'(i4)') layout(2) ! write(chtmp2,'(i4)') NoNy ! call mpp_error( FATAL, 'atmos_aeolus2fms: Requires num_lat_rows/num_pes=int;num_pes='& ! &//chtmp1//';num_lat_rows='//chtmp2 ) !endif call mpp_get_compute_domain(grid_domain, is, ie, js, je) ! --- now we can allocate variables --- allocate(atmos_pelist(npes)) call mpp_get_current_pelist( atmos_pelist, commID=atmos_comm ) !write(*,*) 'current pelist ', atmos_pelist !write(*,'(a,z)') 'current atmos_comm ',atmos_comm deallocate(atmos_pelist) ! get the number of prognostic tracers call get_number_tracers( MODEL_ATMOS, num_prog=nr_tracers) q_ind = get_tracer_index( MODEL_ATMOS, 'sphum' ) if (q_ind==NO_TRACER) & q_ind = get_tracer_index( MODEL_ATMOS, 'mix_rat' ) if (q_ind==NO_TRACER) & call error_mesg('atmosphere_init', 'Cannot find water vapor in ATM tracer table', FATAL) co2_ind = get_tracer_index( MODEL_ATMOS, 'co2' ) allocate (Surf_diff%dtmass (is:ie, js:je) ) allocate (Surf_diff%dflux_t (is:ie, js:je) ) allocate (Surf_diff%dflux_tr (is:ie, js:je, nr_tracers) ) allocate (Surf_diff%delta_t (is:ie, js:je) ) allocate (Surf_diff%delta_tr (is:ie, js:je, nr_tracers) ) allocate (Surf_diff%delta_u (is:ie, js:je) ) allocate (Surf_diff%delta_v (is:ie, js:je) ) Surf_diff%dtmass = 0.0 Surf_diff%dflux_t = 0.0 Surf_diff%dflux_tr = 0.0 Surf_diff%delta_t = 0.0 Surf_diff%delta_tr = 0.0 Surf_diff%delta_u = 0.0 Surf_diff%delta_v = 0.0 ! --- read restart --- !filename = 'INPUT/atmosphere.res.nc' !if (file_exist(filename)) then ! !call read_data(filename, 'tg', tg, grid_domain) ! !call read_data(filename, 'qg', qg, grid_domain) !else ! ! tg-initial atmospheric temperature default value 273 deg Kelvin ! ! qg-initial atmospheric specific humidity ! !tg = tg_init ! !qg = rh*sat_mixing_ratio(tg) !endif ! possibly set up some more tracers ! ... call get_date(Time_init, Time_init_year, Time_init_month, Time_init_day, & Time_init_hours, Time_init_minutes, Time_init_seconds) call get_time(Time, Time_seconds, Time_days) call get_time(Time_step, Time_step_seconds, Time_step_days) call aeolus2fms_interface_check(atmos_comm) call aeolus2fms_init_grid( & NoNx, NoNy, & !NoNz, NoNz_Strato, & lon_0_360, & atmos_comm, pe, mpp_root_pe(), npes, & is, ie, js, je & ) !----- initialize atmos_axes and diagnostics ! If we want to use the FMS diagnostics manager infrastructure, ! we can set it up here. Even if we wanted to handle that from ! the python side only, the flux exchange module of the coupler wants to ! register some diagnostic fields using the atmosphere axes. ! Thus we really need to set those up here. ! The diagnostic axes require _global_ information allocate(lonb(1:NoNx+1), latb(1:NoNy+1), lon(1:NoNx), lat(1:NoNy)) !allocate(level(NoNz+NoNz_Strato), levelb(NoNz+NoNz_Strato+1)) ! Fill arrays with well-known values to ease debugging !do i=1,NoNx ! lonb(i) = i ! lon(i) = i !end do !lonb(NoNx+1) = NoNx+1 !do i=1,NoNy ! latb(i) = i ! lat(i) = i !end do !latb(NoNy+1) = NoNy+1 !do i=1,NoNz+NoNz_Strato ! level(i) = i ! levelb(i) = i !end do !levelb(NoNz+NoNz_Strato+1) = NoNz+NoNz_Strato+1 ! receive global(!) axis labels from the python side call aeolus2fms_get_axes_coords(lon, lonb, lat, latb) write (*,*) 'lon', lon write (*,*) 'lonb', lonb write (*,*) 'lat', lat write (*,*) 'latb', latb !write (*,*) 'level', level !write (*,*) 'levelb', levelb ! ... and also retrieve global(!) axis information from the exchange grid definition, ! to countercheck that python and FMS agree in their grid definitions. allocate(exglon(1:NoNx), exglat(1:NoNy), londiff(1:NoNx), latdiff(1:NoNy)) call get_grid_cell_centers('ATM', tile=1, glon=exglon, glat=exglat) londiff = exglon-lon if (any(abs(londiff) > 1e-13)) then write(*,*) 'Error: longitude axis info non-trivial mismatch between FMS and python' write(*,*) 'FMS', exglon write(*,*) 'python', lon write(*,*) 'FMS-python', exglon-lon endif latdiff = exglat-lat if (any(abs(latdiff) > 1e-13)) then write(*,*) 'Error: latitude axis info non-trivial mismatch between FMS and python' write(*,*) 'FMS', exglat write(*,*) 'python', lat write(*,*) 'FMS-python', exglat-lat endif deallocate(exglon, exglat, londiff, latdiff) id_lonb = diag_axis_init('lonb', lonb, 'degrees_E', 'x', 'longitude edges', set_name=axiset, Domain2=grid_domain) id_latb = diag_axis_init('latb', latb, 'degrees_N', 'y', 'latitude edges', set_name=axiset, Domain2=grid_domain) atmos_axes(1) = diag_axis_init('lon', lon, 'degrees_E', 'x', 'longitude', set_name=axiset, Domain2=grid_domain, edges=id_lonb) atmos_axes(2) = diag_axis_init('lat', lat, 'degrees_N', 'y', 'latitude', set_name=axiset, Domain2=grid_domain, edges=id_latb) !id_levelb = diag_axis_init('levelb', levelb, 'm', 'z', 'level edges', set_name=axiset) !atmos_axes(3) = diag_axis_init('level', level, 'm', 'z', 'level height', set_name=axiset, edges=id_levelb) ! init topography allocate(lon_boundaries(ie-is+2, je-js+2)) allocate(lat_boundaries(ie-is+2, je-js+2)) write(*,*) 'allocate lon_boundaries to sizes ', ie-is+2, je-js+2 call atmosphere_boundary(lon_boundaries, lat_boundaries, global=.false.) !write(*,*) 'lon_boundaries', lon_boundaries*180.0/PI ! output in degrees for better readability !write(*,*) 'lat_boundaries', lat_boundaries*180.0/PI previous_fpe_flags = FOR_SET_FPE(loose_fpe_flags) allocate(HORO(ie-is+1, je-js+1)) got_topo = get_topog_mean(lon_boundaries, lat_boundaries, HORO) if (.not. got_topo) & call error_mesg ('atmosphere_init', 'could not get topography mean', FATAL) allocate(SIGORO(ie-is+1, je-js+1)) got_topo = get_topog_stdev(lon_boundaries, lat_boundaries, SIGORO) if (.not. got_topo) & call error_mesg ('atmosphere_init', 'could not get topography stdev', FATAL) allocate(land_frac(ie-is+1, je-js+1)) ! This is not needed during init !where (HORO > 0) ! land_frac = 1.0 !elsewhere ! land_frac = 0.0 !end where allocate(area(ie-is+1, je-js+1)) call atmosphere_cell_area(area) if (THIRTY_DAY_MONTHS /= get_calendar_type() .and. NOLEAP /= get_calendar_type()) then call mpp_error(FATAL,'atmosphere_init: need integral number of days per year, i.e. calendar_type THIRTY_DAY_MONTHS or NOLEAP') endif DaysPerYear = days_in_year(Time_init) ! prepare for the use of astronomy_mod call astronomy_init() allocate(rad_lat(js:je)) !write(*,*)'size rad_lat ',size(rad_lat),' size lat ',size(lat) !write(*,*) 'bounds rad_lat ', lbound(rad_lat), ubound(rad_lat), ' rad ',lbound(lat),ubound(lat) rad_lat(:) = lat(js:je)*PI/180. !write(*,*)'size rad_lat ',size(rad_lat),' size lat ',size(lat) !write(*,*) 'bounds rad_lat ', lbound(rad_lat), ubound(rad_lat), ' rad ',lbound(lat),ubound(lat) if (fix_astro) then Time_fix_astro = set_date(astro_date(1), astro_date(2), astro_date(3), astro_date(4), astro_date(5), astro_date(6)) call get_time(Time_fix_astro, astro_date(6), fixdayoftheyear) fixdayoftheyear = MODULO(fixdayoftheyear, DaysPerYear)+1 end if call aeolus2fms_init( & Time_init_year, Time_init_month, Time_init_day, Time_init_seconds, & Time_days, Time_seconds, & Time_step_days, Time_step_seconds, & without_topo, & HORO, & SIGORO, & land_frac, & area, & update_land_frac_always, & nr_tracers, & q_ind, & co2_ind & ) ! here we could register some diagnostic fields axes_2d = atmos_axes(1:2) !id_astro_cosz = register_diag_field(mod_name, 'fms_astro_cosz' , (/atmos_axes(2)/), Time, long_name='cosz is cosine of an effective zenith angle that would produce the correct daily solar flux if the sun were fixed at that position for the period of daylight; from FMS astronomy_mod', units='1') !id_astro_solar = register_diag_field(mod_name, 'fms_astro_solar' , (/atmos_axes(2)/), Time, long_name='shortwave flux factor: cosine of zenith angle * daylight fraction / (earth-sun distance squared); from FMS astronomy_mod', units='1') !liq_wat = get_tracer_index(MODEL_ATMOS, 'liq_wat') !ice_wat = get_tracer_index(MODEL_ATMOS, 'ice_wat') id_horo = register_static_field(mod_name, "horo", axes_2d, "Orographic Height from FMS", "m") id_sigoro = register_static_field(mod_name, "sigoro", axes_2d, "Standard Deviation of Orographic Height from FMS", "m") if (id_horo > 0) used = send_data(id_horo, HORO, Time) if (id_sigoro > 0) used = send_data(id_sigoro, SIGORO, Time) id_area = register_static_field(mod_name, "area", axes_2d, "Cell Area", "m**2") if (id_area > 0) used = send_data(id_area, area, Time) !deallocate(lonb, latb, lon, lat) !deallocate(level, levelb) deallocate(HORO) deallocate(SIGORO) deallocate(land_frac) deallocate(area) id_flux_lw_sf_up = register_diag_field(mod_name, 'flux_lwr_sf_up', axes_2d, Time, long_name='lwr upward flux calculated in atmos wrapper', units='W/m2') call diag_manager_end(Time) ! to get static fields onto disk before something crashes !call register_Aeolus2_diag_vars(Time) atmos_stock_water_id = register_diag_field(mod_name, 'waterstock', Time, long_name = 'Total water content of atmosphere', units='') atmos_stock_heat_id = register_diag_field(mod_name, 'heatstock', Time, long_name = 'Total heat content of atmosphere', units='J') !call send_Aeolus2_static_diag_vars(Time) call mpp_clock_end(clockid_init) module_is_initialized = .true. !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler original_fpe_flags = FOR_SET_FPE(original_fpe_flags) write (*,*) 'atmosphere_init finished' end subroutine atmosphere_init ! ! ================================================================================== ! ! ! public routine required for atmospheric components of coupled models ! ! ! atmosphere_down ! This routine calls the dynamical core and the ! "downward pass" of the atmospheric physics. ! It should only be called once per time step and before ! calling atmosphere_up. ! ! ! ! time at the current time level (tau) ! ! ! fraction (0. to 1.) of underlying surface which covered by land ! ! ! surface (skin) temperature (in deg k) ! ! ! surface albedo ! ! ! momentum roughness length (units=m) ! ! ! friction velocity ! ! ! buoyancy scale ! ! ! moisture scale ! ! ! derivative of zonal wind stress w.r.t. the lowest level wind speed ! ! ! derivative of meridional wind stress w.r.t. the lowest level wind speed ! ! ! fraction of oopen sea without ice ! ! ! wind gustiness ! ! ! cosine of the zenith angle ! ! ! cosine of the zenith angle, radiation-weighted, Climber-2 style. ! Ignored with Aeolus2. ! ! ! net shortwave surface flux (down minus up) (in watts/m**2) ! ! ! downward longwave surface flux (in watts/m**2) ! ! ! Surface diffusion terms computed by the vertical diffusion scheme ! In the GFDL-provided example atmosphere implementations, all of the ! fields of this variable are asigned new values here. ! subroutine atmosphere_down ( & ! input parameters & Time, & frac_land, t_surf, albedo, & albedo_vis_dir, albedo_nir_dir, & albedo_vis_dif, albedo_nir_dif, & rough_mom, & u_star, b_star, q_star, & dtau_du, dtau_dv, & ! derivatives of wind stress w.r.t. the lowest level wind speed ! inout parameters & u_flux, v_flux, & ! wind stress frac_open_sea, & ! non-seaice fraction (%) ! output parameters & gust, coszen, coszen_radwt, & flux_sw, & flux_sw_dir, flux_sw_dif, & flux_sw_down_vis_dir, & flux_sw_down_vis_dif, & flux_sw_down_total_dir, & flux_sw_down_total_dif, & flux_sw_vis, & flux_sw_vis_dir, & flux_sw_vis_dif, & flux_lw, & Surf_diff ) type(time_type), intent(in) :: Time real, intent(in), dimension(is:ie,js:je) :: frac_land, t_surf, albedo real, intent(in), dimension(is:ie,js:je) :: albedo_vis_dir, albedo_nir_dir real, intent(in), dimension(is:ie,js:je) :: albedo_vis_dif, albedo_nir_dif real, intent(in), dimension(is:ie,js:je) :: rough_mom real, intent(in), dimension(is:ie,js:je) :: u_star, b_star, q_star, dtau_du, dtau_dv real, intent(inout), dimension(is:ie,js:je) :: u_flux, v_flux real, intent(in), dimension(is:ie,js:je) :: frac_open_sea real, intent(out), dimension(is:ie,js:je) :: gust, coszen real, intent(inout), dimension(is:ie,js:je) :: coszen_radwt ! declare as inout to suppress compiler warning real, intent(out), dimension(is:ie,js:je) :: flux_sw real, intent(out), dimension(is:ie,js:je) :: flux_sw_dir, flux_sw_dif, flux_sw_down_vis_dir real, intent(out), dimension(is:ie,js:je) :: flux_sw_down_vis_dif, flux_sw_down_total_dir real, intent(out), dimension(is:ie,js:je) :: flux_sw_down_total_dif, flux_sw_vis real, intent(out), dimension(is:ie,js:je) :: flux_sw_vis_dir, flux_sw_vis_dif real, intent(out), dimension(is:ie,js:je) :: flux_lw type(surf_diff_type), intent(inout) :: Surf_diff integer :: Time_seconds, Time_days integer :: dayoftheyear real, dimension(js:je) :: cosz, solar, cosz_radwt real, dimension(is:ie,js:je) :: flux_lw_sf_up integer :: j, used call mpp_clock_begin (clockid_down) !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler call clearstatusfpqq() ! do not inherit signalled FP exceptions from FMS original_fpe_flags = FOR_SET_FPE(strict_fpe_flags) call get_time(Time, Time_seconds, Time_days) write (*,*) 'atmosphere_down starting at', Time_days, Time_seconds !type(time_type) :: Time_prev, Time_next ! !Time_prev = Time ! two time-level scheme !Time_next = Time + Time_step_atmos !---- dynamics ----- !call timing_on('fv_dynamics') !call mpp_clock_begin (id_dynam) !call mpp_clock_end (id_dynam) !call timing_off('fv_dynamics') !---- call physics ----- !call timing_on('fv_physics_down') !call mpp_clock_begin (id_phys_down) !call mpp_clock_end (id_phys_down) !call timing_off('fv_physics_down') !call fv_print_chksums( 'Exiting atmosphere_down' ) ! TODO: possible optimization: only re-calculate solar data once per day if (fix_astro) then dayoftheyear = fixdayoftheyear else dayoftheyear = MODULO(Time_days, DaysPerYear)+1 end if !write(*,*)'size rad_lat ',size(rad_lat),' size cosz ',size(cosz),' size solar ', size(solar) if (id_astro_cosz > 0 .or. id_astro_solar > 0) then if (fix_astro) then call daily_mean_solar (rad_lat, Time_fix_astro, cosz, solar) else !write(*,*)'size rad_lat ',size(rad_lat) !write(*,*) 'bounds rad_lat ', lbound(rad_lat), ubound(rad_lat) !write(*,*)'size cosz ',size(cosz) !write(*,*) 'bounds cosz ', lbound(cosz), ubound(cosz) call daily_mean_solar (rad_lat, Time, cosz, solar) end if if (id_astro_cosz > 0) used = send_data(id_astro_cosz, cosz, Time) if (id_astro_solar > 0) used = send_data(id_astro_solar, solar, Time) end if ! the return value solar from FMS-astronomy must be multiplied by solar constant to ! obtain the same meaning as solarm from sinsol solar = solar * solar_constant do j = js,je coszen(:,j) = cosz(j) ! pass back to coupler the non-radiation-weighted coszen enddo !cosz_radwt = 0.0 was already set during initialisation, and is not used with Aeolus2 !write(*,*) 'aeolus2fms_atmosphere_down()', 'delta_t ', size(Surf_diff%delta_t), & ! lbound(Surf_diff%delta_t), ubound(Surf_diff%delta_t), & ! 'delta_tr', size(Surf_diff%delta_tr), & ! lbound(Surf_diff%delta_tr), ubound(Surf_diff%delta_tr), & ! 'dflux_tr', size(Surf_diff%dflux_tr), & ! lbound(Surf_diff%dflux_tr), ubound(Surf_diff%dflux_tr) ! same as in surface_flux(), before interpolation to reduced grid flux_lw_sf_up = stefan*t_surf**4 if (id_flux_lw_sf_up > 0) used = send_data(id_flux_lw_sf_up, flux_lw_sf_up, Time) call aeolus2fms_atmosphere_down ( & ! input parameters & Time_seconds, Time_days, & dayoftheyear, & frac_land, t_surf, albedo, & albedo_vis_dir, albedo_nir_dir, & albedo_vis_dif, albedo_nir_dif, & rough_mom, & u_star, b_star, q_star, & !frac_open_sea, & ! non-seaice fraction (%) dtau_du, dtau_dv, & ! derivatives of wind stress w.r.t. the lowest level wind speed ! inout parameters & u_flux, v_flux, & ! wind stress ! output parameters & gust, & ! coszen, cosz, solar were already calculated above, and thus are now input-parameters for Aeolus coszen, & !cosz_radwt, & solar, & ! here pass the 1-D variables. flux_lw_sf_up, & ! output parameters flux_sw, & flux_sw_dir, flux_sw_dif, & flux_sw_down_vis_dir, & flux_sw_down_vis_dif, & flux_sw_down_total_dir, & flux_sw_down_total_dif, & flux_sw_vis, & flux_sw_vis_dir, & flux_sw_vis_dif, & flux_lw, & Surf_diff%delta_tr, & Surf_diff%dflux_tr, & Surf_diff%delta_t, & Surf_diff%dflux_t, & Surf_diff%dtmass, & Surf_diff%delta_u, & Surf_diff%delta_v & ) !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler original_fpe_flags = FOR_SET_FPE(original_fpe_flags) call mpp_clock_end(clockid_down) end subroutine atmosphere_down ! ================================================================================== ! ! ! public routine required for atmospheric components of coupled models ! ! ! ! atmosphere_up ! This routine calls the "upward pass" of the atmospheric physics, ! spherical diagnostics, and time differencing. The prognostic ! variables are advanced to the next time step. It should only be ! called once per time step and after calling atmosphere_down. ! ! ! ! time at the current time level (tau) ! ! ! fraction (0. to 1.) of underlying surface which covered by land ! ! ! liquid precipitation rate (rain) in kg/m2/s ! ! ! frozen precipitation rate (snow) in kg/m2/s ! ! ! wind gustiness ! ! ! surface diffusion terms computed by the vertical diffusion scheme ! ! ! friction velocity ! ! ! buoyancy scale ! ! ! moisture scale ! subroutine atmosphere_up (Time, frac_land, Surf_diff, lprec, fprec, gust, u_star, b_star, q_star ) type(time_type), intent(in) :: Time real, intent(in), dimension(is:ie,js:je) :: frac_land type(surf_diff_type), intent(in) :: Surf_diff real, intent(out), dimension(is:ie,js:je) :: lprec, fprec, gust ! u_star, b_star, q_star are not used in EBM atmosphere_up real, intent(in), dimension(is:ie,js:je) :: u_star, b_star, q_star ! Surf_diff%delta_t and Surf_diff%delta_tr have just been updated ! in the caller update_atmos_model_up(), and can be used now. ! None of the GFDL-provided atmosphere example models ! uses any of the other Surf_diff fields in the upward routines. integer :: Time_year, Time_month ! integer :: Time_seconds, Time_days ! seconds-of-day and total days since start of time axis, not day-of-year! integer :: used type(time_type) :: Time_next ! to get time stamp of diag output right real :: waterstock, heatstock call mpp_clock_begin(clockid_up) !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler call clearstatusfpqq() ! do not inherit signalled FP exceptions from FMS original_fpe_flags = FOR_SET_FPE(strict_fpe_flags) ! here we need only year and month call get_date(Time, Time_year, Time_month, used, used, used, used) call get_time(Time, Time_seconds, Time_days) !type(time_type) :: !Time_prev, Time_next ! !Time_prev = Time ! two time-level scheme !Time_next = Time + Time_step_atmos !call timing_on('fv_physics_up') !call mpp_clock_begin (id_phys_up) !call mpp_clock_end (id_phys_up) !call timing_off('fv_physics_up') !---- diagnostics for FV dynamics ----- !call timing_on('FV_DIAG') !call mpp_clock_begin(id_fv_diag) !fv_time = Time + Time_step_atmos !call fv_diag(fv_time, nlon, mlat, nlev, beglat, endlat, & ! ncnst, zvir, dt_atmos, .false.) !call timing_off('FV_DIAG') !call mpp_clock_end(id_fv_diag) call aeolus2fms_atmosphere_up ( & ! input parameters Time_year, Time_month, & Time_seconds, Time_days, & ! seconds-of-day and days-since-time-axis-start frac_land, & Surf_diff%delta_t, & Surf_diff%delta_tr, & ! output parameters lprec, fprec, gust, & ! input parameters again u_star, b_star, q_star & ) Time_next = Time + Time_step_atmos !call send_Aeolus2fms_diag_vars(Time_next) if (atmos_stock_water_id > 0) then call get_stock_atm_global(ISTOCK_WATER, waterstock) used = send_data(atmos_stock_water_id, waterstock, Time) end if if (atmos_stock_heat_id > 0) then call get_stock_atm_global(ISTOCK_HEAT, heatstock) used = send_data(atmos_stock_heat_id, heatstock, Time) end if !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler original_fpe_flags = FOR_SET_FPE(original_fpe_flags) call mpp_clock_end(clockid_up) end subroutine atmosphere_up ! ! ===================================================================== ! ! ! ! write out restart file ! Termination routine for atmosphere_mod. ! ! ! public routine required for atmospheric components of coupled models ! write out restart file ! Termination routine for atmosphere_mod. ! ! ! ! ! ! time at the current time level (tau) ! ! ! ! Contains information about the grid cells ! subroutine atmosphere_end (Time, Grid_box) type (time_type), intent(in) :: Time type(grid_box_type), intent(inout) :: Grid_box ! unused integer :: Time_days, Time_seconds if (.NOT. module_is_initialized) return !----- initialize domains for writing global physics data ----- !call set_domain ( fv_domain ) ! Note that in atmos_coupled/atmos_model.F90:atmos_model_local_restart() ! already the fields lprec,fprec,gust,t_bot,q_bot are saved, and restored in ! atmos_model_init() call get_time(Time, Time_seconds, Time_days) call aeolus2fms_finish(Time_seconds, Time_days) ! de-allocate memory ! sigh. the surf_diff_type cannot de-allocated in any atmosphere implementation, ! because it is never passed to any end-subroutine !deallocate (Surf_diff%dtmass) !deallocate (Surf_diff%dflux_t) !deallocate (Surf_diff%dflux_tr) !deallocate (Surf_diff%delta_t) !deallocate (Surf_diff%delta_tr) !deallocate (Surf_diff%delta_u) !deallocate (Surf_diff%delta_v) deallocate(lonb, latb, lon, lat) deallocate(rad_lat) !deallocate(Aeolus2_diag_ids) !deallocate(Aeolus2_diag_types) module_is_initialized = .false. end subroutine atmosphere_end ! !####################################################################### ! ! ! Write out restart data ! subroutine atmosphere_restart(timestamp) ! format yyyymmdd.hhmmss, from date_to_string() character(len=*), intent(in) :: timestamp ! sigh deeply. For the time stamp in the aeolus restart file, we want ! days since 0000-01-01. But all we get is this timestamp string, ! which even cannot be converted back into a time_type, because ! the FMS-provided set_date() expects a different format than ! what comes out of date_to_string(). ! On the other hand, set_date() already has all the calendar-handling built in, ! so we really want to use that and not re-invent the wheel. ! First apporach: lets try some string manipulation in Fortran... ! Convert timestamp to the format 1980-01-01 00:00:00 type(time_type) :: Time integer :: Time_days, Time_seconds character(len=50) :: timebuf timebuf = timestamp(1:4)//'-'//timestamp(5:6)//'-'//timestamp(7:8)& //' '//timestamp(10:11)//':'//timestamp(12:13)//':'//timestamp(14:15) write(*,*) 'atmosphere_restart: timebuf>', timebuf, '<' ! TODO: if necessary, save things into atmosphere.res.nc Time = set_date(timebuf) call get_time(Time, Time_seconds, Time_days) ! Note that in atmos_coupled/atmos_model.F90:atmos_model_local_restart() ! already the fields lprec,fprec,gust,t_bot,q_bot are saved, and restored in ! atmos_model_init() call aeolus2fms_restart(Time_seconds, Time_days) end subroutine atmosphere_restart ! ! ================================================================================== ! ! ! returns the number of longitude and latitude grid points ! for either the local PEs grid (default) or the global grid ! ! ! public routine required for atmospheric components of coupled models ! returns the number of longitude and latitude grid points ! for either the local PEs grid (default) or the global grid ! ! ! ! The number of longitude points in the compute domain. ! ! ! The number of latitude points in the compute domain. ! ! ! Flag that specifies whether the returned compute domain size is ! for the global grid (TRUE) or for the current processor (FALSE). ! subroutine atmosphere_resolution (num_lon, num_lat, global) integer, intent(out) :: num_lon, num_lat logical, intent(in), optional :: global logical :: local local = .TRUE. if (PRESENT(global)) local = .NOT.global if (local) then num_lon = ie - is + 1 num_lat = je - js + 1 else num_lon = NoNx num_lat = NoNy end if end subroutine atmosphere_resolution ! ! ================================================================================== ! ! ! ! returns grid cell areas for the domain ! ! ! ! returns grid cell areas for the domain ! subroutine atmosphere_cell_area (area_out) real, dimension(:,:), intent(out) :: area_out integer :: xsize, ysize xsize = ie-is+1 ysize = je-js+1 if(any(shape(area_out) /= (/xsize,ysize/))) then call mpp_error(FATAL,'atmosphere_cell_area: argument has wrong shape. Its shape is ', & shape(area_out),' It should be ',(/xsize,ysize/)) endif call get_grid_cell_area('ATM', tile=1, cellarea=area_out, domain=grid_domain) !do xsize=is,ie ! do ysize=js,je ! write(*,*) 'i ', xsize, ' j ', ysize, ' area ', area_out(xsize, ysize) ! end do !end do end subroutine atmosphere_cell_area ! ! ================================================================================== ! ! ! returns the longitude and latitude grid box edges ! for either the local PEs grid (default) or the global grid ! ! ! public routine required for atmospheric components of coupled models ! returns the longitude and latitude grid box edges ! for either the local PEs grid (default) or the global grid ! ! ! ! The west-to-east longitude edges of grid boxes (in radians). ! ! ! The south-to-north latitude edges of grid boxes (in radians). ! ! ! Flag that specifies whether the returned grid box edges are ! for the global grid (TRUE) or for the current processor (FALSE). ! subroutine atmosphere_boundary (lon_boundaries, lat_boundaries, global) real, intent(out) :: lon_boundaries(:,:), lat_boundaries(:,:) ! radians !! logical, intent(in), optional :: global logical :: local integer :: i, i1, i2, j1, j2 i1=is; i2=ie+1; j1=js; j2=je+1 local = .TRUE. if( PRESENT(global) )local = .NOT.global if (local) then write(*,*) 'atmosphere_boundary local' else write(*,*) 'atmosphere_boundary global' i1=1; i2=size(lonb); j1=1; j2=size(latb) end if write(*,*) 'shape of lon_boundaries', shape(lon_boundaries) write(*,*) 'shape of lonb', shape(lonb) call flush(6) do i = 1, size(lon_boundaries,2) lon_boundaries(:,i) = lonb(i1:i2) end do do i = 1, size(lon_boundaries,1) lat_boundaries(i,:) = latb(j1:j2) end do lon_boundaries = lon_boundaries/180.0*PI lat_boundaries = lat_boundaries/180.0*PI end subroutine atmosphere_boundary ! ! ================================================================================== ! ! ! ! returns the domain2D variable associated with the coupling grid ! ! ! ! public routine required for atmospheric components of coupled models ! returns the domain2D variable associated with the coupling grid ! note: coupling is done using the mass/temperature grid with no halos ! !OUTPUT ! Domain The domain2D variable describing the grid used for coupling. ! For the B-grid, this corresponds to the temperature grid ! without halos. ! ! ! ! The domain2D variable describing the grid used for coupling. ! For the B-grid, this corresponds to the temperature grid ! without halos. ! subroutine atmosphere_domain (domain) type(domain2D), intent(out) :: domain domain = grid_domain end subroutine atmosphere_domain ! ! ================================================================================== ! ! ! ! returns the axis indices associated with the coupling grid ! ! ! public routine required for atmospheric components of coupled models ! returns the axis indices associated with the coupling grid ! ! ! ! The axis identifiers for the atmospheric grids. ! The size of axes must be least 1 but not greater than 4. ! The axes are returned in the order (/ x, y, p_full, p_half /) ! subroutine get_atmosphere_axes ( axes ) integer, intent(out), dimension(:) :: axes !----- returns the axis indices for the atmospheric (mass) grid ----- if ( size(axes(:)) < 1 .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 ! ! ================================================================================== ! ! ! ! returns temp, sphum, pres, height at the lowest model level ! and surface pressure ! ! ! public routine required for atmospheric components of coupled models ! returns temp, sphum, pres, height at the lowest model level ! and surface pressure ! ! ! ! near surface temperature in degrees Kelvin ! ! ! array of tracers, 1st tracer is ! near surface mixing ratio ! ! ! pressure at which atmos near surface values are assumed to be defined ! ! ! height at which atmos near surface values are assumed to be defined ! ! ! surface pressure ! ! ! sea level pressure ! subroutine get_bottom_mass (t_bot, tr_bot, p_bot, z_bot, p_surf, slp) real, intent(out), dimension(is:ie,js:je) :: t_bot, p_bot, z_bot, p_surf, slp real, intent(out), dimension(:,:,:) :: tr_bot !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler call clearstatusfpqq() ! do not inherit signalled FP exceptions from FMS original_fpe_flags = FOR_SET_FPE(strict_fpe_flags) ! for now, just set some dummy values which prevent the coupler from ! division by zero !t_bot = 274.0 !tr_bot = !p_bot = 1 !z_bot = 1 !p_surf = 1 !slp = 1 call aeolus2fms_get_bottom_mass(t_bot, tr_bot, p_bot, z_bot, p_surf, slp) !! To enable strict floating point exceptions inside this module. !! Taken from the example in the Intel Fortran Compiler Reference Manual. !! Most probably will not work with any other compiler original_fpe_flags = FOR_SET_FPE(original_fpe_flags) end subroutine get_bottom_mass ! ================================================================================== ! ! ! ! returns u and v on the mass grid at the lowest model level ! ! ! public routine required for atmospheric components of coupled models ! returns u and v on the mass grid at the lowest model level ! ! ! ! near surface zonal wind ! ! ! near surface meridional wind ! subroutine get_bottom_wind (u_bot, v_bot) real, intent(out), dimension(is:ie,js:je) :: u_bot, v_bot call aeolus2fms_get_bottom_wind(u_bot, v_bot) end subroutine get_bottom_wind ! ===================================================================== ! ! ! ! Return for this pe the total content of a kind of stock (water, salt, heat) in the atmosphere. ! ! ! Return for this pe the total content of a kind of stock (water, salt, heat) in the atmosphere. ! ! ! ! index to indicate wich kind of stock to get ! ! ! value of requested stock ! subroutine get_stock_pe(idx, val) integer, intent(in) :: idx real, intent(out) :: val call aeolus2fms_get_stock_pe(idx, val) end subroutine get_stock_pe ! ! ===================================================================== ! ! ! ! Return global total content of a kind of stock (water, salt, heat) in the atmosphere. ! ! ! Return global total content of a kind of stock (water, salt, heat) in the atmosphere. ! ! ! ! index to indicate wich kind of stock to get ! ! ! value of requested stock ! subroutine get_stock_atm_global(idx, val) integer, intent(in) :: idx real, intent(out) :: val call aeolus2fms_get_stock_pe(idx, val) call mpp_sum(val) end subroutine get_stock_atm_global ! ! This private subroutine is simply copied from atmos_spectral/tools/spec_mpp.F90 ! because it is not accesible in there from here subroutine print_decomp (npes, layout, Domain) integer, intent(in) :: npes, layout(2) type(domain2D), intent(in) :: Domain integer, dimension(0:npes-1) :: xsize, ysize integer :: i, j, xlist(layout(1)), ylist(layout(2)) type (domain1D) :: Xdom, Ydom call mpp_get_compute_domains ( Domain, xsize=xsize, ysize=ysize ) call mpp_get_domain_components ( Domain, Xdom, Ydom ) call mpp_get_pelist ( Xdom, xlist ) call mpp_get_pelist ( Ydom, ylist ) write (*,100) write (*,110) (xsize(xlist(i)),i=1,layout(1)) write (*,120) (ysize(ylist(j)),j=1,layout(2)) 100 format ('ATMOS MODEL DOMAIN DECOMPOSITION') 110 format (' X-AXIS = ',24i4,/,(11x,24i4)) 120 format (' Y-AXIS = ',24i4,/,(11x,24i4)) end subroutine print_decomp end module atmosphere_mod ! And now one routine that is called back from the python side. ! It is placed outside the module, to make cross-language-calling a little more portable. ! When the python code exits or aborts, the FMS diag_manager buffers are lost. ! Here, we try to flush them before submitting to the crash. subroutine abort_with_flush use diag_manager_mod, only: diag_manager_end use time_manager_mod, only: time_type, set_time implicit none type(time_type) :: t t = set_time(0, 999999) call diag_manager_end(t) return end subroutine abort_with_flush