- Timestamp:
- 2020-05-18T23:58:40+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/r4.0-HEAD_r12713_clem_dan_fixcpl/src/ICE/icedyn_adv_pra.F90
r12832 r12949 82 82 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pa_ip ! melt pond fraction 83 83 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pv_ip ! melt pond volume 84 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pv_il ! melt pond lid thickness84 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pv_il ! melt pond lid thickness 85 85 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_s ! snw heat content 86 86 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_i ! ice heat content … … 92 92 REAL(wp), DIMENSION(jpi,jpj) :: zati1, zati2 93 93 REAL(wp), DIMENSION(jpi,jpj) :: zudy, zvdx 94 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zhi_max, zhs_max, zhip_max 94 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zhi_max, zhs_max, zhip_max, zs_i, zsi_max 95 REAL(wp), DIMENSION(jpi,jpj,nlay_i,jpl) :: ze_i, zei_max 96 REAL(wp), DIMENSION(jpi,jpj,nlay_s,jpl) :: ze_s, zes_max 95 97 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zarea 96 98 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z0ice, z0snw, z0ai, z0smi, z0oi … … 102 104 IF( kt == nit000 .AND. lwp ) WRITE(numout,*) '-- ice_dyn_adv_pra: Prather advection scheme' 103 105 ! 104 ! --- Record max of the surrounding 9-pts ice thick. (for call Hbig) --- ! 106 ! --- Record max of the surrounding 9-pts (for call Hbig) --- ! 107 ! thickness and salinity 108 WHERE( pv_i(:,:,:) >= epsi10 ) ; zs_i(:,:,:) = psv_i(:,:,:) / pv_i(:,:,:) 109 ELSEWHERE ; zs_i(:,:,:) = 0._wp 110 END WHERE 105 111 DO jl = 1, jpl 106 112 DO jj = 2, jpjm1 … … 118 124 & ph_s (ji+1,jj+1,jl), ph_s (ji-1,jj-1,jl), & 119 125 & ph_s (ji+1,jj-1,jl), ph_s (ji-1,jj+1,jl) ) 126 zsi_max (ji,jj,jl) = MAX( epsi20, zs_i (ji,jj,jl), zs_i (ji+1,jj ,jl), zs_i (ji ,jj+1,jl), & 127 & zs_i (ji-1,jj ,jl), zs_i (ji ,jj-1,jl), & 128 & zs_i (ji+1,jj+1,jl), zs_i (ji-1,jj-1,jl), & 129 & zs_i (ji+1,jj-1,jl), zs_i (ji-1,jj+1,jl) ) 120 130 END DO 121 131 END DO 122 132 END DO 123 CALL lbc_lnk_multi( 'icedyn_adv_pra', zhi_max, 'T', 1., zhs_max, 'T', 1., zhip_max, 'T', 1. ) 133 CALL lbc_lnk_multi( 'icedyn_adv_pra', zhi_max, 'T', 1., zhs_max, 'T', 1., zhip_max, 'T', 1., zsi_max, 'T', 1. ) 134 ! 135 ! enthalpies 136 DO jk = 1, nlay_i 137 WHERE( pv_i(:,:,:) >= epsi10 ) ; ze_i(:,:,jk,:) = pe_i(:,:,jk,:) / pv_i(:,:,:) 138 ELSEWHERE ; ze_i(:,:,jk,:) = 0._wp 139 END WHERE 140 END DO 141 DO jk = 1, nlay_s 142 WHERE( pv_s(:,:,:) >= epsi10 ) ; ze_s(:,:,jk,:) = pe_s(:,:,jk,:) / pv_s(:,:,:) 143 ELSEWHERE ; ze_s(:,:,jk,:) = 0._wp 144 END WHERE 145 END DO 146 DO jl = 1, jpl 147 DO jk = 1, nlay_i 148 DO jj = 2, jpjm1 149 DO ji = fs_2, fs_jpim1 150 zei_max(ji,jj,jk,jl) = MAX( epsi20, ze_i(ji,jj,jk,jl), ze_i(ji+1,jj ,jk,jl), ze_i(ji ,jj+1,jk,jl), & 151 & ze_i(ji-1,jj ,jk,jl), ze_i(ji ,jj-1,jk,jl), & 152 & ze_i(ji+1,jj+1,jk,jl), ze_i(ji-1,jj-1,jk,jl), & 153 & ze_i(ji+1,jj-1,jk,jl), ze_i(ji-1,jj+1,jk,jl) ) 154 END DO 155 END DO 156 END DO 157 END DO 158 DO jl = 1, jpl 159 DO jk = 1, nlay_s 160 DO jj = 2, jpjm1 161 DO ji = fs_2, fs_jpim1 162 zes_max(ji,jj,jk,jl) = MAX( epsi20, ze_s(ji,jj,jk,jl), ze_s(ji+1,jj ,jk,jl), ze_s(ji ,jj+1,jk,jl), & 163 & ze_s(ji-1,jj ,jk,jl), ze_s(ji ,jj-1,jk,jl), & 164 & ze_s(ji+1,jj+1,jk,jl), ze_s(ji-1,jj-1,jk,jl), & 165 & ze_s(ji+1,jj-1,jk,jl), ze_s(ji-1,jj+1,jk,jl) ) 166 END DO 167 END DO 168 END DO 169 END DO 170 CALL lbc_lnk( 'icedyn_adv_pra', zei_max, 'T', 1. ) 171 CALL lbc_lnk( 'icedyn_adv_pra', zes_max, 'T', 1. ) 172 ! 124 173 ! 125 174 ! --- If ice drift is too fast, use subtime steps for advection (CFL test for stability) --- ! … … 283 332 ! --- Make sure ice thickness is not too big --- ! 284 333 ! (because ice thickness can be too large where ice concentration is very small) 285 CALL Hbig( zdt, zhi_max, zhs_max, zhip_max, pv_i, pv_s, pa_i, pa_ip, pv_ip, pe_s ) 334 CALL Hbig( zdt, zhi_max, zhs_max, zhip_max, zsi_max, zes_max, zei_max, & 335 & pv_i, pv_s, pa_i, pa_ip, pv_ip, psv_i, pe_s, pe_i ) 286 336 ! 287 337 ! --- Ensure snow load is not too big --- ! … … 635 685 636 686 637 SUBROUTINE Hbig( pdt, phi_max, phs_max, phip_max, pv_i, pv_s, pa_i, pa_ip, pv_ip, pe_s ) 687 SUBROUTINE Hbig( pdt, phi_max, phs_max, phip_max, psi_max, pes_max, pei_max, & 688 & pv_i, pv_s, pa_i, pa_ip, pv_ip, psv_i, pe_s, pe_i ) 638 689 !!------------------------------------------------------------------- 639 690 !! *** ROUTINE Hbig *** … … 649 700 !! ** input : Max thickness of the surrounding 9-points 650 701 !!------------------------------------------------------------------- 651 REAL(wp) , INTENT(in ) :: pdt ! tracer time-step 652 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: phi_max, phs_max, phip_max ! max ice thick from surrounding 9-pts 653 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pv_i, pv_s, pa_i, pa_ip, pv_ip 702 REAL(wp) , INTENT(in ) :: pdt ! tracer time-step 703 REAL(wp), DIMENSION(:,:,:) , INTENT(in ) :: phi_max, phs_max, phip_max, psi_max ! max ice thick from surrounding 9-pts 704 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pes_max 705 REAL(wp), DIMENSION(:,:,:,:), INTENT(in ) :: pei_max 706 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pv_i, pv_s, pa_i, pa_ip, pv_ip, psv_i 654 707 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_s 655 ! 656 INTEGER :: ji, jj, jl ! dummy loop indices 657 REAL(wp) :: z1_dt, zhip, zhi, zhs, zfra 708 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_i 709 ! 710 INTEGER :: ji, jj, jk, jl ! dummy loop indices 711 REAL(wp) :: z1_dt, zhip, zhi, zhs, zsi, zes, zei, zfra 658 712 !!------------------------------------------------------------------- 659 713 ! … … 661 715 ! 662 716 DO jl = 1, jpl 663 664 717 DO jj = 1, jpj 665 718 DO ji = 1, jpi … … 695 748 ENDIF 696 749 ! 750 ! ! -- check s_i -- ! 751 ! if s_i is larger than the surrounding 9 pts => put salt excess in the ocean 752 zsi = psv_i(ji,jj,jl) / pv_i(ji,jj,jl) 753 IF( zsi > psi_max(ji,jj,jl) .AND. pa_i(ji,jj,jl) < 0.15 ) THEN 754 zfra = psi_max(ji,jj,jl) / zsi 755 sfx_res(ji,jj) = sfx_res(ji,jj) + psv_i(ji,jj,jl) * ( 1._wp - zfra ) * rhoi * z1_dt 756 psv_i(ji,jj,jl) = psv_i(ji,jj,jl) * zfra 757 ENDIF 758 ! 697 759 ENDIF 698 760 END DO 699 761 END DO 700 762 END DO 763 ! 764 ! ! -- check e_i/v_i -- ! 765 DO jl = 1, jpl 766 DO jk = 1, nlay_i 767 DO jj = 1, jpj 768 DO ji = 1, jpi 769 IF ( pv_i(ji,jj,jl) > 0._wp ) THEN 770 ! if e_i/v_i is larger than the surrounding 9 pts => put the heat excess in the ocean 771 zei = pe_i(ji,jj,jk,jl) / pv_i(ji,jj,jl) 772 IF( zei > pei_max(ji,jj,jk,jl) .AND. pa_i(ji,jj,jl) < 0.15 ) THEN 773 zfra = pei_max(ji,jj,jk,jl) / zei 774 hfx_res(ji,jj) = hfx_res(ji,jj) - pe_i(ji,jj,jk,jl) * ( 1._wp - zfra ) * z1_dt ! W.m-2 <0 775 pe_i(ji,jj,jk,jl) = pe_i(ji,jj,jk,jl) * zfra 776 ENDIF 777 ENDIF 778 END DO 779 END DO 780 END DO 781 END DO 782 ! ! -- check e_s/v_s -- ! 783 DO jl = 1, jpl 784 DO jk = 1, nlay_s 785 DO jj = 1, jpj 786 DO ji = 1, jpi 787 IF ( pv_s(ji,jj,jl) > 0._wp ) THEN 788 ! if e_s/v_s is larger than the surrounding 9 pts => put the heat excess in the ocean 789 zes = pe_s(ji,jj,jk,jl) / pv_s(ji,jj,jl) 790 IF( zes > pes_max(ji,jj,jk,jl) .AND. pa_i(ji,jj,jl) < 0.15 ) THEN 791 zfra = pes_max(ji,jj,jk,jl) / zes 792 hfx_res(ji,jj) = hfx_res(ji,jj) - pe_s(ji,jj,jk,jl) * ( 1._wp - zfra ) * z1_dt ! W.m-2 <0 793 pe_s(ji,jj,jk,jl) = pe_s(ji,jj,jk,jl) * zfra 794 ENDIF 795 ENDIF 796 END DO 797 END DO 798 END DO 799 END DO 701 800 ! 702 801 END SUBROUTINE Hbig
Note: See TracChangeset
for help on using the changeset viewer.