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.
usrdef_sbc.F90 in NEMO/trunk/src/OCE/USR – NEMO

source: NEMO/trunk/src/OCE/USR/usrdef_sbc.F90 @ 12377

Last change on this file since 12377 was 12377, checked in by acc, 4 years ago

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge --ignore-ancestry \

svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The --ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.

  • Property svn:keywords set to Id
File size: 11.7 KB
Line 
1MODULE usrdef_sbc
2   !!======================================================================
3   !!                     ***  MODULE  usrdef_sbc  ***
4   !!
5   !!                     ===  GYRE configuration  ===
6   !!
7   !! User defined :   surface forcing of a user configuration
8   !!======================================================================
9   !! History :  4.0   ! 2016-03  (S. Flavoni, G. Madec)  user defined interface
10   !!----------------------------------------------------------------------
11
12   !!----------------------------------------------------------------------
13   !!   usrdef_sbc    : user defined surface bounday conditions in GYRE case
14   !!----------------------------------------------------------------------
15   USE oce            ! ocean dynamics and tracers
16   USE dom_oce        ! ocean space and time domain
17   USE sbc_oce        ! Surface boundary condition: ocean fields
18   USE phycst         ! physical constants
19   !
20   USE in_out_manager ! I/O manager
21   USE lib_mpp        ! distribued memory computing library
22   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
23   USE lib_fortran    !
24
25   IMPLICIT NONE
26   PRIVATE
27
28   PUBLIC   usrdef_sbc_oce       ! routine called in sbcmod module
29   PUBLIC   usrdef_sbc_ice_tau   ! routine called by icestp.F90 for ice dynamics
30   PUBLIC   usrdef_sbc_ice_flx   ! routine called by icestp.F90 for ice thermo
31
32   !! * Substitutions
33#  include "do_loop_substitute.h90"
34   !!----------------------------------------------------------------------
35   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
36   !! $Id$
37   !! Software governed by the CeCILL license (see ./LICENSE)
38   !!----------------------------------------------------------------------
39CONTAINS
40
41   SUBROUTINE usrdef_sbc_oce( kt, Kbb )
42      !!---------------------------------------------------------------------
43      !!                    ***  ROUTINE usrdef_sbc  ***
44      !!             
45      !! ** Purpose :   provide at each time-step the GYRE surface boundary
46      !!              condition, i.e. the momentum, heat and freshwater fluxes.
47      !!
48      !! ** Method  :   analytical seasonal cycle for GYRE configuration.
49      !!                CAUTION : never mask the surface stress field !
50      !!
51      !! ** Action  : - set the ocean surface boundary condition, i.e.   
52      !!                   utau, vtau, taum, wndm, qns, qsr, emp, sfx
53      !!
54      !! Reference : Hazeleger, W., and S. Drijfhout, JPO, 30, 677-695, 2000.
55      !!----------------------------------------------------------------------
56      INTEGER, INTENT(in) ::   kt   ! ocean time step
57      INTEGER, INTENT(in) ::   Kbb  ! ocean time index
58      !!
59      INTEGER  ::   ji, jj                 ! dummy loop indices
60      INTEGER  ::   zyear0                 ! initial year
61      INTEGER  ::   zmonth0                ! initial month
62      INTEGER  ::   zday0                  ! initial day
63      INTEGER  ::   zday_year0             ! initial day since january 1st
64      REAL(wp) ::   ztau     , ztau_sais   ! wind intensity and of the seasonal cycle
65      REAL(wp) ::   ztime                  ! time in hour
66      REAL(wp) ::   ztimemax , ztimemin    ! 21th June, and 21th decem. if date0 = 1st january
67      REAL(wp) ::   ztimemax1, ztimemin1   ! 21th June, and 21th decem. if date0 = 1st january
68      REAL(wp) ::   ztimemax2, ztimemin2   ! 21th June, and 21th decem. if date0 = 1st january
69      REAL(wp) ::   ztaun                  ! intensity
70      REAL(wp) ::   zemp_s, zemp_n, zemp_sais, ztstar
71      REAL(wp) ::   zcos_sais1, zcos_sais2, ztrp, zconv, t_star
72      REAL(wp) ::   zsumemp, zsurf
73      REAL(wp) ::   zrhoa  = 1.22         ! Air density kg/m3
74      REAL(wp) ::   zcdrag = 1.5e-3       ! drag coefficient
75      REAL(wp) ::   ztx, zty, zmod, zcoef ! temporary variables
76      REAL(wp) ::   zyydd                 ! number of days in one year
77      !!---------------------------------------------------------------------
78      zyydd = REAL(nyear_len(1),wp)
79
80      ! ---------------------------- !
81      !  heat and freshwater fluxes  !
82      ! ---------------------------- !
83      !same temperature, E-P as in HAZELEGER 2000
84
85      zyear0     =   ndate0 / 10000                             ! initial year
86      zmonth0    = ( ndate0 - zyear0 * 10000 ) / 100            ! initial month
87      zday0      =   ndate0 - zyear0 * 10000 - zmonth0 * 100    ! initial day betwen 1 and 30
88      zday_year0 = ( zmonth0 - 1 ) * 30.+zday0                  ! initial day betwen 1 and 360
89
90      ! current day (in hours) since january the 1st of the current year
91      ztime = REAL( kt ) * rdt / (rmmss * rhhmm)   &       !  total incrementation (in hours)
92         &      - (nyear  - 1) * rjjhh * zyydd             !  minus years since beginning of experiment (in hours)
93
94      ztimemax1 = ((5.*30.)+21.)* 24.                      ! 21th june     at 24h in hours
95      ztimemin1 = ztimemax1 + rjjhh * zyydd / 2            ! 21th december        in hours
96      ztimemax2 = ((6.*30.)+21.)* 24.                      ! 21th july     at 24h in hours
97      ztimemin2 = ztimemax2 - rjjhh * zyydd / 2            ! 21th january         in hours
98      !                                                    ! NB: rjjhh * zyydd / 4 = one seasonal cycle in hours
99
100      ! amplitudes
101      zemp_S    = 0.7       ! intensity of COS in the South
102      zemp_N    = 0.8       ! intensity of COS in the North
103      zemp_sais = 0.1
104      zTstar    = 28.3      ! intemsity from 28.3 a -5 deg
105
106      ! 1/2 period between 21th June and 21th December and between 21th July and 21th January
107      zcos_sais1 = COS( (ztime - ztimemax1) / (ztimemin1 - ztimemax1) * rpi ) 
108      zcos_sais2 = COS( (ztime - ztimemax2) / (ztimemax2 - ztimemin2) * rpi )
109
110      ztrp= - 40.e0        ! retroaction term on heat fluxes (W/m2/K)
111      zconv = 3.16e-5      ! convertion factor: 1 m/yr => 3.16e-5 mm/s
112      DO_2D_11_11
113         ! domain from 15 deg to 50 deg between 27 and 28  degC at 15N, -3
114         ! and 13 degC at 50N 53.5 + or - 11 = 1/4 period :
115         ! 64.5 in summer, 42.5 in winter
116         t_star = zTstar * ( 1. + 1. / 50. * zcos_sais2 )                &
117            &                    * COS( rpi * (gphit(ji,jj) - 5.)               &
118            &                    / ( 53.5 * ( 1 + 11 / 53.5 * zcos_sais2 ) * 2.) )
119         ! 23.5 deg : tropics
120         qsr (ji,jj) =  230 * COS( 3.1415 * ( gphit(ji,jj) - 23.5 * zcos_sais1 ) / ( 0.9 * 180 ) )
121         qns (ji,jj) = ztrp * ( ts(ji,jj,1,jp_tem,Kbb) - t_star ) - qsr(ji,jj)
122         IF( gphit(ji,jj) >= 14.845 .AND. 37.2 >= gphit(ji,jj) ) THEN    ! zero at 37.8 deg, max at 24.6 deg
123            emp  (ji,jj) =   zemp_S * zconv   &
124               &         * SIN( rpi / 2 * (gphit(ji,jj) - 37.2) / (24.6 - 37.2) )  &
125               &         * ( 1 - zemp_sais / zemp_S * zcos_sais1)
126         ELSE
127            emp (ji,jj) =  - zemp_N * zconv   &
128               &         * SIN( rpi / 2 * (gphit(ji,jj) - 37.2) / (46.8 - 37.2) )  &
129               &         * ( 1 - zemp_sais / zemp_N * zcos_sais1 )
130         ENDIF
131      END_2D
132
133      zsumemp = GLOB_SUM( 'usrdef_sbc', emp  (:,:)   ) 
134      zsurf   = GLOB_SUM( 'usrdef_sbc', tmask(:,:,1) ) 
135      zsumemp = zsumemp / zsurf         ! Default GYRE configuration
136
137      ! freshwater (mass flux) and update of qns with heat content of emp
138      emp (:,:) = emp(:,:) - zsumemp * tmask(:,:,1)        ! freshwater flux (=0 in domain average)
139      sfx (:,:) = 0.0_wp                                   ! no salt flux
140      qns (:,:) = qns(:,:) - emp(:,:) * sst_m(:,:) * rcp   ! evap and precip are at SST
141
142
143      ! ---------------------------- !
144      !       momentum fluxes        !
145      ! ---------------------------- !
146      ! same wind as in Wico
147      !test date0 : ndate0 = 010203
148      zyear0  =   ndate0 / 10000
149      zmonth0 = ( ndate0 - zyear0 * 10000 ) / 100
150      zday0   =   ndate0 - zyear0 * 10000 - zmonth0 * 100
151      !Calculates nday_year, day since january 1st
152      zday_year0 = (zmonth0-1)*30.+zday0
153
154      !accumulates days of previous months of this year
155      ! day (in hours) since january the 1st
156      ztime = FLOAT( kt ) * rdt / (rmmss * rhhmm)  &  ! incrementation in hour
157         &     - (nyear - 1) * rjjhh * zyydd          !  - nber of hours the precedent years
158      ztimemax = ((5.*30.)+21.)* 24.               ! 21th june     in hours
159      ztimemin = ztimemax + rjjhh * zyydd / 2      ! 21th december in hours
160      !                                            ! NB: rjjhh * zyydd / 4 = 1 seasonal cycle in hours
161
162      ! mean intensity at 0.105 ; srqt(2) because projected with 45deg angle
163      ztau = 0.105 / SQRT( 2. )
164      ! seasonal oscillation intensity
165      ztau_sais = 0.015
166      ztaun = ztau - ztau_sais * COS( (ztime - ztimemax) / (ztimemin - ztimemax) * rpi )
167      DO_2D_11_11
168        ! domain from 15deg to 50deg and 1/2 period along 14deg
169        ! so 5/4 of half period with seasonal cycle
170        utau(ji,jj) = - ztaun * SIN( rpi * (gphiu(ji,jj) - 15.) / (29.-15.) )
171        vtau(ji,jj) =   ztaun * SIN( rpi * (gphiv(ji,jj) - 15.) / (29.-15.) )
172      END_2D
173
174      ! module of wind stress and wind speed at T-point
175      zcoef = 1. / ( zrhoa * zcdrag ) 
176      DO_2D_00_00
177         ztx = utau(ji-1,jj  ) + utau(ji,jj) 
178         zty = vtau(ji  ,jj-1) + vtau(ji,jj) 
179         zmod = 0.5 * SQRT( ztx * ztx + zty * zty )
180         taum(ji,jj) = zmod
181         wndm(ji,jj) = SQRT( zmod * zcoef )
182      END_2D
183      CALL lbc_lnk_multi( 'usrdef_sbc', taum(:,:), 'T', 1. , wndm(:,:), 'T', 1. )
184
185      ! ---------------------------------- !
186      !  control print at first time-step  !
187      ! ---------------------------------- !
188      IF( kt == nit000 .AND. lwp ) THEN
189         WRITE(numout,*)
190         WRITE(numout,*)'usrdef_sbc_oce : analytical surface fluxes for GYRE configuration'               
191         WRITE(numout,*)'~~~~~~~~~~~ ' 
192         WRITE(numout,*)'           nyear      = ', nyear
193         WRITE(numout,*)'           nmonth     = ', nmonth
194         WRITE(numout,*)'           nday       = ', nday
195         WRITE(numout,*)'           nday_year  = ', nday_year
196         WRITE(numout,*)'           ztime      = ', ztime
197         WRITE(numout,*)'           ztimemax   = ', ztimemax
198         WRITE(numout,*)'           ztimemin   = ', ztimemin
199         WRITE(numout,*)'           ztimemax1  = ', ztimemax1
200         WRITE(numout,*)'           ztimemin1  = ', ztimemin1
201         WRITE(numout,*)'           ztimemax2  = ', ztimemax2
202         WRITE(numout,*)'           ztimemin2  = ', ztimemin2
203         WRITE(numout,*)'           zyear0     = ', zyear0
204         WRITE(numout,*)'           zmonth0    = ', zmonth0
205         WRITE(numout,*)'           zday0      = ', zday0
206         WRITE(numout,*)'           zday_year0 = ', zday_year0
207         WRITE(numout,*)'           zyydd      = ', zyydd
208         WRITE(numout,*)'           zemp_S     = ', zemp_S
209         WRITE(numout,*)'           zemp_N     = ', zemp_N
210         WRITE(numout,*)'           zemp_sais  = ', zemp_sais
211         WRITE(numout,*)'           zTstar     = ', zTstar
212         WRITE(numout,*)'           zsumemp    = ', zsumemp
213         WRITE(numout,*)'           zsurf      = ', zsurf
214         WRITE(numout,*)'           ztrp       = ', ztrp
215         WRITE(numout,*)'           zconv      = ', zconv
216         WRITE(numout,*)'           ndastp     = ', ndastp
217         WRITE(numout,*)'           adatrj     = ', adatrj
218      ENDIF
219      !
220   END SUBROUTINE usrdef_sbc_oce
221
222
223   SUBROUTINE usrdef_sbc_ice_tau( kt )
224      INTEGER, INTENT(in) ::   kt   ! ocean time step
225   END SUBROUTINE usrdef_sbc_ice_tau
226
227
228   SUBROUTINE usrdef_sbc_ice_flx( kt, phs, phi )
229      INTEGER, INTENT(in) ::   kt   ! ocean time step
230      REAL(wp), DIMENSION(:,:,:), INTENT(in)  ::   phs    ! snow thickness
231      REAL(wp), DIMENSION(:,:,:), INTENT(in)  ::   phi    ! ice thickness
232   END SUBROUTINE usrdef_sbc_ice_flx
233
234   !!======================================================================
235END MODULE usrdef_sbc
Note: See TracBrowser for help on using the repository browser.