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.
icbdia.F90 in branches/UKMO/dev_r5518_GO6_under_ice_relax_dr_hook/NEMOGCM/NEMO/OPA_SRC/ICB – NEMO

source: branches/UKMO/dev_r5518_GO6_under_ice_relax_dr_hook/NEMOGCM/NEMO/OPA_SRC/ICB/icbdia.F90 @ 11738

Last change on this file since 11738 was 11738, checked in by marc, 4 years ago

The Dr Hook changes from my perl code.

File size: 35.9 KB
Line 
1MODULE icbdia
2
3   !!======================================================================
4   !!                       ***  MODULE  icbdia  ***
5   !! Icebergs:  initialise variables for iceberg budgets and diagnostics
6   !!======================================================================
7   !! History : 3.3 !  2010-01  (Martin, Adcroft) Original code
8   !!            -  !  2011-03  (Madec)          Part conversion to NEMO form
9   !!            -  !                            Removal of mapping from another grid
10   !!            -  !  2011-04  (Alderson)       Split into separate modules
11   !!            -  !  2011-05  (Alderson)       Budgets are now all here with lots
12   !!            -  !                            of silly routines to call to get values in
13   !!            -  !                            from the right points in the code
14   !!----------------------------------------------------------------------
15   !!----------------------------------------------------------------------
16   !! icb_dia_init : initialise iceberg budgeting
17   !!----------------------------------------------------------------------
18   USE par_oce        ! ocean parameters
19   USE dom_oce        ! ocean domain
20   USE in_out_manager ! nemo IO
21   USE lib_mpp        ! MPP library
22   USE iom            ! I/O library
23   USE icb_oce        ! iceberg variables
24   USE icbutl         ! iceberg utility routines
25
26   USE yomhook, ONLY: lhook, dr_hook
27   USE parkind1, ONLY: jprb, jpim
28
29   IMPLICIT NONE
30   PRIVATE
31
32   PUBLIC   icb_dia_init      ! routine called in icbini.F90 module
33   PUBLIC   icb_dia           ! routine called in icbstp.F90 module
34   PUBLIC   icb_dia_step      ! routine called in icbstp.F90 module
35   PUBLIC   icb_dia_put       ! routine called in icbstp.F90 module
36   PUBLIC   icb_dia_melt      ! routine called in icbthm.F90 module
37   PUBLIC   icb_dia_size      ! routine called in icbthm.F90 module
38   PUBLIC   icb_dia_speed     ! routine called in icbdyn.F90 module
39   PUBLIC   icb_dia_calve     ! routine called in icbclv.F90 module
40   PUBLIC   icb_dia_income    ! routine called in icbclv.F90 module
41
42   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   berg_melt       ! Melting+erosion rate of icebergs     [kg/s/m2]
43   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   buoy_melt       ! Buoyancy component of melting rate   [kg/s/m2]
44   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   eros_melt       ! Erosion component of melting rate    [kg/s/m2]
45   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   conv_melt       ! Convective component of melting rate [kg/s/m2]
46   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   bits_src        ! Mass flux from berg erosion into bergy bits [kg/s/m2]
47   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   bits_melt       ! Melting rate of bergy bits           [kg/s/m2]
48   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   bits_mass       ! Mass distribution of bergy bits      [kg/s/m2]
49   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   virtual_area    ! Virtual surface coverage by icebergs [m2]
50   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE, PUBLIC  ::   berg_mass       ! Mass distribution                    [kg/m2]
51   REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, PUBLIC  ::   real_calving    ! Calving rate into iceberg class at
52   !                                                                          ! calving locations                    [kg/s]
53   
54   REAL(wp), DIMENSION(:,:)  , ALLOCATABLE ::   tmpc                     ! Temporary work space
55   REAL(wp), DIMENSION(:)    , ALLOCATABLE ::   rsumbuf                  ! Temporary work space to reduce mpp exchanges
56   INTEGER , DIMENSION(:)    , ALLOCATABLE ::   nsumbuf                  ! Temporary work space to reduce mpp exchanges
57
58   REAL(wp)                      ::  berg_melt_net
59   REAL(wp)                      ::  bits_src_net
60   REAL(wp)                      ::  bits_melt_net
61   REAL(wp)                      ::  bits_mass_start     , bits_mass_end
62   REAL(wp)                      ::  floating_heat_start , floating_heat_end
63   REAL(wp)                      ::  floating_mass_start , floating_mass_end
64   REAL(wp)                      ::  bergs_mass_start    , bergs_mass_end
65   REAL(wp)                      ::  stored_start        , stored_heat_start
66   REAL(wp)                      ::  stored_end          , stored_heat_end
67   REAL(wp)                      ::  calving_src_net     , calving_out_net
68   REAL(wp)                      ::  calving_src_heat_net, calving_out_heat_net
69   REAL(wp)                      ::  calving_src_heat_used_net
70   REAL(wp)                      ::  calving_rcv_net  , calving_ret_net  , calving_used_net
71   REAL(wp)                      ::  heat_to_bergs_net, heat_to_ocean_net, melt_net
72   REAL(wp)                      ::  calving_to_bergs_net
73
74   INTEGER                       ::  nbergs_start, nbergs_end, nbergs_calved
75   INTEGER                       ::  nbergs_melted
76   INTEGER                       ::  nspeeding_tickets
77   INTEGER , DIMENSION(nclasses) ::  nbergs_calved_by_class
78
79   !!----------------------------------------------------------------------
80   !! NEMO/OPA 3.3 , NEMO Consortium (2011)
81   !! $Id$
82   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
83   !!----------------------------------------------------------------------
84CONTAINS
85
86   SUBROUTINE icb_dia_init( )
87   INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
88   INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
89   REAL(KIND=jprb)               :: zhook_handle
90
91   CHARACTER(LEN=*), PARAMETER :: RoutineName='ICB_DIA_INIT'
92
93   IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
94
95      !!----------------------------------------------------------------------
96      !!----------------------------------------------------------------------
97      !
98      IF( .NOT. ln_bergdia ) RETURN
99
100      ALLOCATE( berg_melt    (jpi,jpj)   )           ;   berg_melt   (:,:)   = 0._wp
101      ALLOCATE( buoy_melt    (jpi,jpj)   )           ;   buoy_melt   (:,:)   = 0._wp
102      ALLOCATE( eros_melt    (jpi,jpj)   )           ;   eros_melt   (:,:)   = 0._wp
103      ALLOCATE( conv_melt    (jpi,jpj)   )           ;   conv_melt   (:,:)   = 0._wp
104      ALLOCATE( bits_src     (jpi,jpj)   )           ;   bits_src    (:,:)   = 0._wp
105      ALLOCATE( bits_melt    (jpi,jpj)   )           ;   bits_melt   (:,:)   = 0._wp
106      ALLOCATE( bits_mass    (jpi,jpj)   )           ;   bits_mass   (:,:)   = 0._wp
107      ALLOCATE( virtual_area (jpi,jpj)   )           ;   virtual_area(:,:)   = 0._wp
108      ALLOCATE( berg_mass    (jpi,jpj)   )           ;   berg_mass   (:,:)   = 0._wp
109      ALLOCATE( real_calving (jpi,jpj,nclasses) )    ;   real_calving(:,:,:) = 0._wp
110      ALLOCATE( tmpc(jpi,jpj) )                      ;   tmpc        (:,:)   = 0._wp
111
112      nbergs_start              = 0
113      nbergs_end                = 0
114      stored_end                = 0._wp
115      nbergs_start              = 0._wp
116      stored_start              = 0._wp
117      nbergs_melted             = 0
118      nbergs_calved             = 0
119      nbergs_calved_by_class(:) = 0
120      nspeeding_tickets         = 0
121      stored_heat_end           = 0._wp
122      floating_heat_end         = 0._wp
123      floating_mass_end         = 0._wp
124      bergs_mass_end            = 0._wp
125      bits_mass_end             = 0._wp
126      stored_heat_start         = 0._wp
127      floating_heat_start       = 0._wp
128      floating_mass_start       = 0._wp
129      bergs_mass_start          = 0._wp
130      bits_mass_start           = 0._wp
131      bits_mass_end             = 0._wp
132      calving_used_net          = 0._wp
133      calving_to_bergs_net      = 0._wp
134      heat_to_bergs_net         = 0._wp
135      heat_to_ocean_net         = 0._wp
136      calving_rcv_net           = 0._wp
137      calving_ret_net           = 0._wp
138      calving_src_net           = 0._wp
139      calving_out_net           = 0._wp
140      calving_src_heat_net      = 0._wp
141      calving_src_heat_used_net = 0._wp
142      calving_out_heat_net      = 0._wp
143      melt_net                  = 0._wp
144      berg_melt_net             = 0._wp
145      bits_melt_net             = 0._wp
146      bits_src_net              = 0._wp
147
148      floating_mass_start       = icb_utl_mass( first_berg )
149      bergs_mass_start          = icb_utl_mass( first_berg, justbergs=.true. )
150      bits_mass_start           = icb_utl_mass( first_berg, justbits=.true. )
151      IF( lk_mpp ) THEN
152         ALLOCATE( rsumbuf(23) )          ; rsumbuf(:) = 0._wp
153         ALLOCATE( nsumbuf(4+nclasses) )  ; nsumbuf(:) = 0
154         rsumbuf(1) = floating_mass_start
155         rsumbuf(2) = bergs_mass_start
156         rsumbuf(3) = bits_mass_start
157         CALL mpp_sum( rsumbuf(1:3), 3 )
158         floating_mass_start = rsumbuf(1)
159         bergs_mass_start = rsumbuf(2)
160         bits_mass_start = rsumbuf(3)
161      ENDIF
162      !
163   IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
164   END SUBROUTINE icb_dia_init
165
166
167   SUBROUTINE icb_dia( ld_budge )
168      !!----------------------------------------------------------------------
169      !! sum all the things we've accumulated so far in the current processor
170      !! in MPP case then add these sums across all processors
171      !! for this we pack variables into buffer so we only need one mpp_sum
172      !!----------------------------------------------------------------------
173      LOGICAL, INTENT(in) ::   ld_budge
174      !
175      INTEGER             ::   ik
176      REAL(wp)            ::   zunused_calving, ztmpsum, zgrdd_berg_mass, zgrdd_bits_mass
177      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
178      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
179      REAL(KIND=jprb)               :: zhook_handle
180
181      CHARACTER(LEN=*), PARAMETER :: RoutineName='ICB_DIA'
182
183      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
184
185      !!----------------------------------------------------------------------
186      !
187      IF( .NOT. ln_bergdia )   RETURN
188
189      zunused_calving      = SUM( berg_grid%calving(:,:) )
190      ztmpsum              = SUM( berg_grid%floating_melt(:,:) * e1e2t(:,:) * tmask_i(:,:) )
191      melt_net             = melt_net + ztmpsum * berg_dt
192      calving_out_net      = calving_out_net + ( zunused_calving + ztmpsum ) * berg_dt
193      ztmpsum              = SUM( berg_melt(:,:) * e1e2t(:,:) * tmask_i(:,:) )
194      berg_melt_net        = berg_melt_net + ztmpsum * berg_dt
195      ztmpsum              = SUM( bits_src(:,:) * e1e2t(:,:) * tmask_i(:,:) )
196      bits_src_net         = bits_src_net + ztmpsum * berg_dt
197      ztmpsum              = SUM( bits_melt(:,:) * e1e2t(:,:) * tmask_i(:,:) )
198      bits_melt_net        = bits_melt_net + ztmpsum * berg_dt
199      ztmpsum              = SUM( src_calving(:,:) * tmask_i(:,:) )
200      calving_ret_net      = calving_ret_net + ztmpsum * berg_dt
201      ztmpsum              = SUM( berg_grid%calving_hflx(:,:) * e1e2t(:,:) * tmask_i(:,:) )
202      calving_out_heat_net = calving_out_heat_net + ztmpsum * berg_dt   ! Units of J
203
204      IF( ld_budge ) THEN
205         stored_end        = SUM( berg_grid%stored_ice(:,:,:) )
206         stored_heat_end   = SUM( berg_grid%stored_heat(:,:) )
207         floating_mass_end = icb_utl_mass( first_berg )
208         bergs_mass_end    = icb_utl_mass( first_berg,justbergs=.true. )
209         bits_mass_end     = icb_utl_mass( first_berg,justbits=.true. )
210         floating_heat_end = icb_utl_heat( first_berg )
211
212         nbergs_end        = icb_utl_count()
213         zgrdd_berg_mass   = SUM( berg_mass (:,:)*e1e2t(:,:)*tmask_i(:,:) )
214         zgrdd_bits_mass   = SUM( bits_mass(:,:)*e1e2t(:,:)*tmask_i(:,:) )
215
216         IF( lk_mpp ) THEN
217            rsumbuf( 1) = stored_end
218            rsumbuf( 2) = stored_heat_end
219            rsumbuf( 3) = floating_mass_end
220            rsumbuf( 4) = bergs_mass_end
221            rsumbuf( 5) = bits_mass_end
222            rsumbuf( 6) = floating_heat_end
223            rsumbuf( 7) = calving_ret_net
224            rsumbuf( 8) = calving_out_net
225            rsumbuf( 9) = calving_rcv_net
226            rsumbuf(10) = calving_src_net
227            rsumbuf(11) = calving_src_heat_net
228            rsumbuf(12) = calving_src_heat_used_net
229            rsumbuf(13) = calving_out_heat_net
230            rsumbuf(14) = calving_used_net
231            rsumbuf(15) = calving_to_bergs_net
232            rsumbuf(16) = heat_to_bergs_net
233            rsumbuf(17) = heat_to_ocean_net
234            rsumbuf(18) = melt_net
235            rsumbuf(19) = berg_melt_net
236            rsumbuf(20) = bits_src_net
237            rsumbuf(21) = bits_melt_net
238            rsumbuf(22) = zgrdd_berg_mass
239            rsumbuf(23) = zgrdd_bits_mass
240
241            CALL mpp_sum( rsumbuf(1:23), 23)
242
243            stored_end                = rsumbuf( 1)
244            stored_heat_end           = rsumbuf( 2)
245            floating_mass_end         = rsumbuf( 3)
246            bergs_mass_end            = rsumbuf( 4)
247            bits_mass_end             = rsumbuf( 5)
248            floating_heat_end         = rsumbuf( 6)
249            calving_ret_net           = rsumbuf( 7)
250            calving_out_net           = rsumbuf( 8)
251            calving_rcv_net           = rsumbuf( 9)
252            calving_src_net           = rsumbuf(10)
253            calving_src_heat_net      = rsumbuf(11)
254            calving_src_heat_used_net = rsumbuf(12)
255            calving_out_heat_net      = rsumbuf(13)
256            calving_used_net          = rsumbuf(14)
257            calving_to_bergs_net      = rsumbuf(15)
258            heat_to_bergs_net         = rsumbuf(16)
259            heat_to_ocean_net         = rsumbuf(17)
260            melt_net                  = rsumbuf(18)
261            berg_melt_net             = rsumbuf(19)
262            bits_src_net              = rsumbuf(20)
263            bits_melt_net             = rsumbuf(21)
264            zgrdd_berg_mass           = rsumbuf(22)
265            zgrdd_bits_mass           = rsumbuf(23)
266
267            nsumbuf(1) = nbergs_end
268            nsumbuf(2) = nbergs_calved
269            nsumbuf(3) = nbergs_melted
270            nsumbuf(4) = nspeeding_tickets
271            DO ik = 1, nclasses
272               nsumbuf(4+ik) = nbergs_calved_by_class(ik)
273            END DO
274
275            CALL mpp_sum( nsumbuf(1:nclasses+4), nclasses+4 )
276
277            nbergs_end        = nsumbuf(1)
278            nbergs_calved     = nsumbuf(2)
279            nbergs_melted     = nsumbuf(3)
280            nspeeding_tickets = nsumbuf(4)
281            DO ik = 1,nclasses
282               nbergs_calved_by_class(ik)= nsumbuf(4+ik)
283            ENDDO
284
285         ENDIF
286
287         CALL report_state( 'stored ice','kg','',stored_start,'',stored_end,'')
288         CALL report_state( 'floating','kg','',floating_mass_start,'',floating_mass_end,'',nbergs_end)
289         CALL report_state( 'icebergs','kg','',bergs_mass_start,'',bergs_mass_end,'')
290         CALL report_state( 'bits','kg','',bits_mass_start,'',bits_mass_end,'')
291         CALL report_istate( 'berg #','',nbergs_start,'',nbergs_end,'')
292         CALL report_ibudget( 'berg #','calved',nbergs_calved, &
293                                       'melted',nbergs_melted, &
294                                       '#',nbergs_start,nbergs_end)
295         CALL report_budget( 'stored mass','kg','calving used',calving_used_net, &
296                                           'bergs',calving_to_bergs_net, &
297                                           'stored mass',stored_start,stored_end)
298         CALL report_budget( 'floating mass','kg','calving used',calving_to_bergs_net, &
299                                             'bergs',melt_net, &
300                                             'stored mass',floating_mass_start,floating_mass_end)
301         CALL report_budget( 'berg mass','kg','calving',calving_to_bergs_net, &
302                                              'melt+eros',berg_melt_net, &
303                                              'berg mass',bergs_mass_start,bergs_mass_end)
304         CALL report_budget( 'bits mass','kg','eros used',bits_src_net, &
305                                              'bergs',bits_melt_net, &
306                                              'stored mass',bits_mass_start,bits_mass_end)
307         CALL report_budget( 'net mass','kg','recvd',calving_rcv_net, &
308                                             'rtrnd',calving_ret_net, &
309                                             'net mass',stored_start+floating_mass_start, &
310                                                        stored_end+floating_mass_end)
311         CALL report_consistant( 'iceberg mass','kg','gridded',zgrdd_berg_mass,'bergs',bergs_mass_end)
312         CALL report_consistant( 'bits mass','kg','gridded',zgrdd_bits_mass,'bits',bits_mass_end)
313         CALL report_state( 'net heat','J','',stored_heat_start+floating_heat_start,'', &
314                                              stored_heat_end+floating_heat_end,'')
315         CALL report_state( 'stored heat','J','',stored_heat_start,'',stored_heat_end,'')
316         CALL report_state( 'floating heat','J','',floating_heat_start,'',floating_heat_end,'')
317         CALL report_budget( 'net heat','J','net heat',calving_src_heat_net, &
318                                            'net heat',calving_out_heat_net, &
319                                            'net heat',stored_heat_start+floating_heat_start, &
320                                                       stored_heat_end+floating_heat_end)
321         CALL report_budget( 'stored heat','J','calving used',calving_src_heat_used_net, &
322                                               'bergs',heat_to_bergs_net, &
323                                               'net heat',stored_heat_start,stored_heat_end)
324         CALL report_budget( 'flting heat','J','calved',heat_to_bergs_net, &
325                                               'melt',heat_to_ocean_net, &
326                                               'net heat',floating_heat_start,floating_heat_end)
327         IF (nn_verbose_level >= 1) THEN
328            CALL report_consistant( 'top interface','kg','from SIS',calving_src_net, &
329                                    'received',calving_rcv_net)
330            CALL report_consistant( 'bot interface','kg','sent',calving_out_net, &
331                                    'returned',calving_ret_net)
332         ENDIF
333         WRITE( numicb, '("calved by class = ",i6,20(",",i6))') (nbergs_calved_by_class(ik),ik=1,nclasses)
334         IF ( nspeeding_tickets > 0 ) WRITE( numicb, '("speeding tickets issued = ",i6)') nspeeding_tickets
335
336         nbergs_start              = nbergs_end
337         stored_start              = stored_end
338         nbergs_melted             = 0
339         nbergs_calved             = 0
340         nbergs_calved_by_class(:) = 0
341         nspeeding_tickets         = 0
342         stored_heat_start         = stored_heat_end
343         floating_heat_start       = floating_heat_end
344         floating_mass_start       = floating_mass_end
345         bergs_mass_start          = bergs_mass_end
346         bits_mass_start           = bits_mass_end
347         calving_used_net          = 0._wp
348         calving_to_bergs_net      = 0._wp
349         heat_to_bergs_net         = 0._wp
350         heat_to_ocean_net         = 0._wp
351         calving_rcv_net           = 0._wp
352         calving_ret_net           = 0._wp
353         calving_src_net           = 0._wp
354         calving_out_net           = 0._wp
355         calving_src_heat_net      = 0._wp
356         calving_src_heat_used_net = 0._wp
357         calving_out_heat_net      = 0._wp
358         melt_net                  = 0._wp
359         berg_melt_net             = 0._wp
360         bits_melt_net             = 0._wp
361         bits_src_net              = 0._wp
362      ENDIF
363      !
364      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
365   END SUBROUTINE icb_dia
366
367
368   SUBROUTINE icb_dia_step
369   INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
370   INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
371   REAL(KIND=jprb)               :: zhook_handle
372
373   CHARACTER(LEN=*), PARAMETER :: RoutineName='ICB_DIA_STEP'
374
375   IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
376
377      !!----------------------------------------------------------------------
378      !! things to reset at the beginning of each timestep
379      !!----------------------------------------------------------------------
380      !
381      IF( .NOT. ln_bergdia ) RETURN
382      berg_melt    (:,:)   = 0._wp
383      buoy_melt    (:,:)   = 0._wp
384      eros_melt    (:,:)   = 0._wp
385      conv_melt    (:,:)   = 0._wp
386      bits_src     (:,:)   = 0._wp
387      bits_melt    (:,:)   = 0._wp
388      bits_mass    (:,:)   = 0._wp
389      berg_mass    (:,:)   = 0._wp
390      virtual_area (:,:)   = 0._wp
391      real_calving (:,:,:) = 0._wp
392      !
393   IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
394   END SUBROUTINE icb_dia_step
395
396
397   SUBROUTINE icb_dia_put
398   INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
399   INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
400   REAL(KIND=jprb)               :: zhook_handle
401
402   CHARACTER(LEN=*), PARAMETER :: RoutineName='ICB_DIA_PUT'
403
404   IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
405
406      !!----------------------------------------------------------------------
407      !!----------------------------------------------------------------------
408      !
409      IF( .NOT. ln_bergdia )   RETURN            !!gm useless iom will control whether it is output or not
410      !
411      CALL iom_put( "berg_total_melt"      , berg_grid%floating_melt(:,:)   )   ! Total melt flux to ocean      [kg/m2/s]
412      CALL iom_put( "berg_total_heat_flux" , berg_grid%calving_hflx(:,:)    )   ! Total iceberg-ocean heat flux [W/m2]
413      CALL iom_put( "berg_melt"        , berg_melt   (:,:)   )   ! Melt rate of icebergs                     [kg/m2/s]
414      CALL iom_put( "berg_buoy_melt"   , buoy_melt   (:,:)   )   ! Buoyancy component of iceberg melt rate   [kg/m2/s]
415      CALL iom_put( "berg_eros_melt"   , eros_melt   (:,:)   )   ! Erosion component of iceberg melt rate    [kg/m2/s]
416      CALL iom_put( "berg_conv_melt"   , conv_melt   (:,:)   )   ! Convective component of iceberg melt rate [kg/m2/s]
417      CALL iom_put( "berg_virtual_area", virtual_area(:,:)   )   ! Virtual coverage by icebergs              [m2]
418      CALL iom_put( "bits_src"         , bits_src    (:,:)   )   ! Mass source of bergy bits                 [kg/m2/s]
419      CALL iom_put( "bits_melt"        , bits_melt   (:,:)   )   ! Melt rate of bergy bits                   [kg/m2/s]
420      CALL iom_put( "bits_mass"        , bits_mass   (:,:)   )   ! Bergy bit density field                   [kg/m2]
421      CALL iom_put( "berg_mass"        , berg_mass   (:,:)   )   ! Iceberg density field                     [kg/m2]
422      CALL iom_put( "berg_real_calving", real_calving(:,:,:) )   ! Calving into iceberg class                [kg/s]
423      !
424   IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
425   END SUBROUTINE icb_dia_put
426
427
428   SUBROUTINE icb_dia_calve( ki, kj, kn, pcalved, pheated )
429      !!----------------------------------------------------------------------
430      !!----------------------------------------------------------------------
431      INTEGER,  INTENT(in)  ::   ki, kj, kn
432      REAL(wp), INTENT(in)  ::   pcalved
433      REAL(wp), INTENT(in)  ::   pheated
434      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
435      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
436      REAL(KIND=jprb)               :: zhook_handle
437
438      CHARACTER(LEN=*), PARAMETER :: RoutineName='ICB_DIA_CALVE'
439
440      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
441
442      !!----------------------------------------------------------------------
443      !
444      IF( .NOT. ln_bergdia ) RETURN
445      real_calving(ki,kj,kn)     = real_calving(ki,kj,kn) + pcalved / berg_dt
446      nbergs_calved              = nbergs_calved              + 1
447      nbergs_calved_by_class(kn) = nbergs_calved_by_class(kn) + 1
448      calving_to_bergs_net       = calving_to_bergs_net + pcalved
449      heat_to_bergs_net          = heat_to_bergs_net    + pheated
450      !
451      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
452   END SUBROUTINE icb_dia_calve
453
454
455   SUBROUTINE icb_dia_income( kt,  pcalving_used, pheat_used )
456      !!----------------------------------------------------------------------
457      !!----------------------------------------------------------------------
458      INTEGER ,                 INTENT(in)  :: kt
459      REAL(wp),                 INTENT(in)  :: pcalving_used
460      REAL(wp), DIMENSION(:,:), INTENT(in)  :: pheat_used
461      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
462      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
463      REAL(KIND=jprb)               :: zhook_handle
464
465      CHARACTER(LEN=*), PARAMETER :: RoutineName='ICB_DIA_INCOME'
466
467      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
468
469      !!----------------------------------------------------------------------
470      !
471      IF( .NOT. ln_bergdia ) RETURN
472      !
473      IF( kt == nit000 ) THEN
474         stored_start = SUM( berg_grid%stored_ice(:,:,:) )
475         IF( lk_mpp ) CALL mpp_sum( stored_start )
476         WRITE(numicb,'(a,es13.6,a)')   'icb_dia_income: initial stored mass=',stored_start,' kg'
477         !
478         stored_heat_start = SUM( berg_grid%stored_heat(:,:) )
479         IF( lk_mpp ) CALL mpp_sum( stored_heat_start )
480         WRITE(numicb,'(a,es13.6,a)')    'icb_dia_income: initial stored heat=',stored_heat_start,' J'
481      ENDIF
482      !
483      calving_rcv_net = calving_rcv_net + SUM( berg_grid%calving(:,:) ) * berg_dt
484      calving_src_net = calving_rcv_net
485      calving_src_heat_net = calving_src_heat_net +  &
486         &                      SUM( berg_grid%calving_hflx(:,:) * e1e2t(:,:) ) * berg_dt   ! Units of J
487      calving_used_net = calving_used_net + pcalving_used * berg_dt
488      calving_src_heat_used_net = calving_src_heat_used_net + SUM( pheat_used(:,:) )
489      !
490      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
491   END SUBROUTINE icb_dia_income
492
493
494   SUBROUTINE icb_dia_size(ki, kj, pWn, pLn, pAbits,   &
495      &                    pmass_scale, pMnew, pnMbits, pz1_e1e2)
496      !!----------------------------------------------------------------------
497      !!----------------------------------------------------------------------
498      INTEGER,  INTENT(in) :: ki, kj
499      REAL(wp), INTENT(in) :: pWn, pLn, pAbits, pmass_scale, pMnew, pnMbits, pz1_e1e2
500      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
501      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
502      REAL(KIND=jprb)               :: zhook_handle
503
504      CHARACTER(LEN=*), PARAMETER :: RoutineName='ICB_DIA_SIZE'
505
506      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
507
508      !!----------------------------------------------------------------------
509      !
510      IF( .NOT. ln_bergdia ) RETURN
511      virtual_area(ki,kj) = virtual_area(ki,kj) + ( pWn * pLn + pAbits ) * pmass_scale      ! m^2
512      berg_mass(ki,kj)    = berg_mass(ki,kj) + pMnew * pz1_e1e2                             ! kg/m2
513      bits_mass(ki,kj)    = bits_mass(ki,kj) + pnMbits * pz1_e1e2                           ! kg/m2
514      !
515      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
516   END SUBROUTINE icb_dia_size
517
518
519   SUBROUTINE icb_dia_speed()
520   INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
521   INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
522   REAL(KIND=jprb)               :: zhook_handle
523
524   CHARACTER(LEN=*), PARAMETER :: RoutineName='ICB_DIA_SPEED'
525
526   IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
527
528      !!----------------------------------------------------------------------
529      !!----------------------------------------------------------------------
530      !
531      IF( .NOT. ln_bergdia ) RETURN
532      nspeeding_tickets = nspeeding_tickets + 1
533      !
534   IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
535   END SUBROUTINE icb_dia_speed
536
537
538   SUBROUTINE icb_dia_melt(ki, kj, pmnew, pheat, pmass_scale,   &
539      &                   pdM, pdMbitsE, pdMbitsM, pdMb, pdMe,   &
540      &                   pdMv, pz1_dt_e1e2 )
541      !!----------------------------------------------------------------------
542      !!----------------------------------------------------------------------
543      INTEGER , INTENT(in) ::   ki, kj
544      REAL(wp), INTENT(in) ::   pmnew, pheat, pmass_scale
545      REAL(wp), INTENT(in) ::   pdM, pdMbitsE, pdMbitsM, pdMb, pdMe, pdMv, pz1_dt_e1e2
546      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
547      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
548      REAL(KIND=jprb)               :: zhook_handle
549
550      CHARACTER(LEN=*), PARAMETER :: RoutineName='ICB_DIA_MELT'
551
552      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
553
554      !!----------------------------------------------------------------------
555      !
556      IF( .NOT. ln_bergdia ) RETURN
557
558      berg_melt (ki,kj) = berg_melt (ki,kj) + pdM      * pz1_dt_e1e2   ! kg/m2/s
559      bits_src  (ki,kj) = bits_src  (ki,kj) + pdMbitsE * pz1_dt_e1e2   ! mass flux into bergy bitskg/m2/s
560      bits_melt (ki,kj) = bits_melt (ki,kj) + pdMbitsM * pz1_dt_e1e2   ! melt rate of bergy bits kg/m2/s
561      buoy_melt (ki,kj) = buoy_melt (ki,kj) + pdMb     * pz1_dt_e1e2   ! kg/m2/s
562      eros_melt (ki,kj) = eros_melt (ki,kj) + pdMe     * pz1_dt_e1e2   ! erosion rate kg/m2/s
563      conv_melt (ki,kj) = conv_melt (ki,kj) + pdMv     * pz1_dt_e1e2   ! kg/m2/s
564      heat_to_ocean_net = heat_to_ocean_net + pheat * pmass_scale * berg_dt         ! J
565      IF( pmnew <= 0._wp ) nbergs_melted = nbergs_melted + 1                        ! Delete the berg if completely melted
566      !
567      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
568   END SUBROUTINE icb_dia_melt
569
570
571   SUBROUTINE report_state( cd_budgetstr, cd_budgetunits, cd_startstr, pstartval, cd_endstr,   &
572      &                     pendval, cd_delstr, kbergs )
573      !!----------------------------------------------------------------------
574      !!----------------------------------------------------------------------
575      CHARACTER*(*), INTENT(in)           :: cd_budgetstr, cd_budgetunits, cd_startstr, cd_endstr, cd_delstr
576      REAL(wp),      INTENT(in)           :: pstartval, pendval
577      INTEGER,       INTENT(in), OPTIONAL :: kbergs
578      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
579      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
580      REAL(KIND=jprb)               :: zhook_handle
581
582      CHARACTER(LEN=*), PARAMETER :: RoutineName='REPORT_STATE'
583
584      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
585
586      !!----------------------------------------------------------------------
587      !
588      IF ( PRESENT(kbergs) ) THEN
589         WRITE(numicb,100) cd_budgetstr // ' state:',                                    &
590                           cd_startstr  // ' start',  pstartval,         cd_budgetunits, &
591                           cd_endstr    // ' end',    pendval,           cd_budgetunits, &
592                           'Delta '     // cd_delstr, pendval-pstartval, cd_budgetunits, &
593                           '# of bergs', kbergs
594      ELSE
595         WRITE(numicb,100) cd_budgetstr // ' state:',                                   &
596                           cd_startstr  // ' start', pstartval,         cd_budgetunits, &
597                           cd_endstr    // ' end',   pendval,           cd_budgetunits, &
598                           cd_delstr    // 'Delta',  pendval-pstartval, cd_budgetunits
599      ENDIF
600      100 FORMAT(a19,3(a18,"=",es14.7,x,a2,:,","),a12,i8)
601      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
602   END SUBROUTINE report_state
603
604
605   SUBROUTINE report_consistant( cd_budgetstr, cd_budgetunits, cd_startstr, pstartval, cd_endstr, pendval)
606      !!----------------------------------------------------------------------
607      !!----------------------------------------------------------------------
608      CHARACTER*(*), INTENT(in) :: cd_budgetstr, cd_budgetunits, cd_startstr, cd_endstr
609      REAL(wp),      INTENT(in) :: pstartval, pendval
610      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
611      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
612      REAL(KIND=jprb)               :: zhook_handle
613
614      CHARACTER(LEN=*), PARAMETER :: RoutineName='REPORT_CONSISTANT'
615
616      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
617
618      !!----------------------------------------------------------------------
619      !
620      WRITE(numicb,200) cd_budgetstr // ' check:',                 &
621                        cd_startstr,    pstartval, cd_budgetunits, &
622                        cd_endstr,      pendval,   cd_budgetunits, &
623                        'error',        (pendval-pstartval)/((pendval+pstartval)+1e-30), 'nd'
624      200 FORMAT(a19,10(a18,"=",es14.7,x,a2,:,","))
625      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
626   END SUBROUTINE report_consistant
627
628
629   SUBROUTINE report_budget( cd_budgetstr, cd_budgetunits, cd_instr, pinval, cd_outstr,   &
630      &                      poutval, cd_delstr, pstartval, pendval)
631      !!----------------------------------------------------------------------
632      !!----------------------------------------------------------------------
633      CHARACTER*(*), INTENT(in) :: cd_budgetstr, cd_budgetunits, cd_instr, cd_outstr, cd_delstr
634      REAL(wp),      INTENT(in) :: pinval, poutval, pstartval, pendval
635      !
636      REAL(wp)                  :: zval
637      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
638      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
639      REAL(KIND=jprb)               :: zhook_handle
640
641      CHARACTER(LEN=*), PARAMETER :: RoutineName='REPORT_BUDGET'
642
643      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
644
645      !!----------------------------------------------------------------------
646      !
647      zval = ( ( pendval - pstartval ) - ( pinval - poutval ) ) /   &
648         &   MAX( 1.e-30, MAX( abs( pendval - pstartval ) , ABS( pinval - poutval ) ) )
649
650      WRITE(numicb,200) cd_budgetstr // ' budget:', &
651         &              cd_instr     // ' in',      pinval,         cd_budgetunits, &
652         &              cd_outstr    // ' out',     poutval,        cd_budgetunits, &
653         &              'Delta '     // cd_delstr,  pinval-poutval, cd_budgetunits, &
654         &              'error',        zval,                       'nd'
655  200 FORMAT(a19,3(a18,"=",es14.7,x,a2,:,","),a8,"=",es10.3,x,a2)
656      !
657      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
658   END SUBROUTINE report_budget
659
660
661   SUBROUTINE report_istate( cd_budgetstr, cd_startstr, pstartval, cd_endstr, pendval, cd_delstr)
662      !!----------------------------------------------------------------------
663      !!----------------------------------------------------------------------
664      CHARACTER*(*), INTENT(in) :: cd_budgetstr, cd_startstr, cd_endstr, cd_delstr
665      INTEGER,       INTENT(in) :: pstartval, pendval
666      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
667      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
668      REAL(KIND=jprb)               :: zhook_handle
669
670      CHARACTER(LEN=*), PARAMETER :: RoutineName='REPORT_ISTATE'
671
672      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
673
674      !
675      WRITE(numicb,100) cd_budgetstr // ' state:',           &
676         &              cd_startstr  // ' start', pstartval, &
677         &              cd_endstr    // ' end',   pendval,   &
678         &              cd_delstr    // 'Delta',  pendval-pstartval
679  100 FORMAT(a19,3(a18,"=",i14,x,:,","))
680      !
681      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
682   END SUBROUTINE report_istate
683
684
685   SUBROUTINE report_ibudget( cd_budgetstr, cd_instr, pinval, cd_outstr, poutval,   &
686      &                       cd_delstr, pstartval, pendval)
687      !!----------------------------------------------------------------------
688      !!----------------------------------------------------------------------
689      CHARACTER*(*), INTENT(in) :: cd_budgetstr, cd_instr, cd_outstr, cd_delstr
690      INTEGER,       INTENT(in) :: pinval, poutval, pstartval, pendval
691      INTEGER(KIND=jpim), PARAMETER :: zhook_in = 0
692      INTEGER(KIND=jpim), PARAMETER :: zhook_out = 1
693      REAL(KIND=jprb)               :: zhook_handle
694
695      CHARACTER(LEN=*), PARAMETER :: RoutineName='REPORT_IBUDGET'
696
697      IF (lhook) CALL dr_hook(RoutineName,zhook_in,zhook_handle)
698
699      !!----------------------------------------------------------------------
700      !
701      WRITE(numicb,200) cd_budgetstr // ' budget:', &
702                        cd_instr     // ' in',      pinval, &
703                        cd_outstr    // ' out',     poutval, &
704                        'Delta '     // cd_delstr,  pinval-poutval, &
705                        'error',                    ( ( pendval - pstartval ) - ( pinval - poutval ) )
706      200 FORMAT(a19,10(a18,"=",i14,x,:,","))
707      IF (lhook) CALL dr_hook(RoutineName,zhook_out,zhook_handle)
708   END SUBROUTINE report_ibudget
709
710   !!======================================================================
711END MODULE icbdia
Note: See TracBrowser for help on using the repository browser.