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.
limflx.F90 in branches/dev_002_LIM/NEMO/LIM_SRC_3 – NEMO

source: branches/dev_002_LIM/NEMO/LIM_SRC_3/limflx.F90 @ 830

Last change on this file since 830 was 825, checked in by ctlod, 16 years ago

dev_002_LIM : add the LIM 3.0 component, see ticketr: #71

File size: 17.0 KB
Line 
1MODULE limflx
2   !!======================================================================
3   !!                       ***  MODULE limflx   ***
4   !!           computation of the flux at the sea ice/ocean interface
5   !!======================================================================
6#if defined key_lim3
7   !!----------------------------------------------------------------------
8   !!   'key_lim3'                                     LIM sea-ice model
9   !!----------------------------------------------------------------------
10   !!   lim_flx  : flux at the ice / ocean interface
11   !! * Modules used
12   USE par_oce
13   USE phycst
14   USE ocfzpt
15   USE ice_oce
16   USE flx_oce
17   USE dom_oce
18   USE ice
19   USE flxblk
20   USE lbclnk
21   USE in_out_manager
22   USE albedo
23   USE limicepoints
24   USE par_ice
25   USE prtctl           ! Print control
26
27   IMPLICIT NONE
28   PRIVATE
29
30   !! * Routine accessibility
31   PUBLIC lim_flx       ! called by lim_step
32
33  !! * Module variables
34     REAL(wp)  ::            &  ! constant values
35         epsi16 = 1e-16   ,  &
36         rzero  = 0.0    ,  &
37         rone   = 1.0
38   !! * Substitutions
39#  include "vectopt_loop_substitute.h90"
40   !!----------------------------------------------------------------------
41   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)
42   !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limflx.F90,v 1.6 2005/03/27 18:34:41 opalod Exp $
43   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
44   !!----------------------------------------------------------------------
45CONTAINS
46
47   SUBROUTINE lim_flx
48      !!-------------------------------------------------------------------
49      !!                ***  ROUTINE lim_flx ***
50      !!               
51      !! 
52      !! ** Purpose : Computes the mass and heat fluxes to the ocean
53      !!         
54      !! ** Action  : - Initialisation of some variables
55      !!              - comput. of the fluxes at the sea ice/ocean interface
56      !!     
57      !! ** Outputs : - fsolar  : solar heat flux at sea ice/ocean interface
58      !!              - fnsolar : non solar heat flux
59      !!              - fsalt   : salt flux at sea ice/ocean interface
60      !!              - fmass   : freshwater flux at sea ice/ocean interface
61      !!
62      !!
63      !! ** References :
64      !!       H. Goosse et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90
65      !!         original    : 00-01 (LIM)
66      !!         addition    : 02-07 (C. Ethe, G. Madec)
67      !!---------------------------------------------------------------------
68
69      !! * Modules used
70      !! * Local variables
71      INTEGER ::   ji, jj, jl         ! dummy loop indices
72
73      INTEGER ::   &
74         ifvt, i1mfr, idfr ,   &  ! some switches
75         iflt, ial, iadv, ifral, ifrdv
76     
77      REAL(wp) ::   &
78         zinda  ,              &  ! switch for testing the values of ice concentration
79!!         zfcm1  ,              &  ! solar  heat fluxes
80!!         zfcm2  ,              &  !  non solar heat fluxes
81           zfold  ,            &
82           zdfseqv,            &    ! increment to the ice salt flux
83#if defined key_lim_fdd   
84         zfons,                &  ! salt exchanges at the ice/ocean interface
85         zpme                     ! freshwater exchanges at the ice/ocean interface
86#else
87         zprs  , zfons,        &  ! salt exchanges at the ice/ocean interface
88         zpmess                   ! freshwater exchanges at the ice/ocean interface
89#endif
90      REAL(wp), DIMENSION(jpi,jpj) ::  &
91         zfcm1  ,              &  ! solar  heat fluxes
92         zfcm2                    !  non solar heat fluxes     
93#if defined key_coupled   
94      REAL(wp), DIMENSION(jpi,jpj) ::  &
95         zalb  ,               &  ! albedo of ice under overcast sky
96         zalcn ,               &  ! albedo of ocean under overcast sky
97         zalbp ,               &  ! albedo of ice under clear sky
98         zaldum                   ! albedo of ocean under clear sky
99#endif
100      !!---------------------------------------------------------------------
101     
102      !---------------------------------!
103      !      Sea ice/ocean interface    !
104      !---------------------------------!
105       
106      !      heat flux at the ocean surface
107      !-------------------------------------------------------
108      ! pfrld is the lead fraction at the previous time step (actually between TRP and THD)
109      ! changed to old_frld and old ht_i
110       
111      DO jj = 1, jpj
112         DO ji = 1, jpi
113            zinda   = 1.0 - MAX( rzero , SIGN( rone , - ( 1.0 - pfrld(ji,jj) ) ) )
114            ifvt    = zinda  *  MAX( rzero , SIGN( rone, -phicif  (ji,jj) ) )  !subscripts are bad here
115            i1mfr   = 1.0 - MAX( rzero , SIGN( rone ,  - ( at_i(ji,jj)       ) ) )
116            idfr    = 1.0 - MAX( rzero , SIGN( rone , ( 1.0 - at_i(ji,jj) ) - pfrld(ji,jj) ) )
117            iflt    = zinda  * (1 - i1mfr) * (1 - ifvt )
118            ial     = ifvt   * i1mfr + ( 1 - ifvt ) * idfr
119            iadv    = ( 1  - i1mfr ) * zinda
120            ifral   = ( 1  - i1mfr * ( 1 - ial ) )   
121            ifrdv   = ( 1  - ifral * ( 1 - ial ) ) * iadv 
122
123       ! switch --- 1.0 ---------------- 0.0 --------------------
124       ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
125       ! zinda   | if pfrld = 1       | if pfrld < 1            |
126       !  -> ifvt| if pfrld old_ht_i
127       ! i1mfr   | if frld = 1        | if frld  < 1            |
128       ! idfr    | if frld <= pfrld    | if frld > pfrld        |
129       ! iflt    |
130       ! ial     |
131       ! iadv    |
132       ! ifral
133       ! ifrdv
134
135            !   computation the solar flux at ocean surface
136            zfcm1(ji,jj)   = pfrld(ji,jj) * qsr_oce(ji,jj)  + ( 1. - pfrld(ji,jj) ) * fstric(ji,jj)
137                ! fstric     Solar flux transmitted trough the ice
138                ! qsr_oce    Net short wave heat flux on free ocean
139! new line
140            fscmbq(ji,jj) = ( 1.0 - pfrld(ji,jj) ) * fstric(ji,jj)
141
142            !  computation the non solar heat flux at ocean surface
143            zfcm2(ji,jj) = - zfcm1(ji,jj)                  &
144               &           + iflt    * ( fscmbq(ji,jj) )   & ! total abl -> fscmbq is given to the ocean
145! fscmbq and ffltbif are obsolete
146!              &           + iflt * ffltbif(ji,jj) !!! only if one category is used
147               &           + ifral   * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) / rdt_ice   &
148               &           + ifrdv   * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) / rdt_ice                     &
149               &           + fhmec(ji,jj)     & ! new contribution due to snow melt in ridging!!
150               &           + fheat_rpo(ji,jj) & ! contribution from ridge formation
151               &           + fheat_res(ji,jj)
152                ! fscmbq  Part of the solar radiation transmitted through the ice and going to the ocean
153                !         computed in limthd_zdf.F90
154                ! ffltbif Total heat content of the ice (brine pockets+ice) / delta_t
155                ! qcmif   Energy needed to bring the ocean surface layer until its freezing (ok)
156                ! qldif   heat balance of the lead (or of the open ocean)
157                ! qfvbq   i think this is wrong!
158                ! ---> Array used to store energy in case of total lateral ablation
159                ! qfvbq latent heat uptake/release after accretion/ablation
160                ! qdtcn Energy from the turbulent oceanic heat flux heat flux coming in the lead
161
162            IF ( num_sal .EQ. 2 ) zfcm2(ji,jj) = zfcm2(ji,jj) + &
163                                  fhbri(ji,jj) ! new contribution due to brine drainage
164
165            ! bottom radiative component is sent to the computation of the
166            ! oceanic heat flux
167            fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj)     
168
169            ! used to compute the oceanic heat flux at the next time step
170            fsolar (ji,jj) = zfcm1(ji,jj)                                       ! solar heat flux
171            fnsolar(ji,jj) = zfcm2(ji,jj) - fdtcn(ji,jj)                        ! non solar heat flux
172!                           ! fdtcn : turbulent oceanic heat flux
173
174
175            IF ( ( ji .EQ. jiindex ) .AND. ( jj .EQ. jjindex) ) THEN
176               WRITE(numout,*) ' lim_flx : heat fluxes '
177               WRITE(numout,*) ' fsolar    : ', fsolar(jiindex,jjindex)
178               WRITE(numout,*) ' zfcm1     : ', zfcm1(jiindex,jjindex)
179               WRITE(numout,*) ' pfrld     : ', pfrld(jiindex,jjindex)
180               WRITE(numout,*) ' qsr_oce   : ', qsr_oce(jiindex,jjindex)
181               WRITE(numout,*) ' fstric    : ', fstric (jiindex,jjindex)
182               WRITE(numout,*)
183               WRITE(numout,*) ' fnsolar   : ', fnsolar(jiindex,jjindex)
184               WRITE(numout,*) ' zfcm2     : ', zfcm2(jiindex,jjindex)
185               WRITE(numout,*) ' zfcm1     : ', zfcm1(jiindex,jjindex)
186               WRITE(numout,*) ' ifral     : ', ifral
187               WRITE(numout,*) ' ial       : ', ial 
188               WRITE(numout,*) ' qcmif     : ', qcmif(jiindex,jjindex)
189               WRITE(numout,*) ' qldif     : ', qldif(jiindex,jjindex)
190               WRITE(numout,*) ' qcmif / dt: ', qcmif(jiindex,jjindex) / rdt_ice
191               WRITE(numout,*) ' qldif / dt: ', qldif(jiindex,jjindex) / rdt_ice
192               WRITE(numout,*) ' ifrdv     : ', ifrdv
193               WRITE(numout,*) ' qfvbq     : ', qfvbq(jiindex,jjindex)
194               WRITE(numout,*) ' qdtcn     : ', qdtcn(jiindex,jjindex)
195               WRITE(numout,*) ' qfvbq / dt: ', qfvbq(jiindex,jjindex) / rdt_ice
196               WRITE(numout,*) ' qdtcn / dt: ', qdtcn(jiindex,jjindex) / rdt_ice
197               WRITE(numout,*) ' '
198               WRITE(numout,*) ' fdtcn     : ', fdtcn(jiindex,jjindex)
199               WRITE(numout,*) ' fhmec     : ', fhmec(jiindex,jjindex)
200               WRITE(numout,*) ' fheat_rpo : ', fheat_rpo(jiindex,jjindex)
201               WRITE(numout,*) ' fhbri     : ', fhbri(jiindex,jjindex)
202               WRITE(numout,*) ' fheat_res : ', fheat_res(jiindex,jjindex)
203            ENDIF
204         END DO
205      END DO
206       
207      !      mass flux at the ocean surface
208      !-------------------------------------------------------
209!     DO jl = 1, jpl
210!        DO jj = 1, jpj
211!           DO ji = 1, jpi
212!              ! this is probably wrong since rdmicif has already been computed
213!              rdmicif(ji,jj) = rdmicif(ji,jj) + rhoic*d_v_i_thd(ji,jj,jl)
214!           END DO
215!        END DO
216!     END DO
217       
218      DO jj = 1, jpj
219         DO ji = 1, jpi
220#if defined key_lim_fdd
221            !  case of realistic freshwater flux (Tartinville et al., 2001) (presently ACTIVATED)
222            !  -------------------------------------------------------------------------------------
223            !  The idea of this approach is that the system that we consider is the ICE-OCEAN system
224            !  Thus  FW  flux  =  External ( E-P+snow melt)
225            !       Salt flux  =  Exchanges in the ice-ocean system then converted into FW
226            !                     Associated to Ice formation AND Ice melting
227            !                     Even if i see Ice melting as a FW and SALT flux
228            !       
229
230            !  computing freshwater exchanges at the ice/ocean interface
231            zpme = - evap(ji,jj) * ( 1.0 - at_i(ji,jj) )     &   !  evaporation over oceanic fraction
232               &   + tprecip(ji,jj)                          &   !  total precipitation
233! old fashioned way               
234!              &   - sprecip(ji,jj) * ( 1. - pfrld(ji,jj) )  &   !  remov. snow precip over ice
235               &   - sprecip(ji,jj) * ( 1. - (pfrld(ji,jj)**betas) )  &   !  remov. snow precip over ice
236               &   - rdmsnif(ji,jj) / rdt_ice                &   !  freshwaterflux due to snow melting
237! new contribution from snow falling when ridging
238               &   + fmmec(ji,jj)
239           
240            !  computing salt exchanges at the ice/ocean interface
241            !  sice should be the same as computed with the ice model
242            zfons =  ( soce - sice ) * ( rdmicif(ji,jj) / rdt_ice ) 
243! SOCE
244            zfons =  ( sss_io(ji,jj) - sice ) * ( rdmicif(ji,jj) / rdt_ice ) 
245           
246            !  salt flux for constant salinity
247            fsalt(ji,jj)      =  zfons / ( sss_io(ji,jj) + epsi16 ) + fsalt_res(ji,jj)
248            zfold             = fsalt(ji,jj)
249           
250            !  salt flux for variable salinity
251            zinda             = 1.0 - MAX( rzero , SIGN( rone , - ( 1.0 - pfrld(ji,jj) ) ) )
252            !correcting brine and salt fluxes
253            fsbri(ji,jj)      =  zinda*fsbri(ji,jj)
254            !  converting the salt fluxes from ice to a freshwater flux from ocean
255            fsalt_res(ji,jj)  =  fsalt_res(ji,jj) / ( sss_io(ji,jj) + epsi16 )
256            fseqv(ji,jj)      =  fseqv(ji,jj)     / ( sss_io(ji,jj) + epsi16 )
257            fsbri(ji,jj)      =  fsbri(ji,jj)     / ( sss_io(ji,jj) + epsi16 )
258            fsalt_rpo(ji,jj)  =  fsalt_rpo(ji,jj) / ( sss_io(ji,jj) + epsi16 )
259
260            !  freshwater mass exchange (positive to the ice, negative for the ocean ?)
261            !  actually it's a salt flux (so it's minus freshwater flux)
262            !  if sea ice grows, zfons is positive, fsalt also
263            !  POSITIVE SALT FLUX FROM THE ICE TO THE OCEAN
264            !  POSITIVE FRESHWATER FLUX FROM THE OCEAN TO THE ICE [kg.m-2.s-1]
265
266            fmass(ji,jj) = - zpme 
267
268#else
269!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
270! ON S'EN TAPE hhhhhhaaaaaaaaaaaaaaaaaahahahahahahahahahahahaha
271!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
272            !  case of freshwater flux equivalent as salt flux
273            !  dilution effect due to evaporation and precipitation
274            zprs  = ( tprecip(ji,jj) - sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) ) * soce 
275!SOCE
276            zprs  = ( tprecip(ji,jj) - sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) ) * sss_io(ji,jj)
277            !  freshwater flux
278            zfons = rdmicif(ji,jj) * ( soce - sice )  &  !  fwf : ice formation and melting
279               &   -  dmgwi(ji,jj) * sice             &  !  fwf : salt flx needed to bring the fresh snow to sea/ice salinity
280               &   + rdmsnif(ji,jj) * soce               !  fwf to ocean due to snow melting
281!SOCE
282            zfons = rdmicif(ji,jj) * ( sss_io(ji,jj) - sice )  &  !  fwf : ice formation and melting
283               &   -  dmgwi(ji,jj) * sice             &  !  fwf : salt flx needed to bring the fresh snow to sea/ice salinity
284               &   + rdmsnif(ji,jj) * sss_io(ji,jj)      !  fwf to ocean due to snow melting
285            !  salt exchanges at the ice/ocean interface
286            zpmess         =  zprs - zfons / rdt_ice - evap(ji,jj) * soce * ( 1.0 - at_i(ji,jj) )
287!SOCE
288            zpmess         =  zprs - zfons / rdt_ice - evap(ji,jj) * sss_io(ji,jj) * ( 1.0 - at_i(ji,jj) )
289
290            fsalt(ji,jj) =  - zpmess
291#endif
292!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
293
294         END DO
295      END DO
296
297      fsalt(:,:) = fseqv(:,:) + fsalt_res(:,:) + fsalt_rpo(:,:)
298      IF (num_sal.eq.2) THEN
299         !In case of variable salinity the salt flux has to be accounted for differently
300         ! Brine drainage has to be added
301         fsalt(:,:) = fsbri(:,:) + fseqv(:,:) + fsalt_res(:,:) + fsalt_rpo(:,:)
302      ENDIF
303
304      !-------------------------------------------------------------------!
305      !  computation of others transmitting variables from ice to ocean   !
306      !------------------------------------------ ------------------------!
307
308      !-----------------------------------------------!
309      !   Storing the transmitted variables           !
310      !-----------------------------------------------!
311
312      ftaux (:,:) = - tio_u(:,:) * rau0   ! taux ( ice: N/m2/rau0, ocean: N/m2 )
313      ftauy (:,:) = - tio_v(:,:) * rau0   ! tauy ( ice: N/m2/rau0, ocean: N/m2 )
314      freeze(:,:) = at_i(:,:)             ! Sea ice cover
315      tn_ice(:,:,:) = t_su(:,:,:)
316
317#if defined key_coupled           
318      zalb  (:,:) = 0.e0
319      zalcn (:,:) = 0.e0
320      zalbp (:,:) = 0.e0
321      zaldum(:,:) = 0.e0
322
323      !------------------------------------------------!
324      !  2) Computation of snow/ice and ocean albedo   !
325      !------------------------------------------------!
326      CALL flx_blk_albedo( zalb, zalcn, zalbp, zaldum )
327
328      alb_ice(:,:) =  0.5 * zalbp(:,:) + 0.5 * zalb (:,:)   ! Ice albedo
329#endif
330
331     IF(ln_ctl) THEN
332        CALL prt_ctl(tab2d_1=fsolar, clinfo1=' lim_flx: fsolar : ', tab2d_2=fnsolar, clinfo2=' fnsolar : ')
333        CALL prt_ctl(tab2d_1=fmass , clinfo1=' lim_flx: fmass  : ', tab2d_2=fsalt  , clinfo2=' fsalt   : ')
334        CALL prt_ctl(tab2d_1=ftaux , clinfo1=' lim_flx: ftaux  : ', tab2d_2=ftauy  , clinfo2=' ftauy   : ')
335        CALL prt_ctl(tab2d_1=freeze, clinfo1=' lim_flx: freeze : ')
336        CALL prt_ctl(tab3d_1=tn_ice, clinfo1=' lim_flx: tn_ice : ', kdim=jpl)
337     ENDIF
338
339    END SUBROUTINE lim_flx
340
341#else
342   !!----------------------------------------------------------------------
343   !!   Default option :        Empty module           NO LIM sea-ice model
344   !!----------------------------------------------------------------------
345CONTAINS
346   SUBROUTINE lim_flx         ! Empty routine
347   END SUBROUTINE lim_flx
348#endif
349
350END MODULE limflx
Note: See TracBrowser for help on using the repository browser.