MODULE dyndmp !!====================================================================== !! *** MODULE dyndmp *** !! Ocean dynamics: internal restoring trend on momentum (U and V current) !!====================================================================== !! History : 3.5 ! 2013-08 (D. Calvert) Original code !!---------------------------------------------------------------------- !!---------------------------------------------------------------------- !! dyn_dmp_alloc : allocate dyndmp arrays !! dyn_dmp_init : namelist read, parameter control and resto coeff. !! dyn_dmp : update the momentum trend with the internal damping !!---------------------------------------------------------------------- USE oce ! ocean: variables USE dom_oce ! ocean: domain variables USE c1d ! 1D vertical configuration USE trdmod ! ocean: trends USE trdmod_oce ! ocean: trends variables USE tradmp ! ocean: internal damping USE zdf_oce ! ocean: vertical physics USE phycst ! physical constants USE dtauvd ! data: U & V current USE zdfmxl ! vertical physics: mixed layer depth USE in_out_manager ! I/O manager USE lib_mpp ! MPP library USE prtctl ! Print control USE wrk_nemo ! Memory allocation USE timing ! Timing IMPLICIT NONE PRIVATE PUBLIC dyn_dmp_init ! routine called by nemogcm.F90 PUBLIC dyn_dmp ! routine called by step_c1d.F90 LOGICAL, PUBLIC :: ln_dyndmp ! Flag for Newtonian damping REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: utrdmp ! damping U current trend (m/s2) REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: vtrdmp ! damping V current trend (m/s2) REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: resto_uv ! restoring coeff. on U & V current !! * Substitutions # include "domzgr_substitute.h90" # include "vectopt_loop_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/OPA 3.3 , NEMO Consortium (2010) !! $Id: dyndmp.F90 3294 2012-01-28 16:44:18Z rblod $ !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) !!---------------------------------------------------------------------- CONTAINS INTEGER FUNCTION dyn_dmp_alloc() !!---------------------------------------------------------------------- !! *** FUNCTION dyn_dmp_alloc *** !!---------------------------------------------------------------------- ALLOCATE( utrdmp(jpi,jpj,jpk), vtrdmp(jpi,jpj,jpk), resto_uv(jpi,jpj,jpk), STAT= dyn_dmp_alloc ) ! IF( lk_mpp ) CALL mpp_sum ( dyn_dmp_alloc ) IF( dyn_dmp_alloc > 0 ) CALL ctl_warn('dyn_dmp_alloc: allocation of arrays failed') ! END FUNCTION dyn_dmp_alloc SUBROUTINE dyn_dmp_init !!---------------------------------------------------------------------- !! *** ROUTINE dyn_dmp_init *** !! !! ** Purpose : Initialization for the Newtonian damping !! !! ** Method : - read the ln_dyndmp parameter from the namc1d_dyndmp namelist !! - allocate damping arrays !! - check the parameters of the namtra_dmp namelist !! - calculate damping coefficient !!---------------------------------------------------------------------- NAMELIST/namc1d_dyndmp/ ln_dyndmp INTEGER :: ios !!---------------------------------------------------------------------- REWIND( numnam_ref ) ! Namelist namc1d_dyndmp in reference namelist : READ ( numnam_ref, namc1d_dyndmp, IOSTAT = ios, ERR = 901) 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namc1d_dyndmp in reference namelist', lwp ) REWIND( numnam_cfg ) ! Namelist namc1d_dyndmp in configuration namelist : Parameters of the run READ ( numnam_cfg, namc1d_dyndmp, IOSTAT = ios, ERR = 902 ) 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namc1d_dyndmp in configuration namelist', lwp ) IF(lwm) WRITE ( numond, namc1d_dyndmp ) IF(lwp) THEN ! control print WRITE(numout,*) WRITE(numout,*) 'dyn_dmp_init : U and V current Newtonian damping' WRITE(numout,*) '~~~~~~~~~~~~' WRITE(numout,*) ' Namelist namc1d_dyndmp : Set damping flag' WRITE(numout,*) ' add a damping term or not ln_dyndmp = ', ln_dyndmp WRITE(numout,*) ' Namelist namtra_dmp : Set damping parameters' WRITE(numout,*) ' horizontal damping option nn_hdmp = ', nn_hdmp WRITE(numout,*) ' mixed layer damping option nn_zdmp = ', nn_zdmp, '(non-C1D zoom: forced to 0)' WRITE(numout,*) ' surface time scale (days) rn_surf = ', rn_surf WRITE(numout,*) ' bottom time scale (days) rn_bot = ', rn_bot WRITE(numout,*) ' depth of transition (meters) rn_dep = ', rn_dep WRITE(numout,*) ' create a damping.coeff file nn_file = ', nn_file WRITE(numout,*) ENDIF IF( ln_dyndmp ) THEN ! !== allocate the data arrays ==! IF( dyn_dmp_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dyn_dmp_init: unable to allocate arrays' ) ! #if ! defined key_c1d SELECT CASE ( nn_hdmp ) !== control print of horizontal option ==! CASE ( -1 ) ; IF(lwp) WRITE(numout,*) ' momentum damping in the Med & Red seas only' CASE ( 1:90 ) ; IF(lwp) WRITE(numout,*) ' momentum damping poleward of', nn_hdmp, ' degrees' CASE DEFAULT WRITE(ctmp1,*) ' bad flag value for nn_hdmp = ', nn_hdmp CALL ctl_stop(ctmp1) END SELECT ! #endif SELECT CASE ( nn_zdmp ) !== control print of vertical option ==! CASE ( 0 ) ; IF(lwp) WRITE(numout,*) ' momentum damping throughout the water column' CASE ( 1 ) ; IF(lwp) WRITE(numout,*) ' no momentum damping in the turbocline (avt > 5 cm2/s)' CASE ( 2 ) ; IF(lwp) WRITE(numout,*) ' no momentum damping in the mixed layer' CASE DEFAULT WRITE(ctmp1,*) ' bad flag value for nn_zdmp = ', nn_zdmp CALL ctl_stop(ctmp1) END SELECT ! IF( .NOT. ln_uvd_dyndmp ) THEN ! force the initialization of U & V current data for damping CALL ctl_warn( 'dyn_dmp_init: U & V current read data not initialized, we force ln_uvd_dyndmp=T' ) CALL dta_uvd_init( ld_dyndmp=ln_dyndmp ) ENDIF ! utrdmp(:,:,:) = 0._wp ! internal damping trends vtrdmp(:,:,:) = 0._wp ! !== Damping coefficients calculation: ==! ! !== use tradmp.F90 subroutines dtacof, dtacof_zoom and cofdis ==! ! !!! NOTE: these need to be altered for use in this module if ! !!! they are to be used outside the C1D context ! !!! (use of U,V grid variables) IF( lzoom .AND. .NOT. lk_c1d ) THEN ; CALL dtacof_zoom( resto_uv ) ELSE ; CALL dtacof( nn_hdmp, rn_surf, rn_bot, rn_dep, nn_file, 'DYN', resto_uv ) ENDIF ! ENDIF ! END SUBROUTINE dyn_dmp_init SUBROUTINE dyn_dmp( kt ) !!---------------------------------------------------------------------- !! *** ROUTINE dyn_dmp *** !! !! ** Purpose : Compute the momentum trends due to a newtonian damping !! of the ocean velocities towards the given data and add it to the !! general momentum trends. !! !! ** Method : Compute Newtonian damping towards u_dta and v_dta !! and add to the general momentum trends: !! ua = ua + resto_uv * (u_dta - ub) !! va = va + resto_uv * (v_dta - vb) !! The trend is computed either throughout the water column !! (nn_zdmp=0), where the vertical mixing is weak (nn_zdmp=1) or !! below the well mixed layer (nn_zdmp=2) !! !! ** Action : - (ua,va) momentum trends updated with the damping trend !!---------------------------------------------------------------------- ! INTEGER, INTENT(in) :: kt ! ocean time-step index !! INTEGER :: ji, jj, jk ! dummy loop indices REAL(wp) :: zua, zva ! local scalars REAL(wp), POINTER, DIMENSION(:,:,:,:) :: zuv_dta ! Read in data !!---------------------------------------------------------------------- ! IF( nn_timing == 1 ) CALL timing_start( 'dyn_dmp' ) ! CALL wrk_alloc( jpi, jpj, jpk, 2, zuv_dta ) ! ! !== read and interpolate U & V current data at kt ==! CALL dta_uvd( kt, zuv_dta ) !!! NOTE: This subroutine must be altered for use outside !!! the C1D context (use of U,V grid variables) ! SELECT CASE ( nn_zdmp ) !== Calculate/add Newtonian damping to the momentum trend ==! ! CASE( 0 ) ! Newtonian damping throughout the water column DO jk = 1, jpkm1 DO jj = 2, jpjm1 DO ji = fs_2, fs_jpim1 ! vector opt. zua = resto_uv(ji,jj,jk) * ( zuv_dta(ji,jj,jk,1) - ub(ji,jj,jk) ) zva = resto_uv(ji,jj,jk) * ( zuv_dta(ji,jj,jk,2) - vb(ji,jj,jk) ) ua(ji,jj,jk) = ua(ji,jj,jk) + zua va(ji,jj,jk) = va(ji,jj,jk) + zva utrdmp(ji,jj,jk) = zua ! save the trends vtrdmp(ji,jj,jk) = zva END DO END DO END DO ! CASE ( 1 ) ! no damping above the turbocline (avt > 5 cm2/s) DO jk = 1, jpkm1 DO jj = 2, jpjm1 DO ji = fs_2, fs_jpim1 ! vector opt. IF( avt(ji,jj,jk) <= 5.e-4_wp ) THEN zua = resto_uv(ji,jj,jk) * ( zuv_dta(ji,jj,jk,1) - ub(ji,jj,jk) ) zva = resto_uv(ji,jj,jk) * ( zuv_dta(ji,jj,jk,2) - vb(ji,jj,jk) ) ELSE zua = 0._wp zva = 0._wp ENDIF ua(ji,jj,jk) = ua(ji,jj,jk) + zua va(ji,jj,jk) = va(ji,jj,jk) + zva utrdmp(ji,jj,jk) = zua ! save the trends vtrdmp(ji,jj,jk) = zva END DO END DO END DO ! CASE ( 2 ) ! no damping in the mixed layer DO jk = 1, jpkm1 DO jj = 2, jpjm1 DO ji = fs_2, fs_jpim1 ! vector opt. IF( fsdept(ji,jj,jk) >= hmlp (ji,jj) ) THEN zua = resto_uv(ji,jj,jk) * ( zuv_dta(ji,jj,jk,1) - ub(ji,jj,jk) ) zva = resto_uv(ji,jj,jk) * ( zuv_dta(ji,jj,jk,2) - vb(ji,jj,jk) ) ELSE zua = 0._wp zva = 0._wp ENDIF ua(ji,jj,jk) = ua(ji,jj,jk) + zua va(ji,jj,jk) = va(ji,jj,jk) + zva utrdmp(ji,jj,jk) = zua ! save the trends vtrdmp(ji,jj,jk) = zva END DO END DO END DO ! END SELECT ! ! ! Trend diagnostic IF( l_trddyn ) CALL trd_mod( utrdmp, vtrdmp, jpdyn_trd_dat, 'DYN', kt ) ! ! ! Control print IF( ln_ctl ) CALL prt_ctl( tab3d_1=ua(:,:,:), clinfo1=' dmp - Ua: ', mask1=umask, & & tab3d_2=va(:,:,:), clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) ! CALL wrk_dealloc( jpi, jpj, jpk, 2, zuv_dta ) ! IF( nn_timing == 1 ) CALL timing_stop( 'dyn_dmp') ! END SUBROUTINE dyn_dmp !!====================================================================== END MODULE dyndmp