- Timestamp:
- 2017-06-28T17:55:50+02:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90
r7753 r8239 117 117 REAL(wp) :: zm1, zm2, zm3, zmassU, zmassV ! ice/snow mass 118 118 REAL(wp) :: zdelta, zp_delf, zds2, zdt, zdt2, zdiv, zdiv2 ! temporary scalars 119 REAL(wp) :: zTauO, zTauB, zTauE, z Cor, zvel! temporary scalars119 REAL(wp) :: zTauO, zTauB, zTauE, zvel ! temporary scalars 120 120 121 121 REAL(wp) :: zsig1, zsig2 ! internal ice stress 122 122 REAL(wp) :: zresm ! Maximal error on ice velocity 123 123 REAL(wp) :: zintb, zintn ! dummy argument 124 REAL(wp) :: zfac_x, zfac_y 124 125 125 126 REAL(wp), POINTER, DIMENSION(:,:) :: z1_e1t0, z1_e2t0 ! scale factors … … 140 141 ! ocean surface (ssh_m) if ice is not embedded 141 142 ! ice top surface if ice is embedded 143 REAL(wp), POINTER, DIMENSION(:,:) :: zCorx, zCory ! Coriolis stress array 144 REAL(wp), POINTER, DIMENSION(:,:) :: ztaux_oi, ztauy_oi ! Ocean-to-ice stress array 145 142 146 REAL(wp), POINTER, DIMENSION(:,:) :: zswitchU, zswitchV ! dummy arrays 143 147 REAL(wp), POINTER, DIMENSION(:,:) :: zmaskU, zmaskV ! mask for ice presence … … 154 158 CALL wrk_alloc( jpi,jpj, zds, zs1, zs2, zs12, zu_ice, zv_ice, zresr, zpice ) 155 159 CALL wrk_alloc( jpi,jpj, zswitchU, zswitchV, zmaskU, zmaskV, zfmask, zwf ) 160 CALL wrk_alloc( jpi,jpj, zCorx, zCory) 161 CALL wrk_alloc( jpi,jpj, ztaux_oi, ztauy_oi) 156 162 157 163 #if defined key_agrif … … 427 433 & + ( u_iceV(ji,jj) - u_oceV(ji,jj) ) * ( u_iceV(ji,jj) - u_oceV(ji,jj) ) ) 428 434 435 ! Ocean-to-Ice stress 436 ztauy_oi(ji,jj) = zTauO * ( v_oce(ji,jj) - v_ice(ji,jj) ) 437 429 438 ! tau_bottom/v_ice 430 439 zvel = MAX( zepsi, SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) ) … … 432 441 433 442 ! Coriolis at V-points (energy conserving formulation) 434 zCor = - 0.25_wp * r1_e2v(ji,jj) * &443 zCory(ji,jj) = - 0.25_wp * r1_e2v(ji,jj) * & 435 444 & ( zmf(ji,jj ) * ( e2u(ji,jj ) * u_ice(ji,jj ) + e2u(ji-1,jj ) * u_ice(ji-1,jj ) ) & 436 445 & + zmf(ji,jj+1) * ( e2u(ji,jj+1) * u_ice(ji,jj+1) + e2u(ji-1,jj+1) * u_ice(ji-1,jj+1) ) ) 437 446 438 447 ! Sum of external forces (explicit solution) = F + tau_ia + Coriolis + spg + tau_io 439 zTauE = zfV(ji,jj) + zTauV_ia(ji,jj) + zCor + zspgV(ji,jj) + zTauO * ( v_oce(ji,jj) - v_ice(ji,jj))448 zTauE = zfV(ji,jj) + zTauV_ia(ji,jj) + zCory(ji,jj) + zspgV(ji,jj) + ztauy_oi(ji,jj) 440 449 441 450 ! landfast switch => 0 = static friction ; 1 = sliding friction … … 451 460 ! Bouillon 2013 452 461 !!v_ice(ji,jj) = ( zmV_t(ji,jj) * ( zbeta * v_ice(ji,jj) + v_ice_b(ji,jj) ) & 453 !! & + zfV(ji,jj) + zCor + zTauV_ia(ji,jj) + zTauO * v_oce(ji,jj) + zspgV(ji,jj) &462 !! & + zfV(ji,jj) + zCory(ji,jj) + zTauV_ia(ji,jj) + zTauO * v_oce(ji,jj) + zspgV(ji,jj) & 454 463 !! & ) / MAX( zmV_t(ji,jj) * ( zbeta + 1._wp ) + zTauO - zTauB, zepsi ) * zswitchV(ji,jj) 455 464 … … 471 480 & + ( v_iceU(ji,jj) - v_oceU(ji,jj) ) * ( v_iceU(ji,jj) - v_oceU(ji,jj) ) ) 472 481 482 ! Ocean-to-Ice stress 483 ztaux_oi(ji,jj) = zTauO * ( u_oce(ji,jj) - u_ice(ji,jj) ) 484 473 485 ! tau_bottom/u_ice 474 486 zvel = MAX( zepsi, SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) ) … … 476 488 477 489 ! Coriolis at U-points (energy conserving formulation) 478 zCor = 0.25_wp * r1_e1u(ji,jj) * &490 zCorx(ji,jj) = 0.25_wp * r1_e1u(ji,jj) * & 479 491 & ( zmf(ji ,jj) * ( e1v(ji ,jj) * v_ice(ji ,jj) + e1v(ji ,jj-1) * v_ice(ji ,jj-1) ) & 480 492 & + zmf(ji+1,jj) * ( e1v(ji+1,jj) * v_ice(ji+1,jj) + e1v(ji+1,jj-1) * v_ice(ji+1,jj-1) ) ) 481 493 482 494 ! Sum of external forces (explicit solution) = F + tau_ia + Coriolis + spg + tau_io 483 zTauE = zfU(ji,jj) + zTauU_ia(ji,jj) + zCor + zspgU(ji,jj) + zTauO * ( u_oce(ji,jj) - u_ice(ji,jj))495 zTauE = zfU(ji,jj) + zTauU_ia(ji,jj) + zCorx(ji,jj) + zspgU(ji,jj) + ztaux_oi(ji,jj) 484 496 485 497 ! landfast switch => 0 = static friction ; 1 = sliding friction … … 495 507 ! Bouillon 2013 496 508 !!u_ice(ji,jj) = ( zmU_t(ji,jj) * ( zbeta * u_ice(ji,jj) + u_ice_b(ji,jj) ) & 497 !! & + zfU(ji,jj) + zCor + zTauU_ia(ji,jj) + zTauO * u_oce(ji,jj) + zspgU(ji,jj) &509 !! & + zfU(ji,jj) + zCorx(ji,jj) + zTauU_ia(ji,jj) + zTauO * u_oce(ji,jj) + zspgU(ji,jj) & 498 510 !! & ) / MAX( zmU_t(ji,jj) * ( zbeta + 1._wp ) + zTauO - zTauB, zepsi ) * zswitchU(ji,jj) 499 511 END DO … … 516 528 & + ( v_iceU(ji,jj) - v_oceU(ji,jj) ) * ( v_iceU(ji,jj) - v_oceU(ji,jj) ) ) 517 529 530 ! Ocean-to-Ice stress 531 ztaux_oi(ji,jj) = zTauO * ( u_oce(ji,jj) - u_ice(ji,jj) ) 532 518 533 ! tau_bottom/u_ice 519 534 zvel = MAX( zepsi, SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) ) … … 521 536 522 537 ! Coriolis at U-points (energy conserving formulation) 523 zCor = 0.25_wp * r1_e1u(ji,jj) * &538 zCorx(ji,jj) = 0.25_wp * r1_e1u(ji,jj) * & 524 539 & ( zmf(ji ,jj) * ( e1v(ji ,jj) * v_ice(ji ,jj) + e1v(ji ,jj-1) * v_ice(ji ,jj-1) ) & 525 540 & + zmf(ji+1,jj) * ( e1v(ji+1,jj) * v_ice(ji+1,jj) + e1v(ji+1,jj-1) * v_ice(ji+1,jj-1) ) ) 526 541 527 542 ! Sum of external forces (explicit solution) = F + tau_ia + Coriolis + spg + tau_io 528 zTauE = zfU(ji,jj) + zTauU_ia(ji,jj) + zCor + zspgU(ji,jj) + zTauO * ( u_oce(ji,jj) - u_ice(ji,jj))543 zTauE = zfU(ji,jj) + zTauU_ia(ji,jj) + zCorx(ji,jj) + zspgU(ji,jj) + ztaux_oi(ji,jj) 529 544 530 545 ! landfast switch => 0 = static friction ; 1 = sliding friction … … 540 555 ! Bouillon 2013 541 556 !!u_ice(ji,jj) = ( zmU_t(ji,jj) * ( zbeta * u_ice(ji,jj) + u_ice_b(ji,jj) ) & 542 !! & + zfU(ji,jj) + zCor + zTauU_ia(ji,jj) + zTauO * u_oce(ji,jj) + zspgU(ji,jj) &557 !! & + zfU(ji,jj) + zCorx(ji,jj) + zTauU_ia(ji,jj) + zTauO * u_oce(ji,jj) + zspgU(ji,jj) & 543 558 !! & ) / MAX( zmU_t(ji,jj) * ( zbeta + 1._wp ) + zTauO - zTauB, zepsi ) * zswitchU(ji,jj) 544 559 END DO … … 559 574 & + ( u_iceV(ji,jj) - u_oceV(ji,jj) ) * ( u_iceV(ji,jj) - u_oceV(ji,jj) ) ) 560 575 576 ! Ocean-to-Ice stress 577 ztauy_oi(ji,jj) = zTauO * ( v_oce(ji,jj) - v_ice(ji,jj) ) 578 561 579 ! tau_bottom/v_ice 562 580 zvel = MAX( zepsi, SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) ) … … 564 582 565 583 ! Coriolis at V-points (energy conserving formulation) 566 zCor = - 0.25_wp * r1_e2v(ji,jj) * &584 zCory(ji,jj) = - 0.25_wp * r1_e2v(ji,jj) * & 567 585 & ( zmf(ji,jj ) * ( e2u(ji,jj ) * u_ice(ji,jj ) + e2u(ji-1,jj ) * u_ice(ji-1,jj ) ) & 568 586 & + zmf(ji,jj+1) * ( e2u(ji,jj+1) * u_ice(ji,jj+1) + e2u(ji-1,jj+1) * u_ice(ji-1,jj+1) ) ) 569 587 570 588 ! Sum of external forces (explicit solution) = F + tau_ia + Coriolis + spg + tau_io 571 zTauE = zfV(ji,jj) + zTauV_ia(ji,jj) + zCor + zspgV(ji,jj) + zTauO * ( v_oce(ji,jj) - v_ice(ji,jj))589 zTauE = zfV(ji,jj) + zTauV_ia(ji,jj) + zCory(ji,jj) + zspgV(ji,jj) + ztauy_oi(ji,jj) 572 590 573 591 ! landfast switch => 0 = static friction (tau_icebfr > zTauE); 1 = sliding friction … … 583 601 ! Bouillon 2013 584 602 !!v_ice(ji,jj) = ( zmV_t(ji,jj) * ( zbeta * v_ice(ji,jj) + v_ice_b(ji,jj) ) & 585 !! & + zfV(ji,jj) + zCor + zTauV_ia(ji,jj) + zTauO * v_oce(ji,jj) + zspgV(ji,jj) &603 !! & + zfV(ji,jj) + zCory(ji,jj) + zTauV_ia(ji,jj) + zTauO * v_oce(ji,jj) + zspgV(ji,jj) & 586 604 !! & ) / MAX( zmV_t(ji,jj) * ( zbeta + 1._wp ) + zTauO - zTauB, zepsi ) * zswitchV(ji,jj) 587 605 END DO … … 662 680 663 681 !------------------------------------------------------------------------------! 664 ! 5) Control prints of residual and charge ellipse 682 ! 5) SIMIP diagnostics 683 !------------------------------------------------------------------------------! 684 685 DO jj = 2, jpjm1 686 DO ji = 2, jpim1 687 rswitch = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) ! 1 if ice, 0 if no ice 688 689 ! Stress tensor invariants (normal and shear stress N/m) 690 diag_sig1(ji,jj) = ( zs1(ji,jj) + zs2(ji,jj) ) * rswitch ! normal stress 691 diag_sig2(ji,jj) = SQRT( ( zs1(ji,jj) - zs2(ji,jj) )**2 + 4*zs12(ji,jj)**2 ) * rswitch ! shear stress 692 693 ! Stress terms of the momentum equation (N/m2) 694 diag_dssh_dx(ji,jj) = zspgU(ji,jj) * rswitch ! sea surface slope stress term 695 diag_dssh_dy(ji,jj) = zspgV(ji,jj) * rswitch 696 697 diag_corstrx(ji,jj) = zCorx(ji,jj) * rswitch ! Coriolis stress term 698 diag_corstry(ji,jj) = zCory(ji,jj) * rswitch 699 700 diag_intstrx(ji,jj) = zfU(ji,jj) * rswitch ! internal stress term 701 diag_intstry(ji,jj) = zfV(ji,jj) * rswitch 702 703 diag_utau_oi(ji,jj) = ztaux_oi(ji,jj) * rswitch ! oceanic stress 704 diag_vtau_oi(ji,jj) = ztauy_oi(ji,jj) * rswitch 705 706 ! 2D ice mass, snow mass, area transport arrays (X, Y) 707 zfac_x = 0.5 * u_ice(ji,jj) * e2u(ji,jj) * rswitch 708 zfac_y = 0.5 * v_ice(ji,jj) * e1v(ji,jj) * rswitch 709 710 diag_xmtrp_ice(ji,jj) = rhoic * zfac_x * ( vt_i(ji+1,jj) + vt_i(ji,jj) ) ! ice mass transport, X-component 711 diag_ymtrp_ice(ji,jj) = rhoic * zfac_y * ( vt_i(ji,jj+1) + vt_i(ji,jj) ) ! '' Y- '' 712 713 diag_xmtrp_snw(ji,jj) = rhosn * zfac_x * ( vt_s(ji+1,jj) + vt_s(ji,jj) ) ! snow mass transport, X-component 714 diag_ymtrp_snw(ji,jj) = rhosn * zfac_y * ( vt_s(ji,jj+1) + vt_s(ji,jj) ) ! '' Y- '' 715 716 diag_xatrp(ji,jj) = zfac_x * ( at_i(ji+1,jj) + at_i(ji,jj) ) ! area transport, X-component 717 diag_yatrp(ji,jj) = zfac_y * ( at_i(ji,jj+1) + at_i(ji,jj) ) ! '' Y- '' 718 719 END DO 720 END DO 721 722 CALL lbc_lnk_multi( diag_sig1 , 'T', 1., diag_sig2 , 'T', 1., & 723 & diag_dssh_dx, 'U', -1., diag_dssh_dy, 'V', -1., & 724 & diag_corstrx, 'U', -1., diag_corstry, 'V', -1., & 725 & diag_intstrx, 'U', -1., diag_intstry, 'V', -1. ) 726 727 CALL lbc_lnk_multi( diag_utau_oi, 'U', -1., diag_vtau_oi, 'V', -1. ) 728 729 ! 730 !------------------------------------------------------------------------------! 731 ! 6) Control prints of residual and charge ellipse 665 732 !------------------------------------------------------------------------------! 666 733 ! … … 703 770 CALL wrk_dealloc( jpi,jpj, zds, zs1, zs2, zs12, zu_ice, zv_ice, zresr, zpice ) 704 771 CALL wrk_dealloc( jpi,jpj, zswitchU, zswitchV, zmaskU, zmaskV, zfmask, zwf ) 772 CALL wrk_dealloc( jpi,jpj, zCorx, zCory ) 773 CALL wrk_dealloc( jpi,jpj, ztaux_oi, ztauy_oi ) 705 774 706 775 END SUBROUTINE lim_rhg
Note: See TracChangeset
for help on using the changeset viewer.