MODULE limsbc_2 !!====================================================================== !! *** MODULE limsbc_2 *** !! computation of the flux at the sea ice/ocean interface !!====================================================================== !! History : 00-01 (H. Goosse) Original code !! 02-07 (C. Ethe, G. Madec) re-writing F90 !! 06-07 (G. Madec) surface module !!---------------------------------------------------------------------- #if defined key_lim2 !!---------------------------------------------------------------------- !! 'key_lim2' LIM 2.0 sea-ice model !!---------------------------------------------------------------------- !!---------------------------------------------------------------------- !! lim_sbc_2 : flux at the ice / ocean interface !!---------------------------------------------------------------------- USE par_oce ! ocean parameters USE dom_oce ! ocean domain USE sbc_ice ! surface boundary condition USE sbc_oce ! surface boundary condition USE phycst ! physical constants USE ice_2 ! LIM sea-ice variables USE lbclnk ! ocean lateral boundary condition USE in_out_manager ! I/O manager USE diaar5, ONLY : lk_diaar5 USE iom ! USE albedo ! albedo parameters USE prtctl ! Print control USE cpl_oasis3, ONLY : lk_cpl IMPLICIT NONE PRIVATE PUBLIC lim_sbc_2 ! called by sbc_ice_lim_2 REAL(wp) :: epsi16 = 1.e-16 ! constant values REAL(wp) :: rzero = 0.e0 REAL(wp) :: rone = 1.e0 REAL(wp), DIMENSION(jpi,jpj) :: soce_r REAL(wp), DIMENSION(jpi,jpj) :: sice_r !! * Substitutions # include "vectopt_loop_substitute.h90" !!---------------------------------------------------------------------- !! LIM 2.0, UCL-LOCEAN-IPSL (2006) !! $Id$ !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) !!---------------------------------------------------------------------- CONTAINS SUBROUTINE lim_sbc_2( kt ) !!------------------------------------------------------------------- !! *** ROUTINE lim_sbc_2 *** !! !! ** Purpose : Update surface ocean boundary condition over areas !! that are at least partially covered by sea-ice !! !! ** Action : - comput. of the momentum, heat and freshwater/salt !! fluxes at the ice-ocean interface. !! - Update !! !! ** Outputs : - qsr : sea heat flux: solar !! - qns : sea heat flux: non solar !! - emp : freshwater budget: volume flux !! - emps : freshwater budget: concentration/dillution !! - utau : sea surface i-stress (ocean referential) !! - vtau : sea surface j-stress (ocean referential) !! - fr_i : ice fraction !! - tn_ice : sea-ice surface temperature !! - alb_ice : sea-ice alberdo (lk_cpl=T) !! !! References : Goosse, H. et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90. !! Tartinville et al. 2001 Ocean Modelling, 3, 95-108. !!--------------------------------------------------------------------- INTEGER :: kt ! number of iteration !! INTEGER :: ji, jj ! dummy loop indices INTEGER :: ifvt, i1mfr, idfr ! some switches INTEGER :: iflt, ial, iadv, ifral, ifrdv REAL(wp) :: zrdtir ! 1. / rdt_ice REAL(wp) :: zqsr , zqns ! solar & non solar heat flux REAL(wp) :: zinda ! switch for testing the values of ice concentration REAL(wp) :: zfons ! salt exchanges at the ice/ocean interface REAL(wp) :: zemp ! freshwater exchanges at the ice/ocean interface REAL(wp) :: zfrldu, zfrldv ! lead fraction at U- & V-points REAL(wp) :: zutau , zvtau ! lead fraction at U- & V-points REAL(wp) :: zu_io , zv_io ! 2 components of the ice-ocean velocity ! interface 2D --> 3D REAL(wp), DIMENSION(jpi,jpj,1) :: zalb ! albedo of ice under overcast sky REAL(wp), DIMENSION(jpi,jpj,1) :: zalbp ! albedo of ice under clear sky REAL(wp) :: zsang, zmod, zztmp, zfm REAL(wp), DIMENSION(jpi,jpj) :: ztio_u, ztio_v ! component of ocean stress below sea-ice at I-point REAL(wp), DIMENSION(jpi,jpj) :: ztiomi ! module of ocean stress below sea-ice at I-point REAL(wp), DIMENSION(jpi,jpj) :: zqnsoce ! save qns before its modification by ice model !!--------------------------------------------------------------------- zrdtir = 1. / rdt_ice IF( kt == nit000 ) THEN IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) 'lim_sbc_2 : LIM 2.0 sea-ice - surface boundary condition' IF(lwp) WRITE(numout,*) '~~~~~~~~~ ' soce_r(:,:) = soce sice_r(:,:) = sice ! IF( cp_cfg == "orca" ) THEN ! ocean/ice salinity in the Baltic sea DO jj = 1, jpj DO ji = 1, jpi IF( glamt(ji,jj) >= 14. .AND. glamt(ji,jj) <= 32. .AND. gphit(ji,jj) >= 54. .AND. gphit(ji,jj) <= 66. ) THEN soce_r(ji,jj) = 4.e0 sice_r(ji,jj) = 2.e0 END IF END DO END DO ! END IF END IF !------------------------------------------! ! heat flux at the ocean surface ! !------------------------------------------! !!gm !!gm CAUTION !!gm re-verifies the non solar expression, especially over open ocen !!gm zqnsoce(:,:) = qns(:,:) DO jj = 1, jpj DO ji = 1, jpi zinda = 1.0 - MAX( rzero , SIGN( rone, - ( 1.0 - pfrld(ji,jj) ) ) ) ifvt = zinda * MAX( rzero , SIGN( rone, - phicif(ji,jj) ) ) i1mfr = 1.0 - MAX( rzero , SIGN( rone, - ( 1.0 - frld(ji,jj) ) ) ) idfr = 1.0 - MAX( rzero , SIGN( rone, frld(ji,jj) - pfrld(ji,jj) ) ) iflt = zinda * (1 - i1mfr) * (1 - ifvt ) ial = ifvt * i1mfr + ( 1 - ifvt ) * idfr iadv = ( 1 - i1mfr ) * zinda ifral = ( 1 - i1mfr * ( 1 - ial ) ) ifrdv = ( 1 - ifral * ( 1 - ial ) ) * iadv !!$ zinda = 1.0 - AINT( pfrld(ji,jj) ) ! = 0. if pure ocean else 1. (at previous time) !!$ !!$ i1mfr = 1.0 - AINT( frld(ji,jj) ) ! = 0. if pure ocean else 1. (at current time) !!$ !!$ IF( phicif(ji,jj) <= 0. ) THEN ; ifvt = zinda ! = 1. if (snow and no ice at previous time) else 0. ??? !!$ ELSE ; ifvt = 0. !!$ ENDIF !!$ !!$ IF( frld(ji,jj) >= pfrld(ji,jj) ) THEN ; idfr = 0. ! = 0. if lead fraction increases from previous to current !!$ ELSE ; idfr = 1. !!$ ENDIF !!$ !!$ iflt = zinda * (1 - i1mfr) * (1 - ifvt ) ! = 1. if ice (not only snow) at previous and pure ocean at current !!$ !!$ ial = ifvt * i1mfr + ( 1 - ifvt ) * idfr !!$! snow no ice ice ice or nothing lead fraction increases !!$! at previous now at previous !!$! -> ice aera increases ??? -> ice aera decreases ??? !!$ !!$ iadv = ( 1 - i1mfr ) * zinda !!$! pure ocean ice at !!$! at current previous !!$! -> = 1. if ice disapear between previous and current !!$ !!$ ifral = ( 1 - i1mfr * ( 1 - ial ) ) !!$! ice at ??? !!$! current !!$! -> ??? !!$ !!$ ifrdv = ( 1 - ifral * ( 1 - ial ) ) * iadv !!$! ice disapear !!$ !!$ ! computation the solar flux at ocean surface #if defined key_coupled zqsr = qsr_tot(ji,jj) + ( fstric(ji,jj) - qsr_ice(ji,jj,1) ) * ( 1.0 - pfrld(ji,jj) ) #else zqsr = pfrld(ji,jj) * qsr(ji,jj) + ( 1. - pfrld(ji,jj) ) * fstric(ji,jj) #endif ! computation the non solar heat flux at ocean surface zqns = - ( 1. - thcm(ji,jj) ) * zqsr & ! part of the solar energy used in leads & + iflt * ( fscmbq(ji,jj) + ffltbif(ji,jj) ) & & + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * zrdtir & & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * zrdtir fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj) ! ??? qsr (ji,jj) = zqsr ! solar heat flux qns (ji,jj) = zqns - fdtcn(ji,jj) ! non solar heat flux END DO END DO CALL iom_put( 'hflx_ice_cea', - fdtcn(:,:) ) CALL iom_put( 'qns_io_cea', qns(:,:) - zqnsoce(:,:) * pfrld(:,:) ) CALL iom_put( 'qsr_io_cea', fstric(:,:) * (1. - pfrld(:,:)) ) !------------------------------------------! ! mass flux at the ocean surface ! !------------------------------------------! !!gm !!gm CAUTION !!gm re-verifies the emp & emps expression, especially the absence of 1-frld on zfm !!gm DO jj = 1, jpj DO ji = 1, jpi #if defined key_coupled zemp = emp_tot(ji,jj) - emp_ice(ji,jj) * ( 1. - pfrld(ji,jj) ) & ! & + rdmsnif(ji,jj) * zrdtir ! freshwaterflux due to snow melting #else !!$ ! computing freshwater exchanges at the ice/ocean interface !!$ zpme = - evap(ji,jj) * frld(ji,jj) & ! evaporation over oceanic fraction !!$ & + tprecip(ji,jj) & ! total precipitation !!$ & - sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) & ! remov. snow precip over ice !!$ & - rdmsnif(ji,jj) / rdt_ice ! freshwaterflux due to snow melting ! computing freshwater exchanges at the ice/ocean interface zemp = + emp(ji,jj) * frld(ji,jj) & ! e-p budget over open ocean fraction & - tprecip(ji,jj) * ( 1. - frld(ji,jj) ) & ! liquid precipitation reaches directly the ocean & + sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) & ! taking into account change in ice cover within the time step & + rdmsnif(ji,jj) * zrdtir ! freshwaterflux due to snow melting ! ! ice-covered fraction: #endif ! computing salt exchanges at the ice/ocean interface zfons = ( soce_r(ji,jj) - sice_r(ji,jj) ) * ( rdmicif(ji,jj) * zrdtir ) ! converting the salt flux from ice to a freshwater flux from ocean zfm = zfons / ( sss_m(ji,jj) + epsi16 ) emps(ji,jj) = zemp + zfm ! surface ocean concentration/dilution effect (use on SSS evolution) emp (ji,jj) = zemp ! surface ocean volume flux (use on sea-surface height evolution) END DO END DO IF( lk_diaar5 ) THEN CALL iom_put( 'isnwmlt_cea' , rdmsnif(:,:) * zrdtir ) CALL iom_put( 'fsal_virt_cea', soce_r(:,:) * rdmicif(:,:) * zrdtir ) CALL iom_put( 'fsal_real_cea', - sice_r(:,:) * rdmicif(:,:) * zrdtir ) ENDIF !------------------------------------------! ! momentum flux at the ocean surface ! !------------------------------------------! IF ( ln_limdyn ) THEN ! Update the stress over ice-over area (only in ice-dynamic case) ! ! otherwise the atmosphere-ocean stress is used everywhere ! ... ice stress over ocean with a ice-ocean rotation angle (at I-point) !CDIR NOVERRCHK DO jj = 1, jpj !CDIR NOVERRCHK DO ji = 1, jpi ! ... change the cosinus angle sign in the south hemisphere zsang = SIGN(1.e0, gphif(ji,jj) ) * sangvg ! ... ice velocity relative to the ocean at I-point zu_io = u_ice(ji,jj) - u_oce(ji,jj) zv_io = v_ice(ji,jj) - v_oce(ji,jj) zmod = SQRT( zu_io * zu_io + zv_io * zv_io ) zztmp = rhoco * zmod ! ... components of ice stress over ocean with a ice-ocean rotation angle (at I-point) ztio_u(ji,jj) = zztmp * ( cangvg * zu_io - zsang * zv_io ) ztio_v(ji,jj) = zztmp * ( cangvg * zv_io + zsang * zu_io ) ! ... module of ice stress over ocean (at I-point) ztiomi(ji,jj) = zztmp * zmod ! END DO END DO DO jj = 2, jpjm1 DO ji = 2, jpim1 ! NO vector opt. ! ... components of ice-ocean stress at U and V-points (from I-point values) zutau = 0.5 * ( ztio_u(ji+1,jj) + ztio_u(ji+1,jj+1) ) zvtau = 0.5 * ( ztio_v(ji,jj+1) + ztio_v(ji+1,jj+1) ) ! ... open-ocean (lead) fraction at U- & V-points (from T-point values) zfrldu = 0.5 * ( frld(ji,jj) + frld(ji+1,jj ) ) zfrldv = 0.5 * ( frld(ji,jj) + frld(ji ,jj+1) ) ! ... update components of surface ocean stress (ice-cover wheighted) utau(ji,jj) = zfrldu * utau(ji,jj) + ( 1. - zfrldu ) * zutau vtau(ji,jj) = zfrldv * vtau(ji,jj) + ( 1. - zfrldv ) * zvtau ! ... module of ice-ocean stress at T-points (from I-point values) zztmp = 0.25 * ( ztiomi(ji,jj) + ztiomi(ji+1,jj) + ztiomi(ji,jj+1) + ztiomi(ji+1,jj+1) ) ! ... update module of surface ocean stress (ice-cover wheighted) taum(ji,jj) = frld(ji,jj) * taum(ji,jj) + ( 1. - frld(ji,jj) ) * zztmp ! END DO END DO ! boundary condition on the stress (utau,vtau,taum) CALL lbc_lnk( utau, 'U', -1. ) CALL lbc_lnk( vtau, 'V', -1. ) CALL lbc_lnk( taum, 'T', 1. ) ENDIF !-----------------------------------------------! ! Coupling variables ! !-----------------------------------------------! IF ( lk_cpl ) THEN ! Ice surface temperature tn_ice(:,:,1) = sist(:,:) ! sea-ice surface temperature ! Computation of snow/ice and ocean albedo CALL albedo_ice( tn_ice, reshape( hicif, (/jpi,jpj,1/) ), reshape( hsnif, (/jpi,jpj,1/) ), zalbp, zalb ) alb_ice(:,:,1) = 0.5 * ( zalbp(:,:,1) + zalb (:,:,1) ) ! Ice albedo (mean clear and overcast skys) CALL iom_put( "icealb_cea", alb_ice(:,:,1) * fr_i(:,:) ) ! ice albedo ENDIF IF(ln_ctl) THEN CALL prt_ctl(tab2d_1=qsr , clinfo1=' lim_sbc: qsr : ', tab2d_2=qns , clinfo2=' qns : ') CALL prt_ctl(tab2d_1=emp , clinfo1=' lim_sbc: emp : ', tab2d_2=emps , clinfo2=' emps : ') CALL prt_ctl(tab2d_1=utau , clinfo1=' lim_sbc: utau : ', mask1=umask, & & tab2d_2=vtau , clinfo2=' vtau : ' , mask2=vmask ) CALL prt_ctl(tab2d_1=fr_i , clinfo1=' lim_sbc: fr_i : ', tab2d_2=tn_ice(:,:,1), clinfo2=' tn_ice : ') ENDIF END SUBROUTINE lim_sbc_2 #else !!---------------------------------------------------------------------- !! Default option : Dummy module NO LIM 2.0 sea-ice model !!---------------------------------------------------------------------- CONTAINS SUBROUTINE lim_sbc_2 ! Dummy routine END SUBROUTINE lim_sbc_2 #endif !!====================================================================== END MODULE limsbc_2