24 #define _GET_VAR1 get_var1_real 26 #define _GET_VAR1 get_var1_double 149 use fms_mod
, only: file_exist, read_data, field_exist, write_version_number
150 use fms_mod
, only: open_namelist_file, check_nml_error, close_file
151 use fms_mod
, only: get_mosaic_tile_file, error_mesg
152 use fms_io_mod
, only: get_tile_string, field_size, free_restart_type
153 use fms_io_mod
, only: restart_file_type, register_restart_field
154 use fms_io_mod
, only: save_restart, restore_state, set_filename_appendix, get_global_att_value
155 use mpp_mod
, only: mpp_error, fatal, note, mpp_pe, mpp_root_pe
156 use mpp_mod
, only: stdlog, input_nml_file
157 use mpp_parameter_mod
, only: agrid_param=>agrid
158 use mpp_domains_mod
, only: mpp_get_tile_id, domain2d, mpp_update_domains, north, east
159 use tracer_manager_mod
, only: get_tracer_names, get_number_tracers, get_tracer_index
160 use tracer_manager_mod
, only: set_tracer_profile
161 use field_manager_mod
, only: model_atmos
163 use constants_mod
, only: pi=>pi_8, omega, grav, kappa, rdgas, rvgas, cp_air
171 use fv_mp_mod, only: is_master, fill_corners, ydir, mp_reduce_min, mp_reduce_max
189 use mpp_domains_mod
, only: mpp_get_data_domain, mpp_get_global_domain, mpp_get_compute_domain
198 real,
parameter::
zvir = rvgas/rdgas - 1.
202 character(len=27),
parameter ::
source_fv3gfs =
'FV3GFS GAUSSIAN NEMSIO FILE' 208 #include<file_version.h> 215 type(domain2d),
intent(inout) :: fv_domain
216 logical,
intent(IN) :: cold_start
217 real,
intent(IN) :: dt_atmos
222 real,
pointer,
dimension(:,:,:) :: grid, agrid
223 real,
pointer,
dimension(:,:) :: fC, f0
225 integer :: is, ie, js, je
226 integer :: isd, ied, jsd, jed, ng
227 integer :: sphum, liq_wat, ice_wat, rainwat, snowwat, graupel
229 integer :: liq_aero, ice_aero
232 integer :: spfo, spfo2, spfo3
247 grid => atm%gridstruct%grid
248 agrid => atm%gridstruct%agrid
250 fc => atm%gridstruct%fC
251 f0 => atm%gridstruct%f0
257 fc(i,j) = 2.*omega*( -1.*cos(grid(i,j,1))*cos(grid(i,j,2))*sin(alpha) + &
258 sin(grid(i,j,2))*cos(alpha) )
264 f0(i,j) = 2.*omega*( -1.*cos(agrid(i,j,1))*cos(agrid(i,j,2))*sin(alpha) + &
265 sin(agrid(i,j,2))*cos(alpha) )
269 call mpp_update_domains( f0, fv_domain )
270 if ( atm%gridstruct%cubed_sphere .and. (.not. atm%gridstruct%bounded_domain))
then 271 call fill_corners(f0, atm%npx, atm%npy, ydir)
275 if ( atm%flagstruct%mountain )
then 278 if (.not. atm%neststruct%nested) atm%phis = 0.
282 if ( atm%flagstruct%ncep_ic )
then 288 if (.not. cold_start)
then 290 if(is_master())
write(*,*)
'All tracers except sphum replaced by FV IC' 293 elseif ( atm%flagstruct%nggps_ic )
then 297 elseif ( atm%flagstruct%ecmwf_ic )
then 298 if( is_master() )
write(*,*)
'Calling get_ecmwf_ic' 309 call prt_maxmin(
'PS', atm%ps, is, ie, js, je, ng, 1, 0.01)
310 call prt_maxmin(
'T', atm%pt, is, ie, js, je, ng, atm%npz, 1.)
311 if (.not.atm%flagstruct%hydrostatic)
call prt_maxmin(
'W', atm%w, is, ie, js, je, ng, atm%npz, 1.)
312 call prt_maxmin(
'SPHUM', atm%q(:,:,:,1), is, ie, js, je, ng, atm%npz, 1.)
313 if ( atm%flagstruct%nggps_ic )
then 314 call prt_maxmin(
'TS', atm%ts, is, ie, js, je, 0, 1, 1.)
316 if ( atm%flagstruct%nggps_ic .or. atm%flagstruct%ecmwf_ic )
then 317 sphum = get_tracer_index(model_atmos,
'sphum')
318 liq_wat = get_tracer_index(model_atmos,
'liq_wat')
319 ice_wat = get_tracer_index(model_atmos,
'ice_wat')
320 rainwat = get_tracer_index(model_atmos,
'rainwat')
321 snowwat = get_tracer_index(model_atmos,
'snowwat')
322 graupel = get_tracer_index(model_atmos,
'graupel')
324 spfo = get_tracer_index(model_atmos,
'spfo')
325 spfo2 = get_tracer_index(model_atmos,
'spfo2')
326 spfo3 = get_tracer_index(model_atmos,
'spfo3')
328 o3mr = get_tracer_index(model_atmos,
'o3mr')
331 liq_aero = get_tracer_index(model_atmos,
'liq_aero')
332 ice_aero = get_tracer_index(model_atmos,
'ice_aero')
336 call prt_maxmin(
'liq_wat', atm%q(:,:,:,liq_wat), is, ie, js, je, ng, atm%npz, 1.)
338 call prt_maxmin(
'ice_wat', atm%q(:,:,:,ice_wat), is, ie, js, je, ng, atm%npz, 1.)
340 call prt_maxmin(
'rainwat', atm%q(:,:,:,rainwat), is, ie, js, je, ng, atm%npz, 1.)
342 call prt_maxmin(
'snowwat', atm%q(:,:,:,snowwat), is, ie, js, je, ng, atm%npz, 1.)
344 call prt_maxmin(
'graupel', atm%q(:,:,:,graupel), is, ie, js, je, ng, atm%npz, 1.)
347 call prt_maxmin(
'SPFO', atm%q(:,:,:,spfo), is, ie, js, je, ng, atm%npz, 1.)
349 call prt_maxmin(
'SPFO2', atm%q(:,:,:,spfo2), is, ie, js, je, ng, atm%npz, 1.)
351 call prt_maxmin(
'SPFO3', atm%q(:,:,:,spfo3), is, ie, js, je, ng, atm%npz, 1.)
354 call prt_maxmin(
'O3MR', atm%q(:,:,:,o3mr), is, ie, js, je, ng, atm%npz, 1.)
358 call prt_maxmin(
'liq_aero',atm%q(:,:,:,liq_aero),is, ie, js, je, ng, atm%npz, 1.)
360 call prt_maxmin(
'ice_aero',atm%q(:,:,:,ice_aero),is, ie, js, je, ng, atm%npz, 1.)
378 type(domain2d),
intent(inout) :: fv_domain
379 integer :: tile_id(1)
380 character(len=64) :: fname
381 character(len=7) :: gn
383 integer :: jbeg, jend
385 real,
allocatable :: g_dat2(:,:,:)
386 real,
allocatable :: pt_coarse(:,:,:)
387 integer isc_p, iec_p, jsc_p, jec_p, isg, ieg, jsg,jeg
389 integer :: is, ie, js, je
390 integer :: isd, ied, jsd, jed, ng
402 if (atm%grid_number > 1)
then 404 write(gn,
'(A5, I2.2)')
".nest", atm%grid_number
409 tile_id = mpp_get_tile_id( fv_domain )
411 call get_tile_string(fname,
'INPUT/fv_core.res'//trim(gn)//
'.tile', tile_id(n),
'.nc' )
412 if (mpp_pe() == mpp_root_pe()) print*,
'external_ic: looking for ', fname
415 if( file_exist(fname) )
then 416 call read_data(fname,
'phis', atm%phis(is:ie,js:je), &
417 domain=fv_domain, tile_count=n)
419 call surfdrv( atm%npx, atm%npy, atm%gridstruct%grid_64, atm%gridstruct%agrid_64, &
420 atm%gridstruct%area_64, atm%gridstruct%dx, atm%gridstruct%dy, &
421 atm%gridstruct%dxa, atm%gridstruct%dya, &
422 atm%gridstruct%dxc, atm%gridstruct%dyc, atm%gridstruct%sin_sg, &
423 atm%phis, atm%flagstruct%stretch_fac, &
424 atm%neststruct%nested, atm%gridstruct%bounded_domain, &
425 atm%neststruct%npx_global, atm%domain, &
426 atm%flagstruct%grid_number, atm%bd )
427 call mpp_error(note,
'terrain datasets generated using USGS data')
432 call mpp_update_domains( atm%phis, atm%domain )
433 ftop =
g_sum(atm%domain, atm%phis(is:ie,js:je), is, ie, js, je, ng, atm%gridstruct%area_64, 1)
435 call prt_maxmin(
'ZS', atm%phis, is, ie, js, je, ng, 1, 1./grav)
436 if(is_master())
write(*,*)
'mean terrain height (m)=', ftop/grav
466 use gfs_restart
, only : gfs_restart_type
472 type(domain2d),
intent(inout) :: fv_domain
473 real,
intent(in) :: dt_atmos
475 real,
dimension(:),
allocatable:: ak, bk
476 real,
dimension(:,:),
allocatable:: wk2, ps, oro_g
477 real,
dimension(:,:,:),
allocatable:: ud, vd, u_w, v_w, u_s, v_s, omga, temp
478 real,
dimension(:,:,:),
allocatable:: zh(:,:,:)
479 real,
dimension(:,:,:,:),
allocatable:: q
480 real,
dimension(:,:),
allocatable :: phis_coarse
481 real rdg, wt, qt, m_fac, pe1
482 integer:: n, npx, npy, npz, itoa, nt, ntprog, ntdiag, ntracers, ntrac, iq
483 integer :: is, ie, js, je
484 integer :: isd, ied, jsd, jed
485 integer :: ios, ierr, unit, id_res
486 type(restart_file_type) :: oro_restart, sfc_restart, gfs_restart
487 character(len=6) :: gn, stile_name
488 character(len=64) :: tracer_name
489 character(len=64) :: fn_gfs_ctl =
'gfs_ctrl.nc' 490 character(len=64) :: fn_gfs_ics =
'gfs_data.nc' 491 character(len=64) :: fn_sfc_ics =
'sfc_data.nc' 492 character(len=64) :: fn_oro_ics =
'oro_data.nc' 495 logical :: filtered_terrain = .true.
496 logical :: gfs_dwinds = .true.
498 logical :: checker_tr = .false.
499 integer :: nt_checker = 0
500 real(kind=R_GRID),
dimension(2):: p1, p2, p3
501 real(kind=R_GRID),
dimension(3):: e1, e2, ex, ey
502 integer:: i,j,k,nts, ks
503 integer:: liq_wat, ice_wat, rainwat, snowwat, graupel, tke, ntclamt
504 namelist /external_ic_nml/ filtered_terrain, levp, gfs_dwinds, &
505 checker_tr, nt_checker
509 call mpp_error(note,
'Using external_IC::get_nggps_ic which is valid only for data which has been & 510 &horizontally interpolated to the current cubed-sphere grid')
511 #ifdef INTERNAL_FILE_NML 512 read (input_nml_file,external_ic_nml,iostat=ios)
513 ierr = check_nml_error(ios,
'external_ic_nml')
515 unit=open_namelist_file()
516 read (unit,external_ic_nml,iostat=ios)
517 ierr = check_nml_error(ios,
'external_ic_nml')
518 call close_file(unit)
522 call write_version_number (
'EXTERNAL_IC_MOD::get_nggps_ic', version )
523 write(unit, nml=external_ic_nml)
526 if (atm%flagstruct%external_eta)
then 527 if (filtered_terrain)
then 528 call mpp_error(note,
'External_IC::get_nggps_ic - use externally-generated, filtered terrain & 529 &and NCEP pressure levels (no vertical remapping)')
530 else if (.not. filtered_terrain)
then 531 call mpp_error(note,
'External_IC::get_nggps_ic - use externally-generated, raw terrain & 532 &and NCEP pressure levels (no vertical remapping)')
535 if (filtered_terrain)
then 536 call mpp_error(note,
'External_IC::get_nggps_ic - use externally-generated, filtered terrain & 537 &and FV3 pressure levels (vertical remapping)')
538 else if (.not. filtered_terrain)
then 539 call mpp_error(note,
'External_IC::get_nggps_ic - use externally-generated, raw terrain & 540 &and FV3 pressure levels (vertical remapping)')
553 write(*,22001)is,ie,js,je,isd,ied,jsd,jed
554 22001
format(
' enter get_nggps_ic is=',i4,
' ie=',i4,
' js=',i4,
' je=',i4,
' isd=',i4,
' ied=',i4,
' jsd=',i4,
' jed=',i4)
555 call get_number_tracers(model_atmos, num_tracers=ntracers, num_prog=ntprog)
556 ntdiag = ntracers-ntprog
559 if (atm%grid_number > 1)
then 560 write(gn,
'(A4, I2.2)')
"nest", atm%grid_number
564 call set_filename_appendix(
'')
567 if (.not. file_exist(
'INPUT/'//trim(fn_gfs_ctl), no_domain=.true.))
then 568 call mpp_error(fatal,
'==> Error in External_ic::get_nggps_ic: file '//trim(fn_gfs_ctl)//
' for NGGPS IC does not exist')
570 call mpp_error(note,
'==> External_ic::get_nggps_ic: using control file '//trim(fn_gfs_ctl)//
' for NGGPS IC')
573 call read_data (
'INPUT/'//trim(fn_gfs_ctl),
'ntrac', ntrac, no_domain=.true.)
580 call get_data_source(
source,atm%flagstruct%regional)
582 call mpp_error(note,
"READING FROM REGRIDDED FV3GFS NEMSIO FILE")
586 allocate (wk2(levp+1,2))
587 allocate (ak(levp+1))
588 allocate (bk(levp+1))
589 call read_data(
'INPUT/'//trim(fn_gfs_ctl),
'vcoord',wk2, no_domain=.true.)
590 ak(1:levp+1) = wk2(1:levp+1,1)
591 bk(1:levp+1) = wk2(1:levp+1,2)
594 if (.not. file_exist(
'INPUT/'//trim(fn_oro_ics), domain=atm%domain))
then 595 call mpp_error(fatal,
'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_oro_ics)//
' for NGGPS IC does not exist')
597 call mpp_error(note,
'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_oro_ics)//
' for NGGPS IC')
599 if (.not. file_exist(
'INPUT/'//trim(fn_sfc_ics), domain=atm%domain))
then 600 call mpp_error(fatal,
'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_sfc_ics)//
' for NGGPS IC does not exist')
602 call mpp_error(note,
'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_sfc_ics)//
' for NGGPS IC')
604 if (.not. file_exist(
'INPUT/'//trim(fn_gfs_ics), domain=atm%domain))
then 605 call mpp_error(fatal,
'==> Error in External_ic::get_nggps_ic: tiled file '//trim(fn_gfs_ics)//
' for NGGPS IC does not exist')
607 call mpp_error(note,
'==> External_ic::get_nggps_ic: using tiled data file '//trim(fn_gfs_ics)//
' for NGGPS IC')
609 allocate (zh(is:ie,js:je,levp+1))
610 allocate (ps(is:ie,js:je))
611 allocate (omga(is:ie,js:je,levp))
612 allocate (q(is:ie,js:je,levp,ntracers))
613 allocate ( u_w(is:ie+1, js:je, 1:levp) )
614 allocate ( v_w(is:ie+1, js:je, 1:levp) )
615 allocate ( u_s(is:ie, js:je+1, 1:levp) )
616 allocate ( v_s(is:ie, js:je+1, 1:levp) )
617 allocate (temp(is:ie,js:je,levp))
620 if (atm%neststruct%nested)
then 621 allocate(phis_coarse(isd:ied,jsd:jed))
624 phis_coarse(i,j) = atm%phis(i,j)
631 id_res = register_restart_field(sfc_restart, fn_sfc_ics,
'tsea', atm%ts, domain=atm%domain)
634 if (filtered_terrain)
then 635 id_res = register_restart_field(oro_restart, fn_oro_ics,
'orog_filt', atm%phis, domain=atm%domain)
636 elseif (.not. filtered_terrain)
then 637 id_res = register_restart_field(oro_restart, fn_oro_ics,
'orog_raw', atm%phis, domain=atm%domain)
640 if ( atm%flagstruct%full_zs_filter)
then 641 allocate (oro_g(isd:ied,jsd:jed))
644 id_res = register_restart_field(oro_restart, fn_oro_ics,
'land_frac', oro_g, domain=atm%domain)
645 call mpp_update_domains(oro_g, atm%domain)
646 if (atm%neststruct%nested)
then 651 if ( atm%flagstruct%fv_land )
then 653 id_res = register_restart_field(oro_restart, fn_oro_ics,
'stddev', atm%sgh, domain=atm%domain)
655 id_res = register_restart_field(oro_restart, fn_oro_ics,
'land_frac', atm%oro, domain=atm%domain)
659 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'ps', ps, domain=atm%domain)
662 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'u_w', u_w, domain=atm%domain,position=east)
664 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'v_w', v_w, domain=atm%domain,position=east)
666 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'u_s', u_s, domain=atm%domain,position=north)
668 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'v_s', v_s, domain=atm%domain,position=north)
671 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'w', omga, domain=atm%domain)
673 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'ZH', zh, domain=atm%domain)
676 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
't', temp, mandatory=.false., &
680 call get_tracer_names(model_atmos, nt, tracer_name)
682 id_res = register_restart_field(gfs_restart, fn_gfs_ics, trim(tracer_name), q(:,:,:,nt), &
683 mandatory=.false.,domain=atm%domain)
688 call get_tracer_names(model_atmos, nt, tracer_name)
690 call set_tracer_profile (model_atmos, nt, atm%q(:,:,:,nt) )
692 do nt = ntprog+1, ntracers
693 call get_tracer_names(model_atmos, nt, tracer_name)
695 call set_tracer_profile (model_atmos, nt, atm%qdiag(:,:,:,nt) )
699 call restore_state (oro_restart)
700 call restore_state (sfc_restart)
701 call restore_state (gfs_restart)
704 call free_restart_type(oro_restart)
705 call free_restart_type(sfc_restart)
706 call free_restart_type(gfs_restart)
709 atm%phis = atm%phis*grav
713 if (atm%flagstruct%external_eta)
then 714 itoa = levp - npz + 1
716 atm%ak(1:npz+1) = ak(itoa:levp+1)
717 atm%bk(1:npz+1) = bk(itoa:levp+1)
721 if(is_master())
write(*,*)
'GFS ak(1)=', ak(1),
' ak(2)=', ak(2)
722 ak(1) = max(1.e-9, ak(1))
730 if (n==1.and.atm%flagstruct%regional)
then 740 call remap_scalar(atm, levp, npz, ntracers, ak, bk, ps, q, zh, omga, temp)
742 allocate ( ud(is:ie, js:je+1, 1:levp) )
743 allocate ( vd(is:ie+1,js:je, 1:levp) )
750 p1(:) = atm%gridstruct%grid(i, j,1:2)
751 p2(:) = atm%gridstruct%grid(i+1,j,1:2)
760 p1(:) = atm%gridstruct%grid(i,j ,1:2)
761 p2(:) = atm%gridstruct%grid(i,j+1,1:2)
779 if (atm%neststruct%nested)
then 780 if (is_master())
write(*,*)
'Blending nested and coarse grid topography' 785 wt = max(0.,min(1.,
real(5 - min(i,j,npx-i,npy-j,5))/5. ))
786 atm%phis(i,j) = (1.-wt)*atm%phis(i,j) + wt*phis_coarse(i,j)
793 if ( atm%flagstruct%full_zs_filter )
then 795 call mpp_update_domains(atm%phis, atm%domain)
797 call fv3_zs_filter( atm%bd, isd, ied, jsd, jed, npx, npy, atm%neststruct%npx_global, &
798 atm%flagstruct%stretch_fac, atm%gridstruct%bounded_domain, atm%domain, &
799 atm%gridstruct%area_64, atm%gridstruct%dxa, atm%gridstruct%dya, &
800 atm%gridstruct%dx, atm%gridstruct%dy, atm%gridstruct%dxc, &
801 atm%gridstruct%dyc, atm%gridstruct%grid_64, atm%gridstruct%agrid_64, &
802 atm%gridstruct%sin_sg, atm%phis, oro_g)
807 if ( atm%flagstruct%n_zs_filter > 0 )
then 809 if ( atm%flagstruct%nord_zs_filter == 2 )
then 811 atm%gridstruct%area_64, atm%gridstruct%dx, atm%gridstruct%dy, &
812 atm%gridstruct%dxc, atm%gridstruct%dyc, atm%gridstruct%sin_sg, &
813 atm%flagstruct%n_zs_filter,
cnst_0p20*atm%gridstruct%da_min, &
814 .false., oro_g, atm%gridstruct%bounded_domain, &
816 if ( is_master() )
write(*,*)
'Warning !!! del-2 terrain filter has been applied ', &
817 atm%flagstruct%n_zs_filter,
' times' 818 else if( atm%flagstruct%nord_zs_filter == 4 )
then 820 atm%gridstruct%dx, atm%gridstruct%dy, &
821 atm%gridstruct%dxc, atm%gridstruct%dyc, atm%gridstruct%sin_sg, &
822 atm%flagstruct%n_zs_filter, .false., oro_g, &
823 atm%gridstruct%bounded_domain, &
825 if ( is_master() )
write(*,*)
'Warning !!! del-4 terrain filter has been applied ', &
826 atm%flagstruct%n_zs_filter,
' times' 831 if ( atm%neststruct%nested .and. ( atm%flagstruct%n_zs_filter > 0 .or. atm%flagstruct%full_zs_filter ) )
then 836 wt = max(0.,min(1.,
real(5 - min(i,j,npx-i,npy-j,5))/5. ))
837 atm%phis(i,j) = (1.-wt)*atm%phis(i,j) + wt*phis_coarse(i,j)
840 deallocate(phis_coarse)
843 call mpp_update_domains( atm%phis, atm%domain, complete=.true. )
844 liq_wat = get_tracer_index(model_atmos,
'liq_wat')
845 ice_wat = get_tracer_index(model_atmos,
'ice_wat')
846 rainwat = get_tracer_index(model_atmos,
'rainwat')
847 snowwat = get_tracer_index(model_atmos,
'snowwat')
848 graupel = get_tracer_index(model_atmos,
'graupel')
849 ntclamt = get_tracer_index(model_atmos,
'cld_amt')
855 if ( atm%flagstruct%nwat == 6 )
then 856 qt = wt*(1. + atm%q(i,j,k,liq_wat) + &
857 atm%q(i,j,k,ice_wat) + &
858 atm%q(i,j,k,rainwat) + &
859 atm%q(i,j,k,snowwat) + &
860 atm%q(i,j,k,graupel))
862 qt = wt*(1. + sum(atm%q(i,j,k,2:atm%flagstruct%nwat)))
865 if (ntclamt > 0) atm%q(i,j,k,ntclamt) = 0.0
876 if ( atm%flagstruct%nwat == 6 )
then 877 qt = wt*(1. + atm%q(i,j,k,liq_wat) + &
878 atm%q(i,j,k,ice_wat) + &
879 atm%q(i,j,k,rainwat) + &
880 atm%q(i,j,k,snowwat) + &
881 atm%q(i,j,k,graupel))
883 qt = wt*(1. + sum(atm%q(i,j,k,2:atm%flagstruct%nwat)))
887 atm%q(i,j,k,iq) = m_fac * atm%q(i,j,k,iq)
890 if (ntclamt > 0) atm%q(i,j,k,ntclamt) = 0.0
898 tke = get_tracer_index(model_atmos,
'sgs_tke')
904 atm%q(i,j,k,tke) = 0.00
912 nts = ntracers - nt_checker+1
914 npz, atm%q(:,:,:,nts:ntracers), &
915 atm%gridstruct%agrid_64(is:ie,js:je,1), &
916 atm%gridstruct%agrid_64(is:ie,js:je,2), 9., 9.)
919 atm%flagstruct%make_nh = .false.
934 type(domain2d),
intent(inout) :: fv_domain
935 integer,
intent(in):: nq
938 real :: ak_HIWPP(65), bk_HIWPP(65)
940 0, 0.00064247, 0.0013779, 0.00221958, 0.00318266, 0.00428434, &
941 0.00554424, 0.00698457, 0.00863058, 0.0105108, 0.01265752, 0.01510711, &
942 0.01790051, 0.02108366, 0.02470788, 0.02883038, 0.0335146, 0.03883052, &
943 0.04485493, 0.05167146, 0.0593705, 0.06804874, 0.0777715, 0.08832537, &
944 0.09936614, 0.1105485, 0.1215294, 0.1319707, 0.1415432, 0.1499307, &
945 0.1568349, 0.1619797, 0.1651174, 0.166116, 0.1650314, 0.1619731, &
946 0.1570889, 0.1505634, 0.1426143, 0.1334867, 0.1234449, 0.1127635, &
947 0.1017171, 0.09057051, 0.07956908, 0.06893117, 0.05884206, 0.04945029, &
948 0.04086614, 0.03316217, 0.02637553, 0.0205115, 0.01554789, 0.01143988, &
949 0.00812489, 0.0055272, 0.00356223, 0.00214015, 0.00116899, 0.00055712, &
950 0.00021516, 5.741e-05, 5.75e-06, 0, 0 /
953 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, &
954 3.697e-05, 0.00043106, 0.00163591, 0.00410671, 0.00829402, 0.01463712, &
955 0.02355588, 0.03544162, 0.05064684, 0.06947458, 0.09216691, 0.1188122, &
956 0.1492688, 0.1832962, 0.2205702, 0.2606854, 0.3031641, 0.3474685, &
957 0.3930182, 0.4392108, 0.4854433, 0.5311348, 0.5757467, 0.6187996, &
958 0.659887, 0.6986829, 0.7349452, 0.7685147, 0.7993097, 0.8273188, &
959 0.8525907, 0.8752236, 0.895355, 0.913151, 0.9287973, 0.9424911, &
960 0.9544341, 0.9648276, 0.9738676, 0.9817423, 0.9886266, 0.9946712, 1 /
962 character(len=128) :: fname
963 real(kind=4),
allocatable:: wk1(:), wk2(:,:), wk3(:,:,:)
964 real,
dimension(:),
allocatable:: lat, lon, ak0, bk0
965 real,
dimension(:,:,:),
allocatable:: ud, vd
966 real,
dimension(:,:,:,:),
allocatable:: qp
967 real(kind=4),
dimension(:,:),
allocatable:: psncep, zsncep, psc
968 real(kind=4),
dimension(:,:,:),
allocatable:: uncep, vncep, tncep, zhncep
969 real(kind=4),
dimension(:,:,:,:),
allocatable:: qncep
970 real,
dimension(:,:),
allocatable:: psc_r8
971 real,
dimension(:,:,:),
allocatable:: pt_c, pt_d, gzc
972 real:: s2c(atm%bd%is:atm%bd%ie,atm%bd%js:atm%bd%je,4)
973 real:: s2c_c(atm%bd%is:atm%bd%ie+1,atm%bd%js:atm%bd%je,4)
974 real:: s2c_d(atm%bd%is:atm%bd%ie,atm%bd%js:atm%bd%je+1,4)
975 integer,
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je):: id1, id2, jdc
976 integer,
dimension(Atm%bd%is:Atm%bd%ie+1,Atm%bd%js:Atm%bd%je):: &
978 integer,
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je+1):: &
980 real :: tmean, utmp, vtmp
981 integer:: i, j, k, im, jm, km, npz, npt
982 integer:: i1, i2, j1, ncid
983 integer:: jbeg, jend, jn
985 logical:: read_ts = .true.
986 logical:: land_ts = .false.
988 integer :: is, ie, js, je
989 integer :: isd, ied, jsd, jed
990 real(kind=R_GRID),
dimension(2):: p1, p2, p3
991 real(kind=R_GRID),
dimension(3):: e1, e2, ex, ey
992 integer :: id_res, ntprog, ntracers, ks, iq, nt
1004 call get_number_tracers(model_atmos, num_tracers=ntracers, num_prog=ntprog)
1005 if(is_master())
write(*,*)
'ntracers = ', ntracers,
'ntprog = ',ntprog
1011 fname = atm%flagstruct%res_latlon_dynamics
1013 if( file_exist(fname) )
then 1019 im = tsize(1); jm = tsize(2); km = tsize(3)
1021 if(is_master())
write(*,*) fname
1022 if(is_master())
write(*,*)
' NCEP IC dimensions:', tsize
1024 allocate ( lon(im) )
1025 allocate ( lat(jm) )
1027 call _get_var1(ncid,
'lon', im, lon )
1028 call _get_var1(ncid,
'lat', jm, lat )
1038 allocate ( ak0(km+1) )
1039 allocate ( bk0(km+1) )
1044 ak0(k) = ak_hiwpp(k)
1045 bk0(k) = bk_hiwpp(k)
1048 call _get_var1(ncid,
'hyai', km+1, ak0, found )
1049 if ( .not. found ) ak0(:) = 0.
1051 call _get_var1(ncid,
'hybi', km+1, bk0 )
1053 if( is_master() )
then 1055 write(*,*) k, ak0(k), bk0(k)
1060 ak0(:) = ak0(:) * 1.e5
1063 if ( bk0(1) < 1.e-9 ) ak0(1) = max(1.e-9, ak0(1))
1066 call mpp_error(fatal,
'==> Error in get_external_ic: Expected file '//trim(fname)//
' for NCEP IC does not exist')
1070 call remap_coef( is, ie, js, je, isd, ied, jsd, jed, &
1071 im, jm, lon, lat, id1, id2, jdc, s2c , atm%gridstruct%agrid)
1074 jbeg = jm-1; jend = 2
1078 jbeg = min(jbeg, j1)
1079 jend = max(jend, j1+1)
1083 if(is_master())
write(*,*)
'jbeg, jend = ', jbeg, jend
1085 allocate ( psncep(im,jbeg:jend) )
1086 allocate ( zsncep(im,jbeg:jend) )
1088 call get_var3_r4( ncid,
'PS', 1,im, jbeg,jend, 1,1, psncep )
1089 if(is_master())
write(*,*)
'done reading psncep' 1090 call get_var3_r4( ncid,
'PHIS', 1,im, jbeg,jend, 1,1, zsncep )
1091 zsncep(:,:) = zsncep(:,:)/grav
1092 if(is_master())
write(*,*)
'done reading zsncep' 1094 allocate ( tncep(1:im,jbeg:jend, 1:km) )
1095 call get_var3_r4( ncid,
'T', 1,im, jbeg,jend, 1,km, tncep )
1096 if(is_master())
write(*,*)
'done reading tncep' 1098 allocate ( wk3(1:im,jbeg:jend, 1:km) )
1099 allocate ( qncep(1:im,jbeg:jend, 1:km,2) )
1100 call get_var3_r4( ncid,
'Q', 1,im, jbeg,jend, 1,km, wk3 )
1101 if(is_master())
write(*,*)
'done reading sphumncep' 1102 qncep(:,:,:,1) = wk3(:,:,:)
1103 call get_var3_r4( ncid,
'CWAT', 1,im, jbeg,jend, 1,km, wk3 )
1104 if(is_master())
write(*,*)
'done reading cwatncep' 1105 qncep(:,:,:,2) = wk3(:,:,:)
1114 tncep(i,j,k) = tncep(i,j,k)/(1.+
zvir*qncep(i,j,k,1))
1121 allocate ( zhncep(1:im,jbeg:jend, km+1) )
1122 jn = jend - jbeg + 1
1124 call compute_zh(im, jn, km, ak0, bk0, psncep, zsncep, tncep, qncep, 2, zhncep )
1128 if(is_master())
write(*,*)
'done compute zhncep' 1131 allocate (psc(is:ie,js:je))
1132 allocate (psc_r8(is:ie,js:je))
1139 psc(i,j) = s2c(i,j,1)*psncep(i1,j1 ) + s2c(i,j,2)*psncep(i2,j1 ) + &
1140 s2c(i,j,3)*psncep(i2,j1+1) + s2c(i,j,4)*psncep(i1,j1+1)
1143 deallocate ( psncep )
1146 allocate (gzc(is:ie,js:je,km+1))
1153 gzc(i,j,k) = s2c(i,j,1)*zhncep(i1,j1 ,k) + s2c(i,j,2)*zhncep(i2,j1 ,k) + &
1154 s2c(i,j,3)*zhncep(i2,j1+1,k) + s2c(i,j,4)*zhncep(i1,j1+1,k)
1158 deallocate ( zhncep )
1160 if(is_master())
write(*,*)
'done interpolate psncep/zhncep into cubic grid psc/gzc!' 1163 allocate ( wk2(im,jm) )
1169 if ( .not. land_ts )
then 1170 allocate ( wk1(im) )
1174 call get_var3_r4( ncid,
'ORO', 1,im, j,j, 1,1, wk1 )
1178 if( abs(wk1(i)-1.) > 0.99 )
then 1179 tmean = tmean + wk2(i,j)
1186 if ( npt /= 0 )
then 1187 tmean= tmean /
real(npt)
1189 if( abs(wk1(i)-1.) <= 0.99 )
then 1192 elseif ( i==im )
then 1197 if ( abs(wk1(i2)-1.)>0.99 )
then 1198 wk2(i,j) = wk2(i2,j)
1199 elseif ( abs(wk1(i1)-1.)>0.99 )
then 1200 wk2(i,j) = wk2(i1,j)
1216 atm%ts(i,j) = s2c(i,j,1)*wk2(i1,j1 ) + s2c(i,j,2)*wk2(i2,j1 ) + &
1217 s2c(i,j,3)*wk2(i2,j1+1) + s2c(i,j,4)*wk2(i1,j1+1)
1220 call prt_maxmin(
'SST_model', atm%ts, is, ie, js, je, 0, 1, 1.)
1224 call ncep2fms(im, jm, lon, lat, wk2)
1225 if( is_master() )
then 1226 write(*,*)
'External_ic_mod: i_sst=', i_sst,
' j_sst=', j_sst
1227 call pmaxmin(
'SST_ncep_fms', sst_ncep, i_sst, j_sst, 1.)
1235 allocate ( qp(is:ie,js:je,km,2) )
1242 qp(i,j,k,1) = s2c(i,j,1)*qncep(i1,j1 ,k,1) + s2c(i,j,2)*qncep(i2,j1 ,k,1) + &
1243 s2c(i,j,3)*qncep(i2,j1+1,k,1) + s2c(i,j,4)*qncep(i1,j1+1,k,1)
1244 qp(i,j,k,2) = s2c(i,j,1)*qncep(i1,j1 ,k,2) + s2c(i,j,2)*qncep(i2,j1 ,k,2) + &
1245 s2c(i,j,3)*qncep(i2,j1+1,k,2) + s2c(i,j,4)*qncep(i1,j1+1,k,2)
1252 psc_r8(:,:) = psc(:,:)
1256 call remap_scalar(atm, km, npz, 2, ak0, bk0, psc_r8, qp, gzc)
1257 call mpp_update_domains(atm%phis, atm%domain)
1258 if(is_master())
write(*,*)
'done remap_scalar' 1264 allocate (pt_c(isd:ied+1,jsd:jed ,2))
1265 allocate (pt_d(isd:ied ,jsd:jed+1,2))
1266 allocate (ud(is:ie , js:je+1, km))
1267 allocate (vd(is:ie+1, js:je , km))
1270 isd, ied, jsd, jed, &
1271 atm%gridstruct%grid, pt_c, pt_d)
1275 call remap_coef( is, ie+1, js, je, isd, ied+1, jsd, jed, &
1276 im, jm, lon, lat, id1_c, id2_c, jdc_c, s2c_c, pt_c)
1279 jbeg = jm-1; jend = 2
1283 jbeg = min(jbeg, j1)
1284 jend = max(jend, j1+1)
1289 allocate ( uncep(1:im,jbeg:jend, 1:km) )
1290 allocate ( vncep(1:im,jbeg:jend, 1:km) )
1292 call get_var3_r4( ncid,
'U', 1,im, jbeg,jend, 1,km, uncep )
1293 if(is_master())
write(*,*)
'first time done reading Uncep' 1294 call get_var3_r4( ncid,
'V', 1,im, jbeg,jend, 1,km, vncep )
1295 if(is_master())
write(*,*)
'first time done reading Vncep' 1305 p1(:) = atm%gridstruct%grid(i,j ,1:2)
1306 p2(:) = atm%gridstruct%grid(i,j+1,1:2)
1310 utmp = s2c_c(i,j,1)*uncep(i1,j1 ,k) + &
1311 s2c_c(i,j,2)*uncep(i2,j1 ,k) + &
1312 s2c_c(i,j,3)*uncep(i2,j1+1,k) + &
1313 s2c_c(i,j,4)*uncep(i1,j1+1,k)
1314 vtmp = s2c_c(i,j,1)*vncep(i1,j1 ,k) + &
1315 s2c_c(i,j,2)*vncep(i2,j1 ,k) + &
1316 s2c_c(i,j,3)*vncep(i2,j1+1,k) + &
1317 s2c_c(i,j,4)*vncep(i1,j1+1,k)
1323 deallocate ( uncep, vncep )
1327 call remap_coef( is, ie, js, je+1, isd, ied, jsd, jed+1, &
1328 im, jm, lon, lat, id1_d, id2_d, jdc_d, s2c_d, pt_d)
1329 deallocate ( pt_c, pt_d )
1332 jbeg = jm-1; jend = 2
1336 jbeg = min(jbeg, j1)
1337 jend = max(jend, j1+1)
1342 allocate ( uncep(1:im,jbeg:jend, 1:km) )
1343 allocate ( vncep(1:im,jbeg:jend, 1:km) )
1345 call get_var3_r4( ncid,
'U', 1,im, jbeg,jend, 1,km, uncep )
1346 if(is_master())
write(*,*)
'second time done reading uec' 1348 call get_var3_r4( ncid,
'V', 1,im, jbeg,jend, 1,km, vncep )
1349 if(is_master())
write(*,*)
'second time done reading vec' 1359 p1(:) = atm%gridstruct%grid(i, j,1:2)
1360 p2(:) = atm%gridstruct%grid(i+1,j,1:2)
1364 utmp = s2c_d(i,j,1)*uncep(i1,j1 ,k) + &
1365 s2c_d(i,j,2)*uncep(i2,j1 ,k) + &
1366 s2c_d(i,j,3)*uncep(i2,j1+1,k) + &
1367 s2c_d(i,j,4)*uncep(i1,j1+1,k)
1368 vtmp = s2c_d(i,j,1)*vncep(i1,j1 ,k) + &
1369 s2c_d(i,j,2)*vncep(i2,j1 ,k) + &
1370 s2c_d(i,j,3)*vncep(i2,j1+1,k) + &
1371 s2c_d(i,j,4)*vncep(i1,j1+1,k)
1376 deallocate ( uncep, vncep )
1378 call remap_dwinds(km, npz, ak0, bk0, psc_r8, ud, vd, atm)
1379 deallocate ( ud, vd )
1395 use gfs_restart
, only : gfs_restart_type
1401 type(domain2d),
intent(inout) :: fv_domain
1403 real :: ak_ec(138), bk_ec(138)
1404 data ak_ec/ 0.000000, 2.000365, 3.102241, 4.666084, 6.827977, 9.746966, &
1405 13.605424, 18.608931, 24.985718, 32.985710, 42.879242, 54.955463, &
1406 69.520576, 86.895882, 107.415741, 131.425507, 159.279404, 191.338562, &
1407 227.968948, 269.539581, 316.420746, 368.982361, 427.592499, 492.616028, &
1408 564.413452, 643.339905, 729.744141, 823.967834, 926.344910, 1037.201172, &
1409 1156.853638, 1285.610352, 1423.770142, 1571.622925, 1729.448975, 1897.519287, &
1410 2076.095947, 2265.431641, 2465.770508, 2677.348145, 2900.391357, 3135.119385, &
1411 3381.743652, 3640.468262, 3911.490479, 4194.930664, 4490.817383, 4799.149414, &
1412 5119.895020, 5452.990723, 5798.344727, 6156.074219, 6526.946777, 6911.870605, &
1413 7311.869141, 7727.412109, 8159.354004, 8608.525391, 9076.400391, 9562.682617, &
1414 10065.978516, 10584.631836, 11116.662109, 11660.067383, 12211.547852, 12766.873047, &
1415 13324.668945, 13881.331055, 14432.139648, 14975.615234, 15508.256836, 16026.115234, &
1416 16527.322266, 17008.789063, 17467.613281, 17901.621094, 18308.433594, 18685.718750, &
1417 19031.289063, 19343.511719, 19620.042969, 19859.390625, 20059.931641, 20219.664063, &
1418 20337.863281, 20412.308594, 20442.078125, 20425.718750, 20361.816406, 20249.511719, &
1419 20087.085938, 19874.025391, 19608.572266, 19290.226563, 18917.460938, 18489.707031, &
1420 18006.925781, 17471.839844, 16888.687500, 16262.046875, 15596.695313, 14898.453125, &
1421 14173.324219, 13427.769531, 12668.257813, 11901.339844, 11133.304688, 10370.175781, &
1422 9617.515625, 8880.453125, 8163.375000, 7470.343750, 6804.421875, 6168.531250, &
1423 5564.382813, 4993.796875, 4457.375000, 3955.960938, 3489.234375, 3057.265625, &
1424 2659.140625, 2294.242188, 1961.500000, 1659.476563, 1387.546875, 1143.250000, &
1425 926.507813, 734.992188, 568.062500, 424.414063, 302.476563, 202.484375, &
1426 122.101563, 62.781250, 22.835938, 3.757813, 0.000000, 0.000000 /
1428 data bk_ec/ 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, &
1429 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, &
1430 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, &
1431 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, &
1432 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, &
1433 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, &
1434 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, &
1435 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, &
1436 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, 0.000000, &
1437 0.000000, 0.000007, 0.000024, 0.000059, 0.000112, 0.000199, &
1438 0.000340, 0.000562, 0.000890, 0.001353, 0.001992, 0.002857, &
1439 0.003971, 0.005378, 0.007133, 0.009261, 0.011806, 0.014816, &
1440 0.018318, 0.022355, 0.026964, 0.032176, 0.038026, 0.044548, &
1441 0.051773, 0.059728, 0.068448, 0.077958, 0.088286, 0.099462, &
1442 0.111505, 0.124448, 0.138313, 0.153125, 0.168910, 0.185689, &
1443 0.203491, 0.222333, 0.242244, 0.263242, 0.285354, 0.308598, &
1444 0.332939, 0.358254, 0.384363, 0.411125, 0.438391, 0.466003, &
1445 0.493800, 0.521619, 0.549301, 0.576692, 0.603648, 0.630036, &
1446 0.655736, 0.680643, 0.704669, 0.727739, 0.749797, 0.770798, &
1447 0.790717, 0.809536, 0.827256, 0.843881, 0.859432, 0.873929, &
1448 0.887408, 0.899900, 0.911448, 0.922096, 0.931881, 0.940860, &
1449 0.949064, 0.956550, 0.963352, 0.969513, 0.975078, 0.980072, &
1450 0.984542, 0.988500, 0.991984, 0.995003, 0.997630, 1.000000 /
1454 real,
dimension(64):: ak_sj, bk_sj
1455 data ak_sj/64.247, 137.790, 221.958, &
1456 318.266, 428.434, 554.424, &
1457 698.457, 863.05803, 1051.07995, &
1458 1265.75194, 1510.71101, 1790.05098, &
1459 2108.36604, 2470.78817, 2883.03811, &
1460 3351.46002, 3883.05187, 4485.49315, &
1461 5167.14603, 5937.04991, 6804.87379, &
1462 7780.84698, 8875.64338, 10100.20534, &
1463 11264.35673, 12190.64366, 12905.42546, &
1464 13430.87867, 13785.88765, 13986.77987, &
1465 14047.96335, 13982.46770, 13802.40331, &
1466 13519.33841, 13144.59486, 12689.45608, &
1467 12165.28766, 11583.57006, 10955.84778, &
1468 10293.60402, 9608.08306, 8910.07678, &
1469 8209.70131, 7516.18560, 6837.69250, &
1470 6181.19473, 5552.39653, 4955.72632, &
1471 4394.37629, 3870.38682, 3384.76586, &
1472 2937.63489, 2528.37666, 2155.78385, &
1473 1818.20722, 1513.68173, 1240.03585, &
1474 994.99144, 776.23591, 581.48797, &
1475 408.53400, 255.26520, 119.70243, 0. /
1477 data bk_sj/0.00000, 0.00000, 0.00000, &
1478 0.00000, 0.00000, 0.00000, &
1479 0.00000, 0.00000, 0.00000, &
1480 0.00000, 0.00000, 0.00000, &
1481 0.00000, 0.00000, 0.00000, &
1482 0.00000, 0.00000, 0.00000, &
1483 0.00000, 0.00000, 0.00000, &
1484 0.00000, 0.00000, 0.00000, &
1485 0.00201, 0.00792, 0.01755, &
1486 0.03079, 0.04751, 0.06761, &
1487 0.09097, 0.11746, 0.14690, &
1488 0.17911, 0.21382, 0.25076, &
1489 0.28960, 0.32994, 0.37140, &
1490 0.41353, 0.45589, 0.49806, &
1491 0.53961, 0.58015, 0.61935, &
1492 0.65692, 0.69261, 0.72625, &
1493 0.75773, 0.78698, 0.81398, &
1494 0.83876, 0.86138, 0.88192, &
1495 0.90050, 0.91722, 0.93223, &
1496 0.94565, 0.95762, 0.96827, &
1497 0.97771, 0.98608, 0.99347, 1./
1499 character(len=128) :: fname
1500 real,
allocatable:: wk2(:,:)
1501 real(kind=4),
allocatable:: wk2_r4(:,:)
1502 real,
dimension(:,:,:),
allocatable:: ud, vd
1503 real,
allocatable:: wc(:,:,:)
1504 real(kind=4),
allocatable:: uec(:,:,:), vec(:,:,:), tec(:,:,:), wec(:,:,:)
1505 real(kind=4),
allocatable:: psec(:,:), zsec(:,:), zhec(:,:,:), qec(:,:,:,:)
1506 real(kind=4),
allocatable:: psc(:,:)
1507 real(kind=4),
allocatable:: sphumec(:,:,:)
1508 real,
allocatable:: psc_r8(:,:), zhc(:,:,:), qc(:,:,:,:)
1509 real,
allocatable:: lat(:), lon(:), ak0(:), bk0(:)
1510 real,
allocatable:: pt_c(:,:,:), pt_d(:,:,:)
1511 real:: s2c(atm%bd%is:atm%bd%ie,atm%bd%js:atm%bd%je,4)
1512 real:: s2c_c(atm%bd%is:atm%bd%ie+1,atm%bd%js:atm%bd%je,4)
1513 real:: s2c_d(atm%bd%is:atm%bd%ie,atm%bd%js:atm%bd%je+1,4)
1514 integer,
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je):: &
1516 integer,
dimension(Atm%bd%is:Atm%bd%ie+1,Atm%bd%js:Atm%bd%je):: &
1518 integer,
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je+1):: &
1521 integer:: i, j, k, n, im, jm, km, npz, npt
1522 integer:: i1, i2, j1, ncid
1523 integer:: jbeg, jend, jn
1525 logical:: read_ts = .true.
1526 logical:: land_ts = .false.
1528 integer :: is, ie, js, je
1529 integer :: isd, ied, jsd, jed
1530 integer :: sphum, liq_wat, ice_wat, rainwat, snowwat, graupel
1532 integer :: spfo, spfo2, spfo3
1536 real:: wt, qt, m_fac
1537 real(kind=8) :: scale_value, offset, ptmp
1538 real(kind=R_GRID),
dimension(2):: p1, p2, p3
1539 real(kind=R_GRID),
dimension(3):: e1, e2, ex, ey
1540 real,
allocatable:: ps_gfs(:,:), zh_gfs(:,:,:)
1542 real,
allocatable:: spfo_gfs(:,:,:), spfo2_gfs(:,:,:), spfo3_gfs(:,:,:)
1544 real,
allocatable:: o3mr_gfs(:,:,:)
1546 real,
allocatable:: ak_gfs(:), bk_gfs(:)
1547 integer :: id_res, ntprog, ntracers, ks, iq, nt
1548 character(len=64) :: tracer_name
1549 integer :: levp_gfs = 64
1550 type(restart_file_type) :: oro_restart, gfs_restart
1551 character(len=64) :: fn_oro_ics =
'oro_data.nc' 1552 character(len=64) :: fn_gfs_ics =
'gfs_data.nc' 1553 character(len=64) :: fn_gfs_ctl =
'gfs_ctrl.nc' 1554 logical :: filtered_terrain = .true.
1555 namelist /external_ic_nml/ filtered_terrain
1567 call get_number_tracers(model_atmos, num_tracers=ntracers, num_prog=ntprog)
1568 if(is_master())
write(*,*)
'ntracers = ', ntracers,
'ntprog = ',ntprog
1570 sphum = get_tracer_index(model_atmos,
'sphum')
1571 liq_wat = get_tracer_index(model_atmos,
'liq_wat')
1572 ice_wat = get_tracer_index(model_atmos,
'ice_wat')
1573 rainwat = get_tracer_index(model_atmos,
'rainwat')
1574 snowwat = get_tracer_index(model_atmos,
'snowwat')
1575 graupel = get_tracer_index(model_atmos,
'graupel')
1577 spfo = get_tracer_index(model_atmos,
'spfo')
1578 spfo2 = get_tracer_index(model_atmos,
'spfo2')
1579 spfo3 = get_tracer_index(model_atmos,
'spfo3')
1581 o3mr = get_tracer_index(model_atmos,
'o3mr')
1584 if (is_master())
then 1585 print *,
'sphum = ', sphum
1586 print *,
'liq_wat = ', liq_wat
1587 if ( atm%flagstruct%nwat .eq. 6 )
then 1588 print *,
'rainwat = ', rainwat
1589 print *,
'iec_wat = ', ice_wat
1590 print *,
'snowwat = ', snowwat
1591 print *,
'graupel = ', graupel
1594 print *,
' spfo3 = ', spfo3
1595 print *,
' spfo = ', spfo
1596 print *,
' spfo2 = ', spfo2
1598 print *,
' o3mr = ', o3mr
1604 if (atm%flagstruct%external_eta)
then 1609 if (filtered_terrain)
then 1610 id_res = register_restart_field(oro_restart, fn_oro_ics,
'orog_filt', atm%phis, domain=atm%domain)
1611 elseif (.not. filtered_terrain)
then 1612 id_res = register_restart_field(oro_restart, fn_oro_ics,
'orog_raw', atm%phis, domain=atm%domain)
1614 call restore_state (oro_restart)
1615 call free_restart_type(oro_restart)
1616 atm%phis = atm%phis*grav
1617 if(is_master())
write(*,*)
'done reading model terrain from oro_data.nc' 1618 call mpp_update_domains( atm%phis, atm%domain )
1622 allocate (spfo3_gfs(is:ie,js:je,levp_gfs))
1623 allocate ( spfo_gfs(is:ie,js:je,levp_gfs))
1624 allocate (spfo2_gfs(is:ie,js:je,levp_gfs))
1626 allocate (o3mr_gfs(is:ie,js:je,levp_gfs))
1628 allocate (ps_gfs(is:ie,js:je))
1629 allocate (zh_gfs(is:ie,js:je,levp_gfs+1))
1632 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'spfo3', spfo3_gfs, &
1633 mandatory=.false.,domain=atm%domain)
1634 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'spfo', spfo_gfs, &
1635 mandatory=.false.,domain=atm%domain)
1636 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'spfo2', spfo2_gfs, &
1637 mandatory=.false.,domain=atm%domain)
1639 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'o3mr', o3mr_gfs, &
1640 mandatory=.false.,domain=atm%domain)
1642 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'ps', ps_gfs, domain=atm%domain)
1643 id_res = register_restart_field(gfs_restart, fn_gfs_ics,
'ZH', zh_gfs, domain=atm%domain)
1644 call restore_state (gfs_restart)
1645 call free_restart_type(gfs_restart)
1649 allocate (wk2(levp_gfs+1,2))
1650 allocate (ak_gfs(levp_gfs+1))
1651 allocate (bk_gfs(levp_gfs+1))
1652 call read_data(
'INPUT/'//trim(fn_gfs_ctl),
'vcoord',wk2, no_domain=.true.)
1653 ak_gfs(1:levp_gfs+1) = wk2(1:levp_gfs+1,1)
1654 bk_gfs(1:levp_gfs+1) = wk2(1:levp_gfs+1,2)
1657 if ( bk_gfs(1) < 1.e-9 ) ak_gfs(1) = max(1.e-9, ak_gfs(1))
1661 if(is_master())
write(*,*)
'Reading spfo3 from GFS_data.nc:' 1662 if(is_master())
write(*,*)
'spfo3 =', iq
1663 call remap_scalar_single(atm, levp_gfs, npz, ak_gfs, bk_gfs, ps_gfs, spfo3_gfs, zh_gfs, iq)
1665 if(is_master())
write(*,*)
'Reading spfo from GFS_data.nc:' 1666 if(is_master())
write(*,*)
'spfo =', iq
1667 call remap_scalar_single(atm, levp_gfs, npz, ak_gfs, bk_gfs, ps_gfs, spfo_gfs, zh_gfs, iq)
1669 if(is_master())
write(*,*)
'Reading spfo2 from GFS_data.nc:' 1670 if(is_master())
write(*,*)
'spfo2 =', iq
1671 call remap_scalar_single(atm, levp_gfs, npz, ak_gfs, bk_gfs, ps_gfs, spfo2_gfs, zh_gfs, iq)
1674 if(is_master())
write(*,*)
'Reading o3mr from GFS_data.nc:' 1675 if(is_master())
write(*,*)
'o3mr =', iq
1676 call remap_scalar_single(atm, levp_gfs, npz, ak_gfs, bk_gfs, ps_gfs, o3mr_gfs, zh_gfs, iq)
1679 deallocate (ak_gfs, bk_gfs)
1680 deallocate (ps_gfs, zh_gfs)
1682 deallocate (spfo3_gfs)
1683 deallocate ( spfo_gfs)
1684 deallocate (spfo2_gfs)
1686 deallocate (o3mr_gfs)
1690 fname = atm%flagstruct%res_latlon_dynamics
1692 if( file_exist(fname) )
then 1695 call get_ncdim1( ncid,
'longitude', tsize(1) )
1696 call get_ncdim1( ncid,
'latitude', tsize(2) )
1699 im = tsize(1); jm = tsize(2); km = tsize(3)
1701 if(is_master())
write(*,*) fname
1702 if(is_master())
write(*,*)
' ECMWF IC dimensions:', tsize
1704 allocate ( lon(im) )
1705 allocate ( lat(jm) )
1707 call _get_var1(ncid,
'longitude', im, lon )
1708 call _get_var1(ncid,
'latitude', jm, lat )
1718 allocate ( ak0(km+1) )
1719 allocate ( bk0(km+1) )
1727 if( is_master() )
then 1729 write(*,*) k, ak0(k), bk0(k)
1734 if ( bk0(1) < 1.e-9 ) ak0(1) = max(1.e-9, ak0(1))
1737 call mpp_error(fatal,
'==> Error in get_external_ic: Expected file '//trim(fname)//
' for NCEP IC does not exist')
1741 call remap_coef( is, ie, js, je, isd, ied, jsd, jed, &
1742 im, jm, lon, lat, id1, id2, jdc, s2c , atm%gridstruct%agrid )
1745 jbeg = jm-1; jend = 2
1749 jbeg = min(jbeg, j1)
1750 jend = max(jend, j1+1)
1754 if(is_master())
write(*,*)
'jbeg, jend = ', jbeg, jend
1756 allocate ( psec(im,jbeg:jend) )
1757 allocate ( zsec(im,jbeg:jend) )
1758 allocate ( wk2_r4(im,jbeg:jend) )
1760 call get_var2_r4( ncid,
'lnsp', 1,im, jbeg,jend, wk2_r4 )
1763 psec(:,:) = exp(wk2_r4(:,:)*scale_value + offset)
1764 if(is_master())
write(*,*)
'done reading psec' 1766 call get_var2_r4( ncid,
'z', 1,im, jbeg,jend, wk2_r4 )
1769 zsec(:,:) = (wk2_r4(:,:)*scale_value + offset)/grav
1770 if(is_master())
write(*,*)
'done reading zsec' 1772 deallocate ( wk2_r4 )
1775 allocate ( tec(1:im,jbeg:jend, 1:km) )
1777 call get_var3_r4( ncid,
't', 1,im, jbeg,jend, 1,km, tec )
1780 tec(:,:,:) = tec(:,:,:)*scale_value + offset
1781 if(is_master())
write(*,*)
'done reading tec' 1784 allocate ( sphumec(1:im,jbeg:jend, 1:km) )
1786 call get_var3_r4( ncid,
'q', 1,im, jbeg,jend, 1,km, sphumec(:,:,:) )
1789 sphumec(:,:,:) = sphumec(:,:,:)*scale_value + offset
1790 if(is_master())
write(*,*)
'done reading sphum ec' 1793 allocate ( qec(1:im,jbeg:jend,1:km,5) )
1796 if (n == sphum)
then 1797 qec(:,:,:,sphum) = sphumec(:,:,:)
1798 deallocate ( sphumec )
1799 else if (n == liq_wat)
then 1800 call get_var3_r4( ncid,
'clwc', 1,im, jbeg,jend, 1,km, qec(:,:,:,liq_wat) )
1803 qec(:,:,:,liq_wat) = qec(:,:,:,liq_wat)*scale_value + offset
1804 if(is_master())
write(*,*)
'done reading clwc ec' 1805 else if (n == rainwat)
then 1806 call get_var3_r4( ncid,
'crwc', 1,im, jbeg,jend, 1,km, qec(:,:,:,rainwat) )
1809 qec(:,:,:,rainwat) = qec(:,:,:,rainwat)*scale_value + offset
1810 if(is_master())
write(*,*)
'done reading crwc ec' 1811 else if (n == ice_wat)
then 1812 call get_var3_r4( ncid,
'ciwc', 1,im, jbeg,jend, 1,km, qec(:,:,:,ice_wat) )
1815 qec(:,:,:,ice_wat) = qec(:,:,:,ice_wat)*scale_value + offset
1816 if(is_master())
write(*,*)
'done reading ciwc ec' 1817 else if (n == snowwat)
then 1818 call get_var3_r4( ncid,
'cswc', 1,im, jbeg,jend, 1,km, qec(:,:,:,snowwat) )
1821 qec(:,:,:,snowwat) = qec(:,:,:,snowwat)*scale_value + offset
1822 if(is_master())
write(*,*)
'done reading cswc ec' 1824 if(is_master())
write(*,*)
'nq is more then 5!' 1831 allocate ( zhec(1:im,jbeg:jend, km+1) )
1832 jn = jend - jbeg + 1
1834 call compute_zh(im, jn, km, ak0, bk0, psec, zsec, tec, qec, 5, zhec )
1835 if(is_master())
write(*,*)
'done compute zhec' 1840 allocate (psc(is:ie,js:je))
1841 allocate (psc_r8(is:ie,js:je))
1846 psec(i,j) = log(psec(i,j))
1856 ptmp = s2c(i,j,1)*psec(i1,j1 ) + s2c(i,j,2)*psec(i2,j1 ) + &
1857 s2c(i,j,3)*psec(i2,j1+1) + s2c(i,j,4)*psec(i1,j1+1)
1858 psc(i,j) = exp(ptmp)
1860 psc(i,j) = s2c(i,j,1)*psec(i1,j1 ) + s2c(i,j,2)*psec(i2,j1 ) + &
1861 s2c(i,j,3)*psec(i2,j1+1) + s2c(i,j,4)*psec(i1,j1+1)
1867 allocate (zhc(is:ie,js:je,km+1))
1876 zhc(i,j,k) = s2c(i,j,1)*zhec(i1,j1 ,k) + s2c(i,j,2)*zhec(i2,j1 ,k) + &
1877 s2c(i,j,3)*zhec(i2,j1+1,k) + s2c(i,j,4)*zhec(i1,j1+1,k)
1883 if(is_master())
write(*,*)
'done interpolate psec/zhec into cubic grid psc/zhc!' 1886 allocate ( qc(is:ie,js:je,km,6) )
1897 qc(i,j,k,n) = s2c(i,j,1)*qec(i1,j1 ,k,n) + s2c(i,j,2)*qec(i2,j1 ,k,n) + &
1898 s2c(i,j,3)*qec(i2,j1+1,k,n) + s2c(i,j,4)*qec(i1,j1+1,k,n)
1904 qc(:,:,:,graupel) = 0.
1907 if(is_master())
write(*,*)
'done interpolate tracers (qec) into cubic (qc)' 1910 allocate ( wec(1:im,jbeg:jend, 1:km) )
1911 allocate ( wc(is:ie,js:je,km))
1913 call get_var3_r4( ncid,
'w', 1,im, jbeg,jend, 1,km, wec )
1916 wec(:,:,:) = wec(:,:,:)*scale_value + offset
1927 wc(i,j,k) = s2c(i,j,1)*wec(i1,j1 ,k) + s2c(i,j,2)*wec(i2,j1 ,k) + &
1928 s2c(i,j,3)*wec(i2,j1+1,k) + s2c(i,j,4)*wec(i1,j1+1,k)
1935 if(is_master())
write(*,*)
'done reading and interpolate vertical wind (w) into cubic' 1938 psc_r8(:,:) = psc(:,:)
1941 call remap_scalar(atm, km, npz, 6, ak0, bk0, psc_r8, qc, zhc, wc)
1942 call mpp_update_domains(atm%phis, atm%domain)
1943 if(is_master())
write(*,*)
'done remap_scalar' 1951 allocate (pt_c(isd:ied+1,jsd:jed ,2))
1952 allocate (pt_d(isd:ied ,jsd:jed+1,2))
1953 allocate (ud(is:ie , js:je+1, km))
1954 allocate (vd(is:ie+1, js:je , km))
1957 isd, ied, jsd, jed, &
1958 atm%gridstruct%grid, pt_c, pt_d)
1962 call remap_coef( is, ie+1, js, je, isd, ied+1, jsd, jed, &
1963 im, jm, lon, lat, id1_c, id2_c, jdc_c, s2c_c, pt_c)
1966 jbeg = jm-1; jend = 2
1970 jbeg = min(jbeg, j1)
1971 jend = max(jend, j1+1)
1976 allocate ( uec(1:im,jbeg:jend, 1:km) )
1977 allocate ( vec(1:im,jbeg:jend, 1:km) )
1979 call get_var3_r4( ncid,
'u', 1,im, jbeg,jend, 1,km, uec )
1985 uec(i,j,k) = uec(i,j,k)*scale_value + offset
1989 if(is_master())
write(*,*)
'first time done reading uec' 1991 call get_var3_r4( ncid,
'v', 1,im, jbeg,jend, 1,km, vec )
1997 vec(i,j,k) = vec(i,j,k)*scale_value + offset
2002 if(is_master())
write(*,*)
'first time done reading vec' 2012 p1(:) = atm%gridstruct%grid(i,j ,1:2)
2013 p2(:) = atm%gridstruct%grid(i,j+1,1:2)
2017 utmp = s2c_c(i,j,1)*uec(i1,j1 ,k) + &
2018 s2c_c(i,j,2)*uec(i2,j1 ,k) + &
2019 s2c_c(i,j,3)*uec(i2,j1+1,k) + &
2020 s2c_c(i,j,4)*uec(i1,j1+1,k)
2021 vtmp = s2c_c(i,j,1)*vec(i1,j1 ,k) + &
2022 s2c_c(i,j,2)*vec(i2,j1 ,k) + &
2023 s2c_c(i,j,3)*vec(i2,j1+1,k) + &
2024 s2c_c(i,j,4)*vec(i1,j1+1,k)
2030 deallocate ( uec, vec )
2034 call remap_coef( is, ie, js, je+1, isd, ied, jsd, jed+1, &
2035 im, jm, lon, lat, id1_d, id2_d, jdc_d, s2c_d, pt_d)
2036 deallocate ( pt_c, pt_d )
2039 jbeg = jm-1; jend = 2
2043 jbeg = min(jbeg, j1)
2044 jend = max(jend, j1+1)
2049 allocate ( uec(1:im,jbeg:jend, 1:km) )
2050 allocate ( vec(1:im,jbeg:jend, 1:km) )
2052 call get_var3_r4( ncid,
'u', 1,im, jbeg,jend, 1,km, uec )
2055 uec(:,:,:) = uec(:,:,:)*scale_value + offset
2056 if(is_master())
write(*,*)
'second time done reading uec' 2058 call get_var3_r4( ncid,
'v', 1,im, jbeg,jend, 1,km, vec )
2061 vec(:,:,:) = vec(:,:,:)*scale_value + offset
2062 if(is_master())
write(*,*)
'second time done reading vec' 2072 p1(:) = atm%gridstruct%grid(i, j,1:2)
2073 p2(:) = atm%gridstruct%grid(i+1,j,1:2)
2077 utmp = s2c_d(i,j,1)*uec(i1,j1 ,k) + &
2078 s2c_d(i,j,2)*uec(i2,j1 ,k) + &
2079 s2c_d(i,j,3)*uec(i2,j1+1,k) + &
2080 s2c_d(i,j,4)*uec(i1,j1+1,k)
2081 vtmp = s2c_d(i,j,1)*vec(i1,j1 ,k) + &
2082 s2c_d(i,j,2)*vec(i2,j1 ,k) + &
2083 s2c_d(i,j,3)*vec(i2,j1+1,k) + &
2084 s2c_d(i,j,4)*vec(i1,j1+1,k)
2089 deallocate ( uec, vec )
2091 call remap_dwinds(km, npz, ak0, bk0, psc_r8, ud, vd, atm)
2092 deallocate ( ud, vd )
2100 wt = atm%delp(i,j,k)
2101 if ( atm%flagstruct%nwat .eq. 2 )
then 2102 qt = wt*(1.+atm%q(i,j,k,liq_wat))
2103 elseif ( atm%flagstruct%nwat .eq. 6 )
then 2104 qt = wt*(1. + atm%q(i,j,k,liq_wat) + &
2105 atm%q(i,j,k,ice_wat) + &
2106 atm%q(i,j,k,rainwat) + &
2107 atm%q(i,j,k,snowwat) + &
2108 atm%q(i,j,k,graupel))
2112 atm%q(i,j,k,iq) = m_fac * atm%q(i,j,k,iq)
2114 atm%delp(i,j,k) = qt
2120 deallocate ( ak0, bk0 )
2122 deallocate ( psc_r8 )
2123 deallocate ( lat, lon )
2125 atm%flagstruct%make_nh = .false.
2130 subroutine get_fv_ic( Atm, fv_domain, nq )
2132 type(domain2d),
intent(inout) :: fv_domain
2133 integer,
intent(in):: nq
2135 character(len=128) :: fname, tracer_name
2136 real,
allocatable:: ps0(:,:), gz0(:,:), u0(:,:,:), v0(:,:,:), t0(:,:,:), dp0(:,:,:), q0(:,:,:,:)
2137 real,
allocatable:: ua(:,:,:), va(:,:,:)
2138 real,
allocatable:: lat(:), lon(:), ak0(:), bk0(:)
2139 integer :: i, j, k, im, jm, km, npz, tr_ind
2150 fname = atm%flagstruct%res_latlon_dynamics
2152 if( file_exist(fname) )
then 2153 call field_size(fname,
'T', tsize, field_found=found)
2154 if(is_master())
write(*,*)
'Using lat-lon FV restart:', fname
2157 im = tsize(1); jm = tsize(2); km = tsize(3)
2158 if(is_master())
write(*,*)
'External IC dimensions:', tsize
2160 call mpp_error(fatal,
'==> Error in get_external_ic: field not found')
2164 allocate ( lon(im) )
2165 allocate ( lat(jm) )
2168 lon(i) = (0.5 +
real(i-1)) * 2.*pi/
real(im)
2172 lat(j) = -0.5*pi +
real(j-1)*pi/
real(jm-1) 2175 allocate ( ak0(1:km+1) )
2176 allocate ( bk0(1:km+1) )
2177 allocate ( ps0(1:im,1:jm) )
2178 allocate ( gz0(1:im,1:jm) )
2179 allocate ( u0(1:im,1:jm,1:km) )
2180 allocate ( v0(1:im,1:jm,1:km) )
2181 allocate ( t0(1:im,1:jm,1:km) )
2182 allocate ( dp0(1:im,1:jm,1:km) )
2184 call read_data (fname,
'ak', ak0)
2185 call read_data (fname,
'bk', bk0)
2186 call read_data (fname,
'Surface_geopotential', gz0)
2187 call read_data (fname,
'U', u0)
2188 call read_data (fname,
'V', v0)
2189 call read_data (fname,
'T', t0)
2190 call read_data (fname,
'DELP', dp0)
2193 if(is_master())
call pmaxmin(
'ZS_data', gz0, im, jm, 1./grav)
2194 if(mpp_pe()==1)
call pmaxmin(
'U_data', u0, im*jm, km, 1.)
2195 if(mpp_pe()==1)
call pmaxmin(
'V_data', v0, im*jm, km, 1.)
2196 if(mpp_pe()==2)
call pmaxmin(
'T_data', t0, im*jm, km, 1.)
2197 if(mpp_pe()==3)
call pmaxmin(
'DEL-P', dp0, im*jm, km, 0.01)
2201 call mpp_error(fatal,
'==> Error in get_external_ic: Expected file '//trim(fname)//
' for dynamics does not exist')
2205 fname = atm%flagstruct%res_latlon_tracers
2207 if( file_exist(fname) )
then 2208 if(is_master())
write(*,*)
'Using lat-lon tracer restart:', fname
2210 allocate ( q0(im,jm,km,atm%ncnst) )
2214 call get_tracer_names(model_atmos, tr_ind, tracer_name)
2215 if (field_exist(fname,tracer_name))
then 2216 call read_data(fname, tracer_name, q0(1:im,1:jm,1:km,tr_ind))
2217 call mpp_error(note,
'==> Have read tracer '//trim(tracer_name)//
' from '//trim(fname))
2222 call mpp_error(fatal,
'==> Error in get_external_ic: Expected file '//trim(fname)//
' for tracers does not exist')
2226 allocate ( ua(im,jm,km) )
2227 allocate ( va(im,jm,km) )
2229 call d2a3d(u0, v0, ua, va, im, jm, km, lon)
2234 if(mpp_pe()==4)
call pmaxmin(
'UA', ua, im*jm, km, 1.)
2235 if(mpp_pe()==4)
call pmaxmin(
'VA', va, im*jm, km, 1.)
2246 ps0(i,j) = ps0(i,j) + dp0(i,j,k)
2251 if (is_master())
call pmaxmin(
'PS_data (mb)', ps0, im, jm, 0.01)
2256 call remap_xyz( im, 1, jm, jm, km, npz, nq, atm%ncnst, lon, lat, ak0, bk0, &
2257 ps0, gz0, ua, va, t0, q0, atm )
2275 subroutine ncep2fms(im, jm, lon, lat, wk)
2277 integer,
intent(in):: im, jm
2278 real,
intent(in):: lon(im), lat(jm)
2279 real(kind=4),
intent(in):: wk(im,jm)
2286 real:: c1, c2, c3, c4
2287 integer i,j, i1, i2, jc, i0, j0, it, jt
2290 rdlon(i) = 1. / (lon(i+1) - lon(i))
2292 rdlon(im) = 1. / (lon(1) + 2.*pi - lon(im))
2295 rdlat(j) = 1. / (lat(j+1) - lat(j))
2302 delx = 360./
real(i_sst)
2303 dely = 180./
real(j_sst)
2308 yc = (-90. + dely * (0.5+
real(j-1))) *
deg2rad 2309 if ( yc<lat(1) )
then 2312 elseif ( yc>lat(jm) )
then 2317 if ( yc>=lat(j0) .and. yc<=lat(j0+1) )
then 2320 b1 = (yc-lat(jc)) * rdlat(jc)
2329 xc = delx * (0.5+
real(i-1)) *
deg2rad 2330 if ( xc>lon(im) )
then 2332 a1 = (xc-lon(im)) * rdlon(im)
2333 elseif ( xc<lon(1) )
then 2335 a1 = (xc+2.*pi-lon(im)) * rdlon(im)
2338 if ( xc>=lon(i0) .and. xc<=lon(i0+1) )
then 2341 a1 = (xc-lon(i1)) * rdlon(i0)
2348 if ( a1<0.0 .or. a1>1.0 .or. b1<0.0 .or. b1>1.0 )
then 2349 write(*,*)
'gid=', mpp_pe(), i,j,a1, b1
2352 c1 = (1.-a1) * (1.-b1)
2357 sst_ncep(i,j) = c1*wk(i1,jc ) + c2*wk(i2,jc ) + &
2358 c3*wk(i2,jc+1) + c4*wk(i1,jc+1)
2366 subroutine remap_coef( is, ie, js, je, isd, ied, jsd, jed, &
2367 im, jm, lon, lat, id1, id2, jdc, s2c, agrid )
2369 integer,
intent(in):: is, ie, js, je, isd, ied, jsd, jed
2370 integer,
intent(in):: im, jm
2371 real,
intent(in):: lon(im), lat(jm)
2372 real,
intent(out):: s2c(is:ie,js:je,4)
2373 integer,
intent(out),
dimension(is:ie,js:je):: id1, id2, jdc
2374 real,
intent(in):: agrid(isd:ied,jsd:jed,2)
2379 integer i,j, i1, i2, jc, i0, j0
2382 rdlon(i) = 1. / (lon(i+1) - lon(i))
2384 rdlon(im) = 1. / (lon(1) + 2.*pi - lon(im))
2387 rdlat(j) = 1. / (lat(j+1) - lat(j))
2395 if ( agrid(i,j,1)>lon(im) )
then 2397 a1 = (agrid(i,j,1)-lon(im)) * rdlon(im)
2398 elseif ( agrid(i,j,1)<lon(1) )
then 2400 a1 = (agrid(i,j,1)+2.*pi-lon(im)) * rdlon(im)
2403 if ( agrid(i,j,1)>=lon(i0) .and. agrid(i,j,1)<=lon(i0+1) )
then 2405 a1 = (agrid(i,j,1)-lon(i1)) * rdlon(i0)
2412 if ( agrid(i,j,2)<lat(1) )
then 2415 elseif ( agrid(i,j,2)>lat(jm) )
then 2420 if ( agrid(i,j,2)>=lat(j0) .and. agrid(i,j,2)<=lat(j0+1) )
then 2422 b1 = (agrid(i,j,2)-lat(jc)) * rdlat(jc)
2429 if ( a1<0.0 .or. a1>1.0 .or. b1<0.0 .or. b1>1.0 )
then 2430 write(*,*)
'gid=', mpp_pe(), i,j,a1, b1
2433 s2c(i,j,1) = (1.-a1) * (1.-b1)
2434 s2c(i,j,2) = a1 * (1.-b1)
2435 s2c(i,j,3) = a1 * b1
2436 s2c(i,j,4) = (1.-a1) * b1
2446 subroutine remap_scalar(Atm, km, npz, ncnst, ak0, bk0, psc, qa, zh, omga, t_in)
2448 integer,
intent(in):: km, npz, ncnst
2449 real,
intent(in):: ak0(km+1), bk0(km+1)
2450 real,
intent(in),
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je):: psc
2451 real,
intent(in),
optional,
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je,km):: omga, t_in
2452 real,
intent(in),
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je,km,ncnst):: qa
2453 real,
intent(in),
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je,km+1):: zh
2455 real,
dimension(Atm%bd%is:Atm%bd%ie,km+1):: pe0
2456 real,
dimension(Atm%bd%is:Atm%bd%ie,npz):: qn1, dp2
2457 real,
dimension(Atm%bd%is:Atm%bd%ie,npz+1):: pe1
2458 real qp(atm%bd%is:atm%bd%ie,km)
2459 real wk(atm%bd%is:atm%bd%ie,atm%bd%js:atm%bd%je)
2460 real,
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je):: z500
2462 real(kind=R_GRID),
dimension(Atm%bd%is:Atm%bd%ie,npz+1):: pn1
2463 real(kind=R_GRID):: gz_fv(npz+1)
2464 real(kind=R_GRID),
dimension(2*km+1):: gz, pn
2465 real(kind=R_GRID),
dimension(Atm%bd%is:Atm%bd%ie,km+1):: pn0
2466 real(kind=R_GRID):: pst
2468 integer i,j,k,l,m, k2,iq
2469 integer sphum, liq_wat, ice_wat, rainwat, snowwat, graupel, cld_amt, liq_aero, ice_aero
2471 integer spfo, spfo2, spfo3
2475 integer :: is, ie, js, je
2482 sphum = get_tracer_index(model_atmos,
'sphum')
2483 liq_wat = get_tracer_index(model_atmos,
'liq_wat')
2484 ice_wat = get_tracer_index(model_atmos,
'ice_wat')
2485 rainwat = get_tracer_index(model_atmos,
'rainwat')
2486 snowwat = get_tracer_index(model_atmos,
'snowwat')
2487 graupel = get_tracer_index(model_atmos,
'graupel')
2488 cld_amt = get_tracer_index(model_atmos,
'cld_amt')
2490 spfo = get_tracer_index(model_atmos,
'spfo')
2491 spfo2 = get_tracer_index(model_atmos,
'spfo2')
2492 spfo3 = get_tracer_index(model_atmos,
'spfo3')
2494 o3mr = get_tracer_index(model_atmos,
'o3mr')
2496 liq_aero = get_tracer_index(model_atmos,
'liq_aero')
2497 ice_aero = get_tracer_index(model_atmos,
'ice_aero')
2501 if (mpp_pe()==1)
then 2502 print *,
'In remap_scalar:' 2503 print *,
'ncnst = ', ncnst
2504 print *,
'nwat = ', atm%flagstruct%nwat
2505 print *,
'sphum = ', sphum
2506 print *,
'clwmr = ', liq_wat
2508 print *,
'spfo3 = ', spfo3
2509 print *,
' spfo = ', spfo
2510 print *,
'spfo2 = ', spfo2
2512 print *,
' o3mr = ', o3mr
2514 print *,
'liq_aero = ', liq_aero
2515 print *,
'ice_aero = ', ice_aero
2516 if ( atm%flagstruct%nwat .eq. 6 )
then 2517 print *,
'rainwat = ', rainwat
2518 print *,
'ice_wat = ', ice_wat
2519 print *,
'snowwat = ', snowwat
2520 print *,
'graupel = ', graupel
2524 if ( sphum/=1 )
then 2525 call mpp_error(fatal,
'SPHUM must be 1st tracer')
2529 atm%phis(is:ie,js:je) = zh(is:ie,js:je,km+1)*grav
2532 if (atm%flagstruct%ecmwf_ic)
then 2533 if (cld_amt .gt. 0) atm%q(i,j,k,cld_amt) = 0.
2543 pe0(i,k) = ak0(k) + bk0(k)*psc(i,j)
2544 pn0(i,k) = log(pe0(i,k))
2551 gz(k) = zh(i,j,k)*grav
2557 gz(k) = 2.*gz(km+1) - gz(l)
2558 pn(k) = 2.*pn(km+1) - pn(l)
2562 if( atm%phis(i,j).le.gz(k) .and. atm%phis(i,j).ge.gz(k+1) )
then 2563 pst = pn(k) + (pn(k+1)-pn(k))*(gz(k)-atm%phis(i,j))/(gz(k)-gz(k+1))
2567 123 atm%ps(i,j) = exp(pst)
2574 if( pst.le.pn(k+1) .and. pst.ge.pn(k) )
then 2575 z500(i,j) = (gz(k+1) + (gz(k)-gz(k+1))*(pn(k+1)-pst)/(pn(k+1)-pn(k)))/grav
2584 pe1(i,1) = atm%ak(1)
2585 pn1(i,1) = log(pe1(i,1))
2589 pe1(i,k) = atm%ak(k) + atm%bk(k)*atm%ps(i,j)
2590 pn1(i,k) = log(pe1(i,k))
2597 dp2(i,k) = pe1(i,k+1) - pe1(i,k)
2598 atm%delp(i,j,k) = dp2(i,k)
2604 if (floor(qa(is,j,1,iq)) > -999)
then 2607 qp(i,k) = qa(i,j,k,iq)
2610 call mappm(km, pe0, qp, npz, pe1, qn1, is,ie, 0, 8, atm%ptop)
2611 if ( iq==sphum )
then 2612 call fillq(ie-is+1, npz, 1, qn1, dp2)
2614 call fillz(ie-is+1, npz, 1, qn1, dp2)
2619 atm%q(i,j,k,iq) = qn1(i,k)
2630 if ( pn1(i,1) .lt. pn0(i,1) )
then 2631 call mpp_error(fatal,
'FV3 top higher than external data')
2636 gz(k) = zh(i,j,k)*grav
2641 gz(k) = 2.*gz(km+1) - gz(l)
2642 pn(k) = 2.*pn(km+1) - pn(l)
2646 gz_fv(npz+1) = atm%phis(i,j)
2652 #ifdef USE_ISOTHERMO 2654 if ( (pn1(i,k).le.pn(l+1)) .and. (pn1(i,k).ge.pn(l)) )
then 2655 gz_fv(k) = gz(l) + (gz(l+1)-gz(l))*(pn1(i,k)-pn(l))/(pn(l+1)-pn(l))
2657 elseif ( pn1(i,k) .gt. pn(km+1) )
then 2659 gz_fv(k) = gz(km+1) + (gz_fv(npz+1)-gz(km+1))*(pn1(i,k)-pn(km+1))/(pn1(i,npz+1)-pn(km+1))
2665 if ( (pn1(i,k).le.pn(l+1)) .and. (pn1(i,k).ge.pn(l)) )
then 2666 gz_fv(k) = gz(l) + (gz(l+1)-gz(l))*(pn1(i,k)-pn(l))/(pn(l+1)-pn(l))
2675 atm%peln(i,k,j) = pn1(i,k)
2684 atm%pt(i,j,k) = (gz_fv(k)-gz_fv(k+1))/( rdgas*(pn1(i,k+1)-pn1(i,k))*
virq(atm%q(i,j,k,:)) )
2686 atm%pt(i,j,k) = (gz_fv(k)-gz_fv(k+1))/( rdgas*(pn1(i,k+1)-pn1(i,k))*(1.+
zvir*atm%q(i,j,k,sphum)) )
2694 qp(i,k) = t_in(i,j,k)
2697 call mappm(km, log(pe0), qp, npz, log(pe1), qn1, is,ie, 2, 4, atm%ptop)
2700 atm%pt(i,j,k) = qn1(i,k)
2704 if ( .not. atm%flagstruct%hydrostatic )
then 2706 atm%delz(i,j,k) = (gz_fv(k+1) - gz_fv(k)) / grav
2717 if ((atm%flagstruct%nwat .eq. 3 .or. atm%flagstruct%nwat .eq. 6) .and. &
2718 (atm%flagstruct%ncep_ic .or. atm%flagstruct%nggps_ic))
then 2722 qn1(i,k) = atm%q(i,j,k,liq_wat)
2723 if (cld_amt .gt. 0) atm%q(i,j,k,cld_amt) = 0.
2725 if ( atm%pt(i,j,k) > 273.16 )
then 2726 atm%q(i,j,k,liq_wat) = qn1(i,k)
2727 atm%q(i,j,k,ice_wat) = 0.
2728 #ifdef ORIG_CLOUDS_PART 2729 else if ( atm%pt(i,j,k) < 258.16 )
then 2730 atm%q(i,j,k,liq_wat) = 0.
2731 atm%q(i,j,k,ice_wat) = qn1(i,k)
2733 atm%q(i,j,k,liq_wat) = qn1(i,k)*((atm%pt(i,j,k)-258.16)/15.)
2734 atm%q(i,j,k,ice_wat) = qn1(i,k) - atm%q(i,j,k,liq_wat)
2737 else if ( atm%pt(i,j,k) < 233.16 )
then 2738 atm%q(i,j,k,liq_wat) = 0.
2739 atm%q(i,j,k,ice_wat) = qn1(i,k)
2742 atm%q(i,j,k,liq_wat) = qn1(i,k)*((atm%pt(i,j,k)-233.16)/40.)
2743 atm%q(i,j,k,ice_wat) = qn1(i,k) - atm%q(i,j,k,liq_wat)
2745 if (atm%pt(i,j,k)<258.16 .and. atm%q(i,j,k-1,ice_wat)>1.e-5 )
then 2746 atm%q(i,j,k,liq_wat) = 0.
2747 atm%q(i,j,k,ice_wat) = qn1(i,k)
2749 atm%q(i,j,k,liq_wat) = qn1(i,k)*((atm%pt(i,j,k)-233.16)/40.)
2750 atm%q(i,j,k,ice_wat) = qn1(i,k) - atm%q(i,j,k,liq_wat)
2755 if (atm%flagstruct%nwat .eq. 6 )
then 2756 atm%q(i,j,k,rainwat) = 0.
2757 atm%q(i,j,k,snowwat) = 0.
2758 atm%q(i,j,k,graupel) = 0.
2760 atm%q(i,j,k,ice_wat), atm%q(i,j,k,snowwat) )
2773 if ( (.not. atm%flagstruct%hydrostatic) .and. (.not. atm%flagstruct%ncep_ic) )
then 2776 qp(i,k) = omga(i,j,k)
2779 call mappm(km, pe0, qp, npz, pe1, qn1, is,ie, -1, 4, atm%ptop)
2783 atm%w(i,j,k) = qn1(i,k)
2789 atm%w(i,j,k) = qn1(i,k)/atm%delp(i,j,k)*atm%delz(i,j,k)
2798 if (.not. atm%flagstruct%hydrostatic)
call p_maxmin(
'delz_model', atm%delz, is, ie, js, je, npz, 1.)
2799 call p_maxmin(
'sphum_model', atm%q(is:ie,js:je,1:npz,sphum), is, ie, js, je, npz, 1.)
2800 call p_maxmin(
'liq_wat_model', atm%q(is:ie,js:je,1:npz,liq_wat), is, ie, js, je, npz, 1.)
2801 if (ice_wat .gt. 0)
call p_maxmin(
'ice_wat_model', atm%q(is:ie,js:je,1:npz,ice_wat), is, ie, js, je, npz, 1.)
2802 call p_maxmin(
'PS_model (mb)', atm%ps(is:ie,js:je), is, ie, js, je, 1, 0.01)
2803 call p_maxmin(
'PT_model', atm%pt(is:ie,js:je,1:npz), is, ie, js, je, npz, 1.)
2804 call pmaxmn(
'ZS_model', atm%phis(is:ie,js:je)/grav, is, ie, js, je, 1, 1., atm%gridstruct%area_64, atm%domain)
2805 call pmaxmn(
'ZS_data', zh(is:ie,js:je,km+1), is, ie, js, je, 1, 1., atm%gridstruct%area_64, atm%domain)
2808 wk(i,j) = atm%phis(i,j)/grav - zh(i,j,km+1)
2816 call pmaxmn(
'ZS_diff (m)', wk, is, ie, js, je, 1, 1., atm%gridstruct%area_64, atm%domain)
2818 if (.not.atm%gridstruct%bounded_domain)
then 2819 call prt_gb_nh_sh(
'DATA_IC Z500', is,ie, js,je, z500, atm%gridstruct%area_64(is:ie,js:je), atm%gridstruct%agrid_64(is:ie,js:je,2))
2820 if ( .not. atm%flagstruct%hydrostatic ) &
2821 call prt_height(
'fv3_IC Z500', is,ie, js,je, 3, npz, 500.e2, atm%phis, atm%delz, atm%peln, &
2822 atm%gridstruct%area_64(is:ie,js:je), atm%gridstruct%agrid_64(is:ie,js:je,2))
2827 wk(i,j) = atm%ps(i,j) - psc(i,j)
2830 call pmaxmn(
'PS_diff (mb)', wk, is, ie, js, je, 1, 0.01, atm%gridstruct%area_64, atm%domain)
2832 if (is_master())
write(*,*)
'done remap_scalar' 2838 integer,
intent(in):: km, npz, iq
2839 real,
intent(in):: ak0(km+1), bk0(km+1)
2840 real,
intent(in),
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je):: psc
2841 real,
intent(in),
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je,km):: qa
2842 real,
intent(in),
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je,km+1):: zh
2844 real,
dimension(Atm%bd%is:Atm%bd%ie,km+1):: pe0
2845 real,
dimension(Atm%bd%is:Atm%bd%ie,npz):: qn1, dp2
2846 real,
dimension(Atm%bd%is:Atm%bd%ie,npz+1):: pe1
2847 real qp(atm%bd%is:atm%bd%ie,km)
2848 real wk(atm%bd%is:atm%bd%ie,atm%bd%js:atm%bd%je)
2850 real(kind=R_GRID),
dimension(Atm%bd%is:Atm%bd%ie,npz+1):: pn1
2851 real(kind=R_GRID):: gz_fv(npz+1)
2852 real(kind=R_GRID),
dimension(2*km+1):: gz, pn
2853 real(kind=R_GRID),
dimension(Atm%bd%is:Atm%bd%ie,km+1):: pn0
2854 real(kind=R_GRID):: pst
2856 integer i,j,k, k2, l
2857 integer :: is, ie, js, je
2858 real,
allocatable:: ps_temp(:,:)
2867 allocate(ps_temp(is:ie,js:je))
2872 pe0(i,k) = ak0(k) + bk0(k)*psc(i,j)
2873 pn0(i,k) = log(pe0(i,k))
2880 gz(k) = zh(i,j,k)*grav
2886 gz(k) = 2.*gz(km+1) - gz(l)
2887 pn(k) = 2.*pn(km+1) - pn(l)
2891 if( atm%phis(i,j).le.gz(k) .and. atm%phis(i,j).ge.gz(k+1) )
then 2892 pst = pn(k) + (pn(k+1)-pn(k))*(gz(k)-atm%phis(i,j))/(gz(k)-gz(k+1))
2896 123 ps_temp(i,j) = exp(pst)
2900 pe1(i,1) = atm%ak(1)
2901 pn1(i,1) = log(pe1(i,1))
2905 pe1(i,k) = atm%ak(k) + atm%bk(k)*ps_temp(i,j)
2906 pn1(i,k) = log(pe1(i,k))
2913 dp2(i,k) = pe1(i,k+1) - pe1(i,k)
2923 call mappm(km, pe0, qp, npz, pe1, qn1, is,ie, 0, 8, atm%ptop)
2925 call fillq(ie-is+1, npz, 1, qn1, dp2)
2927 call fillz(ie-is+1, npz, 1, qn1, dp2)
2932 atm%q(i,j,k,iq) = qn1(i,k)
2937 call p_maxmin(
'o3mr remap', atm%q(is:ie,js:je,1:npz,iq), is, ie, js, je, npz, 1.)
2945 real,
intent(inout):: ql, qr, qi, qs
2946 real,
parameter:: qi0_max = 2.0e-3
2947 real,
parameter:: ql0_max = 2.5e-3
2950 if ( ql > ql0_max )
then 2955 if ( qi > qi0_max )
then 2963 subroutine remap_dwinds(km, npz, ak0, bk0, psc, ud, vd, Atm)
2965 integer,
intent(in):: km, npz
2966 real,
intent(in):: ak0(km+1), bk0(km+1)
2967 real,
intent(in):: psc(atm%bd%is:atm%bd%ie,atm%bd%js:atm%bd%je)
2968 real,
intent(in):: ud(atm%bd%is:atm%bd%ie,atm%bd%js:atm%bd%je+1,km)
2969 real,
intent(in):: vd(atm%bd%is:atm%bd%ie+1,atm%bd%js:atm%bd%je,km)
2971 real,
dimension(Atm%bd%isd:Atm%bd%ied,Atm%bd%jsd:Atm%bd%jed):: psd
2972 real,
dimension(Atm%bd%is:Atm%bd%ie+1, km+1):: pe0
2973 real,
dimension(Atm%bd%is:Atm%bd%ie+1,npz+1):: pe1
2974 real,
dimension(Atm%bd%is:Atm%bd%ie+1,npz):: qn1
2976 integer :: is, ie, js, je
2977 integer :: isd, ied, jsd, jed
2989 if (atm%gridstruct%bounded_domain)
then 2992 psd(i,j) = atm%ps(i,j)
3002 call mpp_update_domains( psd, atm%domain, complete=.false. )
3003 call mpp_update_domains( atm%ps, atm%domain, complete=.true. )
3013 pe0(i,k) = ak0(k) + bk0(k)*0.5*(psd(i,j-1)+psd(i,j))
3018 pe1(i,k) = atm%ak(k) + atm%bk(k)*0.5*(atm%ps(i,j-1)+atm%ps(i,j))
3021 call mappm(km, pe0(is:ie,1:km+1), ud(is:ie,j,1:km), npz, pe1(is:ie,1:npz+1), &
3022 qn1(is:ie,1:npz), is,ie, -1, 8, atm%ptop)
3025 atm%u(i,j,k) = qn1(i,k)
3031 if ( j/=(je+1) )
then 3035 pe0(i,k) = ak0(k) + bk0(k)*0.5*(psd(i-1,j)+psd(i,j))
3040 pe1(i,k) = atm%ak(k) + atm%bk(k)*0.5*(atm%ps(i-1,j)+atm%ps(i,j))
3043 call mappm(km, pe0(is:ie+1,1:km+1), vd(is:ie+1,j,1:km), npz, pe1(is:ie+1,1:npz+1), &
3044 qn1(is:ie+1,1:npz), is,ie+1, -1, 8, atm%ptop)
3047 atm%v(i,j,k) = qn1(i,k)
3055 if (is_master())
write(*,*)
'done remap_dwinds' 3060 subroutine remap_winds(im, jm, km, npz, ak0, bk0, psc, ua, va, Atm)
3062 integer,
intent(in):: im, jm, km, npz
3063 real,
intent(in):: ak0(km+1), bk0(km+1)
3064 real,
intent(in):: psc(atm%bd%is:atm%bd%ie,atm%bd%js:atm%bd%je)
3065 real,
intent(in),
dimension(Atm%bd%is:Atm%bd%ie,Atm%bd%js:Atm%bd%je,km):: ua, va
3067 real,
dimension(Atm%bd%isd:Atm%bd%ied,Atm%bd%jsd:Atm%bd%jed,npz):: ut, vt
3068 real,
dimension(Atm%bd%is:Atm%bd%ie, km+1):: pe0
3069 real,
dimension(Atm%bd%is:Atm%bd%ie,npz+1):: pe1
3070 real,
dimension(Atm%bd%is:Atm%bd%ie,npz):: qn1
3073 integer :: is, ie, js, je
3074 integer :: isd, ied, jsd, jed
3091 pe0(i,k) = ak0(k) + bk0(k)*psc(i,j)
3097 pe1(i,k) = atm%ak(k) + atm%bk(k)*atm%ps(i,j)
3104 call mappm(km, pe0, ua(is:ie,j,1:km), npz, pe1, qn1, is,ie, -1, 8, atm%ptop)
3107 ut(i,j,k) = qn1(i,k)
3113 call mappm(km, pe0, va(is:ie,j,1:km), npz, pe1, qn1, is,ie, -1, 8, atm%ptop)
3116 vt(i,j,k) = qn1(i,k)
3122 call prt_maxmin(
'UT', ut, is, ie, js, je, ng, npz, 1.)
3123 call prt_maxmin(
'VT', vt, is, ie, js, je, ng, npz, 1.)
3124 call prt_maxmin(
'UA_top',ut(:,:,1), is, ie, js, je, ng, 1, 1.)
3129 call cubed_a2d(atm%npx, atm%npy, npz, ut, vt, atm%u, atm%v, atm%gridstruct, atm%domain, atm%bd )
3131 if (is_master())
write(*,*)
'done remap_winds' 3136 subroutine remap_xyz( im, jbeg, jend, jm, km, npz, nq, ncnst, lon, lat, ak0, bk0, ps0, gz0, &
3137 ua, va, ta, qa, Atm )
3140 integer,
intent(in):: im, jm, km, npz, nq, ncnst
3141 integer,
intent(in):: jbeg, jend
3142 real,
intent(in):: lon(im), lat(jm), ak0(km+1), bk0(km+1)
3143 real,
intent(in):: gz0(im,jbeg:jend), ps0(im,jbeg:jend)
3144 real,
intent(in),
dimension(im,jbeg:jend,km):: ua, va, ta
3145 real,
intent(in),
dimension(im,jbeg:jend,km,ncnst):: qa
3147 real,
pointer,
dimension(:,:,:) :: agrid
3150 real,
dimension(Atm%bd%isd:Atm%bd%ied,Atm%bd%jsd:Atm%bd%jed,npz):: ut, vt
3151 real,
dimension(Atm%bd%is:Atm%bd%ie,km):: up, vp, tp
3152 real,
dimension(Atm%bd%is:Atm%bd%ie,km+1):: pe0, pn0
3153 real pt0(km), gz(km+1), pk0(km+1)
3154 real qp(atm%bd%is:atm%bd%ie,km,ncnst)
3155 real,
dimension(Atm%bd%is:Atm%bd%ie,npz):: qn1
3156 real,
dimension(Atm%bd%is:Atm%bd%ie,npz+1):: pe1, pn1
3159 real:: a1, b1, c1, c2, c3, c4
3160 real:: gzc, psc, pst
3164 integer i,j,k, i1, i2, jc, i0, j0, iq
3167 integer :: is, ie, js, je
3168 integer :: isd, ied, jsd, jed, ng
3181 agrid => atm%gridstruct%agrid
3183 sphum = get_tracer_index(model_atmos,
'sphum')
3185 if ( sphum/=1 )
then 3186 call mpp_error(fatal,
'SPHUM must be 1st tracer')
3189 pk0(1) = ak0(1)**kappa
3192 rdlon(i) = 1. / (lon(i+1) - lon(i))
3194 rdlon(im) = 1. / (lon(1) + 2.*pi - lon(im))
3197 rdlat(j) = 1. / (lat(j+1) - lat(j))
3205 pn0(i,1) = log(ak0(1))
3211 if ( agrid(i,j,1)>lon(im) )
then 3213 a1 = (agrid(i,j,1)-lon(im)) * rdlon(im)
3214 elseif ( agrid(i,j,1)<lon(1) )
then 3216 a1 = (agrid(i,j,1)+2.*pi-lon(im)) * rdlon(im)
3219 if ( agrid(i,j,1)>=lon(i0) .and. agrid(i,j,1)<=lon(i0+1) )
then 3221 a1 = (agrid(i,j,1)-lon(i1)) * rdlon(i0)
3229 if ( agrid(i,j,2)<lat(1) )
then 3232 elseif ( agrid(i,j,2)>lat(jm) )
then 3237 if ( agrid(i,j,2)>=lat(j0) .and. agrid(i,j,2)<=lat(j0+1) )
then 3239 b1 = (agrid(i,j,2)-lat(jc)) * rdlat(jc)
3247 if ( a1<0.0 .or. a1>1.0 .or. b1<0.0 .or. b1>1.0 )
then 3248 write(*,*) i,j,a1, b1
3251 c1 = (1.-a1) * (1.-b1)
3257 psc = c1*ps0(i1,jc ) + c2*ps0(i2,jc ) + &
3258 c3*ps0(i2,jc+1) + c4*ps0(i1,jc+1)
3261 gzc = c1*gz0(i1,jc ) + c2*gz0(i2,jc ) + &
3262 c3*gz0(i2,jc+1) + c4*gz0(i1,jc+1)
3268 qp(i,k,iq) = c1*qa(i1,jc, k,iq) + c2*qa(i2,jc, k,iq) + &
3269 c3*qa(i2,jc+1,k,iq) + c4*qa(i1,jc+1,k,iq)
3275 up(i,k) = c1*ua(i1,jc, k) + c2*ua(i2,jc, k) + &
3276 c3*ua(i2,jc+1,k) + c4*ua(i1,jc+1,k)
3277 vp(i,k) = c1*va(i1,jc, k) + c2*va(i2,jc, k) + &
3278 c3*va(i2,jc+1,k) + c4*va(i1,jc+1,k)
3279 tp(i,k) = c1*ta(i1,jc, k) + c2*ta(i2,jc, k) + &
3280 c3*ta(i2,jc+1,k) + c4*ta(i1,jc+1,k)
3283 tp(i,k) = tp(i,k)*
virq(qp(i,k,:))
3285 tp(i,k) = tp(i,k)*(1.+
zvir*qp(i,k,sphum))
3291 pe0(i,k) = ak0(k) + bk0(k)*psc
3292 pn0(i,k) = log(pe0(i,k))
3293 pk0(k) = pe0(i,k)**kappa
3304 gz(k) = gz(k+1) + rdgas*tp(i,k)*(pn0(i,k+1)-pn0(i,k))
3309 pkx = (pk0(km+1)-pk0(km))*(kappa*(pn0(i,km+1)-pn0(i,km)))
3310 pkx = exp( kappax*log(pkx) )
3311 pt0(km) = tp(i,km)/pkx
3313 pt0(km) = tp(i,km)/(pk0(km+1)-pk0(km))*(kappa*(pn0(i,km+1)-pn0(i,km)))
3315 if( atm%phis(i,j)>gzc )
then 3317 if( atm%phis(i,j) < gz(k) .and. &
3318 atm%phis(i,j) >= gz(k+1) )
then 3319 pst = pk0(k) + (pk0(k+1)-pk0(k))*(gz(k)-atm%phis(i,j))/(gz(k)-gz(k+1))
3326 pst = pk0(km+1) + (gzc-atm%phis(i,j))/(cp_air*pt0(km)*pkx)
3328 pst = pk0(km+1) + (gzc-atm%phis(i,j))/(cp_air*pt0(km))
3333 123 atm%ps(i,j) = pst**(1./(kappa*kappax))
3335 123 atm%ps(i,j) = pst**(1./kappa)
3343 pe1(i,1) = atm%ak(1)
3344 pn1(i,1) = log(pe1(i,1))
3348 pe1(i,k) = atm%ak(k) + atm%bk(k)*atm%ps(i,j)
3349 pn1(i,k) = log(pe1(i,k))
3355 atm%delp(i,j,k) = pe1(i,k+1) - pe1(i,k)
3363 call mappm(km, pe0, up, npz, pe1, qn1, is,ie, -1, 9, atm%ptop)
3366 ut(i,j,k) = qn1(i,k)
3372 call mappm(km, pe0, vp, npz, pe1, qn1, is,ie, -1, 9, atm%ptop)
3375 vt(i,j,k) = qn1(i,k)
3385 call mappm(km, pe0, qp(is,1,iq), npz, pe1, qn1, is,ie, 0, 11, atm%ptop)
3388 atm%q(i,j,k,iq) = qn1(i,k)
3397 call mappm(km, pn0, tp, npz, pn1, qn1, is,ie, 1, 9, atm%ptop)
3401 atm%pt(i,j,k) = qn1(i,k)/
virq(atm%q(i,j,k,:))
3403 atm%pt(i,j,k) = qn1(i,k)/(1.+
zvir*atm%q(i,j,k,sphum))
3410 call prt_maxmin(
'PS_model', atm%ps, is, ie, js, je, ng, 1, 0.01)
3411 call prt_maxmin(
'UT', ut, is, ie, js, je, ng, npz, 1.)
3412 call prt_maxmin(
'VT', vt, is, ie, js, je, ng, npz, 1.)
3417 call cubed_a2d(atm%npx, atm%npy, npz, ut, vt, atm%u, atm%v, atm%gridstruct, atm%domain, atm%bd )
3419 if (is_master())
write(*,*)
'done remap_xyz' 3424 subroutine cubed_a2d( npx, npy, npz, ua, va, u, v, gridstruct, fv_domain, bd )
3425 use mpp_domains_mod
, only: mpp_update_domains
3428 integer,
intent(in):: npx, npy, npz
3429 real,
intent(inout),
dimension(bd%isd:bd%ied,bd%jsd:bd%jed,npz):: ua, va
3430 real,
intent(out):: u(bd%isd:bd%ied, bd%jsd:bd%jed+1,npz)
3431 real,
intent(out):: v(bd%isd:bd%ied+1,bd%jsd:bd%jed ,npz)
3433 type(domain2d),
intent(INOUT) :: fv_domain
3435 real v3(3,bd%is-1:bd%ie+1,bd%js-1:bd%je+1)
3436 real ue(3,bd%is-1:bd%ie+1,bd%js:bd%je+1)
3437 real ve(3,bd%is:bd%ie+1,bd%js-1:bd%je+1)
3438 real,
dimension(bd%is:bd%ie):: ut1, ut2, ut3
3439 real,
dimension(bd%js:bd%je):: vt1, vt2, vt3
3440 integer i, j, k, im2, jm2
3442 real(kind=R_GRID),
pointer,
dimension(:,:,:) :: vlon, vlat
3443 real(kind=R_GRID),
pointer,
dimension(:) :: edge_vect_w, edge_vect_e, edge_vect_s, edge_vect_n
3444 real(kind=R_GRID),
pointer,
dimension(:,:,:,:) :: ew, es
3446 integer :: is, ie, js, je
3447 integer :: isd, ied, jsd, jed
3458 vlon => gridstruct%vlon
3459 vlat => gridstruct%vlat
3461 edge_vect_w => gridstruct%edge_vect_w
3462 edge_vect_e => gridstruct%edge_vect_e
3463 edge_vect_s => gridstruct%edge_vect_s
3464 edge_vect_n => gridstruct%edge_vect_n
3469 call mpp_update_domains(ua, fv_domain, complete=.false.)
3470 call mpp_update_domains(va, fv_domain, complete=.true.)
3479 v3(1,i,j) = ua(i,j,k)*vlon(i,j,1) + va(i,j,k)*vlat(i,j,1)
3480 v3(2,i,j) = ua(i,j,k)*vlon(i,j,2) + va(i,j,k)*vlat(i,j,2)
3481 v3(3,i,j) = ua(i,j,k)*vlon(i,j,3) + va(i,j,k)*vlat(i,j,3)
3489 ue(1,i,j) = 0.5*(v3(1,i,j-1) + v3(1,i,j))
3490 ue(2,i,j) = 0.5*(v3(2,i,j-1) + v3(2,i,j))
3491 ue(3,i,j) = 0.5*(v3(3,i,j-1) + v3(3,i,j))
3497 ve(1,i,j) = 0.5*(v3(1,i-1,j) + v3(1,i,j))
3498 ve(2,i,j) = 0.5*(v3(2,i-1,j) + v3(2,i,j))
3499 ve(3,i,j) = 0.5*(v3(3,i-1,j) + v3(3,i,j))
3504 if (.not. gridstruct%bounded_domain)
then 3509 vt1(j) = edge_vect_w(j)*ve(1,i,j-1)+(1.-edge_vect_w(j))*ve(1,i,j)
3510 vt2(j) = edge_vect_w(j)*ve(2,i,j-1)+(1.-edge_vect_w(j))*ve(2,i,j)
3511 vt3(j) = edge_vect_w(j)*ve(3,i,j-1)+(1.-edge_vect_w(j))*ve(3,i,j)
3513 vt1(j) = edge_vect_w(j)*ve(1,i,j+1)+(1.-edge_vect_w(j))*ve(1,i,j)
3514 vt2(j) = edge_vect_w(j)*ve(2,i,j+1)+(1.-edge_vect_w(j))*ve(2,i,j)
3515 vt3(j) = edge_vect_w(j)*ve(3,i,j+1)+(1.-edge_vect_w(j))*ve(3,i,j)
3525 if ( (ie+1)==npx )
then 3529 vt1(j) = edge_vect_e(j)*ve(1,i,j-1)+(1.-edge_vect_e(j))*ve(1,i,j)
3530 vt2(j) = edge_vect_e(j)*ve(2,i,j-1)+(1.-edge_vect_e(j))*ve(2,i,j)
3531 vt3(j) = edge_vect_e(j)*ve(3,i,j-1)+(1.-edge_vect_e(j))*ve(3,i,j)
3533 vt1(j) = edge_vect_e(j)*ve(1,i,j+1)+(1.-edge_vect_e(j))*ve(1,i,j)
3534 vt2(j) = edge_vect_e(j)*ve(2,i,j+1)+(1.-edge_vect_e(j))*ve(2,i,j)
3535 vt3(j) = edge_vect_e(j)*ve(3,i,j+1)+(1.-edge_vect_e(j))*ve(3,i,j)
3550 ut1(i) = edge_vect_s(i)*ue(1,i-1,j)+(1.-edge_vect_s(i))*ue(1,i,j)
3551 ut2(i) = edge_vect_s(i)*ue(2,i-1,j)+(1.-edge_vect_s(i))*ue(2,i,j)
3552 ut3(i) = edge_vect_s(i)*ue(3,i-1,j)+(1.-edge_vect_s(i))*ue(3,i,j)
3554 ut1(i) = edge_vect_s(i)*ue(1,i+1,j)+(1.-edge_vect_s(i))*ue(1,i,j)
3555 ut2(i) = edge_vect_s(i)*ue(2,i+1,j)+(1.-edge_vect_s(i))*ue(2,i,j)
3556 ut3(i) = edge_vect_s(i)*ue(3,i+1,j)+(1.-edge_vect_s(i))*ue(3,i,j)
3566 if ( (je+1)==npy )
then 3570 ut1(i) = edge_vect_n(i)*ue(1,i-1,j)+(1.-edge_vect_n(i))*ue(1,i,j)
3571 ut2(i) = edge_vect_n(i)*ue(2,i-1,j)+(1.-edge_vect_n(i))*ue(2,i,j)
3572 ut3(i) = edge_vect_n(i)*ue(3,i-1,j)+(1.-edge_vect_n(i))*ue(3,i,j)
3574 ut1(i) = edge_vect_n(i)*ue(1,i+1,j)+(1.-edge_vect_n(i))*ue(1,i,j)
3575 ut2(i) = edge_vect_n(i)*ue(2,i+1,j)+(1.-edge_vect_n(i))*ue(2,i,j)
3576 ut3(i) = edge_vect_n(i)*ue(3,i+1,j)+(1.-edge_vect_n(i))*ue(3,i,j)
3590 u(i,j,k) = ue(1,i,j)*es(1,i,j,1) + &
3591 ue(2,i,j)*es(2,i,j,1) + &
3592 ue(3,i,j)*es(3,i,j,1)
3597 v(i,j,k) = ve(1,i,j)*ew(1,i,j,2) + &
3598 ve(2,i,j)*ew(2,i,j,2) + &
3599 ve(3,i,j)*ew(3,i,j,2)
3608 subroutine d2a3d(u, v, ua, va, im, jm, km, lon)
3609 integer,
intent(in):: im, jm, km
3610 real,
intent(in ) :: lon(im)
3611 real,
intent(in ),
dimension(im,jm,km):: u, v
3612 real,
intent(out),
dimension(im,jm,km):: ua, va
3614 real :: coslon(im),sinlon(im)
3624 sinlon(i) = sin(lon(i))
3625 coslon(i) = cos(lon(i))
3631 ua(i,j,k) = 0.5*(u(i,j,k) + u(i,j+1,k))
3637 va(i,j,k) = 0.5*(v(i,j,k) + v(i+1,j,k))
3639 va(im,j,k) = 0.5*(v(im,j,k) + v(1,j,k))
3646 us = us + (ua(i+imh,2,k)-ua(i,2,k))*sinlon(i) &
3647 + (va(i,2,k)-va(i+imh,2,k))*coslon(i)
3648 vs = vs + (ua(i+imh,2,k)-ua(i,2,k))*coslon(i) &
3649 + (va(i+imh,2,k)-va(i,2,k))*sinlon(i)
3654 ua(i,1,k) = -us*sinlon(i) - vs*coslon(i)
3655 va(i,1,k) = us*coslon(i) - vs*sinlon(i)
3656 ua(i+imh,1,k) = -ua(i,1,k)
3657 va(i+imh,1,k) = -va(i,1,k)
3664 un = un + (ua(i+imh,jm-1,k)-ua(i,jm-1,k))*sinlon(i) &
3665 + (va(i+imh,jm-1,k)-va(i,jm-1,k))*coslon(i)
3666 vn = vn + (ua(i,jm-1,k)-ua(i+imh,jm-1,k))*coslon(i) &
3667 + (va(i+imh,jm-1,k)-va(i,jm-1,k))*sinlon(i)
3673 ua(i,jm,k) = -un*sinlon(i) + vn*coslon(i)
3674 va(i,jm,k) = -un*coslon(i) - vn*sinlon(i)
3675 ua(i+imh,jm,k) = -ua(i,jm,k)
3676 va(i+imh,jm,k) = -va(i,jm,k)
3680 end subroutine d2a3d 3683 subroutine pmaxmin( qname, a, im, jm, fac )
3685 integer,
intent(in):: im, jm
3686 character(len=*) :: qname
3690 real qmin(jm), qmax(jm)
3698 pmax = max(pmax, a(i,j))
3699 pmin = min(pmin, a(i,j))
3710 pmax = max(pmax, qmax(j))
3711 pmin = min(pmin, qmin(j))
3714 write(*,*) qname,
' max = ', pmax*fac,
' min = ', pmin*fac
3718 subroutine pmaxmn(qname, q, is, ie, js, je, km, fac, area, domain)
3719 character(len=*),
intent(in):: qname
3720 integer,
intent(in):: is, ie, js, je
3721 integer,
intent(in):: km
3722 real,
intent(in):: q(is:ie, js:je, km)
3723 real,
intent(in):: fac
3724 real(kind=R_GRID),
intent(IN):: area(is-3:ie+3, js-3:je+3)
3725 type(domain2d),
intent(INOUT) :: domain
3727 real qmin, qmax, gmean
3737 if( q(i,j,k) < qmin )
then 3739 elseif( q(i,j,k) > qmax )
then 3746 call mp_reduce_min(qmin)
3747 call mp_reduce_max(qmax)
3749 gmean =
g_sum(domain, q(is,js,km), is, ie, js, je, 3, area, 1, reproduce=.true.)
3750 if(is_master())
write(6,*) qname, qmax*fac, qmin*fac, gmean*fac
3754 subroutine p_maxmin(qname, q, is, ie, js, je, km, fac)
3755 character(len=*),
intent(in):: qname
3756 integer,
intent(in):: is, ie, js, je, km
3757 real,
intent(in):: q(is:ie, js:je, km)
3758 real,
intent(in):: fac
3767 if( q(i,j,k) < qmin )
then 3769 elseif( q(i,j,k) > qmax )
then 3775 call mp_reduce_min(qmin)
3776 call mp_reduce_max(qmax)
3777 if(is_master())
write(6,*) qname, qmax*fac, qmin*fac
3781 subroutine fillq(im, km, nq, q, dp)
3782 integer,
intent(in):: im
3783 integer,
intent(in):: km
3784 integer,
intent(in):: nq
3785 real ,
intent(in):: dp(im,km)
3786 real ,
intent(inout) :: q(im,km,nq)
3788 integer i, k, ic, k1
3795 if( q(i,k,ic) < 0. )
then 3796 q(i,k1,ic) = q(i,k1,ic) + q(i,k,ic)*dp(i,k)/dp(i,k1)
3805 if( q(i,k,ic) < 0. )
then 3806 q(i,k1,ic) = q(i,k1,ic) + q(i,k,ic)*dp(i,k)/dp(i,k1)
3814 end subroutine fillq 3816 subroutine compute_zh(im, jm, levp, ak0, bk0, ps, zs, t, q, nq, zh )
3818 integer,
intent(in):: levp, im,jm, nq
3819 real,
intent(in),
dimension(levp+1):: ak0, bk0
3820 real(kind=4),
intent(in),
dimension(im,jm):: ps, zs
3821 real(kind=4),
intent(in),
dimension(im,jm,levp):: t
3822 real(kind=4),
intent(in),
dimension(im,jm,levp,nq):: q
3823 real(kind=4),
intent(out),
dimension(im,jm,levp+1):: zh
3825 real,
dimension(im,levp+1):: pe0, pn0
3835 pn0(i,1) = log(pe0(i,1))
3836 zh(i,j,levp+1) = zs(i,j)
3841 pe0(i,k) = ak0(k) + bk0(k)*ps(i,j)
3842 pn0(i,k) = log(pe0(i,k))
3849 zh(i,j,k) = zh(i,j,k+1)+(t(i,j,k)*(1.+
zvir*q(i,j,k,1))*(pn0(i,k+1)-pn0(i,k)))*(rdgas/grav)
3858 subroutine get_staggered_grid( is, ie, js, je, isd, ied, jsd, jed, pt_b, pt_c, pt_d)
3859 integer,
intent(in):: is, ie, js, je, isd, ied, jsd, jed
3860 real,
dimension(isd:ied+1,jsd:jed+1,2),
intent(in) :: pt_b
3861 real,
dimension(isd:ied+1,jsd:jed ,2),
intent(out) :: pt_c
3862 real,
dimension(isd:ied ,jsd:jed+1,2),
intent(out) :: pt_d
3864 real(kind=R_GRID),
dimension(2):: p1, p2, p3
3869 p1(:) = pt_b(i, j,1:2)
3870 p2(:) = pt_b(i+1,j,1:2)
3872 pt_d(i,j,1:2) = p3(:)
3878 p1(:) = pt_b(i,j ,1:2)
3879 p2(:) = pt_b(i,j+1,1:2)
3881 pt_c(i,j,1:2) = p3(:)
subroutine, public mid_pt_sphere(p1, p2, pm)
subroutine, public prt_height(qname, is, ie, js, je, ng, km, press, phis, delz, peln, area, lat)
pure real function, public vicpqd(q)
subroutine, public get_var3_r4(ncid, var3_name, is, ie, js, je, ks, ke, var3, time_slice)
The module 'fv_mp_mod' is a single program multiple data (SPMD) parallel decompostion/communication m...
subroutine timing_off(blk_name)
The subroutine 'timing_off' stops a timer.
The type 'fv_grid_type' is made up of grid-dependent information from fv_grid_tools and fv_grid_utils...
character(len=80), public source
integer, parameter, public h_stagger
character(len=27), parameter source_fv3gfs
subroutine, public mappm(km, pe1, q1, kn, pe2, q2, i1, i2, iv, kord, ptop)
The subroutine 'mappm' is a general-purpose routine for remapping one set of vertical levels to anoth...
subroutine, public start_regional_cold_start(Atm, dt_atmos, ak, bk, levp, is, ie, js, je, isd, ied, jsd, jed)
The module 'multi_gases' peforms multi constitutents computations.
subroutine compute_zh(im, jm, levp, ak0, bk0, ps, zs, t, q, nq, zh)
subroutine, public del2_cubed_sphere(npx, npy, q, area, dx, dy, dxc, dyc, sin_sg, nmax, cd, zero_ocean, oro, bounded_domain, domain, bd)
subroutine, public p_var(km, ifirst, ilast, jfirst, jlast, ptop, ptop_min, delp, delz, pt, ps, pe, peln, pk, pkz, cappa, q, ng, nq, area, dry_mass, adjust_dry_mass, mountain, moist_phys, hydrostatic, nwat, domain, adiabatic, make_nh)
the subroutine 'p_var' computes auxiliary pressure variables for a hydrostatic state.
real function, public g_sum(domain, p, ifirst, ilast, jfirst, jlast, ngc, area, mode, reproduce)
The function 'g_sum' is the fast version of 'globalsum'.
real function, public inner_prod(v1, v2)
The module 'fv_io' contains restart facilities for FV core.
subroutine get_staggered_grid(is, ie, js, je, isd, ied, jsd, jed, pt_b, pt_c, pt_d)
subroutine, public get_latlon_vector(pp, elon, elat)
subroutine, public get_unit_vect2(e1, e2, uc)
subroutine remap_xyz(im, jbeg, jend, jm, km, npz, nq, ncnst, lon, lat, ak0, bk0, ps0, gz0, ua, va, ta, qa, Atm)
integer, parameter, public r_grid
pure real function, public virqd(q)
The module fv_nwp_nudge contains routines for nudging to input analyses. note This module is currentl...
real(kind=r_grid), parameter cnst_0p20
subroutine, public get_var1_real(ncid, var1_name, im, var1, var_exist)
subroutine, public surfdrv(npx, npy, grid, agrid, area, dx, dy, dxa, dya, dxc, dyc, sin_sg, phis, stretch_fac, nested, bounded_domain, npx_global, domain, grid_number, bd)
The module 'fv_timing' contains FV3 timers.
subroutine, public remap_dwinds(km, npz, ak0, bk0, psc, ud, vd, Atm)
pure real function, public virq(q)
subroutine, public fv3_zs_filter(bd, isd, ied, jsd, jed, npx, npy, npx_global, stretch_fac, bounded_domain, domain, area, dxa, dya, dx, dy, dxc, dyc, grid, agrid, sin_sg, phis, oro)
The module 'boundary' contains utility routines for grid nesting and boundary conditions.
subroutine pmaxmin(qname, a, im, jm, fac)
subroutine, public get_var1_double(ncid, var1_name, im, var1, var_exist)
The 'get_var' subroutines read in variables from netcdf files.
real, parameter, public ptop_min
subroutine, public fv_io_read_tracers(fv_domain, Atm)
The subroutine 'fv_io_read_tracers' reads in only tracers from restart files.
The module 'fv_mapz' contains the vertical mapping routines .
subroutine, public fillz(im, km, nq, q, dp)
The subroutine 'fillz' is for mass-conservative filling of nonphysical negative values in the tracers...
The module 'fv_arrays' contains the 'fv_atmos_type' and associated datatypes.
interface 'nested_grid_BC' includes subroutines 'nested_grid_BC_2d' and 'nested_grid_BC_3d' that fetc...
The module 'sim_nc' is a netcdf file reader.
subroutine ncep2fms(im, jm, lon, lat, wk)
subroutine, public set_external_eta(ak, bk, ptop, ks)
The subroutine 'set_external_eta' sets 'ptop' (model top) and 'ks' (first level of pure pressure coor...
subroutine, public get_var2_real(ncid, var_name, im, jm, var2)
subroutine pmaxmn(qname, q, is, ie, js, je, km, fac, area, domain)
subroutine, public set_eta(km, ks, ptop, ak, bk, npz_type)
The module 'fv_eta' contains routine to set up the reference (Eulerian) pressure coordinate.
subroutine cubed_a2d(npx, npy, npz, ua, va, u, v, gridstruct, fv_domain, bd)
The subroutine 'cubed_a2d' transforms the wind from the A Grid to the D Grid.
subroutine, public get_var2_r4(ncid, var2_name, is, ie, js, je, var2, time_slice)
The module 'external_ic_mod' contains routines that read in and remap initial conditions.
subroutine d2a3d(u, v, ua, va, im, jm, km, lon)
real, dimension(:,:), allocatable, public sgh_g
subroutine, public open_ncfile(iflnm, ncid)
subroutine remap_coef(is, ie, js, je, isd, ied, jsd, jed, im, jm, lon, lat, id1, id2, jdc, s2c, agrid)
subroutine, public get_data_source(source, regional)
subroutine, public get_var_att_double(ncid, var_name, att_name, value)
subroutine timing_on(blk_name)
The subroutine 'timing_on' starts a timer.
subroutine get_fv_ic(Atm, fv_domain, nq)
subroutine get_ncep_ic(Atm, fv_domain, nq)
The subroutine 'get_ncep_ic' reads in the specified NCEP analysis or reanalysis dataset.
subroutine remap_winds(im, jm, km, npz, ak0, bk0, psc, ua, va, Atm)
integer, parameter, public v_stagger
@ The module 'fv_diagnostics' contains routines to compute diagnosic fields.
The module 'fv_grid_utils' contains routines for setting up and computing grid-related quantities...
subroutine, public prt_gb_nh_sh(qname, is, ie, js, je, a2, area, lat)
real, dimension(:,:), allocatable, public oro_g
subroutine, public get_external_ic(Atm, fv_domain, cold_start, dt_atmos)
subroutine, public prt_maxmin(qname, q, is, ie, js, je, n_g, km, fac)
subroutine, public get_cubed_sphere_terrain(Atm, fv_domain)
subroutine, public checker_tracers(i0, i1, j0, j1, ifirst, ilast, jfirst, jlast, nq, km, q, lon, lat, nx, ny, rn)
subroutine mp_auto_conversion(ql, qr, qi, qs)
subroutine remap_scalar_single(Atm, km, npz, ak0, bk0, psc, qa, zh, iq)
subroutine, public remap_scalar(Atm, km, npz, ncnst, ak0, bk0, psc, qa, zh, omga, t_in)
subroutine get_ecmwf_ic(Atm, fv_domain)
The subroutine 'get_ecmwf_ic' reads in initial conditions from ECMWF analyses (EXPERIMENTAL: contact ...
subroutine, public get_ncdim1(ncid, var1_name, im)
subroutine, public del4_cubed_sphere(npx, npy, q, area, dx, dy, dxc, dyc, sin_sg, nmax, zero_ocean, oro, bounded_domain, domain, bd)
subroutine get_nggps_ic(Atm, fv_domain, dt_atmos)
The subroutine 'get_nggps_ic' reads in data after it has been preprocessed with NCEP/EMC orography ma...
subroutine p_maxmin(qname, q, is, ie, js, je, km, fac)
subroutine, public extrapolation_bc(q, istag, jstag, npx, npy, bd, pd_in, debug_in)
The subroutine 'extrapolation_BC' performs linear extrapolation into the halo region.
subroutine fillq(im, km, nq, q, dp)
integer, parameter, public u_stagger
subroutine, public close_ncfile(ncid)