00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 module ice_step_mod
00023
00024
00025
00026 use ice_atmo
00027 use ice_calendar
00028 use ice_communicate
00029 use ice_diagnostics
00030 use ice_domain
00031 use ice_dyn_evp
00032 use ice_fileunits
00033 use ice_flux
00034 use ice_grid
00035 use ice_history
00036 use ice_restart
00037 use ice_itd
00038 use ice_kinds_mod
00039 use ice_mechred
00040 use ice_ocean
00041 use ice_orbital
00042 use ice_shortwave
00043 use ice_state
00044 use ice_therm_itd
00045 use ice_therm_vertical
00046 use ice_timers
00047 use ice_transport_driver
00048 use ice_transport_remap
00049 use perf_mod, only: t_startf, t_stopf, t_barrierf
00050
00051 implicit none
00052 private
00053 save
00054
00055
00056
00057 public :: step_therm2, step_dynamics, &
00058 prep_radiation, step_radiation
00059
00060
00061
00062
00063
00064 contains
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079 subroutine prep_radiation(dt)
00080
00081
00082
00083
00084
00085 real (kind=dbl_kind), intent(in) ::
00086 dt
00087
00088
00089
00090 integer (kind=int_kind) ::
00091 i,j,n,iblk
00092
00093 if (calc_Tsfc) then
00094
00095
00096 do iblk = 1, nblocks
00097 call prep_radiation_iblk(dt, iblk)
00098 end do
00099
00100
00101 else
00102
00103
00104 do iblk = 1, nblocks
00105 do n = 1, ncat
00106 do j = 1, ny_block
00107 do i = 1, nx_block
00108 fswsfcn(i,j,n,iblk) = c0
00109 fswintn(i,j,n,iblk) = c0
00110 fswthrun(i,j,n,iblk) = c0
00111 enddo
00112 enddo
00113 enddo
00114 Iswabsn(:,:,:,iblk) = c0
00115 Sswabsn(:,:,:,iblk) = c0
00116 enddo
00117
00118 endif
00119
00120 end subroutine prep_radiation
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135 subroutine prep_radiation_iblk (dt, iblk)
00136
00137
00138
00139
00140
00141 real (kind=dbl_kind), intent(in) ::
00142 dt
00143
00144 integer (kind=int_kind), intent(in) ::
00145 iblk
00146
00147
00148
00149 integer (kind=int_kind) ::
00150 i, j, ij ,
00151 ilo,ihi,jlo,jhi,
00152 n ,
00153 il1, il2 ,
00154 sl1, sl2
00155
00156 integer (kind=int_kind) ::
00157 icells
00158
00159 integer (kind=int_kind), dimension(nx_block*ny_block) ::
00160 indxi, indxj
00161
00162
00163 real (kind=dbl_kind), dimension (nx_block,ny_block) ::
00164 fsn
00165 real (kind=dbl_kind), dimension (nx_block,ny_block,nslyr) ::
00166 rhosnwn ,
00167 rsnwn
00168
00169
00170 real (kind=dbl_kind), dimension (nx_block,ny_block) ::
00171 fpn ,
00172 hpn
00173
00174 real (kind=dbl_kind) :: netsw, netsw_old, ar
00175
00176 type (block) ::
00177 this_block
00178
00179 logical (kind=log_kind) ::
00180 l_stop
00181
00182 integer (kind=int_kind) ::
00183 istop, jstop
00184
00185 l_stop = .false.
00186
00187 this_block = get_block(blocks_ice(iblk),iblk)
00188 ilo = this_block%ilo
00189 ihi = this_block%ihi
00190 jlo = this_block%jlo
00191 jhi = this_block%jhi
00192
00193
00194
00195
00196
00197 do j = jlo, jhi
00198 do i = ilo, ihi
00199 if (aice(i,j,iblk) > c0 .and. scale_factor(i,j,iblk) > puny) then
00200 netsw = swvdr(i,j,iblk)*(c1 - alvdr_gbm(i,j,iblk)) &
00201 + swvdf(i,j,iblk)*(c1 - alvdf_gbm(i,j,iblk)) &
00202 + swidr(i,j,iblk)*(c1 - alidr_gbm(i,j,iblk)) &
00203 + swidf(i,j,iblk)*(c1 - alidf_gbm(i,j,iblk))
00204 scale_factor(i,j,iblk) = netsw / scale_factor(i,j,iblk)
00205 else
00206 scale_factor(i,j,iblk) = c1
00207 endif
00208 fswfac(i,j,iblk) = scale_factor(i,j,iblk)
00209 enddo
00210 enddo
00211
00212 do n = 1, ncat
00213
00214
00215
00216
00217
00218 icells = 0
00219 do j = jlo, jhi
00220 do i = ilo, ihi
00221 if (aicen(i,j,n,iblk) > puny) then
00222 icells = icells + 1
00223 indxi(icells) = i
00224 indxj(icells) = j
00225 endif
00226 enddo
00227 enddo
00228
00229
00230
00231
00232
00233 il1 = ilyr1(n)
00234 il2 = ilyrn(n)
00235 sl1 = slyr1(n)
00236 sl2 = slyrn(n)
00237
00238 do ij = 1, icells
00239 i = indxi(ij)
00240 j = indxj(ij)
00241
00242 fswsfcn(i,j,n,iblk) = scale_factor(i,j,iblk)*fswsfcn (i,j,n,iblk)
00243 fswintn(i,j,n,iblk) = scale_factor(i,j,iblk)*fswintn (i,j,n,iblk)
00244 fswthrun(i,j,n,iblk)= scale_factor(i,j,iblk)*fswthrun(i,j,n,iblk)
00245 Sswabsn(i,j,sl1:sl2,iblk) = &
00246 scale_factor(i,j,iblk)*Sswabsn(i,j,sl1:sl2,iblk)
00247 Iswabsn(i,j,il1:il2,iblk) = &
00248 scale_factor(i,j,iblk)*Iswabsn(i,j,il1:il2,iblk)
00249 enddo
00250 enddo
00251
00252 end subroutine prep_radiation_iblk
00253
00254
00255
00256
00257
00258
00259
00260
00261
00262
00263
00264
00265
00266
00267
00268
00269
00270
00271
00272 subroutine step_therm2 (dt)
00273
00274
00275
00276
00277
00278 real (kind=dbl_kind), intent(in) ::
00279 dt
00280
00281
00282
00283 integer (kind=int_kind) ::
00284 iblk,
00285 i, j
00286
00287 call ice_timer_start(timer_column)
00288 call ice_timer_start(timer_thermo)
00289
00290
00291
00292 do iblk = 1, nblocks
00293 call step_therm2_iblk(dt, iblk)
00294 end do
00295
00296
00297
00298
00299
00300
00301 call ice_timer_start(timer_bound)
00302 call bound_state (aicen, trcrn, &
00303 vicen, vsnon, &
00304 eicen, esnon)
00305 call ice_timer_stop(timer_bound)
00306
00307
00308 do iblk = 1, nblocks
00309
00310
00311
00312
00313
00314 call aggregate (nx_block, ny_block, &
00315 aicen(:,:,:,iblk), trcrn(:,:,:,:,iblk), &
00316 vicen(:,:,:,iblk), vsnon(:,:, :,iblk), &
00317 eicen(:,:,:,iblk), esnon(:,:, :,iblk), &
00318 aice (:,:, iblk), trcr (:,:,:, iblk), &
00319 vice (:,:, iblk), vsno (:,:, iblk), &
00320 eice (:,:, iblk), esno (:,:, iblk), &
00321 aice0(:,:, iblk), tmask(:,:, iblk), &
00322 trcr_depend)
00323
00324
00325
00326
00327
00328
00329 do j = 1, ny_block
00330 do i = 1, nx_block
00331 daidtt(i,j,iblk) = (aice(i,j,iblk) - daidtt(i,j,iblk)) / dt
00332 dvidtt(i,j,iblk) = (vice(i,j,iblk) - dvidtt(i,j,iblk)) / dt
00333 enddo
00334 enddo
00335
00336 enddo
00337
00338
00339
00340 call ice_timer_stop(timer_thermo)
00341 call ice_timer_stop(timer_column)
00342
00343 end subroutine step_therm2
00344
00345
00346
00347
00348
00349
00350
00351
00352
00353
00354
00355
00356
00357
00358
00359
00360
00361
00362
00363
00364 subroutine step_therm2_iblk (dt, iblk)
00365
00366
00367
00368
00369
00370
00371 real (kind=dbl_kind), intent(in) ::
00372 dt
00373
00374 integer (kind=int_kind), intent(in) ::
00375 iblk
00376
00377
00378
00379
00380
00381
00382
00383
00384 integer (kind=int_kind) ::
00385 ilo,ihi,jlo,jhi,
00386 i, j, n
00387
00388 integer (kind=int_kind) ::
00389 icells
00390
00391 integer (kind=int_kind), dimension(nx_block*ny_block) ::
00392 indxi, indxj
00393
00394 type (block) ::
00395 this_block
00396
00397 logical (kind=log_kind) ::
00398 l_stop
00399
00400 integer (kind=int_kind) ::
00401 istop, jstop
00402
00403 l_stop = .false.
00404
00405 this_block = get_block(blocks_ice(iblk),iblk)
00406 ilo = this_block%ilo
00407 ihi = this_block%ihi
00408 jlo = this_block%jlo
00409 jhi = this_block%jhi
00410
00411
00412
00413
00414
00415 do j = 1, ny_block
00416 do i = 1, nx_block
00417 fresh (i,j,iblk) = fresh(i,j,iblk) &
00418 + frain(i,j,iblk)*aice(i,j,iblk)
00419 enddo
00420 enddo
00421
00422
00423
00424
00425
00426
00427 call ice_timer_start(timer_catconv)
00428
00429 if (kitd == 1) then
00430
00431
00432
00433 call aggregate_area (nx_block, ny_block, &
00434 aicen(:,:,:,iblk), &
00435 aice (:,:, iblk), aice0(:,:,iblk))
00436
00437
00438
00439
00440
00441 icells = 0
00442 do j = jlo,jhi
00443 do i = ilo,ihi
00444 if (aice(i,j,iblk) > puny) then
00445 icells = icells + 1
00446 indxi(icells) = i
00447 indxj(icells) = j
00448 endif
00449 enddo
00450 enddo
00451
00452 if (icells > 0) then
00453
00454 call linear_itd (nx_block, ny_block, &
00455 icells, indxi, indxj, &
00456 trcr_depend, &
00457 aicen_init(:,:,:,iblk), &
00458 vicen_init(:,:,:,iblk), &
00459 aicen (:,:,:,iblk), &
00460 trcrn (:,:,:,:,iblk), &
00461 vicen (:,:,:,iblk), &
00462 vsnon (:,:,:,iblk), &
00463 eicen (:,:,:,iblk), &
00464 esnon (:,:,:,iblk), &
00465 aice (:,:, iblk), &
00466 aice0 (:,:, iblk), &
00467 l_stop, &
00468 istop, jstop)
00469
00470 if (l_stop) then
00471 write (nu_diag,*) 'istep1, my_task, iblk =', &
00472 istep1, my_task, iblk
00473 write (nu_diag,*) 'Global block:', this_block%block_id
00474 if (istop > 0 .and. jstop > 0) &
00475 write(nu_diag,*) 'Global i and j:', &
00476 this_block%i_glob(istop), &
00477 this_block%j_glob(jstop)
00478 call abort_ice ('ice: Linear ITD error')
00479 endif
00480
00481 endif
00482
00483 endif
00484
00485 call ice_timer_stop(timer_catconv)
00486
00487
00488
00489
00490
00491
00492 icells = 0
00493 do j = 1, ny_block
00494 do i = 1, nx_block
00495 if (tmask(i,j,iblk)) then
00496 icells = icells + 1
00497 indxi(icells) = i
00498 indxj(icells) = j
00499 endif
00500 enddo
00501 enddo
00502
00503 call add_new_ice (nx_block, ny_block, &
00504 icells, &
00505 indxi, indxj, &
00506 tmask (:,:, iblk), dt, &
00507 aicen (:,:,:,iblk), &
00508 trcrn (:,:,:,:,iblk), &
00509 vicen (:,:,:,iblk), &
00510 eicen (:,:,:,iblk), &
00511 aice0 (:,:, iblk), &
00512 aice (:,:, iblk), &
00513 frzmlt (:,:, iblk), &
00514 frazil (:,:, iblk), &
00515 frz_onset(:,:, iblk), yday, &
00516 fresh (:,:, iblk), &
00517 fsalt (:,:, iblk), &
00518 Tf (:,:, iblk), l_stop, &
00519 istop, jstop)
00520
00521 if (l_stop) then
00522 write (nu_diag,*) 'istep1, my_task, iblk =', &
00523 istep1, my_task, iblk
00524 write (nu_diag,*) 'Global block:', this_block%block_id
00525 if (istop > 0 .and. jstop > 0) &
00526 write(nu_diag,*) 'Global i and j:', &
00527 this_block%i_glob(istop), &
00528 this_block%j_glob(jstop)
00529 call abort_ice ('ice: add_new_ice error')
00530 endif
00531
00532
00533
00534
00535 call lateral_melt (nx_block, ny_block, &
00536 ilo, ihi, jlo, jhi, &
00537 dt, &
00538 fresh (:,:, iblk), &
00539 fsalt (:,:, iblk), &
00540 fhocn (:,:, iblk), &
00541 fsoot (:,:,:,iblk), &
00542 rside (:,:, iblk), &
00543 meltl (:,:, iblk), &
00544 aicen (:,:,:,iblk), &
00545 vicen (:,:,:,iblk), &
00546 vsnon (:,:,:,iblk), &
00547 eicen (:,:,:,iblk), &
00548 esnon (:,:,:,iblk), &
00549 trcrn (:,:,:,:,iblk) )
00550
00551
00552
00553
00554
00555
00556
00557
00558
00559
00560
00561
00562
00563
00564
00565
00566
00567
00568
00569
00570
00571
00572
00573 call cleanup_itd (nx_block, ny_block, &
00574 ilo, ihi, jlo, jhi, &
00575 dt, &
00576 aicen (:,:,:,iblk), trcrn (:,:,:,:,iblk), &
00577 vicen (:,:,:,iblk), vsnon (:,:, :,iblk), &
00578 eicen (:,:,:,iblk), esnon (:,:, :,iblk), &
00579 aice0 (:,:, iblk), aice (:,:,iblk), &
00580 trcr_depend, &
00581 fresh (:,:, iblk), fsalt (:,:, iblk), &
00582 fhocn (:,:, iblk), fsoot (:,:,:,iblk), &
00583 tr_aero, &
00584 heat_capacity, l_stop, &
00585 istop, jstop)
00586
00587 if (l_stop) then
00588 write (nu_diag,*) 'istep1, my_task, iblk =', &
00589 istep1, my_task, iblk
00590 write (nu_diag,*) 'Global block:', this_block%block_id
00591 if (istop > 0 .and. jstop > 0) &
00592 write(nu_diag,*) 'Global i and j:', &
00593 this_block%i_glob(istop), &
00594 this_block%j_glob(jstop)
00595 call abort_ice ('ice: ITD cleanup error')
00596 endif
00597
00598 end subroutine step_therm2_iblk
00599
00600
00601
00602
00603
00604
00605
00606
00607
00608
00609
00610
00611
00612
00613
00614
00615
00616
00617
00618
00619
00620 subroutine step_dynamics(dt_dyn,dt_thm)
00621
00622
00623
00624
00625
00626 real (kind=dbl_kind), intent(in) ::
00627 dt_dyn ,
00628 dt_thm
00629
00630
00631
00632 type (block) ::
00633 this_block
00634
00635 integer (kind=int_kind) ::
00636 iblk ,
00637 i,j ,
00638 ilo,ihi,jlo,jhi
00639
00640 integer (kind=int_kind) ::
00641 icells
00642
00643 integer (kind=int_kind), dimension(nx_block*ny_block) ::
00644 indxi, indxj
00645
00646 logical (kind=log_kind) ::
00647 l_stop
00648
00649 integer (kind=int_kind) ::
00650 istop, jstop
00651
00652 call init_history_dyn
00653
00654 dynCnt = dynCnt + 1
00655
00656
00657
00658
00659
00660 call t_barrierf ('cice_dyn_evp_BARRIER',MPI_COMM_ICE)
00661 call t_startf ('cice_dyn_evp')
00662 if (kdyn == 1) call evp (dt_dyn)
00663 call t_stopf ('cice_dyn_evp')
00664
00665
00666
00667
00668
00669 call t_barrierf ('cice_dyn_horz_transport_BARRIER',MPI_COMM_ICE)
00670 call t_startf ('cice_dyn_horz_transport')
00671 if (advection == 'upwind') then
00672 call transport_upwind (dt_dyn)
00673 else
00674 call transport_remap (dt_dyn)
00675 endif
00676 call t_stopf ('cice_dyn_horz_transport')
00677
00678
00679
00680
00681
00682 call ice_timer_start(timer_column)
00683 call ice_timer_start(timer_ridge)
00684 call t_barrierf ('cice_dyn_ridge_BARRIER',MPI_COMM_ICE)
00685 call t_startf ('cice_dyn_ridge')
00686
00687 l_stop = .false.
00688
00689
00690
00691 do iblk = 1, nblocks
00692 this_block = get_block(blocks_ice(iblk), iblk)
00693 ilo = this_block%ilo
00694 ihi = this_block%ihi
00695 jlo = this_block%jlo
00696 jhi = this_block%jhi
00697
00698
00699
00700
00701
00702
00703
00704
00705 icells = 0
00706 do j = jlo, jhi
00707 do i = ilo, ihi
00708
00709
00710 if (tmask(i,j,iblk)) then
00711 icells = icells + 1
00712 indxi(icells) = i
00713 indxj(icells) = j
00714 endif
00715 enddo
00716 enddo
00717
00718 if (icells > 0) then
00719
00720 call ridge_ice (nx_block, ny_block, &
00721 dt_dyn, dt_thm, icells, &
00722 indxi, indxj, &
00723
00724 rdg_conv(:,:, iblk), rdg_shear (:,:, iblk), &
00725 aicen (:,:,:,iblk), trcrn (:,:,:,:,iblk), &
00726 vicen (:,:,:,iblk), vsnon (:,:,:,iblk), &
00727 eicen (:,:,:,iblk), esnon (:,:,:,iblk), &
00728 aice0 (:,:, iblk), &
00729 trcr_depend, l_stop, &
00730 istop, jstop, &
00731 dardg1dt(:,:,iblk), dardg2dt (:,:,iblk), &
00732 dvirdgdt(:,:,iblk), opening (:,:,iblk), &
00733 fresh (:,:,iblk), fhocn (:,:,iblk), &
00734 fsoot (:,:,:,iblk))
00735
00736 if (l_stop) then
00737 write (nu_diag,*) 'istep1, my_task, iblk =', &
00738 istep1, my_task, iblk
00739 write (nu_diag,*) 'Global block:', this_block%block_id
00740 if (istop > 0 .and. jstop > 0) &
00741 write(nu_diag,*) 'Global i and j:', &
00742 this_block%i_glob(istop), &
00743 this_block%j_glob(jstop)
00744 call abort_ice ('ice: Ridging error')
00745 endif
00746
00747 endif
00748
00749 enddo
00750
00751
00752 call ice_timer_stop(timer_ridge)
00753 call t_stopf ('cice_dyn_ridge')
00754
00755 call t_barrierf ('cice_dyn_column_BARRIER',MPI_COMM_ICE)
00756 call t_startf ('cice_dyn_column')
00757
00758
00759
00760 do iblk = 1, nblocks
00761 this_block = get_block(blocks_ice(iblk), iblk)
00762 ilo = this_block%ilo
00763 ihi = this_block%ihi
00764 jlo = this_block%jlo
00765 jhi = this_block%jhi
00766
00767
00768
00769
00770
00771
00772 call cleanup_itd (nx_block, ny_block, &
00773 ilo, ihi, jlo, jhi, &
00774 dt_thm, &
00775 aicen (:,:,:,iblk), trcrn (:,:,:,:,iblk), &
00776 vicen (:,:,:,iblk), vsnon (:,:, :,iblk), &
00777 eicen (:,:,:,iblk), esnon (:,:, :,iblk), &
00778 aice0 (:,:, iblk), aice (:,:,iblk), &
00779 trcr_depend, &
00780 fresh (:,:, iblk), fsalt (:,:, iblk), &
00781 fhocn (:,:, iblk), fsoot (:,:,:,iblk), &
00782 tr_aero, &
00783 heat_capacity, l_stop, &
00784 istop, jstop)
00785
00786 if (l_stop) then
00787 write (nu_diag,*) 'istep1, my_task, iblk =', &
00788 istep1, my_task, iblk
00789 write (nu_diag,*) 'Global block:', this_block%block_id
00790 if (istop > 0 .and. jstop > 0) &
00791 write(nu_diag,*) 'Global i and j:', &
00792 this_block%i_glob(istop), &
00793 this_block%j_glob(jstop)
00794 call abort_ice ('ice: ITD cleanup error')
00795 endif
00796
00797 enddo
00798
00799
00800 call t_stopf ('cice_dyn_column')
00801
00802
00803
00804
00805
00806 call t_barrierf ('cice_dyn_bound_BARRIER',MPI_COMM_ICE)
00807 call t_startf ('cice_dyn_bound')
00808 call ice_timer_start(timer_bound)
00809 call bound_state (aicen, trcrn, &
00810 vicen, vsnon, &
00811 eicen, esnon)
00812 call ice_timer_stop(timer_bound)
00813 call t_stopf ('cice_dyn_bound')
00814
00815 call t_barrierf ('cice_dyn_agg_BARRIER',MPI_COMM_ICE)
00816 call t_startf ('cice_dyn_agg')
00817
00818
00819 do iblk = 1, nblocks
00820
00821
00822
00823
00824
00825 call aggregate (nx_block, ny_block, &
00826 aicen(:,:,:,iblk), trcrn(:,:,:,:,iblk), &
00827 vicen(:,:,:,iblk), vsnon(:,:, :,iblk), &
00828 eicen(:,:,:,iblk), esnon(:,:, :,iblk), &
00829 aice (:,:, iblk), trcr (:,:,:, iblk), &
00830 vice (:,:, iblk), vsno (:,:, iblk), &
00831 eice (:,:, iblk), esno (:,:, iblk), &
00832 aice0(:,:, iblk), tmask(:,:, iblk), &
00833 trcr_depend)
00834
00835
00836
00837
00838
00839 this_block = get_block(blocks_ice(iblk),iblk)
00840 ilo = this_block%ilo
00841 ihi = this_block%ihi
00842 jlo = this_block%jlo
00843 jhi = this_block%jhi
00844
00845 do j = jlo,jhi
00846 do i = ilo,ihi
00847 dvidtd(i,j,iblk) = (vice(i,j,iblk) - dvidtd(i,j,iblk)) /dt_dyn
00848 daidtd(i,j,iblk) = (aice(i,j,iblk) - daidtd(i,j,iblk)) /dt_dyn
00849 enddo
00850 enddo
00851
00852 enddo
00853
00854
00855 call t_stopf ('cice_dyn_agg')
00856 call ice_timer_stop(timer_column)
00857
00858 end subroutine step_dynamics
00859
00860
00861
00862
00863
00864
00865
00866
00867
00868
00869
00870
00871
00872
00873 subroutine step_radiation(dt)
00874
00875
00876
00877
00878
00879 real (kind=dbl_kind), intent(in) ::
00880 dt
00881
00882
00883
00884 integer (kind=int_kind) ::
00885 i,j,n,iblk
00886
00887
00888 alvdr(:,:,:) = c0
00889 alvdf(:,:,:) = c0
00890 alidr(:,:,:) = c0
00891 alidf(:,:,:) = c0
00892 Sswabsn(:,:,:,:) = c0
00893
00894
00895
00896
00897
00898
00899
00900
00901 if (calc_Tsfc) then
00902
00903
00904 do iblk = 1, nblocks
00905 call step_radiation_iblk(dt, iblk)
00906 end do
00907
00908
00909 else
00910
00911
00912 do iblk = 1, nblocks
00913 do n = 1, ncat
00914 do j = 1, ny_block
00915 do i = 1, nx_block
00916 alvdrn(i,j,n,iblk) = c0
00917 alidrn(i,j,n,iblk) = c0
00918 alvdfn(i,j,n,iblk) = c0
00919 alidfn(i,j,n,iblk) = c0
00920 fswsfcn(i,j,n,iblk) = c0
00921 fswintn(i,j,n,iblk) = c0
00922 fswthrun(i,j,n,iblk) = c0
00923 enddo
00924 enddo
00925 enddo
00926 Iswabsn(:,:,:,iblk) = c0
00927 Sswabsn(:,:,:,iblk) = c0
00928 enddo
00929
00930 endif
00931
00932 end subroutine step_radiation
00933
00934
00935
00936
00937
00938
00939
00940
00941
00942
00943
00944
00945
00946
00947 subroutine step_radiation_iblk (dt, iblk)
00948
00949
00950
00951
00952
00953 real (kind=dbl_kind), intent(in) ::
00954 dt
00955
00956 integer (kind=int_kind), intent(in) ::
00957 iblk
00958
00959
00960
00961 integer (kind=int_kind) ::
00962 i, j, ij ,
00963 ilo,ihi,jlo,jhi,
00964 n ,
00965 il1, il2 ,
00966 sl1, sl2
00967
00968 integer (kind=int_kind) ::
00969 icells
00970
00971 integer (kind=int_kind), dimension(nx_block*ny_block) ::
00972 indxi, indxj
00973
00974
00975 real (kind=dbl_kind), dimension (nx_block,ny_block) ::
00976 fsn
00977 real (kind=dbl_kind), dimension (nx_block,ny_block,nslyr) ::
00978 rhosnwn ,
00979 rsnwn
00980
00981
00982 real (kind=dbl_kind), dimension (nx_block,ny_block) ::
00983 fpn ,
00984 hpn
00985
00986 type (block) ::
00987 this_block
00988
00989 logical (kind=log_kind) ::
00990 l_stop
00991
00992 integer (kind=int_kind) ::
00993 istop, jstop
00994
00995 l_stop = .false.
00996
00997 this_block = get_block(blocks_ice(iblk),iblk)
00998 ilo = this_block%ilo
00999 ihi = this_block%ihi
01000 jlo = this_block%jlo
01001 jhi = this_block%jhi
01002
01003
01004
01005
01006
01007
01008
01009
01010 if (trim(shortwave) == 'dEdd') then
01011
01012
01013 icells = 0
01014 do j = 1, ny_block
01015 do i = 1, nx_block
01016 if (tmask(i,j,iblk)) then
01017 icells = icells + 1
01018 indxi(icells) = i
01019 indxj(icells) = j
01020 endif
01021 enddo
01022 enddo
01023
01024 call compute_coszen (nx_block, ny_block, &
01025 icells, &
01026 indxi, indxj, &
01027 tlat (:,:,iblk), tlon(:,:,iblk), &
01028 coszen(:,:,iblk), dt)
01029
01030 else
01031 coszen(:,:,iblk) = p5
01032 endif
01033
01034 do n = 1, ncat
01035
01036
01037
01038
01039
01040 icells = 0
01041 do j = jlo, jhi
01042 do i = ilo, ihi
01043 if (aicen(i,j,n,iblk) > puny) then
01044 icells = icells + 1
01045 indxi(icells) = i
01046 indxj(icells) = j
01047 endif
01048 enddo
01049 enddo
01050
01051
01052
01053
01054
01055 il1 = ilyr1(n)
01056 il2 = ilyrn(n)
01057 sl1 = slyr1(n)
01058 sl2 = slyrn(n)
01059
01060 if (trim(shortwave) == 'dEdd') then
01061
01062
01063
01064
01065
01066 call shortwave_dEdd_set_snow(nx_block, ny_block, &
01067 icells, &
01068 indxi, indxj, &
01069 aicen(:,:,n,iblk), vsnon(:,:,n,iblk), &
01070 trcrn(:,:,nt_Tsfc,n,iblk), fsn, &
01071 rhosnwn, rsnwn)
01072
01073
01074 if (.not. tr_pond) then
01075
01076
01077 call shortwave_dEdd_set_pond(nx_block, ny_block, &
01078 icells, &
01079 indxi, indxj, &
01080 aicen(:,:,n,iblk), &
01081 trcrn(:,:,nt_Tsfc,n,iblk), &
01082 fsn, fpn, &
01083 hpn)
01084
01085 else
01086
01087
01088 fpn(:,:) = apondn(:,:,n,iblk)
01089 hpn(:,:) = hpondn(:,:,n,iblk)
01090
01091 endif
01092
01093 call shortwave_dEdd(nx_block, ny_block, &
01094 icells, &
01095 indxi, indxj, &
01096 coszen(:,:, iblk), &
01097 aicen(:,:,n,iblk), vicen(:,:,n,iblk), &
01098 vsnon(:,:,n,iblk), fsn, &
01099 rhosnwn, rsnwn, &
01100 fpn, hpn, &
01101 trcrn(:,:,:,n,iblk), tarea(:,:,iblk), &
01102 swvdr(:,:, iblk), swvdf(:,:, iblk), &
01103 swidr(:,:, iblk), swidf(:,:, iblk), &
01104 alvdrn(:,:,n,iblk),alvdfn(:,:,n,iblk), &
01105 alidrn(:,:,n,iblk),alidfn(:,:,n,iblk), &
01106 fswsfcn(:,:,n,iblk),fswintn(:,:,n,iblk),&
01107 fswthrun(:,:,n,iblk), &
01108 Sswabsn(:,:,sl1:sl2,iblk), &
01109 Iswabsn(:,:,il1:il2,iblk), &
01110 albicen(:,:,n,iblk), &
01111 albsnon(:,:,n,iblk),albpndn(:,:,n,iblk))
01112
01113 else
01114
01115 Sswabsn(:,:,sl1:sl2,iblk) = c0
01116
01117 call shortwave_ccsm3(nx_block, ny_block, &
01118 icells, &
01119 indxi, indxj, &
01120 aicen(:,:,n,iblk), vicen(:,:,n,iblk), &
01121 vsnon(:,:,n,iblk), &
01122 trcrn(:,:,nt_Tsfc,n,iblk), &
01123 swvdr(:,:, iblk), swvdf(:,:, iblk), &
01124 swidr(:,:, iblk), swidf(:,:, iblk), &
01125 alvdrn(:,:,n,iblk),alidrn(:,:,n,iblk), &
01126 alvdfn(:,:,n,iblk),alidfn(:,:,n,iblk), &
01127 fswsfcn(:,:,n,iblk),fswintn(:,:,n,iblk),&
01128 fswthrun(:,:,n,iblk), &
01129 Iswabsn(:,:,il1:il2,iblk), &
01130 albicen(:,:,n,iblk),albsnon(:,:,n,iblk))
01131
01132 endif
01133
01134 enddo
01135
01136 end subroutine step_radiation_iblk
01137
01138
01139
01140 end module ice_step_mod
01141
01142