New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 5581 for branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90 – NEMO

Ignore:
Timestamp:
2015-07-10T13:28:53+02:00 (9 years ago)
Author:
timgraham
Message:

Merged head of trunk into branch

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90

    r4765 r5581  
    2525   USE par_oce          ! ocean parameters 
    2626   USE phycst           ! physical constants 
    27    USE par_ice          ! ice parameters 
    2827   USE dom_oce          ! ocean domain 
    29    USE dom_ice,    ONLY : tms, area 
    3028   USE ice              ! LIM sea-ice variables 
    3129   USE sbc_ice          ! Surface boundary condition: sea-ice fields 
    3230   USE sbc_oce          ! Surface boundary condition: ocean fields 
    3331   USE sbccpl 
    34    USE cpl_oasis3, ONLY : lk_cpl 
    35    USE oce       , ONLY : iatte, oatte, sshn, sshb, snwice_mass, snwice_mass_b, snwice_fmass 
     32   USE oce       , ONLY : sshn, sshb, snwice_mass, snwice_mass_b, snwice_fmass 
    3633   USE albedo           ! albedo parameters 
    3734   USE lbclnk           ! ocean lateral boundary condition - MPP exchanges 
     
    4138   USE prtctl           ! Print control 
    4239   USE lib_fortran      ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
    43    USE traqsr           ! clem: add penetration of solar flux into the calculation of heat budget 
     40   USE traqsr           ! add penetration of solar flux in the calculation of heat budget 
    4441   USE iom 
    4542   USE domvvl           ! Variable volume 
     43   USE limctl 
     44   USE limcons 
    4645 
    4746   IMPLICIT NONE 
    4847   PRIVATE 
    4948 
    50    PUBLIC   lim_sbc_init   ! called by ice_init 
     49   PUBLIC   lim_sbc_init   ! called by sbc_lim_init 
    5150   PUBLIC   lim_sbc_flx    ! called by sbc_ice_lim 
    5251   PUBLIC   lim_sbc_tau    ! called by sbc_ice_lim 
    53  
    54    REAL(wp)  ::   epsi10 = 1.e-10   ! 
    55    REAL(wp)  ::   epsi20 = 1.e-20   ! 
    5652 
    5753   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   utau_oce, vtau_oce   ! air-ocean surface i- & j-stress     [N/m2] 
     
    9894      !!              - fr_i    : ice fraction 
    9995      !!              - tn_ice  : sea-ice surface temperature 
    100       !!              - alb_ice : sea-ice alberdo (lk_cpl=T) 
     96      !!              - alb_ice : sea-ice albedo (only useful in coupled mode) 
    10197      !! 
    10298      !! References : Goosse, H. et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90. 
    10399      !!              Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 
     100      !!              These refs are now obsolete since everything has been revised 
     101      !!              The ref should be Rousset et al., 2015 
    104102      !!--------------------------------------------------------------------- 
    105       INTEGER, INTENT(in) ::   kt    ! number of iteration 
    106       ! 
    107       INTEGER  ::   ji, jj, jl, jk           ! dummy loop indices 
    108       REAL(wp) ::   zinda, zemp      ! local scalars 
    109       REAL(wp) ::   zf_mass         ! Heat flux associated with mass exchange ice->ocean (W.m-2) 
    110       REAL(wp) ::   zfcm1           ! New solar flux received by the ocean 
    111       REAL(wp), POINTER, DIMENSION(:,:,:) ::   zalb, zalbp     ! 2D/3D workspace 
     103      INTEGER, INTENT(in) ::   kt                                  ! number of iteration 
     104      INTEGER  ::   ji, jj, jl, jk                                 ! dummy loop indices 
     105      REAL(wp) ::   zqmass                                         ! Heat flux associated with mass exchange ice->ocean (W.m-2) 
     106      REAL(wp) ::   zqsr                                           ! New solar flux received by the ocean 
     107      ! 
     108      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zalb_cs, zalb_os     ! 2D/3D workspace 
    112109      !!--------------------------------------------------------------------- 
    113        
    114       IF( lk_cpl )   CALL wrk_alloc( jpi, jpj, jpl, zalb, zalbp ) 
    115110 
    116111      ! make calls for heat fluxes before it is modified 
    117       CALL iom_put( "qsr_oce" , qsr(:,:) * pfrld(:,:) )   !     solar flux at ocean surface 
    118       CALL iom_put( "qns_oce" , qns(:,:) * pfrld(:,:) )   ! non-solar flux at ocean surface 
    119       CALL iom_put( "qsr_ice" , SUM( qsr_ice(:,:,:) * old_a_i(:,:,:), dim=3 ) )  !     solar flux at ice surface 
    120       CALL iom_put( "qns_ice" , SUM( qns_ice(:,:,:) * old_a_i(:,:,:), dim=3 ) )  ! non-solar flux at ice surface 
    121       CALL iom_put( "qtr_ice" , SUM( ftr_ice(:,:,:) * old_a_i(:,:,:), dim=3 ) )  !     solar flux transmitted thru ice 
    122       CALL iom_put( "qt_oce"  , ( qsr(:,:) + qns(:,:) ) * pfrld(:,:) )   
    123       CALL iom_put( "qt_ice"  , SUM( ( qns_ice(:,:,:) + qsr_ice(:,:,:) ) * old_a_i(:,:,:), dim=3 ) ) 
     112      IF( iom_use('qsr_oce') )   CALL iom_put( "qsr_oce" , qsr_oce(:,:) * pfrld(:,:) )                                   !     solar flux at ocean surface 
     113      IF( iom_use('qns_oce') )   CALL iom_put( "qns_oce" , qns_oce(:,:) * pfrld(:,:) + qemp_oce(:,:) )                   ! non-solar flux at ocean surface 
     114      IF( iom_use('qsr_ice') )   CALL iom_put( "qsr_ice" , SUM( qsr_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) )                 !     solar flux at ice surface 
     115      IF( iom_use('qns_ice') )   CALL iom_put( "qns_ice" , SUM( qns_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) + qemp_ice(:,:) ) ! non-solar flux at ice surface 
     116      IF( iom_use('qtr_ice') )   CALL iom_put( "qtr_ice" , SUM( ftr_ice(:,:,:) * a_i_b(:,:,:), dim=3 ) )                 !     solar flux transmitted thru ice 
     117      IF( iom_use('qt_oce' ) )   CALL iom_put( "qt_oce"  , ( qsr_oce(:,:) + qns_oce(:,:) ) * pfrld(:,:) + qemp_oce(:,:) )   
     118      IF( iom_use('qt_ice' ) )   CALL iom_put( "qt_ice"  , SUM( ( qns_ice(:,:,:) + qsr_ice(:,:,:) )   & 
     119         &                                                      * a_i_b(:,:,:), dim=3 ) + qemp_ice(:,:) ) 
     120      IF( iom_use('qemp_oce' ) )   CALL iom_put( "qemp_oce"  , qemp_oce(:,:) )   
     121      IF( iom_use('qemp_ice' ) )   CALL iom_put( "qemp_ice"  , qemp_ice(:,:) )   
    124122 
    125123      ! pfrld is the lead fraction at the previous time step (actually between TRP and THD) 
     
    130128            !      heat flux at the ocean surface      ! 
    131129            !------------------------------------------! 
    132             zinda   = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ( 1._wp - pfrld(ji,jj) ) ) ) ! 1 if ice 
    133  
    134             ! Solar heat flux reaching the ocean = zfcm1 (W.m-2)  
     130            ! Solar heat flux reaching the ocean = zqsr (W.m-2)  
    135131            !--------------------------------------------------- 
    136             IF( lk_cpl ) THEN ! be carfeful: not been tested yet 
    137                ! original line 
    138                zfcm1 = qsr_tot(ji,jj) 
    139                !!!zfcm1 = qsr_tot(ji,jj) + ftr_ice(ji,jj) * ( 1._wp - pfrld(ji,jj) ) / ( 1._wp - zinda + zinda * iatte(ji,jj) ) 
    140                DO jl = 1, jpl 
    141                   zfcm1 = zfcm1 - ( qsr_ice(ji,jj,jl) - ftr_ice(ji,jj,jl) ) * old_a_i(ji,jj,jl) 
    142                END DO 
    143             ELSE 
    144                !!!zfcm1   = pfrld(ji,jj) * qsr(ji,jj)  + & 
    145                !!!     &    ( 1._wp - pfrld(ji,jj) ) * ftr_ice(ji,jj) / ( 1._wp - zinda + zinda * iatte(ji,jj) ) 
    146                zfcm1   = pfrld(ji,jj) * qsr(ji,jj) 
    147                DO jl = 1, jpl 
    148                   zfcm1   = zfcm1 + old_a_i(ji,jj,jl) * ftr_ice(ji,jj,jl) 
    149                END DO 
    150             ENDIF 
     132            zqsr = qsr_tot(ji,jj) 
     133            DO jl = 1, jpl 
     134               zqsr = zqsr - a_i_b(ji,jj,jl) * (  qsr_ice(ji,jj,jl) - ftr_ice(ji,jj,jl) )  
     135            END DO 
    151136 
    152137            ! Total heat flux reaching the ocean = hfx_out (W.m-2)  
    153138            !--------------------------------------------------- 
    154             zf_mass        = hfx_thd(ji,jj) + hfx_dyn(ji,jj) + hfx_res(ji,jj) ! heat flux from snow is 0 (T=0 degC) 
    155             hfx_out(ji,jj) = hfx_out(ji,jj) + zf_mass + zfcm1 
     139            zqmass         = hfx_thd(ji,jj) + hfx_dyn(ji,jj) + hfx_res(ji,jj) ! heat flux from snow is 0 (T=0 degC) 
     140            hfx_out(ji,jj) = hfx_out(ji,jj) + zqmass + zqsr 
     141 
     142            ! Add the residual from heat diffusion equation (W.m-2) 
     143            !------------------------------------------------------- 
     144            hfx_out(ji,jj) = hfx_out(ji,jj) + hfx_err_dif(ji,jj) 
    156145 
    157146            ! New qsr and qns used to compute the oceanic heat flux at the next time step 
    158147            !--------------------------------------------------- 
    159             qsr(ji,jj) = zfcm1                                       
    160             qns(ji,jj) = hfx_out(ji,jj) - zfcm1               
     148            qsr(ji,jj) = zqsr                                       
     149            qns(ji,jj) = hfx_out(ji,jj) - zqsr               
    161150 
    162151            !------------------------------------------! 
     
    171160            !                     Even if i see Ice melting as a FW and SALT flux 
    172161            !         
    173             !  computing freshwater exchanges at the ice/ocean interface 
    174             IF( lk_cpl ) THEN  
    175                zemp = - emp_tot(ji,jj) + emp_ice(ji,jj) * ( 1. - pfrld(ji,jj) )    &   ! 
    176                   &   + wfx_snw(ji,jj) 
    177             ELSE 
    178                zemp =   emp(ji,jj)     *           pfrld(ji,jj)            &   ! evaporation over oceanic fraction 
    179                   &   - tprecip(ji,jj) * ( 1._wp - pfrld(ji,jj) )          &   ! all precipitation reach the ocean 
    180                   &   + sprecip(ji,jj) * ( 1._wp - pfrld(ji,jj)**betas )       ! except solid precip intercepted by sea-ice 
    181             ENDIF 
    182  
    183162            ! mass flux from ice/ocean 
    184163            wfx_ice(ji,jj) = wfx_bog(ji,jj) + wfx_bom(ji,jj) + wfx_sum(ji,jj) + wfx_sni(ji,jj)   & 
     
    186165 
    187166            ! mass flux at the ocean/ice interface 
    188             fmmflx(ji,jj) = - wfx_ice(ji,jj) * rdt_ice                   ! F/M mass flux save at least for biogeochemical model 
    189             emp(ji,jj)    = zemp - wfx_ice(ji,jj) - wfx_snw(ji,jj) - wfx_sub(ji,jj)   ! mass flux + F/M mass flux (always ice/ocean mass exchange) 
     167            fmmflx(ji,jj) = - ( wfx_ice(ji,jj) + wfx_snw(ji,jj) ) * r1_rdtice  ! F/M mass flux save at least for biogeochemical model 
     168            emp(ji,jj)    = emp_oce(ji,jj) - wfx_ice(ji,jj) - wfx_snw(ji,jj)   ! mass flux + F/M mass flux (always ice/ocean mass exchange) 
    190169             
    191170         END DO 
     
    205184         snwice_mass_b(:,:) = snwice_mass(:,:)                   
    206185         ! new mass per unit area 
    207          snwice_mass  (:,:) = tms(:,:) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:)  )  
     186         snwice_mass  (:,:) = tmask(:,:,1) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:)  )  
    208187         ! time evolution of snow+ice mass 
    209188         snwice_fmass (:,:) = ( snwice_mass(:,:) - snwice_mass_b(:,:) ) * r1_rdtice 
     
    216195      tn_ice(:,:,:) = t_su(:,:,:)           ! Ice surface temperature                       
    217196 
    218       !------------------------------------------------! 
    219       !    Computation of snow/ice and ocean albedo    ! 
    220       !------------------------------------------------! 
    221       IF( lk_cpl ) THEN          ! coupled case 
    222          CALL albedo_ice( t_su, ht_i, ht_s, zalbp, zalb )                  ! snow/ice albedo 
    223          alb_ice(:,:,:) =  0.5_wp * zalbp(:,:,:) + 0.5_wp * zalb (:,:,:)   ! Ice albedo (mean clear and overcast skys) 
    224       ENDIF 
    225  
     197      !------------------------------------------------------------------------! 
     198      !    Snow/ice albedo (only if sent to coupler, useless in forced mode)   ! 
     199      !------------------------------------------------------------------------! 
     200      CALL wrk_alloc( jpi, jpj, jpl, zalb_cs, zalb_os )     
     201      CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os )  ! cloud-sky and overcast-sky ice albedos 
     202      alb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) 
     203      CALL wrk_dealloc( jpi, jpj, jpl, zalb_cs, zalb_os ) 
     204 
     205      ! conservation test 
     206      IF( ln_limdiahsb ) CALL lim_cons_final( 'limsbc' ) 
     207 
     208      ! control prints 
     209      IF( ln_icectl )   CALL lim_prt( kt, iiceprt, jiceprt, 3, ' - Final state lim_sbc - ' ) 
    226210 
    227211      IF(ln_ctl) THEN 
     
    231215         CALL prt_ctl( tab3d_1=tn_ice, clinfo1=' lim_sbc: tn_ice : ', kdim=jpl ) 
    232216      ENDIF 
    233       ! 
    234       IF( lk_cpl )   CALL wrk_dealloc( jpi, jpj, jpl, zalb, zalbp ) 
    235       !  
     217 
    236218   END SUBROUTINE lim_sbc_flx 
    237219 
     
    271253      ! 
    272254      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !==  Ice time-step only  ==!   (i.e. surface module time-step) 
    273 !CDIR NOVERRCHK 
    274255         DO jj = 2, jpjm1                             !* update the modulus of stress at ocean surface (T-point) 
    275 !CDIR NOVERRCHK 
    276256            DO ji = fs_2, fs_jpim1 
    277257               !                                               ! 2*(U_ice-U_oce) at T-point 
     
    323303      !! ** input   : Namelist namicedia 
    324304      !!------------------------------------------------------------------- 
    325       REAL(wp) :: zsum, zarea 
    326       ! 
    327305      INTEGER  ::   ji, jj, jk                       ! dummy loop indices 
    328306      REAL(wp) ::   zcoefu, zcoefv, zcoeff          ! local scalar 
     
    344322         END WHERE 
    345323      ENDIF 
    346       ! clem modif 
    347       IF( .NOT. ln_rstart ) THEN 
    348          iatte(:,:) = 1._wp 
    349          oatte(:,:) = 1._wp 
    350       ENDIF 
    351       ! 
    352       ! clem: snwice_mass in the restart file now 
     324      ! 
    353325      IF( .NOT. ln_rstart ) THEN 
    354326         !                                      ! embedded sea ice 
    355327         IF( nn_ice_embd /= 0 ) THEN            ! mass exchanges between ice and ocean (case 1 or 2) set the snow+ice mass 
    356             snwice_mass  (:,:) = tms(:,:) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:)  ) 
     328            snwice_mass  (:,:) = tmask(:,:,1) * ( rhosn * vt_s(:,:) + rhoic * vt_i(:,:)  ) 
    357329            snwice_mass_b(:,:) = snwice_mass(:,:) 
    358330         ELSE 
Note: See TracChangeset for help on using the changeset viewer.