diff --git a/model/fv_control.F90 b/model/fv_control.F90 index fa377579e..45e89576e 100644 --- a/model/fv_control.F90 +++ b/model/fv_control.F90 @@ -996,7 +996,7 @@ subroutine read_namelist_fv_grid_nml ! Read Main namelist read (f_unit,fv_grid_nml,iostat=ios) ierr = check_nml_error(ios,'fv_grid_nml') - rewind (f_unit) + call close_file (f_unit) #endif call write_version_number ( 'FV_CONTROL_MOD', version ) unit = stdlog() diff --git a/model/fv_regional_bc.F90 b/model/fv_regional_bc.F90 index 5c0bd7972..98285b7a3 100644 --- a/model/fv_regional_bc.F90 +++ b/model/fv_regional_bc.F90 @@ -112,15 +112,19 @@ module fv_regional_mod ! integer, parameter :: iend_nest = 1346 ! integer, parameter :: jend_nest = 1290 - integer,parameter :: nvars_core=7 & !<-- # of prognostic variables in restart file - ,nvars_tracers=8 !<-- # of prognostic variables in restart file - + integer,parameter :: nvars_core=7 & !<-- # of prognostic variables in core restart file + ,ndims_core=6 & !<-- # of core restart dimensions + ,ndims_tracers=4 !<-- # of tracer restart dimensions +! real,parameter :: blend_exp1=0.5,blend_exp2=10. !<-- Define the exponential dropoff of weights ! for prescribed external values in the ! blending rows inside the domain boundary. real :: current_time_in_seconds ! - integer,save :: ncid,next_time_to_read_bcs,npz,ntracers + integer,save :: isd_mod,ied_mod,jsd_mod,jed_mod +! + integer,save :: ncid,next_time_to_read_bcs,nfields_tracers & + ,npz,ntracers ! integer,save :: k_split,n_split ! @@ -153,6 +157,8 @@ module fv_regional_mod logical,save :: north_bc,south_bc,east_bc,west_bc & ,begin_regional_restart=.true. + logical,dimension(:),allocatable,save :: blend_this_tracer + character(len=50) :: filename_core='INPUT/fv_core.res.temp.nc' character(len=50) :: filename_core_new='RESTART/fv_core.res.tile1_new.nc' character(len=50) :: filename_tracers='INPUT/fv_tracer.res.temp.nc' @@ -329,6 +335,11 @@ subroutine setup_regional_BC(Atm, dt_atmos & west_bc =.false. ! ! write(0,*)' enter setup_regional_BC isd=',isd,' ied=',ied,' jsd=',jsd,' jed=',jed + isd_mod=isd + ied_mod=ied + jsd_mod=jsd + jed_mod=jed +! !----------------------------------------------------------------------- !*** Which side(s) of the domain does this task lie on if any? !----------------------------------------------------------------------- @@ -412,20 +423,10 @@ subroutine setup_regional_BC(Atm, dt_atmos & !*** Compute the index limits within the boundary region on each !*** side of the domain for both scalars and winds. Since the !*** domain does not move then the computations need to be done -!*** only once. Likewise find and save the locations of the -!*** available tracers in the tracers array. +!*** only once. !----------------------------------------------------------------------- ! call compute_regional_bc_indices(Atm%regional_bc_bounds) -! - sphum_index =get_tracer_index(MODEL_ATMOS, 'sphum') - liq_water_index =get_tracer_index(MODEL_ATMOS, 'liq_wat') - ice_water_index =get_tracer_index(MODEL_ATMOS, 'ice_wat') - rain_water_index=get_tracer_index(MODEL_ATMOS, 'rainwat') - snow_water_index=get_tracer_index(MODEL_ATMOS, 'snowwat') - graupel_index =get_tracer_index(MODEL_ATMOS, 'graupel') - cld_amt_index =get_tracer_index(MODEL_ATMOS, 'cld_amt') - o3mr_index =get_tracer_index(MODEL_ATMOS, 'o3mr') ! !----------------------------------------------------------------------- ! @@ -435,7 +436,8 @@ subroutine setup_regional_BC(Atm, dt_atmos & ! !----------------------------------------------------------------------- ! - ntracers=Atm%ncnst - Atm%flagstruct%dnats !<-- # of advected tracers +! ntracers=Atm%ncnst - Atm%flagstruct%dnats !<-- # of advected tracers + ntracers=Atm%ncnst !<-- Total # of tracers npz=Atm%npz !<-- # of layers in vertical configuration of integration klev_out=npz ! @@ -733,6 +735,24 @@ subroutine setup_regional_BC(Atm, dt_atmos & enddo enddo endif +! + sphum_index = get_tracer_index(MODEL_ATMOS, 'sphum') + liq_water_index = get_tracer_index(MODEL_ATMOS, 'liq_wat') + ice_water_index = get_tracer_index(MODEL_ATMOS, 'ice_wat') + rain_water_index = get_tracer_index(MODEL_ATMOS, 'rainwat') + snow_water_index = get_tracer_index(MODEL_ATMOS, 'snowwat') + graupel_index = get_tracer_index(MODEL_ATMOS, 'graupel') + cld_amt_index = get_tracer_index(MODEL_ATMOS, 'cld_amt') + o3mr_index = get_tracer_index(MODEL_ATMOS, 'o3mr') +! write(0,*)' setup_regional_bc' +! write(0,*)' sphum_index=',sphum_index +! write(0,*)' liq_water_index=',liq_water_index +! write(0,*)' ice_water_index=',ice_water_index +! write(0,*)' rain_water_index=',rain_water_index +! write(0,*)' snow_water_index=',snow_water_index +! write(0,*)' graupel_index=',graupel_index +! write(0,*)' cld_amt_index=',cld_amt_index +! write(0,*)' o3mr_index=',o3mr_index ! !----------------------------------------------------------------------- !*** When nudging of specific humidity is selected then we need a @@ -1067,8 +1087,6 @@ subroutine read_regional_lon_lat i_start_data=2*(isd+nhalo_model)-1 j_start_data=2*(jsd+nhalo_model)-1 ! -! write(0,11110)i_start_data,j_start_data -11110 format(' i_start_data=',i5,' j_start_data=',i5) !--------------- !*** Longitude !--------------- @@ -1248,7 +1266,7 @@ subroutine start_regional_cold_start(Atm, dt_atmos, ak, bk, levp & ,isd, ied, jsd, jed & ,ak, bk ) call regional_bc_t1_to_t0(BC_t1, BC_t0 & ! - ,Atm%npz & !<-- Move BC t1 data + ,Atm%npz & !<-- Move BC t1 data to t0. ,ntracers & ,Atm%regional_bc_bounds ) ! ! @@ -1277,13 +1295,14 @@ subroutine start_regional_cold_start(Atm, dt_atmos, ak, bk, levp & enddo ! !----------------------------------------------------------------------- -!*** If the GSI will need a restart file that includes the -!*** fields' boundary rows then create that file and define -!*** its dimensions and variables. +!*** If the GSI will need restart files that includes the +!*** fields' boundary rows. Those files were already created. +!*** Prepare the objects that hold their variables' names and +!*** values. !----------------------------------------------------------------------- ! if(Atm%flagstruct%write_restart_with_bcs)then - call create_restart_with_bcs(Atm) + call prepare_full_fields(Atm) endif ! !----------------------------------------------------------------------- @@ -1389,14 +1408,14 @@ subroutine start_regional_restart(Atm, dt_atmos & endif ! !----------------------------------------------------------------------- -!*** If the GSI will need a restart file that includes the +!*** If the GSI will need restart files that include the !*** fields' boundary rows after this forecast or forecast -!*** segment completes then create that file and define -!*** its dimensions and variables. +!*** segment completes then prepare the objects that will +!*** hold their variables' names and values. !----------------------------------------------------------------------- ! if(Atm%flagstruct%write_restart_with_bcs)then - call create_restart_with_bcs(Atm) + call prepare_full_fields(Atm) endif ! !----------------------------------------------------------------------- @@ -2356,8 +2375,7 @@ subroutine regional_bc_data(Atm,bc_hour & !*** FV3's modified virtual potential temperature. !----------------------------------------------------------------------- ! - call convert_to_virt_pot_temp(isd,ied,jsd,jed,npz & - ,sphum_index,liq_water_index ) + call convert_to_virt_pot_temp(isd,ied,jsd,jed,npz) ! !----------------------------------------------------------------------- !*** If nudging of the specific humidity has been selected then @@ -2924,17 +2942,14 @@ subroutine compute_cappa(i1,i2,j1,j2,cappa,temp,liq_wat,sphum) implicit none !----------------------------------------------------------------------- ! -!--------------------- -!*** Input variables -!--------------------- +!------------------------ +!*** Argument variables +!------------------------ ! integer,intent(in) :: i1,i2,j1,j2 ! - real,dimension(i1:i2,j1:j2,1:npz) :: cappa,temp,liq_wat,sphum -! -!---------------------- -!*** Output variables -!---------------------- + real,dimension(i1:i2,j1:j2,1:npz),intent(in) :: temp,liq_wat,sphum + real,dimension(i1:i2,j1:j2,1:npz),intent(inout) :: cappa ! !--------------------- !*** Local variables @@ -2967,7 +2982,7 @@ subroutine compute_cappa(i1,i2,j1,j2,cappa,temp,liq_wat,sphum) ql=qd-qs qv=max(0.,sphum(i,j,k)) cvm=(1.-(qv+qd))*cv_air + qv*cv_vap + ql*c_liq + qs*c_ice - ! +! cappa(i,j,k)=rdgas/(rdgas+cvm/(1.+zvir*sphum(i,j,k))) ! enddo @@ -3044,7 +3059,7 @@ subroutine read_regional_bc_file(is_input,ie_input & ! character(len=80) :: var_name !<-- Variable name in the boundary NetCDF file ! - logical :: call_get_var + logical :: call_get_var,is_root_pe logical :: required_local ! !----------------------------------------------------------------------- @@ -3060,6 +3075,8 @@ subroutine read_regional_bc_file(is_input,ie_input & else required_local=.true. endif +! + is_root_pe=(mpp_pe()==mpp_root_pe()) ! !----------------------------------------------------------------------- !*** Loop through the four sides of the domain. @@ -3227,13 +3244,16 @@ subroutine read_regional_bc_file(is_input,ie_input & ! if(call_get_var)then if (present(array_4d)) then !<-- 4-D variable - status=nf90_inq_varid(ncid,trim(var_name),var_id) !<-- Get this variable's ID. + status=nf90_inq_varid(ncid,trim(var_name),var_id) !<-- Get this variable's ID. if (required_local) then call check(status) endif if (status /= nf90_noerr) then if (east_bc.and.is_master()) write(0,*)' WARNING: Tracer ',trim(var_name),' not in input file' array_4d(:,:,:,tlev)=0. !<-- Tracer not in input so set to zero in boundary. +! + blend_this_tracer(tlev)=.false. !<-- Tracer not in input so do not apply blending. +! else call check(nf90_get_var(ncid,var_id & ,array_4d(i_start_array:i_end_array & !<-- Fill this task's domain boundary halo. @@ -3325,6 +3345,11 @@ subroutine allocate_regional_BC_arrays(side & allocate(BC_side%divgd_BC(is_0:ie_0,js_0:je_0,klev)) ; BC_side%divgd_BC=real_snan ! allocate(BC_side%q_BC (is_0:ie_0,js_0:je_0,1:klev,1:ntracers)) ; BC_side%q_BC=real_snan +! + if(.not.allocated(blend_this_tracer))then + allocate(blend_this_tracer(1:ntracers)) + blend_this_tracer=.true. !<-- Start with blending all tracers. + endif ! #ifndef SW_DYNAMICS allocate(BC_side%pt_BC (is_0:ie_0,js_0:je_0,klev)) ; BC_side%pt_BC=real_snan @@ -3415,14 +3440,14 @@ subroutine remap_scalar_nggps_regional_bc(Atm & ! !--------------------------------------------------------------------------------- ! - sphum = get_tracer_index(MODEL_ATMOS, 'sphum') - liq_wat = get_tracer_index(MODEL_ATMOS, 'liq_wat') - ice_wat = get_tracer_index(MODEL_ATMOS, 'ice_wat') - rainwat = get_tracer_index(MODEL_ATMOS, 'rainwat') - snowwat = get_tracer_index(MODEL_ATMOS, 'snowwat') - graupel = get_tracer_index(MODEL_ATMOS, 'graupel') - cld_amt = get_tracer_index(MODEL_ATMOS, 'cld_amt') - o3mr = get_tracer_index(MODEL_ATMOS, 'o3mr') + sphum = sphum_index + liq_wat = liq_water_index + ice_wat = ice_water_index + rainwat = rain_water_index + snowwat = snow_water_index + graupel = graupel_index + cld_amt = cld_amt_index + o3mr = o3mr_index k2 = max(10, km/2) @@ -3587,7 +3612,7 @@ subroutine remap_scalar_nggps_regional_bc(Atm & BC_side%q_BC(i,j,k,iq) = qn1(i,k) enddo enddo - endif ! skip cld_amt in the remap since it is not included in the input + endif enddo !--------------------------------------------------- @@ -4272,7 +4297,7 @@ subroutine regional_boundary_update(array & ! real,dimension(:,:,:),pointer :: bc_t0,bc_t1 !<-- Boundary data at the two bracketing times. ! - logical :: call_interp + logical :: blend,call_interp ! !--------------------------------------------------------------------- !********************************************************************* @@ -4282,9 +4307,11 @@ subroutine regional_boundary_update(array & return endif ! + blend=.true. iq=0 if(present(index4))then iq=index4 + blend=blend_this_tracer(iq) endif ! !--------------------------------------------------------------------- @@ -4486,7 +4513,7 @@ subroutine regional_boundary_update(array & ,fcst_time & ,bc_update_interval & ,i1_blend,i2_blend,j1_blend,j2_blend & - ,i_bc,j_bc, nside, bc_vbl_name ) + ,i_bc,j_bc,nside,bc_vbl_name,blend ) endif ! !--------------------------------------------------------------------- @@ -4616,7 +4643,7 @@ subroutine bc_time_interpolation(array & ,fcst_time & ,bc_update_interval & ,i1_blend,i2_blend,j1_blend,j2_blend & - ,i_bc,j_bc,nside, bc_vbl_name ) + ,i_bc,j_bc,nside,bc_vbl_name,blend ) !--------------------------------------------------------------------- !*** Update the boundary region of the input array at the given @@ -4645,12 +4672,12 @@ subroutine bc_time_interpolation(array & ! real,intent(in) :: fcst_time !<-- Current forecast time (sec) ! - real :: rdenom -! - real,dimension(lbnd1:ubnd1,lbnd2:ubnd2,1:ubnd_z) :: bc_t0 & !<-- Interpolate between these - ,bc_t1 ! two boundary region states. + real,dimension(lbnd1:ubnd1,lbnd2:ubnd2,1:ubnd_z),intent(in) :: bc_t0 & !<-- Interpolate between these + ,bc_t1 ! two boundary region states. ! character(len=*),intent(in) :: bc_vbl_name +! + logical,intent(in) :: blend !<-- Can blending be applied to this variable? ! !--------------------- !*** Output variables @@ -4665,7 +4692,7 @@ subroutine bc_time_interpolation(array & ! integer :: i,j,k ! - real :: blend_value,factor_dist,fraction_interval + real :: blend_value,factor_dist,fraction_interval,rdenom ! !--------------------------------------------------------------------- !********************************************************************* @@ -4691,6 +4718,15 @@ subroutine bc_time_interpolation(array & enddo ! !--------------------------------------------------------------------- +!*** If this tracer is not in the external BC file then it will not +!*** be blended. +!--------------------------------------------------------------------- +! + if(.not.blend)then + return + endif +! +!--------------------------------------------------------------------- !*** Use specified external data to blend with integration values !*** across nrows_blend rows immediately within the domain's !*** boundary rows. The weighting of the external data drops @@ -4990,8 +5026,7 @@ end subroutine regional_bc_t1_to_t0 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ !--------------------------------------------------------------------- ! - subroutine convert_to_virt_pot_temp(isd,ied,jsd,jed,npz & - ,sphum,liq_wat ) + subroutine convert_to_virt_pot_temp(isd,ied,jsd,jed,npz) ! !----------------------------------------------------------------------- !*** Convert the incoming sensible temperature to virtual potential @@ -5000,13 +5035,11 @@ subroutine convert_to_virt_pot_temp(isd,ied,jsd,jed,npz & implicit none !----------------------------------------------------------------------- ! -!--------------------- -!*** Input arguments -!--------------------- +!------------------------ +!*** Argument variables +!------------------------ ! integer,intent(in) :: isd,ied,jsd,jed,npz -! - integer,intent(in) :: liq_wat,sphum ! !--------------------- !*** Local variables @@ -5136,11 +5169,11 @@ subroutine compute_vpt ! do j=j1,j2 do i=i1,i2 - dp1 = zvir*q(i,j,k,sphum) + dp1 = zvir*q(i,j,k,sphum_index) #ifdef USE_COND #ifdef MOIST_CAPPA - cvm=(1.-q(i,j,k,sphum)+q_con(i,j,k))*cv_air & - +q(i,j,k,sphum)*cv_vap+q(i,j,k,liq_wat)*c_liq + cvm=(1.-q(i,j,k,sphum_index)+q_con(i,j,k))*cv_air & + +q(i,j,k,sphum_index)*cv_vap+q(i,j,k,liq_water_index)*c_liq pkz=exp(cappa(i,j,k)*log(rdg*delp(i,j,k)*pt(i,j,k) & *(1.+dp1)*(1.-q_con(i,j,k))/delz(i,j,k))) #else @@ -5718,25 +5751,35 @@ end subroutine dump_field_2d !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ !----------------------------------------------------------------------- ! - subroutine create_restart_with_bcs(Atm) + subroutine prepare_full_fields(Atm) ! !----------------------------------------------------------------------- -!*** Create the netcdf files into which full fields of the -!*** restart variables will be written INCLUDING BOUNDARY ROWS -!*** so the GSI can update both the interior and BCs. +!*** Prepare the objects that will hold the names and values of +!*** the core and tracer fields to be written into the expanded +!*** restart files that include the boundary rows so the GSI +!*** can update both the interior and BCs. !----------------------------------------------------------------------- ! - integer,parameter :: ndims_core=6 & !<-- # of core restart dimensions - ,ndims_tracers=4 !<-- # of tracer restart dimensions +!------------------------ +!*** Argument variables +!------------------------ ! type(fv_atmos_type),target,intent(inout) :: Atm !<-- Atm object for the current domain ! - integer :: dimid,n,na,natts & - ,ncid_core,ncid_core_new & - ,ncid_tracers,ncid_tracers_new & - ,nctype,ndims,ngatts,nn,nv_core,nv_tracers,var_id +!--------------------- +!*** Local variables +!--------------------- +! + integer :: index,istat,n & + ,ncid_core_new & + ,ncid_tracers_new & + ,ndims,nkount,nv_core,nv_tracers & + ,var_id +! + integer :: lbnd1,lbnd2,lbnd3,ubnd1,ubnd2,ubnd3 ! integer,dimension(ndims_core) :: dim_lengths_core +! integer,dimension(ndims_tracers) :: dim_lengths_tracers ! integer,dimension(1:4) :: dimids=(/0,0,0,0/) @@ -5772,7 +5815,9 @@ subroutine create_restart_with_bcs(Atm) !----------------------------------------------------------------------- !*** All tasks are given pointers into the model data that will !*** be written to the new restart file. The following are the -!*** prognostic variables in the core rstart file. +!*** prognostic variables in the core restart file. Note that +!*** we must add the halo region back into DZ since we need the +!*** domain boundary points for all the fields. !----------------------------------------------------------------------- ! allocate(fields_core(1:nvars_core)) @@ -5786,7 +5831,15 @@ subroutine create_restart_with_bcs(Atm) fields_core(3)%ptr=>Atm%w fields_core(3)%name='W' ! - fields_core(4)%ptr=>Atm%delz + lbnd1=lbound(Atm%delz,1) + ubnd1=ubound(Atm%delz,1) + lbnd2=lbound(Atm%delz,2) + ubnd2=ubound(Atm%delz,2) + lbnd3=1 + ubnd3=ubound(Atm%delz,3) + allocate(fields_core(4)%ptr(lbnd1-nhalo_model:ubnd1+nhalo_model & + ,lbnd2-nhalo_model:ubnd2+nhalo_model & + ,lbnd3:ubnd3)) fields_core(4)%name='DZ' ! fields_core(5)%ptr=>Atm%pt @@ -5798,267 +5851,70 @@ subroutine create_restart_with_bcs(Atm) allocate(fields_core(7)%ptr(lbound(Atm%phis,1):ubound(Atm%phis,1) & ,lbound(Atm%phis,2):ubound(Atm%phis,2) & ,1:1)) - fields_core(7)%ptr(:,:,1)=Atm%phis(:,:) !<-- For generality treat the 2-D phis as 3-D + fields_core(7)%ptr(:,:,1)=Atm%phis(:,:) !<-- For generality treat the 2-D phis as 3-D fields_core(7)%name='phis' ! !----------------------------------------------------------------------- -!*** Also point at the proper tracer arrays and set their names. +!*** We need to point at the tracers in the model's tracer array. +!*** Those tracers depend on the physics that was selected so they +!*** cannot be pre-specified like the variables in the core restart +!*** file were. Read them from the expanded tracer restart file +!*** that was created prior to the start for the forecast. !----------------------------------------------------------------------- ! - allocate(fields_tracers(1:nvars_tracers)) + call check(nf90_open(path=filename_tracers_new & !<-- The expanded tracer restart file. + ,mode=nf90_nowrite & !<-- File access. + ,ncid=ncid_tracers_new )) !<-- The expanded tracer restart file's ID ! - lbnd_x_tracers=lbound(Atm%q,1) - ubnd_x_tracers=ubound(Atm%q,1) - lbnd_y_tracers=lbound(Atm%q,2) - ubnd_y_tracers=ubound(Atm%q,2) -! - fields_tracers(1)%ptr=>Atm%q(:,:,:,sphum_index) - fields_tracers(1)%name='sphum' -! - fields_tracers(2)%ptr=>Atm%q(:,:,:,liq_water_index) - fields_tracers(2)%name='liq_wat' -! - fields_tracers(3)%ptr=>Atm%q(:,:,:,ice_water_index) - fields_tracers(3)%name='ice_wat' -! - fields_tracers(4)%ptr=>Atm%q(:,:,:,rain_water_index) - fields_tracers(4)%name='rainwat' -! - fields_tracers(5)%ptr=>Atm%q(:,:,:,snow_water_index) - fields_tracers(5)%name='snowwat' -! - fields_tracers(6)%ptr=>Atm%q(:,:,:,graupel_index) - fields_tracers(6)%name='graupel' -! - fields_tracers(7)%ptr=>Atm%q(:,:,:,o3mr_index) - fields_tracers(7)%name='o3mr' -! - fields_tracers(8)%ptr=>Atm%q(:,:,:,cld_amt_index) - fields_tracers(8)%name='cld_amt' -! -!----------------------------------------------------------------------- -!*** Only one compute task prepares the new restart files. The task -!*** with the highest rank knows the upper bounds of the domain so -!*** it does the work. All other tasks may exit. -!----------------------------------------------------------------------- -! - if(mpp_pe()/=Atm%layout(1)*Atm%layout(2)-1)then - return - endif -! - call check(nf90_create(filename_core_new & - ,cmode=or(nf90_clobber,nf90_64bit_offset) & - ,ncid=ncid_core_new)) -! -!----------------------------------------------------------------------- -!*** Define the output file's dimensions and insert them into the file. -!----------------------------------------------------------------------- -! - dim_lengths_core(1)=Atm%bd%ied+nhalo_model !<-- x - dim_lengths_core(2)=dim_lengths_core(1)+1 !<-- x+1 - dim_lengths_core(3)=Atm%bd%jed+nhalo_model+1 !<-- y+1 - dim_lengths_core(4)=dim_lengths_core(3)-1 !<-- y - dim_lengths_core(5)=Atm%npz !<-- z - dim_lengths_core(6)=nf90_unlimited !<-- time -! - do n=1,ndims_core - call check(nf90_def_dim(ncid_core_new & - ,dim_names_core(n) & - ,dim_lengths_core(n) & - ,dimid)) - enddo -! -!----------------------------------------------------------------------- -!*** The new file's variables must be defined while that file -!*** is still in define mode. Define each of the restart file's -!*** variables in the new file. -!----------------------------------------------------------------------- + call check(nf90_inquire(ncid =ncid_tracers_new & !<-- The expanded tracer restart file's ID. + ,nvariables=nv_tracers )) !<-- The TOTAL number of tracer restart file variables. ! - call check(nf90_open(filename_core,nf90_nowrite,ncid_core)) !<-- The core restart file's ID -! - call check(nf90_inquire(ncid_core,nvariables=nv_core)) !<-- The TOTAL number of core restart file variables -! - do n=1,nv_core - var_id=n - call check(nf90_inquire_variable(ncid_core,var_id,var_name,nctype & !<-- Name and type of this variable - ,ndims,dimids,natts)) !<-- # of dimensions and attributes in this variable -! - call check(nf90_def_var(ncid_core_new,var_name,nctype,dimids(1:ndims),var_id)) !<-- Define the variable in the new file. -! -!----------------------------------------------------------------------- -!*** Copy this variable's attributes to the new core file's -!*** variable. -!----------------------------------------------------------------------- -! - if(natts>0)then - do na=1,natts - call check(nf90_inq_attname(ncid_core,var_id,na,att_name)) !<-- Get the attribute's name and ID from restart. - call check(nf90_copy_att(ncid_core,var_id,att_name,ncid_core_new,var_id)) !<-- Copy to the new file. - enddo + nfields_tracers=nv_tracers-ndims_tracers !<-- # of 3-D tracer fields + allocate(fields_tracers(1:nfields_tracers),stat=istat) + if(istat/=0)then + call mpp_error(FATAL,' Failed to allocate fields_tracers.') + else + if(is_master())then + write(0,33012)nfields_tracers +33012 format(' Allocated fields_tracers(1:',i3,')') endif -! - enddo -! -!----------------------------------------------------------------------- -!*** Find the number of global attributes in the restart file and -!*** copy them to the new file. -!----------------------------------------------------------------------- -! - call check(nf90_inquire(ncid_core,nattributes=ngatts)) -! - do n=1,ngatts - call check(nf90_inq_attname(ncid_core,nf90_global,n,att_name)) - call check(nf90_copy_att(ncid_core,nf90_global,att_name,ncid_core_new,nf90_global)) - enddo -! -!----------------------------------------------------------------------- -! - call check(nf90_enddef(ncid_core_new)) !<-- Put the output file into data mode. -! -!----------------------------------------------------------------------- -!*** Insert the dimension values. -!----------------------------------------------------------------------- -! - do n=1,ndims_core-1 - allocate(dim_values(1:dim_lengths_core(n))) - do nn=1,dim_lengths_core(n) - dim_values(nn)=nn - enddo - var_id=n - call check(nf90_put_var(ncid_core_new,var_id & - ,dim_values(:) & - ,start=(/1/) & - ,count=(/dim_lengths_core(n)/))) - deallocate(dim_values) - enddo -! - write( 0,*)' nf90_unlimited=',nf90_unlimited - var_id=ndims_core !<-- Time is the final dimension; treat it separately. - allocate(dim_values(1)) - dim_values(1)=1 - call check(nf90_put_var(ncid_core_new,var_id & - ,dim_values(:))) - deallocate(dim_values) -! -!----------------------------------------------------------------------- -! - call check(nf90_close(ncid_core_new)) - call check(nf90_close(ncid_core)) -! -!----------------------------------------------------------------------- -!*** The second file to be handled is the tracer restart file. -!----------------------------------------------------------------------- -! - call check(nf90_create(filename_tracers_new & - ,cmode=or(nf90_clobber,nf90_64bit_offset) & - ,ncid=ncid_tracers_new)) -! -!----------------------------------------------------------------------- -!*** Define the output file's dimensions and insert them into the file. -!----------------------------------------------------------------------- -! - dim_lengths_tracers(1)=Atm%bd%ied+nhalo_model !<-- x - dim_lengths_tracers(2)=Atm%bd%jed+nhalo_model !<-- y - dim_lengths_tracers(3)=Atm%npz !<-- z - dim_lengths_tracers(4)=nf90_unlimited !<-- time -! - do n=1,ndims_tracers - call check(nf90_def_dim(ncid_tracers_new & - ,dim_names_tracers(n) & - ,dim_lengths_tracers(n) & - ,dimid)) - enddo -! -!----------------------------------------------------------------------- -!*** The new file's variables must be defined while that file -!*** is still in define mode. Define each of the restart file's -!*** variables in the new file. -!----------------------------------------------------------------------- -! - call check(nf90_open(filename_tracers,nf90_nowrite,ncid_tracers)) !<-- The tracer restart file's ID -! - call check(nf90_inquire(ncid_tracers,nvariables=nv_tracers)) !<-- The TOTAL number of tracer restart file variables + endif + nkount=0 ! do n=1,nv_tracers var_id=n - call check(nf90_inquire_variable(ncid_tracers,var_id,var_name,nctype & !<-- Name and type of this variable - ,ndims,dimids,natts)) !<-- # of dimensions and attributes in this variable -! - call check(nf90_def_var(ncid_tracers_new,var_name,nctype,dimids(1:ndims),var_id)) !<-- Define the variable in the new file. -! -!----------------------------------------------------------------------- -!*** Copy this variable's attributes to the new tracers file's -!*** variable. -!----------------------------------------------------------------------- -! - if(natts>0)then - do na=1,natts - call check(nf90_inq_attname(ncid_tracers,var_id,na,att_name)) !<-- Get the attribute's name and ID from restart. - call check(nf90_copy_att(ncid_tracers,var_id,att_name,ncid_tracers_new,var_id)) !<-- Copy to the new file. - enddo + call check(nf90_inquire_variable(ncid =ncid_tracers_new & !<-- The file's ID. + ,varid=var_id & !<-- The variable's ID. + ,name =var_name )) !<-- The variable's name. +! + if(n>ndims_tracers)then + nkount=nkount+1 + fields_tracers(nkount)%name=trim(var_name) + index=get_tracer_index(MODEL_ATMOS, trim(var_name)) + fields_tracers(nkount)%ptr=>Atm%q(:,:,:, index) endif ! enddo ! !----------------------------------------------------------------------- -!*** Find the number of global attributes in the restart file and -!*** copy them to the new file. -!----------------------------------------------------------------------- -! - call check(nf90_inquire(ncid_tracers,nattributes=ngatts)) -! - do n=1,ngatts - call check(nf90_inq_attname(ncid_tracers,nf90_global,n,att_name)) - call check(nf90_copy_att(ncid_tracers,nf90_global,att_name,ncid_tracers_new,nf90_global)) - enddo -! -!----------------------------------------------------------------------- -! - call check(nf90_enddef(ncid_tracers_new)) !<-- Put the output file into data mode. -! -!----------------------------------------------------------------------- -!*** Insert the dimension values. -!----------------------------------------------------------------------- -! - do n=1,ndims_tracers-1 - allocate(dim_values(1:dim_lengths_tracers(n))) - do nn=1,dim_lengths_tracers(n) - dim_values(nn)=nn - enddo - var_id=n - call check(nf90_put_var(ncid_tracers_new,var_id & - ,dim_values(:) & - ,start=(/1/) & - ,count=(/dim_lengths_tracers(n)/))) - deallocate(dim_values) - enddo -! - var_id=ndims_tracers !<-- Time is the final dimension; treat it separately. - allocate(dim_values(1)) - dim_values(1)=1 - call check(nf90_put_var(ncid_tracers_new,var_id & - ,dim_values(:))) - deallocate(dim_values) -! -!----------------------------------------------------------------------- ! call check(nf90_close(ncid_tracers_new)) - call check(nf90_close(ncid_tracers)) ! !----------------------------------------------------------------------- ! - end subroutine create_restart_with_bcs + end subroutine prepare_full_fields ! !----------------------------------------------------------------------- -!-------------------------------------------------------------------------------------- +!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ +!----------------------------------------------------------------------- ! subroutine write_full_fields(Atm) ! -!-------------------------------------------------------------------------------------- +!----------------------------------------------------------------------- !*** Write out full fields of the primary restart variables !*** INCLUDING BOUNDARY ROWS so the GSI can include BCs in its !*** update. This is done in a restart look-alike file. -!-------------------------------------------------------------------------------------- +!----------------------------------------------------------------------- ! type(fv_atmos_type), intent(inout), target :: Atm ! @@ -6177,30 +6033,35 @@ subroutine write_full_fields(Atm) !*** halo points need to be converted to sensible temperature. !*** Each boundary task will now do the conversion before the !*** values are gathered onto the root task for writing out. +!*** Also since the DZ array no longer has halo points we must +!*** insert values back into the domain boundary rows of the +!*** object holding DZ. !----------------------------------------------------------------------- ! - if(trim(fields_core(nv)%name)=='T')then - call sensible_temp(istart,iend,jstart,jend,nz & - ,Atm & - ,fields_core(nv)%ptr(istart:iend,jstart:jend,:)) + if(trim(fields_core(nv)%name)=='T' & + .or. & + trim(fields_core(nv)%name)=='DZ') then +! + call apply_delz_boundary(istart,iend,jstart,jend,nz & + ,Atm & + ,fields_core(nv)%name & + ,fields_core(nv)%ptr(istart:iend,jstart:jend,:)) endif ! !----------------------------------------------------------------------- -!*** Since we are gathering onto a single task then do so one layer -!*** at a time to avoid potential memory problems for large high -!*** resolution domains. Then that task writes the full data to the -!*** new file. +!*** Gather onto a single task one layer at a time. That task +!*** writes the full data to the new larger restart file. !----------------------------------------------------------------------- ! do k=1,nz - call mpp_gather(istart,iend,jstart,jend & - ,pelist, fields_core(nv)%ptr(istart:iend,jstart:jend,k) & + call mpp_gather(istart,iend,jstart,jend & + ,pelist, fields_core(nv)%ptr(istart:iend,jstart:jend,k) & ,global_field(:,:,k), is_root_pe, halo, halo) ! if(is_root_pe)then - call check(nf90_put_var(ncid_core_new,var_id & - ,global_field(:,:,k) & - ,start=(/1,1,k/) & + call check(nf90_put_var(ncid_core_new,var_id & + ,global_field(:,:,k) & + ,start=(/1,1,k/) & ,count=(/count_i,count_j,1/))) endif enddo @@ -6219,7 +6080,7 @@ subroutine write_full_fields(Atm) ! if(is_root_pe)then call check(nf90_open(filename_tracers_new,nf90_write,ncid_tracers_new)) !<-- Open the new netcdf file - write(0,*)' Opened core restart with BCs: ',trim(filename_tracers_new) + write(0,*)' Opened tracer restart with BCs: ',trim(filename_tracers_new) endif ! !----------------------------------------------------------------------- @@ -6227,7 +6088,7 @@ subroutine write_full_fields(Atm) !*** boundary rows? !----------------------------------------------------------------------- ! - call mpp_get_global_domain (atm%domain, isg, ieg, jsg, jeg, position=CENTER ) + call mpp_get_global_domain (Atm%domain, isg, ieg, jsg, jeg, position=CENTER ) istart_g=isg-halo iend_g =ieg+halo jstart_g=jsg-halo @@ -6250,6 +6111,11 @@ subroutine write_full_fields(Atm) !----------------------------------------------------------------------- !*** The following are local bounds based on the global indexing. !----------------------------------------------------------------------- +! + lbnd_x_tracers=lbound(Atm%q,1) + ubnd_x_tracers=ubound(Atm%q,1) + lbnd_y_tracers=lbound(Atm%q,2) + ubnd_y_tracers=ubound(Atm%q,2) ! istart=lbnd_x_tracers if(istart>1)then @@ -6301,7 +6167,7 @@ subroutine write_full_fields(Atm) !*** Loop through that file's prognostic tracers. !----------------------------------------------------------------------- ! - vbls_tracers: do nv=1,nvars_tracers + vbls_tracers: do nv=1,nfields_tracers ! var_name=fields_tracers(nv)%name if(is_root_pe)then @@ -6309,10 +6175,8 @@ subroutine write_full_fields(Atm) endif ! !----------------------------------------------------------------------- -!*** Since we are gathering onto a single task then do so one layer -!*** at a time to avoid potential memory problems for large high -!*** resolution domains. Then that task writes the full data to the -!*** new file. +!*** Gather onto a single task one layer at a time. That task +!*** writes the full data to the new larger restart file. !----------------------------------------------------------------------- ! do k=1,nz @@ -6321,9 +6185,9 @@ subroutine write_full_fields(Atm) ,global_field(:,:,k), is_root_pe, halo, halo) ! if(is_root_pe)then - call check(nf90_put_var(ncid_tracers_new,var_id & - ,global_field(:,:,k) & - ,start=(/1,1,k/) & + call check(nf90_put_var(ncid_tracers_new,var_id & + ,global_field(:,:,k) & + ,start=(/1,1,k/) & ,count=(/count_i,count_j,1/))) endif enddo @@ -6341,13 +6205,15 @@ end subroutine write_full_fields !--------------------------------------------------------------------- !--------------------------------------------------------------------- ! - subroutine sensible_temp(istart,iend,jstart,jend,nz & - ,Atm & - ,temp) + subroutine apply_delz_boundary(istart,iend,jstart,jend,nz & + ,Atm & + ,name & + ,field) ! !--------------------------------------------------------------------- -!*** Convert the special potential temperature in the domain -!*** halo rows to sensible temperature. +!*** Use the current boundary values of delz to convert the +!*** boundary potential temperature to sensible temperature +!*** and to fill in the boundary rows of the 3D delz array. !--------------------------------------------------------------------- ! !------------------------ @@ -6355,16 +6221,19 @@ subroutine sensible_temp(istart,iend,jstart,jend,nz & !------------------------ ! integer,intent(in) :: istart,iend,jstart,jend +! + character(len=*),intent(in) :: name ! type(fv_atmos_type),intent(inout) :: Atm ! - real,dimension(istart:iend,jstart:jend,1:nz),intent(inout) :: temp + real,dimension(istart:iend,jstart:jend,1:nz),intent(inout) :: field ! !--------------------- !*** Local variables !--------------------- ! integer :: i1,i2,j1,j2,nz + integer :: lbnd1,lbnd2,ubnd1,ubnd2,i,j,k ! real :: rdg ! @@ -6373,6 +6242,28 @@ subroutine sensible_temp(istart,iend,jstart,jend,nz & !--------------------------------------------------------------------- !********************************************************************* !--------------------------------------------------------------------- +! +!fill interior delz points before dealing with boundaries +! +!--------------------------------------------------------------------- +!*** Fill the interior of the full delz array using Atm%delz +!*** which does not have a boundary. +!--------------------------------------------------------------------- +! + if (trim(name)=='DZ') then + lbnd1=lbound(Atm%delz,1) + ubnd1=ubound(Atm%delz,1) + lbnd2=lbound(Atm%delz,2) + ubnd2=ubound(Atm%delz,2) +! + do k=1,nz + do j=lbnd2,ubnd2 + do i=lbnd1,ubnd1 + field(i,j,k)=Atm%delz(i,j,k) + enddo + enddo + enddo + endif ! if(.not.(north_bc.or.south_bc.or.east_bc.or.west_bc))then return !<-- Tasks not on the boundary may exit. @@ -6388,7 +6279,13 @@ subroutine sensible_temp(istart,iend,jstart,jend,nz & j1=jstart j2=jstart+nhalo_model-1 delz_ptr=>delz_auxiliary%north - call compute_halo_t +! + if(trim(name)=='T')then + call compute_halo_t + elseif(trim(name)=='DZ')then + call fill_delz + endif +! endif ! if(south_bc)then @@ -6397,7 +6294,13 @@ subroutine sensible_temp(istart,iend,jstart,jend,nz & j1=jend-nhalo_model+1 j2=jend delz_ptr=>delz_auxiliary%south - call compute_halo_t +! + if(trim(name)=='T')then + call compute_halo_t + elseif(trim(name)=='DZ')then + call fill_delz + endif +! endif ! if(east_bc)then @@ -6411,7 +6314,13 @@ subroutine sensible_temp(istart,iend,jstart,jend,nz & j2=jend-nhalo_model endif delz_ptr=>delz_auxiliary%east - call compute_halo_t +! + if(trim(name)=='T')then + call compute_halo_t + elseif(trim(name)=='DZ')then + call fill_delz + endif +! endif ! if(west_bc)then @@ -6425,7 +6334,13 @@ subroutine sensible_temp(istart,iend,jstart,jend,nz & j2=jend-nhalo_model endif delz_ptr=>delz_auxiliary%west - call compute_halo_t +! + if(trim(name)=='T')then + call compute_halo_t + elseif(trim(name)=='DZ')then + call fill_delz + endif +! endif ! !--------------------------------------------------------------------- @@ -6456,8 +6371,8 @@ subroutine compute_halo_t part1=(1.+dp1)*(1.-Atm%q_con(i,j,k)) part2=rdg*Atm%delp(i,j,k)*(1.+dp1)*(1.-Atm%q_con(i,j,k)) & /delz_ptr(i,j,k) - temp(i,j,k)=exp((log(temp(i,j,k))-log(part1)+cappa*log(part2)) & - /(1.-cappa)) + field(i,j,k)=exp((log(field(i,j,k))-log(part1)+cappa*log(part2)) & + /(1.-cappa)) enddo enddo enddo @@ -6466,7 +6381,34 @@ subroutine compute_halo_t end subroutine compute_halo_t !--------------------------------------------------------------------- ! - end subroutine sensible_temp + subroutine fill_delz +! +!--------------------------------------------------------------------- +! + integer :: i,j,k + integer :: lbnd1,lbnd2,ubnd1,ubnd2 +! +!--------------------------------------------------------------------- +!********************************************************************* +!--------------------------------------------------------------------- +! +!--------------------------------------------------------------------- +!*** Now fill the boundary rows using data from the BC files. +!--------------------------------------------------------------------- +! + do k=1,nz + do j=j1,j2 + do i=i1,i2 + field(i,j,k)=delz_ptr(i,j,k) + enddo + enddo + enddo +! +!--------------------------------------------------------------------- + end subroutine fill_delz +!--------------------------------------------------------------------- +! + end subroutine apply_delz_boundary ! !--------------------------------------------------------------------- !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ diff --git a/tools/external_ic.F90 b/tools/external_ic.F90 index d122f2dd5..432703501 100644 --- a/tools/external_ic.F90 +++ b/tools/external_ic.F90 @@ -577,7 +577,6 @@ subroutine get_nggps_ic (Atm, fv_domain, dt_atmos ) call get_data_source(source,Atm%flagstruct%regional) if (trim(source) == source_fv3gfs) then call mpp_error(NOTE, "READING FROM REGRIDDED FV3GFS NEMSIO FILE") - levp = 65 endif ! !--- read in ak and bk from the gfs control file using fms_io read_data --- @@ -810,7 +809,7 @@ subroutine get_nggps_ic (Atm, fv_domain, dt_atmos ) Atm%gridstruct%dxc, Atm%gridstruct%dyc, Atm%gridstruct%sin_sg, & Atm%flagstruct%n_zs_filter, cnst_0p20*Atm%gridstruct%da_min, & .false., oro_g, Atm%gridstruct%bounded_domain, & - Atm%domain, Atm%bd) + Atm%domain, Atm%bd) if ( is_master() ) write(*,*) 'Warning !!! del-2 terrain filter has been applied ', & Atm%flagstruct%n_zs_filter, ' times' else if( Atm%flagstruct%nord_zs_filter == 4 ) then @@ -818,7 +817,7 @@ subroutine get_nggps_ic (Atm, fv_domain, dt_atmos ) Atm%gridstruct%dx, Atm%gridstruct%dy, & Atm%gridstruct%dxc, Atm%gridstruct%dyc, Atm%gridstruct%sin_sg, & Atm%flagstruct%n_zs_filter, .false., oro_g, & - Atm%gridstruct%bounded_domain, & + Atm%gridstruct%bounded_domain, & Atm%domain, Atm%bd) if ( is_master() ) write(*,*) 'Warning !!! del-4 terrain filter has been applied ', & Atm%flagstruct%n_zs_filter, ' times' diff --git a/tools/fv_eta.F90_65lyrs b/tools/fv_eta.F90_65lyrs deleted file mode 100644 index 3c49530f8..000000000 --- a/tools/fv_eta.F90_65lyrs +++ /dev/null @@ -1,3093 +0,0 @@ -!*********************************************************************** -!* GNU Lesser General Public License -!* -!* This file is part of the FV3 dynamical core. -!* -!* The FV3 dynamical core is free software: you can redistribute it -!* and/or modify it under the terms of the -!* GNU Lesser General Public License as published by the -!* Free Software Foundation, either version 3 of the License, or -!* (at your option) any later version. -!* -!* The FV3 dynamical core is distributed in the hope that it will be -!* useful, but WITHOUT ANYWARRANTY; 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 Lesser General Public -!* License along with the FV3 dynamical core. -!* If not, see . -!*********************************************************************** - -!>@brief The module 'fv_eta' contains routine to set up the reference -!! (Eulerian) pressure coordinate - -module fv_eta_mod - -! -! -! -! -! -! -! -! -! -! -! -! -! -! -! -! -! -!
Module NameFunctions Included
constants_modkappa, grav, cp_air, rdgas
fv_mp_modis_master
mpp_modmpp_error, FATAL
- - use constants_mod, only: kappa, grav, cp_air, rdgas - use fv_mp_mod, only: is_master - use mpp_mod, only: FATAL, mpp_error - implicit none - private - public set_eta, set_external_eta, get_eta_level, compute_dz_var, & - compute_dz_L32, compute_dz_L101, set_hybrid_z, compute_dz, & - gw_1d, sm1_edge, hybrid_z_dz - - contains - -!!!NOTE: USE_VAR_ETA not used in fvGFS -#ifdef USE_VAR_ETA - subroutine set_eta(km, ks, ptop, ak, bk) -! This is the easy to use version of the set_eta - integer, intent(in):: km ! vertical dimension - integer, intent(out):: ks ! number of pure p layers - real:: a60(61),b60(61) -! Thfollowing L63 setting is the same as NCEP GFS's L64 except the top -! 3 layers - data a60/300.0000, 430.00000, 558.00000, & - 700.00000, 863.05803, 1051.07995, & - 1265.75194, 1510.71101, 1790.05098, & - 2108.36604, 2470.78817, 2883.03811, & - 3351.46002, 3883.05187, 4485.49315, & - 5167.14603, 5937.04991, 6804.87379, & - 7780.84698, 8875.64338, 10100.20534, & - 11264.35673, 12190.64366, 12905.42546, & - 13430.87867, 13785.88765, 13986.77987, & - 14047.96335, 13982.46770, 13802.40331, & - 13519.33841, 13144.59486, 12689.45608, & - 12165.28766, 11583.57006, 10955.84778, & - 10293.60402, 9608.08306, 8910.07678, & - 8209.70131, 7516.18560, 6837.69250, & - 6181.19473, 5552.39653, 4955.72632, & - 4394.37629, 3870.38682, 3384.76586, & - 2937.63489, 2528.37666, 2155.78385, & - 1818.20722, 1513.68173, 1240.03585, & - 994.99144, 776.23591, 581.48797, & - 408.53400, 255.26520, 119.70243, 0. / - - data b60/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00201, 0.00792, 0.01755, & - 0.03079, 0.04751, 0.06761, & - 0.09097, 0.11746, 0.14690, & - 0.17911, 0.21382, 0.25076, & - 0.28960, 0.32994, 0.37140, & - 0.41353, 0.45589, 0.49806, & - 0.53961, 0.58015, 0.61935, & - 0.65692, 0.69261, 0.72625, & - 0.75773, 0.78698, 0.81398, & - 0.83876, 0.86138, 0.88192, & - 0.90050, 0.91722, 0.93223, & - 0.94565, 0.95762, 0.96827, & - 0.97771, 0.98608, 0.99347, 1./ - real, intent(out):: ak(km+1) - real, intent(out):: bk(km+1) - real, intent(out):: ptop ! model top (Pa) - real pint, stretch_fac - integer k - real :: s_rate = -1.0 ! dummy value to not use var_les - - pint = 100.E2 - -!- Notes --------------------------------- -! low-top: ptop = 100. ! ~45 km -! mid-top: ptop = 10. ! ~60 km -! hi -top: ptop = 1. ! ~80 km -!----------------------------------------- - select case (km) - -! Optimal number = 8 * N -1 (for 8 openMP threads) -! 16 * M -1 (for 16 openMP threads) - -#ifdef HIWPP -#ifdef SUPER_K - case (20) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (24) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (30) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (40) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (50) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (60) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (80) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 -#else - case (30) ! For Baroclinic Instability Test - ptop = 2.26e2 - pint = 250.E2 - stretch_fac = 1.03 - case (40) - ptop = 50.e2 ! For super cell test - pint = 300.E2 - stretch_fac = 1.03 - case (50) ! Mountain waves? - ptop = 30.e2 - stretch_fac = 1.025 - case (60) ! For Baroclinic Instability Test -#ifdef GFSL60 - ks = 20 - ptop = a60(1) - pint = a60(ks+1) - do k=1,km+1 - ak(k) = a60(k) - bk(k) = b60(k) - enddo -#else - ptop = 3.e2 -! pint = 250.E2 - pint = 300.E2 ! revised for Moist test - stretch_fac = 1.03 -#endif -#endif - case (64) -!!! ptop = 3.e2 - ptop = 2.0e2 - pint = 300.E2 - stretch_fac = 1.03 -#else -! *Very-low top: for idealized super-cell simulation: - case (50) - ptop = 50.e2 - pint = 250.E2 - stretch_fac = 1.03 - case (60) - ptop = 40.e2 - pint = 250.E2 - stretch_fac = 1.03 - case (90) ! super-duper cell - ptop = 40.e2 - stretch_fac = 1.025 -#endif -! Low-top: - case (31) ! N = 4, M=2 - ptop = 100. - stretch_fac = 1.035 - case (32) ! N = 4, M=2 - ptop = 100. - stretch_fac = 1.035 - case (39) ! N = 5 - ptop = 100. - stretch_fac = 1.035 - case (41) - ptop = 100. - stretch_fac = 1.035 - case (47) ! N = 6, M=3 - ptop = 100. - stretch_fac = 1.035 - case (51) - ptop = 100. - stretch_fac = 1.03 - case (52) ! very low top - ptop = 30.e2 ! for special DPM RCE experiments - stretch_fac = 1.03 -! Mid-top: - case (55) ! N = 7 - ptop = 10. - stretch_fac = 1.035 -! Hi-top: - case (63) ! N = 8, M=4 - ptop = 1. - ! c360 or c384 - stretch_fac = 1.035 - case (71) ! N = 9 - ptop = 1. - stretch_fac = 1.03 - case (79) ! N = 10, M=5 - ptop = 1. - stretch_fac = 1.03 - case (127) ! N = 10, M=5 - ptop = 1. - stretch_fac = 1.03 - case (151) - ptop = 75.e2 - pint = 500.E2 - s_rate = 1.01 - case default - ptop = 1. - stretch_fac = 1.03 - end select - -#ifdef MOUNTAIN_WAVES - call mount_waves(km, ak, bk, ptop, ks, pint) -#else - if (s_rate > 0.) then - call var_les(km, ak, bk, ptop, ks, pint, s_rate) - else - if ( km > 79 ) then - call var_hi2(km, ak, bk, ptop, ks, pint, stretch_fac) - elseif (km==5 .or. km==10 ) then -! Equivalent Shallow Water: for NGGPS, variable-resolution testing - ptop = 500.e2 - ks = 0 - do k=1,km+1 - bk(k) = real(k-1) / real (km) - ak(k) = ptop*(1.-bk(k)) - enddo - else -#ifndef GFSL60 - call var_hi(km, ak, bk, ptop, ks, pint, stretch_fac) -#endif - endif -#endif - endif - - ptop = ak(1) - pint = ak(ks+1) - - end subroutine set_eta - - subroutine mount_waves(km, ak, bk, ptop, ks, pint) - integer, intent(in):: km - real, intent(out):: ak(km+1), bk(km+1) - real, intent(out):: ptop, pint - integer, intent(out):: ks -! Local - real, parameter:: p00 = 1.E5 - real, dimension(km+1):: ze, pe1, peln, eta - real, dimension(km):: dz, dlnp - real ztop, t0, dz0, sum1, tmp1 - real ep, es, alpha, beta, gama, s_fac - integer k, k500 - - pint = 300.e2 -! s_fac = 1.05 -! dz0 = 500. - if ( km <= 60 ) then - s_fac = 1.0 - dz0 = 500. - else - s_fac = 1. - dz0 = 250. - endif - -! Basic parameters for HIWPP mountain waves - t0 = 300. -! ztop = 20.0e3; 500-m resolution in halft of the vertical domain -! ztop = real(km-1)*500. -!----------------------- -! Compute temp ptop based on isothermal atm -! ptop = p00*exp(-grav*ztop/(rdgas*t0)) - -! Lowest half has constant resolution - ze(km+1) = 0. - do k=km, km-19, -1 - ze(k) = ze(k+1) + dz0 - enddo - -! Stretching from 10-km and up: - do k=km-20, 3, -1 - dz0 = s_fac * dz0 - ze(k) = ze(k+1) + dz0 - enddo - ze(2) = ze(3) + sqrt(2.)*dz0 - ze(1) = ze(2) + 2.0*dz0 - -! call sm1_edge(1, 1, 1, 1, km, 1, 1, ze, 1) - -! Given z --> p - do k=1,km - dz(k) = ze(k) - ze(k+1) - dlnp(k) = grav*dz(k) / (rdgas*t0) - enddo - - pe1(km+1) = p00 - peln(km+1) = log(p00) - do k=km,1,-1 - peln(k) = peln(k+1) - dlnp(k) - pe1(k) = exp(peln(k)) - enddo - -! Comnpute new ptop - ptop = pe1(1) - -! Pe(k) = ak(k) + bk(k) * PS -! Locate pint and KS - ks = 0 - do k=2,km - if ( pint < pe1(k)) then - ks = k-1 - exit - endif - enddo - - if ( is_master() ) then - write(*,*) 'For (input) PINT=', 0.01*pint, ' KS=', ks, 'pint(computed)=', 0.01*pe1(ks+1) - write(*,*) 'Modified ptop =', ptop, ' ztop=', ze(1)/1000. - do k=1,km - write(*,*) k, 'ze =', ze(k)/1000. - enddo - endif - pint = pe1(ks+1) - -#ifdef NO_UKMO_HB - do k=1,ks+1 - ak(k) = pe1(k) - bk(k) = 0. - enddo - - do k=ks+2,km+1 - bk(k) = (pe1(k) - pint) / (pe1(km+1)-pint) ! bk == sigma - ak(k) = pe1(k) - bk(k) * pe1(km+1) - enddo - bk(km+1) = 1. - ak(km+1) = 0. -#else -! Problematic for non-hydrostatic - do k=1,km+1 - eta(k) = pe1(k) / pe1(km+1) - enddo - ep = eta(ks+1) - es = eta(km) -! es = 1. - alpha = (ep**2-2.*ep*es) / (es-ep)**2 - beta = 2.*ep*es**2 / (es-ep)**2 - gama = -(ep*es)**2 / (es-ep)**2 - -! Pure pressure: - do k=1,ks+1 - ak(k) = eta(k)*1.e5 - bk(k) = 0. - enddo - - do k=ks+2, km - ak(k) = alpha*eta(k) + beta + gama/eta(k) - ak(k) = ak(k)*1.e5 - enddo - ak(km+1) = 0. - - do k=ks+2, km - bk(k) = (pe1(k) - ak(k))/pe1(km+1) - enddo - bk(km+1) = 1. -#endif - - if ( is_master() ) then - tmp1 = ak(ks+1) - do k=ks+1,km - tmp1 = max(tmp1, (ak(k)-ak(k+1))/max(1.E-5, (bk(k+1)-bk(k))) ) - enddo - write(*,*) 'Hybrid Sigma-P: minimum allowable surface pressure (hpa)=', tmp1/100. - endif - - end subroutine mount_waves - -#else - !>This is the version of set_eta used in fvGFS and AM4 - !>@note 01/2018: 'set_eta' is being cleaned up. - subroutine set_eta(km, ks, ptop, ak, bk) - integer, intent(in):: km !< vertical dimension - integer, intent(out):: ks !< number of pure p layers - real, intent(out):: ak(km+1) - real, intent(out):: bk(km+1) - real, intent(out):: ptop !< model top (Pa) -! local - real a24(25),b24(25) !< GFDL AM2L24 - real a26(27),b26(27) !< Jablonowski & Williamson 26-level - real a32(33),b32(33) - real a32w(33),b32w(33) - real a47(48),b47(48) - real a48(49),b48(49) - real a52(53),b52(53) - real a54(55),b54(55) - real a56(57),b56(57) - real a60(61),b60(61) - real a63(64),b63(64) - real a64(65),b64(65) - real a68(69),b68(69) !< cjg: grid with enhanced PBL resolution - real a96(97),b96(97) !< cjg: grid with enhanced PBL resolution - real a100(101),b100(101) - real a104(105),b104(105) - real a125(126),b125(126) - - real:: p0=1000.E2 - real:: pc=200.E2 - - real pt, pint, lnpe, dlnp - real press(km+1), pt1(km) - integer k - -! Definition: press(i,j,k) = ak(k) + bk(k) * ps(i,j) - -!----------------------------------------------- -! GFDL AM2-L24: modified by SJL at the model top -!----------------------------------------------- -! data a24 / 100.0000, 1050.0000, 3474.7942, 7505.5556, 12787.2428, & - data a24 / 100.0000, 903.4465, 3474.7942, 7505.5556, 12787.2428, & - 19111.3683, 21854.9274, 22884.1866, 22776.3058, 21716.1604, & - 20073.2963, 18110.5123, 16004.7832, 13877.6253, 11812.5452, & - 9865.8840, 8073.9726, 6458.0834, 5027.9899, 3784.6085, & - 2722.0086, 1828.9752, 1090.2396, 487.4595, 0.0000 / - - data b24 / 0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0000000, & - 0.0000000, 0.0435679, 0.1102275, 0.1922249, 0.2817656, & - 0.3694997, 0.4532348, 0.5316253, 0.6038733, 0.6695556, & - 0.7285176, 0.7808017, 0.8265992, 0.8662148, 0.9000406, & - 0.9285364, 0.9522140, 0.9716252, 0.9873523, 1.0000000 / - -! Jablonowski & Williamson 26-level setup - data a26 / 219.4067, 489.5209, 988.2418, 1805.2010, 2983.7240, 4462.3340, & - 6160.5870, 7851.2430, 7731.2710, 7590.1310, 7424.0860, & - 7228.7440, 6998.9330, 6728.5740, 6410.5090, 6036.3220, & - 5596.1110, 5078.2250, 4468.9600, 3752.1910, 2908.9490, & - 2084.739, 1334.443, 708.499, 252.1360, 0.0, 0.0 / - - data b26 / 0.0, 0.0, 0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0000000,& - 0.0000000, 0.01505309, 0.03276228, 0.05359622, 0.07810627, & - 0.1069411, 0.1408637, 0.1807720, 0.2277220, 0.2829562, & - 0.3479364, 0.4243822, 0.5143168, 0.6201202, 0.7235355, & - 0.8176768, 0.8962153, 0.9534761, 0.9851122, 1.0000000 / - - -! High-resolution troposphere setup -#ifdef OLD_32 -! Revised Apr 14, 2004: PINT = 245.027 mb - data a32/100.00000, 400.00000, 818.60211, & - 1378.88653, 2091.79519, 2983.64084, & - 4121.78960, 5579.22148, 7419.79300, & - 9704.82578, 12496.33710, 15855.26306, & - 19839.62499, 24502.73262, 28177.10152, & - 29525.28447, 29016.34358, 27131.32792, & - 24406.11225, 21326.04907, 18221.18357, & - 15275.14642, 12581.67796, 10181.42843, & - 8081.89816, 6270.86956, 4725.35001, & - 3417.39199, 2317.75459, 1398.09473, & - 632.49506, 0.00000, 0.00000 / - - data b32/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.01711, & - 0.06479, 0.13730, 0.22693, & - 0.32416, 0.42058, 0.51105, & - 0.59325, 0.66628, 0.73011, & - 0.78516, 0.83217, 0.87197, & - 0.90546, 0.93349, 0.95685, & - 0.97624, 0.99223, 1.00000 / -#else -! SJL June 26, 2012 -! pint= 55.7922 - data a32/100.00000, 400.00000, 818.60211, & - 1378.88653, 2091.79519, 2983.64084, & - 4121.78960, 5579.22148, 6907.19063, & - 7735.78639, 8197.66476, 8377.95525, & - 8331.69594, 8094.72213, 7690.85756, & - 7139.01788, 6464.80251, 5712.35727, & - 4940.05347, 4198.60465, 3516.63294, & - 2905.19863, 2366.73733, 1899.19455, & - 1497.78137, 1156.25252, 867.79199, & - 625.59324, 423.21322, 254.76613, & - 115.06646, 0.00000, 0.00000 / - - data b32/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00513, & - 0.01969, 0.04299, 0.07477, & - 0.11508, 0.16408, 0.22198, & - 0.28865, 0.36281, 0.44112, & - 0.51882, 0.59185, 0.65810, & - 0.71694, 0.76843, 0.81293, & - 0.85100, 0.88331, 0.91055, & - 0.93338, 0.95244, 0.96828, & - 0.98142, 0.99223, 1.00000 / -#endif - -!--------------------- -! Wilson's 32L settings: -!--------------------- -! Top changed to 0.01 mb - data a32w/ 1.00, 26.6378, 84.5529, 228.8592, & - 539.9597, 1131.7087, 2141.8082, 3712.0454, & - 5963.5317, 8974.1873, 12764.5388, 17294.5911, & - 20857.7007, 22221.8651, 22892.7202, 22891.1641, & - 22286.0724, 21176.0846, 19673.0671, 17889.0989, & - 15927.5060, 13877.6239, 11812.5474, 9865.8830, & - 8073.9717, 6458.0824, 5027.9893, 3784.6104, & - 2722.0093, 1828.9741, 1090.2397, 487.4575, & - 0.0000 / - - data b32w/ 0.0000, 0.0000, 0.0000, 0.0000, & - 0.0000, 0.0000, 0.0000, 0.0000, & - 0.0000, 0.0000, 0.0000, 0.0000, & - 0.0159, 0.0586, 0.1117, 0.1734, & - 0.2415, 0.3137, 0.3878, 0.4619, & - 0.5344, 0.6039, 0.6696, 0.7285, & - 0.7808, 0.8266, 0.8662, 0.9000, & - 0.9285, 0.9522, 0.9716, 0.9874, & - 1.0000 / - - -#ifdef OLD_L47 -! QBO setting with ptop = 0.1 mb and p_full=0.17 mb; pint ~ 100 mb - data a47/ 10.00000, 24.45365, 48.76776, & - 85.39458, 133.41983, 191.01402, & - 257.94919, 336.63306, 431.52741, & - 548.18995, 692.78825, 872.16512, & - 1094.18467, 1368.11917, 1704.99489, & - 2117.91945, 2622.42986, 3236.88281, & - 3982.89623, 4885.84733, 5975.43260, & - 7286.29500, 8858.72424, 10739.43477, & - 12982.41110, 15649.68745, 18811.37629, & - 22542.71275, 25724.93857, 27314.36781, & - 27498.59474, 26501.79312, 24605.92991, & - 22130.51655, 19381.30274, 16601.56419, & - 13952.53231, 11522.93244, 9350.82303, & - 7443.47723, 5790.77434, 4373.32696, & - 3167.47008, 2148.51663, 1293.15510, & - 581.62429, 0.00000, 0.00000 / - - data b47/ 0.0000, 0.0000, 0.0000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.01188, 0.04650, & - 0.10170, 0.17401, 0.25832, & - 0.34850, 0.43872, 0.52448, & - 0.60307, 0.67328, 0.73492, & - 0.78834, 0.83418, 0.87320, & - 0.90622, 0.93399, 0.95723, & - 0.97650, 0.99223, 1.00000 / -#else -! Oct 23, 2012 -! QBO setting with ptop = 0.1 mb, pint ~ 60 mb - data a47/ 10.00000, 24.45365, 48.76776, & - 85.39458, 133.41983, 191.01402, & - 257.94919, 336.63306, 431.52741, & - 548.18995, 692.78825, 872.16512, & - 1094.18467, 1368.11917, 1704.99489, & - 2117.91945, 2622.42986, 3236.88281, & - 3982.89623, 4885.84733, 5975.43260, & - 7019.26669, 7796.15848, 8346.60209, & - 8700.31838, 8878.27554, 8894.27179, & - 8756.46404, 8469.60171, 8038.92687, & - 7475.89006, 6803.68067, 6058.68992, & - 5285.28859, 4526.01565, 3813.00206, & - 3164.95553, 2589.26318, 2085.96929, & - 1651.11596, 1278.81205, 962.38875, & - 695.07046, 470.40784, 282.61654, & - 126.92745, 0.00000, 0.00000 / - data b47/ 0.0000, 0.0000, 0.0000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00267, 0.01063, 0.02393, & - 0.04282, 0.06771, 0.09917, & - 0.13786, 0.18444, 0.23925, & - 0.30193, 0.37100, 0.44379, & - 0.51695, 0.58727, 0.65236, & - 0.71094, 0.76262, 0.80757, & - 0.84626, 0.87930, 0.90731, & - 0.93094, 0.95077, 0.96733, & - 0.98105, 0.99223, 1.00000 / -#endif - - data a48/ & - 1.00000, 2.69722, 5.17136, & - 8.89455, 14.24790, 22.07157, & - 33.61283, 50.48096, 74.79993, & - 109.40055, 158.00460, 225.44108, & - 317.89560, 443.19350, 611.11558, & - 833.74392, 1125.83405, 1505.20759, & - 1993.15829, 2614.86254, 3399.78420, & - 4382.06240, 5600.87014, 7100.73115, & - 8931.78242, 11149.97021, 13817.16841, & - 17001.20930, 20775.81856, 23967.33875, & - 25527.64563, 25671.22552, 24609.29622, & - 22640.51220, 20147.13482, 17477.63530, & - 14859.86462, 12414.92533, 10201.44191, & - 8241.50255, 6534.43202, 5066.17865, & - 3815.60705, 2758.60264, 1870.64631, & - 1128.33931, 510.47983, 0.00000, & - 0.00000 / - - data b48/ & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.01253, & - 0.04887, 0.10724, 0.18455, & - 0.27461, 0.36914, 0.46103, & - 0.54623, 0.62305, 0.69099, & - 0.75016, 0.80110, 0.84453, & - 0.88127, 0.91217, 0.93803, & - 0.95958, 0.97747, 0.99223, & - 1.00000 / - -! High PBL resolution with top at 1 mb -! SJL modified May 7, 2013 to ptop ~ 100 mb - data a54/100.00000, 254.83931, 729.54278, & - 1602.41121, 2797.50667, 4100.18977, & - 5334.87140, 6455.24153, 7511.80944, & - 8580.26355, 9714.44293, 10938.62253, & - 12080.36051, 12987.13921, 13692.75084, & - 14224.92180, 14606.55444, 14856.69953, & - 14991.32121, 15023.90075, 14965.91493, & - 14827.21612, 14616.33505, 14340.72252, & - 14006.94280, 13620.82849, 13187.60470, & - 12711.98873, 12198.27003, 11650.37451, & - 11071.91608, 10466.23819, 9836.44706, & - 9185.43852, 8515.96231, 7831.01080, & - 7135.14301, 6436.71659, 5749.00215, & - 5087.67188, 4465.67510, 3889.86419, & - 3361.63433, 2879.51065, 2441.02496, & - 2043.41345, 1683.80513, 1359.31122, & - 1067.09135, 804.40101, 568.62625, & - 357.32525, 168.33263, 0.00000, & - 0.00000 / - - data b54/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00180, 0.00694, 0.01510, & - 0.02601, 0.03942, 0.05515, & - 0.07302, 0.09288, 0.11459, & - 0.13803, 0.16307, 0.18960, & - 0.21753, 0.24675, 0.27716, & - 0.30866, 0.34115, 0.37456, & - 0.40879, 0.44375, 0.47935, & - 0.51551, 0.55215, 0.58916, & - 0.62636, 0.66334, 0.69946, & - 0.73395, 0.76622, 0.79594, & - 0.82309, 0.84780, 0.87020, & - 0.89047, 0.90876, 0.92524, & - 0.94006, 0.95336, 0.96529, & - 0.97596, 0.98551, 0.99400, & - 1.00000 / - - -! The 56-L setup - data a56/ 10.00000, 24.97818, 58.01160, & - 115.21466, 199.29210, 309.39897, & - 445.31785, 610.54747, 812.28518, & - 1059.80882, 1363.07092, 1732.09335, & - 2176.91502, 2707.68972, 3334.70962, & - 4068.31964, 4918.76594, 5896.01890, & - 7009.59166, 8268.36324, 9680.41211, & - 11252.86491, 12991.76409, 14901.95764, & - 16987.01313, 19249.15733, 21689.24182, & - 23845.11055, 25330.63353, 26243.52467, & - 26663.84998, 26657.94696, 26281.61371, & - 25583.05256, 24606.03265, 23393.39510, & - 21990.28845, 20445.82122, 18811.93894, & - 17139.59660, 15473.90350, 13850.50167, & - 12294.49060, 10821.62655, 9440.57746, & - 8155.11214, 6965.72496, 5870.70511, & - 4866.83822, 3949.90019, 3115.03562, & - 2357.07879, 1670.87329, 1051.65120, & - 495.51399, 0.00000, 0.00000 / - - data b56 /0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00462, 0.01769, 0.03821, & - 0.06534, 0.09834, 0.13659, & - 0.17947, 0.22637, 0.27660, & - 0.32929, 0.38343, 0.43791, & - 0.49162, 0.54361, 0.59319, & - 0.63989, 0.68348, 0.72391, & - 0.76121, 0.79545, 0.82679, & - 0.85537, 0.88135, 0.90493, & - 0.92626, 0.94552, 0.96286, & - 0.97840, 0.99223, 1.00000 / - - data a60/ 1.7861000000e-01, 1.0805100000e+00, 3.9647100000e+00, & - 9.7516000000e+00, 1.9816580000e+01, 3.6695950000e+01, & - 6.2550570000e+01, 9.9199620000e+01, 1.4792505000e+02, & - 2.0947487000e+02, 2.8422571000e+02, 3.7241721000e+02, & - 4.7437835000e+02, 5.9070236000e+02, 7.2236063000e+02, & - 8.7076746000e+02, 1.0378138800e+03, 1.2258877300e+03, & - 1.4378924600e+03, 1.6772726600e+03, 1.9480506400e+03, & - 2.2548762700e+03, 2.6030909400e+03, 2.9988059200e+03, & - 3.4489952300e+03, 3.9616028900e+03, 4.5456641600e+03, & - 5.2114401700e+03, 5.9705644000e+03, 6.8361981800e+03, & - 7.8231906000e+03, 8.9482351000e+03, 1.0230010660e+04, & - 1.1689289750e+04, 1.3348986860e+04, 1.5234111060e+04, & - 1.7371573230e+04, 1.9789784580e+04, 2.2005564550e+04, & - 2.3550115120e+04, 2.4468583320e+04, 2.4800548800e+04, & - 2.4582445070e+04, 2.3849999620e+04, 2.2640519740e+04, & - 2.0994737150e+04, 1.8957848730e+04, 1.6579413230e+04, & - 1.4080071030e+04, 1.1753630920e+04, 9.6516996300e+03, & - 7.7938009300e+03, 6.1769062800e+03, 4.7874276000e+03, & - 3.6050497500e+03, 2.6059860700e+03, 1.7668328200e+03, & - 1.0656131200e+03, 4.8226201000e+02, 0.0000000000e+00, & - 0.0000000000e+00 / - - - data b60/ 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 5.0600000000e-03, & - 2.0080000000e-02, 4.4900000000e-02, 7.9360000000e-02, & - 1.2326000000e-01, 1.7634000000e-01, 2.3820000000e-01, & - 3.0827000000e-01, 3.8581000000e-01, 4.6989000000e-01, & - 5.5393000000e-01, 6.2958000000e-01, 6.9642000000e-01, & - 7.5458000000e-01, 8.0463000000e-01, 8.4728000000e-01, & - 8.8335000000e-01, 9.1368000000e-01, 9.3905000000e-01, & - 9.6020000000e-01, 9.7775000000e-01, 9.9223000000e-01, & - 1.0000000000e+00 / - -! This is activated by USE_GFSL63 -! Thfollowing L63 setting is the same as NCEP GFS's L64 except the top -! 3 layers - data a63/64.247, 137.790, 221.958, & - 318.266, 428.434, 554.424, & - 698.457, 863.05803, 1051.07995, & - 1265.75194, 1510.71101, 1790.05098, & - 2108.36604, 2470.78817, 2883.03811, & - 3351.46002, 3883.05187, 4485.49315, & - 5167.14603, 5937.04991, 6804.87379, & - 7780.84698, 8875.64338, 10100.20534, & - 11264.35673, 12190.64366, 12905.42546, & - 13430.87867, 13785.88765, 13986.77987, & - 14047.96335, 13982.46770, 13802.40331, & - 13519.33841, 13144.59486, 12689.45608, & - 12165.28766, 11583.57006, 10955.84778, & - 10293.60402, 9608.08306, 8910.07678, & - 8209.70131, 7516.18560, 6837.69250, & - 6181.19473, 5552.39653, 4955.72632, & - 4394.37629, 3870.38682, 3384.76586, & - 2937.63489, 2528.37666, 2155.78385, & - 1818.20722, 1513.68173, 1240.03585, & - 994.99144, 776.23591, 581.48797, & - 408.53400, 255.26520, 119.70243, 0. / - - data b63/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00201, 0.00792, 0.01755, & - 0.03079, 0.04751, 0.06761, & - 0.09097, 0.11746, 0.14690, & - 0.17911, 0.21382, 0.25076, & - 0.28960, 0.32994, 0.37140, & - 0.41353, 0.45589, 0.49806, & - 0.53961, 0.58015, 0.61935, & - 0.65692, 0.69261, 0.72625, & - 0.75773, 0.78698, 0.81398, & - 0.83876, 0.86138, 0.88192, & - 0.90050, 0.91722, 0.93223, & - 0.94565, 0.95762, 0.96827, & - 0.97771, 0.98608, 0.99347, 1./ -#ifdef GFSL64 - data a64/20.00000, 68.00000, 137.79000, & - 221.95800, 318.26600, 428.43400, & - 554.42400, 698.45700, 863.05803, & - 1051.07995, 1265.75194, 1510.71101, & - 1790.05098, 2108.36604, 2470.78817, & - 2883.03811, 3351.46002, 3883.05187, & - 4485.49315, 5167.14603, 5937.04991, & - 6804.87379, 7780.84698, 8875.64338, & - 9921.40745, 10760.99844, 11417.88354, & - 11911.61193, 12258.61668, 12472.89642, & - 12566.58298, 12550.43517, 12434.26075, & - 12227.27484, 11938.39468, 11576.46910, & - 11150.43640, 10669.41063, 10142.69482, & - 9579.72458, 8989.94947, 8382.67090, & - 7766.85063, 7150.91171, 6542.55077, & - 5948.57894, 5374.81094, 4825.99383, & - 4305.79754, 3816.84622, 3360.78848, & - 2938.39801, 2549.69756, 2194.08449, & - 1870.45732, 1577.34218, 1313.00028, & - 1075.52114, 862.90778, 673.13815, & - 504.22118, 354.22752, 221.32110, & - 103.78014, 0./ - data b64/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00179, 0.00705, 0.01564, & - 0.02749, 0.04251, 0.06064, & - 0.08182, 0.10595, 0.13294, & - 0.16266, 0.19492, 0.22950, & - 0.26615, 0.30455, 0.34435, & - 0.38516, 0.42656, 0.46815, & - 0.50949, 0.55020, 0.58989, & - 0.62825, 0.66498, 0.69987, & - 0.73275, 0.76351, 0.79208, & - 0.81845, 0.84264, 0.86472, & - 0.88478, 0.90290, 0.91923, & - 0.93388, 0.94697, 0.95865, & - 0.96904, 0.97826, 0.98642, & - 0.99363, 1./ -#else - data a64/1.00000, 3.90000, 8.70000, & - 15.42000, 24.00000, 34.50000, & - 47.00000, 61.50000, 78.60000, & - 99.13500, 124.12789, 154.63770, & - 191.69700, 236.49300, 290.38000, & - 354.91000, 431.82303, 523.09300, & - 630.92800, 757.79000, 906.45000, & - 1079.85000, 1281.00000, 1515.00000, & - 1788.00000, 2105.00000, 2470.00000, & - 2889.00000, 3362.00000, 3890.00000, & - 4475.00000, 5120.00000, 5830.00000, & - 6608.00000, 7461.00000, 8395.00000, & - 9424.46289, 10574.46880, 11864.80270, & - 13312.58890, 14937.03710, 16759.70700, & - 18804.78710, 21099.41210, 23674.03710, & - 26562.82810, 29804.11720, 32627.31640, & - 34245.89840, 34722.28910, 34155.19920, & - 32636.50390, 30241.08200, 27101.44920, & - 23362.20700, 19317.05270, 15446.17090, & - 12197.45210, 9496.39941, 7205.66992, & - 5144.64307, 3240.79346, 1518.62134, & - 0.00000, 0.00000 / - - data b64/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00813, & - 0.03224, 0.07128, 0.12445, & - 0.19063, 0.26929, 0.35799, & - 0.45438, 0.55263, 0.64304, & - 0.71703, 0.77754, 0.82827, & - 0.87352, 0.91502, 0.95235, & - 0.98511, 1.00000 / -#endif -!-->cjg - data a68/1.00000, 2.68881, 5.15524, & - 8.86683, 14.20349, 22.00278, & - 33.50807, 50.32362, 74.56680, & - 109.05958, 157.51214, 224.73844, & - 316.90481, 441.81219, 609.21090, & - 831.14537, 1122.32514, 1500.51628, & - 1986.94617, 2606.71274, 3389.18802, & - 4368.40473, 5583.41379, 7078.60015, & - 8903.94455, 11115.21886, 13774.60566, & - 16936.82070, 20340.47045, 23193.71492, & - 24870.36141, 25444.59363, 25252.57081, & - 24544.26211, 23474.29096, 22230.65331, & - 20918.50731, 19589.96280, 18296.26682, & - 17038.02866, 15866.85655, 14763.18943, & - 13736.83624, 12794.11850, 11930.72442, & - 11137.17217, 10404.78946, 9720.03954, & - 9075.54055, 8466.72650, 7887.12346, & - 7333.90490, 6805.43028, 6297.33773, & - 5805.78227, 5327.94995, 4859.88765, & - 4398.63854, 3942.81761, 3491.08449, & - 3043.04531, 2598.71608, 2157.94527, & - 1720.87444, 1287.52805, 858.02944, & - 432.71276, 8.10905, 0.00000 / - - data b68/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00283, 0.01590, & - 0.04412, 0.08487, 0.13284, & - 0.18470, 0.23828, 0.29120, & - 0.34211, 0.39029, 0.43518, & - 0.47677, 0.51536, 0.55091, & - 0.58331, 0.61263, 0.63917, & - 0.66333, 0.68552, 0.70617, & - 0.72555, 0.74383, 0.76117, & - 0.77765, 0.79335, 0.80838, & - 0.82287, 0.83693, 0.85069, & - 0.86423, 0.87760, 0.89082, & - 0.90392, 0.91689, 0.92973, & - 0.94244, 0.95502, 0.96747, & - 0.97979, 0.99200, 1.00000 / - - data a96/1.00000, 2.35408, 4.51347, & - 7.76300, 12.43530, 19.26365, & - 29.33665, 44.05883, 65.28397, & - 95.48274, 137.90344, 196.76073, & - 277.45330, 386.81095, 533.37018, & - 727.67600, 982.60677, 1313.71685, & - 1739.59104, 2282.20281, 2967.26766, & - 3824.58158, 4888.33404, 6197.38450, & - 7795.49158, 9731.48414, 11969.71024, & - 14502.88894, 17304.52434, 20134.76139, & - 22536.63814, 24252.54459, 25230.65591, & - 25585.72044, 25539.91412, 25178.87141, & - 24644.84493, 23978.98781, 23245.49366, & - 22492.11600, 21709.93990, 20949.64473, & - 20225.94258, 19513.31158, 18829.32485, & - 18192.62250, 17589.39396, 17003.45386, & - 16439.01774, 15903.91204, 15396.39758, & - 14908.02140, 14430.65897, 13967.88643, & - 13524.16667, 13098.30227, 12687.56457, & - 12287.08757, 11894.41553, 11511.54106, & - 11139.22483, 10776.01912, 10419.75711, & - 10067.11881, 9716.63489, 9369.61967, & - 9026.69066, 8687.29884, 8350.04978, & - 8013.20925, 7677.12187, 7343.12994, & - 7011.62844, 6681.98102, 6353.09764, & - 6025.10535, 5699.10089, 5375.54503, & - 5053.63074, 4732.62740, 4413.38037, & - 4096.62775, 3781.79777, 3468.45371, & - 3157.19882, 2848.25306, 2541.19150, & - 2236.21942, 1933.50628, 1632.83741, & - 1334.35954, 1038.16655, 744.22318, & - 452.71094, 194.91899, 0.00000, & - 0.00000 / - - data b96/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00193, & - 0.00974, 0.02538, 0.04876, & - 0.07817, 0.11081, 0.14514, & - 0.18007, 0.21486, 0.24866, & - 0.28088, 0.31158, 0.34030, & - 0.36701, 0.39210, 0.41554, & - 0.43733, 0.45774, 0.47707, & - 0.49540, 0.51275, 0.52922, & - 0.54495, 0.56007, 0.57459, & - 0.58850, 0.60186, 0.61471, & - 0.62715, 0.63922, 0.65095, & - 0.66235, 0.67348, 0.68438, & - 0.69510, 0.70570, 0.71616, & - 0.72651, 0.73675, 0.74691, & - 0.75700, 0.76704, 0.77701, & - 0.78690, 0.79672, 0.80649, & - 0.81620, 0.82585, 0.83542, & - 0.84492, 0.85437, 0.86375, & - 0.87305, 0.88229, 0.89146, & - 0.90056, 0.90958, 0.91854, & - 0.92742, 0.93623, 0.94497, & - 0.95364, 0.96223, 0.97074, & - 0.97918, 0.98723, 0.99460, & - 1.00000 / -!<--cjg -! -! Ultra high troposphere resolution - data a100/100.00000, 300.00000, 800.00000, & - 1762.35235, 3106.43596, 4225.71874, & - 4946.40525, 5388.77387, 5708.35540, & - 5993.33124, 6277.73673, 6571.49996, & - 6877.05339, 7195.14327, 7526.24920, & - 7870.82981, 8229.35361, 8602.30193, & - 8990.16936, 9393.46399, 9812.70768, & - 10248.43625, 10701.19980, 11171.56286, & - 11660.10476, 12167.41975, 12694.11735, & - 13240.82253, 13808.17600, 14396.83442, & - 15007.47066, 15640.77407, 16297.45067, & - 16978.22343, 17683.83253, 18415.03554, & - 19172.60771, 19957.34218, 20770.05022, & - 21559.14829, 22274.03147, 22916.87519, & - 23489.70456, 23994.40187, 24432.71365, & - 24806.25734, 25116.52754, 25364.90190, & - 25552.64670, 25680.92203, 25750.78675, & - 25763.20311, 25719.04113, 25619.08274, & - 25464.02630, 25254.49482, 24991.06137, & - 24674.32737, 24305.11235, 23884.79781, & - 23415.77059, 22901.76510, 22347.84738, & - 21759.93950, 21144.07284, 20505.73136, & - 19849.54271, 19179.31412, 18498.23400, & - 17809.06809, 17114.28232, 16416.10343, & - 15716.54833, 15017.44246, 14320.43478, & - 13627.01116, 12938.50682, 12256.11762, & - 11580.91062, 10913.83385, 10255.72526, & - 9607.32122, 8969.26427, 8342.11044, & - 7726.33606, 7122.34405, 6530.46991, & - 5950.98721, 5384.11279, 4830.01153, & - 4288.80090, 3760.55514, 3245.30920, & - 2743.06250, 2253.78294, 1777.41285, & - 1313.88054, 863.12371, 425.13088, & - 0.00000, 0.00000 / - - - data b100/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00052, 0.00209, 0.00468, & - 0.00828, 0.01288, 0.01849, & - 0.02508, 0.03266, 0.04121, & - 0.05075, 0.06126, 0.07275, & - 0.08521, 0.09866, 0.11308, & - 0.12850, 0.14490, 0.16230, & - 0.18070, 0.20009, 0.22042, & - 0.24164, 0.26362, 0.28622, & - 0.30926, 0.33258, 0.35605, & - 0.37958, 0.40308, 0.42651, & - 0.44981, 0.47296, 0.49591, & - 0.51862, 0.54109, 0.56327, & - 0.58514, 0.60668, 0.62789, & - 0.64872, 0.66919, 0.68927, & - 0.70895, 0.72822, 0.74709, & - 0.76554, 0.78357, 0.80117, & - 0.81835, 0.83511, 0.85145, & - 0.86736, 0.88286, 0.89794, & - 0.91261, 0.92687, 0.94073, & - 0.95419, 0.96726, 0.97994, & - 0.99223, 1.00000 / - - data a104/ & - 1.8827062944e-01, 7.7977549145e-01, 2.1950593583e+00, & - 4.9874566624e+00, 9.8041418997e+00, 1.7019717163e+01, & - 2.7216579591e+01, 4.0518628401e+01, 5.6749646818e+01, & - 7.5513868331e+01, 9.6315093333e+01, 1.1866706195e+02, & - 1.4216835396e+02, 1.6653733709e+02, 1.9161605772e+02, & - 2.1735580129e+02, 2.4379516604e+02, 2.7103771847e+02, & - 2.9923284173e+02, 3.2856100952e+02, 3.5922338766e+02, & - 3.9143507908e+02, 4.2542117983e+02, 4.6141487902e+02, & - 4.9965698106e+02, 5.4039638379e+02, 5.8389118154e+02, & - 6.3041016829e+02, 6.8023459505e+02, 7.3366009144e+02, & - 7.9099869949e+02, 8.5258099392e+02, 9.1875827946e+02, & - 9.8990486716e+02, 1.0664204381e+03, 1.1487325074e+03, & - 1.2372990044e+03, 1.3326109855e+03, 1.4351954993e+03, & - 1.5456186222e+03, 1.6644886848e+03, 1.7924597105e+03, & - 1.9302350870e+03, 2.0785714934e+03, 2.2382831070e+03, & - 2.4102461133e+03, 2.5954035462e+03, 2.7947704856e+03, & - 3.0094396408e+03, 3.2405873512e+03, 3.4894800360e+03, & - 3.7574811281e+03, 4.0460585279e+03, 4.3567926151e+03, & - 4.6913848588e+03, 5.0516670674e+03, 5.4396113207e+03, & - 5.8573406270e+03, 6.3071403487e+03, 6.7914704368e+03, & - 7.3129785102e+03, 7.8745138115e+03, 8.4791420557e+03, & - 9.1301611750e+03, 9.8311179338e+03, 1.0585825354e+04, & - 1.1398380836e+04, 1.2273184781e+04, 1.3214959424e+04, & - 1.4228767429e+04, 1.5320029596e+04, 1.6494540743e+04, & - 1.7758482452e+04, 1.9118430825e+04, 2.0422798801e+04, & - 2.1520147587e+04, 2.2416813461e+04, 2.3118184510e+04, & - 2.3628790785e+04, 2.3952411814e+04, 2.4092209011e+04, & - 2.4050892106e+04, 2.3830930156e+04, 2.3434818358e+04, & - 2.2865410898e+04, 2.2126326004e+04, 2.1222420323e+04, & - 2.0160313690e+04, 1.8948920926e+04, 1.7599915822e+04, & - 1.6128019809e+04, 1.4550987232e+04, 1.2889169132e+04, & - 1.1164595563e+04, 9.4227665517e+03, 7.7259097899e+03, & - 6.1538244381e+03, 4.7808126007e+03, 3.5967415552e+03, & - 2.5886394104e+03, 1.7415964865e+03, 1.0393721271e+03, & - 4.6478852032e+02, 7.0308342481e-13, 0.0000000000e+00 / - - - data b104/ & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 1.5648447298e-03, & - 6.2617046389e-03, 1.4104157933e-02, 2.5118187415e-02, & - 3.9340510972e-02, 5.6816335609e-02, 7.7596328431e-02, & - 1.0173255472e-01, 1.2927309709e-01, 1.6025505622e-01, & - 1.9469566981e-01, 2.3258141217e-01, 2.7385520518e-01, & - 3.1840233814e-01, 3.6603639170e-01, 4.1648734767e-01, & - 4.6939496013e-01, 5.2431098738e-01, 5.8071350676e-01, & - 6.3803478105e-01, 6.9495048840e-01, 7.4963750338e-01, & - 7.9975208897e-01, 8.4315257576e-01, 8.8034012292e-01, & - 9.1184389721e-01, 9.3821231526e-01, 9.6000677644e-01, & - 9.7779792223e-01, 9.9216315122e-01, 1.0000000000e+00 / - -! IFS-like L125(top 12 levels removed from IFSL137) - data a125/ 64., & - 86.895882, 107.415741, 131.425507, 159.279404, 191.338562, & - 227.968948, 269.539581, 316.420746, 368.982361, 427.592499, 492.616028, & - 564.413452, 643.339905, 729.744141, 823.967834, 926.344910, 1037.201172, & - 1156.853638, 1285.610352, 1423.770142, 1571.622925, 1729.448975, 1897.519287, & - 2076.095947, 2265.431641, 2465.770508, 2677.348145, 2900.391357, 3135.119385, & - 3381.743652, 3640.468262, 3911.490479, 4194.930664, 4490.817383, 4799.149414, & - 5119.895020, 5452.990723, 5798.344727, 6156.074219, 6526.946777, 6911.870605, & - 7311.869141, 7727.412109, 8159.354004, 8608.525391, 9076.400391, 9562.682617, & - 10065.978516, 10584.631836, 11116.662109, 11660.067383, 12211.547852, 12766.873047, & - 13324.668945, 13881.331055, 14432.139648, 14975.615234, 15508.256836, 16026.115234, & - 16527.322266, 17008.789063, 17467.613281, 17901.621094, 18308.433594, 18685.718750, & - 19031.289063, 19343.511719, 19620.042969, 19859.390625, 20059.931641, 20219.664063, & - 20337.863281, 20412.308594, 20442.078125, 20425.718750, 20361.816406, 20249.511719, & - 20087.085938, 19874.025391, 19608.572266, 19290.226563, 18917.460938, 18489.707031, & - 18006.925781, 17471.839844, 16888.687500, 16262.046875, 15596.695313, 14898.453125, & - 14173.324219, 13427.769531, 12668.257813, 11901.339844, 11133.304688, 10370.175781, & - 9617.515625, 8880.453125, 8163.375000, 7470.343750, 6804.421875, 6168.531250, & - 5564.382813, 4993.796875, 4457.375000, 3955.960938, 3489.234375, 3057.265625, & - 2659.140625, 2294.242188, 1961.500000, 1659.476563, 1387.546875, 1143.250000, & - 926.507813, 734.992188, 568.062500, 424.414063, 302.476563, 202.484375, & - 122.101563, 62.781250, 22.835938, 3.757813, 0.000000, 0.000000 / - - data b125/ 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000007, 0.000024, 0.000059, 0.000112, 0.000199, & - 0.000340, 0.000562, 0.000890, 0.001353, 0.001992, 0.002857, & - 0.003971, 0.005378, 0.007133, 0.009261, 0.011806, 0.014816, & - 0.018318, 0.022355, 0.026964, 0.032176, 0.038026, 0.044548, & - 0.051773, 0.059728, 0.068448, 0.077958, 0.088286, 0.099462, & - 0.111505, 0.124448, 0.138313, 0.153125, 0.168910, 0.185689, & - 0.203491, 0.222333, 0.242244, 0.263242, 0.285354, 0.308598, & - 0.332939, 0.358254, 0.384363, 0.411125, 0.438391, 0.466003, & - 0.493800, 0.521619, 0.549301, 0.576692, 0.603648, 0.630036, & - 0.655736, 0.680643, 0.704669, 0.727739, 0.749797, 0.770798, & - 0.790717, 0.809536, 0.827256, 0.843881, 0.859432, 0.873929, & - 0.887408, 0.899900, 0.911448, 0.922096, 0.931881, 0.940860, & - 0.949064, 0.956550, 0.963352, 0.969513, 0.975078, 0.980072, & - 0.984542, 0.988500, 0.991984, 0.995003, 0.997630, 1.000000 / - - select case (km) - - case (24) - - ks = 5 - do k=1,km+1 - ak(k) = a24(k) - bk(k) = b24(k) - enddo - - case (26) - - ks = 7 - do k=1,km+1 - ak(k) = a26(k) - bk(k) = b26(k) - enddo - - case (32) -#ifdef OLD_32 - ks = 13 ! high-res trop_32 setup -#else - ks = 7 -#endif - do k=1,km+1 - ak(k) = a32(k) - bk(k) = b32(k) - enddo - - case (47) -! ks = 27 ! high-res trop-strat - ks = 20 ! Oct 23, 2012 - do k=1,km+1 - ak(k) = a47(k) - bk(k) = b47(k) - enddo - - case (48) - ks = 28 - do k=1,km+1 - ak(k) = a48(k) - bk(k) = b48(k) - enddo - - case (52) - ks = 35 ! pint = 223 - do k=1,km+1 - ak(k) = a52(k) - bk(k) = b52(k) - enddo - - case (54) - ks = 11 ! pint = 109.4 - do k=1,km+1 - ak(k) = a54(k) - bk(k) = b54(k) - enddo - - case (56) - ks = 26 - do k=1,km+1 - ak(k) = a56(k) - bk(k) = b56(k) - enddo - - case (60) - ks = 37 - do k=1,km+1 - ak(k) = a60(k) - bk(k) = b60(k) - enddo - - - case (64) -#ifdef GFSL64 - ks = 23 -#else - ks = 46 -#endif - do k=1,km+1 - ak(k) = a64(k) - bk(k) = b64(k) - enddo -!-->cjg - case (68) - ks = 27 - do k=1,km+1 - ak(k) = a68(k) - bk(k) = b68(k) - enddo - - case (96) - ks = 27 - do k=1,km+1 - ak(k) = a96(k) - bk(k) = b96(k) - enddo -!<--cjg - - case (100) - ks = 38 - do k=1,km+1 - ak(k) = a100(k) - bk(k) = b100(k) - enddo - - case (104) - ks = 73 - do k=1,km+1 - ak(k) = a104(k) - bk(k) = b104(k) - enddo - -#ifndef TEST_GWAVES - case (10) -!-------------------------------------------------- -! Pure sigma-coordinate with uniform spacing in "z" -!-------------------------------------------------- -! - pt = 2000. ! model top pressure (pascal) -! pt = 100. ! 1 mb - press(1) = pt - press(km+1) = p0 - dlnp = (log(p0) - log(pt)) / real(km) - - lnpe = log(press(km+1)) - do k=km,2,-1 - lnpe = lnpe - dlnp - press(k) = exp(lnpe) - enddo - -! Search KS - ks = 0 - do k=1,km - if(press(k) >= pc) then - ks = k-1 - goto 123 - endif - enddo -123 continue - - if(ks /= 0) then - do k=1,ks - ak(k) = press(k) - bk(k) = 0. - enddo - endif - - pint = press(ks+1) - do k=ks+1,km - ak(k) = pint*(press(km)-press(k))/(press(km)-pint) - bk(k) = (press(k) - ak(k)) / press(km+1) - enddo - ak(km+1) = 0. - bk(km+1) = 1. - -! do k=2,km -! bk(k) = real(k-1) / real(km) -! ak(k) = pt * ( 1. - bk(k) ) -! enddo -#endif - -! The following 4 selections are better for non-hydrostatic -! Low top: - case (31) - ptop = 300. - pint = 100.E2 - call var_dz(km, ak, bk, ptop, ks, pint, 1.035) -#ifndef TEST_GWAVES - case (41) - ptop = 100. - pint = 100.E2 - call var_hi(km, ak, bk, ptop, ks, pint, 1.035) -#endif - case (51) - ptop = 100. - pint = 100.E2 - call var_hi(km, ak, bk, ptop, ks, pint, 1.035) -! Mid-top: - case (55) - ptop = 10. - pint = 100.E2 -! call var_dz(km, ak, bk, ptop, ks, pint, 1.035) - call var_hi(km, ak, bk, ptop, ks, pint, 1.035) -#ifdef USE_GFSL63 -! GFS L64 equivalent setting - case (63) - ks = 23 - ptop = a63(1) - pint = a63(ks+1) - do k=1,km+1 - ak(k) = a63(k) - bk(k) = b63(k) - enddo -#else - case (63) - ptop = 1. ! high top - pint = 100.E2 - call var_hi(km, ak, bk, ptop, ks, pint, 1.035) -#endif -! NGGPS_GFS - case (91) - pint = 100.E2 - ptop = 40. - call var_gfs(km, ak, bk, ptop, ks, pint, 1.029) -! call var_gfs(km, ak, bk, ptop, ks, pint, 1.03) - case (95) -! Mid-top settings: - pint = 100.E2 - ptop = 20. - call var_gfs(km, ak, bk, ptop, ks, pint, 1.028) - case (127) - ptop = 1. - pint = 75.E2 - call var_gfs(km, ak, bk, ptop, ks, pint, 1.028) -! IFS-like L125 - case (125) - ks = 33 - ptop = a125(1) - pint = a125(ks+1) - do k=1,km+1 - ak(k) = a125(k) - bk(k) = b125(k) - enddo - case default - -#ifdef TEST_GWAVES -!-------------------------------------------------- -! Pure sigma-coordinate with uniform spacing in "z" -!-------------------------------------------------- - call gw_1d(km, 1000.E2, ak, bk, ptop, 10.E3, pt1) - ks = 0 - pint = ak(1) -#else - -!---------------------------------------------------------------- -! Sigma-coordinate with uniform spacing in sigma and ptop = 1 mb -!---------------------------------------------------------------- - pt = 100. -! One pressure layer - ks = 1 -! pint = pt + 0.5*1.E5/real(km) ! SJL: 20120327 - pint = pt + 1.E5/real(km) - - ak(1) = pt - bk(1) = 0. - ak(2) = pint - bk(2) = 0. - - do k=3,km+1 - bk(k) = real(k-2) / real(km-1) - ak(k) = pint - bk(k)*pint - enddo - ak(km+1) = 0. - bk(km+1) = 1. -#endif - end select - ptop = ak(1) - pint = ak(ks+1) - - end subroutine set_eta -#endif - -!>@brief The subroutine 'set_external_eta' sets 'ptop' (model top) and -!! 'ks' (first level of pure pressure coordinates given the coefficients -!! 'ak' and 'bk' - subroutine set_external_eta(ak, bk, ptop, ks) - implicit none - real, intent(in) :: ak(:) - real, intent(in) :: bk(:) - real, intent(out) :: ptop !< model top (Pa) - integer, intent(out) :: ks !< number of pure p layers - !--- local variables - integer :: k - real :: eps = 1.d-7 - - ptop = ak(1) - ks = 1 - do k = 1, size(bk(:)) - if (bk(k).lt.eps) ks = k - enddo - !--- change ks to layers from levels - ks = ks - 1 - if (is_master()) write(6,*) ' ptop & ks ', ptop, ks - - end subroutine set_external_eta - - - subroutine var_les(km, ak, bk, ptop, ks, pint, s_rate) - implicit none - integer, intent(in):: km - real, intent(in):: ptop - real, intent(in):: s_rate !< between [1. 1.1] - real, intent(out):: ak(km+1), bk(km+1) - real, intent(inout):: pint - integer, intent(out):: ks -! Local - real, parameter:: p00 = 1.E5 - real, dimension(km+1):: ze, pe1, peln, eta - real, dimension(km):: dz, s_fac, dlnp, pm, dp, dk - real ztop, t0, dz0, sum1, tmp1 - real ep, es, alpha, beta, gama - real, parameter:: akap = 2./7. -!---- Tunable parameters: - real:: k_inc = 10 !< number of layers from bottom up to near const dz region - real:: s0 = 0.8 !< lowest layer stretch factor -!----------------------- - real:: s_inc - integer k - - pe1(1) = ptop - peln(1) = log(pe1(1)) - pe1(km+1) = p00 - peln(km+1) = log(pe1(km+1)) - - t0 = 273. - ztop = rdgas/grav*t0*(peln(km+1) - peln(1)) - - s_inc = (1.-s0) / real(k_inc) - s_fac(km) = s0 - - do k=km-1, km-k_inc, -1 - s_fac(k) = s_fac(k+1) + s_inc - enddo - - s_fac(km-k_inc-1) = 0.5*(s_fac(km-k_inc) + s_rate) - - do k=km-k_inc-2, 5, -1 - s_fac(k) = s_rate * s_fac(k+1) - enddo - - s_fac(4) = 0.5*(1.1+s_rate)*s_fac(5) - s_fac(3) = 1.1 *s_fac(4) - s_fac(2) = 1.1 *s_fac(3) - s_fac(1) = 1.1 *s_fac(2) - - sum1 = 0. - do k=1,km - sum1 = sum1 + s_fac(k) - enddo - - dz0 = ztop / sum1 - - do k=1,km - dz(k) = s_fac(k) * dz0 - enddo - - ze(km+1) = 0. - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo - -! Re-scale dz with the stretched ztop - do k=1,km - dz(k) = dz(k) * (ztop/ze(1)) - enddo - - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo -! ze(1) = ztop - - if ( is_master() ) then - write(*,*) 'var_les: computed model top (m)=', ztop, ' bottom/top dz=', dz(km), dz(1) -! do k=1,km -! write(*,*) k, s_fac(k) -! enddo - endif - - call sm1_edge(1, 1, 1, 1, km, 1, 1, ze, 2) - -! Given z --> p - do k=1,km - dz(k) = ze(k) - ze(k+1) - dlnp(k) = grav*dz(k) / (rdgas*t0) - !write(*,*) k, dz(k) - enddo - do k=2,km - peln(k) = peln(k-1) + dlnp(k-1) - pe1(k) = exp(peln(k)) - enddo - - -! Pe(k) = ak(k) + bk(k) * PS -! Locate pint and KS - ks = 0 - do k=2,km - if ( pint < pe1(k)) then - ks = k-1 - exit - endif - enddo - if ( is_master() ) then - write(*,*) 'For (input) PINT=', 0.01*pint, ' KS=', ks, 'pint(computed)=', 0.01*pe1(ks+1) - endif - pint = pe1(ks+1) - - do k=1,km+1 - eta(k) = pe1(k) / pe1(km+1) - enddo - - ep = eta(ks+1) - es = eta(km) -! es = 1. - alpha = (ep**2-2.*ep*es) / (es-ep)**2 - beta = 2.*ep*es**2 / (es-ep)**2 - gama = -(ep*es)**2 / (es-ep)**2 - -! Pure pressure: - do k=1,ks+1 - ak(k) = eta(k)*1.e5 - bk(k) = 0. - enddo - - do k=ks+2, km - ak(k) = alpha*eta(k) + beta + gama/eta(k) - ak(k) = ak(k)*1.e5 - enddo - ak(km+1) = 0. - - do k=ks+2, km - bk(k) = (pe1(k) - ak(k))/pe1(km+1) - enddo - bk(km+1) = 1. - - if ( is_master() ) then - ! write(*,*) 'KS=', ks, 'PINT (mb)=', pint/100. - ! do k=1,km - ! pm(k) = 0.5*(pe1(k)+pe1(k+1))/100. - ! write(*,*) k, pm(k), dz(k) - ! enddo - tmp1 = ak(ks+1) - do k=ks+1,km - tmp1 = max(tmp1, (ak(k)-ak(k+1))/max(1.E-5, (bk(k+1)-bk(k))) ) - enddo - write(*,*) 'Hybrid Sigma-P: minimum allowable surface pressure (hpa)=', tmp1/100. - write(*,800) (pm(k), k=km,1,-1) - endif - - do k=1,km - dp(k) = (pe1(k+1) - pe1(k))/100. - dk(k) = pe1(k+1)**akap - pe1(k)**akap - enddo - -800 format(1x,5(1x,f9.4)) - - end subroutine var_les - - - - subroutine var_gfs(km, ak, bk, ptop, ks, pint, s_rate) - integer, intent(in):: km - real, intent(in):: ptop - real, intent(in):: s_rate !< between [1. 1.1] - real, intent(out):: ak(km+1), bk(km+1) - real, intent(inout):: pint - integer, intent(out):: ks -! Local - real, parameter:: p00 = 1.E5 - real, dimension(km+1):: ze, pe1, peln, eta - real, dimension(km):: dz, s_fac, dlnp - real ztop, t0, dz0, sum1, tmp1 - real ep, es, alpha, beta, gama -!---- Tunable parameters: - integer:: k_inc = 25 !< number of layers from bottom up to near const dz region - real:: s0 = 0.13 !< lowest layer stretch factor -!----------------------- - real:: s_inc - integer k - - pe1(1) = ptop - peln(1) = log(pe1(1)) - pe1(km+1) = p00 - peln(km+1) = log(pe1(km+1)) - - t0 = 270. - ztop = rdgas/grav*t0*(peln(km+1) - peln(1)) - - s_inc = (1.-s0) / real(k_inc) - s_fac(km) = s0 - - do k=km-1, km-k_inc, -1 - s_fac(k) = s_fac(k+1) + s_inc - enddo - - do k=km-k_inc-1, 9, -1 - s_fac(k) = s_rate * s_fac(k+1) - enddo - s_fac(8) = 0.5*(1.1+s_rate)*s_fac(9) - s_fac(7) = 1.10*s_fac(8) - s_fac(6) = 1.15*s_fac(7) - s_fac(5) = 1.20*s_fac(6) - s_fac(4) = 1.26*s_fac(5) - s_fac(3) = 1.33*s_fac(4) - s_fac(2) = 1.41*s_fac(3) - s_fac(1) = 1.60*s_fac(2) - - sum1 = 0. - do k=1,km - sum1 = sum1 + s_fac(k) - enddo - - dz0 = ztop / sum1 - - do k=1,km - dz(k) = s_fac(k) * dz0 - enddo - - ze(km+1) = 0. - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo - -! Re-scale dz with the stretched ztop - do k=1,km - dz(k) = dz(k) * (ztop/ze(1)) - enddo - - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo -! ze(1) = ztop - - if ( is_master() ) then - write(*,*) 'var_gfs: computed model top (m)=', ztop*0.001, ' bottom/top dz=', dz(km), dz(1) -! do k=1,km -! write(*,*) k, s_fac(k) -! enddo - endif - -! call sm1_edge(1, 1, 1, 1, km, 1, 1, ze, 3) - -! Given z --> p - do k=1,km - dz(k) = ze(k) - ze(k+1) - dlnp(k) = grav*dz(k) / (rdgas*t0) - enddo - do k=2,km - peln(k) = peln(k-1) + dlnp(k-1) - pe1(k) = exp(peln(k)) - enddo - -! Pe(k) = ak(k) + bk(k) * PS -! Locate pint and KS - ks = 0 - do k=2,km - if ( pint < pe1(k)) then - ks = k-1 - exit - endif - enddo - if ( is_master() ) then - write(*,*) 'For (input) PINT=', 0.01*pint, ' KS=', ks, 'pint(computed)=', 0.01*pe1(ks+1) - write(*,*) 'ptop =', ptop - endif - pint = pe1(ks+1) - -#ifdef NO_UKMO_HB - do k=1,ks+1 - ak(k) = pe1(k) - bk(k) = 0. - enddo - - do k=ks+2,km+1 - bk(k) = (pe1(k) - pint) / (pe1(km+1)-pint) ! bk == sigma - ak(k) = pe1(k) - bk(k) * pe1(km+1) - enddo - bk(km+1) = 1. - ak(km+1) = 0. -#else -! Problematic for non-hydrostatic - do k=1,km+1 - eta(k) = pe1(k) / pe1(km+1) - enddo - - ep = eta(ks+1) - es = eta(km) -! es = 1. - alpha = (ep**2-2.*ep*es) / (es-ep)**2 - beta = 2.*ep*es**2 / (es-ep)**2 - gama = -(ep*es)**2 / (es-ep)**2 - -! Pure pressure: - do k=1,ks+1 - ak(k) = eta(k)*1.e5 - bk(k) = 0. - enddo - - do k=ks+2, km - ak(k) = alpha*eta(k) + beta + gama/eta(k) - ak(k) = ak(k)*1.e5 - enddo - ak(km+1) = 0. - - do k=ks+2, km - bk(k) = (pe1(k) - ak(k))/pe1(km+1) - enddo - bk(km+1) = 1. -#endif - - if ( is_master() ) then - write(*,*) 'KS=', ks, 'PINT (mb)=', pint/100. - do k=1,km - write(*,*) k, 0.5*(pe1(k)+pe1(k+1))/100., dz(k) - enddo - tmp1 = ak(ks+1) - do k=ks+1,km - tmp1 = max(tmp1, (ak(k)-ak(k+1))/max(1.E-5, (bk(k+1)-bk(k))) ) - enddo - write(*,*) 'Hybrid Sigma-P: minimum allowable surface pressure (hpa)=', tmp1/100. - endif - - end subroutine var_gfs - - subroutine var_hi(km, ak, bk, ptop, ks, pint, s_rate) - integer, intent(in):: km - real, intent(in):: ptop - real, intent(in):: s_rate !< between [1. 1.1] - real, intent(out):: ak(km+1), bk(km+1) - real, intent(inout):: pint - integer, intent(out):: ks -! Local - real, parameter:: p00 = 1.E5 - real, dimension(km+1):: ze, pe1, peln, eta - real, dimension(km):: dz, s_fac, dlnp - real ztop, t0, dz0, sum1, tmp1 - real ep, es, alpha, beta, gama -!---- Tunable parameters: - integer:: k_inc = 15 ! 1.E-8 ) then - pf(1) = (ph(2) - ph(1)) / log(ph(2)/ph(1)) - else - pf(1) = (ph(2) - ph(1)) * kappa/(kappa+1.) - endif - - do k=2,npz - pf(k) = (ph(k+1) - ph(k)) / log(ph(k+1)/ph(k)) - enddo - - end subroutine get_eta_level - - - - subroutine compute_dz(km, ztop, dz) - - integer, intent(in):: km - real, intent(in):: ztop ! try 50.E3 - real, intent(out):: dz(km) -!------------------------------ - real ze(km+1), dzt(km) - integer k - - -! ztop = 30.E3 - dz(1) = ztop / real(km) - dz(km) = 0.5*dz(1) - - do k=2,km-1 - dz(k) = dz(1) - enddo - -! Top: - dz(1) = 2.*dz(2) - - ze(km+1) = 0. - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo - - if ( is_master() ) then - write(*,*) 'Hybrid_z: dz, zm' - do k=1,km - dzt(k) = 0.5*(ze(k)+ze(k+1)) / 1000. - write(*,*) k, dz(k), dzt(k) - enddo - endif - - end subroutine compute_dz - - subroutine compute_dz_var(km, ztop, dz) - - integer, intent(in):: km - real, intent(in):: ztop ! try 50.E3 - real, intent(out):: dz(km) -!------------------------------ - real, parameter:: s_rate = 1.0 - real ze(km+1) - real s_fac(km) - real sum1, dz0 - integer k - - s_fac(km ) = 0.125 - s_fac(km-1) = 0.20 - s_fac(km-2) = 0.30 - s_fac(km-3) = 0.40 - s_fac(km-4) = 0.50 - s_fac(km-5) = 0.60 - s_fac(km-6) = 0.70 - s_fac(km-7) = 0.80 - s_fac(km-8) = 0.90 - s_fac(km-9) = 1. - - do k=km-10, 9, -1 - s_fac(k) = s_rate * s_fac(k+1) - enddo - - s_fac(8) = 1.05*s_fac(9) - s_fac(7) = 1.1 *s_fac(8) - s_fac(6) = 1.15*s_fac(7) - s_fac(5) = 1.2 *s_fac(6) - s_fac(4) = 1.3 *s_fac(5) - s_fac(3) = 1.4 *s_fac(4) - s_fac(2) = 1.5 *s_fac(3) - s_fac(1) = 1.6 *s_fac(2) - - sum1 = 0. - do k=1,km - sum1 = sum1 + s_fac(k) - enddo - - dz0 = ztop / sum1 - - do k=1,km - dz(k) = s_fac(k) * dz0 - enddo - - ze(1) = ztop - ze(km+1) = 0. - do k=km,2,-1 - ze(k) = ze(k+1) + dz(k) - enddo - -! Re-scale dz with the stretched ztop - do k=1,km - dz(k) = dz(k) * (ztop/ze(1)) - enddo - - do k=km,2,-1 - ze(k) = ze(k+1) + dz(k) - enddo - - call sm1_edge(1, 1, 1, 1, km, 1, 1, ze, 2) - - do k=1,km - dz(k) = ze(k) - ze(k+1) - enddo - - end subroutine compute_dz_var - - subroutine compute_dz_L32(km, ztop, dz) - - integer, intent(in):: km - real, intent(out):: dz(km) - real, intent(out):: ztop ! try 50.E3 -!------------------------------ - real dzt(km) - real ze(km+1) - real dz0, dz1, dz2 - real z0, z1, z2 - integer k, k0, k1, k2, n - -!------------------- - k2 = 8 - z2 = 30.E3 -!------------------- - k1 = 21 - z1 = 10.0E3 -!------------------- - k0 = 2 - z0 = 0. - dz0 = 75. ! meters -!------------------- -! Treat the surface layer as a special layer - ze(1) = z0 - dz(1) = dz0 - - ze(2) = dz(1) - dz0 = 1.5*dz0 - dz(2) = dz0 - - ze(3) = ze(2) + dz(2) - - dz1 = 2.*(z1-ze(3) - k1*dz0) / (k1*(k1-1)) - - do k=k0+1,k0+k1 - dz(k) = dz0 + (k-k0)*dz1 - ze(k+1) = ze(k) + dz(k) - enddo - - dz0 = dz(k1+k0) - dz2 = 2.*(z2-ze(k0+k1+1)-k2*dz0) / (k2*(k2-1)) - - do k=k0+k1+1,k0+k1+k2 - dz(k) = dz0 + (k-k0-k1)*dz2 - ze(k+1) = ze(k) + dz(k) - enddo - - dz(km) = 2.*dz(km-1) - ztop = ze(km) + dz(km) - ze(km+1) = ze(km) + dz(km) - - call zflip (dz, 1, km) - - ze(km+1) = 0. - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo - -! if ( is_master() ) then -! write(*,*) 'Hybrid_z: dz, zm' -! do k=1,km -! dzt(k) = 0.5*(ze(k)+ze(k+1)) / 1000. -! write(*,*) k, dz(k), dzt(k) -! enddo -! endif - - end subroutine compute_dz_L32 - - subroutine compute_dz_L101(km, ztop, dz) - - integer, intent(in):: km ! km==101 - real, intent(out):: dz(km) - real, intent(out):: ztop ! try 30.E3 -!------------------------------ - real ze(km+1) - real dz0, dz1 - real:: stretch_f = 1.16 - integer k, k0, k1 - - k1 = 2 - k0 = 25 - dz0 = 40. ! meters - - ze(km+1) = 0. - - do k=km, k0, -1 - dz(k) = dz0 - ze(k) = ze(k+1) + dz(k) - enddo - - do k=k0+1, k1, -1 - dz(k) = stretch_f * dz(k+1) - ze(k) = ze(k+1) + dz(k) - enddo - - dz(1) = 4.0*dz(2) - ze(1) = ze(2) + dz(1) - ztop = ze(1) - - if ( is_master() ) then - write(*,*) 'Hybrid_z: dz, ze' - do k=1,km - write(*,*) k, 0.001*dz(k), 0.001*ze(k) - enddo -! ztop (km) = 20.2859154 - write(*,*) 'ztop (km) =', ztop * 0.001 - endif - - end subroutine compute_dz_L101 - - subroutine set_hybrid_z(is, ie, js, je, ng, km, ztop, dz, rgrav, hs, ze, dz3) - - integer, intent(in):: is, ie, js, je, ng, km - real, intent(in):: rgrav, ztop - real, intent(in):: dz(km) !< Reference vertical resolution for zs=0 - real, intent(in):: hs(is-ng:ie+ng,js-ng:je+ng) - real, intent(inout):: ze(is:ie,js:je,km+1) - real, optional, intent(out):: dz3(is-ng:ie+ng,js-ng:je+ng,km) -! local - logical:: filter_xy = .false. - real, allocatable:: delz(:,:,:) - integer ntimes - real zint - real:: z1(is:ie,js:je) - real:: z(km+1) - real sig, z_rat - integer ks(is:ie,js:je) - integer i, j, k, ks_min, kint - - z(km+1) = 0. - do k=km,1,-1 - z(k) = z(k+1) + dz(k) - enddo - - do j=js,je - do i=is,ie - ze(i,j, 1) = ztop - ze(i,j,km+1) = hs(i,j) * rgrav - enddo - enddo - - do k=2,km - do j=js,je - do i=is,ie - ze(i,j,k) = z(k) - enddo - enddo - enddo - -! Set interface: -#ifndef USE_VAR_ZINT - zint = 12.0E3 - ntimes = 2 - kint = 2 - do k=2,km - if ( z(k)<=zint ) then - kint = k - exit - endif - enddo - - if ( is_master() ) write(*,*) 'Z_coord interface set at k=',kint, ' ZE=', z(kint) - - do j=js,je - do i=is,ie - z_rat = (ze(i,j,kint)-ze(i,j,km+1)) / (z(kint)-z(km+1)) - do k=km,kint+1,-1 - ze(i,j,k) = ze(i,j,k+1) + dz(k)*z_rat - enddo -!-------------------------------------- -! Apply vertical smoother locally to dz -!-------------------------------------- - call sm1_edge(is, ie, js, je, km, i, j, ze, ntimes) - enddo - enddo -#else -! ZINT is a function of local terrain - ntimes = 4 - do j=js,je - do i=is,ie - z1(i,j) = dim(ze(i,j,km+1), 2500.) + 7500. - enddo - enddo - - ks_min = km - do j=js,je - do i=is,ie - do k=km,2,-1 - if ( z(k)>=z1(i,j) ) then - ks(i,j) = k - ks_min = min(ks_min, k) - go to 555 - endif - enddo -555 continue - enddo - enddo - - do j=js,je - do i=is,ie - kint = ks(i,j) + 1 - z_rat = (ze(i,j,kint)-ze(i,j,km+1)) / (z(kint)-z(km+1)) - do k=km,kint+1,-1 - ze(i,j,k) = ze(i,j,k+1) + dz(k)*z_rat - enddo -!-------------------------------------- -! Apply vertical smoother locally to dz -!-------------------------------------- - call sm1_edge(is, ie, js, je, km, i, j, ze, ntimes) - enddo - enddo -#endif - -#ifdef DEV_ETA - if ( filter_xy ) then - allocate (delz(isd:ied, jsd:jed, km) ) - ntimes = 2 - do k=1,km - do j=js,je - do i=is,ie - delz(i,j,k) = ze(i,j,k+1) - ze(i,j,k) - enddo - enddo - enddo - call del2_cubed(delz, 0.2*da_min, npx, npy, km, ntimes) - do k=km,2,-1 - do j=js,je - do i=is,ie - ze(i,j,k) = ze(i,j,k+1) - delz(i,j,k) - enddo - enddo - enddo - deallocate ( delz ) - endif -#endif - if ( present(dz3) ) then - do k=1,km - do j=js,je - do i=is,ie - dz3(i,j,k) = ze(i,j,k+1) - ze(i,j,k) - enddo - enddo - enddo - endif - - end subroutine set_hybrid_z - - - subroutine sm1_edge(is, ie, js, je, km, i, j, ze, ntimes) - integer, intent(in):: is, ie, js, je, km - integer, intent(in):: ntimes, i, j - real, intent(inout):: ze(is:ie,js:je,km+1) -! local: - real, parameter:: df = 0.25 - real dz(km) - real flux(km+1) - integer k, n, k1, k2 - - k2 = km-1 - do k=1,km - dz(k) = ze(i,j,k+1) - ze(i,j,k) - enddo - - do n=1,ntimes - k1 = 2 + (ntimes-n) - - flux(k1 ) = 0. - flux(k2+1) = 0. - do k=k1+1,k2 - flux(k) = df*(dz(k) - dz(k-1)) - enddo - - do k=k1,k2 - dz(k) = dz(k) - flux(k) + flux(k+1) - enddo - enddo - - do k=km,1,-1 - ze(i,j,k) = ze(i,j,k+1) - dz(k) - enddo - - end subroutine sm1_edge - - - - subroutine gw_1d(km, p0, ak, bk, ptop, ztop, pt1) - integer, intent(in):: km - real, intent(in):: p0, ztop - real, intent(inout):: ptop - real, intent(inout):: ak(km+1), bk(km+1) - real, intent(out):: pt1(km) -! Local - logical:: isothermal - real, dimension(km+1):: ze, pe1, pk1 - real, dimension(km):: dz1 - real t0, n2, s0 - integer k - -! Set up vertical coordinare with constant del-z spacing: - isothermal = .false. - t0 = 300. - - if ( isothermal ) then - n2 = grav**2/(cp_air*t0) - else - n2 = 0.0001 - endif - - s0 = grav*grav / (cp_air*n2) - - ze(km+1) = 0. - do k=km,1,-1 - dz1(k) = ztop / real(km) - ze(k) = ze(k+1) + dz1(k) - enddo - -! Given z --> p - do k=1,km+1 - pe1(k) = p0*( (1.-s0/t0) + s0/t0*exp(-n2*ze(k)/grav) )**(1./kappa) - enddo - - ptop = pe1(1) -! if ( is_master() ) write(*,*) 'GW_1D: computed model top (pa)=', ptop - -! Set up "sigma" coordinate - ak(1) = pe1(1) - bk(1) = 0. - do k=2,km - bk(k) = (pe1(k) - pe1(1)) / (pe1(km+1)-pe1(1)) ! bk == sigma - ak(k) = pe1(1)*(1.-bk(k)) - enddo - ak(km+1) = 0. - bk(km+1) = 1. - - do k=1,km+1 - pk1(k) = pe1(k) ** kappa - enddo - -! Compute volume mean potential temperature with hydrostatic eqn: - do k=1,km - pt1(k) = grav*dz1(k) / ( cp_air*(pk1(k+1)-pk1(k)) ) - enddo - - end subroutine gw_1d - - - - subroutine zflip(q, im, km) - integer, intent(in):: im, km - real, intent(inout):: q(im,km) -!--- - integer i, k - real qtmp - - do i = 1, im - do k = 1, (km+1)/2 - qtmp = q(i,k) - q(i,k) = q(i,km+1-k) - q(i,km+1-k) = qtmp - end do - end do - - end subroutine zflip - -end module fv_eta_mod diff --git a/tools/fv_eta.F90_NAM_lyrs b/tools/fv_eta.F90_NAM_lyrs deleted file mode 100644 index 7697c385c..000000000 --- a/tools/fv_eta.F90_NAM_lyrs +++ /dev/null @@ -1,3079 +0,0 @@ -!*********************************************************************** -!* GNU Lesser General Public License -!* -!* This file is part of the FV3 dynamical core. -!* -!* The FV3 dynamical core is free software: you can redistribute it -!* and/or modify it under the terms of the -!* GNU Lesser General Public License as published by the -!* Free Software Foundation, either version 3 of the License, or -!* (at your option) any later version. -!* -!* The FV3 dynamical core is distributed in the hope that it will be -!* useful, but WITHOUT ANYWARRANTY; 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 Lesser General Public -!* License along with the FV3 dynamical core. -!* If not, see . -!*********************************************************************** - -!>@brief The module 'fv_eta' contains routine to set up the reference -!! (Eulerian) pressure coordinate - -module fv_eta_mod - -! -! -! -! -! -! -! -! -! -! -! -! -! -! -! -! -! -!
Module NameFunctions Included
constants_modkappa, grav, cp_air, rdgas
fv_mp_modis_master
mpp_modmpp_error, FATAL
- - use constants_mod, only: kappa, grav, cp_air, rdgas - use fv_mp_mod, only: is_master - use mpp_mod, only: FATAL, mpp_error - implicit none - private - public set_eta, set_external_eta, get_eta_level, compute_dz_var, & - compute_dz_L32, compute_dz_L101, set_hybrid_z, compute_dz, & - gw_1d, sm1_edge, hybrid_z_dz - - contains - -!!!NOTE: USE_VAR_ETA not used in fvGFS -#ifdef USE_VAR_ETA - subroutine set_eta(km, ks, ptop, ak, bk) -! This is the easy to use version of the set_eta - integer, intent(in):: km ! vertical dimension - integer, intent(out):: ks ! number of pure p layers - real:: a60(61),b60(61) -! Thfollowing L63 setting is the same as NCEP GFS's L64 except the top -! 3 layers - data a60/300.0000, 430.00000, 558.00000, & - 700.00000, 863.05803, 1051.07995, & - 1265.75194, 1510.71101, 1790.05098, & - 2108.36604, 2470.78817, 2883.03811, & - 3351.46002, 3883.05187, 4485.49315, & - 5167.14603, 5937.04991, 6804.87379, & - 7780.84698, 8875.64338, 10100.20534, & - 11264.35673, 12190.64366, 12905.42546, & - 13430.87867, 13785.88765, 13986.77987, & - 14047.96335, 13982.46770, 13802.40331, & - 13519.33841, 13144.59486, 12689.45608, & - 12165.28766, 11583.57006, 10955.84778, & - 10293.60402, 9608.08306, 8910.07678, & - 8209.70131, 7516.18560, 6837.69250, & - 6181.19473, 5552.39653, 4955.72632, & - 4394.37629, 3870.38682, 3384.76586, & - 2937.63489, 2528.37666, 2155.78385, & - 1818.20722, 1513.68173, 1240.03585, & - 994.99144, 776.23591, 581.48797, & - 408.53400, 255.26520, 119.70243, 0. / - - data b60/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00201, 0.00792, 0.01755, & - 0.03079, 0.04751, 0.06761, & - 0.09097, 0.11746, 0.14690, & - 0.17911, 0.21382, 0.25076, & - 0.28960, 0.32994, 0.37140, & - 0.41353, 0.45589, 0.49806, & - 0.53961, 0.58015, 0.61935, & - 0.65692, 0.69261, 0.72625, & - 0.75773, 0.78698, 0.81398, & - 0.83876, 0.86138, 0.88192, & - 0.90050, 0.91722, 0.93223, & - 0.94565, 0.95762, 0.96827, & - 0.97771, 0.98608, 0.99347, 1./ - real, intent(out):: ak(km+1) - real, intent(out):: bk(km+1) - real, intent(out):: ptop ! model top (Pa) - real pint, stretch_fac - integer k - real :: s_rate = -1.0 ! dummy value to not use var_les - - pint = 100.E2 - -!- Notes --------------------------------- -! low-top: ptop = 100. ! ~45 km -! mid-top: ptop = 10. ! ~60 km -! hi -top: ptop = 1. ! ~80 km -!----------------------------------------- - select case (km) - -! Optimal number = 8 * N -1 (for 8 openMP threads) -! 16 * M -1 (for 16 openMP threads) - -#ifdef HIWPP -#ifdef SUPER_K - case (20) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (24) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (30) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (40) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (50) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (60) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 - case (80) - ptop = 56.e2 - pint = ptop - stretch_fac = 1.03 -#else - case (30) ! For Baroclinic Instability Test - ptop = 2.26e2 - pint = 250.E2 - stretch_fac = 1.03 - case (40) - ptop = 50.e2 ! For super cell test - pint = 300.E2 - stretch_fac = 1.03 - case (50) ! Mountain waves? - ptop = 30.e2 - stretch_fac = 1.025 - case (60) ! For Baroclinic Instability Test -#ifdef GFSL60 - ks = 20 - ptop = a60(1) - pint = a60(ks+1) - do k=1,km+1 - ak(k) = a60(k) - bk(k) = b60(k) - enddo -#else - ptop = 3.e2 -! pint = 250.E2 - pint = 300.E2 ! revised for Moist test - stretch_fac = 1.03 -#endif -#endif - case (64) -!!! ptop = 3.e2 - ptop = 2.0e2 - pint = 300.E2 - stretch_fac = 1.03 -#else -! *Very-low top: for idealized super-cell simulation: - case (50) - ptop = 50.e2 - pint = 250.E2 - stretch_fac = 1.03 - case (60) - ptop = 40.e2 - pint = 250.E2 - stretch_fac = 1.03 - case (90) ! super-duper cell - ptop = 40.e2 - stretch_fac = 1.025 -#endif -! Low-top: - case (31) ! N = 4, M=2 - ptop = 100. - stretch_fac = 1.035 - case (32) ! N = 4, M=2 - ptop = 100. - stretch_fac = 1.035 - case (39) ! N = 5 - ptop = 100. - stretch_fac = 1.035 - case (41) - ptop = 100. - stretch_fac = 1.035 - case (47) ! N = 6, M=3 - ptop = 100. - stretch_fac = 1.035 - case (51) - ptop = 100. - stretch_fac = 1.03 - case (52) ! very low top - ptop = 30.e2 ! for special DPM RCE experiments - stretch_fac = 1.03 -! Mid-top: - case (55) ! N = 7 - ptop = 10. - stretch_fac = 1.035 -! Hi-top: - case (63) ! N = 8, M=4 - ptop = 1. - ! c360 or c384 - stretch_fac = 1.035 - case (71) ! N = 9 - ptop = 1. - stretch_fac = 1.03 - case (79) ! N = 10, M=5 - ptop = 1. - stretch_fac = 1.03 - case (127) ! N = 10, M=5 - ptop = 1. - stretch_fac = 1.03 - case (151) - ptop = 75.e2 - pint = 500.E2 - s_rate = 1.01 - case default - ptop = 1. - stretch_fac = 1.03 - end select - -#ifdef MOUNTAIN_WAVES - call mount_waves(km, ak, bk, ptop, ks, pint) -#else - if (s_rate > 0.) then - call var_les(km, ak, bk, ptop, ks, pint, s_rate) - else - if ( km > 79 ) then - call var_hi2(km, ak, bk, ptop, ks, pint, stretch_fac) - elseif (km==5 .or. km==10 ) then -! Equivalent Shallow Water: for NGGPS, variable-resolution testing - ptop = 500.e2 - ks = 0 - do k=1,km+1 - bk(k) = real(k-1) / real (km) - ak(k) = ptop*(1.-bk(k)) - enddo - else -#ifndef GFSL60 - call var_hi(km, ak, bk, ptop, ks, pint, stretch_fac) -#endif - endif -#endif - endif - - ptop = ak(1) - pint = ak(ks+1) - - end subroutine set_eta - - subroutine mount_waves(km, ak, bk, ptop, ks, pint) - integer, intent(in):: km - real, intent(out):: ak(km+1), bk(km+1) - real, intent(out):: ptop, pint - integer, intent(out):: ks -! Local - real, parameter:: p00 = 1.E5 - real, dimension(km+1):: ze, pe1, peln, eta - real, dimension(km):: dz, dlnp - real ztop, t0, dz0, sum1, tmp1 - real ep, es, alpha, beta, gama, s_fac - integer k, k500 - - pint = 300.e2 -! s_fac = 1.05 -! dz0 = 500. - if ( km <= 60 ) then - s_fac = 1.0 - dz0 = 500. - else - s_fac = 1. - dz0 = 250. - endif - -! Basic parameters for HIWPP mountain waves - t0 = 300. -! ztop = 20.0e3; 500-m resolution in halft of the vertical domain -! ztop = real(km-1)*500. -!----------------------- -! Compute temp ptop based on isothermal atm -! ptop = p00*exp(-grav*ztop/(rdgas*t0)) - -! Lowest half has constant resolution - ze(km+1) = 0. - do k=km, km-19, -1 - ze(k) = ze(k+1) + dz0 - enddo - -! Stretching from 10-km and up: - do k=km-20, 3, -1 - dz0 = s_fac * dz0 - ze(k) = ze(k+1) + dz0 - enddo - ze(2) = ze(3) + sqrt(2.)*dz0 - ze(1) = ze(2) + 2.0*dz0 - -! call sm1_edge(1, 1, 1, 1, km, 1, 1, ze, 1) - -! Given z --> p - do k=1,km - dz(k) = ze(k) - ze(k+1) - dlnp(k) = grav*dz(k) / (rdgas*t0) - enddo - - pe1(km+1) = p00 - peln(km+1) = log(p00) - do k=km,1,-1 - peln(k) = peln(k+1) - dlnp(k) - pe1(k) = exp(peln(k)) - enddo - -! Comnpute new ptop - ptop = pe1(1) - -! Pe(k) = ak(k) + bk(k) * PS -! Locate pint and KS - ks = 0 - do k=2,km - if ( pint < pe1(k)) then - ks = k-1 - exit - endif - enddo - - if ( is_master() ) then - write(*,*) 'For (input) PINT=', 0.01*pint, ' KS=', ks, 'pint(computed)=', 0.01*pe1(ks+1) - write(*,*) 'Modified ptop =', ptop, ' ztop=', ze(1)/1000. - do k=1,km - write(*,*) k, 'ze =', ze(k)/1000. - enddo - endif - pint = pe1(ks+1) - -#ifdef NO_UKMO_HB - do k=1,ks+1 - ak(k) = pe1(k) - bk(k) = 0. - enddo - - do k=ks+2,km+1 - bk(k) = (pe1(k) - pint) / (pe1(km+1)-pint) ! bk == sigma - ak(k) = pe1(k) - bk(k) * pe1(km+1) - enddo - bk(km+1) = 1. - ak(km+1) = 0. -#else -! Problematic for non-hydrostatic - do k=1,km+1 - eta(k) = pe1(k) / pe1(km+1) - enddo - ep = eta(ks+1) - es = eta(km) -! es = 1. - alpha = (ep**2-2.*ep*es) / (es-ep)**2 - beta = 2.*ep*es**2 / (es-ep)**2 - gama = -(ep*es)**2 / (es-ep)**2 - -! Pure pressure: - do k=1,ks+1 - ak(k) = eta(k)*1.e5 - bk(k) = 0. - enddo - - do k=ks+2, km - ak(k) = alpha*eta(k) + beta + gama/eta(k) - ak(k) = ak(k)*1.e5 - enddo - ak(km+1) = 0. - - do k=ks+2, km - bk(k) = (pe1(k) - ak(k))/pe1(km+1) - enddo - bk(km+1) = 1. -#endif - - if ( is_master() ) then - tmp1 = ak(ks+1) - do k=ks+1,km - tmp1 = max(tmp1, (ak(k)-ak(k+1))/max(1.E-5, (bk(k+1)-bk(k))) ) - enddo - write(*,*) 'Hybrid Sigma-P: minimum allowable surface pressure (hpa)=', tmp1/100. - endif - - end subroutine mount_waves - -#else - !>This is the version of set_eta used in fvGFS and AM4 - !>@note 01/2018: 'set_eta' is being cleaned up. - subroutine set_eta(km, ks, ptop, ak, bk) - integer, intent(in):: km !< vertical dimension - integer, intent(out):: ks !< number of pure p layers - real, intent(out):: ak(km+1) - real, intent(out):: bk(km+1) - real, intent(out):: ptop !< model top (Pa) -! local - real a24(25),b24(25) !< GFDL AM2L24 - real a26(27),b26(27) !< Jablonowski & Williamson 26-level - real a32(33),b32(33) - real a32w(33),b32w(33) - real a47(48),b47(48) - real a48(49),b48(49) - real a52(53),b52(53) - real a54(55),b54(55) - real a56(57),b56(57) - real a60(61),b60(61) - real a63(64),b63(64) - real a64(65),b64(65) - real a68(69),b68(69) !< cjg: grid with enhanced PBL resolution - real a96(97),b96(97) !< cjg: grid with enhanced PBL resolution - real a100(101),b100(101) - real a104(105),b104(105) - real a125(126),b125(126) - - real:: p0=1000.E2 - real:: pc=200.E2 - - real pt, pint, lnpe, dlnp - real press(km+1), pt1(km) - integer k - -! Definition: press(i,j,k) = ak(k) + bk(k) * ps(i,j) - -!----------------------------------------------- -! GFDL AM2-L24: modified by SJL at the model top -!----------------------------------------------- -! data a24 / 100.0000, 1050.0000, 3474.7942, 7505.5556, 12787.2428, & - data a24 / 100.0000, 903.4465, 3474.7942, 7505.5556, 12787.2428, & - 19111.3683, 21854.9274, 22884.1866, 22776.3058, 21716.1604, & - 20073.2963, 18110.5123, 16004.7832, 13877.6253, 11812.5452, & - 9865.8840, 8073.9726, 6458.0834, 5027.9899, 3784.6085, & - 2722.0086, 1828.9752, 1090.2396, 487.4595, 0.0000 / - - data b24 / 0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0000000, & - 0.0000000, 0.0435679, 0.1102275, 0.1922249, 0.2817656, & - 0.3694997, 0.4532348, 0.5316253, 0.6038733, 0.6695556, & - 0.7285176, 0.7808017, 0.8265992, 0.8662148, 0.9000406, & - 0.9285364, 0.9522140, 0.9716252, 0.9873523, 1.0000000 / - -! Jablonowski & Williamson 26-level setup - data a26 / 219.4067, 489.5209, 988.2418, 1805.2010, 2983.7240, 4462.3340, & - 6160.5870, 7851.2430, 7731.2710, 7590.1310, 7424.0860, & - 7228.7440, 6998.9330, 6728.5740, 6410.5090, 6036.3220, & - 5596.1110, 5078.2250, 4468.9600, 3752.1910, 2908.9490, & - 2084.739, 1334.443, 708.499, 252.1360, 0.0, 0.0 / - - data b26 / 0.0, 0.0, 0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0000000,& - 0.0000000, 0.01505309, 0.03276228, 0.05359622, 0.07810627, & - 0.1069411, 0.1408637, 0.1807720, 0.2277220, 0.2829562, & - 0.3479364, 0.4243822, 0.5143168, 0.6201202, 0.7235355, & - 0.8176768, 0.8962153, 0.9534761, 0.9851122, 1.0000000 / - - -! High-resolution troposphere setup -#ifdef OLD_32 -! Revised Apr 14, 2004: PINT = 245.027 mb - data a32/100.00000, 400.00000, 818.60211, & - 1378.88653, 2091.79519, 2983.64084, & - 4121.78960, 5579.22148, 7419.79300, & - 9704.82578, 12496.33710, 15855.26306, & - 19839.62499, 24502.73262, 28177.10152, & - 29525.28447, 29016.34358, 27131.32792, & - 24406.11225, 21326.04907, 18221.18357, & - 15275.14642, 12581.67796, 10181.42843, & - 8081.89816, 6270.86956, 4725.35001, & - 3417.39199, 2317.75459, 1398.09473, & - 632.49506, 0.00000, 0.00000 / - - data b32/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.01711, & - 0.06479, 0.13730, 0.22693, & - 0.32416, 0.42058, 0.51105, & - 0.59325, 0.66628, 0.73011, & - 0.78516, 0.83217, 0.87197, & - 0.90546, 0.93349, 0.95685, & - 0.97624, 0.99223, 1.00000 / -#else -! SJL June 26, 2012 -! pint= 55.7922 - data a32/100.00000, 400.00000, 818.60211, & - 1378.88653, 2091.79519, 2983.64084, & - 4121.78960, 5579.22148, 6907.19063, & - 7735.78639, 8197.66476, 8377.95525, & - 8331.69594, 8094.72213, 7690.85756, & - 7139.01788, 6464.80251, 5712.35727, & - 4940.05347, 4198.60465, 3516.63294, & - 2905.19863, 2366.73733, 1899.19455, & - 1497.78137, 1156.25252, 867.79199, & - 625.59324, 423.21322, 254.76613, & - 115.06646, 0.00000, 0.00000 / - - data b32/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00513, & - 0.01969, 0.04299, 0.07477, & - 0.11508, 0.16408, 0.22198, & - 0.28865, 0.36281, 0.44112, & - 0.51882, 0.59185, 0.65810, & - 0.71694, 0.76843, 0.81293, & - 0.85100, 0.88331, 0.91055, & - 0.93338, 0.95244, 0.96828, & - 0.98142, 0.99223, 1.00000 / -#endif - -!--------------------- -! Wilson's 32L settings: -!--------------------- -! Top changed to 0.01 mb - data a32w/ 1.00, 26.6378, 84.5529, 228.8592, & - 539.9597, 1131.7087, 2141.8082, 3712.0454, & - 5963.5317, 8974.1873, 12764.5388, 17294.5911, & - 20857.7007, 22221.8651, 22892.7202, 22891.1641, & - 22286.0724, 21176.0846, 19673.0671, 17889.0989, & - 15927.5060, 13877.6239, 11812.5474, 9865.8830, & - 8073.9717, 6458.0824, 5027.9893, 3784.6104, & - 2722.0093, 1828.9741, 1090.2397, 487.4575, & - 0.0000 / - - data b32w/ 0.0000, 0.0000, 0.0000, 0.0000, & - 0.0000, 0.0000, 0.0000, 0.0000, & - 0.0000, 0.0000, 0.0000, 0.0000, & - 0.0159, 0.0586, 0.1117, 0.1734, & - 0.2415, 0.3137, 0.3878, 0.4619, & - 0.5344, 0.6039, 0.6696, 0.7285, & - 0.7808, 0.8266, 0.8662, 0.9000, & - 0.9285, 0.9522, 0.9716, 0.9874, & - 1.0000 / - - -#ifdef OLD_L47 -! QBO setting with ptop = 0.1 mb and p_full=0.17 mb; pint ~ 100 mb - data a47/ 10.00000, 24.45365, 48.76776, & - 85.39458, 133.41983, 191.01402, & - 257.94919, 336.63306, 431.52741, & - 548.18995, 692.78825, 872.16512, & - 1094.18467, 1368.11917, 1704.99489, & - 2117.91945, 2622.42986, 3236.88281, & - 3982.89623, 4885.84733, 5975.43260, & - 7286.29500, 8858.72424, 10739.43477, & - 12982.41110, 15649.68745, 18811.37629, & - 22542.71275, 25724.93857, 27314.36781, & - 27498.59474, 26501.79312, 24605.92991, & - 22130.51655, 19381.30274, 16601.56419, & - 13952.53231, 11522.93244, 9350.82303, & - 7443.47723, 5790.77434, 4373.32696, & - 3167.47008, 2148.51663, 1293.15510, & - 581.62429, 0.00000, 0.00000 / - - data b47/ 0.0000, 0.0000, 0.0000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.01188, 0.04650, & - 0.10170, 0.17401, 0.25832, & - 0.34850, 0.43872, 0.52448, & - 0.60307, 0.67328, 0.73492, & - 0.78834, 0.83418, 0.87320, & - 0.90622, 0.93399, 0.95723, & - 0.97650, 0.99223, 1.00000 / -#else -! Oct 23, 2012 -! QBO setting with ptop = 0.1 mb, pint ~ 60 mb - data a47/ 10.00000, 24.45365, 48.76776, & - 85.39458, 133.41983, 191.01402, & - 257.94919, 336.63306, 431.52741, & - 548.18995, 692.78825, 872.16512, & - 1094.18467, 1368.11917, 1704.99489, & - 2117.91945, 2622.42986, 3236.88281, & - 3982.89623, 4885.84733, 5975.43260, & - 7019.26669, 7796.15848, 8346.60209, & - 8700.31838, 8878.27554, 8894.27179, & - 8756.46404, 8469.60171, 8038.92687, & - 7475.89006, 6803.68067, 6058.68992, & - 5285.28859, 4526.01565, 3813.00206, & - 3164.95553, 2589.26318, 2085.96929, & - 1651.11596, 1278.81205, 962.38875, & - 695.07046, 470.40784, 282.61654, & - 126.92745, 0.00000, 0.00000 / - data b47/ 0.0000, 0.0000, 0.0000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00267, 0.01063, 0.02393, & - 0.04282, 0.06771, 0.09917, & - 0.13786, 0.18444, 0.23925, & - 0.30193, 0.37100, 0.44379, & - 0.51695, 0.58727, 0.65236, & - 0.71094, 0.76262, 0.80757, & - 0.84626, 0.87930, 0.90731, & - 0.93094, 0.95077, 0.96733, & - 0.98105, 0.99223, 1.00000 / -#endif - - data a48/ & - 1.00000, 2.69722, 5.17136, & - 8.89455, 14.24790, 22.07157, & - 33.61283, 50.48096, 74.79993, & - 109.40055, 158.00460, 225.44108, & - 317.89560, 443.19350, 611.11558, & - 833.74392, 1125.83405, 1505.20759, & - 1993.15829, 2614.86254, 3399.78420, & - 4382.06240, 5600.87014, 7100.73115, & - 8931.78242, 11149.97021, 13817.16841, & - 17001.20930, 20775.81856, 23967.33875, & - 25527.64563, 25671.22552, 24609.29622, & - 22640.51220, 20147.13482, 17477.63530, & - 14859.86462, 12414.92533, 10201.44191, & - 8241.50255, 6534.43202, 5066.17865, & - 3815.60705, 2758.60264, 1870.64631, & - 1128.33931, 510.47983, 0.00000, & - 0.00000 / - - data b48/ & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.01253, & - 0.04887, 0.10724, 0.18455, & - 0.27461, 0.36914, 0.46103, & - 0.54623, 0.62305, 0.69099, & - 0.75016, 0.80110, 0.84453, & - 0.88127, 0.91217, 0.93803, & - 0.95958, 0.97747, 0.99223, & - 1.00000 / - -! High PBL resolution with top at 1 mb -! SJL modified May 7, 2013 to ptop ~ 100 mb - data a54/100.00000, 254.83931, 729.54278, & - 1602.41121, 2797.50667, 4100.18977, & - 5334.87140, 6455.24153, 7511.80944, & - 8580.26355, 9714.44293, 10938.62253, & - 12080.36051, 12987.13921, 13692.75084, & - 14224.92180, 14606.55444, 14856.69953, & - 14991.32121, 15023.90075, 14965.91493, & - 14827.21612, 14616.33505, 14340.72252, & - 14006.94280, 13620.82849, 13187.60470, & - 12711.98873, 12198.27003, 11650.37451, & - 11071.91608, 10466.23819, 9836.44706, & - 9185.43852, 8515.96231, 7831.01080, & - 7135.14301, 6436.71659, 5749.00215, & - 5087.67188, 4465.67510, 3889.86419, & - 3361.63433, 2879.51065, 2441.02496, & - 2043.41345, 1683.80513, 1359.31122, & - 1067.09135, 804.40101, 568.62625, & - 357.32525, 168.33263, 0.00000, & - 0.00000 / - - data b54/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00180, 0.00694, 0.01510, & - 0.02601, 0.03942, 0.05515, & - 0.07302, 0.09288, 0.11459, & - 0.13803, 0.16307, 0.18960, & - 0.21753, 0.24675, 0.27716, & - 0.30866, 0.34115, 0.37456, & - 0.40879, 0.44375, 0.47935, & - 0.51551, 0.55215, 0.58916, & - 0.62636, 0.66334, 0.69946, & - 0.73395, 0.76622, 0.79594, & - 0.82309, 0.84780, 0.87020, & - 0.89047, 0.90876, 0.92524, & - 0.94006, 0.95336, 0.96529, & - 0.97596, 0.98551, 0.99400, & - 1.00000 / - - -! The 56-L setup - data a56/ 10.00000, 24.97818, 58.01160, & - 115.21466, 199.29210, 309.39897, & - 445.31785, 610.54747, 812.28518, & - 1059.80882, 1363.07092, 1732.09335, & - 2176.91502, 2707.68972, 3334.70962, & - 4068.31964, 4918.76594, 5896.01890, & - 7009.59166, 8268.36324, 9680.41211, & - 11252.86491, 12991.76409, 14901.95764, & - 16987.01313, 19249.15733, 21689.24182, & - 23845.11055, 25330.63353, 26243.52467, & - 26663.84998, 26657.94696, 26281.61371, & - 25583.05256, 24606.03265, 23393.39510, & - 21990.28845, 20445.82122, 18811.93894, & - 17139.59660, 15473.90350, 13850.50167, & - 12294.49060, 10821.62655, 9440.57746, & - 8155.11214, 6965.72496, 5870.70511, & - 4866.83822, 3949.90019, 3115.03562, & - 2357.07879, 1670.87329, 1051.65120, & - 495.51399, 0.00000, 0.00000 / - - data b56 /0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00462, 0.01769, 0.03821, & - 0.06534, 0.09834, 0.13659, & - 0.17947, 0.22637, 0.27660, & - 0.32929, 0.38343, 0.43791, & - 0.49162, 0.54361, 0.59319, & - 0.63989, 0.68348, 0.72391, & - 0.76121, 0.79545, 0.82679, & - 0.85537, 0.88135, 0.90493, & - 0.92626, 0.94552, 0.96286, & - 0.97840, 0.99223, 1.00000 / - -! NAM levels - data a60/200., 1311.4934, 2424.6044, 3541.7594,& - 4662.9584, 5790.2234, 6932.6534, 8095.3034,& - 9278.1734, 10501.4834, 11755.1234, 13049.2034,& - 14403.9434, 15809.2334, 17315.6234, 18953.4434,& - 20783.3534, 22815.4634, 25059.8834, 27567.1634,& - 30148.42896047, 32193.91776039, 33237.35176644, 33332.15200668,& - 32747.34688095, 31710.06232008, 30381.0344269, 28858.71577772,& - 27218.00439794, 25500.31691133, 23734.52294749, 21947.3406187,& - 20167.06984021, 18396.08144096, 16688.20978135, 15067.73749198,& - 13564.49530178, 12183.34512952, 10928.24869364, 9815.02787644,& - 8821.38325756, 7943.05793658, 7181.90985128, 6500.94645341,& - 5932.84856135, 5420.87683616, 4959.15585353, 4522.15047657,& - 4103.63596619, 3703.72540955, 3322.52525084, 2953.65688391,& - 2597.18532669, 2253.10764634, 1915.10585833, 1583.14516612,& - 1257.18953818, 937.3977544 , 623.60136981, 311.11085215,& - 0. / - - data b60/0., 0., 0., 0., 0.,& - 0. , 0. , 0. , 0. , 0. ,& - 0. , 0. , 0. , 0. , 0. ,& - 0. , 0. , 0. , 0. , 0. ,& - 0.0014653 , 0.01021565, 0.0301554 , 0.06025816, 0.09756877,& - 0.13994493, 0.18550048, 0.23318371, 0.2819159 , 0.33120838,& - 0.38067633, 0.42985641, 0.47816985, 0.52569303, 0.57109611,& - 0.61383996, 0.6532309 , 0.68922093, 0.72177094, 0.75052515,& - 0.77610288, 0.79864598, 0.81813309, 0.83553022, 0.85001773,& - 0.86305395, 0.8747947 , 0.88589325, 0.89650986, 0.9066434 ,& - 0.91629284, 0.92562094, 0.93462705, 0.94331221, 0.95183659,& - 0.96020153, 0.96840839, 0.97645359, 0.98434181, 0.99219119, 1. / - -! This is activated by USE_GFSL63 -! Thfollowing L63 setting is the same as NCEP GFS's L64 except the top -! 3 layers - data a63/64.247, 137.790, 221.958, & - 318.266, 428.434, 554.424, & - 698.457, 863.05803, 1051.07995, & - 1265.75194, 1510.71101, 1790.05098, & - 2108.36604, 2470.78817, 2883.03811, & - 3351.46002, 3883.05187, 4485.49315, & - 5167.14603, 5937.04991, 6804.87379, & - 7780.84698, 8875.64338, 10100.20534, & - 11264.35673, 12190.64366, 12905.42546, & - 13430.87867, 13785.88765, 13986.77987, & - 14047.96335, 13982.46770, 13802.40331, & - 13519.33841, 13144.59486, 12689.45608, & - 12165.28766, 11583.57006, 10955.84778, & - 10293.60402, 9608.08306, 8910.07678, & - 8209.70131, 7516.18560, 6837.69250, & - 6181.19473, 5552.39653, 4955.72632, & - 4394.37629, 3870.38682, 3384.76586, & - 2937.63489, 2528.37666, 2155.78385, & - 1818.20722, 1513.68173, 1240.03585, & - 994.99144, 776.23591, 581.48797, & - 408.53400, 255.26520, 119.70243, 0. / - - data b63/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00201, 0.00792, 0.01755, & - 0.03079, 0.04751, 0.06761, & - 0.09097, 0.11746, 0.14690, & - 0.17911, 0.21382, 0.25076, & - 0.28960, 0.32994, 0.37140, & - 0.41353, 0.45589, 0.49806, & - 0.53961, 0.58015, 0.61935, & - 0.65692, 0.69261, 0.72625, & - 0.75773, 0.78698, 0.81398, & - 0.83876, 0.86138, 0.88192, & - 0.90050, 0.91722, 0.93223, & - 0.94565, 0.95762, 0.96827, & - 0.97771, 0.98608, 0.99347, 1./ -#ifdef GFSL64 - data a64/20.00000, 68.00000, 137.79000, & - 221.95800, 318.26600, 428.43400, & - 554.42400, 698.45700, 863.05803, & - 1051.07995, 1265.75194, 1510.71101, & - 1790.05098, 2108.36604, 2470.78817, & - 2883.03811, 3351.46002, 3883.05187, & - 4485.49315, 5167.14603, 5937.04991, & - 6804.87379, 7780.84698, 8875.64338, & - 9921.40745, 10760.99844, 11417.88354, & - 11911.61193, 12258.61668, 12472.89642, & - 12566.58298, 12550.43517, 12434.26075, & - 12227.27484, 11938.39468, 11576.46910, & - 11150.43640, 10669.41063, 10142.69482, & - 9579.72458, 8989.94947, 8382.67090, & - 7766.85063, 7150.91171, 6542.55077, & - 5948.57894, 5374.81094, 4825.99383, & - 4305.79754, 3816.84622, 3360.78848, & - 2938.39801, 2549.69756, 2194.08449, & - 1870.45732, 1577.34218, 1313.00028, & - 1075.52114, 862.90778, 673.13815, & - 504.22118, 354.22752, 221.32110, & - 103.78014, 0./ - data b64/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00179, 0.00705, 0.01564, & - 0.02749, 0.04251, 0.06064, & - 0.08182, 0.10595, 0.13294, & - 0.16266, 0.19492, 0.22950, & - 0.26615, 0.30455, 0.34435, & - 0.38516, 0.42656, 0.46815, & - 0.50949, 0.55020, 0.58989, & - 0.62825, 0.66498, 0.69987, & - 0.73275, 0.76351, 0.79208, & - 0.81845, 0.84264, 0.86472, & - 0.88478, 0.90290, 0.91923, & - 0.93388, 0.94697, 0.95865, & - 0.96904, 0.97826, 0.98642, & - 0.99363, 1./ -#else - data a64/1.00000, 3.90000, 8.70000, & - 15.42000, 24.00000, 34.50000, & - 47.00000, 61.50000, 78.60000, & - 99.13500, 124.12789, 154.63770, & - 191.69700, 236.49300, 290.38000, & - 354.91000, 431.82303, 523.09300, & - 630.92800, 757.79000, 906.45000, & - 1079.85000, 1281.00000, 1515.00000, & - 1788.00000, 2105.00000, 2470.00000, & - 2889.00000, 3362.00000, 3890.00000, & - 4475.00000, 5120.00000, 5830.00000, & - 6608.00000, 7461.00000, 8395.00000, & - 9424.46289, 10574.46880, 11864.80270, & - 13312.58890, 14937.03710, 16759.70700, & - 18804.78710, 21099.41210, 23674.03710, & - 26562.82810, 29804.11720, 32627.31640, & - 34245.89840, 34722.28910, 34155.19920, & - 32636.50390, 30241.08200, 27101.44920, & - 23362.20700, 19317.05270, 15446.17090, & - 12197.45210, 9496.39941, 7205.66992, & - 5144.64307, 3240.79346, 1518.62134, & - 0.00000, 0.00000 / - - data b64/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00813, & - 0.03224, 0.07128, 0.12445, & - 0.19063, 0.26929, 0.35799, & - 0.45438, 0.55263, 0.64304, & - 0.71703, 0.77754, 0.82827, & - 0.87352, 0.91502, 0.95235, & - 0.98511, 1.00000 / -#endif -!-->cjg - data a68/1.00000, 2.68881, 5.15524, & - 8.86683, 14.20349, 22.00278, & - 33.50807, 50.32362, 74.56680, & - 109.05958, 157.51214, 224.73844, & - 316.90481, 441.81219, 609.21090, & - 831.14537, 1122.32514, 1500.51628, & - 1986.94617, 2606.71274, 3389.18802, & - 4368.40473, 5583.41379, 7078.60015, & - 8903.94455, 11115.21886, 13774.60566, & - 16936.82070, 20340.47045, 23193.71492, & - 24870.36141, 25444.59363, 25252.57081, & - 24544.26211, 23474.29096, 22230.65331, & - 20918.50731, 19589.96280, 18296.26682, & - 17038.02866, 15866.85655, 14763.18943, & - 13736.83624, 12794.11850, 11930.72442, & - 11137.17217, 10404.78946, 9720.03954, & - 9075.54055, 8466.72650, 7887.12346, & - 7333.90490, 6805.43028, 6297.33773, & - 5805.78227, 5327.94995, 4859.88765, & - 4398.63854, 3942.81761, 3491.08449, & - 3043.04531, 2598.71608, 2157.94527, & - 1720.87444, 1287.52805, 858.02944, & - 432.71276, 8.10905, 0.00000 / - - data b68/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00283, 0.01590, & - 0.04412, 0.08487, 0.13284, & - 0.18470, 0.23828, 0.29120, & - 0.34211, 0.39029, 0.43518, & - 0.47677, 0.51536, 0.55091, & - 0.58331, 0.61263, 0.63917, & - 0.66333, 0.68552, 0.70617, & - 0.72555, 0.74383, 0.76117, & - 0.77765, 0.79335, 0.80838, & - 0.82287, 0.83693, 0.85069, & - 0.86423, 0.87760, 0.89082, & - 0.90392, 0.91689, 0.92973, & - 0.94244, 0.95502, 0.96747, & - 0.97979, 0.99200, 1.00000 / - - data a96/1.00000, 2.35408, 4.51347, & - 7.76300, 12.43530, 19.26365, & - 29.33665, 44.05883, 65.28397, & - 95.48274, 137.90344, 196.76073, & - 277.45330, 386.81095, 533.37018, & - 727.67600, 982.60677, 1313.71685, & - 1739.59104, 2282.20281, 2967.26766, & - 3824.58158, 4888.33404, 6197.38450, & - 7795.49158, 9731.48414, 11969.71024, & - 14502.88894, 17304.52434, 20134.76139, & - 22536.63814, 24252.54459, 25230.65591, & - 25585.72044, 25539.91412, 25178.87141, & - 24644.84493, 23978.98781, 23245.49366, & - 22492.11600, 21709.93990, 20949.64473, & - 20225.94258, 19513.31158, 18829.32485, & - 18192.62250, 17589.39396, 17003.45386, & - 16439.01774, 15903.91204, 15396.39758, & - 14908.02140, 14430.65897, 13967.88643, & - 13524.16667, 13098.30227, 12687.56457, & - 12287.08757, 11894.41553, 11511.54106, & - 11139.22483, 10776.01912, 10419.75711, & - 10067.11881, 9716.63489, 9369.61967, & - 9026.69066, 8687.29884, 8350.04978, & - 8013.20925, 7677.12187, 7343.12994, & - 7011.62844, 6681.98102, 6353.09764, & - 6025.10535, 5699.10089, 5375.54503, & - 5053.63074, 4732.62740, 4413.38037, & - 4096.62775, 3781.79777, 3468.45371, & - 3157.19882, 2848.25306, 2541.19150, & - 2236.21942, 1933.50628, 1632.83741, & - 1334.35954, 1038.16655, 744.22318, & - 452.71094, 194.91899, 0.00000, & - 0.00000 / - - data b96/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00193, & - 0.00974, 0.02538, 0.04876, & - 0.07817, 0.11081, 0.14514, & - 0.18007, 0.21486, 0.24866, & - 0.28088, 0.31158, 0.34030, & - 0.36701, 0.39210, 0.41554, & - 0.43733, 0.45774, 0.47707, & - 0.49540, 0.51275, 0.52922, & - 0.54495, 0.56007, 0.57459, & - 0.58850, 0.60186, 0.61471, & - 0.62715, 0.63922, 0.65095, & - 0.66235, 0.67348, 0.68438, & - 0.69510, 0.70570, 0.71616, & - 0.72651, 0.73675, 0.74691, & - 0.75700, 0.76704, 0.77701, & - 0.78690, 0.79672, 0.80649, & - 0.81620, 0.82585, 0.83542, & - 0.84492, 0.85437, 0.86375, & - 0.87305, 0.88229, 0.89146, & - 0.90056, 0.90958, 0.91854, & - 0.92742, 0.93623, 0.94497, & - 0.95364, 0.96223, 0.97074, & - 0.97918, 0.98723, 0.99460, & - 1.00000 / -!<--cjg -! -! Ultra high troposphere resolution - data a100/100.00000, 300.00000, 800.00000, & - 1762.35235, 3106.43596, 4225.71874, & - 4946.40525, 5388.77387, 5708.35540, & - 5993.33124, 6277.73673, 6571.49996, & - 6877.05339, 7195.14327, 7526.24920, & - 7870.82981, 8229.35361, 8602.30193, & - 8990.16936, 9393.46399, 9812.70768, & - 10248.43625, 10701.19980, 11171.56286, & - 11660.10476, 12167.41975, 12694.11735, & - 13240.82253, 13808.17600, 14396.83442, & - 15007.47066, 15640.77407, 16297.45067, & - 16978.22343, 17683.83253, 18415.03554, & - 19172.60771, 19957.34218, 20770.05022, & - 21559.14829, 22274.03147, 22916.87519, & - 23489.70456, 23994.40187, 24432.71365, & - 24806.25734, 25116.52754, 25364.90190, & - 25552.64670, 25680.92203, 25750.78675, & - 25763.20311, 25719.04113, 25619.08274, & - 25464.02630, 25254.49482, 24991.06137, & - 24674.32737, 24305.11235, 23884.79781, & - 23415.77059, 22901.76510, 22347.84738, & - 21759.93950, 21144.07284, 20505.73136, & - 19849.54271, 19179.31412, 18498.23400, & - 17809.06809, 17114.28232, 16416.10343, & - 15716.54833, 15017.44246, 14320.43478, & - 13627.01116, 12938.50682, 12256.11762, & - 11580.91062, 10913.83385, 10255.72526, & - 9607.32122, 8969.26427, 8342.11044, & - 7726.33606, 7122.34405, 6530.46991, & - 5950.98721, 5384.11279, 4830.01153, & - 4288.80090, 3760.55514, 3245.30920, & - 2743.06250, 2253.78294, 1777.41285, & - 1313.88054, 863.12371, 425.13088, & - 0.00000, 0.00000 / - - - data b100/0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00000, 0.00000, 0.00000, & - 0.00052, 0.00209, 0.00468, & - 0.00828, 0.01288, 0.01849, & - 0.02508, 0.03266, 0.04121, & - 0.05075, 0.06126, 0.07275, & - 0.08521, 0.09866, 0.11308, & - 0.12850, 0.14490, 0.16230, & - 0.18070, 0.20009, 0.22042, & - 0.24164, 0.26362, 0.28622, & - 0.30926, 0.33258, 0.35605, & - 0.37958, 0.40308, 0.42651, & - 0.44981, 0.47296, 0.49591, & - 0.51862, 0.54109, 0.56327, & - 0.58514, 0.60668, 0.62789, & - 0.64872, 0.66919, 0.68927, & - 0.70895, 0.72822, 0.74709, & - 0.76554, 0.78357, 0.80117, & - 0.81835, 0.83511, 0.85145, & - 0.86736, 0.88286, 0.89794, & - 0.91261, 0.92687, 0.94073, & - 0.95419, 0.96726, 0.97994, & - 0.99223, 1.00000 / - - data a104/ & - 1.8827062944e-01, 7.7977549145e-01, 2.1950593583e+00, & - 4.9874566624e+00, 9.8041418997e+00, 1.7019717163e+01, & - 2.7216579591e+01, 4.0518628401e+01, 5.6749646818e+01, & - 7.5513868331e+01, 9.6315093333e+01, 1.1866706195e+02, & - 1.4216835396e+02, 1.6653733709e+02, 1.9161605772e+02, & - 2.1735580129e+02, 2.4379516604e+02, 2.7103771847e+02, & - 2.9923284173e+02, 3.2856100952e+02, 3.5922338766e+02, & - 3.9143507908e+02, 4.2542117983e+02, 4.6141487902e+02, & - 4.9965698106e+02, 5.4039638379e+02, 5.8389118154e+02, & - 6.3041016829e+02, 6.8023459505e+02, 7.3366009144e+02, & - 7.9099869949e+02, 8.5258099392e+02, 9.1875827946e+02, & - 9.8990486716e+02, 1.0664204381e+03, 1.1487325074e+03, & - 1.2372990044e+03, 1.3326109855e+03, 1.4351954993e+03, & - 1.5456186222e+03, 1.6644886848e+03, 1.7924597105e+03, & - 1.9302350870e+03, 2.0785714934e+03, 2.2382831070e+03, & - 2.4102461133e+03, 2.5954035462e+03, 2.7947704856e+03, & - 3.0094396408e+03, 3.2405873512e+03, 3.4894800360e+03, & - 3.7574811281e+03, 4.0460585279e+03, 4.3567926151e+03, & - 4.6913848588e+03, 5.0516670674e+03, 5.4396113207e+03, & - 5.8573406270e+03, 6.3071403487e+03, 6.7914704368e+03, & - 7.3129785102e+03, 7.8745138115e+03, 8.4791420557e+03, & - 9.1301611750e+03, 9.8311179338e+03, 1.0585825354e+04, & - 1.1398380836e+04, 1.2273184781e+04, 1.3214959424e+04, & - 1.4228767429e+04, 1.5320029596e+04, 1.6494540743e+04, & - 1.7758482452e+04, 1.9118430825e+04, 2.0422798801e+04, & - 2.1520147587e+04, 2.2416813461e+04, 2.3118184510e+04, & - 2.3628790785e+04, 2.3952411814e+04, 2.4092209011e+04, & - 2.4050892106e+04, 2.3830930156e+04, 2.3434818358e+04, & - 2.2865410898e+04, 2.2126326004e+04, 2.1222420323e+04, & - 2.0160313690e+04, 1.8948920926e+04, 1.7599915822e+04, & - 1.6128019809e+04, 1.4550987232e+04, 1.2889169132e+04, & - 1.1164595563e+04, 9.4227665517e+03, 7.7259097899e+03, & - 6.1538244381e+03, 4.7808126007e+03, 3.5967415552e+03, & - 2.5886394104e+03, 1.7415964865e+03, 1.0393721271e+03, & - 4.6478852032e+02, 7.0308342481e-13, 0.0000000000e+00 / - - - data b104/ & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 0.0000000000e+00, & - 0.0000000000e+00, 0.0000000000e+00, 1.5648447298e-03, & - 6.2617046389e-03, 1.4104157933e-02, 2.5118187415e-02, & - 3.9340510972e-02, 5.6816335609e-02, 7.7596328431e-02, & - 1.0173255472e-01, 1.2927309709e-01, 1.6025505622e-01, & - 1.9469566981e-01, 2.3258141217e-01, 2.7385520518e-01, & - 3.1840233814e-01, 3.6603639170e-01, 4.1648734767e-01, & - 4.6939496013e-01, 5.2431098738e-01, 5.8071350676e-01, & - 6.3803478105e-01, 6.9495048840e-01, 7.4963750338e-01, & - 7.9975208897e-01, 8.4315257576e-01, 8.8034012292e-01, & - 9.1184389721e-01, 9.3821231526e-01, 9.6000677644e-01, & - 9.7779792223e-01, 9.9216315122e-01, 1.0000000000e+00 / - -! IFS-like L125(top 12 levels removed from IFSL137) - data a125/ 64., & - 86.895882, 107.415741, 131.425507, 159.279404, 191.338562, & - 227.968948, 269.539581, 316.420746, 368.982361, 427.592499, 492.616028, & - 564.413452, 643.339905, 729.744141, 823.967834, 926.344910, 1037.201172, & - 1156.853638, 1285.610352, 1423.770142, 1571.622925, 1729.448975, 1897.519287, & - 2076.095947, 2265.431641, 2465.770508, 2677.348145, 2900.391357, 3135.119385, & - 3381.743652, 3640.468262, 3911.490479, 4194.930664, 4490.817383, 4799.149414, & - 5119.895020, 5452.990723, 5798.344727, 6156.074219, 6526.946777, 6911.870605, & - 7311.869141, 7727.412109, 8159.354004, 8608.525391, 9076.400391, 9562.682617, & - 10065.978516, 10584.631836, 11116.662109, 11660.067383, 12211.547852, 12766.873047, & - 13324.668945, 13881.331055, 14432.139648, 14975.615234, 15508.256836, 16026.115234, & - 16527.322266, 17008.789063, 17467.613281, 17901.621094, 18308.433594, 18685.718750, & - 19031.289063, 19343.511719, 19620.042969, 19859.390625, 20059.931641, 20219.664063, & - 20337.863281, 20412.308594, 20442.078125, 20425.718750, 20361.816406, 20249.511719, & - 20087.085938, 19874.025391, 19608.572266, 19290.226563, 18917.460938, 18489.707031, & - 18006.925781, 17471.839844, 16888.687500, 16262.046875, 15596.695313, 14898.453125, & - 14173.324219, 13427.769531, 12668.257813, 11901.339844, 11133.304688, 10370.175781, & - 9617.515625, 8880.453125, 8163.375000, 7470.343750, 6804.421875, 6168.531250, & - 5564.382813, 4993.796875, 4457.375000, 3955.960938, 3489.234375, 3057.265625, & - 2659.140625, 2294.242188, 1961.500000, 1659.476563, 1387.546875, 1143.250000, & - 926.507813, 734.992188, 568.062500, 424.414063, 302.476563, 202.484375, & - 122.101563, 62.781250, 22.835938, 3.757813, 0.000000, 0.000000 / - - data b125/ 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, & - 0.000000, 0.000007, 0.000024, 0.000059, 0.000112, 0.000199, & - 0.000340, 0.000562, 0.000890, 0.001353, 0.001992, 0.002857, & - 0.003971, 0.005378, 0.007133, 0.009261, 0.011806, 0.014816, & - 0.018318, 0.022355, 0.026964, 0.032176, 0.038026, 0.044548, & - 0.051773, 0.059728, 0.068448, 0.077958, 0.088286, 0.099462, & - 0.111505, 0.124448, 0.138313, 0.153125, 0.168910, 0.185689, & - 0.203491, 0.222333, 0.242244, 0.263242, 0.285354, 0.308598, & - 0.332939, 0.358254, 0.384363, 0.411125, 0.438391, 0.466003, & - 0.493800, 0.521619, 0.549301, 0.576692, 0.603648, 0.630036, & - 0.655736, 0.680643, 0.704669, 0.727739, 0.749797, 0.770798, & - 0.790717, 0.809536, 0.827256, 0.843881, 0.859432, 0.873929, & - 0.887408, 0.899900, 0.911448, 0.922096, 0.931881, 0.940860, & - 0.949064, 0.956550, 0.963352, 0.969513, 0.975078, 0.980072, & - 0.984542, 0.988500, 0.991984, 0.995003, 0.997630, 1.000000 / - - select case (km) - - case (24) - - ks = 5 - do k=1,km+1 - ak(k) = a24(k) - bk(k) = b24(k) - enddo - - case (26) - - ks = 7 - do k=1,km+1 - ak(k) = a26(k) - bk(k) = b26(k) - enddo - - case (32) -#ifdef OLD_32 - ks = 13 ! high-res trop_32 setup -#else - ks = 7 -#endif - do k=1,km+1 - ak(k) = a32(k) - bk(k) = b32(k) - enddo - - case (47) -! ks = 27 ! high-res trop-strat - ks = 20 ! Oct 23, 2012 - do k=1,km+1 - ak(k) = a47(k) - bk(k) = b47(k) - enddo - - case (48) - ks = 28 - do k=1,km+1 - ak(k) = a48(k) - bk(k) = b48(k) - enddo - - case (52) - ks = 35 ! pint = 223 - do k=1,km+1 - ak(k) = a52(k) - bk(k) = b52(k) - enddo - - case (54) - ks = 11 ! pint = 109.4 - do k=1,km+1 - ak(k) = a54(k) - bk(k) = b54(k) - enddo - - case (56) - ks = 26 - do k=1,km+1 - ak(k) = a56(k) - bk(k) = b56(k) - enddo - - case (60) - ks = 19 - do k=1,km+1 - ak(k) = a60(k) - bk(k) = b60(k) - enddo - - - case (64) -#ifdef GFSL64 - ks = 23 -#else - ks = 46 -#endif - do k=1,km+1 - ak(k) = a64(k) - bk(k) = b64(k) - enddo -!-->cjg - case (68) - ks = 27 - do k=1,km+1 - ak(k) = a68(k) - bk(k) = b68(k) - enddo - - case (96) - ks = 27 - do k=1,km+1 - ak(k) = a96(k) - bk(k) = b96(k) - enddo -!<--cjg - - case (100) - ks = 38 - do k=1,km+1 - ak(k) = a100(k) - bk(k) = b100(k) - enddo - - case (104) - ks = 73 - do k=1,km+1 - ak(k) = a104(k) - bk(k) = b104(k) - enddo - -#ifndef TEST_GWAVES - case (10) -!-------------------------------------------------- -! Pure sigma-coordinate with uniform spacing in "z" -!-------------------------------------------------- -! - pt = 2000. ! model top pressure (pascal) -! pt = 100. ! 1 mb - press(1) = pt - press(km+1) = p0 - dlnp = (log(p0) - log(pt)) / real(km) - - lnpe = log(press(km+1)) - do k=km,2,-1 - lnpe = lnpe - dlnp - press(k) = exp(lnpe) - enddo - -! Search KS - ks = 0 - do k=1,km - if(press(k) >= pc) then - ks = k-1 - goto 123 - endif - enddo -123 continue - - if(ks /= 0) then - do k=1,ks - ak(k) = press(k) - bk(k) = 0. - enddo - endif - - pint = press(ks+1) - do k=ks+1,km - ak(k) = pint*(press(km)-press(k))/(press(km)-pint) - bk(k) = (press(k) - ak(k)) / press(km+1) - enddo - ak(km+1) = 0. - bk(km+1) = 1. - -! do k=2,km -! bk(k) = real(k-1) / real(km) -! ak(k) = pt * ( 1. - bk(k) ) -! enddo -#endif - -! The following 4 selections are better for non-hydrostatic -! Low top: - case (31) - ptop = 300. - pint = 100.E2 - call var_dz(km, ak, bk, ptop, ks, pint, 1.035) -#ifndef TEST_GWAVES - case (41) - ptop = 100. - pint = 100.E2 - call var_hi(km, ak, bk, ptop, ks, pint, 1.035) -#endif - case (51) - ptop = 100. - pint = 100.E2 - call var_hi(km, ak, bk, ptop, ks, pint, 1.035) -! Mid-top: - case (55) - ptop = 10. - pint = 100.E2 -! call var_dz(km, ak, bk, ptop, ks, pint, 1.035) - call var_hi(km, ak, bk, ptop, ks, pint, 1.035) -#ifdef USE_GFSL63 -! GFS L64 equivalent setting - case (63) - ks = 23 - ptop = a63(1) - pint = a63(ks+1) - do k=1,km+1 - ak(k) = a63(k) - bk(k) = b63(k) - enddo -#else - case (63) - ptop = 1. ! high top - pint = 100.E2 - call var_hi(km, ak, bk, ptop, ks, pint, 1.035) -#endif -! NGGPS_GFS - case (91) - pint = 100.E2 - ptop = 40. - call var_gfs(km, ak, bk, ptop, ks, pint, 1.029) -! call var_gfs(km, ak, bk, ptop, ks, pint, 1.03) - case (95) -! Mid-top settings: - pint = 100.E2 - ptop = 20. - call var_gfs(km, ak, bk, ptop, ks, pint, 1.028) - case (127) - ptop = 1. - pint = 75.E2 - call var_gfs(km, ak, bk, ptop, ks, pint, 1.028) -! IFS-like L125 - case (125) - ks = 33 - ptop = a125(1) - pint = a125(ks+1) - do k=1,km+1 - ak(k) = a125(k) - bk(k) = b125(k) - enddo - case default - -#ifdef TEST_GWAVES -!-------------------------------------------------- -! Pure sigma-coordinate with uniform spacing in "z" -!-------------------------------------------------- - call gw_1d(km, 1000.E2, ak, bk, ptop, 10.E3, pt1) - ks = 0 - pint = ak(1) -#else - -!---------------------------------------------------------------- -! Sigma-coordinate with uniform spacing in sigma and ptop = 1 mb -!---------------------------------------------------------------- - pt = 100. -! One pressure layer - ks = 1 -! pint = pt + 0.5*1.E5/real(km) ! SJL: 20120327 - pint = pt + 1.E5/real(km) - - ak(1) = pt - bk(1) = 0. - ak(2) = pint - bk(2) = 0. - - do k=3,km+1 - bk(k) = real(k-2) / real(km-1) - ak(k) = pint - bk(k)*pint - enddo - ak(km+1) = 0. - bk(km+1) = 1. -#endif - end select - ptop = ak(1) - pint = ak(ks+1) - - end subroutine set_eta -#endif - -!>@brief The subroutine 'set_external_eta' sets 'ptop' (model top) and -!! 'ks' (first level of pure pressure coordinates given the coefficients -!! 'ak' and 'bk' - subroutine set_external_eta(ak, bk, ptop, ks) - implicit none - real, intent(in) :: ak(:) - real, intent(in) :: bk(:) - real, intent(out) :: ptop !< model top (Pa) - integer, intent(out) :: ks !< number of pure p layers - !--- local variables - integer :: k - real :: eps = 1.d-7 - - ptop = ak(1) - ks = 1 - do k = 1, size(bk(:)) - if (bk(k).lt.eps) ks = k - enddo - !--- change ks to layers from levels - ks = ks - 1 - if (is_master()) write(6,*) ' ptop & ks ', ptop, ks - - end subroutine set_external_eta - - - subroutine var_les(km, ak, bk, ptop, ks, pint, s_rate) - implicit none - integer, intent(in):: km - real, intent(in):: ptop - real, intent(in):: s_rate !< between [1. 1.1] - real, intent(out):: ak(km+1), bk(km+1) - real, intent(inout):: pint - integer, intent(out):: ks -! Local - real, parameter:: p00 = 1.E5 - real, dimension(km+1):: ze, pe1, peln, eta - real, dimension(km):: dz, s_fac, dlnp, pm, dp, dk - real ztop, t0, dz0, sum1, tmp1 - real ep, es, alpha, beta, gama - real, parameter:: akap = 2./7. -!---- Tunable parameters: - real:: k_inc = 10 !< number of layers from bottom up to near const dz region - real:: s0 = 0.8 !< lowest layer stretch factor -!----------------------- - real:: s_inc - integer k - - pe1(1) = ptop - peln(1) = log(pe1(1)) - pe1(km+1) = p00 - peln(km+1) = log(pe1(km+1)) - - t0 = 273. - ztop = rdgas/grav*t0*(peln(km+1) - peln(1)) - - s_inc = (1.-s0) / real(k_inc) - s_fac(km) = s0 - - do k=km-1, km-k_inc, -1 - s_fac(k) = s_fac(k+1) + s_inc - enddo - - s_fac(km-k_inc-1) = 0.5*(s_fac(km-k_inc) + s_rate) - - do k=km-k_inc-2, 5, -1 - s_fac(k) = s_rate * s_fac(k+1) - enddo - - s_fac(4) = 0.5*(1.1+s_rate)*s_fac(5) - s_fac(3) = 1.1 *s_fac(4) - s_fac(2) = 1.1 *s_fac(3) - s_fac(1) = 1.1 *s_fac(2) - - sum1 = 0. - do k=1,km - sum1 = sum1 + s_fac(k) - enddo - - dz0 = ztop / sum1 - - do k=1,km - dz(k) = s_fac(k) * dz0 - enddo - - ze(km+1) = 0. - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo - -! Re-scale dz with the stretched ztop - do k=1,km - dz(k) = dz(k) * (ztop/ze(1)) - enddo - - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo -! ze(1) = ztop - - if ( is_master() ) then - write(*,*) 'var_les: computed model top (m)=', ztop, ' bottom/top dz=', dz(km), dz(1) -! do k=1,km -! write(*,*) k, s_fac(k) -! enddo - endif - - call sm1_edge(1, 1, 1, 1, km, 1, 1, ze, 2) - -! Given z --> p - do k=1,km - dz(k) = ze(k) - ze(k+1) - dlnp(k) = grav*dz(k) / (rdgas*t0) - !write(*,*) k, dz(k) - enddo - do k=2,km - peln(k) = peln(k-1) + dlnp(k-1) - pe1(k) = exp(peln(k)) - enddo - - -! Pe(k) = ak(k) + bk(k) * PS -! Locate pint and KS - ks = 0 - do k=2,km - if ( pint < pe1(k)) then - ks = k-1 - exit - endif - enddo - if ( is_master() ) then - write(*,*) 'For (input) PINT=', 0.01*pint, ' KS=', ks, 'pint(computed)=', 0.01*pe1(ks+1) - endif - pint = pe1(ks+1) - - do k=1,km+1 - eta(k) = pe1(k) / pe1(km+1) - enddo - - ep = eta(ks+1) - es = eta(km) -! es = 1. - alpha = (ep**2-2.*ep*es) / (es-ep)**2 - beta = 2.*ep*es**2 / (es-ep)**2 - gama = -(ep*es)**2 / (es-ep)**2 - -! Pure pressure: - do k=1,ks+1 - ak(k) = eta(k)*1.e5 - bk(k) = 0. - enddo - - do k=ks+2, km - ak(k) = alpha*eta(k) + beta + gama/eta(k) - ak(k) = ak(k)*1.e5 - enddo - ak(km+1) = 0. - - do k=ks+2, km - bk(k) = (pe1(k) - ak(k))/pe1(km+1) - enddo - bk(km+1) = 1. - - if ( is_master() ) then - ! write(*,*) 'KS=', ks, 'PINT (mb)=', pint/100. - ! do k=1,km - ! pm(k) = 0.5*(pe1(k)+pe1(k+1))/100. - ! write(*,*) k, pm(k), dz(k) - ! enddo - tmp1 = ak(ks+1) - do k=ks+1,km - tmp1 = max(tmp1, (ak(k)-ak(k+1))/max(1.E-5, (bk(k+1)-bk(k))) ) - enddo - write(*,*) 'Hybrid Sigma-P: minimum allowable surface pressure (hpa)=', tmp1/100. - write(*,800) (pm(k), k=km,1,-1) - endif - - do k=1,km - dp(k) = (pe1(k+1) - pe1(k))/100. - dk(k) = pe1(k+1)**akap - pe1(k)**akap - enddo - -800 format(1x,5(1x,f9.4)) - - end subroutine var_les - - - - subroutine var_gfs(km, ak, bk, ptop, ks, pint, s_rate) - integer, intent(in):: km - real, intent(in):: ptop - real, intent(in):: s_rate !< between [1. 1.1] - real, intent(out):: ak(km+1), bk(km+1) - real, intent(inout):: pint - integer, intent(out):: ks -! Local - real, parameter:: p00 = 1.E5 - real, dimension(km+1):: ze, pe1, peln, eta - real, dimension(km):: dz, s_fac, dlnp - real ztop, t0, dz0, sum1, tmp1 - real ep, es, alpha, beta, gama -!---- Tunable parameters: - integer:: k_inc = 25 !< number of layers from bottom up to near const dz region - real:: s0 = 0.13 !< lowest layer stretch factor -!----------------------- - real:: s_inc - integer k - - pe1(1) = ptop - peln(1) = log(pe1(1)) - pe1(km+1) = p00 - peln(km+1) = log(pe1(km+1)) - - t0 = 270. - ztop = rdgas/grav*t0*(peln(km+1) - peln(1)) - - s_inc = (1.-s0) / real(k_inc) - s_fac(km) = s0 - - do k=km-1, km-k_inc, -1 - s_fac(k) = s_fac(k+1) + s_inc - enddo - - do k=km-k_inc-1, 9, -1 - s_fac(k) = s_rate * s_fac(k+1) - enddo - s_fac(8) = 0.5*(1.1+s_rate)*s_fac(9) - s_fac(7) = 1.10*s_fac(8) - s_fac(6) = 1.15*s_fac(7) - s_fac(5) = 1.20*s_fac(6) - s_fac(4) = 1.26*s_fac(5) - s_fac(3) = 1.33*s_fac(4) - s_fac(2) = 1.41*s_fac(3) - s_fac(1) = 1.60*s_fac(2) - - sum1 = 0. - do k=1,km - sum1 = sum1 + s_fac(k) - enddo - - dz0 = ztop / sum1 - - do k=1,km - dz(k) = s_fac(k) * dz0 - enddo - - ze(km+1) = 0. - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo - -! Re-scale dz with the stretched ztop - do k=1,km - dz(k) = dz(k) * (ztop/ze(1)) - enddo - - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo -! ze(1) = ztop - - if ( is_master() ) then - write(*,*) 'var_gfs: computed model top (m)=', ztop*0.001, ' bottom/top dz=', dz(km), dz(1) -! do k=1,km -! write(*,*) k, s_fac(k) -! enddo - endif - -! call sm1_edge(1, 1, 1, 1, km, 1, 1, ze, 3) - -! Given z --> p - do k=1,km - dz(k) = ze(k) - ze(k+1) - dlnp(k) = grav*dz(k) / (rdgas*t0) - enddo - do k=2,km - peln(k) = peln(k-1) + dlnp(k-1) - pe1(k) = exp(peln(k)) - enddo - -! Pe(k) = ak(k) + bk(k) * PS -! Locate pint and KS - ks = 0 - do k=2,km - if ( pint < pe1(k)) then - ks = k-1 - exit - endif - enddo - if ( is_master() ) then - write(*,*) 'For (input) PINT=', 0.01*pint, ' KS=', ks, 'pint(computed)=', 0.01*pe1(ks+1) - write(*,*) 'ptop =', ptop - endif - pint = pe1(ks+1) - -#ifdef NO_UKMO_HB - do k=1,ks+1 - ak(k) = pe1(k) - bk(k) = 0. - enddo - - do k=ks+2,km+1 - bk(k) = (pe1(k) - pint) / (pe1(km+1)-pint) ! bk == sigma - ak(k) = pe1(k) - bk(k) * pe1(km+1) - enddo - bk(km+1) = 1. - ak(km+1) = 0. -#else -! Problematic for non-hydrostatic - do k=1,km+1 - eta(k) = pe1(k) / pe1(km+1) - enddo - - ep = eta(ks+1) - es = eta(km) -! es = 1. - alpha = (ep**2-2.*ep*es) / (es-ep)**2 - beta = 2.*ep*es**2 / (es-ep)**2 - gama = -(ep*es)**2 / (es-ep)**2 - -! Pure pressure: - do k=1,ks+1 - ak(k) = eta(k)*1.e5 - bk(k) = 0. - enddo - - do k=ks+2, km - ak(k) = alpha*eta(k) + beta + gama/eta(k) - ak(k) = ak(k)*1.e5 - enddo - ak(km+1) = 0. - - do k=ks+2, km - bk(k) = (pe1(k) - ak(k))/pe1(km+1) - enddo - bk(km+1) = 1. -#endif - - if ( is_master() ) then - write(*,*) 'KS=', ks, 'PINT (mb)=', pint/100. - do k=1,km - write(*,*) k, 0.5*(pe1(k)+pe1(k+1))/100., dz(k) - enddo - tmp1 = ak(ks+1) - do k=ks+1,km - tmp1 = max(tmp1, (ak(k)-ak(k+1))/max(1.E-5, (bk(k+1)-bk(k))) ) - enddo - write(*,*) 'Hybrid Sigma-P: minimum allowable surface pressure (hpa)=', tmp1/100. - endif - - end subroutine var_gfs - - subroutine var_hi(km, ak, bk, ptop, ks, pint, s_rate) - integer, intent(in):: km - real, intent(in):: ptop - real, intent(in):: s_rate !< between [1. 1.1] - real, intent(out):: ak(km+1), bk(km+1) - real, intent(inout):: pint - integer, intent(out):: ks -! Local - real, parameter:: p00 = 1.E5 - real, dimension(km+1):: ze, pe1, peln, eta - real, dimension(km):: dz, s_fac, dlnp - real ztop, t0, dz0, sum1, tmp1 - real ep, es, alpha, beta, gama -!---- Tunable parameters: - integer:: k_inc = 15 ! 1.E-8 ) then - pf(1) = (ph(2) - ph(1)) / log(ph(2)/ph(1)) - else - pf(1) = (ph(2) - ph(1)) * kappa/(kappa+1.) - endif - - do k=2,npz - pf(k) = (ph(k+1) - ph(k)) / log(ph(k+1)/ph(k)) - enddo - - end subroutine get_eta_level - - - - subroutine compute_dz(km, ztop, dz) - - integer, intent(in):: km - real, intent(in):: ztop ! try 50.E3 - real, intent(out):: dz(km) -!------------------------------ - real ze(km+1), dzt(km) - integer k - - -! ztop = 30.E3 - dz(1) = ztop / real(km) - dz(km) = 0.5*dz(1) - - do k=2,km-1 - dz(k) = dz(1) - enddo - -! Top: - dz(1) = 2.*dz(2) - - ze(km+1) = 0. - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo - - if ( is_master() ) then - write(*,*) 'Hybrid_z: dz, zm' - do k=1,km - dzt(k) = 0.5*(ze(k)+ze(k+1)) / 1000. - write(*,*) k, dz(k), dzt(k) - enddo - endif - - end subroutine compute_dz - - subroutine compute_dz_var(km, ztop, dz) - - integer, intent(in):: km - real, intent(in):: ztop ! try 50.E3 - real, intent(out):: dz(km) -!------------------------------ - real, parameter:: s_rate = 1.0 - real ze(km+1) - real s_fac(km) - real sum1, dz0 - integer k - - s_fac(km ) = 0.125 - s_fac(km-1) = 0.20 - s_fac(km-2) = 0.30 - s_fac(km-3) = 0.40 - s_fac(km-4) = 0.50 - s_fac(km-5) = 0.60 - s_fac(km-6) = 0.70 - s_fac(km-7) = 0.80 - s_fac(km-8) = 0.90 - s_fac(km-9) = 1. - - do k=km-10, 9, -1 - s_fac(k) = s_rate * s_fac(k+1) - enddo - - s_fac(8) = 1.05*s_fac(9) - s_fac(7) = 1.1 *s_fac(8) - s_fac(6) = 1.15*s_fac(7) - s_fac(5) = 1.2 *s_fac(6) - s_fac(4) = 1.3 *s_fac(5) - s_fac(3) = 1.4 *s_fac(4) - s_fac(2) = 1.5 *s_fac(3) - s_fac(1) = 1.6 *s_fac(2) - - sum1 = 0. - do k=1,km - sum1 = sum1 + s_fac(k) - enddo - - dz0 = ztop / sum1 - - do k=1,km - dz(k) = s_fac(k) * dz0 - enddo - - ze(1) = ztop - ze(km+1) = 0. - do k=km,2,-1 - ze(k) = ze(k+1) + dz(k) - enddo - -! Re-scale dz with the stretched ztop - do k=1,km - dz(k) = dz(k) * (ztop/ze(1)) - enddo - - do k=km,2,-1 - ze(k) = ze(k+1) + dz(k) - enddo - - call sm1_edge(1, 1, 1, 1, km, 1, 1, ze, 2) - - do k=1,km - dz(k) = ze(k) - ze(k+1) - enddo - - end subroutine compute_dz_var - - subroutine compute_dz_L32(km, ztop, dz) - - integer, intent(in):: km - real, intent(out):: dz(km) - real, intent(out):: ztop ! try 50.E3 -!------------------------------ - real dzt(km) - real ze(km+1) - real dz0, dz1, dz2 - real z0, z1, z2 - integer k, k0, k1, k2, n - -!------------------- - k2 = 8 - z2 = 30.E3 -!------------------- - k1 = 21 - z1 = 10.0E3 -!------------------- - k0 = 2 - z0 = 0. - dz0 = 75. ! meters -!------------------- -! Treat the surface layer as a special layer - ze(1) = z0 - dz(1) = dz0 - - ze(2) = dz(1) - dz0 = 1.5*dz0 - dz(2) = dz0 - - ze(3) = ze(2) + dz(2) - - dz1 = 2.*(z1-ze(3) - k1*dz0) / (k1*(k1-1)) - - do k=k0+1,k0+k1 - dz(k) = dz0 + (k-k0)*dz1 - ze(k+1) = ze(k) + dz(k) - enddo - - dz0 = dz(k1+k0) - dz2 = 2.*(z2-ze(k0+k1+1)-k2*dz0) / (k2*(k2-1)) - - do k=k0+k1+1,k0+k1+k2 - dz(k) = dz0 + (k-k0-k1)*dz2 - ze(k+1) = ze(k) + dz(k) - enddo - - dz(km) = 2.*dz(km-1) - ztop = ze(km) + dz(km) - ze(km+1) = ze(km) + dz(km) - - call zflip (dz, 1, km) - - ze(km+1) = 0. - do k=km,1,-1 - ze(k) = ze(k+1) + dz(k) - enddo - -! if ( is_master() ) then -! write(*,*) 'Hybrid_z: dz, zm' -! do k=1,km -! dzt(k) = 0.5*(ze(k)+ze(k+1)) / 1000. -! write(*,*) k, dz(k), dzt(k) -! enddo -! endif - - end subroutine compute_dz_L32 - - subroutine compute_dz_L101(km, ztop, dz) - - integer, intent(in):: km ! km==101 - real, intent(out):: dz(km) - real, intent(out):: ztop ! try 30.E3 -!------------------------------ - real ze(km+1) - real dz0, dz1 - real:: stretch_f = 1.16 - integer k, k0, k1 - - k1 = 2 - k0 = 25 - dz0 = 40. ! meters - - ze(km+1) = 0. - - do k=km, k0, -1 - dz(k) = dz0 - ze(k) = ze(k+1) + dz(k) - enddo - - do k=k0+1, k1, -1 - dz(k) = stretch_f * dz(k+1) - ze(k) = ze(k+1) + dz(k) - enddo - - dz(1) = 4.0*dz(2) - ze(1) = ze(2) + dz(1) - ztop = ze(1) - - if ( is_master() ) then - write(*,*) 'Hybrid_z: dz, ze' - do k=1,km - write(*,*) k, 0.001*dz(k), 0.001*ze(k) - enddo -! ztop (km) = 20.2859154 - write(*,*) 'ztop (km) =', ztop * 0.001 - endif - - end subroutine compute_dz_L101 - - subroutine set_hybrid_z(is, ie, js, je, ng, km, ztop, dz, rgrav, hs, ze, dz3) - - integer, intent(in):: is, ie, js, je, ng, km - real, intent(in):: rgrav, ztop - real, intent(in):: dz(km) !< Reference vertical resolution for zs=0 - real, intent(in):: hs(is-ng:ie+ng,js-ng:je+ng) - real, intent(inout):: ze(is:ie,js:je,km+1) - real, optional, intent(out):: dz3(is-ng:ie+ng,js-ng:je+ng,km) -! local - logical:: filter_xy = .false. - real, allocatable:: delz(:,:,:) - integer ntimes - real zint - real:: z1(is:ie,js:je) - real:: z(km+1) - real sig, z_rat - integer ks(is:ie,js:je) - integer i, j, k, ks_min, kint - - z(km+1) = 0. - do k=km,1,-1 - z(k) = z(k+1) + dz(k) - enddo - - do j=js,je - do i=is,ie - ze(i,j, 1) = ztop - ze(i,j,km+1) = hs(i,j) * rgrav - enddo - enddo - - do k=2,km - do j=js,je - do i=is,ie - ze(i,j,k) = z(k) - enddo - enddo - enddo - -! Set interface: -#ifndef USE_VAR_ZINT - zint = 12.0E3 - ntimes = 2 - kint = 2 - do k=2,km - if ( z(k)<=zint ) then - kint = k - exit - endif - enddo - - if ( is_master() ) write(*,*) 'Z_coord interface set at k=',kint, ' ZE=', z(kint) - - do j=js,je - do i=is,ie - z_rat = (ze(i,j,kint)-ze(i,j,km+1)) / (z(kint)-z(km+1)) - do k=km,kint+1,-1 - ze(i,j,k) = ze(i,j,k+1) + dz(k)*z_rat - enddo -!-------------------------------------- -! Apply vertical smoother locally to dz -!-------------------------------------- - call sm1_edge(is, ie, js, je, km, i, j, ze, ntimes) - enddo - enddo -#else -! ZINT is a function of local terrain - ntimes = 4 - do j=js,je - do i=is,ie - z1(i,j) = dim(ze(i,j,km+1), 2500.) + 7500. - enddo - enddo - - ks_min = km - do j=js,je - do i=is,ie - do k=km,2,-1 - if ( z(k)>=z1(i,j) ) then - ks(i,j) = k - ks_min = min(ks_min, k) - go to 555 - endif - enddo -555 continue - enddo - enddo - - do j=js,je - do i=is,ie - kint = ks(i,j) + 1 - z_rat = (ze(i,j,kint)-ze(i,j,km+1)) / (z(kint)-z(km+1)) - do k=km,kint+1,-1 - ze(i,j,k) = ze(i,j,k+1) + dz(k)*z_rat - enddo -!-------------------------------------- -! Apply vertical smoother locally to dz -!-------------------------------------- - call sm1_edge(is, ie, js, je, km, i, j, ze, ntimes) - enddo - enddo -#endif - -#ifdef DEV_ETA - if ( filter_xy ) then - allocate (delz(isd:ied, jsd:jed, km) ) - ntimes = 2 - do k=1,km - do j=js,je - do i=is,ie - delz(i,j,k) = ze(i,j,k+1) - ze(i,j,k) - enddo - enddo - enddo - call del2_cubed(delz, 0.2*da_min, npx, npy, km, ntimes) - do k=km,2,-1 - do j=js,je - do i=is,ie - ze(i,j,k) = ze(i,j,k+1) - delz(i,j,k) - enddo - enddo - enddo - deallocate ( delz ) - endif -#endif - if ( present(dz3) ) then - do k=1,km - do j=js,je - do i=is,ie - dz3(i,j,k) = ze(i,j,k+1) - ze(i,j,k) - enddo - enddo - enddo - endif - - end subroutine set_hybrid_z - - - subroutine sm1_edge(is, ie, js, je, km, i, j, ze, ntimes) - integer, intent(in):: is, ie, js, je, km - integer, intent(in):: ntimes, i, j - real, intent(inout):: ze(is:ie,js:je,km+1) -! local: - real, parameter:: df = 0.25 - real dz(km) - real flux(km+1) - integer k, n, k1, k2 - - k2 = km-1 - do k=1,km - dz(k) = ze(i,j,k+1) - ze(i,j,k) - enddo - - do n=1,ntimes - k1 = 2 + (ntimes-n) - - flux(k1 ) = 0. - flux(k2+1) = 0. - do k=k1+1,k2 - flux(k) = df*(dz(k) - dz(k-1)) - enddo - - do k=k1,k2 - dz(k) = dz(k) - flux(k) + flux(k+1) - enddo - enddo - - do k=km,1,-1 - ze(i,j,k) = ze(i,j,k+1) - dz(k) - enddo - - end subroutine sm1_edge - - - - subroutine gw_1d(km, p0, ak, bk, ptop, ztop, pt1) - integer, intent(in):: km - real, intent(in):: p0, ztop - real, intent(inout):: ptop - real, intent(inout):: ak(km+1), bk(km+1) - real, intent(out):: pt1(km) -! Local - logical:: isothermal - real, dimension(km+1):: ze, pe1, pk1 - real, dimension(km):: dz1 - real t0, n2, s0 - integer k - -! Set up vertical coordinare with constant del-z spacing: - isothermal = .false. - t0 = 300. - - if ( isothermal ) then - n2 = grav**2/(cp_air*t0) - else - n2 = 0.0001 - endif - - s0 = grav*grav / (cp_air*n2) - - ze(km+1) = 0. - do k=km,1,-1 - dz1(k) = ztop / real(km) - ze(k) = ze(k+1) + dz1(k) - enddo - -! Given z --> p - do k=1,km+1 - pe1(k) = p0*( (1.-s0/t0) + s0/t0*exp(-n2*ze(k)/grav) )**(1./kappa) - enddo - - ptop = pe1(1) -! if ( is_master() ) write(*,*) 'GW_1D: computed model top (pa)=', ptop - -! Set up "sigma" coordinate - ak(1) = pe1(1) - bk(1) = 0. - do k=2,km - bk(k) = (pe1(k) - pe1(1)) / (pe1(km+1)-pe1(1)) ! bk == sigma - ak(k) = pe1(1)*(1.-bk(k)) - enddo - ak(km+1) = 0. - bk(km+1) = 1. - - do k=1,km+1 - pk1(k) = pe1(k) ** kappa - enddo - -! Compute volume mean potential temperature with hydrostatic eqn: - do k=1,km - pt1(k) = grav*dz1(k) / ( cp_air*(pk1(k+1)-pk1(k)) ) - enddo - - end subroutine gw_1d - - - - subroutine zflip(q, im, km) - integer, intent(in):: im, km - real, intent(inout):: q(im,km) -!--- - integer i, k - real qtmp - - do i = 1, im - do k = 1, (km+1)/2 - qtmp = q(i,k) - q(i,k) = q(i,km+1-k) - q(i,km+1-k) = qtmp - end do - end do - - end subroutine zflip - -end module fv_eta_mod diff --git a/tools/test_cases.F90 b/tools/test_cases.F90 index ae3019a5a..9146d230e 100644 --- a/tools/test_cases.F90 +++ b/tools/test_cases.F90 @@ -113,7 +113,7 @@ module test_cases_mod use mpp_mod, only: mpp_error, FATAL, mpp_root_pe, mpp_broadcast, mpp_sum use mpp_mod, only: stdlog, input_nml_file - use fms_mod, only: check_nml_error + use fms_mod, only: check_nml_error, close_file, open_namelist_file use mpp_domains_mod, only: mpp_update_domains, domain2d use mpp_parameter_mod, only: AGRID_PARAM=>AGRID,CGRID_NE_PARAM=>CGRID_NE, & SCALAR_PAIR