Changeset 7590 for branches/UKMO
- Timestamp:
- 2017-01-20T16:51:19+01:00 (7 years ago)
- Location:
- branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 10 edited
-
BDY/bdytides.F90 (modified) (3 diffs)
-
DIA/dia25h.F90 (modified) (9 diffs)
-
DIA/diatmb.F90 (modified) (2 diffs)
-
DIA/diawri.F90 (modified) (3 diffs)
-
IOM/iom_def.F90 (modified) (1 diff)
-
IOM/restart.F90 (modified) (3 diffs)
-
ZDF/zdfmxl.F90 (modified) (6 diffs)
-
nemogcm.F90 (modified) (3 diffs)
-
step.F90 (modified) (5 diffs)
-
step_oce.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90
r6332 r7590 50 50 REAL(wp), POINTER, DIMENSION(:,:,:) :: v !: Tidal constituents : V (after nodal cor.) 51 51 END TYPE TIDES_DATA 52 INTEGER, PUBLIC, PARAMETER :: jptides_max = 15 !: Max number of tidal contituents 53 LOGICAL, PUBLIC :: ln_harm_ana_store !: =T Stores data for harmonic Analysis 54 LOGICAL, PUBLIC :: ln_harm_ana_compute !: =T Compute harmonic Analysis 55 LOGICAL, PUBLIC :: ln_harmana_read !: =T Decide to do the analysis 56 !from scratch or continue previous run 52 57 53 58 !$AGRIF_DO_NOT_TREAT … … 90 95 TYPE(MAP_POINTER), DIMENSION(jpbgrd) :: ibmap_ptr !: array of pointers to nbmap 91 96 !! 92 NAMELIST/nambdy_tide/filtide, ln_bdytide_2ddta, ln_bdytide_conj 97 NAMELIST/nambdy_tide/filtide, ln_bdytide_2ddta, ln_bdytide_conj, ln_harm_ana_store, ln_harm_ana_compute, ln_harmana_read 93 98 !!---------------------------------------------------------------------- 94 99 … … 126 131 IF(lwp) WRITE(numout,*) ' assume complex conjugate : ', ln_bdytide_conj 127 132 IF(lwp) WRITE(numout,*) ' Number of tidal components to read: ', nb_harmo 133 IF(lwp) WRITE(numout,*) ' Use PCOMS harmonic ananalysis or not: ', ln_harm_ana_store 134 IF(lwp) WRITE(numout,*) ' Compute Final harmonic ananalysis or not: ', ln_harm_ana_compute 135 IF(lwp) WRITE(numout,*) ' Read in previous days harmonic data or start afresh: ', ln_harmana_read 128 136 IF(lwp) THEN 129 137 WRITE(numout,*) ' Tidal components: ' -
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/DIA/dia25h.F90
r6332 r7590 20 20 USE zdf_oce, ONLY: en 21 21 #endif 22 USE diatmb 22 23 23 24 IMPLICIT NONE … … 30 31 !! * variables for calculating 25-hourly means 31 32 REAL(wp),SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: tn_25h , sn_25h, rinsitu_t_25h 32 REAL(wp),SAVE, ALLOCATABLE, DIMENSION(:,:) :: sshn_25h 33 REAL(wp),SAVE, ALLOCATABLE, DIMENSION(:,:) :: sshn_25h, insitu_bot_25h, temp_bot_25h 33 34 REAL(wp),SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: un_25h , vn_25h , wn_25h 34 35 REAL(wp),SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: avt_25h , avm_25h … … 63 64 INTEGER :: ios ! Local integer output status for namelist read 64 65 INTEGER :: ierror ! Local integer for memory allocation 66 REAL(wp), DIMENSION(jpi,jpj,3) :: zwtmb ! temporary workspace 65 67 ! 66 68 NAMELIST/nam_dia25h/ ln_dia25h … … 98 100 IF( ierror > 0 ) THEN 99 101 CALL ctl_stop( 'dia_25h: unable to allocate rinsitu_t_25h' ) ; RETURN 102 ENDIF 103 ALLOCATE( insitu_bot_25h(jpi,jpj), STAT=ierror ) 104 IF( ierror > 0 ) THEN 105 CALL ctl_stop( 'dia_tide: unable to allocate insitu_bot_25h' ) ; RETURN 106 ENDIF 107 ALLOCATE( temp_bot_25h(jpi,jpj), STAT=ierror ) 108 IF( ierror > 0 ) THEN 109 CALL ctl_stop( 'dia_tide: unable to allocate temp_bot_25h' ) ; RETURN 100 110 ENDIF 101 111 ALLOCATE( un_25h(jpi,jpj,jpk), STAT=ierror ) … … 143 153 CALL theta2t 144 154 rinsitu_t_25h(:,:,:) = rinsitu_t(:,:,:) 155 CALL dia_calctmb( rinsitu_t(:,:,:),zwtmb ) 156 insitu_bot_25h(:,:) = zwtmb(:,:,3) 157 CALL dia_calctmb( tn_25h(:,:,:),zwtmb ) 158 temp_bot_25h(:,:) = zwtmb(:,:,3) 145 159 sshn_25h(:,:) = sshb(:,:) 146 160 un_25h(:,:,:) = ub(:,:,:) … … 237 251 CALL theta2t 238 252 rinsitu_t_25h(:,:,:) = rinsitu_t_25h(:,:,:) + rinsitu_t(:,:,:) 253 CALL dia_calctmb( rinsitu_t(:,:,:),zwtmb ) 254 insitu_bot_25h(:,:) = insitu_bot_25h(:,:) + zwtmb(:,:,3) 255 zw3d(:,:,:) = tsn(:,:,:,jp_tem) 256 CALL dia_calctmb( zw3d,zwtmb ) 257 temp_bot_25h(:,:) = temp_bot_25h(:,:) + zwtmb(:,:,3) 239 258 sshn_25h(:,:) = sshn_25h(:,:) + sshn (:,:) 240 259 un_25h(:,:,:) = un_25h(:,:,:) + un(:,:,:) … … 268 287 sn_25h(:,:,:) = sn_25h(:,:,:) / 25.0_wp 269 288 rinsitu_t_25h(:,:,:) = rinsitu_t_25h(:,:,:) / 25.0_wp 289 insitu_bot_25h(:,:) = insitu_bot_25h(:,:) / 25.0_wp 290 temp_bot_25h(:,:) = temp_bot_25h(:,:) /25.0_wp 270 291 sshn_25h(:,:) = sshn_25h(:,:) / 25.0_wp 271 292 un_25h(:,:,:) = un_25h(:,:,:) / 25.0_wp … … 289 310 zw3d(:,:,:) = rinsitu_t_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 290 311 CALL iom_put("tempis25h", zw3d) ! in-situ temperature 312 zw2d(:,:) = insitu_bot_25h(:,:)*tmask(:,:,1) + zmdi*(1.0-tmask(:,:,1)) 313 CALL iom_put("tempisbot25h", zw2d) ! bottom in-situ temperature 314 zw2d(:,:) = temp_bot_25h(:,:)*tmask(:,:,1) + zmdi*(1.0-tmask(:,:,1)) 315 CALL iom_put("temperbot25h",zw2d) ! bottom potential temperature 291 316 zw3d(:,:,:) = sn_25h(:,:,:)*tmask(:,:,:) + zmdi*(1.0-tmask(:,:,:)) 292 317 CALL iom_put( "salin25h", zw3d ) ! salinity … … 321 346 CALL theta2t 322 347 rinsitu_t_25h(:,:,:) = rinsitu_t(:,:,:) 348 CALL dia_calctmb( rinsitu_t(:,:,:),zwtmb ) 349 insitu_bot_25h(:,:) = zwtmb(:,:,3) 350 CALL dia_calctmb( tn_25h(:,:,:),zwtmb) 351 temp_bot_25h(:,:) = zwtmb(:,:,3) 323 352 sshn_25h(:,:) = sshn (:,:) 324 353 un_25h(:,:,:) = un(:,:,:) -
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/DIA/diatmb.F90
r6332 r7590 19 19 PUBLIC dia_tmb_init ! routine called by nemogcm.F90 20 20 PUBLIC dia_tmb ! routine called by diawri.F90 21 PUBLIC dia_calctmb ! routine called by dia25h.F90 21 22 22 23 !!---------------------------------------------------------------------- … … 108 109 DO jj = 1,jpj 109 110 DO ji = 1,jpi 110 jk = max(1,mbathy(ji,jj) - 1)111 jk = max(1,mbathy(ji,jj)) 111 112 pouttmb(ji,jj,3) = pinfield(ji,jj,jk)*tmask(ji,jj,jk) + zmdi*(1.0-tmask(ji,jj,jk)) 112 113 END DO -
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90
r6332 r7590 46 46 USE diatmb ! Top,middle,bottom output 47 47 USE dia25h ! 25h Mean output 48 USE diaopfoam ! Diaopfoam output 48 49 USE iom 49 50 USE ioipsl … … 189 190 zztmpy = ( bfrva(ji, jj) * vn(ji,jj ,mbkv(ji,jj )) & 190 191 & + bfrva(ji,jj-1) * vn(ji,jj-1,mbkv(ji,jj-1)) ) 191 z2d(ji,jj) = rau0 * SQRT( zztmpx * zztmpx + zztmpy * zztmpy ) * tmask(ji,jj,1)192 z2d(ji,jj) = 0.5_wp * rau0 * SQRT( zztmpx * zztmpx + zztmpy * zztmpy ) * tmask(ji,jj,1) 192 193 ! 193 194 ENDDO … … 389 390 IF (ln_dia25h) THEN 390 391 CALL dia_25h( kt ) 392 ENDIF 393 IF (ln_diaopfoam) THEN 394 CALL dia_diaopfoam 391 395 ENDIF 392 396 ! -
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/IOM/iom_def.F90
r6331 r7590 43 43 INTEGER, PARAMETER, PUBLIC :: jp_i1 = 204 !: write INTEGER(1) 44 44 45 INTEGER, PARAMETER, PUBLIC :: jpmax_files = 100 !: maximum number of simultaneously opened file45 INTEGER, PARAMETER, PUBLIC :: jpmax_files = 200 !: maximum number of simultaneously opened file 46 46 INTEGER, PARAMETER, PUBLIC :: jpmax_vars = 600 !: maximum number of variables in one file 47 47 INTEGER, PARAMETER, PUBLIC :: jpmax_dims = 4 !: maximum number of dimensions for one variable -
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/IOM/restart.F90
r6332 r7590 25 25 USE trdmxl_oce ! ocean active mixed layer tracers trends variables 26 26 USE divcur ! hor. divergence and curl (div & cur routines) 27 USE diurnal_bulk ! diurnal SST 27 28 28 29 IMPLICIT NONE … … 153 154 CALL iom_rstput( kt, nitrst, numrow, 'rhd' , rhd ) 154 155 #endif 156 157 IF (ln_diurnal) CALL iom_rstput( kt, nitrst, numrow, 'Dsst', x_dsst ) 158 155 159 IF( kt == nitrst ) THEN 156 160 CALL iom_close( numrow ) ! close the restart file (only at last time step) … … 265 269 ENDIF 266 270 #endif 271 272 ! Diurnal SST 273 IF( ln_diurnal ) CALL iom_get( numror, jpdom_autoglo, 'Dsst' , x_dsst ) 274 267 275 ! 268 276 IF( neuler == 0 ) THEN ! Euler restart (neuler=0) -
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90
r6331 r7590 22 22 USE timing ! Timing 23 23 USE trc_oce, ONLY : lk_offline ! offline flag 24 USE lbclnk 25 USE eosbn2 ! Equation of state 24 26 25 27 IMPLICIT NONE … … 27 29 28 30 PUBLIC zdf_mxl ! called by step.F90 31 PUBLIC zdf_mxl_tref ! called by asminc.F90 29 32 PUBLIC zdf_mxl_alloc ! Used in zdf_tke_init 30 33 … … 33 36 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hmlp !: mixed layer depth (rho=rho0+zdcrit) [m] 34 37 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hmlpt !: mixed layer depth at t-points [m] 38 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hmld_tref !: mixed layer depth at t-points - temperature criterion [m] 39 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hmld_kara !: mixed layer depth of Kara et al [m] 35 40 36 41 REAL(wp), PUBLIC :: rho_c = 0.01_wp !: density criterion for mixed layer depth 37 42 REAL(wp) :: avt_c = 5.e-4_wp ! Kz criterion for the turbocline depth 43 44 ! Namelist variables for namzdf_karaml 45 46 LOGICAL :: ln_kara ! Logical switch for calculating kara mixed 47 ! layer 48 LOGICAL :: ln_kara_write25h ! Logical switch for 25 hour outputs 49 INTEGER :: jpmld_type ! mixed layer type 50 REAL(wp) :: ppz_ref ! depth of initial T_ref 51 REAL(wp) :: ppdT_crit ! Critical temp diff 52 REAL(wp) :: ppiso_frac ! Fraction of ppdT_crit used 53 54 !Used for 25h mean 55 LOGICAL, PRIVATE :: kara_25h_init = .TRUE. !Logical used to initalise 25h 56 !outputs. Necissary, because we need to 57 !initalise the kara_25h on the zeroth 58 !timestep (i.e in the nemogcm_init call) 59 REAL, PRIVATE, ALLOCATABLE, DIMENSION(:,:) :: hmld_kara_25h 38 60 39 61 !! * Substitutions … … 52 74 zdf_mxl_alloc = 0 ! set to zero if no array to be allocated 53 75 IF( .NOT. ALLOCATED( nmln ) ) THEN 54 ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), STAT= zdf_mxl_alloc ) 76 ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), & 77 & hmld_tref(jpi,jpj), STAT= zdf_mxl_alloc ) 55 78 ! 56 79 IF( lk_mpp ) CALL mpp_sum ( zdf_mxl_alloc ) … … 123 146 END DO 124 147 ! depth of the mixing and mixed layers 148 149 CALL zdf_mxl_kara( kt ) ! kara MLD 150 125 151 DO jj = 1, jpj 126 152 DO ji = 1, jpi … … 146 172 END SUBROUTINE zdf_mxl 147 173 174 175 SUBROUTINE zdf_mxl_tref() 176 !!---------------------------------------------------------------------- 177 !! *** ROUTINE zdf_mxl_tref *** 178 !! 179 !! ** Purpose : Compute the mixed layer depth with temperature criteria. 180 !! 181 !! ** Method : The temperature-defined mixed layer depth is required 182 !! when assimilating SST in a 2D analysis. 183 !! 184 !! ** Action : hmld_tref 185 !!---------------------------------------------------------------------- 186 ! 187 INTEGER :: ji, jj, jk ! dummy loop indices 188 REAL(wp) :: t_ref ! Reference temperature 189 REAL(wp) :: temp_c = 0.2 ! temperature criterion for mixed layer depth 190 !!---------------------------------------------------------------------- 191 ! 192 ! Initialise array 193 IF( zdf_mxl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'zdf_mxl_tref : unable to allocate arrays' ) 194 195 !For the AMM model assimiation uses a temperature based mixed layer depth 196 !This is defined here 197 DO jj = 1, jpj 198 DO ji = 1, jpi 199 hmld_tref(ji,jj)=fsdept(ji,jj,1 ) 200 IF(ssmask(ji,jj) > 0.)THEN 201 t_ref=tsn(ji,jj,1,jp_tem) 202 DO jk=2,jpk 203 IF(ssmask(ji,jj)==0.)THEN 204 hmld_tref(ji,jj)=fsdept(ji,jj,jk ) 205 EXIT 206 ELSEIF( ABS(tsn(ji,jj,jk,jp_tem)-t_ref) < temp_c)THEN 207 hmld_tref(ji,jj)=fsdept(ji,jj,jk ) 208 ELSE 209 EXIT 210 ENDIF 211 ENDDO 212 ENDIF 213 ENDDO 214 ENDDO 215 216 END SUBROUTINE zdf_mxl_tref 217 218 SUBROUTINE zdf_mxl_kara( kt ) 219 !!---------------------------------------------------------------------------------- 220 !! *** ROUTINE zdf_mxl_kara *** 221 ! 222 ! Calculate mixed layer depth according to the definition of 223 ! Kara et al, 2000, JGR, 105, 16803. 224 ! 225 ! If mld_type=1 the mixed layer depth is calculated as the depth at which the 226 ! density has increased by an amount equivalent to a temperature difference of 227 ! 0.8C at the surface. 228 ! 229 ! For other values of mld_type the mixed layer is calculated as the depth at 230 ! which the temperature differs by 0.8C from the surface temperature. 231 ! 232 ! Original version: David Acreman 233 ! 234 !!----------------------------------------------------------------------------------- 235 236 INTEGER, INTENT( in ) :: kt ! ocean time-step index 237 238 NAMELIST/namzdf_karaml/ ln_kara,jpmld_type, ppz_ref, ppdT_crit, & 239 & ppiso_frac, ln_kara_write25h 240 241 ! Local variables 242 REAL, DIMENSION(jpi,jpk) :: ppzdep ! depth for use in calculating d(rho) 243 REAL(wp), DIMENSION(jpi,jpj,jpts) :: ztsn_2d !Local version of tsn 244 LOGICAL :: ll_found(jpi,jpj) ! Is T_b to be found by interpolation ? 245 LOGICAL :: ll_belowml(jpi,jpj,jpk) ! Flag points below mixed layer when ll_found=F 246 INTEGER :: ji, jj, jk ! loop counter 247 INTEGER :: ik_ref(jpi,jpj) ! index of reference level 248 INTEGER :: ik_iso(jpi,jpj) ! index of last uniform temp level 249 REAL :: zT(jpi,jpj,jpk) ! Temperature or denisty 250 REAL :: zT_ref(jpi,jpj) ! reference temperature 251 REAL :: zT_b ! base temperature 252 REAL :: zdTdz(jpi,jpj,jpk-2) ! gradient of zT 253 REAL :: zmoddT(jpi,jpj,jpk-2) ! Absolute temperature difference 254 REAL :: zdz ! depth difference 255 REAL :: zdT ! temperature difference 256 REAL :: zdelta_T(jpi,jpj) ! difference critereon 257 REAL :: zRHO1(jpi,jpj), zRHO2(jpi,jpj) ! Densities 258 INTEGER, SAVE :: idt, i_steps 259 INTEGER, SAVE :: i_cnt_25h ! Counter for 25 hour means 260 INTEGER :: ios ! Local integer output status for namelist read 261 262 263 264 !!------------------------------------------------------------------------------------- 265 266 IF( kt == nit000 ) THEN 267 ! Default values 268 ln_kara = .FALSE. 269 ln_kara_write25h = .FALSE. 270 jpmld_type = 1 271 ppz_ref = 10.0 272 ppdT_crit = 0.2 273 ppiso_frac = 0.1 274 ! Read namelist 275 REWIND ( numnam_ref ) 276 READ ( numnam_ref, namzdf_karaml, IOSTAT=ios, ERR= 901 ) 277 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_karaml in reference namelist', lwp ) 278 279 REWIND( numnam_cfg ) ! Namelist nam_diadiaop in configuration namelist 3D hourly diagnostics 280 READ ( numnam_cfg, namzdf_karaml, IOSTAT = ios, ERR = 902 ) 281 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_karaml in configuration namelist', lwp ) 282 IF(lwm) WRITE ( numond, namzdf_karaml ) 283 284 285 WRITE(numout,*) '===== Kara mixed layer =====' 286 WRITE(numout,*) 'ln_kara = ', ln_kara 287 WRITE(numout,*) 'jpmld_type = ', jpmld_type 288 WRITE(numout,*) 'ppz_ref = ', ppz_ref 289 WRITE(numout,*) 'ppdT_crit = ', ppdT_crit 290 WRITE(numout,*) 'ppiso_frac = ', ppiso_frac 291 WRITE(numout,*) 'ln_kara_write25h = ', ln_kara_write25h 292 WRITE(numout,*) '============================' 293 294 IF ( .NOT.ln_kara ) THEN 295 WRITE(numout,*) "ln_kara not set; Kara mixed layer not calculated" 296 ELSE 297 IF (.NOT. ALLOCATED(hmld_kara) ) ALLOCATE( hmld_kara(jpi,jpj) ) 298 IF ( ln_kara_write25h .AND. kara_25h_init ) THEN 299 i_cnt_25h = 0 300 IF (.NOT. ALLOCATED(hmld_kara_25h) ) & 301 & ALLOCATE( hmld_kara_25h(jpi,jpj) ) 302 hmld_kara_25h = 0._wp 303 IF( nacc == 1 ) THEN 304 idt = INT(rdtmin) 305 ELSE 306 idt = INT(rdt) 307 ENDIF 308 IF( MOD( 3600,idt) == 0 ) THEN 309 i_steps = 3600 / idt 310 ELSE 311 CALL ctl_stop('STOP', & 312 & 'zdf_mxl_kara: timestep must give MOD(3600,rdt)'// & 313 & ' = 0 otherwise no hourly values are possible') 314 ENDIF 315 ENDIF 316 ENDIF 317 ENDIF 318 319 IF ( ln_kara ) THEN 320 321 !set critical ln_kara 322 ztsn_2d = tsn(:,:,1,:) 323 ztsn_2d(:,:,jp_tem) = ztsn_2d(:,:,jp_tem) + ppdT_crit 324 325 ! Set the mixed layer depth criterion at each grid point 326 ppzdep = 0._wp 327 IF( jpmld_type == 1 ) THEN 328 CALL eos ( tsn(:,:,1,:), & 329 & ppzdep(:,:), zRHO1(:,:) ) 330 CALL eos ( ztsn_2d(:,:,:), & 331 & ppzdep(:,:), zRHO2(:,:) ) 332 zdelta_T(:,:) = abs( zRHO1(:,:) - zRHO2(:,:) ) * rau0 333 ! RHO from eos (2d version) doesn't calculate north or east halo: 334 CALL lbc_lnk( zdelta_T, 'T', 1. ) 335 zT(:,:,:) = rhop(:,:,:) 336 ELSE 337 zdelta_T(:,:) = ppdT_crit 338 zT(:,:,:) = tsn(:,:,:,jp_tem) 339 ENDIF 340 341 ! Calculate the gradient of zT and absolute difference for use later 342 DO jk = 1 ,jpk-2 343 zdTdz(:,:,jk) = ( zT(:,:,jk+1) - zT(:,:,jk) ) / fse3w(:,:,jk+1) 344 zmoddT(:,:,jk) = abs( zT(:,:,jk+1) - zT(:,:,jk) ) 345 END DO 346 347 ! Find density/temperature at the reference level (Kara et al use 10m). 348 ! ik_ref is the index of the box centre immediately above or at the reference level 349 ! Find ppz_ref in the array of model level depths and find the ref 350 ! density/temperature by linear interpolation. 351 ik_ref = -1 352 DO jk = jpkm1, 2, -1 353 WHERE( fsdept(:,:,jk) > ppz_ref ) 354 ik_ref(:,:) = jk - 1 355 zT_ref(:,:) = zT(:,:,jk-1) + & 356 & zdTdz(:,:,jk-1) * ( ppz_ref - fsdept(:,:,jk-1) ) 357 ENDWHERE 358 END DO 359 IF ( ANY( ik_ref < 0 ) .OR. ANY( ik_ref > jpkm1 ) ) THEN 360 CALL ctl_stop( "STOP", & 361 & "zdf_mxl_kara: unable to find reference level for kara ML" ) 362 ELSE 363 ! If the first grid box centre is below the reference level then use the 364 ! top model level to get zT_ref 365 WHERE( fsdept(:,:,1) > ppz_ref ) 366 zT_ref = zT(:,:,1) 367 ik_ref = 1 368 ENDWHERE 369 370 ! Search for a uniform density/temperature region where adjacent levels 371 ! differ by less than ppiso_frac * deltaT. 372 ! ik_iso is the index of the last level in the uniform layer 373 ! ll_found indicates whether the mixed layer depth can be found by interpolation 374 ik_iso(:,:) = ik_ref(:,:) 375 ll_found(:,:) = .false. 376 DO jj = 1, nlcj 377 DO ji = 1, nlci 378 !CDIR NOVECTOR 379 DO jk = ik_ref(ji,jj), mbathy(ji,jj)-1 380 IF( zmoddT(ji,jj,jk) > ( ppiso_frac * zdelta_T(ji,jj) ) ) THEN 381 ik_iso(ji,jj) = jk 382 ll_found(ji,jj) = ( zmoddT(ji,jj,jk) > zdelta_T(ji,jj) ) 383 EXIT 384 ENDIF 385 END DO 386 END DO 387 END DO 388 389 ! Use linear interpolation to find depth of mixed layer base where possible 390 hmld_kara(:,:) = ppz_ref 391 DO jj = 1, jpj 392 DO ji = 1, jpi 393 IF( ll_found(ji,jj) .and. tmask(ji,jj,1) == 1.0 ) THEN 394 zdz = abs( zdelta_T(ji,jj) / zdTdz(ji,jj,ik_iso(ji,jj)) ) 395 hmld_kara(ji,jj) = fsdept(ji,jj,ik_iso(ji,jj)) + zdz 396 ENDIF 397 END DO 398 END DO 399 400 ! If ll_found = .false. then calculate MLD using difference of zdelta_T 401 ! from the reference density/temperature 402 403 ! Prevent this section from working on land points 404 WHERE( tmask(:,:,1) /= 1.0 ) 405 ll_found = .true. 406 ENDWHERE 407 408 DO jk = 1, jpk 409 ll_belowml(:,:,jk) = abs( zT(:,:,jk) - zT_ref(:,:) ) >= & 410 & zdelta_T(:,:) 411 END DO 412 413 ! Set default value where interpolation cannot be used (ll_found=false) 414 DO jj = 1, jpj 415 DO ji = 1, jpi 416 IF( .NOT. ll_found(ji,jj) ) & 417 & hmld_kara(ji,jj) = fsdept(ji,jj,mbathy(ji,jj)) 418 END DO 419 END DO 420 421 DO jj = 1, jpj 422 DO ji = 1, jpi 423 !CDIR NOVECTOR 424 DO jk = ik_ref(ji,jj)+1, mbathy(ji,jj) 425 IF( ll_found(ji,jj) ) EXIT 426 IF( ll_belowml(ji,jj,jk) ) THEN 427 zT_b = zT_ref(ji,jj) + zdelta_T(ji,jj) * & 428 & SIGN(1.0, zdTdz(ji,jj,jk-1) ) 429 zdT = zT_b - zT(ji,jj,jk-1) 430 zdz = zdT / zdTdz(ji,jj,jk-1) 431 hmld_kara(ji,jj) = fsdept(ji,jj,jk-1) + zdz 432 EXIT 433 ENDIF 434 END DO 435 END DO 436 END DO 437 438 hmld_kara(:,:) = hmld_kara(:,:) * tmask(:,:,1) 439 440 IF( ln_kara_write25h ) THEN 441 !Double IF required as i_steps not defined if ln_kara_write25h = 442 ! FALSE 443 IF ( ( MOD( kt, i_steps ) == 0 ) .OR. kara_25h_init ) THEN 444 hmld_kara_25h = hmld_kara_25h + hmld_kara 445 i_cnt_25h = i_cnt_25h + 1 446 IF ( kara_25h_init ) kara_25h_init = .FALSE. 447 ENDIF 448 ENDIF 449 450 #if defined key_iomput 451 IF( kt /= nit000 ) THEN 452 CALL iom_put( "mldkara" , hmld_kara ) 453 IF( ( MOD( i_cnt_25h, 25) == 0 ) .AND. ln_kara_write25h ) & 454 CALL iom_put( "kara25h" , ( hmld_kara_25h / 25._wp ) ) 455 ENDIF 456 #endif 457 458 ENDIF 459 ENDIF 460 461 END SUBROUTINE zdf_mxl_kara 462 148 463 !!====================================================================== 149 464 END MODULE zdfmxl -
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90
r6332 r7590 87 87 USE diatmb ! Top,middle,bottom output 88 88 USE dia25h ! 25h mean output 89 USE diaopfoam ! 3D hourly output 90 USE diurnal_bulk ! diurnal bulk SST 89 91 90 92 IMPLICIT NONE … … 405 407 CALL istate_init ! ocean initial state (Dynamics and tracers) 406 408 409 CALL diurnal_sst_bulk_init ! diurnal sst 410 IF ( ln_diurnal ) CALL diurnal_sst_coolskin_init ! cool skin 411 407 412 IF( lk_tide ) CALL tide_init( nit000 ) ! Initialisation of the tidal harmonics 408 413 … … 479 484 CALL dia_tmb_init ! TMB outputs 480 485 CALL dia_25h_init ! 25h mean outputs 486 CALL dia_diaopfoam_init ! 25h mean outputs 481 487 ! 482 488 END SUBROUTINE nemo_init -
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/step.F90
r6610 r7590 225 225 ENDIF 226 226 227 ! Cool skin 228 IF ( ln_diurnal ) THEN 229 IF ( ln_blk_core ) THEN 230 CALL diurnal_sst_coolskin_step( & 231 qns(:,:)+(rn_abs*qsr(:,:)), taum, rhop(:,:,1), rdt) 232 ELSE 233 CALL diurnal_sst_coolskin_step( & 234 qns, taum, rhop(:,:,1), rdt) 235 ENDIF 236 ENDIF 237 227 238 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 228 239 ! diagnostics and outputs (ua, va, tsa used as workspace) … … 237 248 ! 238 249 IF( ln_crs ) CALL crs_fld( kstp ) ! ocean model: online field coarsening & output 250 251 !Diurnal warm layer model 252 IF ( ln_diurnal ) THEN 253 IF ( ln_blk_core ) THEN 254 IF( kstp == nit000 )THEN 255 CALL diurnal_sst_takaya_step( & 256 & qsr(:,:)-(rn_abs*qsr(:,:)), qns(:,:)+(rn_abs*qsr(:,:)), & 257 & taum, rhop(:,:,1), & 258 & rdt, ld_calcfrac = .TRUE.) 259 ELSE 260 CALL diurnal_sst_takaya_step( & 261 & qsr(:,:)-(rn_abs*qsr(:,:)), qns(:,:)+(rn_abs*qsr(:,:)), & 262 & taum, rhop(:,:,1), rdt ) 263 ENDIF 264 ELSE 265 IF( kstp == nit000 )THEN 266 CALL diurnal_sst_takaya_step( & 267 & qsr, qns, taum, rhop(:,:,1), & 268 & rdt, ld_calcfrac = .TRUE.) 269 ELSE 270 CALL diurnal_sst_takaya_step( & 271 & qsr, qns, taum, rhop(:,:,1), rdt ) 272 ENDIF 273 ENDIF 274 ENDIF 239 275 240 276 #if defined key_top … … 338 374 IF( lk_vvl ) CALL dom_vvl_sf_swp( kstp ) ! swap of vertical scale factors 339 375 ! 340 IF( lrst_oce ) CALL rst_write( kstp ) ! write output ocean restart file341 376 342 377 #if defined key_agrif … … 362 397 CALL dia_wri_state( 'output.abort', kstp ) 363 398 ENDIF 399 IF( ln_harm_ana_store ) CALL harm_ana( kstp ) ! Harmonic analysis of tides 364 400 IF( kstp == nit000 ) THEN 365 401 CALL iom_close( numror ) ! close input ocean restart file … … 367 403 IF( lwm.AND.numoni /= -1 ) CALL FLUSH ( numoni ) ! flush output namelist ice 368 404 ENDIF 405 IF( lrst_oce ) CALL rst_write( kstp ) ! write output ocean restart file 369 406 370 407 !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> -
branches/UKMO/CO6_KD490_amm7_oper/NEMOGCM/NEMO/OPA_SRC/step_oce.F90
r6610 r7590 108 108 USE prtctl ! Print control (prt_ctl routine) 109 109 110 USE harmonic_analysis ! harmonic analysis of tides (harm_ana routine) 111 USE bdytides ! harmonic analysis of tides (harm_ana routine) 110 112 USE diaobs ! Observation operator 113 114 USE diurnal_bulk ! diurnal SST bulk routines (diurnal_sst_takaya routine) 115 USE cool_skin ! diurnal cool skin correction (diurnal_sst_coolskin routine) 111 116 112 117 USE timing ! Timing
Note: See TracChangeset
for help on using the changeset viewer.
