source: branches/2014/dev_r4650_UKMO3_masked_damping/NEMOGCM/TOOLS/DMP_TOOLS/src/zoom.F90 @ 4739

Last change on this file since 4739 was 4739, checked in by timgraham, 6 years ago

Updated C1D/dyndmp.F90 and trcdmp.F90 to read restoration coefficient from a file.
Modified namelist_top_ref to match new options
Bug fixes to DMP_TOOLS tool and addition of custom.F90 to allow users to make modifications. Also changed to use working precision (wp) throughout.

File size: 3.8 KB
RevLine 
[4738]1MODULE zoom
2   
3   USE utils
4
5   CONTAINS
6
7   SUBROUTINE dtacof_zoom( presto, mask)
8      !!----------------------------------------------------------------------
9      !!                  ***  ROUTINE dtacof_zoom  ***
10      !!
11      !! ** Purpose :   Compute the damping coefficient for zoom domain
12      !!
13      !! ** Method  : - set along closed boundary due to zoom a damping over
14      !!                6 points with a max time scale of 5 days.
15      !!              - ORCA arctic/antarctic zoom: set the damping along
16      !!                south/north boundary over a latitude strip.
17      !!
18      !! ** Action  : - resto, the damping coeff. for T and S
19      !!----------------------------------------------------------------------
[4739]20      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout)  ::   presto   ! restoring coeff. (s-1)
21      REAL(wp), DIMENSION(jpi,jpj), INTENT(in)  ::   mask   ! restoring coeff. (s-1)
[4738]22      !
23      INTEGER  ::   ji, jj, jn   ! dummy loop indices
[4739]24      REAL(wp) ::   zlat, zlat0, zlat1, zlat2, z1_5d   ! local scalar
25      REAL(wp), DIMENSION(6)  ::   zfact               ! 1Dworkspace
26
27      !Namelist variables
28      LOGICAL :: lzoom_w, lzoom_e, lzoom_n, lzoom_s 
29      NAMELIST/nam_zoom_dmp/lzoom_n,lzoom_e,lzoom_w,lzoom_s
[4738]30      !!----------------------------------------------------------------------
31      !
32      IF( nn_timing == 1 )  CALL timing_start( 'dtacof_zoom')
33      !
34     
[4739]35      ! Read namelist
36      OPEN( UNIT=numnam, FILE='namelist', FORM='FORMATTED', STATUS='OLD' )
37      READ( numnam, nam_dmp_create )
38      CLOSE( numnam )
[4738]39
[4739]40      zfact(1) =  1._wp
41      zfact(2) =  1._wp
42      zfact(3) = 11._wp / 12._wp
43      zfact(4) =  8._wp / 12._wp
44      zfact(5) =  4._wp / 12._wp
45      zfact(6) =  1._wp / 12._wp
46      zfact(:) = zfact(:) / ( 5._wp * rday )    ! 5 days max restoring time scale
[4738]47
[4739]48      presto(:,:) = 0._wp
[4738]49
50      ! damping along the forced closed boundary over 6 grid-points
51      DO jn = 1, 6
52         IF( lzoom_w )   presto( jn, : )                    = zfact(jn)   ! west  closed
53         IF( lzoom_s )   presto( : , jn )                    = zfact(jn)   ! south closed
54         IF( lzoom_e )   presto( jpi+1-jn , : ) = zfact(jn)   ! east  closed
55         IF( lzoom_n )   presto( : , jpj+1-jn ) = zfact(jn)   ! north closed
56      END DO
57
58      !                                           ! ====================================================
59      IF( cp_cfz == "arctic" .OR. cp_cfz == "antarctic" ) THEN   !  ORCA configuration : arctic or antarctic zoom
60         !                                        ! ====================================================
[4739]61         WRITE(numout,*)
62         IF(cp_cfz == "arctic" ) WRITE(numout,*) '              dtacof_zoom : ORCA    Arctic zoom'
63         IF(cp_cfz == "antarctic" ) WRITE(numout,*) '           dtacof_zoom : ORCA Antarctic zoom'
64         WRITE(numout,*)
[4738]65         !
66         !                          ! Initialization :
[4739]67         presto(:,:) = 0._wp
68         zlat0 = 10._wp                     ! zlat0 : latitude strip where resto decreases
69         zlat1 = 30._wp                     ! zlat1 : resto = 1 before zlat1
[4738]70         zlat2 = zlat1 + zlat0              ! zlat2 : resto decreases from 1 to 0 between zlat1 and zlat2
[4739]71         z1_5d = 1._wp / ( 5._wp * rday )   ! z1_5d : 1 / 5days
[4738]72
73         DO jj = 1, jpj
74            DO ji = 1, jpi
75               zlat = ABS( gphit(ji,jj) )
76               IF( zlat1 <= zlat .AND. zlat <= zlat2 ) THEN
[4739]77                  presto(ji,jj) = 0.5_wp * z1_5d * (  1._wp - COS( rpi*(zlat2-zlat)/zlat0 )  ) 
[4738]78               ELSEIF( zlat < zlat1 ) THEN
79                  presto(ji,jj) = z1_5d
80               ENDIF
81            END DO
82         END DO
83         !
84      ENDIF
85      !                             ! Mask resto array
86      presto(:,:) = presto(:,:) * mask(:,:)
87
88   END SUBROUTINE dtacof_zoom
89
90END MODULE zoom
Note: See TracBrowser for help on using the repository browser.