CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
cu_ntiedtke_pre.F90
1
3
5
6 implicit none
7
8 private
9
10 public :: cu_ntiedtke_pre_run
11
12 contains
13
17 subroutine cu_ntiedtke_pre_run (flag_init, flag_restart, kdt, fhour, dtp, t, q, prevst, prevsq, &
18 forcet, forceq, errmsg, errflg)
19
20 use machine, only: kind_phys
21
22 implicit none
23
24 logical, intent(in) :: flag_init
25 logical, intent(in) :: flag_restart
26 integer, intent(in) :: kdt
27 real(kind_phys), intent(in) :: fhour
28 real(kind_phys), intent(in) :: dtp
29 real(kind_phys), intent(in) :: t(:,:)
30 real(kind_phys), intent(in) :: q(:,:)
31 real(kind_phys), intent(in), optional :: prevst(:,:)
32 real(kind_phys), intent(in), optional :: prevsq(:,:)
33 real(kind_phys), intent(out), optional :: forcet(:,:)
34 real(kind_phys), intent(out), optional :: forceq(:,:)
35 character(len=*), intent(out) :: errmsg
36 integer, intent(out) :: errflg
37
38 ! local variables
39 real(kind=kind_phys) :: dtdyn
40
41 ! Initialize CCPP error handling variables
42 errmsg = ''
43 errflg = 0
44
45 ! For restart runs, can assume that prevst and prevsq
46 ! are read from the restart files beforehand, same
47 ! for conv_act.
48 if(flag_init .and. .not.flag_restart) then
49 forcet(:,:)=0.0
50 forceq(:,:)=0.0
51 else
52 dtdyn=3600.0*(fhour)/kdt
53 if(dtp > dtdyn) then
54 forcet(:,:)=(t(:,:) - prevst(:,:))/dtp
55 forceq(:,:)=(q(:,:) - prevsq(:,:))/dtp
56 else
57 forcet(:,:)=(t(:,:) - prevst(:,:))/dtdyn
58 forceq(:,:)=(q(:,:) - prevsq(:,:))/dtdyn
59 endif
60 endif
61
62 end subroutine cu_ntiedtke_pre_run
63
64end module cu_ntiedtke_pre