source: branches/UKMO/dev_1d_bugfixes_tocommit/NEMOGCM/NEMO/OPA_SRC/TRA/trabbc.F90 @ 13191

Last change on this file since 13191 was 13191, checked in by jwhile, 7 months ago

Updates for 1d runnig

File size: 10.1 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 defined key_traldf_c2d || key_traldf_c3d
95      IF( ln_stopack .AND. nn_spp_geot > 0) THEN
96          qgh_trd1(:,:) = qgh_trd0(:,:)
97          CALL spp_gen(kt, qgh_trd1, nn_spp_geot, rn_geot_sd, jk_spp_geot)
98      ENDIF
99#else
100      IF ( ln_stopack .AND. nn_spp_geot > 0 ) &
101         & CALL ctl_stop( 'tra_bbc: parameter perturbation will only work with '// &
102                          'key_traldf_c2d or key_traldf_c3d')
103#endif
104
105
106      DO jj = 2, jpjm1
107         DO ji = 2, jpim1
108            ik = mbkt(ji,jj)
109            zqgh_trd = qgh_trd1(ji,jj) / fse3t(ji,jj,ik)
110            tsa(ji,jj,ik,jp_tem) = tsa(ji,jj,ik,jp_tem) + zqgh_trd
111         END DO
112      END DO
113      !
114      CALL lbc_lnk( tsa(:,:,:,jp_tem) , 'T', 1. )
115      !
116      IF( l_trdtra ) THEN        ! Save the geothermal heat flux trend for diagnostics
117         ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
118         CALL trd_tra( kt, 'TRA', jp_tem, jptra_bbc, ztrdt )
119         CALL wrk_dealloc( jpi, jpj, jpk, ztrdt )
120      ENDIF
121      !
122      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbc  - Ta: ', mask1=tmask, clinfo3='tra-ta' )
123      !
124      IF( nn_timing == 1 )  CALL timing_stop('tra_bbc')
125      !
126   END SUBROUTINE tra_bbc
127
128
129   SUBROUTINE tra_bbc_init
130      !!----------------------------------------------------------------------
131      !!                  ***  ROUTINE tra_bbc_init  ***
132      !!
133      !! ** Purpose :   Compute once for all the trend associated with geothermal
134      !!              heating that will be applied at each time step at the
135      !!              last ocean level
136      !!
137      !! ** Method  :   Read the nambbc namelist and check the parameters.
138      !!
139      !! ** Input   : - Namlist nambbc
140      !!              - NetCDF file  : geothermal_heating.nc ( if necessary )
141      !!
142      !! ** Action  : - read/fix the geothermal heat qgh_trd0
143      !!----------------------------------------------------------------------
144      USE iom
145      !!
146      INTEGER  ::   ji, jj              ! dummy loop indices
147      INTEGER  ::   inum                ! temporary logical unit
148      INTEGER  ::   ios                 ! Local integer output status for namelist read
149      INTEGER  ::   ierror              ! local integer
150      !
151      TYPE(FLD_N)        ::   sn_qgh    ! informations about the geotherm. field to be read
152      CHARACTER(len=256) ::   cn_dir    ! Root directory for location of ssr files
153      !
154      NAMELIST/nambbc/ln_trabbc, nn_geoflx, rn_geoflx_cst, sn_qgh, cn_dir
155      !!----------------------------------------------------------------------
156
157      REWIND( numnam_ref )              ! Namelist nambbc in reference namelist : Bottom momentum boundary condition
158      READ  ( numnam_ref, nambbc, IOSTAT = ios, ERR = 901)
159901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambbc in reference namelist', lwp )
160
161      REWIND( numnam_cfg )              ! Namelist nambbc in configuration namelist : Bottom momentum boundary condition
162      READ  ( numnam_cfg, nambbc, IOSTAT = ios, ERR = 902 )
163902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambbc in configuration namelist', lwp )
164      IF(lwm) WRITE ( numond, nambbc )
165
166      IF(lwp) THEN                     ! Control print
167         WRITE(numout,*)
168         WRITE(numout,*) 'tra_bbc : Bottom Boundary Condition (bbc), apply a Geothermal heating'
169         WRITE(numout,*) '~~~~~~~   '
170         WRITE(numout,*) '   Namelist nambbc : set bbc parameters'
171         WRITE(numout,*) '      Apply a geothermal heating at ocean bottom   ln_trabbc     = ', ln_trabbc
172         WRITE(numout,*) '      type of geothermal flux                      nn_geoflx     = ', nn_geoflx
173         WRITE(numout,*) '      Constant geothermal flux value               rn_geoflx_cst = ', rn_geoflx_cst
174         WRITE(numout,*)
175      ENDIF
176
177      IF( ln_trabbc ) THEN             !==  geothermal heating  ==!
178         !
179         ALLOCATE( qgh_trd0(jpi,jpj) )    ! allocation
180         ALLOCATE( qgh_trd1(jpi,jpj) )    ! allocation
181         !
182         SELECT CASE ( nn_geoflx )        ! geothermal heat flux / (rauO * Cp)
183         !
184         CASE ( 1 )                          !* constant flux
185            IF(lwp) WRITE(numout,*) '      *** constant heat flux  =   ', rn_geoflx_cst
186            qgh_trd0(:,:) = r1_rau0_rcp * rn_geoflx_cst
187            !
188         CASE ( 2 )                          !* variable geothermal heat flux : read the geothermal fluxes in mW/m2
189            IF(lwp) WRITE(numout,*) '      *** variable geothermal heat flux'
190            !
191            ALLOCATE( sf_qgh(1), STAT=ierror )
192            IF( ierror > 0 ) THEN
193               CALL ctl_stop( 'tra_bbc_init: unable to allocate sf_qgh structure' )   ;
194               RETURN
195            ENDIF
196            ALLOCATE( sf_qgh(1)%fnow(jpi,jpj,1)   )
197            IF( sn_qgh%ln_tint )ALLOCATE( sf_qgh(1)%fdta(jpi,jpj,1,2) )
198            ! fill sf_chl with sn_chl and control print
199            CALL fld_fill( sf_qgh, (/ sn_qgh /), cn_dir, 'tra_bbc_init',   &
200               &          'bottom temperature boundary condition', 'nambbc' )
201
202            CALL fld_read( nit000, 1, sf_qgh )                         ! Read qgh data
203            qgh_trd0(:,:) = r1_rau0_rcp * sf_qgh(1)%fnow(:,:,1) * 1.e-3 ! conversion in W/m2
204            !
205         CASE DEFAULT
206            WRITE(ctmp1,*) '     bad flag value for nn_geoflx = ', nn_geoflx
207            CALL ctl_stop( ctmp1 )
208            !
209         END SELECT
210         qgh_trd1(:,:) = qgh_trd0(:,:)
211         !
212      ELSE
213         IF(lwp) WRITE(numout,*) '      *** no geothermal heat flux'
214      ENDIF
215      !
216   END SUBROUTINE tra_bbc_init
217
218   !!======================================================================
219END MODULE trabbc
Note: See TracBrowser for help on using the repository browser.