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.
icbini.F90 in NEMO/trunk/src/OCE/ICB – NEMO

source: NEMO/trunk/src/OCE/ICB/icbini.F90

Last change on this file was 15372, checked in by davestorkey, 2 years ago

trunk: fix for wind forcing of icebergs #2728

  • Property svn:keywords set to Id
File size: 25.0 KB
Line 
1MODULE icbini
2   !!======================================================================
3   !!                       ***  MODULE  icbini  ***
4   !! Icebergs:  initialise variables for iceberg tracking
5   !!======================================================================
6   !! History :   -   !  2010-01  (T. Martin & A. Adcroft)  Original code
7   !!            3.3  !  2011-03  (G. Madec)  Part conversion to NEMO form ; Removal of mapping from another grid
8   !!             -   !  2011-04  (S. Alderson)  Split into separate modules ; Restore restart routines
9   !!             -   !  2011-05  (S. Alderson)  generate_test_icebergs restored ; new forcing arrays with extra halo ;
10   !!             -   !                          north fold exchange arrays added
11   !!----------------------------------------------------------------------
12   !!----------------------------------------------------------------------
13   !!   icb_init     : initialise icebergs
14   !!   icb_ini_gen  : generate test icebergs
15   !!   icb_nam      : read iceberg namelist
16   !!----------------------------------------------------------------------
17   USE dom_oce        ! ocean domain
18   USE in_out_manager ! IO routines and numout in particular
19   USE lib_mpp        ! mpi library and lk_mpp in particular
20   USE sbc_oce        ! ocean  : surface boundary condition
21   USE sbc_ice        ! sea-ice: surface boundary condition
22   USE iom            ! IOM library
23   USE fldread        ! field read
24   USE lbclnk         ! lateral boundary condition - MPP link
25   !
26   USE icb_oce        ! define iceberg arrays
27   USE icbutl         ! iceberg utility routines
28   USE icbrst         ! iceberg restart routines
29   USE icbtrj         ! iceberg trajectory I/O routines
30   USE icbdia         ! iceberg budget routines
31
32   IMPLICIT NONE
33   PRIVATE
34
35   PUBLIC   icb_init  ! routine called in nemogcm.F90 module
36
37   CHARACTER(len=100)                                 ::   cn_dir = './'   !: Root directory for location of icb files
38   TYPE(FLD_N)                                        ::   sn_icb          !: information about the calving file to be read
39   TYPE(FLD), PUBLIC, ALLOCATABLE     , DIMENSION(:)  ::   sf_icb          !: structure: file information, fields read
40                                                                           !: used in icbini and icbstp
41   !! * Substitutions
42#  include "do_loop_substitute.h90"
43   !!----------------------------------------------------------------------
44   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
45   !! $Id$
46   !! Software governed by the CeCILL license (see ./LICENSE)
47   !!----------------------------------------------------------------------
48CONTAINS
49
50   SUBROUTINE icb_init( pdt, kt )
51      !!----------------------------------------------------------------------
52      !!                  ***  ROUTINE dom_init  ***
53      !!
54      !! ** Purpose :   iceberg initialization.
55      !!
56      !! ** Method  : - read the iceberg namelist
57      !!              - find non-overlapping processor interior since we can only
58      !!                have one instance of a particular iceberg
59      !!              - calculate the destinations for north fold exchanges
60      !!              - setup either test icebergs or calving file
61      !!----------------------------------------------------------------------
62      REAL(wp), INTENT(in) ::   pdt   ! iceberg time-step (rn_Dt*nn_fsbc)
63      INTEGER , INTENT(in) ::   kt    ! time step number
64      !
65      INTEGER ::   ji, jj, jn               ! dummy loop indices
66      INTEGER ::   i1, i2, i3               ! local integers
67      INTEGER ::   ii, inum, ivar           !   -       -
68      INTEGER ::   istat1, istat2, istat3   !   -       -
69      CHARACTER(len=300) ::   cl_sdist      ! local character
70      !!----------------------------------------------------------------------
71      !
72      CALL icb_nam               ! Read and print namelist parameters
73      !
74      IF( .NOT. ln_icebergs )   RETURN
75      !
76      ALLOCATE( utau_icb(jpi,jpj), vtau_icb(jpi,jpj) )
77      !
78      !                          ! allocate gridded fields
79      IF( icb_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'icb_alloc : unable to allocate arrays' )
80      !
81      !                          ! initialised variable with extra haloes to zero
82      ssu_e(:,:) = 0._wp   ;   ssv_e(:,:) = 0._wp   ;
83      ua_e(:,:)  = 0._wp   ;   va_e(:,:)  = 0._wp   ;
84      ff_e(:,:)  = 0._wp   ;   sst_e(:,:) = 0._wp   ;
85      fr_e(:,:)  = 0._wp   ;   sss_e(:,:) = 0._wp   ;
86      !
87      IF ( ln_M2016 ) THEN
88         toce_e(:,:,:) = 0._wp
89         uoce_e(:,:,:) = 0._wp
90         voce_e(:,:,:) = 0._wp
91         e3t_e(:,:,:)  = 0._wp
92      END IF
93      !
94#if defined key_si3
95      hi_e(:,:) = 0._wp   ;
96      ui_e(:,:) = 0._wp   ;   vi_e(:,:) = 0._wp   ;
97#endif
98      ssh_e(:,:) = 0._wp  ; 
99      !
100      !                          ! open ascii output file or files for iceberg status information
101      !                          ! note that we choose to do this on all processors since we cannot
102      !                          ! predict where icebergs will be ahead of time
103      IF( nn_verbose_level > 0) THEN
104         CALL ctl_opn( numicb, 'icebergs.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
105      ENDIF
106
107      ! set parameters (mostly from namelist)
108      !
109      berg_dt         = pdt
110      first_width (:) = SQRT(  rn_initial_mass(:) / ( rn_LoW_ratio * rn_rho_bergs * rn_initial_thickness(:) )  )
111      first_length(:) = rn_LoW_ratio * first_width(:)
112      rho_berg_1_oce  = rn_rho_bergs / pp_rho_seawater  ! scale factor used for convertion thickness to draft
113      !
114      ! deepest level affected by icebergs
115      ! can be tuned but the safest is this
116      ! (with z* and z~ the depth of each level change overtime, so the more robust micbkb is jpk)
117      micbkb = jpk
118
119      berg_grid%calving      (:,:)   = 0._wp
120      berg_grid%calving_hflx (:,:)   = 0._wp
121      berg_grid%stored_heat  (:,:)   = 0._wp
122      berg_grid%floating_melt(:,:)   = 0._wp
123      berg_grid%maxclass     (:,:)   = nclasses
124      berg_grid%stored_ice   (:,:,:) = 0._wp
125      berg_grid%tmp          (:,:)   = 0._wp
126      src_calving            (:,:)   = 0._wp
127      src_calving_hflx       (:,:)   = 0._wp
128
129      !                          ! domain for icebergs
130      IF( lk_mpp .AND. jpni == 1 )   CALL ctl_stop( 'icbinit: having ONE processor in x currently does not work' )
131      ! NB: the issue here is simply that cyclic east-west boundary condition have not been coded in mpp case
132      ! for the north fold we work out which points communicate by asking
133      ! lbc_lnk to pass processor number (valid even in single processor case)
134      ! borrow src_calving arrays for this
135      !
136      ! pack i and j together using a scaling of a power of 10
137      nicbpack = 10000
138      IF( jpiglo >= nicbpack )   CALL ctl_stop( 'icbini: processor index packing failure' )
139      nicbfldproc(:) = -1
140
141      DO_2D( 1, 1, 1, 1 )
142         src_calving_hflx(ji,jj) = narea
143         src_calving     (ji,jj) = nicbpack * mjg(jj) + mig(ji)
144      END_2D
145      CALL lbc_lnk( 'icbini', src_calving_hflx, 'T', 1._wp )
146      CALL lbc_lnk( 'icbini', src_calving     , 'T', 1._wp )
147
148      ! work out interior of processor from exchange array
149      ! first entry with narea for this processor is left hand interior index
150      ! last  entry                               is right hand interior index
151      jj = jpj/2
152      nicbdi = -1
153      nicbei = -1
154      DO ji = 1, jpi
155         i3 = INT( src_calving(ji,jj) )
156         i2 = INT( i3/nicbpack )
157         i1 = i3 - i2*nicbpack
158         i3 = INT( src_calving_hflx(ji,jj) )
159         IF( i1 == mig(ji) .AND. i3 == narea ) THEN
160            IF( nicbdi < 0 ) THEN   ;   nicbdi = ji
161            ELSE                    ;   nicbei = ji
162            ENDIF
163         ENDIF
164      END DO
165      !
166      ! repeat for j direction
167      ji = jpi/2
168      nicbdj = -1
169      nicbej = -1
170      DO jj = 1, jpj
171         i3 = INT( src_calving(ji,jj) )
172         i2 = INT( i3/nicbpack )
173         i1 = i3 - i2*nicbpack
174         i3 = INT( src_calving_hflx(ji,jj) )
175         IF( i2 == mjg(jj) .AND. i3 == narea ) THEN
176            IF( nicbdj < 0 ) THEN   ;   nicbdj = jj
177            ELSE                    ;   nicbej = jj
178            ENDIF
179         ENDIF
180      END DO
181      !   
182      ! special for east-west boundary exchange we save the destination index
183      i1 = MAX( nicbdi-1, 1)
184      i3 = INT( src_calving(i1,jpj/2) )
185      jj = INT( i3/nicbpack )
186      ricb_left = REAL( i3 - nicbpack*jj, wp ) - (nn_hls-1)
187      i1 = MIN( nicbei+1, jpi )
188      i3 = INT( src_calving(i1,jpj/2) )
189      jj = INT( i3/nicbpack )
190      ricb_right = REAL( i3 - nicbpack*jj, wp ) - (nn_hls-1)
191     
192      ! north fold
193      IF( l_IdoNFold ) THEN
194         !
195         ! icebergs in row nicbej+1 get passed across fold
196         nicbfldpts(:)  = INT( src_calving(:,nicbej+1) )
197         nicbflddest(:) = INT( src_calving_hflx(:,nicbej+1) )
198         !
199         ! work out list of unique processors to talk to
200         ! pack them into a fixed size array where empty slots are marked by a -1
201         DO ji = nicbdi, nicbei
202            ii = nicbflddest(ji)
203            IF( ii .GT. 0 ) THEN     ! Needed because land suppression can mean
204                                     ! that unused points are not set in edge haloes
205               DO jn = 1, jpni
206                  ! work along array until we find an empty slot
207                  IF( nicbfldproc(jn) == -1 ) THEN
208                     nicbfldproc(jn) = ii
209                     EXIT                             !!gm EXIT should be avoided: use DO WHILE expression instead
210                  ENDIF
211                  ! before we find an empty slot, we may find processor number is already here so we exit
212                  IF( nicbfldproc(jn) == ii ) EXIT
213               END DO
214            ENDIF
215         END DO
216      ENDIF
217      !
218      IF( nn_verbose_level > 0) THEN
219         WRITE(numicb,*) 'processor ', narea
220         WRITE(numicb,*) 'jpi, jpj   ', jpi, jpj
221         WRITE(numicb,*) 'Nis0, Nie0 ', Nis0, Nie0
222         WRITE(numicb,*) 'Njs0, Nje0 ', Njs0, Nje0
223         WRITE(numicb,*) 'berg i interior ', nicbdi, nicbei
224         WRITE(numicb,*) 'berg j interior ', nicbdj, nicbej
225         WRITE(numicb,*) 'berg left       ', ricb_left
226         WRITE(numicb,*) 'berg right      ', ricb_right
227         jj = jpj/2
228         WRITE(numicb,*) "central j line:"
229         WRITE(numicb,*) "i processor"
230         WRITE(numicb,*) (INT(src_calving_hflx(ji,jj)), ji=1,jpi)
231         WRITE(numicb,*) "i point"
232         WRITE(numicb,*) (INT(src_calving(ji,jj)), ji=1,jpi)
233         ji = jpi/2
234         WRITE(numicb,*) "central i line:"
235         WRITE(numicb,*) "j processor"
236         WRITE(numicb,*) (INT(src_calving_hflx(ji,jj)), jj=1,jpj)
237         WRITE(numicb,*) "j point"
238         WRITE(numicb,*) (INT(src_calving(ji,jj)), jj=1,jpj)
239         IF( l_IdoNFold ) THEN
240            WRITE(numicb,*) 'north fold destination points '
241            WRITE(numicb,*) nicbfldpts
242            WRITE(numicb,*) 'north fold destination procs  '
243            WRITE(numicb,*) nicbflddest
244            WRITE(numicb,*) 'north fold destination proclist  '
245            WRITE(numicb,*) nicbfldproc
246         ENDIF
247         CALL flush(numicb)
248      ENDIF
249     
250      src_calving     (:,:) = 0._wp
251      src_calving_hflx(:,:) = 0._wp
252
253      ! definition of extended surface masked needed by icb_bilin_h
254      tmask_e(:,:) = 0._wp   ;   tmask_e(1:jpi,1:jpj) = tmask(:,:,1)
255      umask_e(:,:) = 0._wp   ;   umask_e(1:jpi,1:jpj) = umask(:,:,1)
256      vmask_e(:,:) = 0._wp   ;   vmask_e(1:jpi,1:jpj) = vmask(:,:,1)
257      CALL lbc_lnk_icb( 'icbini', tmask_e, 'T', +1._wp, 1, 1 )
258      CALL lbc_lnk_icb( 'icbini', umask_e, 'U', +1._wp, 1, 1 )
259      CALL lbc_lnk_icb( 'icbini', vmask_e, 'V', +1._wp, 1, 1 )
260
261      ! definition of extended lat/lon array needed by icb_bilin_h
262      rlon_e(:,:) = 0._wp     ;  rlon_e(1:jpi,1:jpj) = glamt(:,:) 
263      rlat_e(:,:) = 0._wp     ;  rlat_e(1:jpi,1:jpj) = gphit(:,:)
264      CALL lbc_lnk_icb( 'icbini', rlon_e, 'T', +1._wp, 1, 1 )
265      CALL lbc_lnk_icb( 'icbini', rlat_e, 'T', +1._wp, 1, 1 )
266      !
267      ! definnitionn of extennded ff_f array needed by icb_utl_interp
268      ff_e(:,:) = 0._wp       ;  ff_e(1:jpi,1:jpj) = ff_f(:,:)
269      CALL lbc_lnk_icb( 'icbini', ff_e, 'F', +1._wp, 1, 1 )
270
271      ! assign each new iceberg with a unique number constructed from the processor number
272      ! and incremented by the total number of processors
273      num_bergs(:) = 0
274      num_bergs(1) = narea - jpnij
275
276      ! when not generating test icebergs we need to setup calving file
277      IF( nn_test_icebergs < 0 .OR. ln_use_calving ) THEN
278         !
279         ! maximum distribution class array does not change in time so read it once
280         cl_sdist = TRIM( cn_dir )//TRIM( sn_icb%clname )
281         CALL iom_open ( cl_sdist, inum )                              ! open file
282         ivar = iom_varid( inum, 'maxclass', ldstop=.FALSE. )
283         IF( ivar > 0 ) THEN
284            CALL iom_get  ( inum, jpdom_global, 'maxclass', src_calving )   ! read the max distribution array
285            berg_grid%maxclass(:,:) = INT( src_calving )
286            src_calving(:,:) = 0._wp
287         ENDIF
288         CALL iom_close( inum )                                     ! close file
289         !
290         IF( nn_verbose_level > 0) THEN
291            WRITE(numicb,*)
292            WRITE(numicb,*) '          calving read in a file'
293         ENDIF
294         ALLOCATE( sf_icb(1), STAT=istat1 )         ! Create sf_icb structure (calving)
295         ALLOCATE( sf_icb(1)%fnow(jpi,jpj,1), STAT=istat2 )
296         ALLOCATE( sf_icb(1)%fdta(jpi,jpj,1,2), STAT=istat3 )
297         IF( istat1+istat2+istat3 > 0 ) THEN
298            CALL ctl_stop( 'sbc_icb: unable to allocate sf_icb structure' )   ;   RETURN
299         ENDIF
300         !                                          ! fill sf_icb with the namelist (sn_icb) and control print
301         CALL fld_fill( sf_icb, (/ sn_icb /), cn_dir, 'icb_init', 'read calving data', 'namicb' )
302         !
303      ENDIF
304
305      IF( .NOT.ln_rstart ) THEN
306         IF( nn_test_icebergs > 0 )   CALL icb_ini_gen()
307      ELSE
308         IF( nn_test_icebergs > 0 ) THEN
309            CALL icb_ini_gen()
310         ELSE
311            CALL icb_rst_read()
312            l_restarted_bergs = .TRUE.
313         ENDIF
314      ENDIF
315      !
316      IF( nn_sample_rate .GT. 0 ) CALL icb_trj_init( nitend )
317      !
318      CALL icb_dia_init()
319      !
320      IF( nn_verbose_level >= 2 )   CALL icb_utl_print('icb_init, initial status', nit000-1)
321      !
322   END SUBROUTINE icb_init
323
324
325   SUBROUTINE icb_ini_gen()
326      !!----------------------------------------------------------------------
327      !!                  ***  ROUTINE icb_ini_gen  ***
328      !!
329      !! ** Purpose :   iceberg generation
330      !!
331      !! ** Method  : - at each grid point of the test box supplied in the namelist
332      !!                generate an iceberg in one class determined by the value of
333      !!                parameter nn_test_icebergs
334      !!----------------------------------------------------------------------
335      INTEGER                         ::   ji, jj, ibergs
336      TYPE(iceberg)                   ::   localberg ! NOT a pointer but an actual local variable
337      TYPE(point)                     ::   localpt
338      INTEGER                         ::   iyr, imon, iday, ihr, imin, isec
339      INTEGER                         ::   iberg
340      !!----------------------------------------------------------------------
341
342      ! For convenience
343      iberg = nn_test_icebergs
344
345      ! call get_date(Time, iyr, imon, iday, ihr, imin, isec)
346      ! Convert nemo time variables from dom_oce into local versions
347      iyr  = nyear
348      imon = nmonth
349      iday = nday
350      ihr = INT(nsec_day/3600)
351      imin = INT((nsec_day-ihr*3600)/60)
352      isec = nsec_day - ihr*3600 - imin*60
353
354      ! no overlap for icebergs since we want only one instance of each across the whole domain
355      ! so restrict area of interest
356      ! use tmask here because tmask_i has been doctored on one side of the north fold line
357
358      DO jj = nicbdj, nicbej
359         DO ji = nicbdi, nicbei
360            IF( tmask(ji,jj,1) > 0._wp        .AND.                                       &
361                rn_test_box(1) < glamt(ji,jj) .AND. glamt(ji,jj) < rn_test_box(2) .AND.   &
362                rn_test_box(3) < gphit(ji,jj) .AND. gphit(ji,jj) < rn_test_box(4) ) THEN
363               localberg%mass_scaling = rn_mass_scaling(iberg)
364               localpt%xi = REAL( mig(ji) - (nn_hls-1), wp )
365               localpt%yj = REAL( mjg(jj) - (nn_hls-1), wp )
366               CALL icb_utl_interp( localpt%xi, localpt%yj, plat=localpt%lat, plon=localpt%lon )   
367               localpt%mass      = rn_initial_mass     (iberg)
368               localpt%thickness = rn_initial_thickness(iberg)
369               localpt%width  = first_width (iberg)
370               localpt%length = first_length(iberg)
371               localpt%year = iyr
372               localpt%day = REAL(iday,wp)+(REAL(ihr,wp)+REAL(imin,wp)/60._wp)/24._wp
373               localpt%mass_of_bits = 0._wp
374               localpt%heat_density = 0._wp
375               localpt%uvel = 0._wp
376               localpt%vvel = 0._wp
377               localpt%kb   = 1
378               CALL icb_utl_incr()
379               localberg%number(:) = num_bergs(:)
380               call icb_utl_add(localberg, localpt)
381            ENDIF
382         END DO
383      END DO
384      !
385      ibergs = icb_utl_count()
386      CALL mpp_sum('icbini', ibergs)
387      IF( nn_verbose_level > 0) THEN
388         WRITE(numicb,'(a,i6,a)') 'diamonds, icb_ini_gen: ',ibergs,' were generated'
389      ENDIF
390      !
391   END SUBROUTINE icb_ini_gen
392
393
394   SUBROUTINE icb_nam
395      !!----------------------------------------------------------------------
396      !!                     ***  ROUTINE icb_nam  ***
397      !!
398      !! ** Purpose :   read iceberg namelist and print the variables.
399      !!
400      !! ** input   : - namberg namelist
401      !!----------------------------------------------------------------------
402      INTEGER  ::   jn      ! dummy loop indices
403      INTEGER  ::   ios     ! Local integer output status for namelist read
404      REAL(wp) ::   zfact   ! local scalar
405      !
406      NAMELIST/namberg/ ln_icebergs    , ln_bergdia     , nn_sample_rate      , rn_initial_mass      ,   &
407         &              rn_distribution, rn_mass_scaling, rn_initial_thickness, nn_verbose_write     ,   &
408         &              rn_rho_bergs   , rn_LoW_ratio   , nn_verbose_level    , ln_operator_splitting,   &
409         &              rn_bits_erosion_fraction        , rn_sicn_shift       , ln_passive_mode      ,   &
410         &              ln_time_average_weight          , nn_test_icebergs    , rn_test_box          ,   &
411         &              ln_use_calving , rn_speed_limit , cn_dir, sn_icb      , ln_M2016             ,   &
412         &              cn_icbrst_indir, cn_icbrst_in   , cn_icbrst_outdir    , cn_icbrst_out        ,   &
413         &              ln_icb_grd
414      !!----------------------------------------------------------------------
415
416#if defined key_agrif
417      IF(lwp) THEN
418         WRITE(numout,*)
419         WRITE(numout,*) 'icb_nam : AGRIF is not compatible with namelist namberg :  '
420         WRITE(numout,*) '~~~~~~~   definition of rn_initial_mass(nclasses) with nclasses as PARAMETER '
421         WRITE(numout,*)
422         WRITE(numout,*) '   ==>>>   force  NO icebergs used. The namelist namberg is not read'
423      ENDIF
424      ln_icebergs = .false.     
425      RETURN
426#else
427      IF(lwp) THEN
428         WRITE(numout,*)
429         WRITE(numout,*) 'icb_nam : iceberg initialization through namberg namelist read'
430         WRITE(numout,*) '~~~~~~~~ '
431      ENDIF
432#endif   
433      !                             !==  read namelist  ==!
434      READ  ( numnam_ref, namberg, IOSTAT = ios, ERR = 901)
435901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namberg in reference namelist' )
436      READ  ( numnam_cfg, namberg, IOSTAT = ios, ERR = 902 )
437902   IF( ios >  0 ) CALL ctl_nam ( ios , 'namberg in configuration namelist' )
438      IF(lwm) WRITE ( numond, namberg )
439      !
440      IF(lwp) WRITE(numout,*)
441      IF( ln_icebergs ) THEN
442         IF(lwp) WRITE(numout,*) '   ==>>>   icebergs are used'
443      ELSE
444         IF(lwp) WRITE(numout,*) '   ==>>>   No icebergs used'
445         RETURN
446      ENDIF
447      !
448      IF( nn_test_icebergs > nclasses ) THEN
449         IF(lwp) WRITE(numout,*)
450         IF(lwp) WRITE(numout,*) '   ==>>>   Resetting of nn_test_icebergs to ', nclasses
451         nn_test_icebergs = nclasses
452      ENDIF
453      !
454      IF( nn_test_icebergs < 0 .AND. .NOT. ln_use_calving ) THEN
455         IF(lwp) WRITE(numout,*)
456         IF(lwp) WRITE(numout,*) '   ==>>>   Resetting ln_use_calving to .true. since we are not using test icebergs'
457         ln_use_calving = .true.
458      ENDIF
459      !
460      IF(lwp) THEN                  ! control print
461         WRITE(numout,*)
462         WRITE(numout,*) 'icb_nam : iceberg initialization through namberg namelist read'
463         WRITE(numout,*) '~~~~~~~~ '
464         WRITE(numout,*) '   Calculate budgets                                            ln_bergdia       = ', ln_bergdia
465         WRITE(numout,*) '   Period between sampling of position for trajectory storage   nn_sample_rate = ', nn_sample_rate
466         WRITE(numout,*) '   Mass thresholds between iceberg classes (kg)                 rn_initial_mass     ='
467         DO jn = 1, nclasses
468            WRITE(numout,'(a,f15.2)') '                                                                ', rn_initial_mass(jn)
469         ENDDO
470         WRITE(numout,*) '   Fraction of calving to apply to this class (non-dim)         rn_distribution     ='
471         DO jn = 1, nclasses
472            WRITE(numout,'(a,f10.4)') '                                                                ', rn_distribution(jn)
473         END DO
474         WRITE(numout,*) '   Ratio between effective and real iceberg mass (non-dim)      rn_mass_scaling     = '
475         DO jn = 1, nclasses
476            WRITE(numout,'(a,f10.2)') '                                                                ', rn_mass_scaling(jn)
477         END DO
478         WRITE(numout,*) '   Total thickness of newly calved bergs (m)                    rn_initial_thickness = '
479         DO jn = 1, nclasses
480            WRITE(numout,'(a,f10.2)') '                                                                ', rn_initial_thickness(jn)
481         END DO
482         WRITE(numout,*) '   Timesteps between verbose messages                           nn_verbose_write    = ', nn_verbose_write
483
484         WRITE(numout,*) '   Density of icebergs                           rn_rho_bergs  = ', rn_rho_bergs
485         WRITE(numout,*) '   Initial ratio L/W for newly calved icebergs   rn_LoW_ratio  = ', rn_LoW_ratio
486         WRITE(numout,*) '   Turn on more verbose output                          level  = ', nn_verbose_level
487         WRITE(numout,*) '   Use first order operator splitting for thermodynamics    ',   &
488            &                    'use_operator_splitting = ', ln_operator_splitting
489         WRITE(numout,*) '   Fraction of erosion melt flux to divert to bergy bits    ',   &
490            &                    'bits_erosion_fraction = ', rn_bits_erosion_fraction
491
492         WRITE(numout,*) '   Use icb module modification from Merino et al. (2016) : ln_M2016 = ', ln_M2016
493         WRITE(numout,*) '       ground icebergs if icb bottom lvl hit the oce bottom level : ln_icb_grd = ', ln_icb_grd
494
495         WRITE(numout,*) '   Shift of sea-ice concentration in erosion flux modulation ',   &
496            &                    '(0<sicn_shift<1)    rn_sicn_shift  = ', rn_sicn_shift
497         WRITE(numout,*) '   Do not add freshwater flux from icebergs to ocean                ',   &
498            &                    '                  passive_mode            = ', ln_passive_mode
499         WRITE(numout,*) '   Time average the weight on the ocean   time_average_weight       = ', ln_time_average_weight
500         WRITE(numout,*) '   Create icebergs in absence of a restart file   nn_test_icebergs  = ', nn_test_icebergs
501         WRITE(numout,*) '                   in lon/lat box                                   = ', rn_test_box
502         WRITE(numout,*) '   Use calving data even if nn_test_icebergs > 0    ln_use_calving  = ', ln_use_calving
503         WRITE(numout,*) '   CFL speed limit for a berg            speed_limit                = ', rn_speed_limit
504         WRITE(numout,*) '   Writing Iceberg status information to icebergs.stat file        '
505      ENDIF
506      !
507      ! ensure that the sum of berg input distribution is equal to one
508      zfact = SUM( rn_distribution )
509      IF( zfact /= 1._wp .AND. 0_wp /= zfact ) THEN
510         rn_distribution(:) = rn_distribution(:) / zfact
511         IF(lwp) THEN
512            WRITE(numout,*)
513            WRITE(numout,*) '      ==>>> CAUTION:    sum of berg input distribution = ', zfact
514            WRITE(numout,*) '            *******     redistribution has been rescaled'
515            WRITE(numout,*) '                        updated berg distribution is :'
516            DO jn = 1, nclasses
517               WRITE(numout,'(a,f10.4)') '                                   ',rn_distribution(jn)
518            END DO
519         ENDIF
520      ENDIF
521      IF( MINVAL( rn_distribution(:) ) < 0._wp ) THEN
522         CALL ctl_stop( 'icb_nam: a negative rn_distribution value encountered ==>> change your namelist namberg' )
523      ENDIF
524      !
525   END SUBROUTINE icb_nam
526
527   !!======================================================================
528END MODULE icbini
Note: See TracBrowser for help on using the repository browser.