Changeset 2528 for trunk/NEMOGCM/NEMO/OFF_SRC/dtadyn.F90
- Timestamp:
- 2010-12-27T18:33:53+01:00 (13 years ago)
- File:
-
- 1 edited
-
trunk/NEMOGCM/NEMO/OFF_SRC/dtadyn.F90 (modified) (34 diffs, 1 prop)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OFF_SRC/dtadyn.F90
- Property svn:eol-style deleted
r1735 r2528 2 2 !!====================================================================== 3 3 !! *** MODULE dtadyn *** 4 !! OFFLINE : interpolation of the physical fields 5 !!===================================================================== 4 !! Off-line : interpolation of the physical fields 5 !!====================================================================== 6 !! History : OPA ! 1992-01 (M. Imbard) Original code 7 !! 8.0 ! 1998-04 (L.Bopp MA Foujols) slopes for isopyc. 8 !! - ! 1998-05 (L. Bopp) read output of coupled run 9 !! 8.2 ! 2001-01 (M. Levy et M. Benjelloul) add netcdf FORMAT 10 !! NEMO 1.0 ! 2005-03 (O. Aumont and A. El Moussaoui) F90 11 !! - ! 2005-12 (C. Ethe) Adapted for DEGINT 12 !! 3.0 ! 2007-06 (C. Ethe) use of iom module 13 !! - ! 2007-09 (C. Ethe) add swap_dyn_data 14 !! 3.3 ! 2010-11 (C. Ethe) Full reorganization of the off-line: phasing with the on-line 15 !!---------------------------------------------------------------------- 6 16 7 17 !!---------------------------------------------------------------------- … … 9 19 !! dta_dyn : Interpolation of the fields 10 20 !!---------------------------------------------------------------------- 11 !! * Modules used12 21 USE oce ! ocean dynamics and tracers variables 13 USE dom_oce ! ocean space and time domain variables 14 USE zdf_oce ! ocean vertical physics 22 USE c1d ! 1D configuration: lk_c1d 23 USE dom_oce ! ocean domain: variables 24 USE zdf_oce ! ocean vertical physics: variables 25 USE sbc_oce ! surface module: variables 26 USE phycst ! physical constants 27 USE trabbl ! active tracer: bottom boundary layer 28 USE ldfslp ! lateral diffusion: iso-neutral slopes 29 USE ldfeiv ! eddy induced velocity coef. 30 USE ldftra_oce ! ocean tracer lateral physics 31 USE zdfmxl ! vertical physics: mixed layer depth 32 USE eosbn2 ! equation of state - Brunt Vaisala frequency 33 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 34 USE zpshde ! z-coord. with partial steps: horizontal derivatives 15 35 USE in_out_manager ! I/O manager 16 USE phycst ! physical constants 17 USE sbc_oce 18 USE ldfslp 19 USE ldfeiv ! eddy induced velocity coef. (ldf_eiv routine) 20 USE ldftra_oce ! ocean tracer lateral physics 21 USE zdfmxl 22 USE trabbl 23 USE eosbn2 24 USE zdfddm ! vertical physics: double diffusion 25 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 26 USE zpshde 36 USE iom ! I/O library 27 37 USE lib_mpp ! distributed memory computing library 38 USE prtctl ! print control 28 39 29 40 IMPLICIT NONE 30 41 PRIVATE 31 42 32 !! * Routine accessibility 33 PUBLIC dta_dyn_init ! called by opa.F90 34 PUBLIC dta_dyn ! called by step.F90 35 36 LOGICAL , PUBLIC :: & 37 lperdyn = .TRUE. , & ! boolean for periodic fields or not 38 lfirdyn = .TRUE. ! boolean for the first call or not 39 40 INTEGER , PUBLIC :: & 41 ndtadyn = 73 , & ! Number of dat in one year 42 ndtatot = 73 , & ! Number of data in the input field 43 nsptint = 1 , & ! type of spatial interpolation 44 nficdyn = 2 ! number of dynamical fields 45 46 CHARACTER(len=45) :: & 47 cfile_grid_T = 'dyna_grid_T.nc', & !: name of the grid_T file 48 cfile_grid_U = 'dyna_grid_U.nc', & !: name of the grid_U file 49 cfile_grid_V = 'dyna_grid_V.nc', & !: name of the grid_V file 50 cfile_grid_W = 'dyna_grid_W.nc' !: name of the grid_W file 43 PUBLIC dta_dyn_init ! called by opa.F90 44 PUBLIC dta_dyn ! called by step.F90 45 46 LOGICAL, PUBLIC :: lperdyn = .TRUE. !: boolean for periodic fields or not 47 LOGICAL, PUBLIC :: lfirdyn = .TRUE. !: boolean for the first call or not 48 49 INTEGER, PUBLIC :: ndtadyn = 73 !: Number of dat in one year 50 INTEGER, PUBLIC :: ndtatot = 73 !: Number of data in the input field 51 INTEGER, PUBLIC :: nsptint = 1 !: type of spatial interpolation 52 53 CHARACTER(len=45) :: cfile_grid_T = 'dyna_grid_T.nc' ! name of the grid_T file 54 CHARACTER(len=45) :: cfile_grid_U = 'dyna_grid_U.nc' ! name of the grid_U file 55 CHARACTER(len=45) :: cfile_grid_V = 'dyna_grid_V.nc' ! name of the grid_V file 56 CHARACTER(len=45) :: cfile_grid_W = 'dyna_grid_W.nc' ! name of the grid_W file 51 57 52 REAL(wp) :: & 53 rnspdta , & !: number of time step per 2 consecutives data 54 rnspdta2 !: rnspdta * 0.5 55 56 INTEGER :: & 57 ndyn1, ndyn2 , & 58 nlecoff = 0 , & ! switch for the first read 59 numfl_t, numfl_u, & 60 numfl_v, numfl_w 61 62 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: & 63 tdta , & ! temperature at two consecutive times 64 sdta , & ! salinity at two consecutive times 65 udta , & ! zonal velocity at two consecutive times 66 vdta , & ! meridional velocity at two consecutive times 67 wdta , & ! vertical velocity at two consecutive times 68 #if defined key_trc_diatrd 69 hdivdta, & ! horizontal divergence 70 #endif 71 avtdta ! vertical diffusivity coefficient 72 73 REAL(wp), DIMENSION(jpi,jpj,2) :: & 74 hmlddta, & ! mixed layer depth at two consecutive times 75 wspddta, & ! wind speed at two consecutive times 76 frlddta, & ! sea-ice fraction at two consecutive times 77 empdta , & ! E-P at two consecutive times 78 qsrdta ! short wave heat flux at two consecutive times 79 58 REAL(wp) :: rnspdta ! number of time step per 2 consecutives data 59 REAL(wp) :: rnspdta2 ! rnspdta * 0.5 60 61 INTEGER :: ndyn1, ndyn2 ! 62 INTEGER :: nlecoff = 0 ! switch for the first read 63 INTEGER :: numfl_t, numfl_u, numfl_v, numfl_w 64 65 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: tdta ! temperature at two consecutive times 66 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: sdta ! salinity at two consecutive times 67 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: udta ! zonal velocity at two consecutive times 68 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: vdta ! meridional velocity at two consecutive times 69 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: wdta ! vertical velocity at two consecutive times 70 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: avtdta ! vertical diffusivity coefficient 71 72 REAL(wp), DIMENSION(jpi,jpj ,2) :: hmlddta ! mixed layer depth at two consecutive times 73 REAL(wp), DIMENSION(jpi,jpj ,2) :: wspddta ! wind speed at two consecutive times 74 REAL(wp), DIMENSION(jpi,jpj ,2) :: frlddta ! sea-ice fraction at two consecutive times 75 REAL(wp), DIMENSION(jpi,jpj ,2) :: empdta ! E-P at two consecutive times 76 REAL(wp), DIMENSION(jpi,jpj ,2) :: qsrdta ! short wave heat flux at two consecutive times 77 REAL(wp), DIMENSION(jpi,jpj ,2) :: bblxdta ! frequency of bbl in the x direction at 2 consecutive times 78 REAL(wp), DIMENSION(jpi,jpj ,2) :: bblydta ! frequency of bbl in the y direction at 2 consecutive times 79 LOGICAL :: l_offbbl 80 80 #if defined key_ldfslp 81 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: & 82 uslpdta , & ! zonal isopycnal slopes 83 vslpdta , & ! meridional isopycnal slopes 84 wslpidta , & ! zonal diapycnal slopes 85 wslpjdta ! meridional diapycnal slopes 86 #endif 87 88 #if ! defined key_off_degrad && defined key_traldf_c2d 89 REAL(wp), DIMENSION(jpi,jpj,2) :: & 90 ahtwdta ! Lateral diffusivity 91 # if defined key_trcldf_eiv 92 REAL(wp), DIMENSION(jpi,jpj,2) :: & 93 aeiwdta ! G&M coefficient 81 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: uslpdta ! zonal isopycnal slopes 82 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: vslpdta ! meridional isopycnal slopes 83 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: wslpidta ! zonal diapycnal slopes 84 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: wslpjdta ! meridional diapycnal slopes 85 #endif 86 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 87 REAL(wp), DIMENSION(jpi,jpj ,2) :: aeiwdta ! G&M coefficient 88 #endif 89 #if defined key_degrad 90 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: ahtudta, ahtvdta, ahtwdta ! Lateral diffusivity 91 # if defined key_traldf_eiv 92 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: aeiudta, aeivdta, aeiwdta ! G&M coefficient 94 93 # endif 95 #endif96 97 #if defined key_off_degrad98 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: &99 ahtudta, ahtvdta, ahtwdta ! Lateral diffusivity100 # if defined key_trcldf_eiv101 REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: &102 aeiudta, aeivdta, aeiwdta ! G&M coefficient103 # endif104 105 #endif106 107 #if defined key_trcbbl_dif || defined key_trcbbl_adv108 REAL(wp), DIMENSION(jpi,jpj,2) :: &109 bblxdta , & ! frequency of bbl in the x direction at 2 consecutive times110 bblydta ! frequency of bbl in the y direction at 2 consecutive times111 94 #endif 112 95 … … 115 98 # include "vectopt_loop_substitute.h90" 116 99 !!---------------------------------------------------------------------- 117 !! OPA 9.0 , LOCEAN-IPSL (2005)118 !! 119 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt100 !! NEMO/OFF 3.3 , NEMO Consortium (2010) 101 !! $Id$ 102 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 120 103 !!---------------------------------------------------------------------- 121 122 104 CONTAINS 123 105 … … 126 108 !! *** ROUTINE dta_dyn *** 127 109 !! 128 !! ** Purpose : Prepares dynamics and physics fields from an 129 !! OPA9 simulation for an off-line simulation 130 !! for passive tracer 110 !! ** Purpose : Prepares dynamics and physics fields from an NEMO run 111 !! for an off-line simulation of passive tracers 131 112 !! 132 113 !! ** Method : calculates the position of DATA to read READ DATA 133 114 !! (example month changement) computes slopes IF needed 134 115 !! interpolates DATA IF needed 135 !! 136 !! ** History : 137 !! ! original : 92-01 (M. Imbard: sub domain) 138 !! ! addition : 98-04 (L.Bopp MA Foujols: slopes for isopyc.) 139 !! ! addition : 98-05 (L. Bopp read output of coupled run) 140 !! ! addition : 05-03 (O. Aumont and A. El Moussaoui) F90 141 !! ! addition : 05-12 (C. Ethe) Adapted for DEGINT 142 !!---------------------------------------------------------------------- 143 !! * Arguments 144 INTEGER, INTENT( in ) :: kt ! ocean time-step index 145 146 !! * Local declarations 147 INTEGER :: iper, iperm1, iswap, izt 148 149 REAL(wp) :: zpdtan, zpdtpe, zdemi, zt 150 REAL(wp) :: zweigh 151 152 ! 0. Initialization 153 ! ----------------- 154 155 IF( lfirdyn ) THEN 156 ! first time step MUST BE nit000 157 IF( kt /= nit000 ) THEN 158 IF (lwp) THEN 159 WRITE (numout,*) ' kt MUST BE EQUAL to nit000. kt = ',kt ,' nit000 = ',nit000 160 STOP 'dtadyn' 161 ENDIF 162 ENDIF 163 ! Initialize the parameters of the interpolation 164 CALL dta_dyn_init 165 ENDIF 166 167 zt = ( FLOAT (kt) + rnspdta2 ) / rnspdta 168 izt = INT( zt ) 169 zweigh = zt - FLOAT( INT(zt) ) 170 171 IF( lperdyn ) THEN 172 iperm1 = MOD( izt, ndtadyn ) 173 ELSE 174 iperm1 = MOD( izt, ndtatot - 1 ) + 1 116 !!---------------------------------------------------------------------- 117 INTEGER, INTENT(in) :: kt ! ocean time-step index 118 !! 119 INTEGER :: iper, iperm1, iswap, izt ! local integers 120 REAL(wp) :: zt, zweigh ! local scalars 121 !!---------------------------------------------------------------------- 122 123 zt = ( REAL(kt,wp) + rnspdta2 ) / rnspdta 124 izt = INT( zt ) 125 zweigh = zt - REAL( INT(zt), wp ) 126 127 IF( lperdyn ) THEN ; iperm1 = MOD( izt, ndtadyn ) 128 ELSE ; iperm1 = MOD( izt, ndtatot - 1 ) + 1 175 129 ENDIF 176 130 … … 181 135 ELSE 182 136 IF( lfirdyn ) THEN 183 IF (lwp) WRITE (numout,*) & 184 & ' dynamic file is not periodic with or without interpolation & 185 & we take the first value for the previous period iperm1 = 0 ' 137 IF(lwp) WRITE (numout,*) 'dta_dyn: dynamic file is not periodic with or without interpolation & 138 & we take the first value for the previous period iperm1 = 0 ' 186 139 END IF 187 140 END IF … … 194 147 195 148 IF( lfirdyn ) THEN 196 ! store the information of the period read 197 ndyn1 = iperm1 149 ndyn1 = iperm1 ! store the information of the period read 198 150 ndyn2 = iper 199 151 200 IF 201 WRITE (numout,*) ' dynamics data read for the period ndyn1 =', ndyn1,&202 & ' and for the period ndyn2 = ', ndyn2152 IF(lwp) THEN 153 WRITE (numout,*) ' dynamics data read for the period ndyn1 =', ndyn1, & 154 & ' and for the period ndyn2 = ', ndyn2 203 155 WRITE (numout,*) ' time step is : ', kt 204 WRITE (numout,*) ' we have ndtadyn = ', ndtadyn,' records in the dynamic file for one year'156 WRITE (numout,*) ' we have ndtadyn = ', ndtadyn, ' records in the dynamic file for one year' 205 157 END IF 206 158 ! 207 IF( iperm1 /= 0 ) THEN ! data read for the iperm1 period 208 CALL dynrea( kt, iperm1 ) 209 ELSE 210 CALL dynrea( kt, 1 ) 211 ENDIF 159 CALL dynrea( kt, MAX( 1, iperm1) ) ! data read for the iperm1 period 212 160 213 #if defined key_ldfslp 214 ! Computes slopes 215 ! Caution : here tn, sn and avt are used as workspace 216 tn (:,:,:) = tdta (:,:,:,2) 217 sn (:,:,:) = sdta (:,:,:,2) 218 avt(:,:,:) = avtdta(:,:,:,2) 161 IF( lk_ldfslp .AND. .NOT. lk_c1d ) THEN ! Computes slopes (here tsn and avt are used as workspace) 162 tsn (:,:,:,jp_tem) = tdta (:,:,:,2) 163 tsn (:,:,:,jp_sal) = sdta (:,:,:,2) 164 avt(:,:,:) = avtdta(:,:,:,2) 219 165 220 CALL eos( tn, sn, rhd, rhop ) ! Time-filtered in situ density 221 CALL bn2( tn, sn, rn2 ) ! before Brunt-Vaisala frequency 222 IF( ln_zps ) & 223 & CALL zps_hde( kt, tn , sn , rhd, & ! Partial steps: before Horizontal DErivative 224 & gtu, gsu, gru, & ! of t, s, rd at the bottom ocean level 225 & gtv, gsv, grv ) 226 CALL zdf_mxl( kt ) ! mixed layer depth 227 CALL ldf_slp( kt, rhd, rn2 ) 166 CALL eos( tsn, rhd, rhop ) ! Time-filtered in situ density 167 CALL bn2( tsn, rn2 ) ! before Brunt-Vaisala frequency 168 IF( ln_zps ) & 169 & CALL zps_hde( kt, jpts, tsn, gtsu, gtsv, & ! Partial steps: before Horizontal DErivative 170 & rhd, gru , grv ) ! of t, s, rd at the bottom ocean level 171 CALL zdf_mxl( kt ) ! mixed layer depth 172 CALL ldf_slp( kt, rhd, rn2 ) 228 173 229 uslpdta (:,:,:,2) = uslp (:,:,:) 230 vslpdta (:,:,:,2) = vslp (:,:,:) 231 wslpidta(:,:,:,2) = wslpi(:,:,:) 232 wslpjdta(:,:,:,2) = wslpj(:,:,:) 233 #endif 234 235 ! swap from record 2 to 1 236 CALL swap_dyn_data 237 174 uslpdta (:,:,:,2) = uslp (:,:,:) 175 vslpdta (:,:,:,2) = vslp (:,:,:) 176 wslpidta(:,:,:,2) = wslpi(:,:,:) 177 wslpjdta(:,:,:,2) = wslpj(:,:,:) 178 END IF 179 ! 180 CALL swap_dyn_data ! swap from record 2 to 1 181 ! 238 182 iswap = 1 ! indicates swap 239 240 CALL dynrea( kt, iper ) ! data read for the iper period 241 242 #if defined key_ldfslp 243 ! Computes slopes 244 ! Caution : here tn, sn and avt are used as workspace 245 tn (:,:,:) = tdta (:,:,:,2) 246 sn (:,:,:) = sdta (:,:,:,2) 247 avt(:,:,:) = avtdta(:,:,:,2) 248 249 CALL eos( tn, sn, rhd, rhop ) ! Time-filtered in situ density 250 CALL bn2( tn, sn, rn2 ) ! before Brunt-Vaisala frequency 251 IF( ln_zps ) & 252 & CALL zps_hde( kt, tn , sn , rhd, & ! Partial steps: before Horizontal DErivative 253 & gtu, gsu, gru, & ! of t, s, rd at the bottom ocean level 254 & gtv, gsv, grv ) 255 CALL zdf_mxl( kt ) ! mixed layer depth 256 CALL ldf_slp( kt, rhd, rn2 ) 257 258 uslpdta (:,:,:,2) = uslp (:,:,:) 259 vslpdta (:,:,:,2) = vslp (:,:,:) 260 wslpidta(:,:,:,2) = wslpi(:,:,:) 261 wslpjdta(:,:,:,2) = wslpj(:,:,:) 262 #endif 263 ! 264 lfirdyn=.FALSE. ! trace the first call 183 ! 184 CALL dynrea( kt, iper ) ! data read for the iper period 185 ! 186 IF( lk_ldfslp .AND. .NOT. lk_c1d ) THEN ! Computes slopes (here tsn and avt are used as workspace) 187 tsn (:,:,:,jp_tem) = tdta (:,:,:,2) 188 tsn (:,:,:,jp_sal) = sdta (:,:,:,2) 189 avt(:,:,:) = avtdta(:,:,:,2) 190 ! 191 CALL eos( tsn, rhd, rhop ) ! now in situ density 192 CALL bn2( tsn, rn2 ) ! now Brunt-Vaisala frequency 193 IF( ln_zps ) CALL zps_hde( kt, jpts, tsn, gtsu, gtsv, & ! Partial steps: before Horizontal DErivative 194 & rhd, gru , grv ) ! of t, s, rd at the bottom ocean level 195 CALL zdf_mxl( kt ) ! mixed layer depth 196 CALL ldf_slp( kt, rhd, rn2 ) ! slope of iso-neutral surfaces 197 ! 198 uslpdta (:,:,:,2) = uslp (:,:,:) 199 vslpdta (:,:,:,2) = vslp (:,:,:) 200 wslpidta(:,:,:,2) = wslpi(:,:,:) 201 wslpjdta(:,:,:,2) = wslpj(:,:,:) 202 END IF 203 ! 204 lfirdyn = .FALSE. ! trace the first call 265 205 ENDIF 266 206 ! … … 269 209 ! 270 210 IF( iperm1 /= ndyn1 ) THEN 271 272 IF( iperm1 == 0 .) THEN273 IF 211 ! 212 IF( iperm1 == 0 ) THEN 213 IF(lwp) THEN 274 214 WRITE (numout,*) ' dynamic file is not periodic with periodic interpolation' 275 215 WRITE (numout,*) ' we take the last value for the last period ' … … 280 220 ENDIF 281 221 ! 282 ! We have to prepare a new read of data : swap from record 2 to 1 283 ! 284 CALL swap_dyn_data 285 286 iswap = 1 ! indicates swap 287 222 CALL swap_dyn_data ! We have to prepare a new read of data : swap from record 2 to 1 223 ! 224 iswap = 1 ! indicates swap 225 ! 288 226 CALL dynrea( kt, iper ) ! data read for the iper period 289 290 #if defined key_ldfslp 291 ! Computes slopes 292 ! Caution : here tn, sn and avt are used as workspace 293 tn (:,:,:) = tdta (:,:,:,2) 294 sn (:,:,:) = sdta (:,:,:,2) 295 avt(:,:,:) = avtdta(:,:,:,2) 296 297 CALL eos( tn, sn, rhd, rhop ) ! Time-filtered in situ density 298 CALL bn2( tn, sn, rn2 ) ! before Brunt-Vaisala frequency 299 IF( ln_zps ) & 300 & CALL zps_hde( kt, tn , sn , rhd, & ! Partial steps: before Horizontal DErivative 301 & gtu, gsu, gru, & ! of t, s, rd at the bottom ocean level 302 & gtv, gsv, grv ) 303 CALL zdf_mxl( kt ) ! mixed layer depth 304 CALL ldf_slp( kt, rhd, rn2 ) 305 306 uslpdta (:,:,:,2) = uslp (:,:,:) 307 vslpdta (:,:,:,2) = vslp (:,:,:) 308 wslpidta(:,:,:,2) = wslpi(:,:,:) 309 wslpjdta(:,:,:,2) = wslpj(:,:,:) 310 #endif 311 312 ! store the information of the period read 313 ndyn1 = ndyn2 227 ! 228 IF( lk_ldfslp .AND. .NOT. lk_c1d ) THEN 229 ! Computes slopes. Caution : here tsn and avt are used as workspace 230 tsn(:,:,:,jp_tem) = tdta (:,:,:,2) 231 tsn(:,:,:,jp_sal) = sdta (:,:,:,2) 232 avt(:,:,:) = avtdta(:,:,:,2) 233 ! 234 CALL eos( tsn, rhd, rhop ) ! now in situ density 235 CALL bn2( tsn, rn2 ) ! now Brunt-Vaisala frequency 236 IF( ln_zps ) CALL zps_hde( kt, jpts, tsn, gtsu, gtsv, & ! Partial steps: before Horizontal DErivative 237 & rhd, gru , grv ) ! of t, s, rd at the bottom ocean level 238 CALL zdf_mxl( kt ) ! mixed layer depth 239 CALL ldf_slp( kt, rhd, rn2 ) ! slope of iso-neutral surfaces 240 ! 241 uslpdta (:,:,:,2) = uslp (:,:,:) 242 vslpdta (:,:,:,2) = vslp (:,:,:) 243 wslpidta(:,:,:,2) = wslpi(:,:,:) 244 wslpjdta(:,:,:,2) = wslpj(:,:,:) 245 END IF 246 ! 247 ndyn1 = ndyn2 ! store the information of the period read 314 248 ndyn2 = iper 315 316 IF 317 WRITE (numout,*) ' dynamics data read for the period ndyn1 =', ndyn1,&318 & ' and for the period ndyn2 = ', ndyn2249 ! 250 IF(lwp) THEN 251 WRITE (numout,*) ' dynamics data read for the period ndyn1 =', ndyn1, & 252 & ' and for the period ndyn2 = ', ndyn2 319 253 WRITE (numout,*) ' time step is : ', kt 320 254 END IF … … 325 259 !---------------------------------------- 326 260 327 IF( nsptint == 0 ) THEN 328 ! No spatial interpolation, data are probably correct 329 ! We have to initialize data if we have changed the period 330 CALL assign_dyn_data 331 ELSE IF( nsptint == 1 ) THEN 332 ! linear interpolation 261 IF( nsptint == 0 ) THEN ! No space interpolation, data are probably correct 262 ! ! We have to initialize data if we have changed the period 263 CALL assign_dyn_data 264 ELSEIF( nsptint == 1 ) THEN ! linear interpolation 333 265 CALL linear_interp_dyn_data( zweigh ) 334 ELSE 335 ! other interpolation 266 ELSE ! other interpolation 336 267 WRITE (numout,*) ' this kind of interpolation do not exist at the moment : we stop' 337 268 STOP 'dtadyn' 338 269 END IF 339 340 ! In any case, we need rhop 341 CALL eos( tn, sn, rhd, rhop ) 342 343 #if ! defined key_off_degrad && defined key_traldf_c2d 344 ! In case of 2D varying coefficients, we need aeiv and aeiu 345 IF( lk_traldf_eiv ) CALL ldf_eiv( kt ) ! eddy induced velocity coefficient 346 #endif 347 270 ! 271 CALL eos( tsn, rhd, rhop ) ! In any case, we need rhop 272 ! 273 #if ! defined key_degrad && defined key_traldf_c2d 274 ! ! In case of 2D varying coefficients, we need aeiv and aeiu 275 IF( lk_traldf_eiv ) CALL dta_eiv( kt ) ! eddy induced velocity coefficient 276 #endif 277 ! 278 IF( .NOT. l_offbbl ) THEN ! Compute bbl coefficients if needed 279 tsb(:,:,:,:) = tsn(:,:,:,:) 280 CALL bbl( kt, 'TRC') 281 END IF 282 IF(ln_ctl) THEN 283 CALL prt_ctl(tab3d_1=tsn(:,:,:,jp_tem), clinfo1=' tn - : ', mask1=tmask, ovlap=1, kdim=jpk ) 284 CALL prt_ctl(tab3d_1=tsn(:,:,:,jp_sal), clinfo1=' sn - : ', mask1=tmask, ovlap=1, kdim=jpk ) 285 CALL prt_ctl(tab3d_1=un , clinfo1=' un - : ', mask1=tmask, ovlap=1, kdim=jpk ) 286 CALL prt_ctl(tab3d_1=vn , clinfo1=' vn - : ', mask1=tmask, ovlap=1, kdim=jpk ) 287 CALL prt_ctl(tab3d_1=wn , clinfo1=' wn - : ', mask1=tmask, ovlap=1, kdim=jpk ) 288 CALL prt_ctl(tab3d_1=avt , clinfo1=' kz - : ', mask1=tmask, ovlap=1, kdim=jpk ) 289 CALL prt_ctl(tab2d_1=fr_i , clinfo1=' fr_i - : ', mask1=tmask, ovlap=1 ) 290 CALL prt_ctl(tab2d_1=hmld , clinfo1=' hmld - : ', mask1=tmask, ovlap=1 ) 291 CALL prt_ctl(tab2d_1=emps , clinfo1=' emps - : ', mask1=tmask, ovlap=1 ) 292 CALL prt_ctl(tab2d_1=wndm , clinfo1=' wspd - : ', mask1=tmask, ovlap=1 ) 293 CALL prt_ctl(tab2d_1=qsr , clinfo1=' qsr - : ', mask1=tmask, ovlap=1 ) 294 ENDIF 295 ! 348 296 END SUBROUTINE dta_dyn 297 349 298 350 299 SUBROUTINE dynrea( kt, kenr ) … … 354 303 !! ** Purpose : READ dynamics fiels from OPA9 netcdf output 355 304 !! 356 !! ** Method : READ the kenr records of DATA and store in 357 !! in udta(...,2), .... 358 !! 359 !! ** History : additions : M. Levy et M. Benjelloul jan 2001 360 !! (netcdf FORMAT) 361 !! 05-03 (O. Aumont and A. El Moussaoui) F90 362 !! 06-07 : (C. Ethe) use of iom module 363 !!---------------------------------------------------------------------- 364 !! * Modules used 365 USE iom 366 367 !! * Arguments 368 INTEGER, INTENT( in ) :: kt, kenr ! time index 369 !! * Local declarations 305 !! ** Method : READ the kenr records of DATA and store in udta(...,2), .... 306 !!---------------------------------------------------------------------- 307 INTEGER, INTENT(in) :: kt, kenr ! time index 308 !! 370 309 INTEGER :: jkenr 371 372 REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 373 zu, zv, zw, zt, zs, zavt , & ! 3-D dynamical fields 374 zhdiv ! horizontal divergence 375 376 REAL(wp), DIMENSION(jpi,jpj) :: & 377 zemp, zqsr, zmld, zice, zwspd, & 378 ztaux, ztauy 379 #if defined key_trcbbl_dif || defined key_trcbbl_adv 380 REAL(wp), DIMENSION(jpi,jpj) :: zbblx, zbbly 381 #endif 382 383 #if ! defined key_off_degrad && defined key_traldf_c2d 384 REAL(wp), DIMENSION(jpi,jpj) :: zahtw 385 # if defined key_trcldf_eiv 310 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zu, zv, zw, zt, zs, zavt , zhdiv ! 3D workspace 311 REAL(wp), DIMENSION(jpi,jpj) :: zemp, zqsr, zmld, zice, zwspd, ztaux, ztauy ! 2D workspace 312 REAL(wp), DIMENSION(jpi,jpj) :: zbblx, zbbly 313 314 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 386 315 REAL(wp), DIMENSION(jpi,jpj) :: zaeiw 387 # endif 388 #endif 389 390 #if defined key_off_degrad 391 REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 392 zahtu, zahtv, zahtw ! Lateral diffusivity 393 # if defined key_trcldf_eiv 394 REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 395 zaeiu, zaeiv, zaeiw ! G&M coefficient 316 #endif 317 #if defined key_degrad 318 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zahtu, zahtv, zahtw ! Lateral diffusivity 319 # if defined key_traldf_eiv 320 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zaeiu, zaeiv, zaeiw ! G&M coefficient 396 321 # endif 397 322 #endif 398 399 !--------------------------------------------------------------- 323 !!---------------------------------------------------------------------- 324 400 325 ! 0. Initialization 401 326 … … 407 332 IF(lwp) THEN 408 333 WRITE(numout,*) 409 WRITE(numout,*) 'Dynrea : read ingdynamical fields, kenr = ', jkenr410 WRITE(numout,*) ' 411 #if defined key_ off_degrad334 WRITE(numout,*) 'Dynrea : read dynamical fields, kenr = ', jkenr 335 WRITE(numout,*) '~~~~~~~' 336 #if defined key_degrad 412 337 WRITE(numout,*) ' Degraded fields' 413 338 #endif … … 442 367 CALL iom_get( numfl_u, jpdom_data, 'vozocrtx', zu (:,:,:), jkenr ) 443 368 CALL iom_get( numfl_v, jpdom_data, 'vomecrty', zv (:,:,:), jkenr ) 444 445 #if defined key_trcbbl_dif || defined key_trcbbl_adv 446 IF( iom_varid( numfl_u, 'sobblcox', ldstop = .FALSE. ) > 0 .AND. & 447 & iom_varid( numfl_v, 'sobblcoy', ldstop = .FALSE. ) > 0 ) THEN 448 CALL iom_get( numfl_u, jpdom_data, 'sobblcox', zbblx(:,:), jkenr ) 449 CALL iom_get( numfl_v, jpdom_data, 'sobblcoy', zbbly(:,:), jkenr ) 450 ELSE 451 CALL bbl_sign( zt, zs, zbblx, zbbly ) 452 ENDIF 453 #endif 369 IF( lk_trabbl .AND. .NOT. lk_c1d .AND. nn_bbl_ldf == 1 ) THEN 370 IF( iom_varid( numfl_u, 'sobblcox', ldstop = .FALSE. ) > 0 .AND. & 371 & iom_varid( numfl_v, 'sobblcoy', ldstop = .FALSE. ) > 0 ) THEN 372 CALL iom_get( numfl_u, jpdom_data, 'sobblcox', zbblx(:,:), jkenr ) 373 CALL iom_get( numfl_v, jpdom_data, 'sobblcoy', zbbly(:,:), jkenr ) 374 l_offbbl = .TRUE. 375 ENDIF 376 ENDIF 454 377 455 378 ! file grid-W … … 458 381 CALL wzv( zu, zv, zw, zhdiv ) 459 382 460 # if defined key_zdfddm 461 CALL iom_get( numfl_w, jpdom_data, 'voddmavs', zavt (:,:,:), jkenr ) 462 #else 463 CALL iom_get( numfl_w, jpdom_data, 'votkeavt', zavt (:,:,:), jkenr ) 464 #endif 465 466 #if ! defined key_off_degrad && defined key_traldf_c2d 467 CALL iom_get( numfl_w, jpdom_data, 'soleahtw', zahtw (:,: ), jkenr ) 468 # if defined key_trcldf_eiv 383 IF( iom_varid( numfl_w, 'voddmavs', ldstop = .FALSE. ) > 0 ) THEN ! avs exist: it is used 384 CALL iom_get( numfl_w, jpdom_data, 'voddmavs', zavt (:,:,:), jkenr ) 385 ELSE ! no avs: use avt 386 CALL iom_get( numfl_w, jpdom_data, 'votkeavt', zavt (:,:,:), jkenr ) 387 ENDIF 388 389 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 469 390 CALL iom_get( numfl_w, jpdom_data, 'soleaeiw', zaeiw (:,: ), jkenr ) 470 # endif 471 #endif 472 473 #if defined key_off_degrad 391 #endif 392 393 #if defined key_degrad 474 394 CALL iom_get( numfl_u, jpdom_data, 'vozoahtu', zahtu(:,:,:), jkenr ) 475 395 CALL iom_get( numfl_v, jpdom_data, 'vomeahtv', zahtv(:,:,:), jkenr ) 476 396 CALL iom_get( numfl_w, jpdom_data, 'voveahtw', zahtw(:,:,:), jkenr ) 477 # if defined key_tr cldf_eiv397 # if defined key_traldf_eiv 478 398 CALL iom_get( numfl_u, jpdom_data, 'vozoaeiu', zaeiu(:,:,:), jkenr ) 479 399 CALL iom_get( numfl_v, jpdom_data, 'vomeaeiv', zaeiv(:,:,:), jkenr ) … … 486 406 wdta(:,:,:,2) = zw(:,:,:) * tmask(:,:,:) 487 407 488 #if defined key_trc_diatrd489 hdivdta(:,:,:,2) = zhdiv(:,:,:) * tmask(:,:,:)490 #endif491 492 408 tdta(:,:,:,2) = zt (:,:,:) * tmask(:,:,:) 493 409 sdta(:,:,:,2) = zs (:,:,:) * tmask(:,:,:) 494 410 avtdta(:,:,:,2) = zavt(:,:,:) * tmask(:,:,:) 495 411 496 #if ! defined key_off_degrad && defined key_traldf_c2d 497 ahtwdta(:,:,2) = zahtw(:,:) * tmask(:,:,1) 498 #if defined key_trcldf_eiv 412 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 499 413 aeiwdta(:,:,2) = zaeiw(:,:) * tmask(:,:,1) 500 414 #endif 501 #endif 502 503 #if defined key_off_degrad 415 416 #if defined key_degrad 504 417 ahtudta(:,:,:,2) = zahtu(:,:,:) * umask(:,:,:) 505 418 ahtvdta(:,:,:,2) = zahtv(:,:,:) * vmask(:,:,:) 506 419 ahtwdta(:,:,:,2) = zahtw(:,:,:) * tmask(:,:,:) 507 # if defined key_tr cldf_eiv420 # if defined key_traldf_eiv 508 421 aeiudta(:,:,:,2) = zaeiu(:,:,:) * umask(:,:,:) 509 422 aeivdta(:,:,:,2) = zaeiv(:,:,:) * vmask(:,:,:) … … 519 432 qsrdta (:,:,2) = zqsr(:,:) * tmask(:,:,1) 520 433 hmlddta(:,:,2) = zmld(:,:) * tmask(:,:,1) 521 522 #if defined key_trcbbl_dif || defined key_trcbbl_adv 523 bblxdta(:,:,2) = MAX( 0., zbblx(:,:) ) 524 bblydta(:,:,2) = MAX( 0., zbbly(:,:) ) 525 526 WHERE( bblxdta(:,:,2) > 2. ) bblxdta(:,:,2) = 0. 527 WHERE( bblydta(:,:,2) > 2. ) bblydta(:,:,2) = 0. 528 #endif 529 434 435 IF( l_offbbl ) THEN 436 bblxdta(:,:,2) = MAX( 0., zbblx(:,:) ) 437 bblydta(:,:,2) = MAX( 0., zbbly(:,:) ) 438 WHERE( bblxdta(:,:,2) > 2. ) bblxdta(:,:,2) = 0. 439 WHERE( bblydta(:,:,2) > 2. ) bblydta(:,:,2) = 0. 440 ENDIF 441 530 442 IF( kt == nitend ) THEN 531 443 CALL iom_close ( numfl_t ) … … 534 446 CALL iom_close ( numfl_w ) 535 447 ENDIF 536 448 ! 537 449 END SUBROUTINE dynrea 538 450 451 539 452 SUBROUTINE dta_dyn_init 540 453 !!---------------------------------------------------------------------- … … 544 457 !! 545 458 !! ** Method : 546 !! 547 !! History : 548 !! ! original : 92-01 (M. Imbard: sub domain) 549 !! ! 98-04 (L.Bopp MA Foujols: slopes for isopyc.) 550 !! ! 98-05 (L. Bopp read output of coupled run) 551 !! ! 05-03 (O. Aumont and A. El Moussaoui) F90 552 !!---------------------------------------------------------------------- 553 !! * Modules used 554 555 !! * Local declarations 556 459 !!---------------------------------------------------------------------- 557 460 REAL(wp) :: znspyr !: number of time step per year 558 559 NAMELIST/namdyn/ ndtadyn, ndtatot, nsptint, nficdyn,lperdyn, &461 !! 462 NAMELIST/namdyn/ ndtadyn, ndtatot, nsptint, lperdyn, & 560 463 & cfile_grid_T, cfile_grid_U, cfile_grid_V, cfile_grid_W 561 464 !!---------------------------------------------------------------------- … … 564 467 ! ====================================== 565 468 566 ! Read Namelist namdyn : Lateral physics on tracers 567 REWIND( numnam ) 469 REWIND( numnam ) ! Read Namelist namdyn : Lateral physics on tracers 568 470 READ ( numnam, namdyn ) 569 471 570 IF(lwp) THEN 472 IF(lwp) THEN ! control print 571 473 WRITE(numout,*) 572 474 WRITE(numout,*) 'namdyn : offline dynamical selection' … … 577 479 WRITE(numout,*) ' total number of elements in the FILE ndtatot = ' , ndtatot 578 480 WRITE(numout,*) ' type of interpolation nsptint = ' , nsptint 579 WRITE(numout,*) ' number of dynamics FILE nficdyn = ' , nficdyn580 481 WRITE(numout,*) ' loop on the same FILE lperdyn = ' , lperdyn 581 482 WRITE(numout,*) ' ' … … 586 487 WRITE(numout,*) ' ' 587 488 ENDIF 588 489 ! 589 490 znspyr = nyear_len(1) * rday / rdt 590 491 rnspdta = znspyr / FLOAT( ndtadyn ) 591 492 rnspdta2 = rnspdta * 0.5 592 493 ! 494 CALL dta_dyn( nit000 ) 495 ! 593 496 END SUBROUTINE dta_dyn_init 594 497 498 595 499 SUBROUTINE wzv( pu, pv, pw, phdiv ) 596 500 !!---------------------------------------------------------------------- … … 599 503 !! ** Purpose : Compute the now vertical velocity after the array swap 600 504 !! 601 !! ** Method : 602 !! ** Method : - Divergence: 603 !! - compute the now divergence given by : 604 !! * z-coordinate 505 !! ** Method : - compute the now divergence given by : 506 !! * z-coordinate ONLY !!!! 605 507 !! hdiv = 1/(e1t*e2t) [ di(e2u u) + dj(e1v v) ] 606 508 !! - Using the incompressibility hypothesis, the vertical 607 509 !! velocity is computed by integrating the horizontal divergence 608 510 !! from the bottom to the surface. 609 !! The boundary conditions are w=0 at the bottom (no flux) and, 610 !! in regid-lid case, w=0 at the sea surface. 611 !! 612 !! 613 !! History : 614 !! 9.0 ! 02-07 (G. Madec) Vector optimization 615 !!---------------------------------------------------------------------- 616 !! * Arguments 617 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: pu, pv !: horizontal velocities 618 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out ) :: pw !: verticla velocity 619 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( inout ) :: phdiv !: horizontal divergence 620 621 !! * Local declarations 511 !! The boundary conditions are w=0 at the bottom (no flux). 512 !!---------------------------------------------------------------------- 513 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pu, pv !: horizontal velocities 514 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out) :: pw !: verticla velocity 515 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: phdiv !: horizontal divergence 516 !! 622 517 INTEGER :: ji, jj, jk 623 518 REAL(wp) :: zu, zu1, zv, zv1, zet 624 625 519 !!---------------------------------------------------------------------- 520 ! 626 521 ! Computation of vertical velocity using horizontal divergence 627 522 phdiv(:,:,:) = 0. … … 629 524 DO jj = 2, jpjm1 630 525 DO ji = fs_2, fs_jpim1 ! vector opt. 631 #if defined key_zco632 zu = pu(ji ,jj ,jk) * umask(ji ,jj ,jk) * e2u(ji ,jj )633 zu1 = pu(ji-1,jj ,jk) * umask(ji-1,jj ,jk) * e2u(ji-1,jj )634 zv = pv(ji ,jj ,jk) * vmask(ji ,jj ,jk) * e1v(ji ,jj )635 zv1 = pv(ji ,jj-1,jk) * vmask(ji ,jj-1,jk) * e1v(ji ,jj-1)636 zet = 1. / ( e1t(ji,jj) * e2t(ji,jj) )637 #else638 526 zu = pu(ji ,jj ,jk) * umask(ji ,jj ,jk) * e2u(ji ,jj ) * fse3u(ji ,jj ,jk) 639 527 zu1 = pu(ji-1,jj ,jk) * umask(ji-1,jj ,jk) * e2u(ji-1,jj ) * fse3u(ji-1,jj ,jk) … … 641 529 zv1 = pv(ji ,jj-1,jk) * vmask(ji ,jj-1,jk) * e1v(ji ,jj-1) * fse3v(ji ,jj-1,jk) 642 530 zet = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 643 #endif644 531 phdiv(ji,jj,jk) = ( zu - zu1 + zv - zv1 ) * zet 645 532 END DO 646 533 END DO 647 ENDDO 648 649 ! Lateral boundary conditions on phdiv 650 CALL lbc_lnk( phdiv, 'T', 1. ) 651 652 534 END DO 535 CALL lbc_lnk( phdiv, 'T', 1. ) ! Lateral boundary conditions on phdiv 536 ! 653 537 ! computation of vertical velocity from the bottom 654 pw(:,:,jpk) = 0. 538 pw(:,:,jpk) = 0._wp 655 539 DO jk = jpkm1, 1, -1 656 540 pw(:,:,jk) = pw(:,:,jk+1) - fse3t(:,:,jk) * phdiv(:,:,jk) 657 541 END DO 658 542 ! 659 543 END SUBROUTINE wzv 544 545 546 SUBROUTINE dta_eiv( kt ) 547 !!---------------------------------------------------------------------- 548 !! *** ROUTINE dta_eiv *** 549 !! 550 !! ** Purpose : Compute the eddy induced velocity coefficient from the 551 !! growth rate of baroclinic instability. 552 !! 553 !! ** Method : Specific to the offline model. Computes the horizontal 554 !! values from the vertical value 555 !!---------------------------------------------------------------------- 556 INTEGER, INTENT( in ) :: kt ! ocean time-step inedx 557 !! 558 INTEGER :: ji, jj ! dummy loop indices 559 !!---------------------------------------------------------------------- 560 ! 561 IF( kt == nit000 ) THEN 562 IF(lwp) WRITE(numout,*) 563 IF(lwp) WRITE(numout,*) 'dta_eiv : eddy induced velocity coefficients' 564 IF(lwp) WRITE(numout,*) '~~~~~~~' 565 ENDIF 566 ! 567 ! Average the diffusive coefficient at u- v- points 568 DO jj = 2, jpjm1 569 DO ji = fs_2, fs_jpim1 ! vector opt. 570 aeiu(ji,jj) = .5 * ( aeiw(ji,jj) + aeiw(ji+1,jj ) ) 571 aeiv(ji,jj) = .5 * ( aeiw(ji,jj) + aeiw(ji ,jj+1) ) 572 END DO 573 END DO 574 CALL lbc_lnk( aeiu, 'U', 1. ) ; CALL lbc_lnk( aeiv, 'V', 1. ) ! lateral boundary condition 575 ! 576 END SUBROUTINE dta_eiv 577 660 578 661 579 SUBROUTINE tau2wnd( ptaux, ptauy, pwspd ) … … 667 585 !! ** Method : |tau|=rhoa*Cd*|U|^2 668 586 !!--------------------------------------------------------------------- 669 !! * Arguments 670 REAL(wp), DIMENSION(jpi,jpj), INTENT( in ) :: & 671 ptaux, ptauy !: wind stress in i-j direction resp. 672 REAL(wp), DIMENSION(jpi,jpj), INTENT( out ) :: & 673 pwspd !: wind speed 674 REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3 675 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient 676 REAL(wp) :: ztx, zty, ztau, zcoef ! temporary variables 677 INTEGER :: ji, jj ! dummy indices 587 REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: ptaux, ptauy ! wind stress in i-j direction resp. 588 REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pwspd ! wind speed 589 !! 590 REAL(wp) :: zrhoa = 1.22_wp ! Air density kg/m3 591 REAL(wp) :: zcdrag = 1.5e-3_wp ! drag coefficient 592 REAL(wp) :: ztx, zty, ztau, zcoef ! temporary variables 593 INTEGER :: ji, jj ! dummy indices 678 594 !!--------------------------------------------------------------------- 679 595 zcoef = 1. / ( zrhoa * zcdrag ) … … 689 605 END DO 690 606 CALL lbc_lnk( pwspd(:,:), 'T', 1. ) 691 607 ! 692 608 END SUBROUTINE tau2wnd 693 609 694 #if defined key_trcbbl_dif || defined key_trcbbl_adv695 696 SUBROUTINE bbl_sign( ptn, psn, pbblx, pbbly )697 !!----------------------------------------------------------------------698 !! *** ROUTINE bbl_sign ***699 !!700 !! ** Purpose : Compute the sign of local gradient of density multiplied by the slope701 !! along the bottom slope gradient : grad( rho) * grad(h)702 !! Need to compute the diffusive bottom boundary layer703 !!704 !! ** Method : When the product grad( rho) * grad(h) < 0 (where grad705 !! is an along bottom slope gradient) an additional lateral diffu-706 !! sive trend along the bottom slope is added to the general tracer707 !! trend, otherwise nothing is done. See trcbbl.F90708 !!709 !!710 !! History :711 !! 9.0 ! 02-07 (G. Madec) Vector optimization712 !!----------------------------------------------------------------------713 !! * Arguments714 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: &715 ptn , & !: temperature716 psn !: salinity717 REAL(wp), DIMENSION(jpi,jpj), INTENT( out ) :: &718 pbblx , pbbly !: sign of bbl in i-j direction resp.719 720 !! * Local declarations721 INTEGER :: ji, jj ! dummy loop indices722 INTEGER :: ik723 REAL(wp) :: &724 ztx, zsx, zhx, zalbetx, zgdrhox, & ! temporary scalars725 zty, zsy, zhy, zalbety, zgdrhoy726 REAL(wp), DIMENSION(jpi,jpj) :: &727 ztnb, zsnb, zdep728 REAL(wp) :: fsalbt, pft, pfs, pfh ! statement function729 !!----------------------------------------------------------------------730 ! ratio alpha/beta731 ! ================732 ! fsalbt: ratio of thermal over saline expension coefficients733 ! pft : potential temperature in degrees celcius734 ! pfs : salinity anomaly (s-35) in psu735 ! pfh : depth in meters736 737 fsalbt( pft, pfs, pfh ) = &738 ( ( ( -0.255019e-07 * pft + 0.298357e-05 ) * pft &739 - 0.203814e-03 ) * pft &740 + 0.170907e-01 ) * pft &741 + 0.665157e-01 &742 +(-0.678662e-05 * pfs - 0.846960e-04 * pft + 0.378110e-02 ) * pfs &743 + ( ( - 0.302285e-13 * pfh &744 - 0.251520e-11 * pfs &745 + 0.512857e-12 * pft * pft ) * pfh &746 - 0.164759e-06 * pfs &747 +( 0.791325e-08 * pft - 0.933746e-06 ) * pft &748 + 0.380374e-04 ) * pfh749 750 ! 0. 2D fields of bottom temperature and salinity, and bottom slope751 ! -----------------------------------------------------------------752 ! mbathy= number of w-level, minimum value=1 (cf domrea.F90)753 # if defined key_vectopt_loop754 jj = 1755 DO ji = 1, jpij ! vector opt. (forced unrolling)756 # else757 DO jj = 1, jpj758 DO ji = 1, jpi759 # endif760 ik = MAX( mbathy(ji,jj) - 1, 1 ) ! vertical index of the bottom ocean T-level761 ztnb(ji,jj) = ptn(ji,jj,ik) * tmask(ji,jj,1) ! masked T and S at ocean bottom762 zsnb(ji,jj) = psn(ji,jj,ik) * tmask(ji,jj,1)763 zdep(ji,jj) = fsdept(ji,jj,ik) ! depth of the ocean bottom T-level764 # if ! defined key_vectopt_loop765 END DO766 # endif767 END DO768 769 !!----------------------------------------------------------------------770 ! 1. Criteria of additional bottom diffusivity: grad(rho).grad(h)<0771 ! --------------------------------------------772 ! Sign of the local density gradient along the i- and j-slopes773 ! multiplied by the slope of the ocean bottom774 775 SELECT CASE ( neos )776 777 CASE ( 0 ) ! Jackett and McDougall (1994) formulation778 779 # if defined key_vectopt_loop780 jj = 1781 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling)782 # else783 DO jj = 1, jpjm1784 DO ji = 1, jpim1785 # endif786 ! temperature, salinity anomalie and depth787 ztx = 0.5 * ( ztnb(ji,jj) + ztnb(ji+1,jj) )788 zsx = 0.5 * ( zsnb(ji,jj) + zsnb(ji+1,jj) ) - 35.0789 zhx = 0.5 * ( zdep(ji,jj) + zdep(ji+1,jj) )790 !791 zty = 0.5 * ( ztnb(ji,jj+1) + ztnb(ji,jj) )792 zsy = 0.5 * ( zsnb(ji,jj+1) + zsnb(ji,jj) ) - 35.0793 zhy = 0.5 * ( zdep(ji,jj+1) + zdep(ji,jj) )794 ! masked ratio alpha/beta795 zalbetx = fsalbt( ztx, zsx, zhx ) * umask(ji,jj,1)796 zalbety = fsalbt( zty, zsy, zhy ) * vmask(ji,jj,1)797 ! local density gradient along i-bathymetric slope798 zgdrhox = zalbetx * ( ztnb(ji+1,jj) - ztnb(ji,jj) ) &799 - ( zsnb(ji+1,jj) - zsnb(ji,jj) )800 ! local density gradient along j-bathymetric slope801 zgdrhoy = zalbety * ( ztnb(ji,jj+1) - ztnb(ji,jj) ) &802 - ( zsnb(ji,jj+1) - zsnb(ji,jj) )803 ! sign of local i-gradient of density multiplied by the i-slope804 pbblx(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhox * ( zdep(ji+1,jj) - zdep(ji,jj) ) )805 ! sign of local j-gradient of density multiplied by the j-slope806 pbbly(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhoy * ( zdep(ji,jj+1) - zdep(ji,jj) ) )807 # if ! defined key_vectopt_loop808 END DO809 # endif810 END DO811 812 CASE ( 1 ) ! Linear formulation function of temperature only813 !814 # if defined key_vectopt_loop815 jj = 1816 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling)817 # else818 DO jj = 1, jpjm1819 DO ji = 1, jpim1820 # endif821 ! local 'density/temperature' gradient along i-bathymetric slope822 zgdrhox = ztnb(ji+1,jj) - ztnb(ji,jj)823 ! local density gradient along j-bathymetric slope824 zgdrhoy = ztnb(ji,jj+1) - ztnb(ji,jj)825 ! sign of local i-gradient of density multiplied by the i-slope826 pbblx(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhox * ( zdep(ji+1,jj) - zdep(ji,jj) ) )827 ! sign of local j-gradient of density multiplied by the j-slope828 pbbly(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhoy * ( zdep(ji,jj+1) - zdep(ji,jj) ) )829 # if ! defined key_vectopt_loop830 END DO831 # endif832 END DO833 834 CASE ( 2 ) ! Linear formulation function of temperature and salinity835 836 # if defined key_vectopt_loop837 jj = 1838 DO ji = 1, jpij-jpi ! vector opt. (forced unrolling)839 # else840 DO jj = 1, jpjm1841 DO ji = 1, jpim1842 # endif843 ! local density gradient along i-bathymetric slope844 zgdrhox = - ( rbeta*( zsnb(ji+1,jj) - zsnb(ji,jj) ) &845 - ralpha*( ztnb(ji+1,jj) - ztnb(ji,jj) ) )846 ! local density gradient along j-bathymetric slope847 zgdrhoy = - ( rbeta*( zsnb(ji,jj+1) - zsnb(ji,jj) ) &848 - ralpha*( ztnb(ji,jj+1) - ztnb(ji,jj) ) )849 ! sign of local i-gradient of density multiplied by the i-slope850 pbblx(ji,jj) = 0.5 - SIGN( 0.5, - zgdrhox * ( zdep(ji+1,jj) - zdep(ji,jj) ) )851 ! sign of local j-gradient of density multiplied by the j-slope852 pbbly(ji,jj) = 0.5 - SIGN( 0.5, -zgdrhoy * ( zdep(ji,jj+1) - zdep(ji,jj) ) )853 # if ! defined key_vectopt_loop854 END DO855 # endif856 END DO857 858 CASE DEFAULT859 860 WRITE(ctmp1,*) ' bad flag value for neos = ', neos861 CALL ctl_stop(ctmp1)862 863 END SELECT864 865 ! Lateral boundary conditions866 CALL lbc_lnk( pbblx, 'U', 1. )867 CALL lbc_lnk( pbbly, 'V', 1. )868 869 END SUBROUTINE bbl_sign870 871 #endif872 610 873 611 SUBROUTINE swap_dyn_data … … 876 614 !! 877 615 !! ** Purpose : swap array data 878 !! 879 !! History : 880 !! 9.0 ! 07-09 (C. Ethe) 881 !!---------------------------------------------------------------------- 882 883 616 !!---------------------------------------------------------------------- 617 ! 884 618 ! swap from record 2 to 1 885 619 tdta (:,:,:,1) = tdta (:,:,:,2) … … 889 623 vdta (:,:,:,1) = vdta (:,:,:,2) 890 624 wdta (:,:,:,1) = wdta (:,:,:,2) 891 #if defined key_trc_diatrd 892 hdivdta(:,:,:,1) = hdivdta(:,:,:,2) 893 #endif 894 895 #if defined key_ldfslp 625 #if defined key_ldfslp && ! defined key_c1d 896 626 uslpdta (:,:,:,1) = uslpdta (:,:,:,2) 897 627 vslpdta (:,:,:,1) = vslpdta (:,:,:,2) … … 904 634 empdta (:,:,1) = empdta (:,:,2) 905 635 qsrdta (:,:,1) = qsrdta (:,:,2) 906 907 #if ! defined key_off_degrad && defined key_traldf_c2d 908 ahtwdta(:,:,1) = ahtwdta(:,:,2) 909 # if defined key_trcldf_eiv 636 IF( l_offbbl ) THEN 637 bblxdta(:,:,1) = bblxdta(:,:,2) 638 bblydta(:,:,1) = bblydta(:,:,2) 639 ENDIF 640 641 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 910 642 aeiwdta(:,:,1) = aeiwdta(:,:,2) 911 # endif 912 #endif 913 914 #if defined key_off_degrad 643 #endif 644 645 #if defined key_degrad 915 646 ahtudta(:,:,:,1) = ahtudta(:,:,:,2) 916 647 ahtvdta(:,:,:,1) = ahtvdta(:,:,:,2) 917 648 ahtwdta(:,:,:,1) = ahtwdta(:,:,:,2) 918 # if defined key_tr cldf_eiv649 # if defined key_traldf_eiv 919 650 aeiudta(:,:,:,1) = aeiudta(:,:,:,2) 920 651 aeivdta(:,:,:,1) = aeivdta(:,:,:,2) … … 922 653 # endif 923 654 #endif 924 925 #if defined key_trcbbl_dif || defined key_trcbbl_adv 926 bblxdta(:,:,1) = bblxdta(:,:,2) 927 bblydta(:,:,1) = bblydta(:,:,2) 928 #endif 929 655 ! 930 656 END SUBROUTINE swap_dyn_data 657 931 658 932 659 SUBROUTINE assign_dyn_data … … 939 666 !!---------------------------------------------------------------------- 940 667 941 t n (:,:,:) = tdta (:,:,:,2)942 sn (:,:,:) = sdta (:,:,:,2)943 avt(:,:,:) = avtdta(:,:,:,2)668 tsn(:,:,:,jp_tem) = tdta (:,:,:,2) 669 tsn(:,:,:,jp_sal) = sdta (:,:,:,2) 670 avt(:,:,:) = avtdta(:,:,:,2) 944 671 945 672 un (:,:,:) = udta (:,:,:,2) 946 673 vn (:,:,:) = vdta (:,:,:,2) 947 674 wn (:,:,:) = wdta (:,:,:,2) 948 949 #if defined key_trc_diatrd 950 hdivn(:,:,:) = hdivdta(:,:,:,2) 951 #endif 952 953 #if defined key_zdfddm 954 avs(:,:,:) = avtdta (:,:,:,2) 955 #endif 956 957 958 #if defined key_ldfslp 675 676 #if defined key_ldfslp && ! defined key_c1d 959 677 uslp (:,:,:) = uslpdta (:,:,:,2) 960 678 vslp (:,:,:) = vslpdta (:,:,:,2) … … 969 687 emps(:,:) = emp(:,:) 970 688 qsr (:,:) = qsrdta (:,:,2) 971 972 #if ! defined key_off_degrad && defined key_traldf_c2d 973 ahtw(:,:) = ahtwdta(:,:,2) 974 # if defined key_trcldf_eiv 689 IF( l_offbbl ) THEN 690 ahu_bbl(:,:) = bblxdta(:,:,2) 691 ahv_bbl(:,:) = bblydta(:,:,2) 692 ENDIF 693 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 975 694 aeiw(:,:) = aeiwdta(:,:,2) 976 # endif 977 #endif 978 979 #if defined key_off_degrad 695 #endif 696 697 #if defined key_degrad 980 698 ahtu(:,:,:) = ahtudta(:,:,:,2) 981 699 ahtv(:,:,:) = ahtvdta(:,:,:,2) 982 700 ahtw(:,:,:) = ahtwdta(:,:,:,2) 983 # if defined key_tr cldf_eiv701 # if defined key_traldf_eiv 984 702 aeiu(:,:,:) = aeiudta(:,:,:,2) 985 703 aeiv(:,:,:) = aeivdta(:,:,:,2) 986 704 aeiw(:,:,:) = aeiwdta(:,:,:,2) 987 705 # endif 988 989 #endif 990 991 #if defined key_trcbbl_dif || defined key_trcbbl_adv 992 bblx(:,:) = bblxdta(:,:,2) 993 bbly(:,:) = bblydta(:,:,2) 994 #endif 995 706 #endif 707 ! 996 708 END SUBROUTINE assign_dyn_data 997 709 710 998 711 SUBROUTINE linear_interp_dyn_data( pweigh ) 999 712 !!---------------------------------------------------------------------- 1000 !! 713 !! *** ROUTINE linear_interp_dyn_data *** 1001 714 !! 1002 715 !! ** Purpose : linear interpolation of data 1003 !! 1004 !!---------------------------------------------------------------------- 1005 !! * Argument 1006 REAL(wp), INTENT( in ) :: pweigh ! weigh 1007 1008 !! * Local declarations 716 !!---------------------------------------------------------------------- 717 REAL(wp), INTENT(in) :: pweigh ! weigh 718 !! 1009 719 REAL(wp) :: zweighm1 1010 720 !!---------------------------------------------------------------------- … … 1012 722 zweighm1 = 1. - pweigh 1013 723 1014 t n (:,:,:) = zweighm1 * tdta (:,:,:,1) + pweigh * tdta (:,:,:,2)1015 sn (:,:,:) = zweighm1 * sdta (:,:,:,1) + pweigh * sdta (:,:,:,2)1016 avt(:,:,:) = zweighm1 * avtdta(:,:,:,1) + pweigh * avtdta(:,:,:,2)724 tsn(:,:,:,jp_tem) = zweighm1 * tdta (:,:,:,1) + pweigh * tdta (:,:,:,2) 725 tsn(:,:,:,jp_sal) = zweighm1 * sdta (:,:,:,1) + pweigh * sdta (:,:,:,2) 726 avt(:,:,:) = zweighm1 * avtdta(:,:,:,1) + pweigh * avtdta(:,:,:,2) 1017 727 1018 728 un (:,:,:) = zweighm1 * udta (:,:,:,1) + pweigh * udta (:,:,:,2) 1019 729 vn (:,:,:) = zweighm1 * vdta (:,:,:,1) + pweigh * vdta (:,:,:,2) 1020 730 wn (:,:,:) = zweighm1 * wdta (:,:,:,1) + pweigh * wdta (:,:,:,2) 1021 1022 #if defined key_trc_diatrd 1023 hdivn(:,:,:) = zweighm1 * hdivdta(:,:,:,1) + pweigh * hdivdta(:,:,:,2) 1024 #endif 1025 1026 #if defined key_zdfddm 1027 avs(:,:,:) = zweighm1 * avtdta (:,:,:,1) + pweigh * avtdta (:,:,:,2) 1028 #endif 1029 1030 1031 #if defined key_ldfslp 731 732 #if defined key_ldfslp && ! defined key_c1d 1032 733 uslp (:,:,:) = zweighm1 * uslpdta (:,:,:,1) + pweigh * uslpdta (:,:,:,2) 1033 734 vslp (:,:,:) = zweighm1 * vslpdta (:,:,:,1) + pweigh * vslpdta (:,:,:,2) … … 1042 743 emps(:,:) = emp(:,:) 1043 744 qsr (:,:) = zweighm1 * qsrdta (:,:,1) + pweigh * qsrdta (:,:,2) 1044 1045 #if ! defined key_off_degrad && defined key_traldf_c2d 1046 ahtw(:,:) = zweighm1 * ahtwdta(:,:,1) + pweigh * ahtwdta(:,:,2) 1047 # if defined key_trcldf_eiv 745 IF( l_offbbl ) THEN 746 ahu_bbl(:,:) = zweighm1 * bblxdta(:,:,1) + pweigh * bblxdta(:,:,2) 747 ahv_bbl(:,:) = zweighm1 * bblydta(:,:,1) + pweigh * bblydta(:,:,2) 748 ENDIF 749 750 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 1048 751 aeiw(:,:) = zweighm1 * aeiwdta(:,:,1) + pweigh * aeiwdta(:,:,2) 1049 # endif 1050 #endif 1051 1052 #if defined key_off_degrad 752 #endif 753 754 #if defined key_degrad 1053 755 ahtu(:,:,:) = zweighm1 * ahtudta(:,:,:,1) + pweigh * ahtudta(:,:,:,2) 1054 756 ahtv(:,:,:) = zweighm1 * ahtvdta(:,:,:,1) + pweigh * ahtvdta(:,:,:,2) 1055 757 ahtw(:,:,:) = zweighm1 * ahtwdta(:,:,:,1) + pweigh * ahtwdta(:,:,:,2) 1056 # if defined key_tr cldf_eiv758 # if defined key_traldf_eiv 1057 759 aeiu(:,:,:) = zweighm1 * aeiudta(:,:,:,1) + pweigh * aeiudta(:,:,:,2) 1058 760 aeiv(:,:,:) = zweighm1 * aeivdta(:,:,:,1) + pweigh * aeivdta(:,:,:,2) … … 1060 762 # endif 1061 763 #endif 1062 1063 #if defined key_trcbbl_dif || defined key_trcbbl_adv 1064 bblx(:,:) = zweighm1 * bblxdta(:,:,1) + pweigh * bblxdta(:,:,2) 1065 bbly(:,:) = zweighm1 * bblydta(:,:,1) + pweigh * bblydta(:,:,2) 1066 #endif 1067 764 ! 1068 765 END SUBROUTINE linear_interp_dyn_data 1069 766 767 !!====================================================================== 1070 768 END MODULE dtadyn
Note: See TracChangeset
for help on using the changeset viewer.
