diff --git a/atmos_model.F90 b/atmos_model.F90 index 29133471a..3ccc1a3ab 100644 --- a/atmos_model.F90 +++ b/atmos_model.F90 @@ -1572,7 +1572,7 @@ subroutine assign_importdata(jdat, rc) type(ESMF_Grid) :: grid type(ESMF_Field) :: dbgField character(19) :: currtimestring - real (kind=GFS_kind_phys), parameter :: z0ice=1.1 ! (in cm) + real (kind=GFS_kind_phys), parameter :: z0ice=1.0 ! (in cm) ! ! real(kind=GFS_kind_phys), parameter :: himax = 8.0 !< maximum ice thickness allowed diff --git a/ccpp/physics b/ccpp/physics index 9880f4407..fc331b814 160000 --- a/ccpp/physics +++ b/ccpp/physics @@ -1 +1 @@ -Subproject commit 9880f4407c9132a61dd11741c6a114ae295d8fba +Subproject commit fc331b8141294bec4190cdc39a23cf093be063c7 diff --git a/cpl/module_block_data.F90 b/cpl/module_block_data.F90 index 02ef0ebc8..9d2cc9192 100644 --- a/cpl/module_block_data.F90 +++ b/cpl/module_block_data.F90 @@ -71,6 +71,7 @@ subroutine block_copy_1d_i4_to_2d_r8(destin_ptr, source_ptr, block, block_index, if (associated(destin_ptr) .and. associated(source_ptr)) then factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -105,6 +106,7 @@ subroutine block_copy_1d_to_2d_r8(destin_ptr, source_ptr, block, block_index, sc if (associated(destin_ptr) .and. associated(source_ptr)) then factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -144,6 +146,7 @@ subroutine block_copy_1dslice_to_2d_r8(destin_ptr, source_ptr, slice, block, blo if (slice > 0 .and. slice <= size(source_ptr, dim=2)) then factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -182,6 +185,7 @@ subroutine block_copy_2d_to_3d_r8(destin_ptr, source_ptr, block, block_index, sc factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor do k = 1, size(source_ptr, dim=2) +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -219,6 +223,7 @@ subroutine block_copy_2d_to_2d_r8(destin_ptr, source_ptr, block, block_index, sc if (associated(destin_ptr) .and. associated(source_ptr)) then factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -253,6 +258,7 @@ subroutine block_array_copy_2d_to_2d_r8(destin_ptr, source_arr, block, block_ind if (associated(destin_ptr)) then factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -290,6 +296,7 @@ subroutine block_copy_3d_to_3d_r8(destin_ptr, source_ptr, block, block_index, sc factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor do k = 1, size(source_ptr, dim=3) +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -326,6 +333,7 @@ subroutine block_array_copy_3d_to_3d_r8(destin_ptr, source_arr, block, block_ind factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor do k = 1, size(source_arr, dim=3) +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -367,6 +375,7 @@ subroutine block_copy_3dslice_to_3d_r8(destin_ptr, source_ptr, slice, block, blo factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor do k = 1, size(source_ptr, dim=3) +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -407,6 +416,7 @@ subroutine block_array_copy_3dslice_to_3d_r8(destin_ptr, source_arr, slice, bloc factor = 1._kind_phys if (present(scale_factor)) factor = scale_factor do k = 1, size(source_arr, dim=3) +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -441,6 +451,7 @@ subroutine block_fill_2d_r8(destin_ptr, fill_value, block, block_index, rc) ! -- begin localrc = ESMF_RC_PTR_NOTALLOC if (associated(destin_ptr)) then +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -474,6 +485,7 @@ subroutine block_fill_3d_r8(destin_ptr, fill_value, block, block_index, rc) localrc = ESMF_RC_PTR_NOTALLOC if (associated(destin_ptr)) then do k = 1, size(destin_ptr, dim=3) +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) @@ -586,6 +598,7 @@ subroutine block_combine_frac_1d_to_2d_r8(destin_ptr, fract1_ptr, fract2_ptr, bl localrc = ESMF_RC_PTR_NOTALLOC if (associated(destin_ptr) .and. & associated(fract1_ptr) .and. associated(fract2_ptr)) then +!$omp parallel do private(ix,ib,jb,i,j) do ix = 1, block%blksz(block_index) ib = block%index(block_index)%ii(ix) jb = block%index(block_index)%jj(ix) diff --git a/io/FV3GFS_io.F90 b/io/FV3GFS_io.F90 index 3827ccb68..b45190f3d 100644 --- a/io/FV3GFS_io.F90 +++ b/io/FV3GFS_io.F90 @@ -1462,7 +1462,7 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain, warm_sta do nb = 1, Atm_block%nblks do ix = 1, Atm_block%blksz(nb) if (Sfcprop(nb)%landfrac(ix) > zero) then - tem = one / Sfcprop(nb)%landfrac(ix) + tem = one / (Sfcprop(nb)%fice(ix)*(one-Sfcprop(nb)%landfrac(ix))+Sfcprop(nb)%landfrac(ix)) Sfcprop(nb)%snodl(ix) = Sfcprop(nb)%snowd(ix) * tem else Sfcprop(nb)%snodl(ix) = zero @@ -1477,7 +1477,7 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain, warm_sta do nb = 1, Atm_block%nblks do ix = 1, Atm_block%blksz(nb) if (Sfcprop(nb)%landfrac(ix) > zero) then - tem = one / Sfcprop(nb)%landfrac(ix) + tem = one / (Sfcprop(nb)%fice(ix)*(one-Sfcprop(nb)%landfrac(ix))+Sfcprop(nb)%landfrac(ix)) Sfcprop(nb)%weasdl(ix) = Sfcprop(nb)%weasd(ix) * tem else Sfcprop(nb)%weasdl(ix) = zero @@ -1501,7 +1501,9 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain, warm_sta !$omp parallel do default(shared) private(nb, ix) do nb = 1, Atm_block%nblks do ix = 1, Atm_block%blksz(nb) - Sfcprop(nb)%zorlw(ix) = Sfcprop(nb)%zorl(ix) !--- compute zorlw from existing variables + if (Sfcprop(nb)%landfrac(ix) < one .and. Sfcprop(nb)%fice(ix) < one) then + Sfcprop(nb)%zorlw(ix) = min(Sfcprop(nb)%zorl(ix), 0.317) + endif enddo enddo endif @@ -1521,7 +1523,9 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain, warm_sta !$omp parallel do default(shared) private(nb, ix) do nb = 1, Atm_block%nblks do ix = 1, Atm_block%blksz(nb) - Sfcprop(nb)%zorli(ix) = Sfcprop(nb)%zorl(ix) !--- compute zorli from existing variables + if (Sfcprop(nb)%fice(ix)*(one-Sfcprop(nb)%landfrac(ix)) > zero) then + Sfcprop(nb)%zorli(ix) = one + endif enddo enddo endif @@ -1547,6 +1551,36 @@ subroutine sfc_prop_restart_read (Sfcprop, Atm_block, Model, fv_domain, warm_sta enddo endif + if (sfc_var2(i,j,47) < -9990.0_r8) then + if (Model%me == Model%master ) call mpp_error(NOTE, 'gfs_driver::surface_props_input - computing snodi') +!$omp parallel do default(shared) private(nb, ix, tem) + do nb = 1, Atm_block%nblks + do ix = 1, Atm_block%blksz(nb) + if (Sfcprop(nb)%fice(ix) > zero) then + tem = one / (Sfcprop(nb)%fice(ix)*(one-Sfcprop(nb)%landfrac(ix))+Sfcprop(nb)%landfrac(ix)) + Sfcprop(nb)%snodi(ix) = min(Sfcprop(nb)%snowd(ix) * tem, 3.0) + else + Sfcprop(nb)%snodi(ix) = zero + endif + enddo + enddo + endif + + if (sfc_var2(i,j,48) < -9990.0_r8) then + if (Model%me == Model%master ) call mpp_error(NOTE, 'gfs_driver::surface_props_input - computing weasdi') +!$omp parallel do default(shared) private(nb, ix, tem) + do nb = 1, Atm_block%nblks + do ix = 1, Atm_block%blksz(nb) + if (Sfcprop(nb)%fice(ix) > zero) then + tem = one / (Sfcprop(nb)%fice(ix)*(one-Sfcprop(nb)%landfrac(ix))+Sfcprop(nb)%landfrac(ix)) + Sfcprop(nb)%weasdi(ix) = Sfcprop(nb)%weasd(ix)*tem + else + Sfcprop(nb)%weasdi(ix) = zero + endif + enddo + enddo + endif + if (Model%use_cice_alb) then if (sfc_var2(i,j,49) < -9990.0_r8) then !$omp parallel do default(shared) private(nb, ix) @@ -3058,7 +3092,7 @@ subroutine fv_phys_bundle_setup(Diag, axes, phys_bundle, fcst_grid, quilting, nb ! implicit none ! - type(GFS_externaldiag_type),intent(in) :: Diag(:) + type(GFS_externaldiag_type),intent(in) :: Diag(:) integer, intent(in) :: axes(:) type(ESMF_FieldBundle),intent(inout) :: phys_bundle(:) type(ESMF_Grid),intent(inout) :: fcst_grid @@ -3099,7 +3133,7 @@ subroutine fv_phys_bundle_setup(Diag, axes, phys_bundle, fcst_grid, quilting, nb !------------------------------------------------------------ ! allocate(bdl_intplmethod(nbdlphys), outputfile(nbdlphys)) - if(mpp_pe()==mpp_root_pe())print *,'in fv_phys bundle,nbdl=',nbdlphys + if(mpp_pe()==mpp_root_pe()) print *,'in fv_phys bundle,nbdl=',nbdlphys do ibdl = 1, nbdlphys loutputfile = .false. call ESMF_FieldBundleGet(phys_bundle(ibdl), name=physbdl_name,rc=rc) @@ -3178,14 +3212,14 @@ subroutine fv_phys_bundle_setup(Diag, axes, phys_bundle, fcst_grid, quilting, nb allocate(udimList(udimCount)) call ESMF_AttributeGet(fcst_grid, convention="NetCDF", purpose="FV3", & name="vertical_dim_labels", valueList=udimList, rc=rc) -! if(mpp_pe()==mpp_root_pe())print *,'in fv3gfsio, vertical +! if(mpp_pe()==mpp_root_pe()) print *,'in fv3gfsio, vertical ! list=',udimList(1:udimCount),'rc=',rc if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return else - if(mpp_pe()==mpp_root_pe())print *,'in fv_dyn bundle,axis_name_vert=',axis_name_vert + if(mpp_pe()==mpp_root_pe()) print *,'in fv_dyn bundle,axis_name_vert=',axis_name_vert call ESMF_AttributeAdd(fcst_grid, convention="NetCDF", purpose="FV3", & attrList=(/"vertical_dim_labels"/), rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, line=__LINE__, file=__FILE__)) return @@ -3207,13 +3241,13 @@ subroutine fv_phys_bundle_setup(Diag, axes, phys_bundle, fcst_grid, quilting, nb direction, edges, Domain, DomainU, axis_data, & num_attributes=num_attributes, attributes=attributes) ! - edgesS='' + edgesS = '' do i = 1,num_axes_phys if(axes(i) == edges) edgesS=axis_name(i) enddo ! Add vertical dimension Attributes to Grid if( id>2 ) then -! if(mpp_pe()==mpp_root_pe())print *,' in dyn add grid, axis_name=', & +! if(mpp_pe()==mpp_root_pe()) print *,' in dyn add grid, axis_name=', & ! trim(axis_name(id)),'axis_data=',axis_data if(trim(edgesS)/='') then call ESMF_AttributeAdd(fcst_grid, convention="NetCDF", purpose="FV3", & @@ -3415,62 +3449,62 @@ subroutine add_field_to_phybundle(var_name,long_name,units,cell_methods, axes,ph ! !*** add field attributes call ESMF_AttributeAdd(field, convention="NetCDF", purpose="FV3", & - attrList=(/"long_name"/), rc=rc) + attrList=(/"long_name"/), rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeSet(field, convention="NetCDF", purpose="FV3", & - name='long_name',value=trim(long_name),rc=rc) + name='long_name',value=trim(long_name),rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeAdd(field, convention="NetCDF", purpose="FV3", & - attrList=(/"units"/), rc=rc) + attrList=(/"units"/), rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeSet(field, convention="NetCDF", purpose="FV3", & - name='units',value=trim(units),rc=rc) + name='units',value=trim(units),rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeAdd(field, convention="NetCDF", purpose="FV3", & - attrList=(/"missing_value"/), rc=rc) + attrList=(/"missing_value"/), rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeSet(field, convention="NetCDF", purpose="FV3", & - name='missing_value',value=missing_value,rc=rc) + name='missing_value',value=missing_value,rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeAdd(field, convention="NetCDF", purpose="FV3", & - attrList=(/"_FillValue"/), rc=rc) + attrList=(/"_FillValue"/), rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeSet(field, convention="NetCDF", purpose="FV3", & - name='_FillValue',value=missing_value,rc=rc) + name='_FillValue',value=missing_value,rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeAdd(field, convention="NetCDF", purpose="FV3", & - attrList=(/"cell_methods"/), rc=rc) + attrList=(/"cell_methods"/), rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeSet(field, convention="NetCDF", purpose="FV3", & - name='cell_methods',value=trim(cell_methods),rc=rc) + name='cell_methods',value=trim(cell_methods),rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) ! call ESMF_AttributeAdd(field, convention="NetCDF", purpose="FV3", & - attrList=(/"output_file"/), rc=rc) + attrList=(/"output_file"/), rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT) call ESMF_AttributeSet(field, convention="NetCDF", purpose="FV3", & - name='output_file',value=trim(output_file),rc=rc) + name='output_file',value=trim(output_file),rc=rc) if (ESMF_LogFoundError(rcToCheck=rc, msg=ESMF_LOGERR_PASSTHRU, & line=__LINE__, file=__FILE__)) call ESMF_Finalize(endflag=ESMF_END_ABORT)