Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

User/wfc/climate nudging #47

Merged
merged 5 commits into from
Oct 9, 2020
Merged

User/wfc/climate nudging #47

Show file tree
Hide file tree
Changes from 4 commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
c6d411d
Add PV diagnostics for 350K and 550K surfaces
wfcooke May 22, 2020
25b3383
Added namelist options for ECDA model.
wfcooke May 22, 2020
cce65ac
Update to correct some OpenMP assignments.
wfcooke Aug 12, 2020
db7b60f
Corrected declaration of variable in OMP region
wfcooke Aug 27, 2020
4dacd7e
Adding missing value for underground points in get_height_given_pressure
wfcooke Oct 6, 2020
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 3 additions & 0 deletions model/fv_arrays.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -60,6 +60,9 @@ module fv_arrays_mod
id_dbz, id_maxdbz, id_basedbz, id_dbz4km, id_dbztop, id_dbz_m10C, &
id_ctz, id_w1km, id_wmaxup, id_wmaxdn, id_cape, id_cin

! Selected theta-level fields from 3D variables:
integer :: id_pv350K, id_pv550K

! Selected p-level fields from 3D variables:
integer :: id_vort200, id_vort500, id_w500, id_w700
integer :: id_vort850, id_w850, id_x850, id_srh25, &
Expand Down
58 changes: 45 additions & 13 deletions tools/fv_diagnostics.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -708,6 +708,10 @@ subroutine fv_diag_init(Atm, axes, Time, npx, npy, npz, p_ref)
!--------------------
idiag%id_pv = register_diag_field ( trim(field), 'pv', axes(1:3), Time, &
'potential vorticity', '1/s', missing_value=missing_value )
idiag%id_pv350K = register_diag_field ( trim(field), 'pv350K', axes(1:2), Time, &
'350-K potential vorticity; needs x350 scaling', '(K m**2) / (kg s)', missing_value=missing_value)
idiag%id_pv550K = register_diag_field ( trim(field), 'pv550K', axes(1:2), Time, &
'550-K potential vorticity; needs x550 scaling', '(K m**2) / (kg s)', missing_value=missing_value)

! -------------------
! Vertical flux correlation terms (good for averages)
Expand Down Expand Up @@ -1577,8 +1581,8 @@ subroutine fv_diag(Atm, zvir, Time, print_freq)
endif

if ( idiag%id_vort200>0 .or. idiag%id_vort500>0 .or. idiag%id_vort850>0 .or. idiag%id_vorts>0 &
.or. idiag%id_vort>0 .or. idiag%id_pv>0 .or. idiag%id_rh>0 .or. idiag%id_x850>0 .or. &
idiag%id_uh03>0 .or. idiag%id_uh25>0) then
.or. idiag%id_vort>0 .or. idiag%id_pv>0 .or. idiag%id_pv350K>0 .or. idiag%id_pv550K>0 &
.or. idiag%id_rh>0 .or. idiag%id_x850>0 .or. idiag%id_uh03>0 .or. idiag%id_uh25>0) then
call get_vorticity(isc, iec, jsc, jec, isd, ied, jsd, jed, npz, Atm(n)%u, Atm(n)%v, wk, &
Atm(n)%gridstruct%dx, Atm(n)%gridstruct%dy, Atm(n)%gridstruct%rarea)

Expand Down Expand Up @@ -1732,11 +1736,36 @@ subroutine fv_diag(Atm, zvir, Time, print_freq)
endif


