Changeset 2443 for branches/nemo_v3_3_beta
- Timestamp:
- 2010-11-29T08:52:36+01:00 (13 years ago)
- Location:
- branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90
r2392 r2443 6 6 !!====================================================================== 7 7 !! History : 1.0 ! 2005-10 (A. Beckmann, G. Madec) reactivate s-coordinate 8 !!---------------------------------------------------------------------- 9 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 10 !! $Id$ 11 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 8 !! 3.3 ! 2010-11 (G. Madec) add mbk. arrays associated to the deepest ocean level 12 9 !!---------------------------------------------------------------------- 13 10 USE par_oce ! ocean parameters … … 172 169 !! masks, bathymetry 173 170 !! --------------------------------------------------------------------- 174 INTEGER , PUBLIC, DIMENSION(jpi,jpj) :: mbathy !: number of ocean level (=0, 1, ... , jpk-1) 175 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: bathy !: ocean depth (meters) 176 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: tmask_i !: interior domain T-point mask 177 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: bmask !: land/ocean mask of barotropic stream function 178 179 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: tmask, umask, vmask, fmask !: land/ocean mask at T-, U-, V- and F-points 180 181 REAL(wp), PUBLIC, DIMENSION(jpiglo) :: tpol, fpol !: north fold mask (nperio= 3 or 4) 171 INTEGER , PUBLIC, DIMENSION(jpi,jpj) :: mbathy !: number of ocean level (=0, 1, ... , jpk-1) 172 INTEGER , PUBLIC, DIMENSION(jpi,jpj) :: mbkt !: vertical index of the bottom last T- ocean level 173 INTEGER , PUBLIC, DIMENSION(jpi,jpj) :: mbku, mbkv !: vertical index of the bottom last U- and W- ocean level 174 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: bathy !: ocean depth (meters) 175 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: tmask_i !: interior domain T-point mask 176 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: bmask !: land/ocean mask of barotropic stream function 177 178 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: tmask, umask, vmask, fmask !: land/ocean mask at T-, U-, V- and F-pts 179 180 REAL(wp), PUBLIC, DIMENSION(jpiglo) :: tpol, fpol !: north fold mask (jperio= 3 or 4) 182 181 183 182 #if defined key_noslip_accurate … … 233 232 END FUNCTION Agrif_CFixed 234 233 #endif 234 !!---------------------------------------------------------------------- 235 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 236 !! $Id$ 237 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 235 238 !!====================================================================== 236 239 END MODULE dom_oce -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90
r2436 r2443 14 14 !! 3.0 ! 2008-06 (G. Madec) insertion of domzgr_zps.h90 & conding style 15 15 !! 3.2 ! 2009-07 (R. Benshila) Suppression of rigid-lid option 16 !! 3.3 ! 2010-11 (G. Madec) add mbk. arrays associated to the deepest ocean level 16 17 !!---------------------------------------------------------------------- 17 18 … … 497 498 IF( .NOT.lk_c1d ) CALL zgr_bat_ctl ! Bathymetry check ! 498 499 ! ! =================== ! 500 ! 501 ! ! ========================= ! 502 CALL zgr_bot_level ! level of ocean bottom ! 503 ! ! ========================= ! 499 504 END SUBROUTINE zgr_bat 500 505 … … 682 687 ! 683 688 END SUBROUTINE zgr_bat_ctl 689 690 691 SUBROUTINE zgr_bot_level 692 !!---------------------------------------------------------------------- 693 !! *** ROUTINE zgr_bot_level *** 694 !! 695 !! ** Purpose : defines the vertical index of ocean bottom (mbk. arrays) 696 !! 697 !! ** Method : computes from mbathy with a minimum value of 1 over land 698 !! 699 !! ** Action : mbkt, mbku, mbkv : vertical indices of the deeptest 700 !! ocean level at t-, u- & v-points 701 !! (min value = 1 over land) 702 !!---------------------------------------------------------------------- 703 INTEGER :: ji, jj ! dummy loop indices 704 REAL(wp), DIMENSION(jpi,jpj) :: zmbk ! 2D workspace 705 !!---------------------------------------------------------------------- 706 ! 707 IF(lwp) WRITE(numout,*) 708 IF(lwp) WRITE(numout,*) ' zgr_bot_level : ocean bottom k-index of T-, U-, V- and W-levels ' 709 IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~' 710 ! 711 mbkt(:,:) = MAX( mbathy(:,:) , 1 ) ! bottom k-index of T-level (=1 over land) 712 ! ! bottom k-index of W-level = mbkt+1 713 DO jj = 1, jpjm1 ! bottom k-index of u- (v-) level 714 DO ji = 1, jpim1 715 mbku(ji,jj) = MIN( mbkt(ji+1,jj ) , mbkt(ji,jj) ) 716 mbkv(ji,jj) = MIN( mbkt(ji ,jj+1) , mbkt(ji,jj) ) 717 END DO 718 END DO 719 ! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk 720 zmbk(:,:) = REAL( mbku(:,:), wp ) ; CALL lbc_lnk(zmbk,'U',1.) ; mbku (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 721 zmbk(:,:) = REAL( mbkv(:,:), wp ) ; CALL lbc_lnk(zmbk,'V',1.) ; mbkv (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 722 ! 723 END SUBROUTINE zgr_bot_level 684 724 685 725 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/TRA/trabbc.F90
r2337 r2443 8 8 !! - ! 2002-11 (A. Bozec) tra_bbc_init: original code 9 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) 10 11 !!---------------------------------------------------------------------- 11 12 13 !!---------------------------------------------------------------------- 12 14 !! tra_bbc : update the tracer trend at ocean bottom 13 15 !! tra_bbc_init : initialization of geothermal heat flux trend 14 16 !!---------------------------------------------------------------------- 15 USE oce ! ocean dynamics and active tracers16 USE dom_oce ! ocean space and time domain17 USE oce ! ocean variables 18 USE dom_oce ! domain: ocean 17 19 USE phycst ! physical constants 18 USE trdmod_oce ! ocean trends19 USE trdtra ! ocean trends20 USE trdmod_oce ! trends: ocean variables 21 USE trdtra ! trends: active tracers 20 22 USE in_out_manager ! I/O manager 21 23 USE prtctl ! Print control … … 32 34 REAL(wp) :: rn_geoflx_cst = 86.4e-3_wp ! Constant value of geothermal heat flux 33 35 34 INTEGER , DIMENSION(:,:), ALLOCATABLE :: nbotlevt ! ocean bottom level index at T-pt35 36 REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE :: qgh_trd0 ! geothermal heating trend 36 37 … … 57 58 !! ocean bottom can be computed once and is added to the temperature 58 59 !! trend juste above the bottom at each time step: 59 !! ta = ta + Qsf / (rau0 rcp e3T) for k= mb athy -160 !! ta = ta + Qsf / (rau0 rcp e3T) for k= mbkt 60 61 !! Where Qsf is the geothermal heat flux. 61 62 !! … … 66 67 !! Emile-Geay and Madec, 2009, Ocean Science. 67 68 !!---------------------------------------------------------------------- 68 !!69 69 INTEGER, INTENT( in ) :: kt ! ocean time-step index 70 70 !! 71 71 INTEGER :: ji, jj, ik ! dummy loop indices 72 REAL(wp) :: zqgh_trd ! geothermal heat flux trend72 REAL(wp) :: zqgh_trd ! geothermal heat flux trend 73 73 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ztrdt 74 74 !!---------------------------------------------------------------------- 75 75 ! 76 76 IF( l_trdtra ) THEN ! Save ta and sa trends 77 77 ALLOCATE( ztrdt(jpi,jpj,jpk) ) ; ztrdt(:,:,:) = tsa(:,:,:,jp_tem) 78 78 ENDIF 79 79 ! 80 80 ! ! Add the geothermal heat flux trend on temperature 81 81 #if defined key_vectopt_loop … … 86 86 DO ji = 2, jpim1 87 87 #endif 88 ik = nbotlevt(ji,jj)88 ik = mbkt(ji,jj) 89 89 zqgh_trd = qgh_trd0(ji,jj) / fse3t(ji,jj,ik) 90 90 tsa(ji,jj,ik,jp_tem) = tsa(ji,jj,ik,jp_tem) + zqgh_trd 91 91 END DO 92 92 END DO 93 93 ! 94 94 IF( l_trdtra ) THEN ! Save the geothermal heat flux trend for diagnostics 95 95 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:) … … 117 117 !! 118 118 !! ** Action : - read/fix the geothermal heat qgh_trd0 119 !! - compute the bottom ocean level nbotlevt120 119 !!---------------------------------------------------------------------- 121 120 USE iom … … 127 126 !!---------------------------------------------------------------------- 128 127 129 REWIND ( numnam )! Read Namelist nambbc : bottom momentum boundary condition130 READ 128 REWIND( numnam ) ! Read Namelist nambbc : bottom momentum boundary condition 129 READ ( numnam, nambbc ) 131 130 132 131 IF(lwp) THEN ! Control print … … 143 142 IF( ln_trabbc ) THEN !== geothermal heating ==! 144 143 ! 145 ALLOCATE( nbotlevt(jpi,jpj) ) ! allocation 146 ALLOCATE( qgh_trd0(jpi,jpj) ) 147 ! 148 DO jj = 1, jpj ! level of the ocean bottom at T-point 149 DO ji = 1, jpi 150 nbotlevt(ji,jj) = MAX( mbathy(ji,jj)-1, 1 ) 151 END DO 152 END DO 144 ALLOCATE( qgh_trd0(jpi,jpj) ) ! allocation 153 145 ! 154 146 SELECT CASE ( nn_geoflx ) ! geothermal heat flux / (rauO * Cp) … … 172 164 ! 173 165 ELSE 174 166 IF(lwp) WRITE(numout,*) ' *** no geothermal heat flux' 175 167 ENDIF 176 168 ! -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/TRA/trabbl.F90
r2287 r2443 4 4 !! Ocean physics : advective and/or diffusive bottom boundary layer scheme 5 5 !!============================================================================== 6 !! History : OPA ! 1996-06 (L. Mortier) Original code 7 !! 8.0 ! 1997-11 (G. Madec) Optimization 8 !! NEMO 1.0 ! 2002-08 (G. Madec) free form + modules 9 !! - ! 2004-01 (A. de Miranda, G. Madec, J.M. Molines ) add advective bbl 10 !! 3.3 ! 2009-11 (G. Madec) merge trabbl and trabbl_adv + style + optimization 11 !! - ! 2010-04 (G. Madec) Campin & Goosse advective bbl 12 !! - ! 2010-06 (C. Ethe, G. Madec) merge TRA-TRC 6 !! History : OPA ! 1996-06 (L. Mortier) Original code 7 !! 8.0 ! 1997-11 (G. Madec) Optimization 8 !! NEMO 1.0 ! 2002-08 (G. Madec) free form + modules 9 !! - ! 2004-01 (A. de Miranda, G. Madec, J.M. Molines ) add advective bbl 10 !! 3.3 ! 2009-11 (G. Madec) merge trabbl and trabbl_adv + style + optimization 11 !! - ! 2010-04 (G. Madec) Campin & Goosse advective bbl 12 !! - ! 2010-06 (C. Ethe, G. Madec) merge TRA-TRC 13 !! - ! 2010-11 (G. Madec) add mbk. arrays associated to the deepest ocean level 13 14 !!---------------------------------------------------------------------- 14 15 #if defined key_trabbl || defined key_esopa … … 24 25 USE oce ! ocean dynamics and active tracers 25 26 USE dom_oce ! ocean space and time domain 26 USE phycst ! 27 USE phycst ! physical constant 27 28 USE eosbn2 ! equation of state 28 USE trdmod_oce ! ocean space and time domain29 USE trdtra ! ocean active tracers trends29 USE trdmod_oce ! trends: ocean variables 30 USE trdtra ! trends: active tracers 30 31 USE iom ! IOM server 31 32 USE in_out_manager ! I/O manager … … 38 39 PUBLIC tra_bbl ! routine called by step.F90 39 40 PUBLIC tra_bbl_init ! routine called by opa.F90 40 PUBLIC tra_bbl_dif ! routine called by tr a_bbl and trc_bbl41 PUBLIC tra_bbl_dif ! routine called by trcbbl.F90 41 42 PUBLIC tra_bbl_adv ! - - - - 42 PUBLIC bbl ! - - - -43 PUBLIC bbl ! routine called by trcbbl.F90 and dtadyn.F90 43 44 44 45 # if defined key_trabbl … … 48 49 # endif 49 50 50 ! !!* Namelist nambbl *51 INTEGER , PUBLIC :: nn_bbl_ldf = 0 !: =1 : diffusive bbl or not (=0)52 INTEGER , PUBLIC :: nn_bbl_adv = 0 !: =1/2 : advective bbl or not (=0)53 ! ! =1 : advective bbl using the bottom ocean velocity54 ! ! =2 : - - using utr_bbl proportional to grad(rho)55 REAL(wp), PUBLIC :: rn_ahtbbl = 1.e +3!: along slope bbl diffusive coefficient [m2/s]56 REAL(wp), PUBLIC :: rn_gambbl = 10. e0!: lateral coeff. for bottom boundary layer scheme [s]51 ! !!* Namelist nambbl * 52 INTEGER , PUBLIC :: nn_bbl_ldf = 0 !: =1 : diffusive bbl or not (=0) 53 INTEGER , PUBLIC :: nn_bbl_adv = 0 !: =1/2 : advective bbl or not (=0) 54 ! ! =1 : advective bbl using the bottom ocean velocity 55 ! ! =2 : - - using utr_bbl proportional to grad(rho) 56 REAL(wp), PUBLIC :: rn_ahtbbl = 1.e3_wp !: along slope bbl diffusive coefficient [m2/s] 57 REAL(wp), PUBLIC :: rn_gambbl = 10.0_wp !: lateral coeff. for bottom boundary layer scheme [s] 57 58 58 59 REAL(wp), DIMENSION(jpi,jpj), PUBLIC :: utr_bbl , vtr_bbl ! u- (v-) transport in the bottom boundary layer 59 60 60 INTEGER , DIMENSION(jpi,jpj) :: mbkt ! vertical index of the bottom ocean T-level61 INTEGER , DIMENSION(jpi,jpj) :: mbku , mbkv ! vertical index of the (upper) bottom ocean U/V-level62 61 INTEGER , DIMENSION(jpi,jpj) :: mbku_d , mbkv_d ! vertical index of the "lower" bottom ocean U/V-level 63 62 INTEGER , DIMENSION(jpi,jpj) :: mgrhu , mgrhv ! = +/-1, sign of grad(H) in u-(v-)direction … … 74 73 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 75 74 !! $Id$ 76 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 77 !!---------------------------------------------------------------------- 78 75 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 76 !!---------------------------------------------------------------------- 79 77 CONTAINS 80 81 78 82 79 SUBROUTINE tra_bbl( kt ) … … 85 82 !! 86 83 !! ** Purpose : Compute the before tracer (t & s) trend associated 87 !! with the bottom boundary layer and add it to the general trend88 !! of tracer equations.84 !! with the bottom boundary layer and add it to the general 85 !! trend of tracer equations. 89 86 !! 90 87 !! ** Method : Depending on namtra_bbl namelist parameters the bbl … … 218 215 !! transport proportional to the along-slope density gradient 219 216 !! 220 !!221 217 !! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591. 222 218 !! Campin, J.-M., and H. Goosse, 1999, Tellus, 412-430. 223 !!224 219 !!---------------------------------------------------------------------- 225 220 INTEGER , INTENT(in ) :: kjpt ! number of tracers … … 233 228 REAL(wp) :: zu_bbl, zv_bbl ! - - 234 229 !!---------------------------------------------------------------------- 235 230 ! 236 231 ! ! =========== 237 232 DO jn = 1, kjpt ! tracer loop … … 547 542 e1e2t_r(:,:) = 1.0 / ( e1t(:,:) * e2t(:,:) ) 548 543 549 ! !* vertical index of bottom t-, u- and v-points 550 DO jj = 1, jpj ! bottom k-index of T-level 551 DO ji = 1, jpi 552 mbkt(ji,jj) = MAX( mbathy(ji,jj) - 1, 1 ) 544 ! !* vertical index of "deep" bottom u- and v-points 545 DO jj = 1, jpjm1 ! (the "shelf" bottom k-indices are mbku and mbkv) 546 DO ji = 1, jpim1 547 mbku_d(ji,jj) = MAX( mbkt(ji+1,jj ) , mbkt(ji,jj) ) ! >= 1 as mbkt=1 over land 548 mbkv_d(ji,jj) = MAX( mbkt(ji ,jj+1) , mbkt(ji,jj) ) 553 549 END DO 554 550 END DO 555 DO jj = 1, jpjm1 ! bottom k-index of u- (v-) level (shelf and deep) 556 DO ji = 1, jpim1 557 mbku (ji,jj) = MAX( MIN( mbathy(ji+1,jj ), mbathy(ji,jj) ) - 1, 1 ) ! "shelf" 558 mbkv (ji,jj) = MAX( MIN( mbathy(ji ,jj+1), mbathy(ji,jj) ) - 1, 1 ) 559 mbku_d(ji,jj) = MAX( MAX( mbathy(ji+1,jj ), mbathy(ji,jj) ) - 1, 1 ) ! "deep" 560 mbkv_d(ji,jj) = MAX( MAX( mbathy(ji ,jj+1), mbathy(ji,jj) ) - 1, 1 ) 561 END DO 562 END DO 563 zmbk(:,:) = FLOAT( mbku (:,:) ) ; CALL lbc_lnk(zmbk,'U',1.) ; mbku (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 564 zmbk(:,:) = FLOAT( mbkv (:,:) ) ; CALL lbc_lnk(zmbk,'V',1.) ; mbkv (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 565 zmbk(:,:) = FLOAT( mbku_d(:,:) ) ; CALL lbc_lnk(zmbk,'U',1.) ; mbku_d(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 566 zmbk(:,:) = FLOAT( mbkv_d(:,:) ) ; CALL lbc_lnk(zmbk,'V',1.) ; mbkv_d(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 551 ! converte into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk 552 zmbk(:,:) = REAL( mbku_d(:,:), wp ) ; CALL lbc_lnk(zmbk,'U',1.) ; mbku_d(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 553 zmbk(:,:) = REAL( mbkv_d(:,:), wp ) ; CALL lbc_lnk(zmbk,'V',1.) ; mbkv_d(:,:) = MAX( INT( zmbk(:,:) ), 1 ) 567 554 568 555 !* sign of grad(H) at u- and v-points
Note: See TracChangeset
for help on using the changeset viewer.