CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
sfc_nst_pre.f90
1
3
5
6 use machine , only : kind_phys
8 use module_nst_parameters , only : zero, one
9
10 implicit none
11
12contains
13
24 subroutine sfc_nst_pre_run &
25 (im, wet, tgice, tsfco, tsurf_wat, &
26 tseal, xt, xz, dt_cool, z_c, tref, cplflx, &
27 oceanfrac, nthreads, errmsg, errflg)
28
29 ! --- inputs:
30 integer, intent(in) :: im, nthreads
31 logical, dimension(:), intent(in) :: wet
32 real (kind=kind_phys), intent(in) :: tgice
33 real (kind=kind_phys), dimension(:), intent(in) :: tsfco, oceanfrac
34 real (kind=kind_phys), dimension(:), intent(in), optional :: xt, xz, dt_cool, z_c
35 logical, intent(in) :: cplflx
36
37 ! --- input/outputs:
38 real (kind=kind_phys), dimension(:), intent(inout) :: tsurf_wat, tseal
39 real (kind=kind_phys), dimension(:), intent(inout), optional :: tref
40
41 ! --- outputs:
42 character(len=*), intent(out) :: errmsg
43 integer, intent(out) :: errflg
44
45 ! --- locals
46 integer :: i
47 real(kind=kind_phys), parameter :: omz1 = 2.0_kind_phys
48 real(kind=kind_phys) :: tem2, dnsst
49 real(kind=kind_phys), dimension(im) :: dtzm, z_c_0
50
51 ! Initialize CCPP error handling variables
52 errmsg = ''
53 errflg = 0
54
55 do i=1,im
56 if (wet(i) .and. oceanfrac(i) > 0.0) then
57 ! tem = (oro(i)-oro_uf(i)) * rlapse
58 ! DH* 20190927 simplyfing this code because tem is zero
59 !tem = zero
60 !tseal(i) = tsfco(i) + tem
61 tseal(i) = tsfco(i)
62 !tsurf_wat(i) = tsurf_wat(i) + tem
63 ! *DH
64 endif
65 enddo
66 !
67 ! update tsfc & tref with T1 from OGCM & NSST Profile if coupled
68 !
69 if (cplflx) then
70 z_c_0 = zero
71 call get_dtzm_2d (xt, xz, dt_cool, z_c_0, wet, zero, omz1, im, 1, nthreads, dtzm)
72 do i=1,im
73 if (wet(i) .and. oceanfrac(i) > zero ) then
74 ! dnsst = tsfc_wat(i) - tref(i) ! retrive/get difference of Ts and Tf
75 tref(i) = max(tgice, tsfco(i) - dtzm(i)) ! update Tf with T1 and NSST T-Profile
76 ! tsfc_wat(i) = max(271.2,tref(i) + dnsst) ! get Ts updated due to Tf update
77 ! tseal(i) = tsfc_wat(i)
78 if (abs(xz(i)) > zero) then
79 tem2 = one / xz(i)
80 else
81 tem2 = zero
82 endif
83 tseal(i) = tref(i) + (xt(i)+xt(i)) * tem2 - dt_cool(i)
84 tsurf_wat(i) = tseal(i)
85 endif
86 enddo
87 endif
88
89 return
90 end subroutine sfc_nst_pre_run
91end module sfc_nst_pre
subroutine, public get_dtzm_2d(xt, xz, dt_cool, zc, wet, z1, z2, nx, ny, nth, dtm)