if ( idiag%id_pv > 0 ) then
! Note: this is expensive computation.
if ( idiag%id_pv > 0 .or. idiag%id_pv350K >0 .or. idiag%id_pv550K >0) then
if (allocated(a3)) deallocate(a3)
allocate ( a3(isc:iec,jsc:jec,npz+1) )
! Modified pv_entropy to get potential temperature at layer interfaces (last variable)
! The values are needed for interpolate_z
! Note: this is expensive computation.
call pv_entropy(isc, iec, jsc, jec, ngc, npz, wk, &
Atm(n)%gridstruct%f0, Atm(n)%pt, Atm(n)%pkz, Atm(n)%delp, grav)
used = send_data ( idiag%id_pv, wk, Time )
Atm(n)%gridstruct%f0, Atm(n)%pt, Atm(n)%pkz, Atm(n)%delp, grav,a3)
if ( idiag%id_pv > 0) then
used = send_data ( idiag%id_pv, wk, Time )
endif
if( idiag%id_pv350K > 0 .or. idiag%id_pv550K >0 ) then
!"pot temp" from pv_entropy is only semi-finished; needs p0^kappa (pk0)
do k=1,npz+1
do j=jsc,jec
do i=isc,iec
a3(i,j,k) = a3(i,j,k) * pk0
enddo
enddo
enddo
!wk as input, which stores pv from pv_entropy;
!use z interpolation, both potential temp and z increase monotonically with height
!interpolate_vertical will apply log operation to 350, and also assumes a different vertical order
call interpolate_z(isc, iec, jsc, jec, npz, 350., a3, wk, a2)
used = send_data( idiag%id_pv350K, a2, Time)
!interpolate_vertical will apply log operation to 350, and also assumes a different vertical order
call interpolate_z(isc, iec, jsc, jec, npz, 550., a3, wk, a2)
used = send_data( idiag%id_pv550K, a2, Time)
endif
deallocate ( a3 )
if (prt_minmax) call prt_maxmin('PV', wk, isc, iec, jsc, jec, 0, 1, 1.)
endif

Expand Down Expand Up @@ -4920,7 +4949,7 @@ end subroutine updraft_helicity



subroutine pv_entropy(is, ie, js, je, ng, km, vort, f_d, pt, pkz, delp, grav)
subroutine pv_entropy(is, ie, js, je, ng, km, vort, f_d, pt, pkz, delp, grav, te)

! !INPUT PARAMETERS:
integer, intent(in):: is, ie, js, je, ng, km
Expand All @@ -4931,7 +4960,9 @@ subroutine pv_entropy(is, ie, js, je, ng, km, vort, f_d, pt, pkz, delp, grav)
real, intent(in):: f_d(is-ng:ie+ng,js-ng:je+ng)

! vort is relative vorticity as input. Becomes PV on output
real, intent(inout):: vort(is:ie,js:je,km)
real, intent(inout):: vort(is:ie,js:je,km)
! output potential temperature at the interface so it can be used for diagnostics
real, intent(out):: te(is:ie,js:je,km+1)

! !DESCRIPTION:
! EPV = 1/r * (vort+f_d) * d(S)/dz; where S is a conservative scalar
Expand All @@ -4942,7 +4973,7 @@ subroutine pv_entropy(is, ie, js, je, ng, km, vort, f_d, pt, pkz, delp, grav)
! z-surface is not that different from the hybrid sigma-p coordinate.
! See page 39, Pedlosky 1979: Geophysical Fluid Dynamics
!
! The follwoing simplified form is strictly correct only if vort is computed on
! The following simplified form is strictly correct only if vort is computed on
! constant z surfaces. In addition hydrostatic approximation is made.
! EPV = - GRAV * (vort+f_d) / del(p) * del(pt) / pt
! where del() is the vertical difference operator.
Expand All @@ -4953,7 +4984,7 @@ subroutine pv_entropy(is, ie, js, je, ng, km, vort, f_d, pt, pkz, delp, grav)
!---------------------------------------------------------------------
!BOC
real w3d(is:ie,js:je,km)
real te(is:ie,js:je,km+1), t2(is:ie,km), delp2(is:ie,km)
real t2(is:ie,km), delp2(is:ie,km)
real te2(is:ie,km+1)
integer i, j, k

Expand All @@ -4966,8 +4997,8 @@ subroutine pv_entropy(is, ie, js, je, ng, km, vort, f_d, pt, pkz, delp, grav)
enddo
#else
! Compute PT at layer edges.
!$OMP parallel do default(none) shared(is,ie,js,je,km,pt,pkz,w3d,delp,te2,te) &
!$OMP private(t2, delp2)
!$OMP parallel do default(none) shared(is,ie,js,je,km,pt,pkz,w3d,delp,te) &
!$OMP private(t2, te2, delp2)
do j=js,je
do k=1,km
do i=is,ie
Expand All @@ -4986,7 +5017,8 @@ subroutine pv_entropy(is, ie, js, je, ng, km, vort, f_d, pt, pkz, delp, grav)
enddo
enddo

!$OMP parallel do default(none) shared(is,ie,js,je,km,vort,f_d,te,w3d,delp,grav)
!$OMP parallel do default(none) shared(is,ie,js,je,km,vort,te,w3d,delp,grav) &
!$OMP private(f_d)
do k=1,km
do j=js,je
do i=is,ie
Expand Down
Loading