00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 module ice_diagnostics
00022
00023
00024
00025 use ice_kinds_mod
00026 use ice_communicate, only: my_task, master_task
00027 use ice_constants
00028 use ice_calendar, only: diagfreq, istep1, istep
00029 use ice_fileunits
00030 use ice_domain_size
00031
00032
00033
00034 implicit none
00035 save
00036
00037
00038 character (len=char_len) :: diag_file
00039
00040
00041
00042 logical (kind=log_kind) ::
00043 print_points ,
00044 print_global
00045
00046 integer (kind=int_kind), parameter ::
00047 npnt = 2
00048
00049
00050 logical (kind=log_kind), parameter ::
00051 check_umax = .false.
00052
00053 real (kind=dbl_kind), parameter ::
00054 umax_stab = 1.0_dbl_kind ,
00055 aice_extmin = 0.15_dbl_kind
00056
00057 real (kind=dbl_kind), dimension(npnt) ::
00058 latpnt ,
00059 lonpnt
00060
00061 integer (kind=int_kind) ::
00062 iindx ,
00063 jindx ,
00064 bindx
00065
00066
00067 real (kind=dbl_kind), dimension(npnt) ::
00068 pdhi ,
00069 pdhs ,
00070 pde ,
00071 plat, plon
00072
00073 integer (kind=int_kind), dimension(npnt) ::
00074 piloc, pjloc, pbloc, pmloc
00075
00076
00077 real (kind=dbl_kind) ::
00078 totmn ,
00079 totms ,
00080 totmin ,
00081 totmis ,
00082 toten ,
00083 totes
00084 real (kind=dbl_kind), dimension(n_aeromx) ::
00085 totaeron ,
00086 totaeros
00087
00088
00089
00090 character (char_len) :: plabel
00091 integer (kind=int_kind), parameter ::
00092 check_step = 999999999,
00093 iblkp = 1,
00094 ip = 3,
00095 jp = 5,
00096 mtask = 0
00097
00098
00099
00100 contains
00101
00102
00103
00104
00105
00106
00107
00108
00109 subroutine runtime_diags (dt)
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123 use ice_broadcast
00124 use ice_global_reductions
00125 use ice_blocks
00126 use ice_domain
00127
00128 use ice_flux
00129 use ice_state
00130 use ice_grid, only: lmask_n, lmask_s, tarean, tareas, grid_type
00131 use ice_therm_vertical, only: calc_Tsfc
00132
00133 #ifdef CCSMCOUPLED
00134 use ice_prescribed_mod, only : prescribed_ice
00135 #endif
00136
00137
00138
00139 real (kind=dbl_kind), intent(in) ::
00140 dt
00141
00142
00143
00144 integer (kind=int_kind) ::
00145 i, j, k, n, ii,jj, iblk
00146
00147
00148 real (kind=dbl_kind) ::
00149 umaxn, hmaxn, shmaxn, arean, snwmxn, extentn,
00150 umaxs, hmaxs, shmaxs, areas, snwmxs, extents,
00151 etotn, mtotn, micen, msnwn, pmaxn, ketotn,
00152 etots, mtots, mices, msnws, pmaxs, ketots,
00153 urmsn, albtotn, arean_alb,
00154 urmss, albtots, areas_alb
00155
00156
00157 real (kind=dbl_kind) ::
00158 rnn, snn, frzn, hnetn, fhocnn, fhatmn, fhfrzn,
00159 rns, sns, frzs, hnets, fhocns, fhatms, fhfrzs,
00160 sfsaltn, sfreshn, evpn, fluxn , delmxn, delmin,
00161 sfsalts, sfreshs, evps, fluxs , delmxs, delmis,
00162 delein, werrn, herrn, msltn, delmsltn, serrn,
00163 deleis, werrs, herrs, mslts, delmslts, serrs,
00164 ftmp,faeron,faeros,fsootn,fsoots
00165
00166
00167 integer (kind=int_kind) ::
00168 kaero, naero
00169 real (kind=dbl_kind) ::
00170 aeromx1n, aeromx1s, aeromx2n, aeromx2s,
00171 aeromx3n, aeromx3s, aoermx4,
00172 aerototn, aerotots
00173
00174
00175 real (kind=dbl_kind), dimension(npnt) ::
00176 paice, pTair, pQa, pfsnow, pfrain, pfsw, pflw,
00177 pTsfc, pevap, pfswabs, pflwout, pflat, pfsens,
00178 pfsurf, pfcondtop, psst, pTf, hiavg, hsavg, pfhocn,
00179 pmeltt, pmeltb, pmeltl, psnoice, pfrazil, pcongel
00180
00181 real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) ::
00182 work1, work2
00183
00184
00185
00186
00187
00188
00189
00190 arean = global_sum(aice, distrb_info, field_loc_center, tarean)
00191 areas = global_sum(aice, distrb_info, field_loc_center, tareas)
00192 arean = arean * m2_to_km2
00193 areas = areas * m2_to_km2
00194
00195
00196 work1(:,:,:) = c0
00197
00198
00199 do iblk = 1, nblocks
00200 do j = 1, ny_block
00201 do i = 1, nx_block
00202 if (aice(i,j,iblk) >= aice_extmin) work1(i,j,iblk) = c1
00203 enddo
00204 enddo
00205 enddo
00206
00207 extentn = global_sum(work1, distrb_info, field_loc_center, &
00208 tarean)
00209 extents = global_sum(work1, distrb_info, field_loc_center, &
00210 tareas)
00211 extentn = extentn * m2_to_km2
00212 extents = extents * m2_to_km2
00213
00214
00215 shmaxn = global_sum(vice, distrb_info, field_loc_center, tarean)
00216 shmaxs = global_sum(vice, distrb_info, field_loc_center, tareas)
00217
00218
00219 snwmxn = global_sum(vsno, distrb_info, field_loc_center, tarean)
00220 snwmxs = global_sum(vsno, distrb_info, field_loc_center, tareas)
00221
00222
00223
00224 do iblk = 1, nblocks
00225 do j = 1, ny_block
00226 do i = 1, nx_block
00227 work1(i,j,iblk) = p5 &
00228 * (rhos*vsno(i,j,iblk) + rhoi*vice(i,j,iblk)) &
00229 * (uvel(i,j,iblk)**2 + vvel(i,j,iblk)**2)
00230 enddo
00231 enddo
00232 enddo
00233
00234 ketotn = global_sum(work1, distrb_info, field_loc_center, tarean)
00235 ketots = global_sum(work1, distrb_info, field_loc_center, tareas)
00236
00237
00238 urmsn = c2*ketotn/(rhoi*shmaxn + rhos*snwmxn + puny)
00239 if (urmsn > puny) then
00240 urmsn = sqrt(urmsn)
00241 else
00242 urmsn = c0
00243 endif
00244
00245 urmss = c2*ketots/(rhoi*shmaxs + rhos*snwmxs + puny)
00246 if (urmss > puny) then
00247 urmss = sqrt(urmss)
00248 else
00249 urmss = c0
00250 endif
00251
00252
00253
00254
00255
00256 do iblk = 1, nblocks
00257 do j = 1, ny_block
00258 do i = 1, nx_block
00259 work1(i,j,iblk) = alvdr(i,j,iblk)*awtvdr &
00260 + alidr(i,j,iblk)*awtidr &
00261 + alvdf(i,j,iblk)*awtvdf &
00262 + alidf(i,j,iblk)*awtidf
00263 enddo
00264 enddo
00265 enddo
00266
00267
00268
00269 do iblk = 1, nblocks
00270 do j = 1, ny_block
00271 do i = 1, nx_block
00272 if (coszen(i,j,iblk) > puny) then
00273 work2(i,j,iblk) = tarean(i,j,iblk)
00274 else
00275 work2(i,j,iblk) = c0
00276 endif
00277 enddo
00278 enddo
00279 enddo
00280
00281
00282 arean_alb = global_sum(aice, distrb_info, field_loc_center, work2)
00283
00284 albtotn = global_sum_prod(aice, work1, distrb_info, &
00285 field_loc_center, work2)
00286
00287 if (arean_alb > c0) then
00288 albtotn = albtotn / arean_alb
00289 else
00290 albtotn = c0
00291 endif
00292
00293
00294 do iblk = 1, nblocks
00295 do j = 1, ny_block
00296 do i = 1, nx_block
00297 if (coszen(i,j,iblk) > puny) then
00298 work2(i,j,iblk) = tareas(i,j,iblk)
00299 else
00300 work2(i,j,iblk) = c0
00301 endif
00302 enddo
00303 enddo
00304 enddo
00305
00306
00307 areas_alb = global_sum(aice, distrb_info, field_loc_center, work2)
00308
00309 albtots = global_sum_prod(aice, work1, distrb_info, &
00310 field_loc_center, work2)
00311
00312 if (areas_alb > c0) then
00313 albtots = albtots / areas_alb
00314 else
00315 albtots = c0
00316 endif
00317
00318
00319 hmaxn = global_maxval(vice, distrb_info, lmask_n)
00320 hmaxs = global_maxval(vice, distrb_info, lmask_s)
00321
00322
00323 if (tr_aero) then
00324
00325 do naero=1,n_aero
00326 faeron = global_sum_prod(faero(:,:,naero,:), aice_init, distrb_info, &
00327 field_loc_center, tarean)
00328 faeros = global_sum_prod(faero(:,:,naero,:), aice_init, distrb_info, &
00329 field_loc_center, tareas)
00330 faeron = faeron*dt
00331 faeros = faeros*dt
00332
00333 fsootn = global_sum_prod(fsoot(:,:,naero,:), aice, distrb_info, &
00334 field_loc_center, tarean)
00335 fsoots = global_sum_prod(fsoot(:,:,naero,:), aice, distrb_info, &
00336 field_loc_center, tareas)
00337 fsootn = fsootn*dt
00338 fsoots = fsoots*dt
00339
00340 do iblk = 1, nblocks
00341 do j = 1, ny_block
00342 do i = 1, nx_block
00343 work1(i,j,iblk) = trcr(i,j,nt_aero +4*(naero-1),iblk) *vsno(i,j,iblk) &
00344 + trcr(i,j,nt_aero+1+4*(naero-1),iblk)*vsno(i,j,iblk) &
00345 + trcr(i,j,nt_aero+2+4*(naero-1),iblk)*vice(i,j,iblk) &
00346 + trcr(i,j,nt_aero+3+4*(naero-1),iblk)*vice(i,j,iblk)
00347 enddo
00348 enddo
00349 enddo
00350 aerototn= global_sum(work1, distrb_info, field_loc_center, tarean)
00351 aerotots= global_sum(work1, distrb_info, field_loc_center, tareas)
00352 aeromx1n = global_maxval(work1, distrb_info, lmask_n)
00353 aeromx1s = global_maxval(work1, distrb_info, lmask_s)
00354 if (my_task == master_task) then
00355 write(nu_diag,*) 'aero: ',naero,' faero : ',&
00356 faeron, faeros
00357 write(nu_diag,*) 'aero: ',naero,' fsoot : ',&
00358 fsootn, fsoots
00359 write(nu_diag,*) 'aero: ',naero,' faero-fsoot : ',&
00360 faeron-fsootn, faeros-fsoots
00361 write(nu_diag,*) 'aero: ',naero,' aerotot : ',&
00362 aerototn, aerotots
00363 write(nu_diag,*) 'aero: ',naero,' aerotot change: ',&
00364 aerototn-totaeron(naero), aerotots-totaeros(naero)
00365 write(nu_diag,*) 'aero: ',naero,' aeromax agg: ',&
00366 aeromx1n,aeromx1s
00367 endif
00368
00369
00370
00371
00372
00373
00374
00375
00376
00377
00378
00379
00380
00381
00382
00383
00384
00385
00386
00387
00388
00389
00390
00391
00392
00393
00394
00395
00396
00397
00398
00399
00400
00401
00402
00403
00404
00405
00406 enddo
00407 endif
00408
00409
00410
00411 do iblk = 1, nblocks
00412 do j = 1, ny_block
00413 do i = 1, nx_block
00414 work1(i,j,iblk) = sqrt(uvel(i,j,iblk)**2 &
00415 + vvel(i,j,iblk)**2)
00416 enddo
00417 enddo
00418 enddo
00419
00420
00421 umaxn = global_maxval(work1, distrb_info, lmask_n)
00422 umaxs = global_maxval(work1, distrb_info, lmask_s)
00423
00424
00425
00426
00427 if (check_umax) then
00428 if (umaxn > umax_stab) then
00429 do iblk = 1, nblocks
00430 do j = 1, ny_block
00431 do i = 1, nx_block
00432 if (abs(work1(i,j,iblk) - umaxn) < puny) then
00433 write(nu_diag,*) ' '
00434 write(nu_diag,*) 'Warning, large ice speed'
00435 write(nu_diag,*) 'my_task, iblk, i, j, umaxn:', &
00436 my_task, iblk, i, j, umaxn
00437 endif
00438 enddo
00439 enddo
00440 enddo
00441 elseif (umaxs > umax_stab) then
00442 do iblk = 1, nblocks
00443 do j = 1, ny_block
00444 do i = 1, nx_block
00445 if (abs(work1(i,j,iblk) - umaxs) < puny) then
00446 write(nu_diag,*) ' '
00447 write(nu_diag,*) 'Warning, large ice speed'
00448 write(nu_diag,*) 'my_task, iblk, i, j, umaxs:', &
00449 my_task, iblk, i, j, umaxs
00450 endif
00451 enddo
00452 enddo
00453 enddo
00454 endif
00455 endif
00456
00457
00458
00459 pmaxn = global_maxval(strength, distrb_info, lmask_n)
00460 pmaxs = global_maxval(strength, distrb_info, lmask_s)
00461
00462 pmaxn = pmaxn / c1000
00463 pmaxs = pmaxs / c1000
00464
00465 if (print_global) then
00466
00467
00468
00469 do iblk = 1, nblocks
00470 do j = 1, ny_block
00471 do i = 1, nx_block
00472 work1(i,j,iblk) = esno(i,j,iblk) + eice(i,j,iblk)
00473 enddo
00474 enddo
00475 enddo
00476
00477
00478 etotn = global_sum(work1, distrb_info, &
00479 field_loc_center, tarean)
00480 etots = global_sum(work1, distrb_info, &
00481 field_loc_center, tareas)
00482
00483
00484
00485
00486
00487
00488
00489
00490
00491
00492 evpn = global_sum_prod(evap, aice, distrb_info, &
00493 field_loc_center, tarean)
00494 evps = global_sum_prod(evap, aice, distrb_info, &
00495 field_loc_center, tareas)
00496 evpn = evpn*dt
00497 evps = evps*dt
00498
00499
00500 sfsaltn = global_sum(fsalt_gbm, distrb_info, &
00501 field_loc_center, tarean)
00502 sfsalts = global_sum(fsalt_gbm, distrb_info, &
00503 field_loc_center, tareas)
00504 sfsaltn = sfsaltn*dt
00505 sfsalts = sfsalts*dt
00506
00507
00508 sfreshn = global_sum(fresh_gbm, distrb_info, &
00509 field_loc_center, tarean)
00510 sfreshs = global_sum(fresh_gbm, distrb_info, &
00511 field_loc_center, tareas)
00512 sfreshn = sfreshn*dt
00513 sfreshs = sfreshs*dt
00514
00515
00516
00517 fhocnn = global_sum(fhocn_gbm, distrb_info, &
00518 field_loc_center, tarean)
00519 fhocns = global_sum(fhocn_gbm, distrb_info, &
00520 field_loc_center, tareas)
00521
00522
00523
00524
00525
00526
00527
00528
00529
00530
00531
00532
00533
00534 if (calc_Tsfc) then
00535
00536
00537 do iblk = 1, nblocks
00538 do j = 1, ny_block
00539 do i = 1, nx_block
00540 work1(i,j,iblk) = &
00541 (fswabs(i,j,iblk) - fswthru(i,j,iblk) &
00542 + flwout(i,j,iblk) &
00543 + fsens (i,j,iblk)) * aice(i,j,iblk) &
00544 + flw (i,j,iblk) * aice_init(i,j,iblk)
00545 enddo
00546 enddo
00547 enddo
00548
00549
00550 else
00551
00552
00553 do iblk = 1, nblocks
00554 do j = 1, ny_block
00555 do i = 1, nx_block
00556 work1(i,j,iblk) = &
00557 (fsurf(i,j,iblk) - flat(i,j,iblk)) &
00558 * aice(i,j,iblk)
00559 enddo
00560 enddo
00561 enddo
00562
00563
00564 endif
00565
00566 fhatmn = global_sum(work1, distrb_info, &
00567 field_loc_center, tarean)
00568 fhatms = global_sum(work1, distrb_info, &
00569 field_loc_center, tareas)
00570
00571
00572
00573 do iblk = 1, nblocks
00574 do j = 1, ny_block
00575 do i = 1, nx_block
00576 work1(i,j,iblk) = max(c0,frzmlt(i,j,iblk))
00577 enddo
00578 enddo
00579 enddo
00580
00581 fhfrzn = global_sum(work1, distrb_info, &
00582 field_loc_center, tarean)
00583 fhfrzs = global_sum(work1, distrb_info, &
00584 field_loc_center, tareas)
00585
00586
00587 rnn = global_sum_prod(frain, aice_init, distrb_info, &
00588 field_loc_center, tarean)
00589 rns = global_sum_prod(frain, aice_init, distrb_info, &
00590 field_loc_center, tareas)
00591 rnn = rnn*dt
00592 rns = rns*dt
00593
00594
00595 snn = global_sum_prod(fsnow, aice_init, distrb_info, &
00596 field_loc_center, tarean)
00597 sns = global_sum_prod(fsnow, aice_init, distrb_info, &
00598 field_loc_center, tareas)
00599 snn = snn*dt
00600 sns = sns*dt
00601
00602
00603
00604 work1(:,:,:) = frazil(:,:,:)*rhoi/dt
00605 frzn = global_sum(work1, distrb_info, &
00606 field_loc_center, tarean)
00607 frzs = global_sum(work1, distrb_info, field_loc_center, &
00608 tareas)
00609 frzn = frzn*dt
00610 frzs = frzs*dt
00611
00612
00613 micen = rhoi*shmaxn
00614 msnwn = rhos*snwmxn
00615 mices = rhoi*shmaxs
00616 msnws = rhos*snwmxs
00617
00618 mtotn = micen + msnwn
00619 mtots = mices + msnws
00620
00621
00622 delmin = mtotn - totmn
00623 delmis = mtots - totms
00624
00625
00626 delmxn = micen - totmin
00627 delmxs = mices - totmis
00628 if (.not. update_ocn_f) then
00629
00630 delmxn = delmxn - frzn
00631 delmxs = delmxs - frzs
00632 endif
00633
00634
00635 fluxn = c0
00636 fluxs = c0
00637 if( arean > c0) then
00638
00639 fluxn = rnn + snn + evpn - sfreshn
00640 if (.not. update_ocn_f) then
00641 fluxn = fluxn + frzn
00642 endif
00643 endif
00644 if( areas > c0) then
00645
00646 fluxs = rns + sns + evps - sfreshs
00647 if (.not. update_ocn_f) then
00648 fluxs = fluxs + frzs
00649 endif
00650 endif
00651
00652 werrn = (fluxn-delmin)/(mtotn+c1)
00653 werrs = (fluxs-delmis)/(mtots+c1)
00654
00655
00656 delein = etotn - toten
00657 deleis = etots - totes
00658
00659 fhatmn = fhatmn + ( - snn * Lfresh + evpn * Lvap ) / dt
00660 fhatms = fhatms + ( - sns * Lfresh + evps * Lvap ) / dt
00661
00662 hnetn = (fhatmn - fhocnn - fhfrzn) * dt
00663 hnets = (fhatms - fhocns - fhfrzs) * dt
00664
00665 herrn = (hnetn - delein) / (etotn - c1)
00666 herrs = (hnets - deleis) / (etots - c1)
00667
00668
00669 msltn = micen*ice_ref_salinity*p001
00670 mslts = mices*ice_ref_salinity*p001
00671
00672
00673 delmsltn = delmxn*ice_ref_salinity*p001
00674 delmslts = delmxs*ice_ref_salinity*p001
00675
00676
00677 serrn = (sfsaltn + delmsltn) / (msltn + c1)
00678 serrs = (sfsalts + delmslts) / (mslts + c1)
00679
00680 endif
00681
00682 if (print_points) then
00683
00684
00685
00686
00687
00688
00689 do n = 1, npnt
00690 if (my_task == pmloc(n)) then
00691 i = piloc(n)
00692 j = pjloc(n)
00693 iblk = pbloc(n)
00694
00695 pTair(n) = Tair(i,j,iblk) - Tffresh
00696 pQa(n) = Qa(i,j,iblk)
00697 pfsnow(n) = fsnow(i,j,iblk)*dt/rhos
00698 pfrain(n) = frain(i,j,iblk)*dt/rhow
00699 pfsw(n) = fsw(i,j,iblk)
00700 pflw(n) = flw(i,j,iblk)
00701 paice(n) = aice(i,j,iblk)
00702
00703 hiavg(n) = c0
00704 hsavg(n) = c0
00705 if (paice(n) /= c0) then
00706 hiavg(n) = vice(i,j,iblk)/paice(n)
00707 hsavg(n) = vsno(i,j,iblk)/paice(n)
00708 endif
00709 pTsfc(n) = trcr(i,j,nt_Tsfc,iblk)
00710 pevap(n) = evap(i,j,iblk)*dt/rhoi
00711 pfswabs(n) = fswabs(i,j,iblk)
00712 pflwout(n) = flwout(i,j,iblk)
00713 pflat(n) = flat(i,j,iblk)
00714 pfsens(n) = fsens(i,j,iblk)
00715 pfsurf(n) = fsurf(i,j,iblk)
00716 pfcondtop(n) = fcondtop(i,j,iblk)
00717 pmeltt(n) = meltt(i,j,iblk)
00718 pmeltb(n) = meltb(i,j,iblk)
00719 pmeltl(n) = meltl(i,j,iblk)
00720 psnoice(n) = snoice(i,j,iblk)
00721 pfrazil(n) = frazil(i,j,iblk)
00722 pcongel(n) = congel(i,j,iblk)
00723 pdhi(n) = vice(i,j,iblk) - pdhi(n)
00724 pdhs(n) = vsno(i,j,iblk) - pdhs(n)
00725 pde(n) = -(eice(i,j,iblk) &
00726 + esno(i,j,iblk) - pde(n)) / dt
00727 psst(n) = sst(i,j,iblk)
00728 pTf(n) = Tf(i,j,iblk)
00729 pfhocn(n) = -fhocn(i,j,iblk)
00730
00731 endif
00732
00733 call broadcast_scalar(pTair (n), pmloc(n))
00734 call broadcast_scalar(pQa (n), pmloc(n))
00735 call broadcast_scalar(pfsnow (n), pmloc(n))
00736 call broadcast_scalar(pfrain (n), pmloc(n))
00737 call broadcast_scalar(pfsw (n), pmloc(n))
00738 call broadcast_scalar(pflw (n), pmloc(n))
00739 call broadcast_scalar(paice (n), pmloc(n))
00740 call broadcast_scalar(hsavg (n), pmloc(n))
00741 call broadcast_scalar(hiavg (n), pmloc(n))
00742 call broadcast_scalar(pTsfc (n), pmloc(n))
00743 call broadcast_scalar(pevap (n), pmloc(n))
00744 call broadcast_scalar(pfswabs (n), pmloc(n))
00745 call broadcast_scalar(pflwout (n), pmloc(n))
00746 call broadcast_scalar(pflat (n), pmloc(n))
00747 call broadcast_scalar(pfsens (n), pmloc(n))
00748 call broadcast_scalar(pfsurf (n), pmloc(n))
00749 call broadcast_scalar(pfcondtop(n), pmloc(n))
00750 call broadcast_scalar(pmeltt (n), pmloc(n))
00751 call broadcast_scalar(pmeltb (n), pmloc(n))
00752 call broadcast_scalar(pmeltl (n), pmloc(n))
00753 call broadcast_scalar(psnoice (n), pmloc(n))
00754 call broadcast_scalar(pfrazil (n), pmloc(n))
00755 call broadcast_scalar(pcongel (n), pmloc(n))
00756 call broadcast_scalar(pdhi (n), pmloc(n))
00757 call broadcast_scalar(pdhs (n), pmloc(n))
00758 call broadcast_scalar(pde (n), pmloc(n))
00759 call broadcast_scalar(psst (n), pmloc(n))
00760 call broadcast_scalar(pTf (n), pmloc(n))
00761 call broadcast_scalar(pfhocn (n), pmloc(n))
00762
00763 enddo
00764 endif
00765
00766
00767
00768
00769
00770 if (my_task == master_task) then
00771 if (grid_type == 'panarctic') then
00772 write (nu_diag,799) 'Arctic diagnostics'
00773 write (nu_diag,801) 'total ice area (km^2) = ',arean
00774 write (nu_diag,801) 'total ice extent(km^2) = ',extentn
00775 write (nu_diag,801) 'total ice volume (m^3) = ',shmaxn
00776 write (nu_diag,801) 'total snw volume (m^3) = ',snwmxn
00777 write (nu_diag,801) 'tot kinetic energy (J) = ',ketotn
00778 write (nu_diag,800) 'rms ice speed (m/s) = ',urmsn
00779 write (nu_diag,800) 'average albedo = ',albtotn
00780 write (nu_diag,800) 'max ice volume (m) = ',hmaxn
00781 write (nu_diag,800) 'max ice speed (m/s) = ',umaxn
00782 write (nu_diag,900) 'max strength (kN/m) = ',pmaxn
00783
00784 if (print_global) then
00785
00786 #ifdef CCSMCOUPLED
00787 if (prescribed_ice) then
00788 write (nu_diag,*) '----------------------------'
00789 write (nu_diag,*) 'This is the prescribed ice option.'
00790 write (nu_diag,*) 'Heat and water will not be conserved.'
00791 else
00792 #endif
00793
00794 write (nu_diag,*) '----------------------------'
00795 write (nu_diag,801) 'arwt rain h2o kg in dt = ',rnn
00796 write (nu_diag,801) 'arwt snow h2o kg in dt = ',snn
00797 write (nu_diag,801) 'arwt evap h2o kg in dt = ',evpn
00798 write (nu_diag,801) 'arwt frzl h2o kg in dt = ',frzn
00799 write (nu_diag,801) 'arwt frsh h2o kg in dt = ',sfreshn
00800
00801 write (nu_diag,801) 'arwt ice mass (kg) = ',micen
00802 write (nu_diag,801) 'arwt snw mass (kg) = ',msnwn
00803
00804 write (nu_diag,801) 'arwt tot mass (kg) = ',mtotn
00805 write (nu_diag,801) 'arwt tot mass chng(kg) = ',delmin
00806 write (nu_diag,801) 'arwt water flux = ',fluxn
00807 if (update_ocn_f) then
00808 write (nu_diag,*) '(=rain+snow+evap-fresh) '
00809 else
00810 write (nu_diag,*) '(=rain+snow+evap+frzl-fresh) '
00811 endif
00812 write (nu_diag,801) 'water flux error = ',werrn
00813 #ifdef CCSMCOUPLED
00814 endif
00815 #endif
00816 write (nu_diag,*) '----------------------------'
00817 write (nu_diag,801) 'arwt atm heat flux (W) = ',fhatmn
00818 write (nu_diag,801) 'arwt ocn heat flux (W) = ',fhocnn
00819 write (nu_diag,801) 'arwt frzl heat flux(W) = ',fhfrzn
00820 write (nu_diag,801) 'arwt tot energy (J) = ',etotn
00821 write (nu_diag,801) 'arwt net heat (J) = ',hnetn
00822 write (nu_diag,801) 'arwt tot energy chng(J)= ',delein
00823 write (nu_diag,801) 'arwt heat error = ',herrn
00824
00825 write (nu_diag,*) '----------------------------'
00826 write (nu_diag,801) 'arwt salt mass (kg) = ',msltn
00827 write (nu_diag,801) 'arwt salt mass chng(kg)= ',delmsltn
00828 write (nu_diag,801) 'arwt salt flx in dt(kg)= ',sfsaltn
00829 write (nu_diag,801) 'arwt salt flx error = ',serrn
00830 write (nu_diag,*) '----------------------------'
00831
00832 endif
00833
00834 else
00835
00836 write(nu_diag,899) 'Arctic','Antarctic'
00837
00838 write(nu_diag,901) 'total ice area (km^2) = ',arean, areas
00839 write(nu_diag,901) 'total ice extent(km^2) = ',extentn,extents
00840 write(nu_diag,901) 'total ice volume (m^3) = ',shmaxn, shmaxs
00841 write(nu_diag,901) 'total snw volume (m^3) = ',snwmxn, snwmxs
00842 write(nu_diag,901) 'tot kinetic energy (J) = ',ketotn, ketots
00843 write(nu_diag,900) 'rms ice speed (m/s) = ',urmsn, urmss
00844 write(nu_diag,900) 'average albedo = ',albtotn,albtots
00845 write(nu_diag,900) 'max ice volume (m) = ',hmaxn, hmaxs
00846 write(nu_diag,900) 'max ice speed (m/s) = ',umaxn, umaxs
00847 write(nu_diag,900) 'max strength (kN/m) = ',pmaxn, pmaxs
00848
00849 if (print_global) then
00850
00851 write(nu_diag,*) '----------------------------'
00852 write(nu_diag,901) 'arwt rain h2o kg in dt = ',rnn,rns
00853 write(nu_diag,901) 'arwt snow h2o kg in dt = ',snn,sns
00854 write(nu_diag,901) 'arwt evap h2o kg in dt = ',evpn,evps
00855 write(nu_diag,901) 'arwt frzl h2o kg in dt = ',frzn,frzs
00856 write(nu_diag,901) 'arwt frsh h2o kg in dt = ',sfreshn,sfreshs
00857
00858 write(nu_diag,901) 'arwt ice mass (kg) = ',micen,mices
00859 write(nu_diag,901) 'arwt snw mass (kg) = ',msnwn,msnws
00860
00861 write(nu_diag,901) 'arwt tot mass (kg) = ',mtotn,mtots
00862 write(nu_diag,901) 'arwt tot mass chng(kg) = ',delmin,delmis
00863 write(nu_diag,901) 'arwt water flux = ',fluxn,fluxs
00864 if (update_ocn_f) then
00865 write (nu_diag,*) '(=rain+snow+evap-fresh) '
00866 else
00867 write (nu_diag,*) '(=rain+snow+evap+frzl-fresh) '
00868 endif
00869 write(nu_diag,901) 'water flux error = ',werrn,werrs
00870
00871 write(nu_diag,*) '----------------------------'
00872 write(nu_diag,901) 'arwt atm heat flux (W) = ',fhatmn,fhatms
00873 write(nu_diag,901) 'arwt ocn heat flux (W) = ',fhocnn,fhocns
00874 write(nu_diag,901) 'arwt frzl heat flux(W) = ',fhfrzn,fhfrzs
00875 write(nu_diag,901) 'arwt tot energy (J) = ',etotn,etots
00876 write(nu_diag,901) 'arwt net heat (J) = ',hnetn,hnets
00877 write(nu_diag,901) 'arwt tot energy chng(J)= ',delein,deleis
00878 write(nu_diag,901) 'arwt heat error = ',herrn,herrs
00879
00880 write(nu_diag,*) '----------------------------'
00881 write(nu_diag,901) 'arwt salt mass (kg) = ',msltn,mslts
00882 write(nu_diag,901) 'arwt salt mass chng(kg)= ',delmsltn, &
00883 delmslts
00884 write(nu_diag,901) 'arwt salt flx in dt(kg)= ',sfsaltn, &
00885 sfsalts
00886 write(nu_diag,901) 'arwt salt flx error = ',serrn,serrs
00887 write(nu_diag,*) '----------------------------'
00888
00889 endif
00890 endif
00891
00892 call flush_fileunit(nu_diag)
00893
00894
00895
00896
00897
00898 if (print_points) then
00899
00900 write(nu_diag,*) ' '
00901 write(nu_diag,902) ' Lat, Long ',plat(1),plon(1), &
00902 plat(2),plon(2)
00903 write(nu_diag,903) ' my_task, iblk, i, j ', &
00904 pmloc(1),pbloc(1),piloc(1),pjloc(1), &
00905 pmloc(2),pbloc(2),piloc(2),pjloc(2)
00906 write(nu_diag,*) '----------atm----------'
00907 write(nu_diag,900) 'air temperature (C) = ',pTair(1),pTair(2)
00908 write(nu_diag,900) 'specific humidity = ',pQa(1),pQa(2)
00909 write(nu_diag,900) 'snowfall (m) = ',pfsnow(1), &
00910 pfsnow(2)
00911 write(nu_diag,900) 'rainfall (m) = ',pfrain(1), &
00912 pfrain(2)
00913 if (.not.calc_Tsfc) then
00914 write(nu_diag,900) 'total surface heat flux= ',pfsurf(1),pfsurf(2)
00915 write(nu_diag,900) 'top sfc conductive flux= ',pfcondtop(1), &
00916 pfcondtop(2)
00917 write(nu_diag,900) 'latent heat flx = ',pflat(1),pflat(2)
00918 else
00919 write(nu_diag,900) 'shortwave radiation sum= ',pfsw(1),pfsw(2)
00920 write(nu_diag,900) 'longwave radiation = ',pflw(1),pflw(2)
00921 endif
00922 write(nu_diag,*) '----------ice----------'
00923 write(nu_diag,900) 'area fraction = ',paice(1),paice(2)
00924 write(nu_diag,900) 'avg ice thickness (m) = ',hiavg(1),hiavg(2)
00925 write(nu_diag,900) 'avg snow depth (m) = ',hsavg(1),hsavg(2)
00926 if (calc_Tsfc) then
00927 write(nu_diag,900) 'surface temperature(C) = ',pTsfc(1),pTsfc(2)
00928 write(nu_diag,900) 'absorbed shortwave flx = ',pfswabs(1), &
00929 pfswabs(2)
00930 write(nu_diag,900) 'outward longwave flx = ',pflwout(1), &
00931 pflwout(2)
00932 write(nu_diag,900) 'sensible heat flx = ',pfsens(1), &
00933 pfsens(2)
00934 write(nu_diag,900) 'latent heat flx = ',pflat(1),pflat(2)
00935 endif
00936 write(nu_diag,900) 'subl/cond (m ice) = ',pevap(1),pevap(2)
00937 write(nu_diag,900) 'top melt (m) = ',pmeltt(1) &
00938 ,pmeltt(2)
00939 write(nu_diag,900) 'bottom melt (m) = ',pmeltb(1) &
00940 ,pmeltb(2)
00941 write(nu_diag,900) 'lateral melt (m) = ',pmeltl(1) &
00942 ,pmeltl(2)
00943 write(nu_diag,900) 'new ice (m) = ',pfrazil(1), &
00944 pfrazil(2)
00945 write(nu_diag,900) 'congelation (m) = ',pcongel(1), &
00946 pcongel(2)
00947 write(nu_diag,900) 'snow-ice (m) = ',psnoice(1), &
00948 psnoice(2)
00949 write(nu_diag,900) 'effective dhi (m) = ',pdhi(1),pdhi(2)
00950 write(nu_diag,900) 'effective dhs (m) = ',pdhs(1),pdhs(2)
00951 write(nu_diag,900) 'intnl enrgy chng(W/m^2)= ',pde (1),pde (2)
00952 write(nu_diag,*) '----------ocn----------'
00953 write(nu_diag,900) 'sst (C) = ',psst(1),psst(2)
00954 write(nu_diag,900) 'freezing temp (C) = ',pTf(1),pTf(2)
00955 write(nu_diag,900) 'heat used (W/m^2) = ',pfhocn(1), &
00956 pfhocn(2)
00957
00958 endif
00959 endif
00960
00961 799 format (27x,a24)
00962 800 format (a25,2x,f24.17)
00963 801 format (a25,2x,1pe24.17)
00964 899 format (27x,a24,2x,a24)
00965 900 format (a25,2x,f24.17,2x,f24.17)
00966 901 format (a25,2x,1pe24.17,2x,1pe24.17)
00967 902 format (a25,10x,f6.1,1x,f6.1,9x,f6.1,1x,f6.1)
00968 903 format (a25,5x,i4,1x,i4,1x,i4,1x,i4,7x,i4,1x,i4,1x,i4,1x,i4)
00969
00970 end subroutine runtime_diags
00971
00972
00973
00974
00975
00976
00977
00978
00979 subroutine init_mass_diags
00980
00981
00982
00983
00984
00985
00986
00987
00988
00989
00990
00991 use ice_global_reductions
00992 use ice_grid
00993 use ice_state
00994 use ice_broadcast
00995
00996
00997
00998
00999
01000 integer (kind=int_kind) :: n, k, ii, jj, i, j, iblk
01001 integer (kind=int_kind) :: naero
01002
01003 real (kind=dbl_kind) ::
01004 shmaxn, snwmxn, shmaxs, snwmxs
01005
01006 real (kind=dbl_kind), dimension(nx_block,ny_block,max_blocks) ::
01007 work1, work2
01008
01009
01010
01011 shmaxn = global_sum(vice, distrb_info, field_loc_center, tarean)
01012 shmaxs = global_sum(vice, distrb_info, field_loc_center, tareas)
01013
01014
01015 snwmxn = global_sum(vsno, distrb_info, field_loc_center, tarean)
01016 snwmxs = global_sum(vsno, distrb_info, field_loc_center, tareas)
01017
01018
01019 totmin = rhoi*shmaxn
01020 totmis = rhoi*shmaxs
01021
01022
01023 totmn = totmin + rhos*snwmxn
01024 totms = totmis + rhos*snwmxs
01025
01026
01027
01028
01029 do iblk = 1, nblocks
01030 do j=1,ny_block
01031 do i=1,nx_block
01032 work1(i,j,iblk) = esno(i,j,iblk) + eice(i,j,iblk)
01033 enddo
01034 enddo
01035 enddo
01036
01037
01038 toten = global_sum(work1, distrb_info, field_loc_center, tarean)
01039 totes = global_sum(work1, distrb_info, field_loc_center, tareas)
01040
01041 if (tr_aero) then
01042 do naero=1,n_aero
01043 do iblk = 1, nblocks
01044 do j = 1, ny_block
01045 do i = 1, nx_block
01046 work1(i,j,iblk) = trcr(i,j,nt_aero +4*(naero-1),iblk)*vsno(i,j,iblk) &
01047 + trcr(i,j,nt_aero+1+4*(naero-1),iblk)*vsno(i,j,iblk) &
01048 + trcr(i,j,nt_aero+2+4*(naero-1),iblk)*vice(i,j,iblk) &
01049 + trcr(i,j,nt_aero+3+4*(naero-1),iblk)*vice(i,j,iblk)
01050 enddo
01051 enddo
01052 enddo
01053 totaeron(naero)= global_sum(work1, distrb_info, field_loc_center, tarean)
01054 totaeros(naero)= global_sum(work1, distrb_info, field_loc_center, tareas)
01055 enddo
01056 endif
01057
01058 if (print_points) then
01059
01060 do n = 1, npnt
01061
01062 if (my_task == pmloc(n)) then
01063 i = piloc(n)
01064 j = pjloc(n)
01065 iblk = pbloc(n)
01066
01067 pdhi(n) = vice(i,j,iblk)
01068 pdhs(n) = vsno(i,j,iblk)
01069 pde(n) = esno(i,j,iblk) + eice(i,j,iblk)
01070 endif
01071
01072 enddo
01073
01074 endif
01075
01076 end subroutine init_mass_diags
01077
01078
01079
01080
01081
01082
01083
01084
01085 subroutine init_diags
01086
01087
01088
01089
01090
01091
01092
01093
01094
01095
01096
01097 use ice_grid
01098 use ice_blocks
01099 use ice_broadcast
01100 use ice_global_reductions
01101 use ice_gather_scatter
01102
01103
01104
01105
01106
01107 real (kind=dbl_kind) ::
01108 latdis ,
01109 londis ,
01110 totdis ,
01111 mindis ,
01112 mindis_g
01113
01114 integer (kind=int_kind) ::
01115 n ,
01116 i,j ,
01117 iblk ,
01118 ilo,ihi,jlo,jhi
01119
01120 character (char_len) :: label(npnt)
01121
01122 type (block) ::
01123 this_block
01124
01125 if (print_points) then
01126
01127 if (my_task==master_task) then
01128 write(nu_diag,*) ' '
01129 write(nu_diag,*) ' Find indices of diagnostic points '
01130 endif
01131
01132
01133 label(1)(1:40) = 'Near North Pole pack ice '
01134 label(2)(1:40) = 'Weddell Sea '
01135
01136 piloc(:) = 0
01137 pjloc(:) = 0
01138 pbloc(:) = 0
01139 pmloc(:) = -999
01140 plat(:) = -999._dbl_kind
01141 plon(:) = -999._dbl_kind
01142
01143
01144 do n = 1, npnt
01145 if (lonpnt(n) > c180) lonpnt(n) = lonpnt(n) - c360
01146
01147 iindx = 0
01148 jindx = 0
01149 bindx = 0
01150 mindis = 540.0_dbl_kind
01151
01152
01153
01154
01155 do iblk = 1, nblocks
01156 this_block = get_block(blocks_ice(iblk),iblk)
01157 ilo = this_block%ilo
01158 ihi = this_block%ihi
01159 jlo = this_block%jlo
01160 jhi = this_block%jhi
01161
01162 do j = jlo, jhi
01163 do i = ilo, ihi
01164 if (hm(i,j,iblk) > p5) then
01165 latdis = abs(latpnt(n)-TLAT(i,j,iblk)*rad_to_deg)
01166 londis = abs(lonpnt(n)-TLON(i,j,iblk)*rad_to_deg) &
01167 * cos(TLAT(i,j,iblk))
01168 totdis = sqrt(latdis**2 + londis**2)
01169 if (totdis < mindis) then
01170 mindis = totdis
01171 jindx = j
01172 iindx = i
01173 bindx = iblk
01174 endif
01175 endif
01176 enddo
01177 enddo
01178 enddo
01179
01180
01181
01182 mindis_g = global_minval(mindis, distrb_info)
01183
01184
01185 if (abs(mindis_g - mindis) < puny) then
01186 piloc(n) = iindx
01187 pjloc(n) = jindx
01188 pbloc(n) = bindx
01189 pmloc(n) = my_task
01190 plat(n) = TLAT(iindx,jindx,bindx)*rad_to_deg
01191 plon(n) = TLON(iindx,jindx,bindx)*rad_to_deg
01192 endif
01193
01194
01195 piloc(n) = global_maxval(piloc(n), distrb_info)
01196 pjloc(n) = global_maxval(pjloc(n), distrb_info)
01197 pbloc(n) = global_maxval(pbloc(n), distrb_info)
01198 pmloc(n) = global_maxval(pmloc(n), distrb_info)
01199 plat(n) = global_maxval(plat(n), distrb_info)
01200 plon(n) = global_maxval(plon(n), distrb_info)
01201
01202
01203 if (my_task==master_task) then
01204 write(nu_diag,*) ' '
01205 write(nu_diag,100) n,latpnt(n),lonpnt(n),plat(n),plon(n), &
01206 piloc(n), pjloc(n), pbloc(n), pmloc(n)
01207 endif
01208 100 format(' found point',i4/ &
01209 ' lat lon TLAT TLON i j block task'/ &
01210 4(f6.1,1x),1x,4(i4,2x) )
01211
01212 enddo
01213 endif
01214
01215 end subroutine init_diags
01216
01217
01218
01219
01220
01221
01222
01223
01224 subroutine print_state(plabel,i,j,iblk)
01225
01226
01227
01228
01229
01230
01231
01232
01233
01234
01235
01236
01237
01238
01239
01240
01241
01242
01243
01244
01245
01246
01247
01248
01249
01250 use ice_state
01251 use ice_itd
01252 use ice_flux
01253
01254
01255
01256 character (len=20), intent(in) :: plabel
01257
01258 integer (kind=int_kind), intent(in) ::
01259 i, j ,
01260 iblk
01261
01262
01263
01264 real (kind=dbl_kind) ::
01265 eidebug, esdebug,
01266 qi, qs, Tsnow
01267
01268 integer (kind=int_kind) :: n, k
01269
01270 write(nu_diag,*) plabel
01271 write(nu_diag,*) 'istep1, my_task, i, j, iblk:', &
01272 istep1, my_task, i, j, iblk
01273 write(nu_diag,*) ' '
01274 write(nu_diag,*) 'aice0', aice0(i,j,iblk)
01275 do n = 1, ncat
01276 write(nu_diag,*) ' '
01277 write(nu_diag,*) 'n =',n
01278 write(nu_diag,*) 'aicen', aicen(i,j,n,iblk)
01279 write(nu_diag,*) 'vicen', vicen(i,j,n,iblk)
01280 write(nu_diag,*) 'vsnon', vsnon(i,j,n,iblk)
01281 if (aicen(i,j,n,iblk) > puny) then
01282 write(nu_diag,*) 'hin', vicen(i,j,n,iblk)/aicen(i,j,n,iblk)
01283 write(nu_diag,*) 'hsn', vsnon(i,j,n,iblk)/aicen(i,j,n,iblk)
01284 endif
01285 write(nu_diag,*) 'Tsfcn',trcrn(i,j,nt_Tsfc,n,iblk)
01286 write(nu_diag,*) ' '
01287 enddo
01288
01289 eidebug = c0
01290 do n = 1,ncat
01291 do k = 1,nilyr
01292 write(nu_diag,*) 'eicen, cat ',n,' layer ',k, &
01293 eicen(i,j,ilyr1(n)+k-1,iblk)
01294 eidebug = eidebug + eicen(i,j,ilyr1(n)+k-1,iblk)
01295 if (aicen(i,j,n,iblk) > puny) then
01296 qi = eicen(i,j,ilyr1(n)+k-1,iblk) / &
01297 (vicen(i,j,n,iblk)/real(nilyr,kind=dbl_kind))
01298 write(nu_diag,*) 'qi/rhoi', qi/rhoi
01299 endif
01300 enddo
01301 write(nu_diag,*) ' '
01302 enddo
01303 write(nu_diag,*) 'eice(i,j)',eidebug
01304 write(nu_diag,*) ' '
01305
01306 esdebug = c0
01307 do n = 1,ncat
01308 if (vsnon(i,j,n,iblk) > puny) then
01309 do k = 1,nslyr
01310 write(nu_diag,*) 'esnon, cat ',n,' layer ',k, &
01311 esnon(i,j,slyr1(n)+k-1,iblk)
01312 esdebug = esdebug + esnon(i,j,slyr1(n)+k-1,iblk)
01313 qs = esnon(i,j,slyr1(n)+k-1,iblk) / &
01314 (vsnon(i,j,n,iblk)/real(nslyr,kind=dbl_kind))
01315 Tsnow = (Lfresh + qs/rhos) / cp_ice
01316 write(nu_diag,*) 'qs/rhos', qs/rhos
01317 write(nu_diag,*) 'Tsnow', Tsnow
01318 enddo
01319 write(nu_diag,*) ' '
01320 endif
01321 enddo
01322 write(nu_diag,*) 'esno(i,j)',esdebug
01323 write(nu_diag,*) ' '
01324
01325 write(nu_diag,*) 'uvel(i,j)',uvel(i,j,iblk)
01326 write(nu_diag,*) 'vvel(i,j)',vvel(i,j,iblk)
01327
01328 write(nu_diag,*) ' '
01329 write(nu_diag,*) 'atm states and fluxes'
01330 write(nu_diag,*) ' uatm = ',uatm (i,j,iblk)
01331 write(nu_diag,*) ' vatm = ',vatm (i,j,iblk)
01332 write(nu_diag,*) ' potT = ',potT (i,j,iblk)
01333 write(nu_diag,*) ' Tair = ',Tair (i,j,iblk)
01334 write(nu_diag,*) ' Qa = ',Qa (i,j,iblk)
01335 write(nu_diag,*) ' rhoa = ',rhoa (i,j,iblk)
01336 write(nu_diag,*) ' swvdr = ',swvdr(i,j,iblk)
01337 write(nu_diag,*) ' swvdf = ',swvdf(i,j,iblk)
01338 write(nu_diag,*) ' swidr = ',swidr(i,j,iblk)
01339 write(nu_diag,*) ' swidf = ',swidf(i,j,iblk)
01340 write(nu_diag,*) ' flw = ',flw (i,j,iblk)
01341 write(nu_diag,*) ' frain = ',frain(i,j,iblk)
01342 write(nu_diag,*) ' fsnow = ',fsnow(i,j,iblk)
01343 write(nu_diag,*) ' '
01344 write(nu_diag,*) 'ocn states and fluxes'
01345 write(nu_diag,*) ' frzmlt = ',frzmlt (i,j,iblk)
01346 write(nu_diag,*) ' sst = ',sst (i,j,iblk)
01347 write(nu_diag,*) ' sss = ',sss (i,j,iblk)
01348 write(nu_diag,*) ' Tf = ',Tf (i,j,iblk)
01349 write(nu_diag,*) ' uocn = ',uocn (i,j,iblk)
01350 write(nu_diag,*) ' vocn = ',vocn (i,j,iblk)
01351 write(nu_diag,*) ' strtltx = ',strtltx(i,j,iblk)
01352 write(nu_diag,*) ' strtlty = ',strtlty(i,j,iblk)
01353 write(nu_diag,*) ' '
01354 write(nu_diag,*) 'srf states and fluxes'
01355 write(nu_diag,*) ' Tref = ',Tref (i,j,iblk)
01356 write(nu_diag,*) ' Qref = ',Qref (i,j,iblk)
01357 write(nu_diag,*) ' fsens = ',fsens (i,j,iblk)
01358 write(nu_diag,*) ' flat = ',flat (i,j,iblk)
01359 write(nu_diag,*) ' evap = ',evap (i,j,iblk)
01360 write(nu_diag,*) ' flwout = ',flwout(i,j,iblk)
01361 write(nu_diag,*) ' '
01362
01363 end subroutine print_state
01364
01365
01366
01367 end module ice_diagnostics
01368
01369
01370
01371
01372
01373
01374