source: branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/LIM_SRC_3/limdiahsb.F90 @ 4346

Last change on this file since 4346 was 4346, checked in by clem, 8 years ago

clean the glob_sum diags for LIM3

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