source: NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/SBC/sbctide.F90 @ 12453

Last change on this file since 12453 was 12453, checked in by jcastill, 12 months ago

First implementation of the branch - compiling after merge

File size: 7.3 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   USE bdytides
19
20   IMPLICIT NONE
21   PUBLIC
22
23   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   pot_astro   !
24
25   !!----------------------------------------------------------------------
26   !!   tidal potential
27   !!----------------------------------------------------------------------
28   !!   sbc_tide            :
29   !!   tide_init_potential :
30   !!----------------------------------------------------------------------
31
32   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   amp_pot, phi_pot
33   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   amp_load, phi_load
34 
35   !!----------------------------------------------------------------------
36   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
37   !! $Id$
38   !! Software governed by the CeCILL license (see ./LICENSE)
39   !!----------------------------------------------------------------------
40CONTAINS
41
42   SUBROUTINE sbc_tide( kt )
43      !!----------------------------------------------------------------------
44      !!                 ***  ROUTINE sbc_tide  ***
45      !!----------------------------------------------------------------------     
46      INTEGER, INTENT( in ) ::   kt     ! ocean time-step
47      INTEGER               ::   jk     ! dummy loop index
48      INTEGER               ::   nsec_day_orig     ! Temporary variable
49      !!----------------------------------------------------------------------
50     
51      IF( nsec_day == NINT(0.5_wp * rdt) .OR. kt == nit000 ) THEN      ! start a new day
52         !
53         IF( kt == nit000 )THEN
54            ALLOCATE( amp_pot(jpi,jpj,nb_harmo),                      &
55               &      phi_pot(jpi,jpj,nb_harmo), pot_astro(jpi,jpj)   )
56            IF( ln_read_load )THEN
57               ALLOCATE( amp_load(jpi,jpj,nb_harmo), phi_load(jpi,jpj,nb_harmo) )
58               CALL tide_init_load
59            ENDIF
60         ENDIF
61         !
62         IF( ln_read_load )THEN
63            amp_pot(:,:,:) = amp_load(:,:,:)
64            phi_pot(:,:,:) = phi_load(:,:,:)
65         ELSE
66            amp_pot(:,:,:) = 0._wp
67            phi_pot(:,:,:) = 0._wp
68         ENDIF
69         pot_astro(:,:) = 0._wp
70         !
71         ! If the run does not start from midnight then need to initialise tides
72         ! at the start of the current day (only occurs when kt==nit000)
73         ! Temporarily set nsec_day to beginning of day.
74         nsec_day_orig = nsec_day
75         IF ( nsec_day /= NINT(0.5_wp * rdt) ) THEN
76            kt_tide = kt - (nsec_day - 0.5_wp * rdt)/rdt
77            nsec_day = NINT(0.5_wp * rdt)
78         ELSE
79            kt_tide = kt 
80         ENDIF
81         CALL tide_harmo( omega_tide, v0tide, utide, ftide, ntide, nb_harmo )
82         !
83         !
84         IF(lwp) THEN
85            WRITE(numout,*)
86            WRITE(numout,*) 'sbc_tide : Update of the components and (re)Init. the potential at kt=', kt
87            WRITE(numout,*) '~~~~~~~~ '
88            DO jk = 1, nb_harmo
89               WRITE(numout,*) Wave(ntide(jk))%cname_tide, utide(jk), ftide(jk), v0tide(jk), omega_tide(jk)
90            END DO
91         ENDIF
92         !
93         IF( ln_tide_pot )   CALL tide_init_potential
94         !
95         ! Reset nsec_day
96         nsec_day = nsec_day_orig 
97      ENDIF
98      !
99   END SUBROUTINE sbc_tide
100
101
102   SUBROUTINE tide_init_potential
103      !!----------------------------------------------------------------------
104      !!                 ***  ROUTINE tide_init_potential  ***
105      !!----------------------------------------------------------------------
106      INTEGER  ::   ji, jj, jk   ! dummy loop indices
107      REAL(wp) ::   zcons, ztmp1, ztmp2, zlat, zlon, ztmp, zamp, zcs   ! local scalar
108      !!----------------------------------------------------------------------
109
110      DO jk = 1, nb_harmo
111         ! love number now provides in tide namelist
112         zcons = dn_love_number * Wave(ntide(jk))%equitide * ftide(jk)
113         DO ji = 1, jpi
114            DO jj = 1, jpj
115               ztmp1 =  ftide(jk) * amp_pot(ji,jj,jk) * COS( phi_pot(ji,jj,jk) + v0tide(jk) + utide(jk) )
116               ztmp2 = -ftide(jk) * amp_pot(ji,jj,jk) * SIN( phi_pot(ji,jj,jk) + v0tide(jk) + utide(jk) )
117               zlat = gphit(ji,jj)*rad !! latitude en radian
118               zlon = glamt(ji,jj)*rad !! longitude en radian
119               ztmp = v0tide(jk) + utide(jk) + Wave(ntide(jk))%nutide * zlon
120               ! le potentiel est composé des effets des astres:
121               IF    ( Wave(ntide(jk))%nutide == 1 )  THEN  ;  zcs = zcons * SIN( 2._wp*zlat )
122               ELSEIF( Wave(ntide(jk))%nutide == 2 )  THEN  ;  zcs = zcons * COS( zlat )**2
123               ! Add tide potential for long period tides
124               ELSEIF( Wave(ntide(jk))%nutide == 0 )  THEN  ;  zcs = zcons * (0.5_wp-1.5_wp*SIN(zlat)**2._wp)
125               ELSE                                         ;  zcs = 0._wp
126               ENDIF
127               ztmp1 = ztmp1 + zcs * COS( ztmp )
128               ztmp2 = ztmp2 - zcs * SIN( ztmp )
129               zamp = SQRT( ztmp1*ztmp1 + ztmp2*ztmp2 )
130               amp_pot(ji,jj,jk) = zamp
131               phi_pot(ji,jj,jk) = ATAN2( -ztmp2 / MAX( 1.e-10_wp , zamp ) ,   &
132                  &                        ztmp1 / MAX( 1.e-10_wp,  zamp )   )
133            END DO
134         END DO
135      END DO
136      !
137   END SUBROUTINE tide_init_potential
138
139   SUBROUTINE tide_init_load
140      !!----------------------------------------------------------------------
141      !!                 ***  ROUTINE tide_init_load  ***
142      !!----------------------------------------------------------------------
143      INTEGER :: inum                 ! Logical unit of input file
144      INTEGER :: ji, jj, itide        ! dummy loop indices
145      REAL(wp), DIMENSION(jpi,jpj) ::   ztr, zti   !: workspace to read in tidal harmonics data
146      !!----------------------------------------------------------------------
147      IF(lwp) THEN
148         WRITE(numout,*)
149         WRITE(numout,*) 'tide_init_load : Initialization of load potential from file'
150         WRITE(numout,*) '~~~~~~~~~~~~~~ '
151      ENDIF
152      !
153      CALL iom_open ( cn_tide_load , inum )
154      !
155      DO itide = 1, nb_harmo
156         CALL iom_get  ( inum, jpdom_data,TRIM(Wave(ntide(itide))%cname_tide)//'_z1', ztr(:,:) )
157         CALL iom_get  ( inum, jpdom_data,TRIM(Wave(ntide(itide))%cname_tide)//'_z2', zti(:,:) )
158         !
159         DO ji=1,jpi
160            DO jj=1,jpj
161               amp_load(ji,jj,itide) =  SQRT( ztr(ji,jj)**2. + zti(ji,jj)**2. )
162               phi_load(ji,jj,itide) = ATAN2(-zti(ji,jj), ztr(ji,jj) )
163            END DO
164         END DO
165         !
166      END DO
167      CALL iom_close( inum )
168      !
169   END SUBROUTINE tide_init_load
170
171  !!======================================================================
172END MODULE sbctide
Note: See TracBrowser for help on using the repository browser.