Line data Source code
1 : !=======================================================================
2 :
3 : ! Diagnostic information output during run
4 : !
5 : ! authors: Elizabeth C. Hunke, LANL
6 :
7 : module icedrv_diagnostics
8 :
9 : use icedrv_kinds
10 : use icedrv_constants, only: nu_diag, nu_diag_out
11 : use icedrv_domain_size, only: nx
12 : use icedrv_domain_size, only: ncat, nfsd, n_iso
13 : use icepack_intfc, only: c0
14 : use icepack_intfc, only: icepack_warnings_flush, icepack_warnings_aborted
15 : use icepack_intfc, only: icepack_query_parameters
16 : use icepack_intfc, only: icepack_query_tracer_flags, icepack_query_tracer_indices
17 : use icedrv_system, only: icedrv_system_abort
18 :
19 : implicit none
20 : private
21 : public :: runtime_diags, &
22 : init_mass_diags, &
23 : icedrv_diagnostics_debug, &
24 : print_state
25 :
26 : ! diagnostic output file
27 : character (len=char_len), public :: diag_file
28 :
29 : ! point print data
30 :
31 : logical (kind=log_kind), public :: &
32 : print_points ! if true, print point data
33 :
34 : integer (kind=int_kind), parameter, public :: &
35 : npnt = 2 ! total number of points to be printed
36 :
37 : character (len=char_len), dimension(nx), public :: nx_names
38 :
39 : ! for water and heat budgets
40 : real (kind=dbl_kind), dimension(nx) :: &
41 : pdhi , & ! change in mean ice thickness (m)
42 : pdhs , & ! change in mean snow thickness (m)
43 : pde ! change in ice and snow energy (W m-2)
44 :
45 : real (kind=dbl_kind), dimension(nx,n_iso) :: &
46 : pdiso ! change in mean isotope concentration
47 :
48 : !=======================================================================
49 :
50 : contains
51 :
52 : !=======================================================================
53 :
54 : ! Writes diagnostic info (max, min, global sums, etc) to standard out
55 : !
56 : ! authors: Elizabeth C. Hunke, LANL
57 : ! Bruce P. Briegleb, NCAR
58 : ! Cecilia M. Bitz, UW
59 :
60 37392 : subroutine runtime_diags (dt)
61 :
62 : use icedrv_arrays_column, only: floe_rad_c
63 : use icedrv_flux, only: evap, fsnow, frazil
64 : use icedrv_flux, only: fswabs, flw, flwout, fsens, fsurf, flat
65 : use icedrv_flux, only: frain, fiso_evap, fiso_ocn, fiso_atm
66 : use icedrv_flux, only: Tair, Qa, fsw, fcondtop
67 : use icedrv_flux, only: meltt, meltb, meltl, snoice
68 : use icedrv_flux, only: dsnow, congel, sst, sss, Tf, fhocn
69 : use icedrv_state, only: aice, vice, vsno, trcr, trcrn, aicen
70 :
71 : real (kind=dbl_kind), intent(in) :: &
72 : dt ! time step
73 :
74 : ! local variables
75 :
76 : integer (kind=int_kind) :: &
77 : n, nc, k
78 :
79 : logical (kind=log_kind) :: &
80 : calc_Tsfc, tr_fsd, tr_iso
81 :
82 : ! fields at diagnostic points
83 : real (kind=dbl_kind) :: &
84 9880 : pTair, pfsnow, pfrain, &
85 14820 : paice, hiavg, hsavg, hbravg, psalt, pTsfc, &
86 9880 : pevap, pfhocn, fsdavg
87 :
88 : real (kind=dbl_kind), dimension (nx) :: &
89 64220 : work1, work2, work3
90 :
91 : real (kind=dbl_kind) :: &
92 4940 : Tffresh, rhos, rhow, rhoi
93 :
94 : logical (kind=log_kind) :: tr_brine
95 : integer (kind=int_kind) :: nt_fbri, nt_Tsfc, nt_fsd, nt_isosno, nt_isoice
96 :
97 : character(len=*), parameter :: subname='(runtime_diags)'
98 :
99 : !-----------------------------------------------------------------
100 : ! query Icepack values
101 : !-----------------------------------------------------------------
102 :
103 37392 : call icepack_query_parameters(calc_Tsfc_out=calc_Tsfc)
104 : call icepack_query_tracer_flags(tr_brine_out=tr_brine, &
105 37392 : tr_fsd_out=tr_fsd,tr_iso_out=tr_iso)
106 : call icepack_query_tracer_indices(nt_fbri_out=nt_fbri, nt_Tsfc_out=nt_Tsfc,&
107 37392 : nt_fsd_out=nt_fsd, nt_isosno_out=nt_isosno, nt_isoice_out=nt_isoice)
108 : call icepack_query_parameters(Tffresh_out=Tffresh, rhos_out=rhos, &
109 37392 : rhow_out=rhow, rhoi_out=rhoi)
110 37392 : call icepack_warnings_flush(nu_diag)
111 37392 : if (icepack_warnings_aborted()) call icedrv_system_abort(string=subname, &
112 0 : file=__FILE__,line= __LINE__)
113 :
114 : !-----------------------------------------------------------------
115 : ! NOTE these are computed for the last timestep only (not avg)
116 : !-----------------------------------------------------------------
117 :
118 37392 : call total_energy (work1)
119 37392 : call total_salt (work2)
120 :
121 186960 : do n = 1, nx
122 149568 : pTair = Tair(n) - Tffresh ! air temperature
123 149568 : pfsnow = fsnow(n)*dt/rhos ! snowfall
124 149568 : pfrain = frain(n)*dt/rhow ! rainfall
125 :
126 149568 : paice = aice(n) ! ice area
127 149568 : hiavg = c0 ! avg snow/ice thickness
128 149568 : fsdavg = c0 ! FSD rep radius
129 149568 : hsavg = c0
130 149568 : hbravg = c0 ! avg brine thickness
131 149568 : psalt = c0
132 149568 : if (paice /= c0) then
133 110929 : hiavg = vice(n)/paice
134 110929 : hsavg = vsno(n)/paice
135 110929 : if (tr_brine) hbravg = trcr(n,nt_fbri)* hiavg
136 110929 : if (tr_fsd) then
137 23760 : do nc = 1, ncat
138 201135 : do k = 1, nfsd
139 : fsdavg = fsdavg &
140 136350 : + trcrn(n,nt_fsd+k-1,nc) * floe_rad_c(k) &
141 265350 : * aicen(n,nc) / paice
142 : end do
143 : end do
144 : end if
145 :
146 : endif
147 149568 : if (vice(n) /= c0) psalt = work2(n)/vice(n)
148 149568 : pTsfc = trcr(n,nt_Tsfc) ! ice/snow sfc temperature
149 149568 : pevap = evap(n)*dt/rhoi ! sublimation/condensation
150 149568 : pdhi(n) = vice(n) - pdhi(n) ! ice thickness change
151 149568 : pdhs(n) = vsno(n) - pdhs(n) ! snow thickness change
152 149568 : pde(n) =-(work1(n)- pde(n))/dt ! ice/snow energy change
153 149568 : pfhocn = -fhocn(n) ! ocean heat used by ice
154 :
155 149568 : work3(:) = c0
156 :
157 157248 : do k = 1, n_iso
158 0 : work3 (n) = (trcr(n,nt_isosno+k-1)*vsno(n) &
159 7680 : +trcr(n,nt_isoice+k-1)*vice(n))
160 157248 : pdiso(n,k) = work3(n) - pdiso(n,k)
161 : enddo
162 :
163 : !-----------------------------------------------------------------
164 : ! start spewing
165 : !-----------------------------------------------------------------
166 :
167 149568 : write(nu_diag_out+n-1,899) nx_names(n)
168 :
169 149568 : write(nu_diag_out+n-1,*) ' '
170 149568 : write(nu_diag_out+n-1,*) '----------atm----------'
171 149568 : write(nu_diag_out+n-1,900) 'air temperature (C) = ',pTair
172 149568 : write(nu_diag_out+n-1,900) 'specific humidity = ',Qa(n)
173 149568 : write(nu_diag_out+n-1,900) 'snowfall (m) = ',pfsnow
174 149568 : write(nu_diag_out+n-1,900) 'rainfall (m) = ',pfrain
175 149568 : if (.not.calc_Tsfc) then
176 0 : write(nu_diag_out+n-1,900) 'total surface heat flux= ', fsurf(n)
177 0 : write(nu_diag_out+n-1,900) 'top sfc conductive flux= ',fcondtop(n)
178 0 : write(nu_diag_out+n-1,900) 'latent heat flux = ',flat(n)
179 : else
180 149568 : write(nu_diag_out+n-1,900) 'shortwave radiation sum= ',fsw(n)
181 149568 : write(nu_diag_out+n-1,900) 'longwave radiation = ',flw(n)
182 : endif
183 149568 : write(nu_diag_out+n-1,*) '----------ice----------'
184 149568 : write(nu_diag_out+n-1,900) 'area fraction = ',aice(n)! ice area
185 149568 : write(nu_diag_out+n-1,900) 'avg ice thickness (m) = ',hiavg
186 149568 : write(nu_diag_out+n-1,900) 'avg snow depth (m) = ',hsavg
187 149568 : write(nu_diag_out+n-1,900) 'avg salinity (ppt) = ',psalt
188 149568 : write(nu_diag_out+n-1,900) 'avg brine thickness (m)= ',hbravg
189 149568 : if (tr_fsd) &
190 5480 : write(nu_diag_out+n-1,900) 'avg fsd rep radius (m) = ',fsdavg
191 :
192 :
193 149568 : if (calc_Tsfc) then
194 149568 : write(nu_diag_out+n-1,900) 'surface temperature(C) = ',pTsfc ! ice/snow
195 149568 : write(nu_diag_out+n-1,900) 'absorbed shortwave flx = ',fswabs(n)
196 149568 : write(nu_diag_out+n-1,900) 'outward longwave flx = ',flwout(n)
197 149568 : write(nu_diag_out+n-1,900) 'sensible heat flx = ',fsens(n)
198 149568 : write(nu_diag_out+n-1,900) 'latent heat flx = ',flat(n)
199 : endif
200 149568 : write(nu_diag_out+n-1,900) 'subl/cond (m ice) = ',pevap ! sublimation/condensation
201 149568 : write(nu_diag_out+n-1,900) 'top melt (m) = ',meltt(n)
202 149568 : write(nu_diag_out+n-1,900) 'bottom melt (m) = ',meltb(n)
203 149568 : write(nu_diag_out+n-1,900) 'lateral melt (m) = ',meltl(n)
204 149568 : write(nu_diag_out+n-1,900) 'new ice (m) = ',frazil(n) ! frazil
205 149568 : write(nu_diag_out+n-1,900) 'congelation (m) = ',congel(n)
206 149568 : write(nu_diag_out+n-1,900) 'snow-ice (m) = ',snoice(n)
207 149568 : write(nu_diag_out+n-1,900) 'snow change (m) = ',dsnow(n)
208 149568 : write(nu_diag_out+n-1,900) 'effective dhi (m) = ',pdhi(n) ! ice thickness change
209 149568 : write(nu_diag_out+n-1,900) 'effective dhs (m) = ',pdhs(n) ! snow thickness change
210 149568 : write(nu_diag_out+n-1,900) 'intnl enrgy chng(W/m^2)= ',pde (n) ! ice/snow energy change
211 149568 : write(nu_diag_out+n-1,*) '----------ocn----------'
212 149568 : write(nu_diag_out+n-1,900) 'sst (C) = ',sst(n) ! sea surface temperature
213 149568 : write(nu_diag_out+n-1,900) 'sss (ppt) = ',sss(n) ! sea surface salinity
214 149568 : write(nu_diag_out+n-1,900) 'freezing temp (C) = ',Tf(n) ! freezing temperature
215 149568 : write(nu_diag_out+n-1,900) 'heat used (W/m^2) = ',pfhocn ! ocean heat used by ice
216 :
217 186960 : if (tr_iso) then
218 10240 : do k = 1, n_iso
219 7680 : write(nu_diag_out+n-1,901) 'isotopic precip = ',fiso_atm(n,k)*dt,k
220 7680 : write(nu_diag_out+n-1,901) 'isotopic evap/cond = ',fiso_evap(n,k)*dt,k
221 7680 : write(nu_diag_out+n-1,901) 'isotopic loss to ocn = ',fiso_ocn(n,k)*dt,k
222 7680 : write(nu_diag_out+n-1,901) 'isotopic gain/loss = ',(fiso_atm(n,k)-fiso_ocn(n,k)+fiso_evap(n,k))*dt,k
223 10240 : write(nu_diag_out+n-1,901) 'isotopic conc chg = ',pdiso(n,k),k
224 : enddo
225 : endif
226 : end do
227 : 899 format (43x,a24)
228 : 900 format (a25,2x,f24.17)
229 : 901 format (a25,2x,f24.17,i6)
230 :
231 37392 : end subroutine runtime_diags
232 :
233 : !=======================================================================
234 :
235 : ! Computes global combined ice and snow mass sum
236 : !
237 : ! author: Elizabeth C. Hunke, LANL
238 :
239 342540 : subroutine init_mass_diags
240 :
241 : use icedrv_state, only: vice, vsno, trcr
242 :
243 : integer (kind=int_kind) :: i, k, nt_isosno, nt_isoice
244 :
245 592860 : real (kind=dbl_kind), dimension (nx) :: work1
246 :
247 : character(len=*), parameter :: subname='(init_mass_diags)'
248 :
249 342540 : call icepack_query_tracer_indices(nt_isosno_out=nt_isosno)
250 342540 : call icepack_query_tracer_indices(nt_isoice_out=nt_isoice)
251 :
252 342540 : call total_energy (work1)
253 1712700 : do i = 1, nx
254 1370160 : pdhi(i) = vice (i)
255 1370160 : pdhs(i) = vsno (i)
256 1370160 : pde (i) = work1(i)
257 1897164 : do k = 1, n_iso
258 0 : pdiso(i,k) = (trcr(i,nt_isosno+k-1)*vsno(i) &
259 1554624 : +trcr(i,nt_isoice+k-1)*vice(i))
260 : enddo
261 : enddo
262 :
263 342540 : end subroutine init_mass_diags
264 :
265 : !=======================================================================
266 :
267 : ! Computes total energy of ice and snow in a grid cell.
268 : !
269 : ! authors: E. C. Hunke, LANL
270 :
271 379932 : subroutine total_energy (work)
272 :
273 : use icedrv_domain_size, only: ncat, nilyr, nslyr
274 : use icedrv_state, only: vicen, vsnon, trcrn
275 :
276 : real (kind=dbl_kind), dimension (nx), intent(out) :: &
277 : work ! total energy
278 :
279 : ! local variables
280 :
281 : integer (kind=int_kind) :: &
282 : i, k, n
283 :
284 : integer (kind=int_kind) :: nt_qice, nt_qsno
285 :
286 : character(len=*), parameter :: subname='(total_energy)'
287 :
288 : !-----------------------------------------------------------------
289 : ! query Icepack values
290 : !-----------------------------------------------------------------
291 :
292 379932 : call icepack_query_tracer_indices(nt_qice_out=nt_qice, nt_qsno_out=nt_qsno)
293 379932 : call icepack_warnings_flush(nu_diag)
294 379932 : if (icepack_warnings_aborted()) call icedrv_system_abort(string=subname, &
295 0 : file=__FILE__,line= __LINE__)
296 :
297 : !-----------------------------------------------------------------
298 : ! Initialize
299 : !-----------------------------------------------------------------
300 :
301 379932 : work(:) = c0
302 :
303 : !-----------------------------------------------------------------
304 : ! Aggregate
305 : !-----------------------------------------------------------------
306 :
307 2179044 : do n = 1, ncat
308 13487964 : do k = 1, nilyr
309 60243372 : do i = 1, nx
310 14955680 : work(i) = work(i) &
311 29911360 : + trcrn(i,nt_qice+k-1,n) &
312 73399940 : * vicen(i,n) / real(nilyr,kind=dbl_kind)
313 : enddo ! i
314 : enddo ! k
315 :
316 4480896 : do k = 1, nslyr
317 13308372 : do i = 1, nx
318 3054240 : work(i) = work(i) &
319 6108480 : + trcrn(i,nt_qsno+k-1,n) &
320 14563500 : * vsnon(i,n) / real(nslyr,kind=dbl_kind)
321 : enddo ! i
322 : enddo ! k
323 : enddo ! n
324 :
325 379932 : end subroutine total_energy
326 :
327 : !=======================================================================
328 :
329 : ! Computes bulk salinity of ice and snow in a grid cell.
330 : ! author: E. C. Hunke, LANL
331 :
332 37392 : subroutine total_salt (work)
333 :
334 : use icedrv_domain_size, only: ncat, nilyr
335 : use icedrv_state, only: vicen, trcrn
336 :
337 : real (kind=dbl_kind), dimension (nx), &
338 : intent(out) :: &
339 : work ! total salt
340 :
341 : ! local variables
342 :
343 : integer (kind=int_kind) :: &
344 : i, k, n
345 :
346 : integer (kind=int_kind) :: nt_sice
347 :
348 : character(len=*), parameter :: subname='(total_salt)'
349 :
350 : !-----------------------------------------------------------------
351 : ! query Icepack values
352 : !-----------------------------------------------------------------
353 :
354 37392 : call icepack_query_tracer_indices(nt_sice_out=nt_sice)
355 37392 : call icepack_warnings_flush(nu_diag)
356 37392 : if (icepack_warnings_aborted()) call icedrv_system_abort(string=subname, &
357 0 : file=__FILE__,line= __LINE__)
358 :
359 : !-----------------------------------------------------------------
360 : ! Initialize
361 : !-----------------------------------------------------------------
362 :
363 37392 : work(:) = c0
364 :
365 : !-----------------------------------------------------------------
366 : ! Aggregate
367 : !-----------------------------------------------------------------
368 :
369 220332 : do n = 1, ncat
370 1464732 : do k = 1, nilyr
371 6404940 : do i = 1, nx
372 598160 : work(i) = work(i) &
373 1196320 : + trcrn(i,nt_sice+k-1,n) &
374 6820160 : * vicen(i,n) / real(nilyr,kind=dbl_kind)
375 : enddo ! i
376 : enddo ! k
377 : enddo ! n
378 :
379 37392 : end subroutine total_salt
380 :
381 : !=======================================================================
382 : !
383 : ! Wrapper for the print_state debugging routine.
384 : ! Useful for debugging in the main driver (see ice.F_debug)
385 : !
386 : ! author Elizabeth C. Hunke, LANL
387 : !
388 0 : subroutine icedrv_diagnostics_debug (plabeld)
389 :
390 : use icedrv_calendar, only: istep1
391 :
392 : character (*), intent(in) :: plabeld
393 :
394 : character(len=*), parameter :: &
395 : subname='(icedrv_diagnostics_debug)'
396 :
397 : ! printing info for routine print_state
398 :
399 : integer (kind=int_kind), parameter :: &
400 : check_step = 1439, & ! begin printing at istep1=check_step
401 : ip = 3 ! i index
402 :
403 0 : if (istep1 >= check_step) then
404 0 : call print_state(plabeld,ip)
405 : endif
406 :
407 0 : end subroutine icedrv_diagnostics_debug
408 :
409 : !=======================================================================
410 :
411 : ! This routine is useful for debugging.
412 : ! Calls to it should be inserted in the form (after thermo, for example)
413 : ! plabel = 'post thermo'
414 : ! if (istep1 >= check_step) call print_state(plabel,ip)
415 : ! 'use ice_diagnostics' may need to be inserted also
416 : ! author: Elizabeth C. Hunke, LANL
417 :
418 0 : subroutine print_state(plabel,i)
419 :
420 : use icedrv_calendar, only: istep1, time
421 : use icedrv_domain_size, only: nilyr, nslyr
422 : use icedrv_state, only: aice0, aicen, vicen, vsnon, uvel, vvel, trcrn
423 : use icedrv_flux, only: uatm, vatm, potT, Tair, Qa, flw, frain, fsnow
424 : use icedrv_flux, only: fsens, flat, evap, flwout
425 : use icedrv_flux, only: swvdr, swvdf, swidr, swidf, rhoa
426 : use icedrv_flux, only: frzmlt, sst, sss, Tf, Tref, Qref, Uref
427 : use icedrv_flux, only: uocn, vocn
428 : use icedrv_flux, only: fsw, fswabs, fswint_ai, fswthru, scale_factor
429 : use icedrv_flux, only: alvdr_ai, alvdf_ai, alidf_ai, alidr_ai
430 :
431 : character (*), intent(in) :: plabel
432 :
433 : integer (kind=int_kind), intent(in) :: &
434 : i ! horizontal index
435 :
436 : ! local variables
437 :
438 : real (kind=dbl_kind) :: &
439 0 : eidebug, esdebug, &
440 0 : qi, qs, Tsnow, &
441 0 : puny, Lfresh, cp_ice, &
442 0 : rhoi, rhos
443 :
444 : integer (kind=int_kind) :: n, k
445 :
446 : integer (kind=int_kind) :: nt_Tsfc, nt_qice, nt_qsno, nt_fsd
447 :
448 : logical (kind=log_kind) :: tr_fsd
449 :
450 : character(len=*), parameter :: subname='(print_state)'
451 :
452 : !-----------------------------------------------------------------
453 : ! query Icepack values
454 : !-----------------------------------------------------------------
455 :
456 0 : call icepack_query_tracer_flags(tr_fsd_out=tr_fsd)
457 : call icepack_query_tracer_indices(nt_Tsfc_out=nt_Tsfc, nt_qice_out=nt_qice, &
458 0 : nt_qsno_out=nt_qsno,nt_fsd_out=nt_fsd)
459 : call icepack_query_parameters(puny_out=puny, Lfresh_out=Lfresh, cp_ice_out=cp_ice, &
460 0 : rhoi_out=rhoi, rhos_out=rhos)
461 0 : call icepack_warnings_flush(nu_diag)
462 0 : if (icepack_warnings_aborted()) call icedrv_system_abort(string=subname, &
463 0 : file=__FILE__,line= __LINE__)
464 :
465 : !-----------------------------------------------------------------
466 : ! write diagnostics
467 : !-----------------------------------------------------------------
468 :
469 0 : write(nu_diag,*) trim(plabel)
470 0 : write(nu_diag,*) 'istep1, i, time', &
471 0 : istep1, i, time
472 0 : write(nu_diag,*) ' '
473 0 : write(nu_diag,*) 'aice0', aice0(i)
474 0 : do n = 1, ncat
475 0 : write(nu_diag,*) ' '
476 0 : write(nu_diag,*) 'n =',n
477 0 : write(nu_diag,*) 'aicen', aicen(i,n)
478 0 : write(nu_diag,*) 'vicen', vicen(i,n)
479 0 : write(nu_diag,*) 'vsnon', vsnon(i,n)
480 0 : if (aicen(i,n) > puny) then
481 0 : write(nu_diag,*) 'hin', vicen(i,n)/aicen(i,n)
482 0 : write(nu_diag,*) 'hsn', vsnon(i,n)/aicen(i,n)
483 : endif
484 0 : write(nu_diag,*) 'Tsfcn',trcrn(i,nt_Tsfc,n)
485 0 : if (tr_fsd) write(nu_diag,*) 'afsdn',trcrn(i,nt_fsd,n) ! fsd cat 1
486 0 : write(nu_diag,*) ' '
487 : enddo ! n
488 :
489 0 : eidebug = c0
490 0 : do n = 1,ncat
491 0 : do k = 1,nilyr
492 0 : qi = trcrn(i,nt_qice+k-1,n)
493 0 : write(nu_diag,*) 'qice, cat ',n,' layer ',k, qi
494 0 : eidebug = eidebug + qi
495 0 : if (aicen(i,n) > puny) then
496 0 : write(nu_diag,*) 'qi/rhoi', qi/rhoi
497 : endif
498 : enddo
499 0 : write(nu_diag,*) ' '
500 : enddo
501 0 : write(nu_diag,*) 'qice(i)',eidebug
502 0 : write(nu_diag,*) ' '
503 :
504 0 : esdebug = c0
505 0 : do n = 1,ncat
506 0 : if (vsnon(i,n) > puny) then
507 0 : do k = 1,nslyr
508 0 : qs = trcrn(i,nt_qsno+k-1,n)
509 0 : write(nu_diag,*) 'qsnow, cat ',n,' layer ',k, qs
510 0 : esdebug = esdebug + qs
511 0 : Tsnow = (Lfresh + qs/rhos) / cp_ice
512 0 : write(nu_diag,*) 'qs/rhos', qs/rhos
513 0 : write(nu_diag,*) 'Tsnow', Tsnow
514 : enddo
515 0 : write(nu_diag,*) ' '
516 : endif
517 : enddo
518 0 : write(nu_diag,*) 'qsnow(i)',esdebug
519 0 : write(nu_diag,*) ' '
520 :
521 0 : write(nu_diag,*) 'uvel(i)',uvel(i)
522 0 : write(nu_diag,*) 'vvel(i)',vvel(i)
523 :
524 0 : write(nu_diag,*) ' '
525 0 : write(nu_diag,*) 'atm states and fluxes'
526 0 : write(nu_diag,*) ' uatm = ',uatm (i)
527 0 : write(nu_diag,*) ' vatm = ',vatm (i)
528 0 : write(nu_diag,*) ' potT = ',potT (i)
529 0 : write(nu_diag,*) ' Tair = ',Tair (i)
530 0 : write(nu_diag,*) ' Qa = ',Qa (i)
531 0 : write(nu_diag,*) ' rhoa = ',rhoa (i)
532 0 : write(nu_diag,*) ' swvdr = ',swvdr(i)
533 0 : write(nu_diag,*) ' swvdf = ',swvdf(i)
534 0 : write(nu_diag,*) ' swidr = ',swidr(i)
535 0 : write(nu_diag,*) ' swidf = ',swidf(i)
536 0 : write(nu_diag,*) ' flw = ',flw (i)
537 0 : write(nu_diag,*) ' frain = ',frain(i)
538 0 : write(nu_diag,*) ' fsnow = ',fsnow(i)
539 0 : write(nu_diag,*) ' '
540 0 : write(nu_diag,*) 'ocn states and fluxes'
541 0 : write(nu_diag,*) ' frzmlt = ',frzmlt (i)
542 0 : write(nu_diag,*) ' sst = ',sst (i)
543 0 : write(nu_diag,*) ' sss = ',sss (i)
544 0 : write(nu_diag,*) ' Tf = ',Tf (i)
545 0 : write(nu_diag,*) ' uocn = ',uocn (i)
546 0 : write(nu_diag,*) ' vocn = ',vocn (i)
547 0 : write(nu_diag,*) ' '
548 0 : write(nu_diag,*) 'srf states and fluxes'
549 0 : write(nu_diag,*) ' Tref = ',Tref (i)
550 0 : write(nu_diag,*) ' Qref = ',Qref (i)
551 0 : write(nu_diag,*) ' Uref = ',Uref (i)
552 0 : write(nu_diag,*) ' fsens = ',fsens (i)
553 0 : write(nu_diag,*) ' flat = ',flat (i)
554 0 : write(nu_diag,*) ' evap = ',evap (i)
555 0 : write(nu_diag,*) ' flwout = ',flwout(i)
556 0 : write(nu_diag,*) ' '
557 0 : write(nu_diag,*) 'shortwave'
558 0 : write(nu_diag,*) ' fsw = ',fsw (i)
559 0 : write(nu_diag,*) ' fswabs = ',fswabs (i)
560 0 : write(nu_diag,*) ' fswint_ai = ',fswint_ai (i)
561 0 : write(nu_diag,*) ' fswthru = ',fswthru (i)
562 0 : write(nu_diag,*) ' scale_factor = ',scale_factor(i)
563 0 : write(nu_diag,*) ' alvdr = ',alvdr_ai (i)
564 0 : write(nu_diag,*) ' alvdf = ',alvdf_ai (i)
565 0 : write(nu_diag,*) ' alidr = ',alidr_ai (i)
566 0 : write(nu_diag,*) ' alidf = ',alidf_ai (i)
567 0 : write(nu_diag,*) ' '
568 :
569 0 : call icepack_warnings_flush(nu_diag)
570 :
571 0 : end subroutine print_state
572 :
573 : !=======================================================================
574 :
575 : end module icedrv_diagnostics
576 :
577 : !=======================================================================
|