Changeset 13613
- Timestamp:
- 2020-10-15T15:42:25+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/ICE/icedyn_adv_pra.F90
r13609 r13613 288 288 CALL adv_x( zdt , zudy , 0._wp , zarea , z0vl , sxvl , sxxvl , syvl , syyvl , sxyvl ) 289 289 ENDIF 290 ENDIF 291 ! 290 ENDIF 291 ! 292 ENDIF 293 294 ! --- Lateral boundary conditions --- ! 295 ! caution: for gradients (sx and sy) the sign changes 296 CALL lbc_lnk_multi( 'icedyn_adv_pra', z0ice , 'T', 1._wp, sxice , 'T', -1._wp, syice , 'T', -1._wp & ! ice volume 297 & , sxxice, 'T', 1._wp, syyice, 'T', 1._wp, sxyice, 'T', 1._wp & 298 & , z0snw , 'T', 1._wp, sxsn , 'T', -1._wp, sysn , 'T', -1._wp & ! snw volume 299 & , sxxsn , 'T', 1._wp, syysn , 'T', 1._wp, sxysn , 'T', 1._wp ) 300 CALL lbc_lnk_multi( 'icedyn_adv_pra', z0smi , 'T', 1._wp, sxsal , 'T', -1._wp, sysal , 'T', -1._wp & ! ice salinity 301 & , sxxsal, 'T', 1._wp, syysal, 'T', 1._wp, sxysal, 'T', 1._wp & 302 & , z0ai , 'T', 1._wp, sxa , 'T', -1._wp, sya , 'T', -1._wp & ! ice concentration 303 & , sxxa , 'T', 1._wp, syya , 'T', 1._wp, sxya , 'T', 1._wp ) 304 CALL lbc_lnk_multi( 'icedyn_adv_pra', z0oi , 'T', 1._wp, sxage , 'T', -1._wp, syage , 'T', -1._wp & ! ice age 305 & , sxxage, 'T', 1._wp, syyage, 'T', 1._wp, sxyage, 'T', 1._wp ) 306 CALL lbc_lnk_multi( 'icedyn_adv_pra', z0es , 'T', 1._wp, sxc0 , 'T', -1._wp, syc0 , 'T', -1._wp & ! snw enthalpy 307 & , sxxc0 , 'T', 1._wp, syyc0 , 'T', 1._wp, sxyc0 , 'T', 1._wp ) 308 CALL lbc_lnk_multi( 'icedyn_adv_pra', z0ei , 'T', 1._wp, sxe , 'T', -1._wp, sye , 'T', -1._wp & ! ice enthalpy 309 & , sxxe , 'T', 1._wp, syye , 'T', 1._wp, sxye , 'T', 1._wp ) 310 IF ( ln_pnd_LEV ) THEN 311 CALL lbc_lnk_multi( 'icedyn_adv_pra', z0ap , 'T', 1._wp, sxap , 'T', -1._wp, syap , 'T', -1._wp & ! melt pond fraction 312 & , sxxap, 'T', 1._wp, syyap, 'T', 1._wp, sxyap, 'T', 1._wp & 313 & , z0vp , 'T', 1._wp, sxvp , 'T', -1._wp, syvp , 'T', -1._wp & ! melt pond volume 314 & , sxxvp, 'T', 1._wp, syyvp, 'T', 1._wp, sxyvp, 'T', 1._wp ) 315 IF ( ln_pnd_lids ) THEN 316 CALL lbc_lnk_multi( 'icedyn_adv_pra', z0vl ,'T', 1._wp, sxvl ,'T', -1._wp, syvl ,'T', -1._wp & ! melt pond lid volume 317 & , sxxvl,'T', 1._wp, syyvl,'T', 1._wp, sxyvl,'T', 1._wp ) 318 ENDIF 292 319 ENDIF 293 320 … … 368 395 !! 369 396 INTEGER :: ji, jj, jl, jcat ! dummy loop indices 397 INTEGER :: jj0 ! dummy loop indices 370 398 REAL(wp) :: zs1max, zslpmax, ztemp ! local scalars 371 399 REAL(wp) :: zs1new, zalf , zalfq , zbt ! - - … … 375 403 REAL(wp), DIMENSION(jpi,jpj) :: zalg, zalg1, zalg1q ! - - 376 404 !----------------------------------------------------------------------- 405 ! in order to avoid lbc_lnk (communications): 406 ! jj loop must be 1:jpj if adv_x is called first 407 ! and 2:jpj-1 if adv_x is called second 408 jj0 = NINT(pcrh) 377 409 ! 378 410 jcat = SIZE( ps0 , 3 ) ! size of input arrays … … 381 413 ! 382 414 ! Limitation of moments. 383 DO_2D( 0, 0, 1, 1 ) 384 ! Initialize volumes of boxes (=area if adv_x first called, =psm otherwise) 385 psm (ji,jj,jl) = MAX( pcrh * e1e2t(ji,jj) + ( 1.0 - pcrh ) * psm(ji,jj,jl) , epsi20 ) 386 ! 387 zslpmax = MAX( 0._wp, ps0(ji,jj,jl) ) 388 zs1max = 1.5 * zslpmax 389 zs1new = MIN( zs1max, MAX( -zs1max, psx(ji,jj,jl) ) ) 390 zs2new = MIN( 2.0 * zslpmax - 0.3334 * ABS( zs1new ), & 391 & MAX( ABS( zs1new ) - zslpmax, psxx(ji,jj,jl) ) ) 392 rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1) ! Case of empty boxes & Apply mask 393 394 ps0 (ji,jj,jl) = zslpmax 395 psx (ji,jj,jl) = zs1new * rswitch 396 psxx(ji,jj,jl) = zs2new * rswitch 397 psy (ji,jj,jl) = psy (ji,jj,jl) * rswitch 398 psyy(ji,jj,jl) = psyy(ji,jj,jl) * rswitch 399 psxy(ji,jj,jl) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj,jl) ) ) * rswitch 400 END_2D 401 402 ! Calculate fluxes and moments between boxes i<-->i+1 403 DO_2D( 0, 0, 1, 1 ) ! Flux from i to i+1 WHEN u GT 0 404 zbet(ji,jj) = MAX( 0._wp, SIGN( 1._wp, put(ji,jj) ) ) 405 zalf = MAX( 0._wp, put(ji,jj) ) * pdt / psm(ji,jj,jl) 406 zalfq = zalf * zalf 407 zalf1 = 1.0 - zalf 408 zalf1q = zalf1 * zalf1 409 ! 410 zfm (ji,jj) = zalf * psm (ji,jj,jl) 411 zf0 (ji,jj) = zalf * ( ps0 (ji,jj,jl) + zalf1 * ( psx(ji,jj,jl) + (zalf1 - zalf) * psxx(ji,jj,jl) ) ) 412 zfx (ji,jj) = zalfq * ( psx (ji,jj,jl) + 3.0 * zalf1 * psxx(ji,jj,jl) ) 413 zfxx(ji,jj) = zalf * psxx(ji,jj,jl) * zalfq 414 zfy (ji,jj) = zalf * ( psy (ji,jj,jl) + zalf1 * psxy(ji,jj,jl) ) 415 zfxy(ji,jj) = zalfq * psxy(ji,jj,jl) 416 zfyy(ji,jj) = zalf * psyy(ji,jj,jl) 417 418 ! Readjust moments remaining in the box. 419 psm (ji,jj,jl) = psm (ji,jj,jl) - zfm(ji,jj) 420 ps0 (ji,jj,jl) = ps0 (ji,jj,jl) - zf0(ji,jj) 421 psx (ji,jj,jl) = zalf1q * ( psx(ji,jj,jl) - 3.0 * zalf * psxx(ji,jj,jl) ) 422 psxx(ji,jj,jl) = zalf1 * zalf1q * psxx(ji,jj,jl) 423 psy (ji,jj,jl) = psy (ji,jj,jl) - zfy(ji,jj) 424 psyy(ji,jj,jl) = psyy(ji,jj,jl) - zfyy(ji,jj) 425 psxy(ji,jj,jl) = zalf1q * psxy(ji,jj,jl) 426 END_2D 427 428 DO_2D( 0, 0, 1, 0 ) ! Flux from i+1 to i when u LT 0. 429 zalf = MAX( 0._wp, -put(ji,jj) ) * pdt / psm(ji+1,jj,jl) 430 zalg (ji,jj) = zalf 431 zalfq = zalf * zalf 432 zalf1 = 1.0 - zalf 433 zalg1 (ji,jj) = zalf1 434 zalf1q = zalf1 * zalf1 435 zalg1q(ji,jj) = zalf1q 436 ! 437 zfm (ji,jj) = zfm (ji,jj) + zalf * psm (ji+1,jj,jl) 438 zf0 (ji,jj) = zf0 (ji,jj) + zalf * ( ps0 (ji+1,jj,jl) & 439 & - zalf1 * ( psx(ji+1,jj,jl) - (zalf1 - zalf ) * psxx(ji+1,jj,jl) ) ) 440 zfx (ji,jj) = zfx (ji,jj) + zalfq * ( psx (ji+1,jj,jl) - 3.0 * zalf1 * psxx(ji+1,jj,jl) ) 441 zfxx (ji,jj) = zfxx(ji,jj) + zalf * psxx(ji+1,jj,jl) * zalfq 442 zfy (ji,jj) = zfy (ji,jj) + zalf * ( psy (ji+1,jj,jl) - zalf1 * psxy(ji+1,jj,jl) ) 443 zfxy (ji,jj) = zfxy(ji,jj) + zalfq * psxy(ji+1,jj,jl) 444 zfyy (ji,jj) = zfyy(ji,jj) + zalf * psyy(ji+1,jj,jl) 445 END_2D 446 447 DO_2D( 0, 0, 0, 0 ) ! Readjust moments remaining in the box. 448 zbt = zbet(ji-1,jj) 449 zbt1 = 1.0 - zbet(ji-1,jj) 450 ! 451 psm (ji,jj,jl) = zbt * psm(ji,jj,jl) + zbt1 * ( psm(ji,jj,jl) - zfm(ji-1,jj) ) 452 ps0 (ji,jj,jl) = zbt * ps0(ji,jj,jl) + zbt1 * ( ps0(ji,jj,jl) - zf0(ji-1,jj) ) 453 psx (ji,jj,jl) = zalg1q(ji-1,jj) * ( psx(ji,jj,jl) + 3.0 * zalg(ji-1,jj) * psxx(ji,jj,jl) ) 454 psxx(ji,jj,jl) = zalg1 (ji-1,jj) * zalg1q(ji-1,jj) * psxx(ji,jj,jl) 455 psy (ji,jj,jl) = zbt * psy (ji,jj,jl) + zbt1 * ( psy (ji,jj,jl) - zfy (ji-1,jj) ) 456 psyy(ji,jj,jl) = zbt * psyy(ji,jj,jl) + zbt1 * ( psyy(ji,jj,jl) - zfyy(ji-1,jj) ) 457 psxy(ji,jj,jl) = zalg1q(ji-1,jj) * psxy(ji,jj,jl) 458 END_2D 459 460 ! Put the temporary moments into appropriate neighboring boxes. 461 DO_2D( 0, 0, 0, 0 ) ! Flux from i to i+1 IF u GT 0. 462 zbt = zbet(ji-1,jj) 463 zbt1 = 1.0 - zbet(ji-1,jj) 464 psm(ji,jj,jl) = zbt * ( psm(ji,jj,jl) + zfm(ji-1,jj) ) + zbt1 * psm(ji,jj,jl) 465 zalf = zbt * zfm(ji-1,jj) / psm(ji,jj,jl) 466 zalf1 = 1.0 - zalf 467 ztemp = zalf * ps0(ji,jj,jl) - zalf1 * zf0(ji-1,jj) 468 ! 469 ps0 (ji,jj,jl) = zbt * ( ps0(ji,jj,jl) + zf0(ji-1,jj) ) + zbt1 * ps0(ji,jj,jl) 470 psx (ji,jj,jl) = zbt * ( zalf * zfx(ji-1,jj) + zalf1 * psx(ji,jj,jl) + 3.0 * ztemp ) + zbt1 * psx(ji,jj,jl) 471 psxx(ji,jj,jl) = zbt * ( zalf * zalf * zfxx(ji-1,jj) + zalf1 * zalf1 * psxx(ji,jj,jl) & 472 & + 5.0 * ( zalf * zalf1 * ( psx (ji,jj,jl) - zfx(ji-1,jj) ) - ( zalf1 - zalf ) * ztemp ) ) & 473 & + zbt1 * psxx(ji,jj,jl) 474 psxy(ji,jj,jl) = zbt * ( zalf * zfxy(ji-1,jj) + zalf1 * psxy(ji,jj,jl) & 475 & + 3.0 * (- zalf1*zfy(ji-1,jj) + zalf * psy(ji,jj,jl) ) ) & 476 & + zbt1 * psxy(ji,jj,jl) 477 psy (ji,jj,jl) = zbt * ( psy (ji,jj,jl) + zfy (ji-1,jj) ) + zbt1 * psy (ji,jj,jl) 478 psyy(ji,jj,jl) = zbt * ( psyy(ji,jj,jl) + zfyy(ji-1,jj) ) + zbt1 * psyy(ji,jj,jl) 479 END_2D 480 481 DO_2D( 0, 0, 0, 0 ) ! Flux from i+1 to i IF u LT 0. 482 zbt = zbet(ji,jj) 483 zbt1 = 1.0 - zbet(ji,jj) 484 psm(ji,jj,jl) = zbt * psm(ji,jj,jl) + zbt1 * ( psm(ji,jj,jl) + zfm(ji,jj) ) 485 zalf = zbt1 * zfm(ji,jj) / psm(ji,jj,jl) 486 zalf1 = 1.0 - zalf 487 ztemp = - zalf * ps0(ji,jj,jl) + zalf1 * zf0(ji,jj) 488 ! 489 ps0 (ji,jj,jl) = zbt * ps0 (ji,jj,jl) + zbt1 * ( ps0(ji,jj,jl) + zf0(ji,jj) ) 490 psx (ji,jj,jl) = zbt * psx (ji,jj,jl) + zbt1 * ( zalf * zfx(ji,jj) + zalf1 * psx(ji,jj,jl) + 3.0 * ztemp ) 491 psxx(ji,jj,jl) = zbt * psxx(ji,jj,jl) + zbt1 * ( zalf * zalf * zfxx(ji,jj) + zalf1 * zalf1 * psxx(ji,jj,jl) & 492 & + 5.0 * ( zalf * zalf1 * ( - psx(ji,jj,jl) + zfx(ji,jj) ) & 493 & + ( zalf1 - zalf ) * ztemp ) ) 494 psxy(ji,jj,jl) = zbt * psxy(ji,jj,jl) + zbt1 * ( zalf * zfxy(ji,jj) + zalf1 * psxy(ji,jj,jl) & 495 & + 3.0 * ( zalf1 * zfy(ji,jj) - zalf * psy(ji,jj,jl) ) ) 496 psy (ji,jj,jl) = zbt * psy (ji,jj,jl) + zbt1 * ( psy (ji,jj,jl) + zfy (ji,jj) ) 497 psyy(ji,jj,jl) = zbt * psyy(ji,jj,jl) + zbt1 * ( psyy(ji,jj,jl) + zfyy(ji,jj) ) 498 END_2D 499 415 DO jj = Njs0 - jj0, Nje0 + jj0 416 417 DO ji = Nis0 - 1, Nie0 + 1 418 419 ! Initialize volumes of boxes (=area if adv_x first called, =psm otherwise) 420 psm (ji,jj,jl) = MAX( pcrh * e1e2t(ji,jj) + ( 1.0 - pcrh ) * psm(ji,jj,jl) , epsi20 ) 421 ! 422 zslpmax = MAX( 0._wp, ps0(ji,jj,jl) ) 423 zs1max = 1.5 * zslpmax 424 zs1new = MIN( zs1max, MAX( -zs1max, psx(ji,jj,jl) ) ) 425 zs2new = MIN( 2.0 * zslpmax - 0.3334 * ABS( zs1new ), MAX( ABS( zs1new ) - zslpmax, psxx(ji,jj,jl) ) ) 426 rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1) ! Case of empty boxes & Apply mask 427 428 ps0 (ji,jj,jl) = zslpmax 429 psx (ji,jj,jl) = zs1new * rswitch 430 psxx(ji,jj,jl) = zs2new * rswitch 431 psy (ji,jj,jl) = psy (ji,jj,jl) * rswitch 432 psyy(ji,jj,jl) = psyy(ji,jj,jl) * rswitch 433 psxy(ji,jj,jl) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj,jl) ) ) * rswitch 434 435 ! Calculate fluxes and moments between boxes i<-->i+1 436 ! ! Flux from i to i+1 WHEN u GT 0 437 zbet(ji,jj) = MAX( 0._wp, SIGN( 1._wp, put(ji,jj) ) ) 438 zalf = MAX( 0._wp, put(ji,jj) ) * pdt / psm(ji,jj,jl) 439 zalfq = zalf * zalf 440 zalf1 = 1.0 - zalf 441 zalf1q = zalf1 * zalf1 442 ! 443 zfm (ji,jj) = zalf * psm (ji,jj,jl) 444 zf0 (ji,jj) = zalf * ( ps0 (ji,jj,jl) + zalf1 * ( psx(ji,jj,jl) + (zalf1 - zalf) * psxx(ji,jj,jl) ) ) 445 zfx (ji,jj) = zalfq * ( psx (ji,jj,jl) + 3.0 * zalf1 * psxx(ji,jj,jl) ) 446 zfxx(ji,jj) = zalf * psxx(ji,jj,jl) * zalfq 447 zfy (ji,jj) = zalf * ( psy (ji,jj,jl) + zalf1 * psxy(ji,jj,jl) ) 448 zfxy(ji,jj) = zalfq * psxy(ji,jj,jl) 449 zfyy(ji,jj) = zalf * psyy(ji,jj,jl) 450 451 ! ! Readjust moments remaining in the box. 452 psm (ji,jj,jl) = psm (ji,jj,jl) - zfm(ji,jj) 453 ps0 (ji,jj,jl) = ps0 (ji,jj,jl) - zf0(ji,jj) 454 psx (ji,jj,jl) = zalf1q * ( psx(ji,jj,jl) - 3.0 * zalf * psxx(ji,jj,jl) ) 455 psxx(ji,jj,jl) = zalf1 * zalf1q * psxx(ji,jj,jl) 456 psy (ji,jj,jl) = psy (ji,jj,jl) - zfy(ji,jj) 457 psyy(ji,jj,jl) = psyy(ji,jj,jl) - zfyy(ji,jj) 458 psxy(ji,jj,jl) = zalf1q * psxy(ji,jj,jl) 459 END DO 460 461 DO ji = Nis0 - 1, Nie0 462 ! ! Flux from i+1 to i when u LT 0. 463 zalf = MAX( 0._wp, -put(ji,jj) ) * pdt / psm(ji+1,jj,jl) 464 zalg (ji,jj) = zalf 465 zalfq = zalf * zalf 466 zalf1 = 1.0 - zalf 467 zalg1 (ji,jj) = zalf1 468 zalf1q = zalf1 * zalf1 469 zalg1q(ji,jj) = zalf1q 470 ! 471 zfm (ji,jj) = zfm (ji,jj) + zalf * psm (ji+1,jj,jl) 472 zf0 (ji,jj) = zf0 (ji,jj) + zalf * ( ps0 (ji+1,jj,jl) & 473 & - zalf1 * ( psx(ji+1,jj,jl) - (zalf1 - zalf ) * psxx(ji+1,jj,jl) ) ) 474 zfx (ji,jj) = zfx (ji,jj) + zalfq * ( psx (ji+1,jj,jl) - 3.0 * zalf1 * psxx(ji+1,jj,jl) ) 475 zfxx (ji,jj) = zfxx(ji,jj) + zalf * psxx(ji+1,jj,jl) * zalfq 476 zfy (ji,jj) = zfy (ji,jj) + zalf * ( psy (ji+1,jj,jl) - zalf1 * psxy(ji+1,jj,jl) ) 477 zfxy (ji,jj) = zfxy(ji,jj) + zalfq * psxy(ji+1,jj,jl) 478 zfyy (ji,jj) = zfyy(ji,jj) + zalf * psyy(ji+1,jj,jl) 479 END DO 480 481 DO ji = Nis0, Nie0 482 ! ! Readjust moments remaining in the box. 483 zbt = zbet(ji-1,jj) 484 zbt1 = 1.0 - zbet(ji-1,jj) 485 ! 486 psm (ji,jj,jl) = zbt * psm(ji,jj,jl) + zbt1 * ( psm(ji,jj,jl) - zfm(ji-1,jj) ) 487 ps0 (ji,jj,jl) = zbt * ps0(ji,jj,jl) + zbt1 * ( ps0(ji,jj,jl) - zf0(ji-1,jj) ) 488 psx (ji,jj,jl) = zalg1q(ji-1,jj) * ( psx(ji,jj,jl) + 3.0 * zalg(ji-1,jj) * psxx(ji,jj,jl) ) 489 psxx(ji,jj,jl) = zalg1 (ji-1,jj) * zalg1q(ji-1,jj) * psxx(ji,jj,jl) 490 psy (ji,jj,jl) = zbt * psy (ji,jj,jl) + zbt1 * ( psy (ji,jj,jl) - zfy (ji-1,jj) ) 491 psyy(ji,jj,jl) = zbt * psyy(ji,jj,jl) + zbt1 * ( psyy(ji,jj,jl) - zfyy(ji-1,jj) ) 492 psxy(ji,jj,jl) = zalg1q(ji-1,jj) * psxy(ji,jj,jl) 493 494 ! Put the temporary moments into appropriate neighboring boxes. 495 ! ! Flux from i to i+1 IF u GT 0. 496 zbt = zbet(ji-1,jj) 497 zbt1 = 1.0 - zbet(ji-1,jj) 498 psm(ji,jj,jl) = zbt * ( psm(ji,jj,jl) + zfm(ji-1,jj) ) + zbt1 * psm(ji,jj,jl) 499 zalf = zbt * zfm(ji-1,jj) / psm(ji,jj,jl) 500 zalf1 = 1.0 - zalf 501 ztemp = zalf * ps0(ji,jj,jl) - zalf1 * zf0(ji-1,jj) 502 ! 503 ps0 (ji,jj,jl) = zbt * ( ps0(ji,jj,jl) + zf0(ji-1,jj) ) + zbt1 * ps0(ji,jj,jl) 504 psx (ji,jj,jl) = zbt * ( zalf * zfx(ji-1,jj) + zalf1 * psx(ji,jj,jl) + 3.0 * ztemp ) + zbt1 * psx(ji,jj,jl) 505 psxx(ji,jj,jl) = zbt * ( zalf * zalf * zfxx(ji-1,jj) + zalf1 * zalf1 * psxx(ji,jj,jl) & 506 & + 5.0 * ( zalf * zalf1 * ( psx (ji,jj,jl) - zfx(ji-1,jj) ) - ( zalf1 - zalf ) * ztemp ) ) & 507 & + zbt1 * psxx(ji,jj,jl) 508 psxy(ji,jj,jl) = zbt * ( zalf * zfxy(ji-1,jj) + zalf1 * psxy(ji,jj,jl) & 509 & + 3.0 * (- zalf1*zfy(ji-1,jj) + zalf * psy(ji,jj,jl) ) ) & 510 & + zbt1 * psxy(ji,jj,jl) 511 psy (ji,jj,jl) = zbt * ( psy (ji,jj,jl) + zfy (ji-1,jj) ) + zbt1 * psy (ji,jj,jl) 512 psyy(ji,jj,jl) = zbt * ( psyy(ji,jj,jl) + zfyy(ji-1,jj) ) + zbt1 * psyy(ji,jj,jl) 513 514 ! ! Flux from i+1 to i IF u LT 0. 515 zbt = zbet(ji,jj) 516 zbt1 = 1.0 - zbet(ji,jj) 517 psm(ji,jj,jl) = zbt * psm(ji,jj,jl) + zbt1 * ( psm(ji,jj,jl) + zfm(ji,jj) ) 518 zalf = zbt1 * zfm(ji,jj) / psm(ji,jj,jl) 519 zalf1 = 1.0 - zalf 520 ztemp = - zalf * ps0(ji,jj,jl) + zalf1 * zf0(ji,jj) 521 ! 522 ps0 (ji,jj,jl) = zbt * ps0 (ji,jj,jl) + zbt1 * ( ps0(ji,jj,jl) + zf0(ji,jj) ) 523 psx (ji,jj,jl) = zbt * psx (ji,jj,jl) + zbt1 * ( zalf * zfx(ji,jj) + zalf1 * psx(ji,jj,jl) + 3.0 * ztemp ) 524 psxx(ji,jj,jl) = zbt * psxx(ji,jj,jl) + zbt1 * ( zalf * zalf * zfxx(ji,jj) + zalf1 * zalf1 * psxx(ji,jj,jl) & 525 & + 5.0 * ( zalf * zalf1 * ( - psx(ji,jj,jl) + zfx(ji,jj) ) & 526 & + ( zalf1 - zalf ) * ztemp ) ) 527 psxy(ji,jj,jl) = zbt * psxy(ji,jj,jl) + zbt1 * ( zalf * zfxy(ji,jj) + zalf1 * psxy(ji,jj,jl) & 528 & + 3.0 * ( zalf1 * zfy(ji,jj) - zalf * psy(ji,jj,jl) ) ) 529 psy (ji,jj,jl) = zbt * psy (ji,jj,jl) + zbt1 * ( psy (ji,jj,jl) + zfy (ji,jj) ) 530 psyy(ji,jj,jl) = zbt * psyy(ji,jj,jl) + zbt1 * ( psyy(ji,jj,jl) + zfyy(ji,jj) ) 531 END DO 532 ! 533 END DO 534 ! 500 535 END DO 501 502 !-- Lateral boundary conditions 503 CALL lbc_lnk_multi( 'icedyn_adv_pra', psm(:,:,1:jcat) , 'T', 1.0_wp, ps0 , 'T', 1.0_wp & 504 & , psx , 'T', -1.0_wp, psy , 'T', -1.0_wp & ! caution gradient ==> the sign changes 505 & , psxx , 'T', 1.0_wp, psyy, 'T', 1.0_wp , psxy, 'T', 1.0_wp ) 506 ! 536 ! 507 537 END SUBROUTINE adv_x 508 538 … … 525 555 !! 526 556 INTEGER :: ji, jj, jl, jcat ! dummy loop indices 557 INTEGER :: ji0 ! dummy loop indices 527 558 REAL(wp) :: zs1max, zslpmax, ztemp ! temporary scalars 528 559 REAL(wp) :: zs1new, zalf , zalfq , zbt ! - - … … 532 563 REAL(wp), DIMENSION(jpi,jpj) :: zalg, zalg1, zalg1q ! - - 533 564 !--------------------------------------------------------------------- 565 ! in order to avoid lbc_lnk (communications): 566 ! ji loop must be 1:jpi if adv_y is called first 567 ! and 2:jpi-1 if adv_y is called second 568 ji0 = NINT(pcrh) 534 569 ! 535 570 jcat = SIZE( ps0 , 3 ) ! size of input arrays … … 538 573 ! 539 574 ! Limitation of moments. 540 DO_2D( 1, 1, 0,0 )575 DO_2D( 1, 1, ji0, ji0 ) 541 576 ! Initialize volumes of boxes (=area if adv_x first called, =psm otherwise) 542 577 psm(ji,jj,jl) = MAX( pcrh * e1e2t(ji,jj) + ( 1.0 - pcrh ) * psm(ji,jj,jl) , epsi20 ) … … 545 580 zs1max = 1.5 * zslpmax 546 581 zs1new = MIN( zs1max, MAX( -zs1max, psy(ji,jj,jl) ) ) 547 zs2new = MIN( ( 2.0 * zslpmax - 0.3334 * ABS( zs1new ) ), & 548 & MAX( ABS( zs1new )-zslpmax, psyy(ji,jj,jl) ) ) 582 zs2new = MIN( ( 2.0 * zslpmax - 0.3334 * ABS( zs1new ) ), MAX( ABS( zs1new )-zslpmax, psyy(ji,jj,jl) ) ) 549 583 rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1) ! Case of empty boxes & Apply mask 550 584 ! … … 555 589 psyy(ji,jj,jl) = zs2new * rswitch 556 590 psxy(ji,jj,jl) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj,jl) ) ) * rswitch 557 END_2D558 591 559 ! Calculate fluxes and moments between boxes j<-->j+1560 DO_2D( 1, 1, 0, 0 ) ! Flux from j to j+1 WHEN v GT 0592 ! Calculate fluxes and moments between boxes j<-->j+1 593 ! ! Flux from j to j+1 WHEN v GT 0 561 594 zbet(ji,jj) = MAX( 0._wp, SIGN( 1._wp, pvt(ji,jj) ) ) 562 595 zalf = MAX( 0._wp, pvt(ji,jj) ) * pdt / psm(ji,jj,jl) … … 573 606 zfxx(ji,jj) = zalf * psxx(ji,jj,jl) 574 607 ! 575 ! Readjust moments remaining in the box.608 ! ! Readjust moments remaining in the box. 576 609 psm (ji,jj,jl) = psm (ji,jj,jl) - zfm(ji,jj) 577 610 ps0 (ji,jj,jl) = ps0 (ji,jj,jl) - zf0(ji,jj) … … 583 616 END_2D 584 617 ! 585 DO_2D( 1, 0, 0, 0 ) ! Flux from j+1 to j when v LT 0. 618 DO_2D( 1, 0, ji0, ji0 ) 619 ! ! Flux from j+1 to j when v LT 0. 586 620 zalf = MAX( 0._wp, -pvt(ji,jj) ) * pdt / psm(ji,jj+1,jl) 587 621 zalg (ji,jj) = zalf … … 602 636 END_2D 603 637 604 ! Readjust moments remaining in the box.605 DO_2D( 0, 0, 0, 0 )638 DO_2D( 0, 0, ji0, ji0 ) 639 ! ! Readjust moments remaining in the box. 606 640 zbt = zbet(ji,jj-1) 607 641 zbt1 = ( 1.0 - zbet(ji,jj-1) ) … … 614 648 psxx(ji,jj,jl) = zbt * psxx(ji,jj,jl) + zbt1 * ( psxx(ji,jj,jl) - zfxx(ji,jj-1) ) 615 649 psxy(ji,jj,jl) = zalg1q(ji,jj-1) * psxy(ji,jj,jl) 616 END_2D 617 618 ! Put the temporary moments into appropriate neighboring boxes. 619 DO_2D( 0, 0, 0, 0 ) ! Flux from j to j+1 IF v GT 0. 650 651 ! Put the temporary moments into appropriate neighboring boxes. 652 ! ! Flux from j to j+1 IF v GT 0. 620 653 zbt = zbet(ji,jj-1) 621 654 zbt1 = 1.0 - zbet(ji,jj-1) … … 636 669 psx (ji,jj,jl) = zbt * ( psx (ji,jj,jl) + zfx (ji,jj-1) ) + zbt1 * psx (ji,jj,jl) 637 670 psxx(ji,jj,jl) = zbt * ( psxx(ji,jj,jl) + zfxx(ji,jj-1) ) + zbt1 * psxx(ji,jj,jl) 638 END_2D 639 640 DO_2D( 0, 0, 0, 0 ) ! Flux from j+1 to j IF v LT 0. 671 672 ! ! Flux from j+1 to j IF v LT 0. 641 673 zbt = zbet(ji,jj) 642 674 zbt1 = 1.0 - zbet(ji,jj) … … 656 688 psxx(ji,jj,jl) = zbt * psxx(ji,jj,jl) + zbt1 * ( psxx(ji,jj,jl) + zfxx(ji,jj) ) 657 689 END_2D 658 690 ! 659 691 END DO 660 661 !-- Lateral boundary conditions662 CALL lbc_lnk_multi( 'icedyn_adv_pra', psm(:,:,1:jcat) , 'T', 1.0_wp, ps0 , 'T', 1.0_wp &663 & , psx , 'T', -1.0_wp, psy , 'T', -1.0_wp & ! caution gradient ==> the sign changes664 & , psxx , 'T', 1.0_wp, psyy, 'T', 1.0_wp , psxy, 'T', 1.0_wp )665 692 ! 666 693 END SUBROUTINE adv_y … … 890 917 ! 891 918 ! ! ice thickness 892 CALL iom_get( numrir, jpdom_auto, 'sxice' , sxice )893 CALL iom_get( numrir, jpdom_auto, 'syice' , syice )919 CALL iom_get( numrir, jpdom_auto, 'sxice' , sxice , psgn = -1._wp ) 920 CALL iom_get( numrir, jpdom_auto, 'syice' , syice , psgn = -1._wp ) 894 921 CALL iom_get( numrir, jpdom_auto, 'sxxice', sxxice ) 895 922 CALL iom_get( numrir, jpdom_auto, 'syyice', syyice ) 896 923 CALL iom_get( numrir, jpdom_auto, 'sxyice', sxyice ) 897 924 ! ! snow thickness 898 CALL iom_get( numrir, jpdom_auto, 'sxsn' , sxsn )899 CALL iom_get( numrir, jpdom_auto, 'sysn' , sysn )925 CALL iom_get( numrir, jpdom_auto, 'sxsn' , sxsn , psgn = -1._wp ) 926 CALL iom_get( numrir, jpdom_auto, 'sysn' , sysn , psgn = -1._wp ) 900 927 CALL iom_get( numrir, jpdom_auto, 'sxxsn' , sxxsn ) 901 928 CALL iom_get( numrir, jpdom_auto, 'syysn' , syysn ) 902 929 CALL iom_get( numrir, jpdom_auto, 'sxysn' , sxysn ) 903 930 ! ! ice concentration 904 CALL iom_get( numrir, jpdom_auto, 'sxa' , sxa )905 CALL iom_get( numrir, jpdom_auto, 'sya' , sya )931 CALL iom_get( numrir, jpdom_auto, 'sxa' , sxa , psgn = -1._wp ) 932 CALL iom_get( numrir, jpdom_auto, 'sya' , sya , psgn = -1._wp ) 906 933 CALL iom_get( numrir, jpdom_auto, 'sxxa' , sxxa ) 907 934 CALL iom_get( numrir, jpdom_auto, 'syya' , syya ) 908 935 CALL iom_get( numrir, jpdom_auto, 'sxya' , sxya ) 909 936 ! ! ice salinity 910 CALL iom_get( numrir, jpdom_auto, 'sxsal' , sxsal )911 CALL iom_get( numrir, jpdom_auto, 'sysal' , sysal )937 CALL iom_get( numrir, jpdom_auto, 'sxsal' , sxsal , psgn = -1._wp ) 938 CALL iom_get( numrir, jpdom_auto, 'sysal' , sysal , psgn = -1._wp ) 912 939 CALL iom_get( numrir, jpdom_auto, 'sxxsal', sxxsal ) 913 940 CALL iom_get( numrir, jpdom_auto, 'syysal', syysal ) 914 941 CALL iom_get( numrir, jpdom_auto, 'sxysal', sxysal ) 915 942 ! ! ice age 916 CALL iom_get( numrir, jpdom_auto, 'sxage' , sxage )917 CALL iom_get( numrir, jpdom_auto, 'syage' , syage )943 CALL iom_get( numrir, jpdom_auto, 'sxage' , sxage , psgn = -1._wp ) 944 CALL iom_get( numrir, jpdom_auto, 'syage' , syage , psgn = -1._wp ) 918 945 CALL iom_get( numrir, jpdom_auto, 'sxxage', sxxage ) 919 946 CALL iom_get( numrir, jpdom_auto, 'syyage', syyage ) … … 922 949 DO jk = 1, nlay_s 923 950 WRITE(zchar1,'(I2.2)') jk 924 znam = 'sxc0'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d ) ; sxc0 (:,:,jk,:) = z3d(:,:,:)925 znam = 'syc0'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d ) ; syc0 (:,:,jk,:) = z3d(:,:,:)951 znam = 'sxc0'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d, psgn = -1._wp ) ; sxc0 (:,:,jk,:) = z3d(:,:,:) 952 znam = 'syc0'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d, psgn = -1._wp ) ; syc0 (:,:,jk,:) = z3d(:,:,:) 926 953 znam = 'sxxc0'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d ) ; sxxc0(:,:,jk,:) = z3d(:,:,:) 927 954 znam = 'syyc0'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d ) ; syyc0(:,:,jk,:) = z3d(:,:,:) … … 931 958 DO jk = 1, nlay_i 932 959 WRITE(zchar1,'(I2.2)') jk 933 znam = 'sxe'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d ) ; sxe (:,:,jk,:) = z3d(:,:,:)934 znam = 'sye'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d ) ; sye (:,:,jk,:) = z3d(:,:,:)960 znam = 'sxe'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d, psgn = -1._wp ) ; sxe (:,:,jk,:) = z3d(:,:,:) 961 znam = 'sye'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d, psgn = -1._wp ) ; sye (:,:,jk,:) = z3d(:,:,:) 935 962 znam = 'sxxe'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d ) ; sxxe(:,:,jk,:) = z3d(:,:,:) 936 963 znam = 'syye'//'_l'//zchar1 ; CALL iom_get( numrir, jpdom_auto, znam , z3d ) ; syye(:,:,jk,:) = z3d(:,:,:) … … 940 967 IF( ln_pnd_LEV ) THEN ! melt pond fraction 941 968 IF( iom_varid( numror, 'sxap', ldstop = .FALSE. ) > 0 ) THEN 942 CALL iom_get( numrir, jpdom_auto, 'sxap' , sxap )943 CALL iom_get( numrir, jpdom_auto, 'syap' , syap )969 CALL iom_get( numrir, jpdom_auto, 'sxap' , sxap , psgn = -1._wp ) 970 CALL iom_get( numrir, jpdom_auto, 'syap' , syap , psgn = -1._wp ) 944 971 CALL iom_get( numrir, jpdom_auto, 'sxxap', sxxap ) 945 972 CALL iom_get( numrir, jpdom_auto, 'syyap', syyap ) 946 973 CALL iom_get( numrir, jpdom_auto, 'sxyap', sxyap ) 947 974 ! ! melt pond volume 948 CALL iom_get( numrir, jpdom_auto, 'sxvp' , sxvp )949 CALL iom_get( numrir, jpdom_auto, 'syvp' , syvp )975 CALL iom_get( numrir, jpdom_auto, 'sxvp' , sxvp , psgn = -1._wp ) 976 CALL iom_get( numrir, jpdom_auto, 'syvp' , syvp , psgn = -1._wp ) 950 977 CALL iom_get( numrir, jpdom_auto, 'sxxvp', sxxvp ) 951 978 CALL iom_get( numrir, jpdom_auto, 'syyvp', syyvp ) … … 958 985 IF ( ln_pnd_lids ) THEN ! melt pond lid volume 959 986 IF( iom_varid( numror, 'sxvl', ldstop = .FALSE. ) > 0 ) THEN 960 CALL iom_get( numrir, jpdom_auto, 'sxvl' , sxvl )961 CALL iom_get( numrir, jpdom_auto, 'syvl' , syvl )987 CALL iom_get( numrir, jpdom_auto, 'sxvl' , sxvl , psgn = -1._wp ) 988 CALL iom_get( numrir, jpdom_auto, 'syvl' , syvl , psgn = -1._wp ) 962 989 CALL iom_get( numrir, jpdom_auto, 'sxxvl', sxxvl ) 963 990 CALL iom_get( numrir, jpdom_auto, 'syyvl', syyvl )
Note: See TracChangeset
for help on using the changeset viewer.