Changeset 2416
- Timestamp:
- 2010-11-21T12:42:54+01:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/TRA/tradmp.F90
r2287 r2416 17 17 #if defined key_tradmp || defined key_esopa 18 18 !!---------------------------------------------------------------------- 19 !! key_tradmpinternal damping19 !! 'key_tradmp' internal damping 20 20 !!---------------------------------------------------------------------- 21 !! tra_dmp : update the tracer trend with the internal damping22 !! tra_dmp_init : initialization, namlist read, parameters control23 !! dtacof_zoom : restoring coefficient for zoom domain24 !! dtacof : restoring coefficient for global domain25 !! cofdis : compute the distance to the coastline21 !! tra_dmp : update the tracer trend with the internal damping 22 !! tra_dmp_init : initialization, namlist read, parameters control 23 !! dtacof_zoom : restoring coefficient for zoom domain 24 !! dtacof : restoring coefficient for global domain 25 !! cofdis : compute the distance to the coastline 26 26 !!---------------------------------------------------------------------- 27 USE oce ! ocean dynamics and tracersvariables28 USE dom_oce ! ocean space and timedomain variables29 USE trdmod_oce ! ocean space and time domainvariables30 USE trdtra ! ocean space and time domain variables31 USE zdf_oce ! oceanvertical physics32 USE phycst ! Define parameters for the routines33 USE dtatem ! temperature data34 USE dtasal ! salinity data35 USE zdfmxl !mixed layer depth36 USE in_out_manager 37 USE lib_mpp ! distribued memory computing38 USE prtctl 27 USE oce ! ocean: variables 28 USE dom_oce ! ocean: domain variables 29 USE trdmod_oce ! ocean: trend variables 30 USE trdtra ! active tracers: trends 31 USE zdf_oce ! ocean: vertical physics 32 USE phycst ! physical constants 33 USE dtatem ! data: temperature 34 USE dtasal ! data: salinity 35 USE zdfmxl ! vertical physics: mixed layer depth 36 USE in_out_manager ! I/O manager 37 USE lib_mpp ! MPP library 38 USE prtctl ! Print control 39 39 40 40 IMPLICIT NONE … … 43 43 PUBLIC tra_dmp ! routine called by step.F90 44 44 PUBLIC tra_dmp_init ! routine called by opa.F90 45 PUBLIC dtacof ! routine called by tradmp.F90 and trcdmp.F9046 PUBLIC dtacof_zoom ! routine called by tradmp.F90 and trcdmp.F9045 PUBLIC dtacof ! routine called by in both tradmp.F90 and trcdmp.F90 46 PUBLIC dtacof_zoom ! routine called by in both tradmp.F90 and trcdmp.F90 47 47 48 48 #if ! defined key_agrif … … 52 52 #endif 53 53 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: strdmp !: damping salinity trend (psu/s) 54 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: ttrdmp !: damping temperature trend (Ce ntigrade/s)54 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: ttrdmp !: damping temperature trend (Celcius/s) 55 55 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: resto !: restoring coeff. on T and S (s-1) 56 56 57 ! !!* Namelist namtra_dmp : T & S newtonian damping *58 INTEGER :: nn_hdmp = -1 ! = 0/-1/'latitude' for damping over T and S59 INTEGER :: nn_zdmp = 0 ! = 0/1/2 flag for damping in the mixed layer60 REAL(wp) :: rn_surf = 50. ! surface time scale for internal damping [days]61 REAL(wp) :: rn_bot = 360. ! bottom time scale for internal damping [days]62 REAL(wp) :: rn_dep = 800. ! depth of transition between rn_surf and rn_bot [meters]63 INTEGER :: nn_file = 2 ! = 1 create a damping.coeff NetCDF file57 ! !!* Namelist namtra_dmp : T & S newtonian damping * 58 INTEGER :: nn_hdmp = -1 ! = 0/-1/'latitude' for damping over T and S 59 INTEGER :: nn_zdmp = 0 ! = 0/1/2 flag for damping in the mixed layer 60 REAL(wp) :: rn_surf = 50._wp ! surface time scale for internal damping [days] 61 REAL(wp) :: rn_bot = 360._wp ! bottom time scale for internal damping [days] 62 REAL(wp) :: rn_dep = 800._wp ! depth of transition between rn_surf and rn_bot [meters] 63 INTEGER :: nn_file = 2 ! = 1 create a damping.coeff NetCDF file 64 64 65 65 !! * Substitutions … … 69 69 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 70 70 !! $Id$ 71 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)71 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 72 72 !!---------------------------------------------------------------------- 73 74 73 CONTAINS 75 74 … … 94 93 INTEGER, INTENT(in) :: kt ! ocean time-step index 95 94 !! 96 REAL(wp) :: zta, zsa ! temporary scalars 97 INTEGER :: ji, jj, jk ! dummy loop indices 98 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: ztrdt, ztrds 99 !!---------------------------------------------------------------------- 100 101 IF( l_trdtra ) THEN !* Save ta and sa trends 102 ALLOCATE( ztrdt(jpi,jpj,jpk) ) ; ztrdt(:,:,:) = tsa(:,:,:,jp_tem) 103 ALLOCATE( ztrds(jpi,jpj,jpk) ) ; ztrds(:,:,:) = tsa(:,:,:,jp_sal) 104 ENDIF 105 106 SELECT CASE ( nn_zdmp ) 95 INTEGER :: ji, jj, jk ! dummy loop indices 96 REAL(wp) :: zta, zsa ! local scalars 97 !!---------------------------------------------------------------------- 98 ! 99 SELECT CASE ( nn_zdmp ) !== type of damping ==! 107 100 ! 108 101 CASE( 0 ) !== newtonian damping throughout the water column ==! … … 114 107 tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) + zta 115 108 tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) + zsa 116 ! save the salinity trend (used in asmtrj) 117 strdmp(ji,jj,jk) = zsa 109 strdmp(ji,jj,jk) = zsa ! save the salinity trend (used in asmtrj) 118 110 ttrdmp(ji,jj,jk) = zta 119 111 END DO … … 125 117 DO jj = 2, jpjm1 126 118 DO ji = fs_2, fs_jpim1 ! vector opt. 127 IF( avt(ji,jj,jk) <= 5.e-4 ) THEN119 IF( avt(ji,jj,jk) <= 5.e-4_wp ) THEN 128 120 zta = resto(ji,jj,jk) * ( t_dta(ji,jj,jk) - tsb(ji,jj,jk,jp_tem) ) 129 121 zsa = resto(ji,jj,jk) * ( s_dta(ji,jj,jk) - tsb(ji,jj,jk,jp_sal) ) 130 122 ELSE 131 zta = 0. e0132 zsa = 0. e0123 zta = 0._wp 124 zsa = 0._wp 133 125 ENDIF 134 126 tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) + zta 135 127 tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) + zsa 136 ! save the salinity trend (used in asmtrj) 137 strdmp(ji,jj,jk) = zsa 128 strdmp(ji,jj,jk) = zsa ! save the salinity trend (used in asmtrj) 138 129 ttrdmp(ji,jj,jk) = zta 139 130 END DO … … 149 140 zsa = resto(ji,jj,jk) * ( s_dta(ji,jj,jk) - tsb(ji,jj,jk,jp_sal) ) 150 141 ELSE 151 zta = 0. e0152 zsa = 0. e0142 zta = 0._wp 143 zsa = 0._wp 153 144 ENDIF 154 ! add the trends to the general tracer trends155 145 tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) + zta 156 146 tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) + zsa 157 ! save the salinity trend (used in asmtrj) 158 strdmp(ji,jj,jk) = zsa 147 strdmp(ji,jj,jk) = zsa ! save the salinity trend (used in asmtrj) 159 148 ttrdmp(ji,jj,jk) = zta 160 149 END DO … … 163 152 ! 164 153 END SELECT 165 154 ! 166 155 IF( l_trdtra ) THEN ! trend diagnostic 167 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:) 168 ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:) 169 CALL trd_tra( kt, 'TRA', jp_tem, jptra_trd_dmp, ztrdt ) 170 CALL trd_tra( kt, 'TRA', jp_sal, jptra_trd_dmp, ztrds ) 171 DEALLOCATE( ztrdt ) ; DEALLOCATE( ztrds ) 156 CALL trd_tra( kt, 'TRA', jp_tem, jptra_trd_dmp, ttrdmp ) 157 CALL trd_tra( kt, 'TRA', jp_sal, jptra_trd_dmp, strdmp ) 172 158 ENDIF 173 159 ! ! Control print … … 191 177 REWIND ( numnam ) ! Read Namelist namtra_dmp : temperature and salinity damping term 192 178 READ ( numnam, namtra_dmp ) 193 IF( lzoom ) nn_zdmp = 0 ! restoring to climatology at closed north or south boundaries 179 180 IF( lzoom ) nn_zdmp = 0 ! restoring to climatology at closed north or south boundaries 194 181 195 182 IF(lwp) THEN ! Namelist print … … 226 213 & CALL ctl_stop( 'no temperature and/or salinity data define key_dtatem and key_dtasal' ) 227 214 228 strdmp(:,:,:) = 0. e0! internal damping salinity trend (used in asmtrj)229 ttrdmp(:,:,:) = 0. e0215 strdmp(:,:,:) = 0._wp ! internal damping salinity trend (used in asmtrj) 216 ttrdmp(:,:,:) = 0._wp 230 217 ! ! Damping coefficients initialization 231 218 IF( lzoom ) THEN ; CALL dtacof_zoom( resto ) … … 250 237 !! ** Action : - resto, the damping coeff. for T and S 251 238 !!---------------------------------------------------------------------- 252 !! * Arguments 253 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: presto !: restoring coeff. (s-1) 254 ! 255 INTEGER :: ji, jj, jk, jn ! dummy loop indices 256 REAL(wp) :: zlat, zlat0, zlat1, zlat2 ! temporary scalar 257 REAL(wp), DIMENSION(6) :: zfact ! temporary workspace 258 !!---------------------------------------------------------------------- 259 260 zfact(1) = 1. 261 zfact(2) = 1. 262 zfact(3) = 11./12. 263 zfact(4) = 8./12. 264 zfact(5) = 4./12. 265 zfact(6) = 1./12. 266 zfact(:) = zfact(:) / ( 5. * rday ) ! 5 days max restoring time scale 267 268 presto(:,:,:) = 0.e0 239 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: presto ! restoring coeff. (s-1) 240 ! 241 INTEGER :: ji, jj, jk, jn ! dummy loop indices 242 REAL(wp) :: zlat, zlat0, zlat1, zlat2, z1_5d ! local scalar 243 REAL(wp), DIMENSION(6) :: zfact ! 1Dworkspace 244 !!---------------------------------------------------------------------- 245 246 zfact(1) = 1._wp 247 zfact(2) = 1._wp 248 zfact(3) = 11._wp / 12._wp 249 zfact(4) = 8._wp / 12._wp 250 zfact(5) = 4._wp / 12._wp 251 zfact(6) = 1._wp / 12._wp 252 zfact(:) = zfact(:) / ( 5._wp * rday ) ! 5 days max restoring time scale 253 254 presto(:,:,:) = 0._wp 269 255 270 256 ! damping along the forced closed boundary over 6 grid-points … … 285 271 ! 286 272 ! ! Initialization : 287 presto(:,:,:) = 0.e0 288 zlat0 = 10. ! zlat0 : latitude strip where resto decreases 289 zlat1 = 30. ! zlat1 : resto = 1 before zlat1 290 zlat2 = zlat1 + zlat0 ! zlat2 : resto decreases from 1 to 0 between zlat1 and zlat2 273 presto(:,:,:) = 0._wp 274 zlat0 = 10._wp ! zlat0 : latitude strip where resto decreases 275 zlat1 = 30._wp ! zlat1 : resto = 1 before zlat1 276 zlat2 = zlat1 + zlat0 ! zlat2 : resto decreases from 1 to 0 between zlat1 and zlat2 277 z1_5d = 1._wp / ( 5._wp * rday ) ! z1_5d : 1 / 5days 291 278 292 279 DO jk = 2, jpkm1 ! Compute arrays resto ; value for internal damping : 5 days … … 295 282 zlat = ABS( gphit(ji,jj) ) 296 283 IF( zlat1 <= zlat .AND. zlat <= zlat2 ) THEN 297 presto(ji,jj,jk) = 0.5 * ( 1./(5.*rday) ) * ( 1. - cos(rpi*(zlat2-zlat)/zlat0))284 presto(ji,jj,jk) = 0.5_wp * z1_5d * ( 1._wp - COS( rpi*(zlat2-zlat)/zlat0 ) ) 298 285 ELSEIF( zlat < zlat1 ) THEN 299 presto(ji,jj,jk) = 1./(5.*rday)286 presto(ji,jj,jk) = z1_5d 300 287 ENDIF 301 288 END DO … … 325 312 USE iom 326 313 USE ioipsl 327 !! * Arguments328 INTEGER , INTENT(in ) :: kn_hdmp ! :damping option329 REAL(wp) , INTENT(in ) :: pn_surf ! :surface time scale (days)330 REAL(wp) , INTENT(in ) :: pn_bot ! :bottom time scale (days)331 REAL(wp) , INTENT(in ) :: pn_dep ! :depth of transition (meters)332 INTEGER , INTENT(in ) :: kn_file ! :save the damping coef on a file or not314 !! 315 INTEGER , INTENT(in ) :: kn_hdmp ! damping option 316 REAL(wp) , INTENT(in ) :: pn_surf ! surface time scale (days) 317 REAL(wp) , INTENT(in ) :: pn_bot ! bottom time scale (days) 318 REAL(wp) , INTENT(in ) :: pn_dep ! depth of transition (meters) 319 INTEGER , INTENT(in ) :: kn_file ! save the damping coef on a file or not 333 320 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 334 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: presto !: restoring coeff. (s-1) 335 ! 336 INTEGER :: ji, jj, jk ! dummy loop indices 337 INTEGER :: ii0, ii1, ij0, ij1 ! - - 338 INTEGER :: inum0 ! logical unit for file restoring damping term 339 INTEGER :: icot ! logical unit for file distance to the coast 340 REAL(wp) :: zinfl, zlon ! temporary scalars 341 REAL(wp) :: zlat, zlat0, zlat1, zlat2 ! - - 342 REAL(wp) :: zsdmp, zbdmp ! - - 321 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: presto ! restoring coeff. (s-1) 322 ! 323 INTEGER :: ji, jj, jk ! dummy loop indices 324 INTEGER :: ii0, ii1, ij0, ij1 ! local integers 325 INTEGER :: inum0, icot ! - - 326 REAL(wp) :: zinfl, zlon ! local scalars 327 REAL(wp) :: zlat, zlat0, zlat1, zlat2 ! - - 328 REAL(wp) :: zsdmp, zbdmp ! - - 343 329 REAL(wp), DIMENSION(jpk) :: zhfac 344 330 REAL(wp), DIMENSION(jpi,jpj) :: zmrs … … 347 333 !!---------------------------------------------------------------------- 348 334 349 ! ====================================350 ! ORCA configuration : global domain351 ! ====================================352 335 ! ! ==================== 336 ! ! ORCA configuration : global domain 337 ! ! ==================== 338 ! 353 339 IF(lwp) WRITE(numout,*) 354 340 IF(lwp) WRITE(numout,*) ' dtacof : Global domain of ORCA' 355 341 IF(lwp) WRITE(numout,*) ' ------------------------------' 356 357 ! ... Initialization : 358 presto(:,:,:) = 0.e0 342 ! 343 presto(:,:,:) = 0._wp 359 344 ! 360 345 IF( kn_hdmp > 0 ) THEN ! Damping poleward of 'nn_hdmp' degrees ! 361 346 ! !-----------------------------------------! 362 347 IF(lwp) WRITE(numout,*) 363 IF(lwp) WRITE(numout,*) ' Damping poleward of ', kn_hdmp, ' deg.'348 IF(lwp) WRITE(numout,*) ' Damping poleward of ', kn_hdmp, ' deg.' 364 349 ! 365 350 CALL iom_open ( 'dist.coast.nc', icot, ldstop = .FALSE. ) … … 373 358 374 359 ! ! Compute arrays resto 375 zinfl = 1000.e3 376 zlat0 = 10. 360 zinfl = 1000.e3_wp ! distance of influence for damping term 361 zlat0 = 10._wp ! latitude strip where resto decreases 377 362 zlat1 = REAL( kn_hdmp ) ! resto = 0 between -zlat1 and zlat1 378 363 zlat2 = zlat1 + zlat0 ! resto increases from 0 to 1 between |zlat1| and |zlat2| … … 382 367 zlat = ABS( gphit(ji,jj) ) 383 368 IF ( zlat1 <= zlat .AND. zlat <= zlat2 ) THEN 384 presto(ji,jj,1) = 0.5 * ( 1. - cos(rpi*(zlat-zlat1)/zlat0 ))369 presto(ji,jj,1) = 0.5_wp * ( 1._wp - COS( rpi*(zlat-zlat1)/zlat0 ) ) 385 370 ELSEIF ( zlat > zlat2 ) THEN 386 presto(ji,jj,1) = 1. 371 presto(ji,jj,1) = 1._wp 387 372 ENDIF 388 373 END DO … … 393 378 DO ji = 1, jpi 394 379 zlat = gphit(ji,jj) 395 zlon = MOD( glamt(ji,jj), 360. )396 IF ( zlat1 < zlat .AND. zlat < zlat2 .AND. 45. < zlon .AND. zlon < 100.) THEN397 presto(ji,jj,1) = 0. e0380 zlon = MOD( glamt(ji,jj), 360._wp ) 381 IF ( zlat1 < zlat .AND. zlat < zlat2 .AND. 45._wp < zlon .AND. zlon < 100._wp ) THEN 382 presto(ji,jj,1) = 0._wp 398 383 ENDIF 399 384 END DO … … 401 386 ENDIF 402 387 403 zsdmp = 1. /(pn_surf * rday)404 zbdmp = 1. /(pn_bot * rday)388 zsdmp = 1._wp / ( pn_surf * rday ) 389 zbdmp = 1._wp / ( pn_bot * rday ) 405 390 DO jk = 2, jpkm1 406 391 DO jj = 1, jpj … … 408 393 zdct(ji,jj,jk) = MIN( zinfl, zdct(ji,jj,jk) ) 409 394 ! ... Decrease the value in the vicinity of the coast 410 presto(ji,jj,jk) = presto(ji,jj,1 ) * 0.5 * ( 1. - COS( rpi*zdct(ji,jj,jk)/zinfl))395 presto(ji,jj,jk) = presto(ji,jj,1 ) * 0.5_wp * ( 1._wp - COS( rpi*zdct(ji,jj,jk)/zinfl) ) 411 396 ! ... Vertical variation from zsdmp (sea surface) to zbdmp (bottom) 412 presto(ji,jj,jk) = presto(ji,jj,jk) * ( zbdmp + (zsdmp-zbdmp)*EXP(-fsdept(ji,jj,jk)/pn_dep))397 presto(ji,jj,jk) = presto(ji,jj,jk) * ( zbdmp + (zsdmp-zbdmp) * EXP(-fsdept(ji,jj,jk)/pn_dep) ) 413 398 END DO 414 399 END DO … … 417 402 ENDIF 418 403 419 404 ! ! ========================= 405 ! ! Med and Red Sea damping (ORCA configuration only) 406 ! ! ========================= 420 407 IF( cp_cfg == "orca" .AND. ( kn_hdmp > 0 .OR. kn_hdmp == -1 ) ) THEN 421 422 ! ! =========================423 ! ! Med and Red Sea damping424 ! ! =========================425 408 IF(lwp)WRITE(numout,*) 426 409 IF(lwp)WRITE(numout,*) ' ORCA configuration: Damping in Med and Red Seas' 427 428 429 zmrs(:,:) = 0.e0 ! damping term on the Med or Red Sea 430 410 ! 411 zmrs(:,:) = 0._wp 412 ! 431 413 SELECT CASE ( jp_cfg ) 432 414 ! ! ======================= 433 415 CASE ( 4 ) ! ORCA_R4 configuration 434 416 ! ! ======================= 435 ! Mediterranean Sea436 ij0 = 50 ; ij1 = 56 437 ii0 = 81 ; ii1 = 91 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0417 ij0 = 50 ; ij1 = 56 ! Mediterranean Sea 418 419 ii0 = 81 ; ii1 = 91 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 438 420 ij0 = 50 ; ij1 = 55 439 ii0 = 75 ; ii1 = 80 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0421 ii0 = 75 ; ii1 = 80 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 440 422 ij0 = 52 ; ij1 = 53 441 ii0 = 70 ; ii1 = 74 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0423 ii0 = 70 ; ii1 = 74 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 442 424 ! Smooth transition from 0 at surface to 1./rday at the 18th level in Med and Red Sea 443 425 DO jk = 1, 17 444 zhfac (jk) = 0.5 *( 1.- COS( rpi*(jk-1)/16. )) / rday426 zhfac (jk) = 0.5_wp * ( 1._wp - COS( rpi * REAL(jk-1,wp) / 16._wp ) ) / rday 445 427 END DO 446 428 DO jk = 18, jpkm1 447 zhfac (jk) = 1. /rday429 zhfac (jk) = 1._wp / rday 448 430 END DO 449 431 ! ! ======================= 450 432 CASE ( 2 ) ! ORCA_R2 configuration 451 433 ! ! ======================= 452 ! Mediterranean Sea 453 ij0 = 96 ; ij1 = 110 454 ii0 = 157 ; ii1 = 181 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1.e0 434 ij0 = 96 ; ij1 = 110 ! Mediterranean Sea 435 ii0 = 157 ; ii1 = 181 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 455 436 ij0 = 100 ; ij1 = 110 456 ii0 = 144 ; ii1 = 156 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0437 ii0 = 144 ; ii1 = 156 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 457 438 ij0 = 100 ; ij1 = 103 458 ii0 = 139 ; ii1 = 143 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0459 ! Decrease before Gibraltar Strait460 ij0 = 101 ; ij1 = 102 461 ii0 = 139 ; ii1 = 141 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0. e0462 ii0 = 142 ; ii1 = 142 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0 / 90.e0463 ii0 = 143 ; ii1 = 143 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.40 e0464 ii0 = 144 ; ii1 = 144 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.75 e0465 ! Red Sea466 ij0 = 87 ; ij1 = 96 467 ii0 = 147 ; ii1 = 163 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0468 ! Decrease before Bab el Mandeb Strait469 ij0 = 91 ; ij1 = 91 470 ii0 = 153 ; ii1 = 160 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.80 e0439 ii0 = 139 ; ii1 = 143 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 440 ! 441 ij0 = 101 ; ij1 = 102 ! Decrease before Gibraltar Strait 442 ii0 = 139 ; ii1 = 141 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0._wp 443 ii0 = 142 ; ii1 = 142 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp / 90._wp 444 ii0 = 143 ; ii1 = 143 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.40_wp 445 ii0 = 144 ; ii1 = 144 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.75_wp 446 ! 447 ij0 = 87 ; ij1 = 96 ! Red Sea 448 ii0 = 147 ; ii1 = 163 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 449 ! 450 ij0 = 91 ; ij1 = 91 ! Decrease before Bab el Mandeb Strait 451 ii0 = 153 ; ii1 = 160 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.80_wp 471 452 ij0 = 90 ; ij1 = 90 472 ii0 = 153 ; ii1 = 160 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.40 e0453 ii0 = 153 ; ii1 = 160 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.40_wp 473 454 ij0 = 89 ; ij1 = 89 474 ii0 = 158 ; ii1 = 160 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0 / 90.e0455 ii0 = 158 ; ii1 = 160 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp / 90._wp 475 456 ij0 = 88 ; ij1 = 88 476 ii0 = 160 ; ii1 = 163 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0. e0457 ii0 = 160 ; ii1 = 163 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0._wp 477 458 ! Smooth transition from 0 at surface to 1./rday at the 18th level in Med and Red Sea 478 459 DO jk = 1, 17 479 zhfac (jk) = 0.5 *( 1.- COS( rpi*(jk-1)/16. )) / rday460 zhfac (jk) = 0.5_wp * ( 1._wp - COS( rpi * REAL(jk-1,wp) / 16._wp ) ) / rday 480 461 END DO 481 462 DO jk = 18, jpkm1 482 zhfac (jk) = 1. /rday463 zhfac (jk) = 1._wp / rday 483 464 END DO 484 465 ! ! ======================= 485 466 CASE ( 05 ) ! ORCA_R05 configuration 486 467 ! ! ======================= 487 ! Mediterranean Sea 488 ii0 = 568 ; ii1 = 574 489 ij0 = 324 ; ij1 = 333 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1.e0 468 ii0 = 568 ; ii1 = 574 ! Mediterranean Sea 469 ij0 = 324 ; ij1 = 333 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 490 470 ii0 = 575 ; ii1 = 658 491 ij0 = 314 ; ij1 = 366 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0492 ! Black Sea (remaining part493 ii0 = 641 ; ii1 = 651 494 ij0 = 367 ; ij1 = 372 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0495 ! Decrease before Gibraltar Strait496 ij0 = 324 ; ij1 = 333 497 ii0 = 565 ; ii1 = 565 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0 / 90.e0498 ii0 = 566 ; ii1 = 566 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.40 499 ii0 = 567 ; ii1 = 567 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.75 500 ! Red Sea501 ii0 = 641 ; ii1 = 665 502 ij0 = 270 ; ij1 = 310 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1. e0503 ! Decrease before Bab el Mandeb Strait504 ii0 = 666 ; ii1 = 675 471 ij0 = 314 ; ij1 = 366 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 472 ! 473 ii0 = 641 ; ii1 = 651 ! Black Sea (remaining part 474 ij0 = 367 ; ij1 = 372 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 475 ! 476 ij0 = 324 ; ij1 = 333 ! Decrease before Gibraltar Strait 477 ii0 = 565 ; ii1 = 565 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp / 90._wp 478 ii0 = 566 ; ii1 = 566 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.40_wp 479 ii0 = 567 ; ii1 = 567 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.75_wp 480 ! 481 ii0 = 641 ; ii1 = 665 ! Red Sea 482 ij0 = 270 ; ij1 = 310 ; zmrs( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 1._wp 483 ! 484 ii0 = 666 ; ii1 = 675 ! Decrease before Bab el Mandeb Strait 505 485 ij0 = 270 ; ij1 = 290 506 486 DO ji = mi0(ii0), mi1(ii1) 507 zmrs( ji , mj0(ij0):mj1(ij1) ) = 0.1 * ABS( FLOAT(ji - mi1(ii1)) )487 zmrs( ji , mj0(ij0):mj1(ij1) ) = 0.1_wp * ABS( FLOAT(ji - mi1(ii1)) ) 508 488 END DO 509 zsdmp = 1. /(pn_surf * rday)510 zbdmp = 1. /(pn_bot * rday)489 zsdmp = 1._wp / ( pn_surf * rday ) 490 zbdmp = 1._wp / ( pn_bot * rday ) 511 491 DO jk = 1, jpk 512 zhfac (jk) = ( zbdmp + (zsdmp-zbdmp) * EXP(-fsdept(1,1,jk)/pn_dep))492 zhfac(jk) = ( zbdmp + (zsdmp-zbdmp) * EXP( -fsdept(1,1,jk)/pn_dep ) ) 513 493 END DO 514 494 ! ! ======================== … … 520 500 521 501 DO jk = 1, jpkm1 522 presto(:,:,jk) = zmrs(:,:) * zhfac(jk) + ( 1. - zmrs(:,:) ) * presto(:,:,jk)502 presto(:,:,jk) = zmrs(:,:) * zhfac(jk) + ( 1._wp - zmrs(:,:) ) * presto(:,:,jk) 523 503 END DO 524 504 525 505 ! Mask resto array and set to 0 first and last levels 526 506 presto(:,:, : ) = presto(:,:,:) * tmask(:,:,:) 527 presto(:,:, 1 ) = 0. e0528 presto(:,:,jpk) = 0. e0507 presto(:,:, 1 ) = 0._wp 508 presto(:,:,jpk) = 0._wp 529 509 ! !--------------------! 530 510 ELSE ! No damping ! … … 594 574 & ' create the "dist.coast.nc" file using IDL' ) 595 575 596 pdct(:,:,:) = 0. e0597 zxt(:,:) = cos( rad * gphit(:,:) ) * cos( rad * glamt(:,:) )598 zyt(:,:) = cos( rad * gphit(:,:) ) * sin( rad * glamt(:,:) )599 zzt(:,:) = sin( rad * gphit(:,:) )576 pdct(:,:,:) = 0._wp 577 zxt(:,:) = COS( rad * gphit(:,:) ) * COS( rad * glamt(:,:) ) 578 zyt(:,:) = COS( rad * gphit(:,:) ) * SIN( rad * glamt(:,:) ) 579 zzt(:,:) = SIN( rad * gphit(:,:) ) 600 580 601 581 … … 610 590 zmask(ji,jj) = ( tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) & 611 591 & + tmask(ji,jj ,jk) + tmask(ji+1,jj ,jk) ) 612 llcotu(ji,jj) = ( tmask(ji,jj, jk) + tmask(ji+1,jj ,jk) == 1. )613 llcotv(ji,jj) = ( tmask(ji,jj ,jk) + tmask(ji ,jj+1,jk) == 1. )614 llcotf(ji,jj) = ( zmask(ji,jj) > 0. ) .AND. ( zmask(ji,jj) < 4.)592 llcotu(ji,jj) = ( tmask(ji,jj, jk) + tmask(ji+1,jj ,jk) == 1._wp ) 593 llcotv(ji,jj) = ( tmask(ji,jj ,jk) + tmask(ji ,jj+1,jk) == 1._wp ) 594 llcotf(ji,jj) = ( zmask(ji,jj) > 0._wp ) .AND. ( zmask(ji,jj) < 4._wp ) 615 595 END DO 616 596 END DO … … 678 658 ! Compute cartesian coordinates of coastline points 679 659 ! and the number of coastline points 680 681 660 icoast = 0 682 661 DO jj = 1, jpj … … 704 683 705 684 ! Distance for the T-points 706 707 685 DO jj = 1, jpj 708 686 DO ji = 1, jpi 709 IF( tmask(ji,jj,jk) == 0. ) THEN710 pdct(ji,jj,jk) = 0. 687 IF( tmask(ji,jj,jk) == 0._wp ) THEN 688 pdct(ji,jj,jk) = 0._wp 711 689 ELSE 712 690 DO jl = 1, icoast … … 727 705 ! ---------------------------------------------------------- 728 706 clname = 'dist.coast' 729 itime = 0730 CALL ymds2ju( 0 , 1 , 1 , 0. e0, zdate0 )707 itime = 0 708 CALL ymds2ju( 0 , 1 , 1 , 0._wp , zdate0 ) 731 709 CALL restini( 'NONE', jpi , jpj , glamt, gphit , & 732 710 & jpk , gdept_0, clname, itime, zdate0, & … … 734 712 CALL restput( icot, 'Tcoast', jpi, jpj, jpk, 0, pdct ) 735 713 CALL restclo( icot ) 736 714 ! 737 715 END SUBROUTINE cofdis 738 716
Note: See TracChangeset
for help on using the changeset viewer.