- Timestamp:
- 2020-05-15T18:26:45+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_v2/src/ICE/icethd_zdf_bl99.F90
r11715 r12938 101 101 REAL(wp) :: zfac ! dummy factor 102 102 ! 103 REAL(wp), DIMENSION(jpij) :: isnow ! switch for presence (1) or absence (0) of snow103 REAL(wp), DIMENSION(jpij) :: isnow ! fraction of sea ice that is snow covered (for thermodynamic use only) 104 104 REAL(wp), DIMENSION(jpij) :: ztsub ! surface temperature at previous iteration 105 105 REAL(wp), DIMENSION(jpij) :: zh_i, z1_h_i ! ice layer thickness … … 147 147 !------------------ 148 148 DO ji = 1, npti 149 isnow(ji) = 1._wp - MAX( 0._wp , SIGN(1._wp, - h_s_1d(ji) ) ) ! is there snow or not 149 150 ! If the snow thickness drops below zhs_min then reduce the snow fraction instead 151 IF( h_s_1d(ji) < zhs_min ) THEN 152 isnow(ji) = h_s_1d(ji) / zhs_min 153 zh_s(ji) = zhs_min * r1_nlay_s 154 ELSE 155 isnow(ji) = 1.0_wp 156 zh_s(ji) = h_s_1d(ji) * r1_nlay_s 157 END IF 158 150 159 ! layer thickness 151 160 zh_i(ji) = h_i_1d(ji) * r1_nlay_i 152 zh_s(ji) = h_s_1d(ji) * r1_nlay_s 161 153 162 END DO 154 163 ! … … 156 165 ELSEWHERE ; z1_h_i(1:npti) = 0._wp 157 166 END WHERE 158 !159 WHERE( zh_s(1:npti) > 0._wp ) zh_s(1:npti) = MAX( zhs_min * r1_nlay_s, zh_s(1:npti) )160 167 ! 161 168 WHERE( zh_s(1:npti) > 0._wp ) ; z1_h_s(1:npti) = 1._wp / zh_s(1:npti) … … 769 776 ! 770 777 ! --- calculate conduction fluxes (positive downward) 771 778 ! bottom ice conduction flux 772 779 DO ji = 1, npti 773 ! ! surface ice conduction flux 774 qcn_ice_top_1d(ji) = - isnow(ji) * zkappa_s(ji,0) * zg1s * ( t_s_1d(ji,1) - t_su_1d(ji) ) & 775 & - ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 * ( t_i_1d(ji,1) - t_su_1d(ji) ) 776 ! ! bottom ice conduction flux 777 qcn_ice_bot_1d(ji) = - zkappa_i(ji,nlay_i) * zg1 * ( t_bo_1d(ji ) - t_i_1d (ji,nlay_i) ) 780 qcn_ice_bot_1d(ji) = - zkappa_i(ji,nlay_i) * zg1 * ( t_bo_1d(ji ) - t_i_1d (ji,nlay_i) ) 778 781 END DO 779 782 ! surface ice conduction flux 783 IF( k_cnd == np_cnd_OFF .OR. k_cnd == np_cnd_EMU ) THEN 784 ! 785 DO ji = 1, npti 786 qcn_ice_top_1d(ji) = - isnow(ji) * zkappa_s(ji,0) * zg1s * ( t_s_1d(ji,1) - t_su_1d(ji) ) & 787 & - ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 * ( t_i_1d(ji,1) - t_su_1d(ji) ) 788 END DO 789 ! 790 ELSEIF( k_cnd == np_cnd_ON ) THEN 791 ! 792 DO ji = 1, npti 793 qcn_ice_top_1d(ji) = qcn_ice_1d(ji) 794 END DO 795 ! 796 ENDIF 797 ! surface ice temperature 798 IF( k_cnd == np_cnd_ON .AND. ln_cndemulate ) THEN 799 ! 800 DO ji = 1, npti 801 t_su_1d(ji) = ( qcn_ice_top_1d(ji) & ! calculate surface temperature 802 & + isnow(ji) * zkappa_s(ji,0) * zg1s * t_s_1d(ji,1) & 803 & + ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 * t_i_1d(ji,1) & 804 & ) / MAX( epsi10, isnow(ji) * zkappa_s(ji,0) * zg1s + ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 ) 805 t_su_1d(ji) = MAX( MIN( t_su_1d(ji), rt0 ), rt0 - 100._wp ) ! cap t_su 806 END DO 807 ! 808 ENDIF 780 809 ! 781 810 ! --- Diagnose the heat loss due to changing non-solar / conduction flux --- ! … … 787 816 END DO 788 817 ! 789 ELSEIF( k_cnd == np_cnd_ON ) THEN790 !791 DO ji = 1, npti792 hfx_err_dif_1d(ji) = hfx_err_dif_1d(ji) - ( qcn_ice_top_1d(ji) - qcn_ice_1d(ji) ) * a_i_1d(ji)793 END DO794 !795 818 ENDIF 796 797 819 ! 798 820 ! --- Diagnose the heat loss due to non-fully converged temperature solution (should not be above 10-4 W-m2) --- ! … … 839 861 ! effective conductivity and 1st layer temperature (needed by Met Office) 840 862 DO ji = 1, npti 841 IF( h_ s_1d(ji) > 0.1_wp ) THEN842 cnd_ice_1d(ji) = 2._wp * zkappa_ s(ji,0)863 IF( h_i_1d(ji) > 0.1_wp ) THEN 864 cnd_ice_1d(ji) = 2._wp * zkappa_i(ji,0) 843 865 ELSE 844 IF( h_i_1d(ji) > 0.1_wp ) THEN 845 cnd_ice_1d(ji) = 2._wp * zkappa_i(ji,0) 846 ELSE 847 cnd_ice_1d(ji) = 2._wp * ztcond_i(ji,0) * 10._wp 848 ENDIF 866 cnd_ice_1d(ji) = 2._wp * ztcond_i(ji,0) * 10._wp 849 867 ENDIF 850 868 t1_ice_1d(ji) = isnow(ji) * t_s_1d(ji,1) + ( 1._wp - isnow(ji) ) * t_i_1d(ji,1) … … 856 874 t_i_1d (1:npti,:) = ztiold (1:npti,:) 857 875 qcn_ice_1d(1:npti) = qcn_ice_top_1d(1:npti) 858 859 !!clem860 ! remettre t_su_1d, qns_ice_1d et dqns_ice_1d comme avant puisqu'on devrait faire comme si on avant conduction = input861 !clem862 876 ENDIF 863 877 !
Note: See TracChangeset
for help on using the changeset viewer.