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.
sbcice_lim_2.F90 in branches/2011/dev_r2802_LOCEAN10_agrif_lim/NEMOGCM/NEMO/OPA_SRC/SBC – NEMO

source: branches/2011/dev_r2802_LOCEAN10_agrif_lim/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim_2.F90 @ 2804

Last change on this file since 2804 was 2804, checked in by rblod, 13 years ago

dev_r2802_LOCEAN10_agrif_lim: first implementation see ticket #848

  • Property svn:keywords set to Id
File size: 13.5 KB
RevLine 
[888]1MODULE sbcice_lim_2
2   !!======================================================================
3   !!                       ***  MODULE  sbcice_lim_2  ***
[2528]4   !! Surface module :  update surface ocean boundary condition over ice covered area using LIM sea-ice model
5   !! Sea-Ice model  :  LIM-2 Sea ice model time-stepping
[888]6   !!======================================================================
[1218]7   !! History :  1.0   !  06-2006  (G. Madec)  from icestp_2.F90
8   !!            3.0   !  08-2008  (S. Masson, E. .... ) coupled interface
[2528]9   !!            3.3   !  05-2009  (G.Garric) addition of the lim2_evp case
[888]10   !!----------------------------------------------------------------------
11#if defined key_lim2
12   !!----------------------------------------------------------------------
[2528]13   !!   'key_lim2' :                                    LIM-2 sea-ice model
[888]14   !!----------------------------------------------------------------------
[2528]15   !!   sbc_ice_lim_2   : sea-ice model time-stepping and update ocean sbc over ice-covered area
[888]16   !!----------------------------------------------------------------------
[2528]17   USE oce              ! ocean dynamics and tracers
18   USE dom_oce          ! ocean space and time domain
[888]19   USE ice_2
[1270]20   USE par_ice_2
[888]21   USE iceini_2
22   USE dom_ice_2
23
[2528]24   USE sbc_oce          ! Surface boundary condition: ocean fields
25   USE sbc_ice          ! Surface boundary condition: ice   fields
26   USE sbcblk_core      ! Surface boundary condition: CORE bulk
27   USE sbcblk_clio      ! Surface boundary condition: CLIO bulk
28   USE sbccpl           ! Surface boundary condition: coupled interface
[888]29   USE albedo
30
[2528]31   USE phycst           ! Define parameters for the routines
32   USE eosbn2           ! equation of state
[888]33   USE limdyn_2
34   USE limtrp_2
35   USE limdmp_2
36   USE limthd_2
[2528]37   USE limsbc_2         ! sea surface boundary condition
[888]38   USE limdia_2
39   USE limwri_2
40   USE limrst_2
41
[2528]42   USE c1d              ! 1D vertical configuration
[888]43
[2528]44   USE lbclnk           ! lateral boundary condition - MPP link
45   USE lib_mpp          ! MPP library
46   USE iom              ! I/O manager library
47   USE in_out_manager   ! I/O manager
48   USE prtctl           ! Print control
49
[2804]50#if defined key_agrif
51   USE agrif_ice
52   USE agrif_lim2_update
53#endif
54
[888]55   IMPLICIT NONE
56   PRIVATE
57
58   PUBLIC sbc_ice_lim_2 ! routine called by sbcmod.F90
59
60   !! * Substitutions
61#  include "domzgr_substitute.h90"
62#  include "vectopt_loop_substitute.h90"
63   !!----------------------------------------------------------------------
[2528]64   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
[1146]65   !! $Id$
[2528]66   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
[888]67   !!----------------------------------------------------------------------
68CONTAINS
69
[1218]70   SUBROUTINE sbc_ice_lim_2( kt, ksbc )
[888]71      !!---------------------------------------------------------------------
72      !!                  ***  ROUTINE sbc_ice_lim_2  ***
73      !!                   
74      !! ** Purpose :   update the ocean surface boundary condition via the
75      !!                Louvain la Neuve Sea Ice Model time stepping
76      !!
77      !! ** Method  :   ice model time stepping
78      !!              - call the ice dynamics routine
79      !!              - call the ice advection/diffusion routine
80      !!              - call the ice thermodynamics routine
81      !!              - call the routine that computes mass and
82      !!                heat fluxes at the ice/ocean interface
83      !!              - save the outputs
84      !!              - save the outputs for restart when necessary
85      !!
86      !! ** Action  : - time evolution of the LIM sea-ice model
87      !!              - update all sbc variables below sea-ice:
[1695]88      !!                utau, vtau, taum, wndm, qns , qsr, emp , emps
[888]89      !!---------------------------------------------------------------------
[2715]90      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released
91      USE wrk_nemo, ONLY:   wrk_3d_1 , wrk_3d_2 , wrk_3d_3   ! 3D workspace
92      !!
[888]93      INTEGER, INTENT(in) ::   kt      ! ocean time step
[1218]94      INTEGER, INTENT(in) ::   ksbc    ! type of sbc ( =3 CLIO bulk ; =4 CORE bulk ; =5 coupled )
[888]95      !!
96      INTEGER  ::   ji, jj   ! dummy loop indices
[2715]97      ! Pointers into workspaces contained in the wrk_nemo module
98      REAL(wp), DIMENSION(:,:,:), POINTER :: zalb_ice_os   ! albedo of the ice under overcast sky
99      REAL(wp), DIMENSION(:,:,:), POINTER :: zalb_ice_cs   ! albedo of ice under clear sky
100      REAL(wp), DIMENSION(:,:,:), POINTER :: zsist         ! surface ice temperature (K)
[888]101      !!----------------------------------------------------------------------
102
[2715]103      IF( wrk_in_use(3, 1,2,3) ) THEN
104         CALL ctl_stop('sbc_ice_lim_2: requested workspace arrays are unavailable')   ;   RETURN
105      ENDIF
106      ! Use pointers to access only sub-arrays of workspaces
107      zalb_ice_os => wrk_3d_1(:,:,1:1)
108      zalb_ice_cs => wrk_3d_2(:,:,1:1)
109      zsist       => wrk_3d_3(:,:,1:1)
110
[888]111      IF( kt == nit000 ) THEN
112         IF(lwp) WRITE(numout,*)
113         IF(lwp) WRITE(numout,*) 'sbc_ice_lim_2 : update ocean surface boudary condition' 
114         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   via Louvain la Neuve Ice Model (LIM) time stepping'
[2528]115         !
[888]116         CALL ice_init_2
[2804]117#if defined key_agrif
118         IF( .NOT. Agrif_Root() ) CALL Agrif_InitValues_cont_lim2  ! AGRIF: set the meshes
119# endif
[888]120      ENDIF
121
[2528]122      !                                        !----------------------!
123      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !  Ice time-step only  !
124         !                                     !----------------------!
[2804]125#if defined key_agrif
126         IF( .NOT. Agrif_Root() ) childfreq = MOD(childfreq,Agrif_rhot()) + 1
127#endif
128         
[2528]129         !  Bulk Formulea !
130         !----------------!
[888]131         ! ... mean surface ocean current at ice dynamics point
[2528]132         SELECT CASE( cp_ice_msh )
133         CASE( 'I' )                  !== B-grid ice dynamics :   I-point (i.e. F-point with sea-ice indexation)
134            DO jj = 2, jpj
135               DO ji = 2, jpi   ! NO vector opt. possible
136                  u_oce(ji,jj) = 0.5_wp * ( ssu_m(ji-1,jj  ) + ssu_m(ji-1,jj-1) ) * tmu(ji,jj)
137                  v_oce(ji,jj) = 0.5_wp * ( ssv_m(ji  ,jj-1) + ssv_m(ji-1,jj-1) ) * tmu(ji,jj)
138               END DO
[888]139            END DO
[2528]140            CALL lbc_lnk( u_oce, 'I', -1. )   ! I-point (i.e. F-point with ice indices)
141            CALL lbc_lnk( v_oce, 'I', -1. )   ! I-point (i.e. F-point with ice indices)
142            !
143         CASE( 'C' )                  !== C-grid ice dynamics :   U & V-points (same as ocean)
144            u_oce(:,:) = ssu_m(:,:)                     ! mean surface ocean current at ice velocity point
145            v_oce(:,:) = ssv_m(:,:)
146            !
147         END SELECT
[888]148
149         ! ... masked sea surface freezing temperature [Kelvin] (set to rt0 over land)
150         tfu(:,:) = tfreez( sss_m ) +  rt0 
151
[1109]152         zsist (:,:,1) = sist (:,:) + rt0 * ( 1. - tmask(:,:,1) )
[888]153
[1218]154         ! ... ice albedo (clear sky and overcast sky)
[2715]155         CALL albedo_ice( zsist, reshape( hicif, (/jpi,jpj,1/) ), &
156                                 reshape( hsnif, (/jpi,jpj,1/) ), &
157                          zalb_ice_cs, zalb_ice_os )
[888]158
159         ! ... Sea-ice surface boundary conditions output from bulk formulae :
[1469]160         !     - utau_ice   ! surface ice stress i-component (I-point)   [N/m2]
161         !     - vtau_ice   ! surface ice stress j-component (I-point)   [N/m2]
[888]162         !     - qns_ice    ! non solar heat flux over ice   (T-point)   [W/m2]
163         !     - qsr_ice    !     solar heat flux over ice   (T-point)   [W/m2]
164         !     - qla_ice    ! latent    heat flux over ice   (T-point)   [W/m2]
165         !     - dqns_ice   ! non solar heat sensistivity    (T-point)   [W/m2]
166         !     - dqla_ice   ! latent    heat sensistivity    (T-point)   [W/m2]
167         !     - tprecip    ! total precipitation            (T-point)   [Kg/m2/s]
168         !     - sprecip    ! solid precipitation            (T-point)   [Kg/m2/s]
169         !     - fr1_i0     ! 1sr fraction of qsr penetration in ice     [%]
170         !     - fr2_i0     ! 2nd fraction of qsr penetration in ice     [%]
171         !
[1218]172         SELECT CASE( ksbc )
[888]173         CASE( 3 )           ! CLIO bulk formulation
[1469]174            CALL blk_ice_clio( zsist, zalb_ice_cs, zalb_ice_os,                         &
175               &                      utau_ice   , vtau_ice   , qns_ice    , qsr_ice,   &
176               &                      qla_ice    , dqns_ice   , dqla_ice   ,            &
177               &                      tprecip    , sprecip    ,                         &
[2528]178               &                      fr1_i0     , fr2_i0     , cp_ice_msh , jpl  )
[914]179
[888]180         CASE( 4 )           ! CORE bulk formulation
[1470]181            CALL blk_ice_core( zsist, u_ice      , v_ice      , zalb_ice_cs,            &
[1469]182               &                      utau_ice   , vtau_ice   , qns_ice    , qsr_ice,   &
183               &                      qla_ice    , dqns_ice   , dqla_ice   ,            &
184               &                      tprecip    , sprecip    ,                         &
[2528]185               &                      fr1_i0     , fr2_i0     , cp_ice_msh , jpl  )
[1218]186         CASE( 5 )           ! Coupled formulation : atmosphere-ice stress only (fluxes provided after ice dynamics)
[1469]187            CALL sbc_cpl_ice_tau( utau_ice , vtau_ice )
[888]188         END SELECT
189
[1482]190         CALL iom_put( 'utau_ice', utau_ice )     ! Wind stress over ice along i-axis at I-point
191         CALL iom_put( 'vtau_ice', vtau_ice )     ! Wind stress over ice along j-axis at I-point
192
[888]193         IF(ln_ctl) THEN         ! print mean trends (used for debugging)
194            CALL prt_ctl_info( 'Ice Forcings ' )
[1469]195            CALL prt_ctl( tab2d_1=tprecip ,clinfo1=' sbc_ice_lim: precip  : ', tab2d_2=sprecip , clinfo2=' Snow    : ' )
196            CALL prt_ctl( tab2d_1=utau_ice,clinfo1=' sbc_ice_lim: utau_ice: ', tab2d_2=vtau_ice, clinfo2=' vtau_ice: ' )
197            CALL prt_ctl( tab2d_1=sst_m   ,clinfo1=' sbc_ice_lim: sst     : ', tab2d_2=sss_m   , clinfo2=' sss     : ' )
[1470]198            CALL prt_ctl( tab2d_1=u_oce   ,clinfo1=' sbc_ice_lim: u_io    : ', tab2d_2=v_oce   , clinfo2=' v_io    : ' )
[1469]199            CALL prt_ctl( tab2d_1=hsnif   ,clinfo1=' sbc_ice_lim: hsnif  1: ', tab2d_2=hicif   , clinfo2=' hicif   : ' )
200            CALL prt_ctl( tab2d_1=frld    ,clinfo1=' sbc_ice_lim: frld   1: ', tab2d_2=sist    , clinfo2=' sist    : ' )
[888]201         ENDIF
202
203         ! ---------------- !
204         !  Ice model step  !
205         ! ---------------- !
[2528]206         numit = numit + nn_fsbc                           ! Ice model time step
[1239]207
[2528]208                           CALL lim_rst_opn_2  ( kt )  ! Open Ice restart file
209         IF( .NOT. lk_c1d ) THEN                       ! Ice dynamics & transport (except in 1D case)
210                           CALL lim_dyn_2      ( kt )      ! Ice dynamics    ( rheology/dynamics )
[2804]211
[2528]212                           CALL lim_trp_2      ( kt )      ! Ice transport   ( Advection/diffusion )
213           IF( ln_limdmp ) CALL lim_dmp_2      ( kt )      ! Ice damping
214         END IF
[1218]215#if defined key_coupled
[2528]216         !                                             ! Ice surface fluxes in coupled mode
217         IF( ksbc == 5 )   CALL sbc_cpl_ice_flx( reshape( frld, (/jpi,jpj,1/) ),                 &
[2715]218            &                                             qns_tot, qns_ice, qsr_tot , qsr_ice,   &
219            &                                             emp_tot, emp_ice, dqns_ice, sprecip,   &
220            !                                optional arguments, used only in 'mixed oce-ice' case
221            &                                             palbi = zalb_ice_cs, psst = sst_m, pist = zsist )
[1218]222#endif
[2528]223                           CALL lim_thd_2      ( kt )      ! Ice thermodynamics
224                           CALL lim_sbc_flx_2  ( kt )      ! update surface ocean mass, heat & salt fluxes
[1481]225
[1055]226         IF( ( MOD( kt+nn_fsbc-1, ninfo ) == 0 .OR. ntmoy == 1 ) .AND. .NOT. lk_mpp )   &
[2528]227            &              CALL lim_dia_2      ( kt )      ! Ice Diagnostics
[1482]228# if ! defined key_iomput
[2528]229                           CALL lim_wri_2      ( kt )      ! Ice outputs
[1482]230# endif
[2528]231         IF( lrst_ice  )   CALL lim_rst_write_2( kt )      ! Ice restart file
[888]232         !
[2528]233      ENDIF                                    ! End sea-ice time step only
[888]234      !
[2528]235      !                                        !--------------------------!
236      !                                        !  at all ocean time step  !
237      !                                        !--------------------------!
[2804]238#if defined key_agrif && defined key_lim2
239      IF( .NOT. Agrif_Root() )   CALL agrif_update_lim2( kt )
240#endif
[2528]241      !                                               
242      !                                              ! Update surface ocean stresses (only in ice-dynamic case)
243      !                                                   ! otherwise the atm.-ocean stresses are used everywhere
244      IF( ln_limdyn    )   CALL lim_sbc_tau_2( kt, ub(:,:,1), vb(:,:,1) )  ! using before instantaneous surf. currents
245      !
[2715]246      IF( wrk_not_released(3, 1,2,3) )   CALL ctl_stop('sbc_ice_lim_2: failed to release workspace arrays')
247      !
[888]248   END SUBROUTINE sbc_ice_lim_2
249
250#else
251   !!----------------------------------------------------------------------
252   !!   Default option           Dummy module      NO LIM 2.0 sea-ice model
253   !!----------------------------------------------------------------------
254CONTAINS
[1218]255   SUBROUTINE sbc_ice_lim_2 ( kt, ksbc )     ! Dummy routine
[2715]256      INTEGER, INTENT(in) ::   kt, ksbc   
[1218]257      WRITE(*,*) 'sbc_ice_lim_2: You should not have seen this print! error?', kt, ksbc
[888]258   END SUBROUTINE sbc_ice_lim_2
259#endif
260
261   !!======================================================================
262END MODULE sbcice_lim_2
Note: See TracBrowser for help on using the repository browser.