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.
icedyn_adv.F90 in NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap/src/ICE – NEMO

source: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps_rewrite_time_filterswap/src/ICE/icedyn_adv.F90 @ 11021

Last change on this file since 11021 was 10413, checked in by clem, 5 years ago

merge dev_r9947_SI3_advection with the trunk. All sette tests passed. There is probably a conservation issue with the new advection scheme that I should solve but the routine icedyn_adv_umx.F90 is going to change anyway in the next couple of weeks. At worst, the old routine can be plugged without harm
  • Property svn:keywords set to Id
File size: 10.2 KB
Line 
1MODULE icedyn_adv
2   !!======================================================================
3   !!                       ***  MODULE icedyn_adv   ***
4   !!   sea-ice: advection
5   !!======================================================================
6   !! History :  4.0  !  2018     (many people)      SI3 [aka Sea Ice cube]
7   !!----------------------------------------------------------------------
8#if defined key_si3
9   !!----------------------------------------------------------------------
10   !!   'key_si3'                                       SI3 sea-ice model
11   !!----------------------------------------------------------------------
12   !!   ice_dyn_adv   : advection of sea ice variables
13   !!----------------------------------------------------------------------
14   USE phycst         ! physical constant
15   USE dom_oce        ! ocean domain
16   USE sbc_oce , ONLY : nn_fsbc   ! frequency of sea-ice call
17   USE ice            ! sea-ice: variables
18   USE icevar         ! sea-ice: operations
19   USE icedyn_adv_pra ! sea-ice: advection scheme (Prather)
20   USE icedyn_adv_umx ! sea-ice: advection scheme (ultimate-macho)
21   USE icectl         ! sea-ice: control prints
22   !
23   USE in_out_manager ! I/O manager
24   USE iom            ! I/O manager library
25   USE lib_mpp        ! MPP library
26   USE lib_fortran    ! fortran utilities (glob_sum + no signed zero)
27   USE timing         ! Timing
28   USE prtctl         ! Print control
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC   ice_dyn_adv        ! called by icestp
34   PUBLIC   ice_dyn_adv_init   ! called by icedyn
35
36   INTEGER ::              nice_adv   ! choice of the type of advection scheme
37   !                                        ! associated indices:
38   INTEGER, PARAMETER ::   np_advPRA = 1   ! Prather scheme
39   INTEGER, PARAMETER ::   np_advUMx = 2   ! Ultimate-Macho scheme
40   !
41   ! ** namelist (namdyn_adv) **
42   INTEGER         ::   nn_UMx       ! order of the UMx advection scheme   
43   !
44   !! * Substitution
45#  include "vectopt_loop_substitute.h90"
46   !!----------------------------------------------------------------------
47   !! NEMO/ICE 4.0 , NEMO Consortium (2018)
48   !! $Id$
49   !! Software governed by the CeCILL licence     (./LICENSE)
50   !!----------------------------------------------------------------------
51CONTAINS
52
53   SUBROUTINE ice_dyn_adv( kt ) 
54      !!----------------------------------------------------------------------
55      !!                   ***  ROUTINE ice_dyn_adv ***
56      !!                   
57      !! ** purpose : advection of sea ice
58      !!
59      !! ** method  : One can choose between
60      !!     a) an Ultimate-Macho scheme (with order defined by nn_UMx) => ln_adv_UMx
61      !!     b) and a second order Prather scheme => ln_adv_Pra
62      !!
63      !! ** action :
64      !!----------------------------------------------------------------------
65      INTEGER, INTENT(in) ::   kt   ! number of iteration
66      !
67      INTEGER ::   jl   ! dummy loop indice
68      REAL(wp), DIMENSION(jpi,jpj) ::   zmask  ! fraction of time step with v_i < 0
69      !!---------------------------------------------------------------------
70      !
71      IF( ln_timing )   CALL timing_start('icedyn_adv')
72      !
73      IF( kt == nit000 .AND. lwp ) THEN
74         WRITE(numout,*)
75         WRITE(numout,*) 'ice_dyn_adv: sea-ice advection'
76         WRITE(numout,*) '~~~~~~~~~~~'
77      ENDIF
78     
79      ! conservation test
80      IF( ln_icediachk )   CALL ice_cons_hsm(0, 'icedyn_adv', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft)
81                     
82      !----------
83      ! Advection
84      !----------
85      SELECT CASE( nice_adv )
86      !                                !-----------------------!
87      CASE( np_advUMx )                ! ULTIMATE-MACHO scheme !
88         !                             !-----------------------!
89         CALL ice_dyn_adv_umx( nn_UMx, kt, u_ice, v_ice, ato_i, v_i, v_s, sv_i, oa_i, a_i, a_ip, v_ip, e_s, e_i )
90      !                                !-----------------------!
91      CASE( np_advPRA )                ! PRATHER scheme        !
92         !                             !-----------------------!
93         CALL ice_dyn_adv_pra(         kt, u_ice, v_ice, ato_i, v_i, v_s, sv_i, oa_i, a_i, a_ip, v_ip, e_s, e_i )
94      END SELECT
95
96      !----------------------------
97      ! Debug the advection schemes
98      !----------------------------
99      ! clem: At least one advection scheme above is not strictly positive => UMx
100      !       In Prather, I am not sure if the fields are bounded by 0 or not (it seems yes)
101      !       In UMx    , advected fields are not perfectly bounded and negative values can appear.
102      !                   These values are usually very small but in some occasions they can also be non-negligible
103      !                   Therefore one needs to bound the advected fields by 0 (though this is not a clean fix)
104      !
105      ! record the negative values resulting from UMx
106      zmask(:,:) = 0._wp ! keep the init to 0 here
107      DO jl = 1, jpl
108         WHERE( v_i(:,:,jl) < 0._wp )   zmask(:,:) = 1._wp
109      END DO
110      IF( iom_use('iceneg_pres') )   CALL iom_put("iceneg_pres", zmask                                      )  ! fraction of time step with v_i < 0
111      IF( iom_use('iceneg_volu') )   CALL iom_put("iceneg_volu", SUM(MIN( v_i, 0. ), dim=3 )                )  ! negative ice volume (only)
112      IF( iom_use('iceneg_hfx' ) )   CALL iom_put("iceneg_hfx" , ( SUM(SUM( MIN( e_i(:,:,1:nlay_i,:), 0. )  &  ! negative ice heat content (only)
113         &                                                                  , dim=4 ), dim=3 ) )* r1_rdtice )  ! -- eq. heat flux [W/m2]
114      !
115      ! ==> conservation is ensured by calling this routine below,
116      !     however the global ice volume is then changed by advection (but errors are small)
117      CALL ice_var_zapneg( ato_i, v_i, v_s, sv_i, oa_i, a_i, a_ip, v_ip, e_s, e_i )
118
119      !------------
120      ! diagnostics
121      !------------
122      diag_trp_ei(:,:) = SUM(SUM( e_i (:,:,1:nlay_i,:) - e_i_b (:,:,1:nlay_i,:), dim=4 ), dim=3 ) * r1_rdtice
123      diag_trp_es(:,:) = SUM(SUM( e_s (:,:,1:nlay_s,:) - e_s_b (:,:,1:nlay_s,:), dim=4 ), dim=3 ) * r1_rdtice
124      diag_trp_sv(:,:) = SUM(     sv_i(:,:,:)          - sv_i_b(:,:,:)                  , dim=3 ) * r1_rdtice
125      diag_trp_vi(:,:) = SUM(     v_i (:,:,:)          - v_i_b (:,:,:)                  , dim=3 ) * r1_rdtice
126      diag_trp_vs(:,:) = SUM(     v_s (:,:,:)          - v_s_b (:,:,:)                  , dim=3 ) * r1_rdtice
127      IF( iom_use('icemtrp') )   CALL iom_put( "icemtrp" , diag_trp_vi * rhoi          )   ! ice mass transport
128      IF( iom_use('snwmtrp') )   CALL iom_put( "snwmtrp" , diag_trp_vs * rhos          )   ! snw mass transport
129      IF( iom_use('salmtrp') )   CALL iom_put( "salmtrp" , diag_trp_sv * rhoi * 1.e-03 )   ! salt mass transport (kg/m2/s)
130      IF( iom_use('dihctrp') )   CALL iom_put( "dihctrp" , -diag_trp_ei                )   ! advected ice heat content (W/m2)
131      IF( iom_use('dshctrp') )   CALL iom_put( "dshctrp" , -diag_trp_es                )   ! advected snw heat content (W/m2)
132
133      ! controls
134      IF( ln_icediachk )   CALL ice_cons_hsm(1, 'icedyn_adv', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) ! conservation
135      IF( ln_icectl    )   CALL ice_prt     (kt, iiceprt, jiceprt,-1, ' - ice dyn & trp - ')                           ! prints
136      IF( ln_timing    )   CALL timing_stop ('icedyn_adv')                                                             ! timing
137      !
138   END SUBROUTINE ice_dyn_adv
139
140
141   SUBROUTINE ice_dyn_adv_init
142      !!-------------------------------------------------------------------
143      !!                  ***  ROUTINE ice_dyn_adv_init  ***
144      !!
145      !! ** Purpose :   Physical constants and parameters linked to the ice
146      !!                dynamics
147      !!
148      !! ** Method  :   Read the namdyn_adv namelist and check the ice-dynamic
149      !!                parameter values called at the first timestep (nit000)
150      !!
151      !! ** input   :   Namelist namdyn_adv
152      !!-------------------------------------------------------------------
153      INTEGER ::   ios, ioptio   ! Local integer output status for namelist read
154      !!
155      NAMELIST/namdyn_adv/ ln_adv_Pra, ln_adv_UMx, nn_UMx
156      !!-------------------------------------------------------------------
157      !
158      REWIND( numnam_ice_ref )         ! Namelist namdyn_adv in reference namelist : Ice dynamics
159      READ  ( numnam_ice_ref, namdyn_adv, IOSTAT = ios, ERR = 901)
160901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namdyn_adv in reference namelist', lwp )
161      REWIND( numnam_ice_cfg )         ! Namelist namdyn_adv in configuration namelist : Ice dynamics
162      READ  ( numnam_ice_cfg, namdyn_adv, IOSTAT = ios, ERR = 902 )
163902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namdyn_adv in configuration namelist', lwp )
164      IF(lwm) WRITE( numoni, namdyn_adv )
165      !
166      IF(lwp) THEN                     ! control print
167         WRITE(numout,*)
168         WRITE(numout,*) 'ice_dyn_adv_init: ice parameters for ice dynamics '
169         WRITE(numout,*) '~~~~~~~~~~~~~~~~'
170         WRITE(numout,*) '   Namelist namdyn_adv:'
171         WRITE(numout,*) '      type of advection scheme (Prather)             ln_adv_Pra = ', ln_adv_Pra 
172         WRITE(numout,*) '      type of advection scheme (Ulimate-Macho)       ln_adv_UMx = ', ln_adv_UMx 
173         WRITE(numout,*) '         order of the Ultimate-Macho scheme          nn_UMx     = ', nn_UMx
174      ENDIF
175      !
176      !                             !== set the choice of ice advection ==!
177      ioptio = 0 
178      IF( ln_adv_Pra ) THEN   ;   ioptio = ioptio + 1   ;   nice_adv = np_advPRA    ;   ENDIF
179      IF( ln_adv_UMx ) THEN   ;   ioptio = ioptio + 1   ;   nice_adv = np_advUMx    ;   ENDIF
180      IF( ioptio /= 1 )   CALL ctl_stop( 'ice_dyn_adv_init: choose one and only one ice adv. scheme (ln_adv_Pra or ln_adv_UMx)' )
181      !
182      IF( ln_adv_Pra )   CALL adv_pra_init  !* read or initialize all required files
183      !
184   END SUBROUTINE ice_dyn_adv_init
185
186#else
187   !!----------------------------------------------------------------------
188   !!   Default option         Empty Module           NO SI3 sea-ice model
189   !!----------------------------------------------------------------------
190#endif
191
192   !!======================================================================
193END MODULE icedyn_adv
194
Note: See TracBrowser for help on using the repository browser.