CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
gwdc_post.f
1
4
5 module gwdc_post
6
7 contains
8
12 subroutine gwdc_post_run( &
13 & im, levs, lssav, ldiag3d, dtf, dtp, con_cp, &
14 & tauctx, taucty, gwdcu, gwdcv, &
15 & dugwd, dvgwd, dtend, dtidx, index_of_x_wind, index_of_y_wind, &
16 & index_of_process_nonorographic_gwd, gu0, gv0, gt0, &
17 & errmsg, errflg)
18
19 use machine, only : kind_phys
20 implicit none
21
22 integer, intent(in) :: im, levs
23 logical, intent(in) :: lssav, ldiag3d
24 real(kind=kind_phys), intent(in) :: dtf, dtp, con_cp
25 real(kind=kind_phys), intent(in) :: &
26 & tauctx(:), taucty(:), gwdcu(:,:), gwdcv(:,:)
27
28 real(kind=kind_phys), intent(inout) :: dugwd(:), dvgwd(:), &
29 & gu0(:,:), gv0(:,:), gt0(:,:)
30 real(kind=kind_phys), intent(inout), optional :: dtend(:,:,:)
31 integer, intent(in) :: dtidx(:,:)
32 integer, intent(in) :: index_of_process_nonorographic_gwd
33 integer, intent(in) :: index_of_x_wind, index_of_y_wind
34
35 character(len=*), intent(out) :: errmsg
36 integer, intent(out) :: errflg
37
38 integer :: i, k, idtend
39 real(kind=kind_phys) :: eng0, eng1
40
41 ! Initialize CCPP error handling variables
42 errmsg = ''
43 errflg = 0
44
45! --- ... write out cloud top stress and wind tendencies
46
47 if (lssav) then
48 dugwd(:) = dugwd(:) + tauctx(:)*dtf
49 dvgwd(:) = dvgwd(:) + taucty(:)*dtf
50 endif ! end if_lssav
51
52 if (ldiag3d) then
53 idtend = dtidx(index_of_x_wind,index_of_process_nonorographic_g&
54 & wd)
55 if(idtend>=1) then
56 dtend(:,:,idtend) = dtend(:,:,idtend) + gwdcu(:,:) * dtf
57 endif
58 idtend = dtidx(index_of_y_wind,index_of_process_nonorographic_g&
59 & wd)
60 if(idtend>=1) then
61 dtend(:,:,idtend) = dtend(:,:,idtend) + gwdcv(:,:) * dtf
62 endif
63 endif
64
65! --- ... update the wind components with gwdc tendencies
66
67 do k = 1, levs
68 do i = 1, im
69 eng0 = 0.5*(gu0(i,k)*gu0(i,k) + gv0(i,k)*gv0(i,k))
70 gu0(i,k) = gu0(i,k) + gwdcu(i,k) * dtp
71 gv0(i,k) = gv0(i,k) + gwdcv(i,k) * dtp
72 eng1 = 0.5*(gu0(i,k)*gu0(i,k) + gv0(i,k)*gv0(i,k))
73 gt0(i,k) = gt0(i,k) + (eng0-eng1)/(dtp*con_cp)
74 enddo
75! if (lprnt) write(7000,*)' gu0=',gu0(ipr,k),' gwdcu=',
76! &gwdcu(ipr,k), ' gv0=', gv0(ipr,k),' gwdcv=',gwdcv(ipr,k)
77! &,' k=',k
78 enddo
79
80 end subroutine gwdc_post_run
81
82 end module gwdc_post
Definition gwdc_post.f:5