New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
dynzad.F90 in NEMO/branches/2020/dev_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final/src/OCE/DYN – NEMO

source: NEMO/branches/2020/dev_r14116_HPC-04_mcastril_Mixed_Precision_implementation_final/src/OCE/DYN/dynzad.F90 @ 14219

Last change on this file since 14219 was 14219, checked in by mcastril, 4 years ago

Add Mixed Precision support by Oriol Tintó

  • Property svn:keywords set to Id
File size: 5.9 KB
Line 
1MODULE dynzad
2   !!======================================================================
3   !!                       ***  MODULE  dynzad  ***
4   !! Ocean dynamics : vertical advection trend
5   !!======================================================================
6   !! History :  OPA  ! 1991-01  (G. Madec) Original code
7   !!   NEMO     0.5  ! 2002-07  (G. Madec) Free form, F90
8   !!----------------------------------------------------------------------
9
10   !!----------------------------------------------------------------------
11   !!   dyn_zad       : vertical advection momentum trend
12   !!----------------------------------------------------------------------
13   USE oce            ! ocean dynamics and tracers
14   USE dom_oce        ! ocean space and time domain
15   USE sbc_oce        ! surface boundary condition: ocean
16   USE trd_oce        ! trends: ocean variables
17   USE trddyn         ! trend manager: dynamics
18   USE sbcwave, ONLY: wsd   ! Surface Waves (add vertical Stokes-drift)
19   !
20   USE in_out_manager ! I/O manager
21   USE lib_mpp        ! MPP library
22   USE prtctl         ! Print control
23   USE timing         ! Timing
24
25   IMPLICIT NONE
26   PRIVATE
27
28   PUBLIC   dyn_zad       ! routine called by dynadv.F90
29
30   !! * Substitutions
31#  include "do_loop_substitute.h90"
32#  include "domzgr_substitute.h90"
33#  include "single_precision_substitute.h90"
34   !!----------------------------------------------------------------------
35   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
36   !! $Id$
37   !! Software governed by the CeCILL license (see ./LICENSE)
38   !!----------------------------------------------------------------------
39CONTAINS
40
41   SUBROUTINE dyn_zad ( kt, Kmm, puu, pvv, Krhs )
42      !!----------------------------------------------------------------------
43      !!                  ***  ROUTINE dynzad  ***
44      !!
45      !! ** Purpose :   Compute the now vertical momentum advection trend and
46      !!      add it to the general trend of momentum equation.
47      !!
48      !! ** Method  :   The now vertical advection of momentum is given by:
49      !!         w dz(u) = u(rhs) + 1/(e1e2u*e3u) mk+1[ mi(e1e2t*ww) dk(u) ]
50      !!         w dz(v) = v(rhs) + 1/(e1e2v*e3v) mk+1[ mj(e1e2t*ww) dk(v) ]
51      !!      Add this trend to the general trend (puu(:,:,:,Krhs),pvv(:,:,:,Krhs)):
52      !!         (u(rhs),v(rhs)) = (u(rhs),v(rhs)) + w dz(u,v)
53      !!
54      !! ** Action  : - Update (puu(:,:,:,Krhs),pvv(:,:,:,Krhs)) with the vert. momentum adv. trends
55      !!              - Send the trends to trddyn for diagnostics (l_trddyn=T)
56      !!----------------------------------------------------------------------
57      INTEGER                             , INTENT( in )  ::  kt               ! ocean time-step inedx
58      INTEGER                             , INTENT( in )  ::  Kmm, Krhs        ! ocean time level indices
59      REAL(dp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) ::  puu, pvv         ! ocean velocities and RHS of momentum equation
60      !
61      INTEGER  ::   ji, jj, jk   ! dummy loop indices
62      REAL(wp) ::   zua, zva     ! local scalars
63      REAL(wp), DIMENSION(jpi,jpj)     ::   zww
64      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zwuw, zwvw
65      REAL(dp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrdu, ztrdv
66      !!----------------------------------------------------------------------
67      !
68      IF( ln_timing )   CALL timing_start('dyn_zad')
69      !
70      IF( kt == nit000 ) THEN
71         IF(lwp) WRITE(numout,*)
72         IF(lwp) WRITE(numout,*) 'dyn_zad : 2nd order vertical advection scheme'
73      ENDIF
74
75      IF( l_trddyn )   THEN           ! Save puu(:,:,:,Krhs) and pvv(:,:,:,Krhs) trends
76         ALLOCATE( ztrdu(jpi,jpj,jpk) , ztrdv(jpi,jpj,jpk) ) 
77         ztrdu(:,:,:) = puu(:,:,:,Krhs) 
78         ztrdv(:,:,:) = pvv(:,:,:,Krhs) 
79      ENDIF
80
81      DO jk = 2, jpkm1                ! Vertical momentum advection at level w and u- and v- vertical
82         DO_2D( 0, 1, 0, 1 )              ! vertical fluxes
83          IF( ln_vortex_force ) THEN
84            zww(ji,jj) = 0.25_wp * e1e2t(ji,jj) * ( ww(ji,jj,jk) + wsd(ji,jj,jk) )
85          ELSE
86            zww(ji,jj) = 0.25_wp * e1e2t(ji,jj) * ww(ji,jj,jk)
87          ENDIF
88         END_2D
89         DO_2D( 0, 0, 0, 0 )              ! vertical momentum advection at w-point
90            zwuw(ji,jj,jk) = ( zww(ji+1,jj  ) + zww(ji,jj) ) * ( puu(ji,jj,jk-1,Kmm) - puu(ji,jj,jk,Kmm) )
91            zwvw(ji,jj,jk) = ( zww(ji  ,jj+1) + zww(ji,jj) ) * ( pvv(ji,jj,jk-1,Kmm) - pvv(ji,jj,jk,Kmm) )
92         END_2D
93      END DO
94      !
95      ! Surface and bottom advective fluxes set to zero
96      DO_2D( 0, 0, 0, 0 )
97         zwuw(ji,jj, 1 ) = 0._wp
98         zwvw(ji,jj, 1 ) = 0._wp
99         zwuw(ji,jj,jpk) = 0._wp
100         zwvw(ji,jj,jpk) = 0._wp
101      END_2D
102      !
103      DO_3D( 0, 0, 0, 0, 1, jpkm1 )   ! Vertical momentum advection at u- and v-points
104         puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zwuw(ji,jj,jk) + zwuw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj)   &
105            &                                      / e3u(ji,jj,jk,Kmm)
106         pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zwvw(ji,jj,jk) + zwvw(ji,jj,jk+1) ) * r1_e1e2v(ji,jj)   &
107            &                                      / e3v(ji,jj,jk,Kmm)
108      END_3D
109
110      IF( l_trddyn ) THEN             ! save the vertical advection trends for diagnostic
111         ztrdu(:,:,:) = puu(:,:,:,Krhs) - ztrdu(:,:,:)
112         ztrdv(:,:,:) = pvv(:,:,:,Krhs) - ztrdv(:,:,:)
113         CALL trd_dyn( ztrdu, ztrdv, jpdyn_zad, kt, Kmm )
114         DEALLOCATE( ztrdu, ztrdv )
115      ENDIF
116      !                               ! Control print
117      IF(sn_cfctl%l_prtctl)   CALL prt_ctl( tab3d_1=CASTWP(puu(:,:,:,Krhs)), clinfo1=' zad  - Ua: ', mask1=umask,   &
118         &                                  tab3d_2=CASTWP(pvv(:,:,:,Krhs)), clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
119      !
120      IF( ln_timing )   CALL timing_stop('dyn_zad')
121      !
122   END SUBROUTINE dyn_zad
123
124   !!======================================================================
125END MODULE dynzad
Note: See TracBrowser for help on using the repository browser.