source: branches/dev_r4028_CNRS_LIM3_MV2014/NEMOGCM/NEMO/LIM_SRC_3/limdiahsb.F90 @ 4626

Last change on this file since 4626 was 4626, checked in by gm, 7 years ago

dev_r4028_CNRS_LIM3_MV2014 : minor corrections on LIM3

File size: 13.3 KB
Line 
1MODULE limdiahsb
2   !!======================================================================
3   !!                       ***  MODULE limdia_hsb   ***
4   !!  LIM-3 sea ice model :   diagnostics of ice model
5   !!======================================================================
6   !! History :  3.4  ! 2012-10  (C. Rousset)  original code
7   !!----------------------------------------------------------------------
8#if defined key_lim3
9   !!----------------------------------------------------------------------
10   !!   'key_lim3'                                       LIM3 sea-ice model
11   !!----------------------------------------------------------------------
12   !!   lim_dia_hsb        : computation and output of the time evolution of keys variables
13   !!   lim_dia_hsb_init   : initialization and namelist read
14   !!----------------------------------------------------------------------
15   USE ice             ! LIM-3: sea-ice variable
16   USE par_ice         ! LIM-3: ice parameters
17   USE dom_ice         ! LIM-3: sea-ice domain
18   USE dom_oce         ! ocean domain
19   USE sbc_oce         ! surface boundary condition: ocean fields
20   USE daymod          ! model calendar
21   USE phycst          ! physical constant
22   USE in_out_manager  ! I/O manager
23   USE lib_mpp         ! MPP library
24   USE timing          ! preformance summary
25   USE iom             ! I/O manager
26   USE lib_fortran     ! glob_sum
27   USE limrst          ! ice restart
28
29   IMPLICIT NONE
30   PRIVATE
31
32   PUBLIC   lim_diahsb        ! routine called by ice_step.F90
33   !!PUBLIC   lim_diahsb_init   ! routine called by ice_init.F90
34   !!PUBLIC   lim_diahsb_rst   ! routine called by ice_init.F90
35
36   REAL(wp) ::   frc_sal, frc_vol   ! global forcing trends
37   REAL(wp) ::   bg_grme            ! global ice growth+melt trends
38   REAL(wp) ::   epsi06 = 1.e-6_wp  ! small number
39   REAL(wp) ::   epsi03 = 1.e-3_wp  ! small number
40
41
42   !! * Substitutions
43#  include "vectopt_loop_substitute.h90"
44
45   !!----------------------------------------------------------------------
46   !! NEMO/OPA 3.6 , NEMO Consortium (2014)
47   !! $Id: limdiahsb.F90 3294 2012-10-18 16:44:18Z rblod $
48   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
49   !!----------------------------------------------------------------------
50CONTAINS
51
52   SUBROUTINE lim_diahsb
53      !!---------------------------------------------------------------------------
54      !!                  ***  ROUTINE lim_diahsb  ***
55      !!     
56      !! ** Purpose: Compute the ice global heat content, salt content and volume conservation
57      !!
58      !!---------------------------------------------------------------------------
59      REAL(wp)   ::   zbg_ivo, zbg_svo, zbg_are, zbg_sal ,zbg_tem ,zbg_ihc ,zbg_shc
60      REAL(wp)   ::   zbg_sfx, zbg_sfx_bri, zbg_sfx_thd, zbg_sfx_res, zbg_sfx_mec 
61      REAL(wp)   ::   zbg_emp, zbg_emp_bog, zbg_emp_lag, zbg_emp_sig, zbg_emp_dyg, zbg_emp_bom, zbg_emp_sum, zbg_emp_res 
62      REAL(wp)   ::   z_frc_vol, z_frc_sal, z_bg_grme 
63      REAL(wp)   ::   z1_area, zcoef                   
64      REAL(wp)   ::   zinda, zindb
65      !!---------------------------------------------------------------------------
66      IF( nn_timing == 1 )   CALL timing_start('lim_diahsb')
67
68      IF( numit == nstart ) CALL lim_diahsb_init 
69
70      z1_area = 1._wp / MAX( glob_sum( area(:,:) * tms(:,:) ), epsi06 )
71
72      zinda = MAX( 0._wp , SIGN( 1._wp , glob_sum( area(:,:) * tms(:,:) ) - epsi06 ) )
73      ! ----------------------- !
74      ! 1 -  Content variations !
75      ! ----------------------- !
76      zbg_ivo = glob_sum( vt_i(:,:) * area(:,:) * tms(:,:) ) ! volume ice
77      zbg_svo = glob_sum( vt_s(:,:) * area(:,:) * tms(:,:) ) ! volume snow
78      zbg_are = glob_sum( at_i(:,:) * area(:,:) * tms(:,:) ) ! area
79      zbg_sal = glob_sum( SUM( smv_i(:,:,:), dim=3 )      * area(:,:) * tms(:,:) )  ! mean salt content
80      zbg_tem = glob_sum( ( tm_i(:,:) - rtt ) * vt_i(:,:) * area(:,:) * tms(:,:) )  ! mean temp content
81
82      !zbg_ihc = glob_sum( et_i(:,:) * area(:,:) * tms(:,:) ) / MAX( zbg_ivo,epsi06 ) ! ice heat content
83      !zbg_shc = glob_sum( et_s(:,:) * area(:,:) * tms(:,:) ) / MAX( zbg_svo,epsi06 ) ! snow heat content
84
85      zbg_ihc = glob_sum( et_i(:,:) * 1.e-11 ) ! ice heat content  [10^9*1.e-11 J]
86      zbg_shc = glob_sum( et_s(:,:) * 1.e-11 ) ! snow heat content [10^9*1.e-11 J]
87
88      zcoef = zinda * z1_area * r1_rau0 * rday
89      zbg_emp     =  glob_sum(         emp(:,:) * area(:,:) * tms(:,:) ) * zcoef
90     
91      zcoef = zinda * z1_area * rhoic * r1_rau0 * rday
92      zbg_emp_bog = glob_sum( diag_bot_gr(:,:) * area(:,:) * tms(:,:) ) * zcoef
93      zbg_emp_lag = glob_sum( diag_lat_gr(:,:) * area(:,:) * tms(:,:) ) * zcoef
94      zbg_emp_sig = glob_sum( diag_sni_gr(:,:) * area(:,:) * tms(:,:) ) * zcoef
95      zbg_emp_dyg = glob_sum( diag_dyn_gr(:,:) * area(:,:) * tms(:,:) ) * zcoef
96      zbg_emp_bom = glob_sum( diag_bot_me(:,:) * area(:,:) * tms(:,:) ) * zcoef
97      zbg_emp_sum = glob_sum( diag_sur_me(:,:) * area(:,:) * tms(:,:) ) * zcoef
98      zbg_emp_res = glob_sum( diag_res_pr(:,:) * area(:,:) * tms(:,:) ) * zcoef
99
100      zcoef = zinda * z1_area * r1_rau0 * rday
101      zbg_sfx     = glob_sum(     sfx(:,:) * area(:,:) * tms(:,:) ) * zcoef
102      zbg_sfx_bri = glob_sum( sfx_bri(:,:) * area(:,:) * tms(:,:) ) * zcoef
103      zbg_sfx_thd = glob_sum( sfx_thd(:,:) * area(:,:) * tms(:,:) ) * zcoef
104      zbg_sfx_res = glob_sum( sfx_res(:,:) * area(:,:) * tms(:,:) ) * zcoef
105      zbg_sfx_mec = glob_sum( sfx_mec(:,:) * area(:,:) * tms(:,:) ) * zcoef
106     
107      ! --------------------------------------------- !
108      ! 2 - Trends due to forcing and ice growth/melt !
109      ! --------------------------------------------- !
110      z_frc_vol = r1_rau0 * glob_sum( - emp(:,:) * area(:,:) * tms(:,:) ) ! volume fluxes
111      z_frc_sal = r1_rau0 * glob_sum(   sfx(:,:) * area(:,:) * tms(:,:) ) ! salt fluxes
112      z_bg_grme = glob_sum(  (  diag_bot_gr(:,:) + diag_lat_gr(:,:) + diag_sni_gr(:,:) + diag_dyn_gr(:,:)     &
113         &                    + diag_bot_me(:,:) + diag_sur_me(:,:) + diag_res_pr(:,:)                    )   &
114         &                 * area(:,:) * tms(:,:) ) ! volume fluxes
115      !
116      frc_vol  = frc_vol + z_frc_vol * rdt_ice
117      frc_sal  = frc_sal + z_frc_sal * rdt_ice
118      bg_grme  = bg_grme + z_bg_grme * rdt_ice
119     
120      ! difference
121      !frc_vol = zbg_ivo - frc_vol
122      !frc_sal = zbg_sal - frc_sal
123     
124      ! ----------------------- !
125      ! 3 - Diagnostics writing !
126      ! ----------------------- !
127      zindb = MAX( 0.d0 , SIGN( 1._wp , zbg_ivo - epsi06 ) )
128      !
129      CALL iom_put( 'ibgvoltot' , zbg_ivo * rhoic * r1_rau0 * 1.e-9        )   ! ice volume (km3 equivalent liquid)         
130      CALL iom_put( 'sbgvoltot' , zbg_svo * rhosn * r1_rau0 * 1.e-9        )   ! snw volume (km3 equivalent liquid)       
131      CALL iom_put( 'ibgarea'   , zbg_are * 1.e-6                          )   ! ice area   (km2)
132      CALL iom_put( 'ibgsaline' , zindb * zbg_sal / MAX( zbg_ivo, epsi06 ) )   ! ice saline (psu)
133      CALL iom_put( 'ibgtemper' , zindb * zbg_tem / MAX( zbg_ivo, epsi06 ) )   ! ice temper (C)
134      CALL iom_put( 'ibgheatco' , zbg_ihc                                  )   ! ice heat content (1.e20 J)       
135      CALL iom_put( 'sbgheatco' , zbg_shc                                  )   ! snw heat content (1.e20 J)
136      CALL iom_put( 'ibgsaltco' , zbg_sal * rhoic * r1_rau0 * 1.e-9        )   ! ice salt content (psu*km3 equivalent liquid)       
137
138      CALL iom_put( 'ibgemp'    , zbg_emp                                  )   ! volume flux emp (m/day liquid)
139      CALL iom_put( 'ibgempbog' , zbg_emp_bog                              )   ! volume flux bottom growth     -(m/day equivalent liquid)
140      CALL iom_put( 'ibgemplag' , zbg_emp_lag                              )   ! volume flux open water growth -
141      CALL iom_put( 'ibgempsig' , zbg_emp_sig                              )   ! volume flux snow ice growth   -
142      CALL iom_put( 'ibgempdyg' , zbg_emp_dyg                              )   ! volume flux dynamic growth    -
143      CALL iom_put( 'ibgempbom' , zbg_emp_bom                              )   ! volume flux bottom melt       -
144      CALL iom_put( 'ibgempsum' , zbg_emp_sum                              )   ! volume flux surface melt      -
145      CALL iom_put( 'ibgempres' , zbg_emp_res                              )   ! volume flux resultant         -
146         
147      CALL iom_put( 'ibgsfx'    , zbg_sfx                                  )   ! salt flux         -(psu*m/day equivalent liquid)       
148      CALL iom_put( 'ibgsfxbri' , zbg_sfx_bri                              )   ! salt flux brines  -     
149      CALL iom_put( 'ibgsfxthd' , zbg_sfx_thd                              )   ! salt flux thermo  -   
150      CALL iom_put( 'ibgsfxmec' , zbg_sfx_mec                              )   ! salt flux dynamic -   
151      CALL iom_put( 'ibgsfxres' , zbg_sfx_res                              )   ! salt flux result  -   
152
153      CALL iom_put( 'ibgfrcvol' , frc_vol * 1.e-9                          )   ! vol - forcing     (km3 equivalent liquid)
154      CALL iom_put( 'ibgfrcsfx' , frc_sal * 1.e-9                          )   ! sal - forcing     (psu*km3 equivalent liquid)   
155      CALL iom_put( 'ibggrme'   , bg_grme * rhoic * r1_rau0 * 1.e-9        )   ! vol growth + melt (km3 equivalent liquid)         
156      !
157      IF( lrst_ice )   CALL lim_diahsb_rst( numit, 'WRITE' )
158      !
159      IF( nn_timing == 1 )   CALL timing_stop('lim_diahsb')
160      !
161   END SUBROUTINE lim_diahsb
162
163
164   SUBROUTINE lim_diahsb_init
165      !!---------------------------------------------------------------------------
166      !!                  ***  ROUTINE lim_diahsb_init  ***
167      !!     
168      !! ** Purpose: Initialization for the heat salt volume budgets
169      !!
170      !! ** Method : Compute initial heat content, salt content and volume
171      !!
172      !! ** Action : - Compute initial heat content, salt content and volume
173      !!             - Initialize forcing trends
174      !!             - Compute coefficients for conversion
175      !!---------------------------------------------------------------------------
176      INTEGER            ::   jk       ! dummy loop indice
177      INTEGER            ::   ierror   ! local integer
178      !!
179      !!NAMELIST/namicehsb/ blabla
180      !!----------------------------------------------------------------------
181      !
182      !!REWIND ( numnam_ice )              ! Read Namelist namicehsb
183      !!READ   ( numnam_ice, namicehsb )
184      !
185      IF(lwp) THEN                   ! Control print
186         WRITE(numout,*)
187         WRITE(numout,*) 'lim_diahsb_init : check the heat and salt budgets'
188         WRITE(numout,*) '~~~~~~~~~~~~'
189      ENDIF
190
191      ! ---------------------------------- !
192      ! 2 - initial conservation variables !
193      ! ---------------------------------- !
194      CALL lim_diahsb_rst( nstart, 'READ' )  !* read or initialize all required files
195      !
196   END SUBROUTINE lim_diahsb_init
197
198
199   SUBROUTINE lim_diahsb_rst( kt, cdrw )
200     !!---------------------------------------------------------------------
201     !!                   ***  ROUTINE limdia_rst  ***
202     !!                     
203     !! ** Purpose :   Read or write DIA file in restart file
204     !!
205     !! ** Method  :   use of IOM library
206     !!----------------------------------------------------------------------
207     INTEGER         , INTENT(in) ::   kt     ! ice time-step
208     CHARACTER(len=*), INTENT(in) ::   cdrw   ! "READ"/"WRITE" flag
209     !
210     INTEGER ::   id1, id2, id3   ! local integers
211     !!----------------------------------------------------------------------
212     !
213     IF( TRIM(cdrw) == 'READ' ) THEN        ! Read/initialise
214        IF( ln_rstart ) THEN                   !* Read the restart file
215           !id1 = iom_varid( numrir, 'frc_vol'  , ldstop = .TRUE. )
216           !
217           IF(lwp) WRITE(numout,*) '~~~~~~~'
218           IF(lwp) WRITE(numout,*) ' lim_diahsb_rst at it= ', kt,' date= ', ndastp
219           IF(lwp) WRITE(numout,*) '~~~~~~~'
220           CALL iom_get( numrir, 'frc_vol', frc_vol )
221           CALL iom_get( numrir, 'frc_sal', frc_sal )
222           CALL iom_get( numrir, 'bg_grme', bg_grme )
223        ELSE
224           IF(lwp) WRITE(numout,*) '~~~~~~~'
225           IF(lwp) WRITE(numout,*) ' lim_diahsb at initial state '
226           IF(lwp) WRITE(numout,*) '~~~~~~~'
227           frc_vol  = 0._wp                                           
228           frc_sal  = 0._wp                                                 
229           bg_grme  = 0._wp                                       
230       ENDIF   
231
232     ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN   ! Create restart file
233        !                                   ! -------------------
234        IF(lwp) WRITE(numout,*) '~~~~~~~'
235        IF(lwp) WRITE(numout,*) ' lim_diahsb_rst at it= ', kt,' date= ', ndastp
236        IF(lwp) WRITE(numout,*) '~~~~~~~'
237        CALL iom_rstput( kt, nitrst, numriw, 'frc_vol'   , frc_vol     )
238        CALL iom_rstput( kt, nitrst, numriw, 'frc_sal'   , frc_sal     )
239        CALL iom_rstput( kt, nitrst, numriw, 'bg_grme'   , bg_grme     )
240        !
241     ENDIF
242     !
243   END SUBROUTINE lim_diahsb_rst
244 
245#else
246   !!----------------------------------------------------------------------
247   !!   Default option :         Empty module          NO LIM sea-ice model
248   !!----------------------------------------------------------------------
249CONTAINS
250   SUBROUTINE lim_diahsb          ! Empty routine
251   END SUBROUTINE lim_diahsb
252#endif
253   !!======================================================================
254END MODULE limdiahsb
Note: See TracBrowser for help on using the repository browser.