- Timestamp:
- 2010-05-06T10:40:07+02:00 (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limsbc_2.F90
r1858 r1859 4 4 !! computation of the flux at the sea ice/ocean interface 5 5 !!====================================================================== 6 !! History : LIM ! 2000-01 (H. Goosse) Original code 7 !! 2.0 ! 2002-07 (C. Ethe, G. Madec) re-writing F90 8 !! - ! 2006-07 (G. Madec) surface module 9 !! 2.1 ! 2010-05 (Y. Aksenov, M. Vancoppenolle, G. Madec) add heat content exchanges 6 !! History : 1.0 ! 2000-01 (H. Goosse) Original code 7 !! 2.0 ! 2002-07 (C. Ethe, G. Madec) re-writing F90 8 !! - ! 2006-07 (G. Madec) surface module 9 !! - ! 2008-07 (C. Talandier,G. Madec) 2D fields for soce and sice 10 !! 2.1 ! 2010-05 (Y. Aksenov G. Madec) salt flux + heat associated with emp 10 11 !!---------------------------------------------------------------------- 11 12 #if defined key_lim2 … … 17 18 USE par_oce ! ocean parameters 18 19 USE dom_oce ! ocean domain 19 USE sbc_ice ! surface boundary condition20 USE sbc_oce ! surface boundary condition20 USE sbc_ice ! ice surface boundary condition 21 USE sbc_oce ! ocean surface boundary condition 21 22 USE phycst ! physical constants 22 USE ice_2 ! LIM sea-ice variables 23 USE albedo ! albedo parameters 24 USE ice_2 ! LIM-2 sea-ice variables 23 25 24 26 USE lbclnk ! ocean lateral boundary condition 25 27 USE in_out_manager ! I/O manager 28 USE iom ! 29 USE prtctl ! Print control 26 30 USE diaar5, ONLY : lk_diaar5 27 USE iom !28 USE albedo ! albedo parameters29 USE prtctl ! Print control30 31 USE cpl_oasis3, ONLY : lk_cpl 31 32 … … 33 34 PRIVATE 34 35 35 PUBLIC lim_sbc_2! called by sbc_ice_lim_236 PUBLIC lim_sbc_2 ! called by sbc_ice_lim_2 36 37 37 38 REAL(wp) :: epsi16 = 1.e-16 ! constant values 38 39 REAL(wp) :: rzero = 0.e0 39 40 REAL(wp) :: rone = 1.e0 40 REAL(wp), DIMENSION(jpi,jpj) :: soce_r 41 REAL(wp), DIMENSION(jpi,jpj) :: sice_r 41 REAL(wp), DIMENSION(jpi,jpj) :: soce_r, sice_r ! ocean and ice 2D constant salinity fields (used if lk_vvl=F) 42 42 43 43 !! * Substitutions 44 44 # include "vectopt_loop_substitute.h90" 45 45 !!---------------------------------------------------------------------- 46 !! LIM 2.0, UCL-LOCEAN-IPSL (2006)46 !! NEMO/LIM 3.3, UCL-LOCEAN-IPSL (2010) 47 47 !! $Id$ 48 48 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 62 62 !! - Update 63 63 !! 64 !! ** Outputs : - qsr : s ea heat flux: solar65 !! - qns : sea heat flux: non solar66 !! - emp : freshwater budget: volume flux67 !! - emps : freshwater budget: concentration/dillution64 !! ** Outputs : - qsr : solar heat flux 65 !! - qns : non-solar heat flux including heat content of mass flux 66 !! - emp : mass flux 67 !! - emps : salt flux due to Freezing/Melting 68 68 !! - utau : sea surface i-stress (ocean referential) 69 69 !! - vtau : sea surface j-stress (ocean referential) … … 75 75 !! Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 76 76 !!--------------------------------------------------------------------- 77 INTEGER :: kt ! number of iteration77 INTEGER, INTENT(in) :: kt ! number of iteration 78 78 !! 79 INTEGER :: ji, jj ! dummy loop indices 80 INTEGER :: ifvt, i1mfr, idfr ! some switches 81 INTEGER :: iflt, ial, iadv, ifral, ifrdv 82 INTEGER :: ii0, ii1, ij0, ij1 ! temporary integers 83 REAL(wp) :: zrdtir ! 1. / rdt_ice 84 REAL(wp) :: zqsr , zqns ! solar & non solar heat flux 85 REAL(wp) :: zinda ! switch for testing the values of ice concentration 86 REAL(wp) :: zfons ! salt exchanges at the ice/ocean interface 87 REAL(wp) :: zemp ! freshwater exchanges at the ice/ocean interface 88 REAL(wp) :: zfrldu, zfrldv ! lead fraction at U- & V-points 89 REAL(wp) :: zutau , zvtau ! lead fraction at U- & V-points 90 REAL(wp) :: zu_io , zv_io ! 2 components of the ice-ocean velocity 91 ! interface 2D --> 3D 92 REAL(wp), DIMENSION(jpi,jpj,1) :: zalb ! albedo of ice under overcast sky 93 REAL(wp), DIMENSION(jpi,jpj,1) :: zalbp ! albedo of ice under clear sky 79 INTEGER :: ji, jj ! dummy loop indices 80 INTEGER :: ifvt, idfr , iadv, i1mfr ! local integers 81 INTEGER :: iflt, ifrdv, ial , ifral ! - - 82 INTEGER :: ii0, ii1, ij0, ij1 ! - - 83 REAL(wp) :: zqsr, zqns, zqhc, zemp ! local scalars 84 REAL(wp) :: zinda, zswitch, zcd ! - - 85 REAL(wp) :: zfrldu, zutau, zu_io ! - - 86 REAL(wp) :: zfrldv, zvtau, zv_io ! - - 87 REAL(wp) :: zemp_snw, zfmm, zfsalt ! - - 94 88 REAL(wp) :: zsang, zmod, zztmp, zfm 95 REAL(wp), DIMENSION(jpi,jpj) :: ztio_u, ztio_v ! component of ocean stress below sea-ice at I-point 96 REAL(wp), DIMENSION(jpi,jpj) :: ztiomi ! module of ocean stress below sea-ice at I-point 97 REAL(wp), DIMENSION(jpi,jpj) :: zqnsoce ! save qns before its modification by ice model 98 89 REAL(wp), DIMENSION(jpi,jpj,1) :: zalb, zalbp ! 3D workspace 90 REAL(wp), DIMENSION(jpi,jpj) :: ztio_u, ztio_v ! 2D workspace 91 REAL(wp), DIMENSION(jpi,jpj) :: ztiomi, zqnsoce ! - - 99 92 !!--------------------------------------------------------------------- 100 101 zrdtir = 1. / rdt_ice 102 93 103 94 IF( kt == nit000 ) THEN 104 95 IF(lwp) WRITE(numout,*) 105 96 IF(lwp) WRITE(numout,*) 'lim_sbc_2 : LIM 2.0 sea-ice - surface boundary condition' 106 97 IF(lwp) WRITE(numout,*) '~~~~~~~~~ ' 107 108 soce_r(:,:) = soce 109 sice_r(:,:) = sice 110 ! 111 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN 112 ! ! ======================= 113 ! ! ORCA_R2 configuration 114 ! ! ======================= 115 ii0 = 145 ; ii1 = 180 ! Baltic Sea 98 ! ! 2D fields for constant ice and ocean salinities 99 soce_r(:,:) = soce ! in order to use different value in the Baltic sea 100 sice_r(:,:) = sice ! which is much less salty than polar regions 101 ! 102 IF( cp_cfg == "orca" ) THEN ! ORCA configuration 103 IF( jp_cfg == 2 ) THEN ! ORCA_R2 configuration 104 ii0 = 145 ; ii1 = 180 ! Baltic Sea 116 105 ij0 = 113 ; ij1 = 130 ; soce_r(mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 117 sice_r(mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 2.e0 106 sice_r(mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 0.50 107 !!gm add here the R1 R05 and R025 cases 108 !! ELSEIF( jp_cfg == 1 ) THEN ! ORCA_R1 configuration 109 !! ELSEIF( jp_cfg == 05 ) THEN ! ORCA_R05 configuration 110 !! ELSEIF( jp_cfg == 025 ) THEN ! ORCA_R025 configuration 111 !! 112 !!gm or better introduce the baltic change as a function of lat/lon of the baltic sea 113 !!end gm 114 ENDIF 118 115 ENDIF 119 116 ! 120 117 ENDIF 121 118 122 !------------------------------------------! 123 ! heat flux at the ocean surface ! 124 !------------------------------------------! 125 126 !!gm 127 !!gm CAUTION 128 !!gm re-verifies the non solar expression, especially over open ocen 129 !!gm 130 zqnsoce(:,:) = qns(:,:) 119 zqnsoce(:,:) = qns(:,:) ! save non-solar flux prior to its modification by ice-ocean fluxes (diag.) 120 ! 121 zswitch = 1 ! standard levitating sea-ice : salt exchange only 122 ! 123 !!gm ice embedment 124 ! SELECT CASE( nn_ice_embd ) ! levitating/embedded sea-ice 125 ! CASE( 0 ) ; zswitch = 1 ! standard levitating sea-ice : salt exchange only 126 ! CASE( 1, 2 ) ; zswitch = 0 ! other levitating sea-ice or embedded sea-ice : salt and volume fluxes 127 ! END SELECT ! 128 !!gm end embedment 129 131 130 DO jj = 1, jpj 132 131 DO ji = 1, jpi 132 ! !------------------------------------------! 133 ! ! heat flux at the ocean surface ! 134 ! !------------------------------------------! 133 135 zinda = 1.0 - MAX( rzero , SIGN( rone, - ( 1.0 - pfrld(ji,jj) ) ) ) 134 136 ifvt = zinda * MAX( rzero , SIGN( rone, - phicif(ji,jj) ) ) … … 138 140 ial = ifvt * i1mfr + ( 1 - ifvt ) * idfr 139 141 iadv = ( 1 - i1mfr ) * zinda 140 ifral = ( 1 - i1mfr * ( 1 - ial ) ) 141 ifrdv = ( 1 - ifral * ( 1 - ial ) ) * iadv 142 143 !!$ zinda = 1.0 - AINT( pfrld(ji,jj) ) ! = 0. if pure ocean else 1. (at previous time) 144 !!$ 145 !!$ i1mfr = 1.0 - AINT( frld(ji,jj) ) ! = 0. if pure ocean else 1. (at current time) 146 !!$ 147 !!$ IF( phicif(ji,jj) <= 0. ) THEN ; ifvt = zinda ! = 1. if (snow and no ice at previous time) else 0. ??? 148 !!$ ELSE ; ifvt = 0. 149 !!$ ENDIF 150 !!$ 151 !!$ IF( frld(ji,jj) >= pfrld(ji,jj) ) THEN ; idfr = 0. ! = 0. if lead fraction increases from previous to current 152 !!$ ELSE ; idfr = 1. 153 !!$ ENDIF 154 !!$ 155 !!$ iflt = zinda * (1 - i1mfr) * (1 - ifvt ) ! = 1. if ice (not only snow) at previous and pure ocean at current 156 !!$ 142 ifral = ( 1 - i1mfr * ( 1 - ial ) ) 143 ifrdv = ( 1 - ifral * ( 1 - ial ) ) * iadv 144 145 !!gm attempt to understand and comment the tricky flags used here.... 146 ! 147 !gm zinda = 1.0 - AINT( pfrld(ji,jj) ) ! = 0. free-ice ocean else 1. (after ice adv, but before ice thermo) 148 !gm i1mfr = 1.0 - AINT( frld(ji,jj) ) ! = 0. free-ice ocean else 1. (after ice thermo) 149 ! 150 !gm IF( phicif(ji,jj) <= 0. ) THEN ; ifvt = zinda ! = 1. if (snow and no ice at previous time) else 0. ??? 151 !gm ELSE ; ifvt = 0. ! correspond to a overmelting of snow in surface ablation 152 !gm ENDIF ! 153 ! 154 !gm IF( frld(ji,jj) >= pfrld(ji,jj) ) THEN ; idfr = 0. ! = 0. if lead fraction increases due to ice thermo 155 !gm ELSE ; idfr = 1. 156 !gm ENDIF 157 ! 158 !!$ iflt = zinda * (1 - i1mfr) * (1 - ifvt ) ! = 1. if ice (not only snow) at previous and pure ocean at current 159 ! 157 160 !!$ ial = ifvt * i1mfr + ( 1 - ifvt ) * idfr 158 161 !!$! snow no ice ice ice or nothing lead fraction increases 159 162 !!$! at previous now at previous 160 163 !!$! -> ice aera increases ??? -> ice aera decreases ??? 161 ! !$164 ! 162 165 !!$ iadv = ( 1 - i1mfr ) * zinda 163 166 !!$! pure ocean ice at 164 167 !!$! at current previous 165 168 !!$! -> = 1. if ice disapear between previous and current 166 ! !$169 ! 167 170 !!$ ifral = ( 1 - i1mfr * ( 1 - ial ) ) 168 171 !!$! ice at ??? 169 172 !!$! current 170 173 !!$! -> ??? 171 ! !$174 ! 172 175 !!$ ifrdv = ( 1 - ifral * ( 1 - ial ) ) * iadv 173 176 !!$! ice disapear 174 !!$ 175 !!$ 176 177 ! computation the solar flux at ocean surface 177 ! 178 ! 179 ! - computation the solar flux at ocean surface 178 180 #if defined key_coupled 179 181 zqsr = qsr_tot(ji,jj) + ( fstric(ji,jj) - qsr_ice(ji,jj,1) ) * ( 1.0 - pfrld(ji,jj) ) … … 181 183 zqsr = pfrld(ji,jj) * qsr(ji,jj) + ( 1. - pfrld(ji,jj) ) * fstric(ji,jj) 182 184 #endif 183 ! computation the non solar heat flux at ocean surface 185 ! 186 ! - computation the non solar heat flux at ocean surface 184 187 zqns = - ( 1. - thcm(ji,jj) ) * zqsr & ! part of the solar energy used in leads 185 & + iflt * ( fscmbq(ji,jj) + ffltbif(ji,jj) ) & 186 & + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * zrdtir & 187 & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * zrdtir 188 189 fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj) ! ??? 190 191 qsr (ji,jj) = zqsr ! solar heat flux 192 qns (ji,jj) = zqns - fdtcn(ji,jj) ! non solar heat flux 188 & + iflt * ( fscmbq(ji,jj) + ffltbif(ji,jj) ) & 189 & + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * r1_rdt_ice & 190 & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * r1_rdt_ice 191 192 ! - store residual heat flux (put in the ocean at the next time-step) 193 fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj) ! ??? 194 ! 195 ! - heat content of mass exchanged between ocean and sea-ice 196 zqhc = ( rdq_snw(ji,jj) + rdq_ice(ji,jj) ) * r1_rdt_ice ! heat flux due to sown & ice heat content exchanges 197 ! 198 qsr(ji,jj) = zqsr ! solar heat flux 199 qns(ji,jj) = zqns - fdtcn(ji,jj) + zqhc ! non solar heat flux 200 201 ! !------------------------------------------! 202 ! ! mass flux at the ocean surface ! 203 ! !------------------------------------------! 204 ! 205 ! mass flux at the ocean-atmosphere interface (open ocean fraction = leads area) 206 #if defined key_coupled 207 ! ! coupled mode: 208 zemp = + emp_tot(ji,jj) & ! net mass flux over the grid cell (ice+ocean area) 209 & - emp_ice(ji,jj) * ( 1. - pfrld(ji,jj) ) ! minus the mass flux intercepted by sea-ice 210 #else 211 ! ! forced mode: 212 zemp = + emp(ji,jj) * frld(ji,jj) & ! mass flux over open ocean fraction 213 & - tprecip(ji,jj) * ( 1. - frld(ji,jj) ) & ! liquid precip. over ice reaches directly the ocean 214 & + sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) & ! snow is intercepted by sea-ice (previous frld) 215 #endif 216 ! 217 ! mass flux at the ocean/ice interface (sea ice fraction) 218 zemp_snw = rdm_snw(ji,jj) * r1_rdt_ice ! snow melting = pure water that enters the ocean 219 zfmm = rdm_ice(ji,jj) * r1_rdt_ice ! Freezing minus Melting (F-M) 220 221 ! salt flux at the ice/ocean interface (sea ice fraction) [PSU*kg/m2/s] 222 zfsalt = - sice_r(ji,jj) * zfmm ! F-M salt exchange 223 zcd = soce_r(ji,jj) * zfmm ! concentration/dilution term due to F-M 224 ! 225 ! salt flux only : add concentration dilution term in salt flux and no F-M term in volume flux 226 ! salt and mass fluxes : non concentartion dilution term in salt flux and add F-M term in volume flux 227 emps(ji,jj) = zfsalt + zswitch * zcd ! salt flux (+ C/D if no ice/ocean mass exchange) 228 emp (ji,jj) = zemp + zemp_snw + ( 1.- zswitch) * zfmm ! mass flux (- F/M mass flux if no ice/ocean mass exchange) 229 ! 193 230 END DO 194 231 END DO 232 195 233 196 234 CALL iom_put( 'hflx_ice_cea', - fdtcn(:,:) ) 197 235 CALL iom_put( 'qns_io_cea', qns(:,:) - zqnsoce(:,:) * pfrld(:,:) ) 198 236 CALL iom_put( 'qsr_io_cea', fstric(:,:) * (1. - pfrld(:,:)) ) 199 200 !------------------------------------------!201 ! mass flux at the ocean surface !202 !------------------------------------------!203 204 !!gm205 !!gm CAUTION206 !!gm re-verifies the emp & emps expression, especially the absence of 1-frld on zfm207 !!gm208 DO jj = 1, jpj209 DO ji = 1, jpi210 211 #if defined key_coupled212 zemp = emp_tot(ji,jj) - emp_ice(ji,jj) * ( 1. - pfrld(ji,jj) ) & !213 & + rdm_snw(ji,jj) * zrdtir ! freshwaterflux due to snow melting214 #else215 !!$ ! computing freshwater exchanges at the ice/ocean interface216 !!$ zpme = - evap(ji,jj) * frld(ji,jj) & ! evaporation over oceanic fraction217 !!$ & + tprecip(ji,jj) & ! total precipitation218 !!$ & - sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) & ! remov. snow precip over ice219 !!$ & - rdm_snw(ji,jj) / rdt_ice ! freshwaterflux due to snow melting220 ! computing freshwater exchanges at the ice/ocean interface221 zemp = + emp(ji,jj) * frld(ji,jj) & ! e-p budget over open ocean fraction222 & - tprecip(ji,jj) * ( 1. - frld(ji,jj) ) & ! liquid precipitation reaches directly the ocean223 & + sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) & ! taking into account change in ice cover within the time step224 & + rdm_snw(ji,jj) * zrdtir ! freshwaterflux due to snow melting225 ! ! ice-covered fraction:226 #endif227 228 ! computing salt exchanges at the ice/ocean interface229 zfons = ( soce_r(ji,jj) - sice_r(ji,jj) ) * ( rdm_ice(ji,jj) * zrdtir )230 231 ! converting the salt flux from ice to a freshwater flux from ocean232 zfm = zfons / ( sss_m(ji,jj) + epsi16 )233 234 emps(ji,jj) = zemp + zfm ! surface ocean concentration/dilution effect (use on SSS evolution)235 emp (ji,jj) = zemp ! surface ocean volume flux (use on sea-surface height evolution)236 237 END DO238 END DO239 237 240 238 IF( lk_diaar5 ) THEN … … 250 248 IF ( ln_limdyn ) THEN ! Update the stress over ice-over area (only in ice-dynamic case) 251 249 ! ! otherwise the atmosphere-ocean stress is used everywhere 252 250 ! 253 251 ! ... ice stress over ocean with a ice-ocean rotation angle (at I-point) 254 252 !CDIR NOVERRCHK … … 290 288 END DO 291 289 END DO 292 293 ! boundary condition on the stress (utau,vtau,taum) 294 CALL lbc_lnk( utau, 'U', -1. ) 295 CALL lbc_lnk( vtau, 'V', -1. ) 290 CALL lbc_lnk( utau, 'U', -1. ) ; CALL lbc_lnk( vtau, 'V', -1. ) ! lateral boundary conditions 296 291 CALL lbc_lnk( taum, 'T', 1. ) 297 292 ! 298 293 ENDIF 299 294 300 295 !-----------------------------------------------! 301 ! Coupling variables!296 ! Storing the transmitted variables ! 302 297 !-----------------------------------------------! 303 304 IF ( lk_cpl ) THEN 305 ! Ice surface temperature 298 !!gm where this is done ????? ==>>> limthd_2 not logic ??? 299 !!gm fr_i(:,:) = 1.0 - frld(:,:) ! sea-ice fraction 300 !!gm 301 302 IF ( lk_cpl ) THEN ! coupled mode : 306 303 tn_ice(:,:,1) = sist(:,:) ! sea-ice surface temperature 307 ! Computation ofsnow/ice and ocean albedo304 ! ! snow/ice and ocean albedo 308 305 CALL albedo_ice( tn_ice, reshape( hicif, (/jpi,jpj,1/) ), reshape( hsnif, (/jpi,jpj,1/) ), zalbp, zalb ) 309 306 alb_ice(:,:,1) = 0.5 * ( zalbp(:,:,1) + zalb (:,:,1) ) ! Ice albedo (mean clear and overcast skys) 307 ! 310 308 CALL iom_put( "icealb_cea", alb_ice(:,:,1) * fr_i(:,:) ) ! ice albedo 311 309 ENDIF … … 318 316 CALL prt_ctl(tab2d_1=fr_i , clinfo1=' lim_sbc: fr_i : ', tab2d_2=tn_ice(:,:,1), clinfo2=' tn_ice : ') 319 317 ENDIF 320 321 318 ! 319 END SUBROUTINE lim_sbc_2 322 320 323 321 #else
Note: See TracChangeset
for help on using the changeset viewer.