source: branches/UKMO/dev_r5518_STOPACK/NEMOGCM/NEMO/OPA_SRC/TRA/trabbc.F90 @ 11384

Last change on this file since 11384 was 11384, checked in by mattmartin, 2 years ago

Included Andrea Storto's STOPACK code into NEMO3.6 branch.

  • Property svn:keywords set to Id
File size: 9.9 KB
Line 
1MODULE trabbc
2   !!==============================================================================
3   !!                       ***  MODULE  trabbc  ***
4   !! Ocean active tracers:  bottom boundary condition (geothermal heat flux)
5   !!==============================================================================
6   !! History :  OPA  ! 1999-10 (G. Madec)  original code
7   !!   NEMO     1.0  ! 2002-08 (G. Madec)  free form + modules
8   !!             -   ! 2002-11 (A. Bozec)  tra_bbc_init: original code
9   !!            3.3  ! 2010-10 (G. Madec)  dynamical allocation + suppression of key_trabbc
10   !!             -   ! 2010-11 (G. Madec)  use mbkt array (deepest ocean t-level)
11   !!----------------------------------------------------------------------
12
13   !!----------------------------------------------------------------------
14   !!   tra_bbc      : update the tracer trend at ocean bottom
15   !!   tra_bbc_init : initialization of geothermal heat flux trend
16   !!----------------------------------------------------------------------
17   USE oce             ! ocean variables
18   USE dom_oce         ! domain: ocean
19   USE phycst          ! physical constants
20   USE trd_oce         ! trends: ocean variables
21   USE trdtra          ! trends manager: tracers
22   USE in_out_manager  ! I/O manager
23   USE iom             ! I/O manager
24   USE fldread         ! read input fields
25   USE lbclnk            ! ocean lateral boundary conditions (or mpp link)
26   USE lib_mpp           ! distributed memory computing library
27   USE prtctl          ! Print control
28   USE wrk_nemo        ! Memory Allocation
29   USE timing          ! Timing
30   USE stopack
31
32   IMPLICIT NONE
33   PRIVATE
34
35   PUBLIC tra_bbc          ! routine called by step.F90
36   PUBLIC tra_bbc_init     ! routine called by opa.F90
37
38   !                                 !!* Namelist nambbc: bottom boundary condition *
39   LOGICAL, PUBLIC ::   ln_trabbc     !: Geothermal heat flux flag
40   INTEGER         ::   nn_geoflx     !  Geothermal flux (=1:constant flux, =2:read in file )
41   REAL(wp)        ::   rn_geoflx_cst !  Constant value of geothermal heat flux
42
43   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   qgh_trd0   ! geothermal heating trend
44   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   qgh_trd1   ! geothermal heating trend
45   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_qgh              ! structure of input qgh (file informations, fields read)
46 
47   !! * Substitutions
48#  include "domzgr_substitute.h90"
49   !!----------------------------------------------------------------------
50   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
51   !! $Id$
52   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
53   !!----------------------------------------------------------------------
54CONTAINS
55
56   SUBROUTINE tra_bbc( kt )
57      !!----------------------------------------------------------------------
58      !!                  ***  ROUTINE tra_bbc  ***
59      !!
60      !! ** Purpose :   Compute the bottom boundary contition on temperature
61      !!              associated with geothermal heating and add it to the
62      !!              general trend of temperature equations.
63      !!
64      !! ** Method  :   The geothermal heat flux set to its constant value of
65      !!              86.4 mW/m2 (Stein and Stein 1992, Huang 1999).
66      !!       The temperature trend associated to this heat flux through the
67      !!       ocean bottom can be computed once and is added to the temperature
68      !!       trend juste above the bottom at each time step:
69      !!            ta = ta + Qsf / (rau0 rcp e3T) for k= mbkt
70      !!       Where Qsf is the geothermal heat flux.
71      !!
72      !! ** Action  : - update the temperature trends (ta) with the trend of
73      !!                the ocean bottom boundary condition
74      !!
75      !! References : Stein, C. A., and S. Stein, 1992, Nature, 359, 123-129.
76      !!              Emile-Geay and Madec, 2009, Ocean Science.
77      !!----------------------------------------------------------------------
78      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
79      !!
80      INTEGER  ::   ji, jj, ik    ! dummy loop indices
81      REAL(wp) ::   zqgh_trd      ! geothermal heat flux trend
82      REAL(wp), POINTER, DIMENSION(:,:,:) ::   ztrdt
83      !!----------------------------------------------------------------------
84      !
85      IF( nn_timing == 1 )  CALL timing_start('tra_bbc')
86      !
87      IF( l_trdtra )   THEN         ! Save ta and sa trends
88         CALL wrk_alloc( jpi, jpj, jpk, ztrdt )
89         ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
90      ENDIF
91      !
92      !                             !  Add the geothermal heat flux trend on temperature
93
94      IF( nn_spp_geot .GT. 0) THEN
95          qgh_trd1 = qgh_trd0
96          CALL spp_gen(kt, qgh_trd0, nn_spp_geot, rn_geot_sd, jk_spp_geot)
97      ENDIF
98      DO jj = 2, jpjm1
99         DO ji = 2, jpim1
100            ik = mbkt(ji,jj)
101            zqgh_trd = qgh_trd1(ji,jj) / fse3t(ji,jj,ik)
102            tsa(ji,jj,ik,jp_tem) = tsa(ji,jj,ik,jp_tem) + zqgh_trd
103         END DO
104      END DO
105      !
106      CALL lbc_lnk( tsa(:,:,:,jp_tem) , 'T', 1. )
107      !
108      IF( l_trdtra ) THEN        ! Save the geothermal heat flux trend for diagnostics
109         ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
110         CALL trd_tra( kt, 'TRA', jp_tem, jptra_bbc, ztrdt )
111         CALL wrk_dealloc( jpi, jpj, jpk, ztrdt )
112      ENDIF
113      !
114      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbc  - Ta: ', mask1=tmask, clinfo3='tra-ta' )
115      !
116      IF( nn_timing == 1 )  CALL timing_stop('tra_bbc')
117      !
118   END SUBROUTINE tra_bbc
119
120
121   SUBROUTINE tra_bbc_init
122      !!----------------------------------------------------------------------
123      !!                  ***  ROUTINE tra_bbc_init  ***
124      !!
125      !! ** Purpose :   Compute once for all the trend associated with geothermal
126      !!              heating that will be applied at each time step at the
127      !!              last ocean level
128      !!
129      !! ** Method  :   Read the nambbc namelist and check the parameters.
130      !!
131      !! ** Input   : - Namlist nambbc
132      !!              - NetCDF file  : geothermal_heating.nc ( if necessary )
133      !!
134      !! ** Action  : - read/fix the geothermal heat qgh_trd0
135      !!----------------------------------------------------------------------
136      USE iom
137      !!
138      INTEGER  ::   ji, jj              ! dummy loop indices
139      INTEGER  ::   inum                ! temporary logical unit
140      INTEGER  ::   ios                 ! Local integer output status for namelist read
141      INTEGER  ::   ierror              ! local integer
142      !
143      TYPE(FLD_N)        ::   sn_qgh    ! informations about the geotherm. field to be read
144      CHARACTER(len=256) ::   cn_dir    ! Root directory for location of ssr files
145      !
146      NAMELIST/nambbc/ln_trabbc, nn_geoflx, rn_geoflx_cst, sn_qgh, cn_dir 
147      !!----------------------------------------------------------------------
148
149      REWIND( numnam_ref )              ! Namelist nambbc in reference namelist : Bottom momentum boundary condition
150      READ  ( numnam_ref, nambbc, IOSTAT = ios, ERR = 901)
151901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambbc in reference namelist', lwp )
152
153      REWIND( numnam_cfg )              ! Namelist nambbc in configuration namelist : Bottom momentum boundary condition
154      READ  ( numnam_cfg, nambbc, IOSTAT = ios, ERR = 902 )
155902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambbc in configuration namelist', lwp )
156      IF(lwm) WRITE ( numond, nambbc )
157
158      IF(lwp) THEN                     ! Control print
159         WRITE(numout,*)
160         WRITE(numout,*) 'tra_bbc : Bottom Boundary Condition (bbc), apply a Geothermal heating'
161         WRITE(numout,*) '~~~~~~~   '
162         WRITE(numout,*) '   Namelist nambbc : set bbc parameters'
163         WRITE(numout,*) '      Apply a geothermal heating at ocean bottom   ln_trabbc     = ', ln_trabbc
164         WRITE(numout,*) '      type of geothermal flux                      nn_geoflx     = ', nn_geoflx
165         WRITE(numout,*) '      Constant geothermal flux value               rn_geoflx_cst = ', rn_geoflx_cst
166         WRITE(numout,*)
167      ENDIF
168
169      IF( ln_trabbc ) THEN             !==  geothermal heating  ==!
170         !
171         ALLOCATE( qgh_trd0(jpi,jpj) )    ! allocation
172         ALLOCATE( qgh_trd1(jpi,jpj) )    ! allocation
173         !
174         SELECT CASE ( nn_geoflx )        ! geothermal heat flux / (rauO * Cp)
175         !
176         CASE ( 1 )                          !* constant flux
177            IF(lwp) WRITE(numout,*) '      *** constant heat flux  =   ', rn_geoflx_cst
178            qgh_trd0(:,:) = r1_rau0_rcp * rn_geoflx_cst
179            !
180         CASE ( 2 )                          !* variable geothermal heat flux : read the geothermal fluxes in mW/m2
181            IF(lwp) WRITE(numout,*) '      *** variable geothermal heat flux'
182            !
183            ALLOCATE( sf_qgh(1), STAT=ierror )
184            IF( ierror > 0 ) THEN
185               CALL ctl_stop( 'tra_bbc_init: unable to allocate sf_qgh structure' )   ;
186               RETURN
187            ENDIF
188            ALLOCATE( sf_qgh(1)%fnow(jpi,jpj,1)   )
189            IF( sn_qgh%ln_tint )ALLOCATE( sf_qgh(1)%fdta(jpi,jpj,1,2) )
190            ! fill sf_chl with sn_chl and control print
191            CALL fld_fill( sf_qgh, (/ sn_qgh /), cn_dir, 'tra_bbc_init',   &
192               &          'bottom temperature boundary condition', 'nambbc' )
193
194            CALL fld_read( nit000, 1, sf_qgh )                         ! Read qgh data
195            qgh_trd0(:,:) = r1_rau0_rcp * sf_qgh(1)%fnow(:,:,1) * 1.e-3 ! conversion in W/m2
196            !
197         CASE DEFAULT
198            WRITE(ctmp1,*) '     bad flag value for nn_geoflx = ', nn_geoflx
199            CALL ctl_stop( ctmp1 )
200            !
201         END SELECT
202         qgh_trd1 = qgh_trd0
203         !
204      ELSE
205         IF(lwp) WRITE(numout,*) '      *** no geothermal heat flux'
206      ENDIF
207      !
208   END SUBROUTINE tra_bbc_init
209
210   !!======================================================================
211END MODULE trabbc
Note: See TracBrowser for help on using the repository browser.