Changeset 3294 for trunk/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90
- Timestamp:
- 2012-01-28T17:44:18+01:00 (12 years ago)
- File:
-
- 1 edited
-
trunk/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90 (modified) (5 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90
r2528 r3294 8 8 !! 3.0 ! 2008-04 (NEMO team) add in the reference version 9 9 !! 3.3 ! 2010-09 (D.Storkey and E.O'Dea) bug fixes 10 !! 3.4 ! 2011 (D. Storkey) rewrite in preparation for OBC-BDY merge 10 11 !!---------------------------------------------------------------------- 11 12 #if defined key_bdy 12 13 !!---------------------------------------------------------------------- 13 !! 'key_bdy' UnstructuredOpen Boundary Condition14 !! 'key_bdy' Open Boundary Condition 14 15 !!---------------------------------------------------------------------- 15 16 !! PUBLIC 16 !! tide_init : read of namelist 17 !! tide_data : read in and initialisation of tidal constituents at boundary 17 !! tide_init : read of namelist and initialisation of tidal harmonics data 18 18 !! tide_update : calculation of tidal forcing at each timestep 19 19 !! PRIVATE … … 24 24 !! vset :/ 25 25 !!---------------------------------------------------------------------- 26 USE timing ! Timing 26 27 USE oce ! ocean dynamics and tracers 27 28 USE dom_oce ! ocean space and time domain … … 33 34 USE bdy_oce ! ocean open boundary conditions 34 35 USE daymod ! calendar 36 USE fldread, ONLY: fld_map 35 37 36 38 IMPLICIT NONE 37 39 PRIVATE 38 40 39 PUBLIC tide_init ! routine called in bdyini 40 PUBLIC tide_data ! routine called in bdyini 41 PUBLIC tide_init ! routine called in nemo_init 41 42 PUBLIC tide_update ! routine called in bdydyn 42 43 43 LOGICAL, PUBLIC :: ln_tide_date !: =T correct tide phases and amplitude for model start date 44 INTEGER, PUBLIC, PARAMETER :: jptides_max = 15 !: Max number of tidal contituents 45 INTEGER, PUBLIC :: ntide !: Actual number of tidal constituents 46 47 CHARACTER(len=80), PUBLIC :: filtide !: Filename root for tidal input files 48 CHARACTER(len= 4), PUBLIC, DIMENSION(jptides_max) :: tide_cpt !: Names of tidal components used. 49 50 INTEGER , PUBLIC, DIMENSION(jptides_max) :: nindx !: ??? 51 REAL(wp), PUBLIC, DIMENSION(jptides_max) :: tide_speed !: Phase speed of tidal constituent (deg/hr) 52 53 REAL(wp), DIMENSION(jpbdim,jptides_max) :: ssh1, ssh2 ! Tidal constituents : SSH 54 REAL(wp), DIMENSION(jpbdim,jptides_max) :: u1 , u2 ! Tidal constituents : U 55 REAL(wp), DIMENSION(jpbdim,jptides_max) :: v1 , v2 ! Tidal constituents : V 44 TYPE, PUBLIC :: TIDES_DATA !: Storage for external tidal harmonics data 45 INTEGER :: ncpt !: Actual number of tidal components 46 REAL(wp), POINTER, DIMENSION(:) :: speed !: Phase speed of tidal constituent (deg/hr) 47 REAL(wp), POINTER, DIMENSION(:,:,:) :: ssh !: Tidal constituents : SSH 48 REAL(wp), POINTER, DIMENSION(:,:,:) :: u !: Tidal constituents : U 49 REAL(wp), POINTER, DIMENSION(:,:,:) :: v !: Tidal constituents : V 50 END TYPE TIDES_DATA 51 52 INTEGER, PUBLIC, PARAMETER :: jptides_max = 15 !: Max number of tidal contituents 53 54 TYPE(TIDES_DATA), PUBLIC, DIMENSION(jp_bdy), TARGET :: tides !: External tidal harmonics data 56 55 57 56 !!---------------------------------------------------------------------- … … 66 65 !! *** SUBROUTINE tide_init *** 67 66 !! 68 !! ** Purpose : - Read in namelist for tides 69 !! 70 !!---------------------------------------------------------------------- 71 INTEGER :: itide ! dummy loop index 67 !! ** Purpose : - Read in namelist for tides and initialise external 68 !! tidal harmonics data 69 !! 70 !!---------------------------------------------------------------------- 71 !! namelist variables 72 !!------------------- 73 LOGICAL :: ln_tide_date !: =T correct tide phases and amplitude for model start date 74 CHARACTER(len=80) :: filtide !: Filename root for tidal input files 75 CHARACTER(len= 4), DIMENSION(jptides_max) :: tide_cpt !: Names of tidal components used. 76 REAL(wp), DIMENSION(jptides_max) :: tide_speed !: Phase speed of tidal constituent (deg/hr) 77 !! 78 INTEGER, DIMENSION(jptides_max) :: nindx !: index of pre-set tidal components 79 INTEGER :: ib_bdy, itide, ib !: dummy loop indices 80 INTEGER :: inum, igrd 81 INTEGER, DIMENSION(3) :: ilen0 !: length of boundary data (from OBC arrays) 82 CHARACTER(len=80) :: clfile !: full file name for tidal input file 83 REAL(wp) :: z_arg, z_atde, z_btde, z1t, z2t 84 REAL(wp),DIMENSION(jptides_max) :: z_vplu, z_ftc 85 REAL(wp),ALLOCATABLE, DIMENSION(:,:,:) :: dta_read !: work space to read in tidal harmonics data 86 !! 87 TYPE(TIDES_DATA), POINTER :: td !: local short cut 72 88 !! 73 89 NAMELIST/nambdy_tide/ln_tide_date, filtide, tide_cpt, tide_speed 74 90 !!---------------------------------------------------------------------- 91 92 IF( nn_timing == 1 ) CALL timing_start('tide_init') 75 93 76 94 IF(lwp) WRITE(numout,*) … … 78 96 IF(lwp) WRITE(numout,*) '~~~~~~~~~' 79 97 80 ! Namelist nambdy_tide : tidal harmonic forcing at open boundaries 81 ln_tide_date = .false. 82 filtide(:) = '' 83 tide_speed(:) = 0.0 84 tide_cpt(:) = '' 85 REWIND( numnam ) ! Read namelist parameters 86 READ ( numnam, nambdy_tide ) 87 ! ! Count number of components specified 88 ntide = jptides_max 89 DO itide = 1, jptides_max 90 IF( tide_cpt(itide) == '' ) THEN 91 ntide = itide-1 92 exit 93 ENDIF 94 END DO 95 96 ! ! find constituents in standard list 97 DO itide = 1, ntide 98 nindx(itide) = 0 99 IF( TRIM( tide_cpt(itide) ) == 'Q1' ) nindx(itide) = 1 100 IF( TRIM( tide_cpt(itide) ) == 'O1' ) nindx(itide) = 2 101 IF( TRIM( tide_cpt(itide) ) == 'P1' ) nindx(itide) = 3 102 IF( TRIM( tide_cpt(itide) ) == 'S1' ) nindx(itide) = 4 103 IF( TRIM( tide_cpt(itide) ) == 'K1' ) nindx(itide) = 5 104 IF( TRIM( tide_cpt(itide) ) == '2N2' ) nindx(itide) = 6 105 IF( TRIM( tide_cpt(itide) ) == 'MU2' ) nindx(itide) = 7 106 IF( TRIM( tide_cpt(itide) ) == 'N2' ) nindx(itide) = 8 107 IF( TRIM( tide_cpt(itide) ) == 'NU2' ) nindx(itide) = 9 108 IF( TRIM( tide_cpt(itide) ) == 'M2' ) nindx(itide) = 10 109 IF( TRIM( tide_cpt(itide) ) == 'L2' ) nindx(itide) = 11 110 IF( TRIM( tide_cpt(itide) ) == 'T2' ) nindx(itide) = 12 111 IF( TRIM( tide_cpt(itide) ) == 'S2' ) nindx(itide) = 13 112 IF( TRIM( tide_cpt(itide) ) == 'K2' ) nindx(itide) = 14 113 IF( TRIM( tide_cpt(itide) ) == 'M4' ) nindx(itide) = 15 114 IF( nindx(itide) == 0 .AND. lwp ) THEN 115 WRITE(ctmp1,*) 'constitunent', itide,':', tide_cpt(itide), 'not in standard list' 116 CALL ctl_warn( ctmp1 ) 117 ENDIF 118 END DO 119 ! ! Parameter control and print 120 IF( ntide < 1 ) THEN 121 CALL ctl_stop( ' Did not find any tidal components in namelist nambdy_tide' ) 122 ELSE 123 IF(lwp) WRITE(numout,*) ' Namelist nambdy_tide : tidal harmonic forcing at open boundaries' 124 IF(lwp) WRITE(numout,*) ' tidal components specified ', ntide 125 IF(lwp) WRITE(numout,*) ' ', tide_cpt(1:ntide) 126 IF(lwp) WRITE(numout,*) ' associated phase speeds (deg/hr) : ' 127 IF(lwp) WRITE(numout,*) ' ', tide_speed(1:ntide) 128 ENDIF 129 130 ! Initialisation of tidal harmonics arrays 131 sshtide(:) = 0.e0 132 utide (:) = 0.e0 133 vtide (:) = 0.e0 134 ! 98 REWIND(numnam) 99 DO ib_bdy = 1, nb_bdy 100 IF( nn_dyn2d_dta(ib_bdy) .ge. 2 ) THEN 101 102 td => tides(ib_bdy) 103 104 ! Namelist nambdy_tide : tidal harmonic forcing at open boundaries 105 ln_tide_date = .false. 106 filtide(:) = '' 107 tide_speed(:) = 0.0 108 tide_cpt(:) = '' 109 110 ! Don't REWIND here - may need to read more than one of these namelists. 111 READ ( numnam, nambdy_tide ) 112 ! ! Count number of components specified 113 td%ncpt = 0 114 DO itide = 1, jptides_max 115 IF( tide_cpt(itide) /= '' ) THEN 116 td%ncpt = td%ncpt + 1 117 ENDIF 118 END DO 119 120 ! Fill in phase speeds from namelist 121 ALLOCATE( td%speed(td%ncpt) ) 122 td%speed = tide_speed(1:td%ncpt) 123 124 ! Find constituents in standard list 125 DO itide = 1, td%ncpt 126 nindx(itide) = 0 127 IF( TRIM( tide_cpt(itide) ) == 'Q1' ) nindx(itide) = 1 128 IF( TRIM( tide_cpt(itide) ) == 'O1' ) nindx(itide) = 2 129 IF( TRIM( tide_cpt(itide) ) == 'P1' ) nindx(itide) = 3 130 IF( TRIM( tide_cpt(itide) ) == 'S1' ) nindx(itide) = 4 131 IF( TRIM( tide_cpt(itide) ) == 'K1' ) nindx(itide) = 5 132 IF( TRIM( tide_cpt(itide) ) == '2N2' ) nindx(itide) = 6 133 IF( TRIM( tide_cpt(itide) ) == 'MU2' ) nindx(itide) = 7 134 IF( TRIM( tide_cpt(itide) ) == 'N2' ) nindx(itide) = 8 135 IF( TRIM( tide_cpt(itide) ) == 'NU2' ) nindx(itide) = 9 136 IF( TRIM( tide_cpt(itide) ) == 'M2' ) nindx(itide) = 10 137 IF( TRIM( tide_cpt(itide) ) == 'L2' ) nindx(itide) = 11 138 IF( TRIM( tide_cpt(itide) ) == 'T2' ) nindx(itide) = 12 139 IF( TRIM( tide_cpt(itide) ) == 'S2' ) nindx(itide) = 13 140 IF( TRIM( tide_cpt(itide) ) == 'K2' ) nindx(itide) = 14 141 IF( TRIM( tide_cpt(itide) ) == 'M4' ) nindx(itide) = 15 142 IF( nindx(itide) == 0 .AND. lwp ) THEN 143 WRITE(ctmp1,*) 'constitunent', itide,':', tide_cpt(itide), 'not in standard list' 144 CALL ctl_warn( ctmp1 ) 145 ENDIF 146 END DO 147 ! ! Parameter control and print 148 IF( td%ncpt < 1 ) THEN 149 CALL ctl_stop( ' Did not find any tidal components in namelist nambdy_tide' ) 150 ELSE 151 IF(lwp) WRITE(numout,*) ' Namelist nambdy_tide : tidal harmonic forcing at open boundaries' 152 IF(lwp) WRITE(numout,*) ' tidal components specified ', td%ncpt 153 IF(lwp) WRITE(numout,*) ' ', tide_cpt(1:td%ncpt) 154 IF(lwp) WRITE(numout,*) ' associated phase speeds (deg/hr) : ' 155 IF(lwp) WRITE(numout,*) ' ', tide_speed(1:td%ncpt) 156 ENDIF 157 158 159 ! Allocate space for tidal harmonics data - 160 ! get size from OBC data arrays 161 ! --------------------------------------- 162 163 ilen0(1) = SIZE( dta_bdy(ib_bdy)%ssh ) 164 ALLOCATE( td%ssh( ilen0(1), td%ncpt, 2 ) ) 165 166 ilen0(2) = SIZE( dta_bdy(ib_bdy)%u2d ) 167 ALLOCATE( td%u( ilen0(2), td%ncpt, 2 ) ) 168 169 ilen0(3) = SIZE( dta_bdy(ib_bdy)%v2d ) 170 ALLOCATE( td%v( ilen0(3), td%ncpt, 2 ) ) 171 172 ALLOCATE( dta_read( MAXVAL(ilen0), 1, 1 ) ) 173 174 175 ! Open files and read in tidal forcing data 176 ! ----------------------------------------- 177 178 DO itide = 1, td%ncpt 179 ! ! SSH fields 180 clfile = TRIM(filtide)//TRIM(tide_cpt(itide))//'_grid_T.nc' 181 IF(lwp) WRITE(numout,*) 'Reading data from file ', clfile 182 CALL iom_open( clfile, inum ) 183 CALL fld_map( inum, 'z1' , dta_read(1:ilen0(1),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,1) ) 184 td%ssh(:,itide,1) = dta_read(1:ilen0(1),1,1) 185 CALL fld_map( inum, 'z2' , dta_read(1:ilen0(1),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,1) ) 186 td%ssh(:,itide,2) = dta_read(1:ilen0(1),1,1) 187 CALL iom_close( inum ) 188 ! ! U fields 189 clfile = TRIM(filtide)//TRIM(tide_cpt(itide))//'_grid_U.nc' 190 IF(lwp) WRITE(numout,*) 'Reading data from file ', clfile 191 CALL iom_open( clfile, inum ) 192 CALL fld_map( inum, 'u1' , dta_read(1:ilen0(2),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,2) ) 193 td%u(:,itide,1) = dta_read(1:ilen0(2),1,1) 194 CALL fld_map( inum, 'u2' , dta_read(1:ilen0(2),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,2) ) 195 td%u(:,itide,2) = dta_read(1:ilen0(2),1,1) 196 CALL iom_close( inum ) 197 ! ! V fields 198 clfile = TRIM(filtide)//TRIM(tide_cpt(itide))//'_grid_V.nc' 199 IF(lwp) WRITE(numout,*) 'Reading data from file ', clfile 200 CALL iom_open( clfile, inum ) 201 CALL fld_map( inum, 'v1' , dta_read(1:ilen0(3),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,3) ) 202 td%v(:,itide,1) = dta_read(1:ilen0(3),1,1) 203 CALL fld_map( inum, 'v2' , dta_read(1:ilen0(3),1:1,1:1) , 1, idx_bdy(ib_bdy)%nbmap(:,3) ) 204 td%v(:,itide,2) = dta_read(1:ilen0(3),1,1) 205 CALL iom_close( inum ) 206 ! 207 END DO ! end loop on tidal components 208 209 IF( ln_tide_date ) THEN ! correct for date factors 210 211 !! used nmonth, nyear and nday from daymod.... 212 ! Calculate date corrects for 15 standard consituents 213 ! This is the initialisation step, so nday, nmonth, nyear are the 214 ! initial date/time of the integration. 215 print *, nday,nmonth,nyear 216 nyear = int(ndate0 / 10000 ) ! initial year 217 nmonth = int((ndate0 - nyear * 10000 ) / 100 ) ! initial month 218 nday = int(ndate0 - nyear * 10000 - nmonth * 100) 219 220 CALL uvset( 0, nday, nmonth, nyear, z_ftc, z_vplu ) 221 222 IF(lwp) WRITE(numout,*) 'Correcting tide for date:', nday, nmonth, nyear 223 224 DO itide = 1, td%ncpt ! loop on tidal components 225 ! 226 IF( nindx(itide) /= 0 ) THEN 227 !!gm use rpi and rad global variable 228 z_arg = 3.14159265d0 * z_vplu(nindx(itide)) / 180.0d0 229 z_atde=z_ftc(nindx(itide))*cos(z_arg) 230 z_btde=z_ftc(nindx(itide))*sin(z_arg) 231 IF(lwp) WRITE(numout,'(2i5,8f10.6)') itide, nindx(itide), td%speed(itide), & 232 & z_ftc(nindx(itide)), z_vplu(nindx(itide)) 233 ELSE 234 z_atde = 1.0_wp 235 z_btde = 0.0_wp 236 ENDIF 237 ! ! elevation 238 igrd = 1 239 DO ib = 1, ilen0(igrd) 240 z1t = z_atde * td%ssh(ib,itide,1) + z_btde * td%ssh(ib,itide,2) 241 z2t = z_atde * td%ssh(ib,itide,2) - z_btde * td%ssh(ib,itide,1) 242 td%ssh(ib,itide,1) = z1t 243 td%ssh(ib,itide,2) = z2t 244 END DO 245 ! ! u 246 igrd = 2 247 DO ib = 1, ilen0(igrd) 248 z1t = z_atde * td%u(ib,itide,1) + z_btde * td%u(ib,itide,2) 249 z2t = z_atde * td%u(ib,itide,2) - z_btde * td%u(ib,itide,1) 250 td%u(ib,itide,1) = z1t 251 td%u(ib,itide,2) = z2t 252 END DO 253 ! ! v 254 igrd = 3 255 DO ib = 1, ilen0(igrd) 256 z1t = z_atde * td%v(ib,itide,1) + z_btde * td%v(ib,itide,2) 257 z2t = z_atde * td%v(ib,itide,2) - z_btde * td%v(ib,itide,1) 258 td%v(ib,itide,1) = z1t 259 td%v(ib,itide,2) = z2t 260 END DO 261 ! 262 END DO ! end loop on tidal components 263 ! 264 ENDIF ! date correction 265 ! 266 ENDIF ! nn_dyn2d_dta(ib_bdy) .ge. 2 267 ! 268 END DO ! loop on ib_bdy 269 270 IF( nn_timing == 1 ) CALL timing_stop('tide_init') 271 135 272 END SUBROUTINE tide_init 136 273 137 274 138 SUBROUTINE tide_data 139 !!---------------------------------------------------------------------- 140 !! *** SUBROUTINE tide_data *** 141 !! 142 !! ** Purpose : - Read in tidal harmonics data and adjust for the start 143 !! time of the model run. 144 !! 145 !!---------------------------------------------------------------------- 146 INTEGER :: itide, igrd, ib ! dummy loop indices 147 CHARACTER(len=80) :: clfile ! full file name for tidal input file 148 INTEGER :: ipi, ipj, inum, idvar ! temporary integers (netcdf read) 149 INTEGER, DIMENSION(6) :: lendta=0 ! length of data in the file (note may be different from nblendta!) 150 REAL(wp) :: z_arg, z_atde, z_btde, z1t, z2t 151 REAL(wp), DIMENSION(jpbdta,1) :: zdta ! temporary array for data fields 152 REAL(wp), DIMENSION(jptides_max) :: z_vplu, z_ftc 153 !!------------------------------------------------------------------------------ 154 155 ! Open files and read in tidal forcing data 156 ! ----------------------------------------- 157 158 ipj = 1 159 160 DO itide = 1, ntide 161 ! ! SSH fields 162 clfile = TRIM(filtide)//TRIM(tide_cpt(itide))//'_grid_T.nc' 163 IF(lwp) WRITE(numout,*) 'Reading data from file ', clfile 164 CALL iom_open( clfile, inum ) 165 igrd = 4 166 IF( nblendta(igrd) <= 0 ) THEN 167 idvar = iom_varid( inum,'z1' ) 168 IF(lwp) WRITE(numout,*) 'iom_file(1)%ndims(idvar) : ',iom_file%ndims(idvar) 169 nblendta(igrd) = iom_file(inum)%dimsz(1,idvar) 170 WRITE(numout,*) 'Dim size for z1 is ', nblendta(igrd) 171 ENDIF 172 ipi = nblendta(igrd) 173 CALL iom_get( inum, jpdom_unknown, 'z1', zdta(1:ipi,1:ipj) ) 174 DO ib = 1, nblenrim(igrd) 175 ssh1(ib,itide) = zdta(nbmap(ib,igrd),1) 176 END DO 177 CALL iom_get( inum, jpdom_unknown, 'z2', zdta(1:ipi,1:ipj) ) 178 DO ib = 1, nblenrim(igrd) 179 ssh2(ib,itide) = zdta(nbmap(ib,igrd),1) 180 END DO 181 CALL iom_close( inum ) 182 ! 183 ! ! U fields 184 clfile = TRIM(filtide)//TRIM(tide_cpt(itide))//'_grid_U.nc' 185 IF(lwp) WRITE(numout,*) 'Reading data from file ', clfile 186 CALL iom_open( clfile, inum ) 187 igrd = 5 188 IF( lendta(igrd) <= 0 ) THEN 189 idvar = iom_varid( inum,'u1' ) 190 lendta(igrd) = iom_file(inum)%dimsz(1,idvar) 191 WRITE(numout,*) 'Dim size for u1 is ',lendta(igrd) 192 ENDIF 193 ipi = lendta(igrd) 194 CALL iom_get( inum, jpdom_unknown, 'u1', zdta(1:ipi,1:ipj) ) 195 DO ib = 1, nblenrim(igrd) 196 u1(ib,itide) = zdta(nbmap(ib,igrd),1) 197 END DO 198 CALL iom_get( inum, jpdom_unknown, 'u2', zdta(1:ipi,1:ipj) ) 199 DO ib = 1, nblenrim(igrd) 200 u2(ib,itide) = zdta(nbmap(ib,igrd),1) 201 END DO 202 CALL iom_close( inum ) 203 ! 204 ! ! V fields 205 clfile = TRIM(filtide)//TRIM(tide_cpt(itide))//'_grid_V.nc' 206 if(lwp) write(numout,*) 'Reading data from file ', clfile 207 CALL iom_open( clfile, inum ) 208 igrd = 6 209 IF( lendta(igrd) <= 0 ) THEN 210 idvar = iom_varid( inum,'v1' ) 211 lendta(igrd) = iom_file(inum)%dimsz(1,idvar) 212 WRITE(numout,*) 'Dim size for v1 is ', lendta(igrd) 213 ENDIF 214 ipi = lendta(igrd) 215 CALL iom_get( inum, jpdom_unknown, 'v1', zdta(1:ipi,1:ipj) ) 216 DO ib = 1, nblenrim(igrd) 217 v1(ib,itide) = zdta(nbmap(ib,igrd),1) 218 END DO 219 CALL iom_get( inum, jpdom_unknown, 'v2', zdta(1:ipi,1:ipj) ) 220 DO ib=1, nblenrim(igrd) 221 v2(ib,itide) = zdta(nbmap(ib,igrd),1) 222 END DO 223 CALL iom_close( inum ) 224 ! 225 END DO ! end loop on tidal components 226 227 IF( ln_tide_date ) THEN ! correct for date factors 228 229 !! used nmonth, nyear and nday from daymod.... 230 ! Calculate date corrects for 15 standard consituents 231 ! This is the initialisation step, so nday, nmonth, nyear are the 232 ! initial date/time of the integration. 233 print *, nday,nmonth,nyear 234 nyear = int(ndate0 / 10000 ) ! initial year 235 nmonth = int((ndate0 - nyear * 10000 ) / 100 ) ! initial month 236 nday = int(ndate0 - nyear * 10000 - nmonth * 100) 237 238 CALL uvset( 0, nday, nmonth, nyear, z_ftc, z_vplu ) 239 240 IF(lwp) WRITE(numout,*) 'Correcting tide for date:', nday, nmonth, nyear 241 242 DO itide = 1, ntide ! loop on tidal components 243 ! 244 IF( nindx(itide) /= 0 ) THEN 245 !!gm use rpi and rad global variable 246 z_arg = 3.14159265d0 * z_vplu(nindx(itide)) / 180.0d0 247 z_atde=z_ftc(nindx(itide))*cos(z_arg) 248 z_btde=z_ftc(nindx(itide))*sin(z_arg) 249 IF(lwp) WRITE(numout,'(2i5,8f10.6)') itide, nindx(itide), tide_speed(itide), & 250 & z_ftc(nindx(itide)), z_vplu(nindx(itide)) 251 ELSE 252 z_atde = 1.0_wp 253 z_btde = 0.0_wp 254 ENDIF 255 ! ! elevation 256 igrd = 4 257 DO ib = 1, nblenrim(igrd) 258 z1t = z_atde * ssh1(ib,itide) + z_btde * ssh2(ib,itide) 259 z2t = z_atde * ssh2(ib,itide) - z_btde * ssh1(ib,itide) 260 ssh1(ib,itide) = z1t 261 ssh2(ib,itide) = z2t 262 END DO 263 ! ! u 264 igrd = 5 265 DO ib = 1, nblenrim(igrd) 266 z1t = z_atde * u1(ib,itide) + z_btde * u2(ib,itide) 267 z2t = z_atde * u2(ib,itide) - z_btde * u1(ib,itide) 268 u1(ib,itide) = z1t 269 u2(ib,itide) = z2t 270 END DO 271 ! ! v 272 igrd = 6 273 DO ib = 1, nblenrim(igrd) 274 z1t = z_atde * v1(ib,itide) + z_btde * v2(ib,itide) 275 z2t = z_atde * v2(ib,itide) - z_btde * v1(ib,itide) 276 v1(ib,itide) = z1t 277 v2(ib,itide) = z2t 278 END DO 279 ! 280 END DO ! end loop on tidal components 281 ! 282 ENDIF ! date correction 283 ! 284 END SUBROUTINE tide_data 285 286 287 SUBROUTINE tide_update ( kt, jit ) 275 SUBROUTINE tide_update ( kt, idx, dta, td, jit, time_offset ) 288 276 !!---------------------------------------------------------------------- 289 277 !! *** SUBROUTINE tide_update *** 290 278 !! 291 !! ** Purpose : - Add tidal forcing to ssh bdy, ubtbdy and vbtbdyarrays.279 !! ** Purpose : - Add tidal forcing to ssh, u2d and v2d OBC data arrays. 292 280 !! 293 281 !!---------------------------------------------------------------------- 294 INTEGER, INTENT( in ) :: kt ! Main timestep counter282 INTEGER, INTENT( in ) :: kt ! Main timestep counter 295 283 !!gm doctor jit ==> kit 296 INTEGER, INTENT( in ) :: jit ! Barotropic timestep counter (for timesplitting option) 297 !! 298 INTEGER :: itide, igrd, ib ! dummy loop indices 299 REAL(wp) :: z_arg, z_sarg ! 284 TYPE(OBC_INDEX), INTENT( in ) :: idx ! OBC indices 285 TYPE(OBC_DATA), INTENT(inout) :: dta ! OBC external data 286 TYPE(TIDES_DATA),INTENT( in ) :: td ! tidal harmonics data 287 INTEGER,INTENT(in),OPTIONAL :: jit ! Barotropic timestep counter (for timesplitting option) 288 INTEGER,INTENT( in ), OPTIONAL :: time_offset ! time offset in units of timesteps. NB. if jit 289 ! is present then units = subcycle timesteps. 290 ! time_offset = 0 => get data at "now" time level 291 ! time_offset = -1 => get data at "before" time level 292 ! time_offset = +1 => get data at "after" time level 293 ! etc. 294 !! 295 INTEGER :: itide, igrd, ib ! dummy loop indices 296 INTEGER :: time_add ! time offset in units of timesteps 297 REAL(wp) :: z_arg, z_sarg 300 298 REAL(wp), DIMENSION(jptides_max) :: z_sist, z_cost 301 299 !!---------------------------------------------------------------------- 302 300 301 IF( nn_timing == 1 ) CALL timing_start('tide_update') 302 303 time_add = 0 304 IF( PRESENT(time_offset) ) THEN 305 time_add = time_offset 306 ENDIF 307 303 308 ! Note tide phase speeds are in deg/hour, so we need to convert the 304 309 ! elapsed time in seconds to hours by dividing by 3600.0 305 IF( jit == 0 ) THEN 306 z_arg = kt * rdt * rad / 3600.0 307 ELSE ! we are in a barotropic subcycle (for timesplitting option) 308 ! z_arg = ( (kt-1) * rdt + jit * rdt / REAL(nn_baro,lwp) ) * rad / 3600.0 309 z_arg = ( (kt-1) * rdt + jit * rdt / REAL(nn_baro,wp) ) * rad / 3600.0 310 IF( PRESENT(jit) ) THEN 311 z_arg = ( (kt-1) * rdt + (jit+time_add) * rdt / REAL(nn_baro,wp) ) * rad / 3600.0 312 ELSE 313 z_arg = (kt+time_add) * rdt * rad / 3600.0 310 314 ENDIF 311 315 312 DO itide = 1, ntide313 z_sarg = z_arg * t ide_speed(itide)316 DO itide = 1, td%ncpt 317 z_sarg = z_arg * td%speed(itide) 314 318 z_cost(itide) = COS( z_sarg ) 315 319 z_sist(itide) = SIN( z_sarg ) 316 320 END DO 317 321 318 ! summing of tidal constituents into BDY arrays 319 sshtide(:) = 0.0 320 utide (:) = 0.0 321 vtide (:) = 0.0 322 ! 323 DO itide = 1, ntide 324 igrd=4 ! SSH on tracer grid. 325 DO ib = 1, nblenrim(igrd) 326 sshtide(ib) =sshtide(ib)+ ssh1(ib,itide)*z_cost(itide) + ssh2(ib,itide)*z_sist(itide) 327 ! if(lwp) write(numout,*) 'z',ib,itide,sshtide(ib), ssh1(ib,itide),ssh2(ib,itide) 322 DO itide = 1, td%ncpt 323 igrd=1 ! SSH on tracer grid. 324 DO ib = 1, idx%nblenrim(igrd) 325 dta%ssh(ib) = dta%ssh(ib) + td%ssh(ib,itide,1)*z_cost(itide) + td%ssh(ib,itide,2)*z_sist(itide) 326 ! if(lwp) write(numout,*) 'z', ib, itide, dta%ssh(ib), td%ssh(ib,itide,1),td%ssh(ib,itide,2) 328 327 END DO 329 igrd= 5! U grid330 DO ib=1, nblenrim(igrd)331 utide(ib) = utide(ib)+ u1(ib,itide)*z_cost(itide) + u2(ib,itide)*z_sist(itide)332 ! if(lwp) write(numout,*) 'u',ib,itide,utide(ib), u1(ib,itide),u2(ib,itide)328 igrd=2 ! U grid 329 DO ib=1, idx%nblenrim(igrd) 330 dta%u2d(ib) = dta%u2d(ib) + td%u(ib,itide,1)*z_cost(itide) + td%u(ib,itide,2)*z_sist(itide) 331 ! if(lwp) write(numout,*) 'u',ib,itide,utide(ib), td%u(ib,itide,1),td%u(ib,itide,2) 333 332 END DO 334 igrd= 6! V grid335 DO ib=1, nblenrim(igrd)336 vtide(ib) = vtide(ib)+ v1(ib,itide)*z_cost(itide) + v2(ib,itide)*z_sist(itide)337 ! if(lwp) write(numout,*) 'v',ib,itide,vtide(ib), v1(ib,itide),v2(ib,itide)333 igrd=3 ! V grid 334 DO ib=1, idx%nblenrim(igrd) 335 dta%v2d(ib) = dta%v2d(ib) + td%v(ib,itide,1)*z_cost(itide) + td%v(ib,itide,2)*z_sist(itide) 336 ! if(lwp) write(numout,*) 'v',ib,itide,vtide(ib), td%v(ib,itide,1),td%v(ib,itide,2) 338 337 END DO 339 338 END DO 339 ! 340 IF( nn_timing == 1 ) CALL timing_stop('tide_update') 340 341 ! 341 342 END SUBROUTINE tide_update
Note: See TracChangeset
for help on using the changeset viewer.
