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.
trabbc.F90 in trunk/NEMOGCM/NEMO/OPA_SRC/TRA – NEMO

source: trunk/NEMOGCM/NEMO/OPA_SRC/TRA/trabbc.F90 @ 7698

Last change on this file since 7698 was 7698, checked in by mocavero, 7 years ago

update trunk with OpenMP parallelization
  • Property svn:keywords set to Id
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   !
23   USE in_out_manager ! I/O manager
24   USE iom            ! xIOS
25   USE fldread        ! read input fields
26   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
27   USE lib_mpp        ! distributed memory computing library
28   USE prtctl         ! Print control
29   USE wrk_nemo       ! Memory Allocation
30   USE timing         ! Timing
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 , ALLOCATABLE, DIMENSION(:,:) ::   qgh_trd0   ! geothermal heating trend
44
45   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_qgh   ! structure of input qgh (file informations, fields read)
46 
47   !!----------------------------------------------------------------------
48   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
49   !! $Id$
50   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
51   !!----------------------------------------------------------------------
52CONTAINS
53
54   SUBROUTINE tra_bbc( kt )
55      !!----------------------------------------------------------------------
56      !!                  ***  ROUTINE tra_bbc  ***
57      !!
58      !! ** Purpose :   Compute the bottom boundary contition on temperature
59      !!              associated with geothermal heating and add it to the
60      !!              general trend of temperature equations.
61      !!
62      !! ** Method  :   The geothermal heat flux set to its constant value of
63      !!              86.4 mW/m2 (Stein and Stein 1992, Huang 1999).
64      !!       The temperature trend associated to this heat flux through the
65      !!       ocean bottom can be computed once and is added to the temperature
66      !!       trend juste above the bottom at each time step:
67      !!            ta = ta + Qsf / (rau0 rcp e3T) for k= mbkt
68      !!       Where Qsf is the geothermal heat flux.
69      !!
70      !! ** Action  : - update the temperature trends with geothermal heating trend
71      !!              - send the trend for further diagnostics (ln_trdtra=T)
72      !!
73      !! References : Stein, C. A., and S. Stein, 1992, Nature, 359, 123-129.
74      !!              Emile-Geay and Madec, 2009, Ocean Science.
75      !!----------------------------------------------------------------------
76      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
77      !
78      INTEGER  ::   ji, jj, jk    ! dummy loop indices
79      REAL(wp), POINTER, DIMENSION(:,:,:) ::   ztrdt
80      !!----------------------------------------------------------------------
81      !
82      IF( nn_timing == 1 )  CALL timing_start('tra_bbc')
83      !
84      IF( l_trdtra )   THEN         ! Save the input temperature trend
85         CALL wrk_alloc( jpi,jpj,jpk,   ztrdt )
86!$OMP PARALLEL DO schedule(static) private(jk, jj, ji)
87         DO jk = 1, jpk
88            DO jj = 1, jpj
89               DO ji = 1, jpi
90                  ztrdt(ji,jj,jk) = tsa(ji,jj,jk,jp_tem)
91               END DO
92            END DO
93         END DO
94      ENDIF
95      !                             !  Add the geothermal trend on temperature
96!$OMP PARALLEL DO schedule(static) private(jj, ji)
97      DO jj = 2, jpjm1
98         DO ji = 2, jpim1
99            tsa(ji,jj,mbkt(ji,jj),jp_tem) = tsa(ji,jj,mbkt(ji,jj),jp_tem) + qgh_trd0(ji,jj) / e3t_n(ji,jj,mbkt(ji,jj))
100         END DO
101      END DO
102      !
103      CALL lbc_lnk( tsa(:,:,:,jp_tem) , 'T', 1. )
104      !
105      IF( l_trdtra ) THEN        ! Send the trend for diagnostics
106!$OMP PARALLEL DO schedule(static) private(jk, jj, ji)
107         DO jk = 1, jpk
108            DO jj = 1, jpj
109               DO ji = 1, jpi
110                  ztrdt(ji,jj,jk) = tsa(ji,jj,jk,jp_tem) - ztrdt(ji,jj,jk)
111               END DO
112            END DO
113         END DO
114         CALL trd_tra( kt, 'TRA', jp_tem, jptra_bbc, ztrdt )
115         CALL wrk_dealloc( jpi,jpj,jpk,   ztrdt )
116      ENDIF
117      !
118      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbc  - Ta: ', mask1=tmask, clinfo3='tra-ta' )
119      !
120      IF( nn_timing == 1 )  CALL timing_stop('tra_bbc')
121      !
122   END SUBROUTINE tra_bbc
123
124
125   SUBROUTINE tra_bbc_init
126      !!----------------------------------------------------------------------
127      !!                  ***  ROUTINE tra_bbc_init  ***
128      !!
129      !! ** Purpose :   Compute once for all the trend associated with geothermal
130      !!              heating that will be applied at each time step at the
131      !!              last ocean level
132      !!
133      !! ** Method  :   Read the nambbc namelist and check the parameters.
134      !!
135      !! ** Input   : - Namlist nambbc
136      !!              - NetCDF file  : geothermal_heating.nc ( if necessary )
137      !!
138      !! ** Action  : - read/fix the geothermal heat qgh_trd0
139      !!----------------------------------------------------------------------
140      INTEGER  ::   ji, jj              ! dummy loop indices
141      INTEGER  ::   inum                ! temporary logical unit
142      INTEGER  ::   ios                 ! Local integer output status for namelist read
143      INTEGER  ::   ierror              ! local integer
144      !
145      TYPE(FLD_N)        ::   sn_qgh    ! informations about the geotherm. field to be read
146      CHARACTER(len=256) ::   cn_dir    ! Root directory for location of ssr files
147      !
148      NAMELIST/nambbc/ln_trabbc, nn_geoflx, rn_geoflx_cst, sn_qgh, cn_dir 
149      !!----------------------------------------------------------------------
150      !
151      REWIND( numnam_ref )              ! Namelist nambbc in reference namelist : Bottom momentum boundary condition
152      READ  ( numnam_ref, nambbc, IOSTAT = ios, ERR = 901)
153901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambbc in reference namelist', lwp )
154      !
155      REWIND( numnam_cfg )              ! Namelist nambbc in configuration namelist : Bottom momentum boundary condition
156      READ  ( numnam_cfg, nambbc, IOSTAT = ios, ERR = 902 )
157902   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambbc in configuration namelist', lwp )
158      IF(lwm) WRITE ( numond, nambbc )
159      !
160      IF(lwp) THEN                     ! Control print
161         WRITE(numout,*)
162         WRITE(numout,*) 'tra_bbc : Bottom Boundary Condition (bbc), apply a Geothermal heating'
163         WRITE(numout,*) '~~~~~~~   '
164         WRITE(numout,*) '   Namelist nambbc : set bbc parameters'
165         WRITE(numout,*) '      Apply a geothermal heating at ocean bottom   ln_trabbc     = ', ln_trabbc
166         WRITE(numout,*) '      type of geothermal flux                      nn_geoflx     = ', nn_geoflx
167         WRITE(numout,*) '      Constant geothermal flux value               rn_geoflx_cst = ', rn_geoflx_cst
168         WRITE(numout,*)
169      ENDIF
170      !
171      IF( ln_trabbc ) THEN             !==  geothermal heating  ==!
172         !
173         ALLOCATE( qgh_trd0(jpi,jpj) )    ! allocation
174         !
175         SELECT CASE ( nn_geoflx )        ! geothermal heat flux / (rauO * Cp)
176         !
177         CASE ( 1 )                          !* constant flux
178            IF(lwp) WRITE(numout,*) '      *** constant heat flux  =   ', rn_geoflx_cst
179!$OMP PARALLEL DO schedule(static) private(jj, ji)
180            DO jj = 1, jpj
181               DO ji = 1, jpi
182                  qgh_trd0(ji,jj) = r1_rau0_rcp * rn_geoflx_cst
183               END DO
184            END DO
185            !
186         CASE ( 2 )                          !* variable geothermal heat flux : read the geothermal fluxes in mW/m2
187            IF(lwp) WRITE(numout,*) '      *** variable geothermal heat flux'
188            !
189            ALLOCATE( sf_qgh(1), STAT=ierror )
190            IF( ierror > 0 ) THEN
191               CALL ctl_stop( 'tra_bbc_init: unable to allocate sf_qgh structure' )   ;
192               RETURN
193            ENDIF
194            ALLOCATE( sf_qgh(1)%fnow(jpi,jpj,1)   )
195            IF( sn_qgh%ln_tint )ALLOCATE( sf_qgh(1)%fdta(jpi,jpj,1,2) )
196            ! fill sf_chl with sn_chl and control print
197            CALL fld_fill( sf_qgh, (/ sn_qgh /), cn_dir, 'tra_bbc_init',   &
198               &          'bottom temperature boundary condition', 'nambbc', no_print )
199
200            CALL fld_read( nit000, 1, sf_qgh )                         ! Read qgh data
201!$OMP PARALLEL DO schedule(static) private(jj, ji)
202            DO jj = 1, jpj
203               DO ji = 1, jpi
204                  qgh_trd0(ji,jj) = r1_rau0_rcp * sf_qgh(1)%fnow(ji,jj,1) * 1.e-3 ! conversion in W/m2
205               END DO
206            END DO
207            !
208         CASE DEFAULT
209            WRITE(ctmp1,*) '     bad flag value for nn_geoflx = ', nn_geoflx
210            CALL ctl_stop( ctmp1 )
211         END SELECT
212         !
213      ELSE
214         IF(lwp) WRITE(numout,*) '      *** no geothermal heat flux'
215      ENDIF
216      !
217   END SUBROUTINE tra_bbc_init
218
219   !!======================================================================
220END MODULE trabbc
Note: See TracBrowser for help on using the repository browser.