CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
get_phi_fv3.F90
1
3
5
6 use machine, only: kind_phys
7 use physcons, only: con_fvirt
8
9!--- public declarations
10 public get_phi_fv3_run
11
12!--- local variables
13 real(kind=kind_phys), parameter :: zero = 0.0_kind_phys
14 real(kind=kind_phys), parameter :: half = 0.5_kind_phys
15 real(kind=kind_phys), parameter :: one = 1.0_kind_phys
16
17contains
18
19!! \section arg_table_get_phi_fv3_run Argument Table
20!! \htmlinclude get_phi_fv3_run.html
21!!
22 subroutine get_phi_fv3_run(ix, levs, con_fvirt, gt0, gq01, del_gz, phii, phil, errmsg, errflg)
23
24 implicit none
25
26 ! Interface variables
27 integer, intent(in) :: ix, levs
28 real(kind=kind_phys), intent(in) :: con_fvirt
29 real(kind=kind_phys), dimension(:,:), intent(in) :: gt0
30 real(kind=kind_phys), dimension(:,:), intent(in) :: gq01
31 real(kind=kind_phys), dimension(:,:), intent(inout) :: del_gz
32 real(kind=kind_phys), dimension(:,:), intent(out) :: phii
33 real(kind=kind_phys), dimension(:,:), intent(out) :: phil
34 character(len=*), intent(out) :: errmsg
35 integer, intent(out) :: errflg
36
37 ! Local variables
38 integer :: i, k
39
40 ! Initialize CCPP error handling variables
41 errmsg = ''
42 errflg = 0
43
44! SJL: Adjust the height hydrostatically in a way consistent with FV3 discretization
45 do i=1,ix
46 phii(i,1) = zero
47 enddo
48 do k=1,levs
49 do i=1,ix
50 del_gz(i,k) = del_gz(i,k)*gt0(i,k) * &
51 & (one + con_fvirt*max(zero,gq01(i,k)))
52 phii(i,k+1) = phii(i,k) + del_gz(i,k)
53 phil(i,k) = half*(phii(i,k) + phii(i,k+1))
54 enddo
55 enddo
56
57 end subroutine get_phi_fv3_run
58
59end module get_phi_fv3
This module contains some of the most frequently used math and physics constants for GCM models.
Definition physcons.F90:36