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.
sbctide.F90 in NEMO/branches/UKMO/NEMO_4.0.4_CO9_shelf_climate/src/OCE/SBC – NEMO

source: NEMO/branches/UKMO/NEMO_4.0.4_CO9_shelf_climate/src/OCE/SBC/sbctide.F90 @ 15547

Last change on this file since 15547 was 15547, checked in by hadjt, 3 years ago

Moving ftide, utide, v0tide output to sbctides, adding tide_t output to diaharm_fast

File size: 10.0 KB
Line 
1MODULE sbctide
2   !!======================================================================
3   !!                       ***  MODULE  sbctide  ***
4   !! Initialization of tidal forcing
5   !!======================================================================
6   !! History :  9.0  !  2007  (O. Le Galloudec)  Original code
7   !!----------------------------------------------------------------------
8   USE oce            ! ocean dynamics and tracers variables
9   USE dom_oce        ! ocean space and time domain
10   USE phycst         ! physical constant
11   USE daymod         ! calandar
12   USE tideini        !
13   !
14   USE in_out_manager ! I/O units
15   USE iom            ! xIOs server
16   USE ioipsl         ! NetCDF IPSL library
17   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
18
19   USE bdytides ! davbyr - Access to love number
20
21   USE tide_mod
22
23   IMPLICIT NONE
24   PUBLIC
25
26   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   pot_astro   !
27
28   !!----------------------------------------------------------------------
29   !!   tidal potential
30   !!----------------------------------------------------------------------
31   !!   sbc_tide            :
32   !!   tide_init_potential :
33   !!----------------------------------------------------------------------
34
35   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   amp_pot, phi_pot
36   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   amp_load, phi_load
37 
38   !!----------------------------------------------------------------------
39   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
40   !! $Id$
41   !! Software governed by the CeCILL license (see ./LICENSE)
42   !!----------------------------------------------------------------------
43CONTAINS
44
45   SUBROUTINE sbc_tide( kt )
46      !!----------------------------------------------------------------------
47      !!                 ***  ROUTINE sbc_tide  ***
48      !!----------------------------------------------------------------------     
49      INTEGER, INTENT( in ) ::   kt     ! ocean time-step
50      INTEGER               ::   jk     ! dummy loop index
51      INTEGER               ::   nsec_day_orig     ! Temporary variable
52      CHARACTER (len=40) :: tmp_name
53      !!----------------------------------------------------------------------
54     
55      IF( nsec_day == NINT(0.5_wp * rdt) .OR. kt == nit000 ) THEN      ! start a new day
56         !
57         IF( kt == nit000 )THEN
58            ALLOCATE( amp_pot(jpi,jpj,nb_harmo),                      &
59               &      phi_pot(jpi,jpj,nb_harmo), pot_astro(jpi,jpj)   )
60            IF( ln_read_load )THEN
61               ALLOCATE( amp_load(jpi,jpj,nb_harmo), phi_load(jpi,jpj,nb_harmo) )
62               CALL tide_init_load
63            ENDIF
64         ENDIF
65         !
66         IF( ln_read_load )THEN
67            amp_pot(:,:,:) = amp_load(:,:,:)
68            phi_pot(:,:,:) = phi_load(:,:,:)
69         ELSE
70            amp_pot(:,:,:) = 0._wp
71            phi_pot(:,:,:) = 0._wp
72         ENDIF
73         pot_astro(:,:) = 0._wp
74         !
75         ! If the run does not start from midnight then need to initialise tides
76         ! at the start of the current day (only occurs when kt==nit000)
77         ! Temporarily set nsec_day to beginning of day.
78         nsec_day_orig = nsec_day
79         IF ( nsec_day /= NINT(0.5_wp * rdt) ) THEN
80            kt_tide = kt - (nsec_day - 0.5_wp * rdt)/rdt
81            nsec_day = NINT(0.5_wp * rdt)
82         ELSE
83            kt_tide = kt 
84         ENDIF
85         CALL tide_harmo( omega_tide, v0tide, utide, ftide, ntide, nb_harmo )
86         !
87         !
88         IF(lwp) THEN
89            WRITE(numout,*)
90            WRITE(numout,*) 'sbc_tide : Update of the components and (re)Init. the potential at kt=', kt
91            WRITE(numout,*) '~~~~~~~~ '
92            DO jk = 1, nb_harmo
93               WRITE(numout,*) Wave(ntide(jk))%cname_tide, utide(jk), ftide(jk), v0tide(jk), omega_tide(jk)
94            END DO
95         ENDIF
96         !
97         IF( ln_tide_pot )   CALL tide_init_potential
98         !
99         ! Reset nsec_day
100         nsec_day = nsec_day_orig 
101      ENDIF
102
103
104        DO jk = 1, nb_harmo
105
106
107
108
109            tmp_name=TRIM(Wave(ntide(jk))%cname_tide)//'_utide'
110            IF( iom_use(TRIM(tmp_name)) )  THEN
111            !    IF(lwp) WRITE(numout,*) "harm_ana_out: iom_put: ",TRIM(tmp_name),'; shape = ', SHAPE(anau(jh) )
112                CALL iom_put( TRIM(tmp_name), utide(jk) )
113            !ELSE
114            !    IF(lwp) WRITE(numout,*) "harm_ana_out: not requested: ",TRIM(tmp_name)
115            ENDIF   
116
117            tmp_name=TRIM(Wave(ntide(jk))%cname_tide)//'_v0tide'
118            IF( iom_use(TRIM(tmp_name)) )  THEN
119            !    IF(lwp) WRITE(numout,*) "harm_ana_out: iom_put: ",TRIM(tmp_name),'; shape = ', SHAPE(anau(jh) )
120                CALL iom_put( TRIM(tmp_name), v0tide(jk) )
121            !ELSE
122            !    IF(lwp) WRITE(numout,*) "harm_ana_out: not requested: ",TRIM(tmp_name)
123            ENDIF   
124
125            tmp_name=TRIM(Wave(ntide(jk))%cname_tide)//'_v0tide_origin'
126            IF( iom_use(TRIM(tmp_name)) )  THEN
127            !    IF(lwp) WRITE(numout,*) "harm_ana_out: iom_put: ",TRIM(tmp_name),'; shape = ', SHAPE(anau(jh) )
128                CALL iom_put( TRIM(tmp_name), v0linearintercept(jk) )
129            !ELSE
130            !    IF(lwp) WRITE(numout,*) "harm_ana_out: not requested: ",TRIM(tmp_name)
131            ENDIF   
132
133
134            tmp_name=TRIM(Wave(ntide(jk))%cname_tide)//'_ftide'
135            IF( iom_use(TRIM(tmp_name)) )  THEN
136            !    IF(lwp) WRITE(numout,*) "harm_ana_out: iom_put: ",TRIM(tmp_name),'; shape = ', SHAPE(anau(jh) )
137                CALL iom_put( TRIM(tmp_name), ftide(jk) )
138            !ELSE
139            !    IF(lwp) WRITE(numout,*) "harm_ana_out: not requested: ",TRIM(tmp_name)
140            ENDIF   
141
142
143            tmp_name=TRIM(Wave(ntide(jk))%cname_tide)//'_per'
144            IF( iom_use(TRIM(tmp_name)) )  THEN
145            !    IF(lwp) WRITE(numout,*) "harm_ana_out: iom_put: ",TRIM(tmp_name),'; shape = ', SHAPE(anau(jh) )
146                CALL iom_put( TRIM(tmp_name), 2*rpi/(3600.0_wp*omega_tide(jk)) )
147            !ELSE
148            !    IF(lwp) WRITE(numout,*) "harm_ana_out: not requested: ",TRIM(tmp_name)
149            ENDIF   
150
151
152            tmp_name=TRIM(Wave(ntide(jk))%cname_tide)//'_freq'
153            IF( iom_use(TRIM(tmp_name)) )  THEN
154            !    IF(lwp) WRITE(numout,*) "harm_ana_out: iom_put: ",TRIM(tmp_name),'; shape = ', SHAPE(anau(jh) )
155                CALL iom_put( TRIM(tmp_name), omega_tide(jk) )
156            !ELSE
157            !    IF(lwp) WRITE(numout,*) "harm_ana_out: not requested: ",TRIM(tmp_name)
158            ENDIF   
159
160
161        END DO
162      !
163   END SUBROUTINE sbc_tide
164
165
166
167   SUBROUTINE tide_init_potential
168      !!----------------------------------------------------------------------
169      !!                 ***  ROUTINE tide_init_potential  ***
170      !!----------------------------------------------------------------------
171      INTEGER  ::   ji, jj, jk   ! dummy loop indices
172      REAL(wp) ::   zcons, ztmp1, ztmp2, zlat, zlon, ztmp, zamp, zcs   ! local scalar
173      !!----------------------------------------------------------------------
174
175      DO jk = 1, nb_harmo
176         ! davbyr - Insert variable Love number where once was 0.7
177         zcons = dn_love_number * Wave(ntide(jk))%equitide * ftide(jk)
178         ! END davbyr
179         DO ji = 1, jpi
180            DO jj = 1, jpj
181               ztmp1 =  ftide(jk) * amp_pot(ji,jj,jk) * COS( phi_pot(ji,jj,jk) + v0tide(jk) + utide(jk) )
182               ztmp2 = -ftide(jk) * amp_pot(ji,jj,jk) * SIN( phi_pot(ji,jj,jk) + v0tide(jk) + utide(jk) )
183               zlat = gphit(ji,jj)*rad !! latitude en radian
184               zlon = glamt(ji,jj)*rad !! longitude en radian
185               ztmp = v0tide(jk) + utide(jk) + Wave(ntide(jk))%nutide * zlon
186               ! le potentiel est composé des effets des astres:
187               IF    ( Wave(ntide(jk))%nutide == 1 )  THEN  ;  zcs = zcons * SIN( 2._wp*zlat )
188               ELSEIF( Wave(ntide(jk))%nutide == 2 )  THEN  ;  zcs = zcons * COS( zlat )**2
189               ! davbyr - Include long period tidal forcing
190               ELSEIF( Wave(ntide(jk))%nutide == 0 )  THEN  ;  zcs = zcons * (0.5_wp-1.5_wp*SIN(zlat)**2._wp)
191               ! END - davbyr
192               ELSE                                         ;  zcs = 0._wp
193               ENDIF
194               ztmp1 = ztmp1 + zcs * COS( ztmp )
195               ztmp2 = ztmp2 - zcs * SIN( ztmp )
196               zamp = SQRT( ztmp1*ztmp1 + ztmp2*ztmp2 )
197               amp_pot(ji,jj,jk) = zamp
198               phi_pot(ji,jj,jk) = ATAN2( -ztmp2 / MAX( 1.e-10_wp , zamp ) ,   &
199                  &                        ztmp1 / MAX( 1.e-10_wp,  zamp )   )
200            END DO
201         END DO
202      END DO
203      !
204   END SUBROUTINE tide_init_potential
205
206   SUBROUTINE tide_init_load
207      !!----------------------------------------------------------------------
208      !!                 ***  ROUTINE tide_init_load  ***
209      !!----------------------------------------------------------------------
210      INTEGER :: inum                 ! Logical unit of input file
211      INTEGER :: ji, jj, itide        ! dummy loop indices
212      REAL(wp), DIMENSION(jpi,jpj) ::   ztr, zti   !: workspace to read in tidal harmonics data
213      !!----------------------------------------------------------------------
214      IF(lwp) THEN
215         WRITE(numout,*)
216         WRITE(numout,*) 'tide_init_load : Initialization of load potential from file'
217         WRITE(numout,*) '~~~~~~~~~~~~~~ '
218      ENDIF
219      !
220      CALL iom_open ( cn_tide_load , inum )
221      !
222      DO itide = 1, nb_harmo
223         CALL iom_get  ( inum, jpdom_data,TRIM(Wave(ntide(itide))%cname_tide)//'_z1', ztr(:,:) )
224         CALL iom_get  ( inum, jpdom_data,TRIM(Wave(ntide(itide))%cname_tide)//'_z2', zti(:,:) )
225         !
226         DO ji=1,jpi
227            DO jj=1,jpj
228               amp_load(ji,jj,itide) =  SQRT( ztr(ji,jj)**2. + zti(ji,jj)**2. )
229               phi_load(ji,jj,itide) = ATAN2(-zti(ji,jj), ztr(ji,jj) )
230            END DO
231         END DO
232         !
233      END DO
234      CALL iom_close( inum )
235      !
236   END SUBROUTINE tide_init_load
237
238  !!======================================================================
239END MODULE sbctide
Note: See TracBrowser for help on using the repository browser.