1 | MODULE icbini |
---|
2 | |
---|
3 | !!====================================================================== |
---|
4 | !! *** MODULE icbini *** |
---|
5 | !! Ocean physics: initialise variables for iceberg tracking |
---|
6 | !!====================================================================== |
---|
7 | !! History : 3.3.1 ! 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 | !! - ! Restore restart routines |
---|
12 | !! - ! 2011-05 (Alderson) generate_test_icebergs restored |
---|
13 | !! - ! 2011-05 (Alderson) new forcing arrays with extra halo |
---|
14 | !! - ! 2011-05 (Alderson) north fold exchange arrays added |
---|
15 | !!---------------------------------------------------------------------- |
---|
16 | !!---------------------------------------------------------------------- |
---|
17 | !! icb_init : initialise |
---|
18 | !! icb_gen : generate test icebergs |
---|
19 | !! icb_nam : read iceberg namelist |
---|
20 | !!---------------------------------------------------------------------- |
---|
21 | USE dom_oce ! ocean domain |
---|
22 | USE in_out_manager ! IO routines and numout in particular |
---|
23 | USE lib_mpp ! mpi library and lk_mpp in particular |
---|
24 | USE sbc_oce |
---|
25 | USE iom |
---|
26 | USE fldread |
---|
27 | USE lbclnk |
---|
28 | |
---|
29 | USE icb_oce ! define iceberg arrays |
---|
30 | USE icbutl ! iceberg utility routines |
---|
31 | USE icbrst ! iceberg restart routines |
---|
32 | USE icbtrj ! iceberg trajectory I/O routines |
---|
33 | USE icbdia ! iceberg budget routines |
---|
34 | |
---|
35 | IMPLICIT NONE |
---|
36 | PRIVATE |
---|
37 | |
---|
38 | CHARACTER(len=100), PRIVATE :: cn_dir = './' ! Root directory for location of icb files |
---|
39 | TYPE(FLD_N), PRIVATE :: sn_icb ! information about the calving file to be read |
---|
40 | |
---|
41 | PUBLIC icb_init ! routine called in nemogcm.F90 module |
---|
42 | PUBLIC icb_gen ! routine called in icbclv.F90 module |
---|
43 | PRIVATE icb_nam |
---|
44 | PRIVATE icb_alloc |
---|
45 | |
---|
46 | CONTAINS |
---|
47 | |
---|
48 | INTEGER FUNCTION icb_alloc() |
---|
49 | !!---------------------------------------------------------------------- |
---|
50 | !! *** ROUTINE icb_alloc *** |
---|
51 | !!---------------------------------------------------------------------- |
---|
52 | ! |
---|
53 | INTEGER :: ill |
---|
54 | ! |
---|
55 | ! open ascii output file or files for iceberg status information |
---|
56 | ! note that we choose to do this on all processors since we cannot predict where |
---|
57 | ! icebergs will be ahead of time |
---|
58 | ! |
---|
59 | CALL ctl_opn( numicb, 'icebergs.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', & |
---|
60 | -1, numout, lwp, narea ) |
---|
61 | |
---|
62 | icb_alloc = 0 |
---|
63 | ALLOCATE(berg_grid, STAT=ill) |
---|
64 | icb_alloc = icb_alloc + ill |
---|
65 | ! |
---|
66 | ALLOCATE( berg_grid%calving (jpi,jpj) , STAT=ill) |
---|
67 | icb_alloc = icb_alloc + ill |
---|
68 | ALLOCATE( berg_grid%calving_hflx (jpi,jpj) , STAT=ill) |
---|
69 | icb_alloc = icb_alloc + ill |
---|
70 | ALLOCATE( berg_grid%stored_heat (jpi,jpj) , STAT=ill) |
---|
71 | icb_alloc = icb_alloc + ill |
---|
72 | ALLOCATE( berg_grid%floating_melt(jpi,jpj) , STAT=ill) |
---|
73 | icb_alloc = icb_alloc + ill |
---|
74 | ALLOCATE( berg_grid%maxclass(jpi,jpj) , STAT=ill) |
---|
75 | icb_alloc = icb_alloc + ill |
---|
76 | ! |
---|
77 | ALLOCATE( berg_grid%stored_ice (jpi,jpj,nclasses) , STAT=ill) |
---|
78 | icb_alloc = icb_alloc + ill |
---|
79 | ! |
---|
80 | ALLOCATE( berg_grid%tmp (jpi,jpj) , STAT=ill) |
---|
81 | icb_alloc = icb_alloc + ill |
---|
82 | ! |
---|
83 | ! expanded arrays for bilinear interpolation |
---|
84 | ALLOCATE( uo_e(0:jpi+1,0:jpj+1) , STAT=ill) |
---|
85 | icb_alloc = icb_alloc + ill |
---|
86 | ALLOCATE( vo_e(0:jpi+1,0:jpj+1) , STAT=ill) |
---|
87 | icb_alloc = icb_alloc + ill |
---|
88 | ALLOCATE( ff_e(0:jpi+1,0:jpj+1) , STAT=ill) |
---|
89 | icb_alloc = icb_alloc + ill |
---|
90 | ALLOCATE( ua_e(0:jpi+1,0:jpj+1) , STAT=ill) |
---|
91 | icb_alloc = icb_alloc + ill |
---|
92 | ALLOCATE( va_e(0:jpi+1,0:jpj+1) , STAT=ill) |
---|
93 | icb_alloc = icb_alloc + ill |
---|
94 | #if defined key_lim2 || defined key_lim3 |
---|
95 | ALLOCATE( ui_e(0:jpi+1,0:jpj+1) , STAT=ill) |
---|
96 | icb_alloc = icb_alloc + ill |
---|
97 | ALLOCATE( vi_e(0:jpi+1,0:jpj+1) , STAT=ill) |
---|
98 | icb_alloc = icb_alloc + ill |
---|
99 | #endif |
---|
100 | ALLOCATE( ssh_e(0:jpi+1,0:jpj+1) , STAT=ill) |
---|
101 | icb_alloc = icb_alloc + ill |
---|
102 | ALLOCATE( first_width (nclasses) , STAT=ill) |
---|
103 | icb_alloc = icb_alloc + ill |
---|
104 | ALLOCATE( first_length (nclasses) , STAT=ill) |
---|
105 | icb_alloc = icb_alloc + ill |
---|
106 | ALLOCATE( src_calving(jpi,jpj) , STAT=ill) |
---|
107 | icb_alloc = icb_alloc + ill |
---|
108 | ALLOCATE( src_calving_hflx(jpi,jpj) , STAT=ill) |
---|
109 | icb_alloc = icb_alloc + ill |
---|
110 | |
---|
111 | ALLOCATE( nicbfldpts(jpi) , STAT=ill) |
---|
112 | icb_alloc = icb_alloc + ill |
---|
113 | ALLOCATE( nicbflddest(jpi) , STAT=ill) |
---|
114 | icb_alloc = icb_alloc + ill |
---|
115 | ALLOCATE( nicbfldproc(jpni) , STAT=ill) |
---|
116 | icb_alloc = icb_alloc + ill |
---|
117 | |
---|
118 | IF( lk_mpp ) CALL mpp_sum ( icb_alloc ) |
---|
119 | IF( icb_alloc > 0 ) CALL ctl_warn('icb_alloc: allocation of arrays failed') |
---|
120 | |
---|
121 | END FUNCTION icb_alloc |
---|
122 | |
---|
123 | SUBROUTINE icb_init( pdt, kt ) |
---|
124 | !!---------------------------------------------------------------------- |
---|
125 | !! *** ROUTINE dom_init *** |
---|
126 | !! |
---|
127 | !! ** Purpose : iceberg initialization. |
---|
128 | !! |
---|
129 | !! ** Method : - blah blah |
---|
130 | !!---------------------------------------------------------------------- |
---|
131 | REAL(wp) , INTENT(in) :: pdt ! iceberg time-step (rdt*nn_fsbc) |
---|
132 | INTEGER , INTENT(in) :: kt ! time step number |
---|
133 | ! |
---|
134 | INTEGER :: ji, jj, jn ! loop counters |
---|
135 | INTEGER :: i1, i2, i3 ! dummy integers |
---|
136 | INTEGER :: ii, inum, ivar |
---|
137 | INTEGER :: istat1, istat2, istat3 |
---|
138 | CHARACTER(len=300) :: cl_sdist |
---|
139 | !!---------------------------------------------------------------------- |
---|
140 | |
---|
141 | CALL icb_nam ! Read and print namelist parameters |
---|
142 | IF( .NOT. ln_icebergs ) RETURN |
---|
143 | |
---|
144 | ! ! allocate gridded fields |
---|
145 | IF( icb_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'icb_alloc : unable to allocate arrays' ) |
---|
146 | |
---|
147 | ! set parameters (mostly from namelist) |
---|
148 | ! |
---|
149 | berg_dt = pdt |
---|
150 | first_width (:) = SQRT( rn_initial_mass(:) / ( rn_LoW_ratio * rn_rho_bergs * rn_initial_thickness(:) ) ) |
---|
151 | first_length(:) = rn_LoW_ratio * first_width(:) |
---|
152 | |
---|
153 | berg_grid%calving (:,:) = 0._wp |
---|
154 | berg_grid%calving_hflx (:,:) = 0._wp |
---|
155 | berg_grid%stored_heat (:,:) = 0._wp |
---|
156 | berg_grid%floating_melt(:,:) = 0._wp |
---|
157 | berg_grid%maxclass (:,:) = nclasses |
---|
158 | berg_grid%stored_ice (:,:,:) = 0._wp |
---|
159 | berg_grid%tmp (:,:) = 0._wp |
---|
160 | src_calving (:,:) = 0._wp |
---|
161 | src_calving_hflx (:,:) = 0._wp |
---|
162 | |
---|
163 | ! domain for icebergs |
---|
164 | IF( lk_mpp .AND. jpni == 1 ) & |
---|
165 | CALL ctl_stop('icbinit: having ONE processor in x currently does not work') |
---|
166 | |
---|
167 | ! for the north fold we work out which points communicate by asking |
---|
168 | ! lbc_lnk to pass processor number (valid even in single processor case) |
---|
169 | ! borrow src_calving arrays for this |
---|
170 | ! |
---|
171 | ! pack i and j together using a scaling of a power of 10 |
---|
172 | nicbpack = 10000 |
---|
173 | IF( jpiglo >= nicbpack ) CALL ctl_stop('icbini: processor index packing failure') |
---|
174 | nicbfldproc(:) = -1 |
---|
175 | |
---|
176 | DO jj = 1, jpj |
---|
177 | DO ji = 1, jpi |
---|
178 | src_calving_hflx(ji,jj) = narea |
---|
179 | src_calving(ji,jj) = nicbpack*(njmpp+jj-1) + nimpp+ji-1 |
---|
180 | ENDDO |
---|
181 | ENDDO |
---|
182 | CALL lbc_lnk( src_calving_hflx, 'T', 1._wp ) |
---|
183 | CALL lbc_lnk( src_calving , 'T', 1._wp ) |
---|
184 | |
---|
185 | ! work out interior of processor from exchange array |
---|
186 | ! first entry with narea for this processor is left hand interior index |
---|
187 | ! last entry is right hand interior index |
---|
188 | jj = jpj/2 |
---|
189 | nicbdi = -1 |
---|
190 | nicbei = -1 |
---|
191 | DO ji = 1,jpi |
---|
192 | i3 = INT( src_calving(ji,jj) ) |
---|
193 | i2 = INT( i3/nicbpack ) |
---|
194 | i1 = i3 - i2*nicbpack |
---|
195 | IF( i1 == nimpp+ji-1 ) THEN |
---|
196 | IF( nicbdi < 0 ) THEN |
---|
197 | nicbdi = ji |
---|
198 | ELSE |
---|
199 | nicbei = ji |
---|
200 | ENDIF |
---|
201 | ENDIF |
---|
202 | END DO |
---|
203 | |
---|
204 | ! repeat for j direction |
---|
205 | ji = jpi/2 |
---|
206 | nicbdj = -1 |
---|
207 | nicbej = -1 |
---|
208 | DO jj = 1,jpj |
---|
209 | i3 = INT( src_calving(ji,jj) ) |
---|
210 | i2 = INT( i3/nicbpack ) |
---|
211 | i1 = i3 - i2*nicbpack |
---|
212 | IF( i2 == njmpp+jj-1 ) THEN |
---|
213 | IF( nicbdj < 0 ) THEN |
---|
214 | nicbdj = jj |
---|
215 | ELSE |
---|
216 | nicbej = jj |
---|
217 | ENDIF |
---|
218 | ENDIF |
---|
219 | END DO |
---|
220 | |
---|
221 | ! special for east-west boundary exchange we save the destination index |
---|
222 | i1 = MAX( nicbdi-1, 1) |
---|
223 | i3 = INT( src_calving(i1,jpj/2) ) |
---|
224 | jj = INT( i3/nicbpack ) |
---|
225 | ricb_left = REAL( i3 - nicbpack*jj, wp ) |
---|
226 | i1 = MIN( nicbei+1, jpi ) |
---|
227 | i3 = INT( src_calving(i1,jpj/2) ) |
---|
228 | jj = INT( i3/nicbpack ) |
---|
229 | ricb_right = REAL( i3 - nicbpack*jj, wp ) |
---|
230 | |
---|
231 | ! north fold |
---|
232 | IF( npolj > 0 ) THEN |
---|
233 | ! |
---|
234 | ! icebergs in row nicbej+1 get passed across fold |
---|
235 | nicbfldpts(:) = INT( src_calving(:,nicbej+1) ) |
---|
236 | nicbflddest(:) = INT( src_calving_hflx(:,nicbej+1) ) |
---|
237 | ! |
---|
238 | ! work out list of unique processors to talk to |
---|
239 | DO ji = nicbdi, nicbei |
---|
240 | ii = nicbflddest(ji) |
---|
241 | DO jn = 1, jpni |
---|
242 | IF( nicbfldproc(jn) == -1 ) THEN |
---|
243 | nicbfldproc(jn) = ii |
---|
244 | EXIT |
---|
245 | ENDIF |
---|
246 | IF( nicbfldproc(jn) == ii ) EXIT |
---|
247 | ENDDO |
---|
248 | ENDDO |
---|
249 | ENDIF |
---|
250 | ! |
---|
251 | IF( nn_verbose_level > 0) THEN |
---|
252 | WRITE(numicb,*) 'processor ', narea |
---|
253 | WRITE(numicb,*) 'jpi, jpj ', jpi, jpj |
---|
254 | WRITE(numicb,*) 'nldi, nlei ', nldi, nlei |
---|
255 | WRITE(numicb,*) 'nldj, nlej ', nldj, nlej |
---|
256 | WRITE(numicb,*) 'berg i interior ', nicbdi, nicbei |
---|
257 | WRITE(numicb,*) 'berg j interior ', nicbdj, nicbej |
---|
258 | WRITE(numicb,*) 'berg left ', ricb_left |
---|
259 | WRITE(numicb,*) 'berg right ', ricb_right |
---|
260 | jj = jpj/2 |
---|
261 | WRITE(numicb,*) "central j line:" |
---|
262 | WRITE(numicb,*) "i processor" |
---|
263 | WRITE(numicb,*) (INT(src_calving_hflx(ji,jj)), ji=1,jpi) |
---|
264 | WRITE(numicb,*) "i point" |
---|
265 | WRITE(numicb,*) (INT(src_calving(ji,jj)), ji=1,jpi) |
---|
266 | ji = jpi/2 |
---|
267 | WRITE(numicb,*) "central i line:" |
---|
268 | WRITE(numicb,*) "j processor" |
---|
269 | WRITE(numicb,*) (INT(src_calving_hflx(ji,jj)), jj=1,jpj) |
---|
270 | WRITE(numicb,*) "j point" |
---|
271 | WRITE(numicb,*) (INT(src_calving(ji,jj)), jj=1,jpj) |
---|
272 | IF( npolj > 0 ) THEN |
---|
273 | WRITE(numicb,*) 'north fold destination points ' |
---|
274 | WRITE(numicb,*) nicbfldpts |
---|
275 | WRITE(numicb,*) 'north fold destination procs ' |
---|
276 | WRITE(numicb,*) nicbflddest |
---|
277 | ENDIF |
---|
278 | CALL flush(numicb) |
---|
279 | ENDIF |
---|
280 | |
---|
281 | src_calving(:,:) = 0._wp |
---|
282 | src_calving_hflx(:,:) = 0._wp |
---|
283 | |
---|
284 | ! assign each new iceberg with a unique number constructed from the processor number |
---|
285 | ! and incremented by the total number of processors |
---|
286 | num_bergs(:) = 0 |
---|
287 | num_bergs(1) = narea - jpnij |
---|
288 | |
---|
289 | ! when not generating test icebergs we need to setup calving file |
---|
290 | IF( nn_test_icebergs < 0 ) THEN |
---|
291 | ! |
---|
292 | ! maximum distribution class array does not change in time so read it once |
---|
293 | cl_sdist = TRIM( cn_dir )//TRIM( sn_icb%clname ) |
---|
294 | CALL iom_open ( cl_sdist, inum ) ! open file |
---|
295 | ivar = iom_varid( inum, 'maxclass', ldstop=.FALSE. ) |
---|
296 | IF( ivar > 0 ) THEN |
---|
297 | CALL iom_get ( inum, jpdom_data, 'maxclass', src_calving ) ! read the max distribution array |
---|
298 | berg_grid%maxclass(:,:) = INT( src_calving ) |
---|
299 | src_calving(:,:) = 0._wp |
---|
300 | ENDIF |
---|
301 | CALL iom_close( inum ) ! close file |
---|
302 | ! |
---|
303 | WRITE(numicb,*) |
---|
304 | WRITE(numicb,*) ' calving read in a file' |
---|
305 | ALLOCATE( sf_icb(1), STAT=istat1 ) ! Create sf_icb structure (calving) |
---|
306 | ALLOCATE( sf_icb(1)%fnow(jpi,jpj,1), STAT=istat2 ) |
---|
307 | ALLOCATE( sf_icb(1)%fdta(jpi,jpj,1,2), STAT=istat3 ) |
---|
308 | IF( istat1+istat2+istat3 > 0 ) THEN |
---|
309 | CALL ctl_stop( 'sbc_icb: unable to allocate sf_icb structure' ) ; RETURN |
---|
310 | ENDIF |
---|
311 | ! ! fill sf_icb with the namelist (sn_icb) and control print |
---|
312 | CALL fld_fill( sf_icb, (/ sn_icb /), cn_dir, 'icb_init', 'read calving data', 'namicb' ) |
---|
313 | ! |
---|
314 | ENDIF |
---|
315 | |
---|
316 | IF( .NOT.ln_rstart ) THEN |
---|
317 | IF( nn_test_icebergs > 0 ) CALL icb_gen() |
---|
318 | ELSE |
---|
319 | IF( nn_test_icebergs > 0 ) THEN |
---|
320 | CALL icb_gen() |
---|
321 | ELSE |
---|
322 | ! |
---|
323 | CALL icebergs_read_restart() |
---|
324 | l_restarted_bergs = .TRUE. |
---|
325 | ENDIF |
---|
326 | ENDIF |
---|
327 | ! |
---|
328 | IF( nn_sample_rate .GT. 0 ) CALL traj_init( nitend ) |
---|
329 | ! |
---|
330 | CALL icb_budget_init() |
---|
331 | ! |
---|
332 | IF( nn_verbose_level >= 2 ) CALL print_bergs('icb_init, initial status', nit000-1) |
---|
333 | ! |
---|
334 | END SUBROUTINE icb_init |
---|
335 | |
---|
336 | SUBROUTINE icb_gen() |
---|
337 | |
---|
338 | ! Local variables |
---|
339 | INTEGER :: ji, jj, ibergs |
---|
340 | TYPE(iceberg) :: localberg ! NOT a pointer but an actual local variable |
---|
341 | TYPE(point) :: localpt |
---|
342 | INTEGER :: iyr, imon, iday, ihr, imin, isec |
---|
343 | INTEGER :: iberg |
---|
344 | |
---|
345 | ! For convenience |
---|
346 | iberg = nn_test_icebergs |
---|
347 | |
---|
348 | ! call get_date(Time, iyr, imon, iday, ihr, imin, isec) |
---|
349 | ! Convert nemo time variables from dom_oce into local versions |
---|
350 | iyr = nyear |
---|
351 | imon = nmonth |
---|
352 | iday = nday |
---|
353 | ihr = INT(nsec_day/3600) |
---|
354 | imin = INT((nsec_day-ihr*3600)/60) |
---|
355 | isec = nsec_day - ihr*3600 - imin*60 |
---|
356 | |
---|
357 | ! no overlap for icebergs since we want only one instance of each across the whole domain |
---|
358 | ! so restrict area of interest |
---|
359 | ! use tmask here because tmask_i has been doctored on one side of the north fold line |
---|
360 | |
---|
361 | DO jj = nicbdj,nicbej |
---|
362 | DO ji = nicbdi,nicbei |
---|
363 | IF (tmask(ji,jj,1) .GT. 0._wp .AND. & |
---|
364 | gphit(ji,jj) .GT. rn_test_box(3) .AND. gphit(ji,jj) .LT. rn_test_box(4) .AND. & |
---|
365 | glamt(ji,jj) .GT. rn_test_box(1) .AND. glamt(ji,jj) .LT. rn_test_box(2)) THEN |
---|
366 | localberg%mass_scaling = rn_mass_scaling(iberg) |
---|
367 | localpt%xi = REAL( nimpp+ji-1, wp ) |
---|
368 | localpt%yj = REAL( njmpp+jj-1, wp ) |
---|
369 | localpt%lon = bilin(glamt, localpt%xi, localpt%yj, 'T', 0, 0 ) |
---|
370 | localpt%lat = bilin(gphit, localpt%xi, localpt%yj, 'T', 0, 0 ) |
---|
371 | localpt%mass = rn_initial_mass(iberg) |
---|
372 | localpt%thickness = rn_initial_thickness(iberg) |
---|
373 | localpt%width = first_width(iberg) |
---|
374 | localpt%length = first_length(iberg) |
---|
375 | localpt%year = iyr |
---|
376 | localpt%day = FLOAT(iday)+(FLOAT(ihr)+FLOAT(imin)/60._wp)/24._wp |
---|
377 | localpt%mass_of_bits = 0._wp |
---|
378 | localpt%heat_density = 0._wp |
---|
379 | localpt%uvel = 0._wp |
---|
380 | localpt%vvel = 0._wp |
---|
381 | CALL increment_kounter() |
---|
382 | localberg%number(:) = num_bergs(:) |
---|
383 | call add_new_berg_to_list(localberg, localpt) |
---|
384 | ENDIF |
---|
385 | ENDDO |
---|
386 | ENDDO |
---|
387 | |
---|
388 | ibergs = count_bergs() |
---|
389 | IF( lk_mpp ) CALL mpp_sum(ibergs) |
---|
390 | WRITE(numicb,'(a,i6,a)') 'diamonds, icb_gen: ',ibergs,' were generated' |
---|
391 | |
---|
392 | END SUBROUTINE icb_gen |
---|
393 | |
---|
394 | SUBROUTINE icb_nam |
---|
395 | !!---------------------------------------------------------------------- |
---|
396 | !! *** ROUTINE icb_nam *** |
---|
397 | !! |
---|
398 | !! ** Purpose : read domaine namelists and print the variables. |
---|
399 | !! |
---|
400 | !! ** input : - namberg namelist |
---|
401 | !!---------------------------------------------------------------------- |
---|
402 | |
---|
403 | NAMELIST/namberg/ ln_icebergs, ln_bergdia, nn_sample_rate, rn_initial_mass, & |
---|
404 | & rn_distribution, rn_mass_scaling, rn_initial_thickness, nn_verbose_write, & |
---|
405 | & rn_rho_bergs, rn_LoW_ratio, nn_verbose_level, ln_operator_splitting, & |
---|
406 | & rn_bits_erosion_fraction, rn_sicn_shift, ln_passive_mode, & |
---|
407 | & ln_time_average_weight, nn_test_icebergs, rn_test_box, rn_speed_limit, & |
---|
408 | & cn_dir, sn_icb |
---|
409 | INTEGER :: jn ! loop counter |
---|
410 | REAL(wp) :: zfact |
---|
411 | !!---------------------------------------------------------------------- |
---|
412 | |
---|
413 | ! (NB: frequency positive => hours, negative => months) |
---|
414 | ! ! file ! frequency ! variable ! time intep ! clim ! 'yearly' or ! weights ! rotation ! |
---|
415 | ! ! name ! (hours) ! name ! (T/F) ! (T/F) ! 'monthly' ! filename ! pairs ! |
---|
416 | sn_icb = FLD_N( 'calving' , -1 , 'calving' , .TRUE. , .TRUE. , 'yearly' , '' , '' ) |
---|
417 | |
---|
418 | REWIND( numnam ) ! Namelist namrun : iceberg parameters |
---|
419 | READ ( numnam, namberg ) |
---|
420 | |
---|
421 | IF( .NOT. ln_icebergs .AND. lwp ) THEN |
---|
422 | WRITE(numout,*) |
---|
423 | WRITE(numout,*) 'icb_nam : iceberg initialization through namberg namelist read' |
---|
424 | WRITE(numout,*) '~~~~~~~~~~~~~ ' |
---|
425 | WRITE(numout,*) 'NO icebergs used' |
---|
426 | RETURN |
---|
427 | ENDIF |
---|
428 | |
---|
429 | IF( nn_test_icebergs > nclasses ) THEN |
---|
430 | IF( lwp ) WRITE(numout,*) 'Resetting nn_test_icebergs to ',nclasses |
---|
431 | nn_test_icebergs = nclasses |
---|
432 | ENDIF |
---|
433 | |
---|
434 | zfact = SUM( rn_distribution ) |
---|
435 | IF( zfact < 1._wp ) THEN |
---|
436 | IF( zfact <= 0._wp ) THEN |
---|
437 | CALL ctl_stop('icb_init: sum of berg distribution equal to zero') |
---|
438 | ELSE |
---|
439 | rn_distribution(:) = rn_distribution(:) / zfact |
---|
440 | CALL ctl_warn('icb_init: sum of berg input distribution not equal to one and so RESCALED') |
---|
441 | ENDIF |
---|
442 | ENDIF |
---|
443 | |
---|
444 | IF(lwp) THEN ! control print |
---|
445 | WRITE(numout,*) |
---|
446 | WRITE(numout,*) 'icb_nam : iceberg initialization through namberg namelist read' |
---|
447 | WRITE(numout,*) '~~~~~~~~~~~~~ ' |
---|
448 | WRITE(numout,*) ' Calculate budgets ln_bergdia = ', ln_bergdia |
---|
449 | WRITE(numout,*) ' Period between sampling of position for trajectory storage nn_sample_rate = ', nn_sample_rate |
---|
450 | WRITE(numout,*) ' Mass thresholds between iceberg classes (kg) rn_initial_mass =' |
---|
451 | DO jn=1,nclasses |
---|
452 | WRITE(numout,'(a,f15.2)') ' ',rn_initial_mass(jn) |
---|
453 | ENDDO |
---|
454 | WRITE(numout,*) ' Fraction of calving to apply to this class (non-dim) rn_distribution =' |
---|
455 | DO jn=1,nclasses |
---|
456 | WRITE(numout,'(a,f10.2)') ' ',rn_distribution(jn) |
---|
457 | ENDDO |
---|
458 | WRITE(numout,*) ' Ratio between effective and real iceberg mass (non-dim) rn_mass_scaling = ' |
---|
459 | DO jn=1,nclasses |
---|
460 | WRITE(numout,'(a,f10.2)') ' ',rn_mass_scaling(jn) |
---|
461 | ENDDO |
---|
462 | WRITE(numout,*) ' Total thickness of newly calved bergs (m) rn_initial_thickness = ' |
---|
463 | DO jn=1,nclasses |
---|
464 | WRITE(numout,'(a,f10.2)') ' ',rn_initial_thickness(jn) |
---|
465 | ENDDO |
---|
466 | WRITE(numout,*) ' Timesteps between verbose messages nn_verbose_write = ', nn_verbose_write |
---|
467 | |
---|
468 | WRITE(numout,*) ' Density of icebergs rn_rho_bergs = ', rn_rho_bergs |
---|
469 | WRITE(numout,*) ' Initial ratio L/W for newly calved icebergs rn_LoW_ratio = ', rn_LoW_ratio |
---|
470 | WRITE(numout,*) ' Turn on more verbose output level = ', nn_verbose_level |
---|
471 | WRITE(numout,*) ' Use first order operator splitting for thermodynamics ', & |
---|
472 | & 'use_operator_splitting = ', ln_operator_splitting |
---|
473 | WRITE(numout,*) ' Fraction of erosion melt flux to divert to bergy bits ', & |
---|
474 | & 'bits_erosion_fraction = ', rn_bits_erosion_fraction |
---|
475 | |
---|
476 | WRITE(numout,*) ' Shift of sea-ice concentration in erosion flux modulation ', & |
---|
477 | & '(0<sicn_shift<1) sicn_shift = ', rn_sicn_shift |
---|
478 | WRITE(numout,*) ' Do not add freshwater flux from icebergs to ocean ', & |
---|
479 | & ' passive_mode = ', ln_passive_mode |
---|
480 | WRITE(numout,*) ' Time average the weight on the ocean time_average_weight = ', ln_time_average_weight |
---|
481 | WRITE(numout,*) ' Create icebergs in absence of a restart file nn_test_icebergs = ', nn_test_icebergs |
---|
482 | WRITE(numout,*) ' in lon/lat box = ', rn_test_box |
---|
483 | WRITE(numout,*) ' CFL speed limit for a berg speed_limit = ', rn_speed_limit |
---|
484 | WRITE(numout,*) ' Writing Iceberg status information to icebergs.stat file ' |
---|
485 | ENDIF |
---|
486 | ! |
---|
487 | END SUBROUTINE icb_nam |
---|
488 | |
---|
489 | !!====================================================================== |
---|
490 | |
---|
491 | END MODULE icbini |
---|