source: codes/icosagcm/trunk/src/kernels/caldyn_vert_NH.k90 @ 580

Last change on this file since 580 was 580, checked in by dubos, 7 years ago

trunk : upgrading to devel
File size: 2.7 KB
Line 
1   !--------------------------------------------------------------------------
2   !---------------------------- caldyn_vert_NH ----------------------------------
3   IF (ll_begin==1) THEN
4      !DIR$ SIMD
5      DO ij=ij_begin_ext, ij_end_ext
6         w_ij = ( W(ij,1)+.5*W(ij,1 +1) )/mass(ij,1)
7         wflux_ij = .5*(wflux(ij,1)+wflux(ij,1 +1))
8         W_etadot(ij,1) = wflux_ij*w_ij
9         eta_dot(ij,1) = wflux_ij / mass(ij,1)
10         wcov(ij,1) = w_ij*(geopot(ij,1 +1)-geopot(ij,1))
11      END DO
12   END IF
13   DO l = ll_beginp1, ll_endm1
14      !DIR$ SIMD
15      DO ij=ij_begin_ext, ij_end_ext
16         w_ij = .5*( W(ij,l)+W(ij,l+1) )/mass(ij,l)
17         wflux_ij = .5*(wflux(ij,l)+wflux(ij,l+1))
18         W_etadot(ij,l) = wflux_ij*w_ij
19         eta_dot(ij,l) = wflux_ij / mass(ij,l)
20         wcov(ij,l) = w_ij*(geopot(ij,l+1)-geopot(ij,l))
21      END DO
22   END DO
23   IF(ll_end==llm) THEN
24      !DIR$ SIMD
25      DO ij=ij_begin_ext, ij_end_ext
26         w_ij = ( .5*W(ij,llm)+W(ij,llm+1) )/mass(ij,llm)
27         wflux_ij = .5*(wflux(ij,llm)+wflux(ij,llm+1))
28         W_etadot(ij,llm) = wflux_ij*w_ij
29         eta_dot(ij,llm) = wflux_ij / mass(ij,llm)
30         wcov(ij,llm) = w_ij*(geopot(ij,llm+1)-geopot(ij,llm))
31      END DO
32   END IF
33   ! add NH term to du
34   DO l = ll_begin, ll_end
35      !DIR$ SIMD
36      DO ij=ij_begin, ij_end
37         du(ij+u_right,l) = du(ij+u_right,l) - .5*(wcov(ij,l)+wcov(ij+t_right,l))*ne_right*(eta_dot(ij+t_right,l)-eta_dot(ij,l))
38         du(ij+u_lup,l) = du(ij+u_lup,l) - .5*(wcov(ij,l)+wcov(ij+t_lup,l))*ne_lup*(eta_dot(ij+t_lup,l)-eta_dot(ij,l))
39         du(ij+u_ldown,l) = du(ij+u_ldown,l) - .5*(wcov(ij,l)+wcov(ij+t_ldown,l))*ne_ldown*(eta_dot(ij+t_ldown,l)-eta_dot(ij,l))
40      END DO
41   END DO
42   ! add NH terms to dW, dPhi
43   ! FIXME : TODO top and bottom
44   DO l = ll_beginp1, ll_end
45      !DIR$ SIMD
46      DO ij=ij_begin, ij_end
47         dPhi(ij,l)=dPhi(ij,l)-wflux(ij,l)*(geopot(ij,l+1)-geopot(ij,l-1))/(mass(ij,l-1)+mass(ij,l))
48      END DO
49   END DO
50   ! We need a barrier here because we compute W_etadot above and do a vertical difference below
51   !$OMP BARRIER
52   IF (ll_begin==1) THEN
53      !DIR$ SIMD
54      DO ij=ij_begin, ij_end
55         dW(ij,1) = dW(ij,1) - W_etadot(ij,1)
56      END DO
57   END IF
58   DO l = ll_beginp1, ll_end
59      !DIR$ SIMD
60      DO ij=ij_begin, ij_end
61         dW(ij,l) = dW(ij,l) + W_etadot(ij,l-1) - W_etadot(ij,l)
62      END DO
63   END DO
64   IF(ll_endp1==llm+1) THEN
65      !DIR$ SIMD
66      DO ij=ij_begin, ij_end
67         dW(ij,llm+1) = dW(ij,llm+1) + W_etadot(ij,llm+1 -1)
68      END DO
69   END IF
70   !---------------------------- caldyn_vert_NH ----------------------------------
71   !--------------------------------------------------------------------------
Note: See TracBrowser for help on using the repository browser.