00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044 module ice_mechred
00045
00046
00047
00048 use ice_kinds_mod
00049 use ice_domain_size
00050 use ice_constants
00051 use ice_fileunits
00052 use ice_itd, only: hin_max, ilyr1, slyr1, column_sum, &
00053 column_conservation_check, compute_tracers
00054
00055
00056
00057 implicit none
00058 save
00059
00060
00061
00062
00063
00064 integer (kind=int_kind) ::
00065 kstrength ,
00066
00067 krdg_partic , &
00068
00069 krdg_redist
00070
00071
00072 real (kind=dbl_kind), parameter ::
00073 Cf = 17._dbl_kind ,
00074 Cs = p25 ,
00075 Cp = p5*gravit*(rhow-rhoi)*rhoi/rhow,
00076 fsnowrdg = p5 ,
00077 Gstar = p15 ,
00078
00079 astar = p05 , &
00080
00081 maxraft= c1 , &
00082
00083 Hstar = c25 , &
00084
00085
00086 mu_rdg = c4 , &
00087
00088 Pstar = 2.75e4_dbl_kind, &
00089
00090 Cstar = c20
00091
00092
00093 logical (kind=log_kind), parameter ::
00094 l_conservation_check = .true.
00095
00096
00097
00098
00099
00100 contains
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119 subroutine ridge_ice (nx_block, ny_block, &
00120 dt_dyn, dt_thm, icells, &
00121 indxi, indxj, &
00122 rdg_conv, rdg_shear, &
00123 aicen, trcrn, &
00124 vicen, vsnon, &
00125 eicen, esnon, &
00126 aice0, &
00127 trcr_depend, l_stop, &
00128 istop, jstop, &
00129 dardg1dt, dardg2dt, &
00130 dvirdgdt, opening, &
00131 fresh, fhocn, &
00132 fsoot)
00133
00134
00135
00136
00137
00138 integer (kind=int_kind), intent(in) ::
00139 nx_block, ny_block,
00140 icells
00141
00142 integer (kind=int_kind), dimension (nx_block*ny_block),
00143 intent(in) ::
00144 indxi, indxj
00145
00146 real (kind=dbl_kind), intent(in) ::
00147 dt_dyn ,
00148 dt_thm
00149
00150 real (kind=dbl_kind), dimension(nx_block,ny_block), intent(in) ::
00151 rdg_conv,
00152 rdg_shear
00153
00154 real (kind=dbl_kind), dimension (nx_block,ny_block,ncat),
00155 intent(inout) ::
00156 aicen ,
00157 vicen ,
00158 vsnon
00159
00160 real (kind=dbl_kind), dimension (nx_block,ny_block,ntrcr,ncat),
00161 intent(inout) ::
00162 trcrn
00163
00164 real (kind=dbl_kind), dimension (nx_block,ny_block,ntilyr),
00165 intent(inout) ::
00166 eicen
00167
00168 real (kind=dbl_kind), dimension (nx_block,ny_block,ntslyr),
00169 intent(inout) ::
00170 esnon
00171
00172 real (kind=dbl_kind), dimension (nx_block,ny_block),
00173 intent(inout) ::
00174 aice0
00175
00176 integer (kind=int_kind), dimension(ntrcr), intent(in) ::
00177 trcr_depend
00178
00179 logical (kind=log_kind), intent(out) ::
00180 l_stop
00181
00182 integer (kind=int_kind), intent(out) ::
00183 istop, jstop
00184
00185
00186 real (kind=dbl_kind), dimension(nx_block,ny_block),
00187 intent(inout), optional ::
00188 dardg1dt ,
00189 dardg2dt ,
00190 dvirdgdt ,
00191 opening ,
00192 fresh ,
00193 fhocn
00194
00195 real (kind=dbl_kind), dimension(nx_block,ny_block,n_aeromx),
00196 intent(inout), optional ::
00197 fsoot
00198
00199
00200
00201 real (kind=dbl_kind), dimension (icells) ::
00202 asum ,
00203 aksum ,
00204 msnow_mlt ,
00205 esnow_mlt ,
00206 closing_net,
00207
00208 divu_adv , &
00209 opning , &
00210
00211
00212 ardg1 , &
00213 ardg2 , &
00214 virdg , &
00215 aopen
00216
00217 real (kind=dbl_kind), dimension (icells,n_aeromx) ::
00218 msoot
00219
00220 real (kind=dbl_kind), dimension (icells,0:ncat) ::
00221 apartic
00222
00223
00224 real (kind=dbl_kind), dimension (icells,ncat) ::
00225 hrmin ,
00226 hrmax ,
00227 hrexp ,
00228 krdg
00229
00230 real (kind=dbl_kind), dimension (icells) ::
00231 vice_init, vice_final,
00232 vsno_init, vsno_final,
00233 eice_init, eice_final,
00234 esno_init, esno_final
00235
00236 integer (kind=int_kind), parameter ::
00237 nitermax = 20
00238
00239 integer (kind=int_kind) ::
00240 i,j ,
00241 n ,
00242 niter ,
00243 ij
00244
00245 real (kind=dbl_kind) ::
00246 dti
00247
00248 logical (kind=log_kind) ::
00249 iterate_ridging,
00250 asum_error
00251
00252 character (len=char_len) ::
00253 fieldid
00254
00255
00256
00257
00258
00259 l_stop = .false.
00260 istop = 0
00261 jstop = 0
00262
00263 do ij = 1, icells
00264 msnow_mlt(ij) = c0
00265 esnow_mlt(ij) = c0
00266 msoot (ij,:) = c0
00267 ardg1 (ij) = c0
00268 ardg2 (ij) = c0
00269 virdg (ij) = c0
00270
00271 enddo
00272
00273
00274
00275
00276 call asum_ridging (nx_block, ny_block, &
00277 icells, indxi, indxj, &
00278 aicen, aice0, &
00279 asum)
00280
00281
00282
00283
00284 call ridge_prep (nx_block, ny_block, &
00285 icells, indxi, indxj, &
00286 dt_dyn, &
00287 rdg_conv, rdg_shear, &
00288 asum, closing_net, &
00289 divu_adv, opning)
00290
00291
00292
00293
00294
00295 if (l_conservation_check) then
00296
00297 call column_sum (nx_block, ny_block, &
00298 icells, indxi, indxj, &
00299 ncat, &
00300 vicen, vice_init)
00301
00302 call column_sum (nx_block, ny_block, &
00303 icells, indxi, indxj, &
00304 ncat, &
00305 vsnon, vsno_init)
00306
00307 call column_sum (nx_block, ny_block, &
00308 icells, indxi, indxj, &
00309 ntilyr, &
00310 eicen, eice_init)
00311
00312 call column_sum (nx_block, ny_block, &
00313 icells, indxi, indxj, &
00314 ntslyr, &
00315 esnon, esno_init)
00316
00317 endif
00318
00319 do niter = 1, nitermax
00320
00321
00322
00323
00324
00325
00326 call ridge_itd (nx_block, ny_block, &
00327 icells, indxi, indxj, &
00328 aicen, vicen, &
00329 aice0, &
00330 aksum, apartic, &
00331 hrmin, hrmax, &
00332 hrexp, krdg)
00333
00334
00335
00336
00337
00338 call ridge_shift (nx_block, ny_block, &
00339 icells, indxi, indxj, &
00340 dt_dyn, &
00341 aicen, trcrn, &
00342 vicen, vsnon, &
00343 eicen, esnon, &
00344 aice0, trcr_depend, &
00345 aksum, apartic, &
00346 hrmin, hrmax, &
00347 hrexp, krdg, &
00348 closing_net, opning, &
00349 ardg1, ardg2, &
00350 virdg, aopen, &
00351 msnow_mlt, esnow_mlt, &
00352 msoot, &
00353 l_stop, &
00354 istop, jstop)
00355
00356 if (l_stop) return
00357
00358
00359
00360
00361
00362
00363
00364 call asum_ridging (nx_block, ny_block, &
00365 icells, indxi, indxj, &
00366 aicen, aice0, &
00367 asum)
00368
00369 call ridge_check (nx_block, ny_block, &
00370 icells, indxi, indxj, &
00371 dt_dyn, &
00372 asum, closing_net, &
00373 divu_adv, opning, &
00374 iterate_ridging)
00375
00376
00377
00378
00379
00380 if (iterate_ridging) then
00381 write(nu_diag,*) 'Repeat ridging, niter =', niter
00382 else
00383 exit
00384 endif
00385
00386 if (niter == nitermax) then
00387 write(nu_diag,*) ' '
00388 write(nu_diag,*) 'Exceeded max number of ridging iterations'
00389 write(nu_diag,*) 'max =',nitermax
00390 l_stop = .true.
00391 return
00392 endif
00393
00394 enddo
00395
00396
00397
00398
00399
00400
00401 if (l_conservation_check) then
00402
00403 call column_sum (nx_block, ny_block, &
00404 icells, indxi, indxj, &
00405 ncat, &
00406 vicen, vice_final)
00407
00408 call column_sum (nx_block, ny_block, &
00409 icells, indxi, indxj, &
00410 ncat, &
00411 vsnon, vsno_final)
00412
00413 call column_sum (nx_block, ny_block, &
00414 icells, indxi, indxj, &
00415 ntilyr, &
00416 eicen, eice_final)
00417
00418 call column_sum (nx_block, ny_block, &
00419 icells, indxi, indxj, &
00420 ntslyr, &
00421 esnon, esno_final)
00422
00423 do ij = 1, icells
00424 vsno_final(ij) = vsno_final(ij) + msnow_mlt(ij)/rhos
00425 esno_final(ij) = esno_final(ij) + esnow_mlt(ij)
00426 enddo
00427
00428 fieldid = 'vice, ridging'
00429 call column_conservation_check (nx_block, ny_block, &
00430 icells, indxi, indxj, &
00431 fieldid, &
00432 vice_init, vice_final, &
00433 puny, l_stop, &
00434 istop, jstop)
00435 if (l_stop) return
00436
00437 fieldid = 'vsno, ridging'
00438 call column_conservation_check (nx_block, ny_block, &
00439 icells, indxi, indxj, &
00440 fieldid, &
00441 vsno_init, vsno_final, &
00442 puny, l_stop, &
00443 istop, jstop)
00444 if (l_stop) return
00445
00446 fieldid = 'eice, ridging'
00447 call column_conservation_check (nx_block, ny_block, &
00448 icells, indxi, indxj, &
00449 fieldid, &
00450 eice_init, eice_final, &
00451 puny*Lfresh*rhoi, &
00452 l_stop, &
00453 istop, jstop)
00454 if (l_stop) return
00455
00456 fieldid = 'esno, ridging'
00457 call column_conservation_check (nx_block, ny_block, &
00458 icells, indxi, indxj, &
00459 fieldid, &
00460 esno_init, esno_final, &
00461 puny*Lfresh*rhos, &
00462 l_stop, &
00463 istop, jstop)
00464 if (l_stop) return
00465
00466 endif
00467
00468
00469
00470
00471
00472 dti = c1/dt_dyn
00473
00474 if (present(dardg1dt)) then
00475 do ij = 1, icells
00476 i = indxi(ij)
00477 j = indxj(ij)
00478 dardg1dt(i,j) = ardg1(ij)*dti
00479 enddo
00480 endif
00481 if (present(dardg2dt)) then
00482 do ij = 1, icells
00483 i = indxi(ij)
00484 j = indxj(ij)
00485 dardg2dt(i,j) = ardg2(ij)*dti
00486 enddo
00487 endif
00488 if (present(dvirdgdt)) then
00489 do ij = 1, icells
00490 i = indxi(ij)
00491 j = indxj(ij)
00492 dvirdgdt(i,j) = virdg(ij)*dti
00493 enddo
00494 endif
00495 if (present(opening)) then
00496 do ij = 1, icells
00497 i = indxi(ij)
00498 j = indxj(ij)
00499 opening(i,j) = aopen(ij)*dti
00500 enddo
00501 endif
00502
00503
00504
00505
00506
00507 dti = c1/dt_thm
00508
00509 if (present(fsoot)) then
00510 do ij = 1, icells
00511 i = indxi(ij)
00512 j = indxj(ij)
00513 fsoot(i,j,:) = fsoot(i,j,:) + msoot(ij,:)*dti
00514 enddo
00515 endif
00516 if (present(fresh)) then
00517 do ij = 1, icells
00518 i = indxi(ij)
00519 j = indxj(ij)
00520 fresh(i,j) = fresh(i,j) + msnow_mlt(ij)*dti
00521 enddo
00522 endif
00523 if (present(fhocn)) then
00524 do ij = 1, icells
00525 i = indxi(ij)
00526 j = indxj(ij)
00527 fhocn(i,j) = fhocn(i,j) + esnow_mlt(ij)*dti
00528 enddo
00529 endif
00530
00531
00532
00533
00534
00535 do ij = 1, icells
00536 i = indxi(ij)
00537 j = indxj(ij)
00538 if (abs(asum(ij) - c1) > puny) then
00539 l_stop = .true.
00540 istop = i
00541 jstop = j
00542
00543 write(nu_diag,*) ' '
00544 write(nu_diag,*) 'Ridging error: total area > 1'
00545 write(nu_diag,*) 'i, j, area:', i, j, asum(ij)
00546 write(nu_diag,*) 'n, aicen:'
00547 write(nu_diag,*) 0, aice0(i,j)
00548 do n = 1, ncat
00549 write(nu_diag,*) n, aicen(i,j,n)
00550 enddo
00551 return
00552 endif
00553 enddo
00554
00555 end subroutine ridge_ice
00556
00557
00558
00559
00560
00561
00562
00563
00564
00565
00566
00567
00568
00569
00570
00571
00572
00573
00574
00575
00576 subroutine asum_ridging (nx_block, ny_block, &
00577 icells, indxi, indxj, &
00578 aicen, aice0, &
00579 asum)
00580
00581
00582
00583
00584
00585 integer (kind=int_kind), intent(in) ::
00586 nx_block, ny_block,
00587 icells
00588
00589 integer (kind=int_kind), dimension (nx_block*ny_block),
00590 intent(in) ::
00591 indxi, indxj
00592
00593
00594 real (kind=dbl_kind), dimension (nx_block,ny_block,ncat),
00595 intent(in) ::
00596 aicen
00597
00598 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(in) ::
00599 aice0
00600
00601 real (kind=dbl_kind), dimension (icells), intent(out)::
00602 asum
00603
00604
00605
00606 integer (kind=int_kind) ::
00607 i, j, n,
00608 ij
00609
00610
00611
00612
00613 do ij = 1, icells
00614 i = indxi(ij)
00615 j = indxj(ij)
00616 asum(ij) = aice0(i,j)
00617 enddo
00618
00619
00620
00621
00622
00623 do n = 1, ncat
00624
00625
00626
00627 do ij = 1, icells
00628 i = indxi(ij)
00629 j = indxj(ij)
00630 asum(ij) = asum(ij) + aicen(i,j,n)
00631 enddo
00632 enddo
00633
00634 end subroutine asum_ridging
00635
00636
00637
00638
00639
00640
00641
00642
00643
00644
00645
00646
00647
00648
00649
00650 subroutine ridge_prep (nx_block, ny_block, &
00651 icells, indxi, indxj, &
00652 dt, &
00653 rdg_conv, rdg_shear, &
00654 asum, closing_net, &
00655 divu_adv, opning)
00656
00657
00658
00659
00660
00661 integer (kind=int_kind), intent(in) ::
00662 nx_block, ny_block,
00663 icells
00664
00665 integer (kind=int_kind), dimension (nx_block*ny_block),
00666 intent(in) ::
00667 indxi, indxj
00668
00669
00670 real (kind=dbl_kind), intent(in) ::
00671 dt
00672
00673 real (kind=dbl_kind), dimension(nx_block,ny_block), intent(in) ::
00674 rdg_conv,
00675 rdg_shear
00676
00677 real (kind=dbl_kind), dimension(icells),
00678 intent(inout)::
00679 asum
00680
00681 real (kind=dbl_kind), dimension(icells),
00682 intent(out)::
00683 closing_net,
00684 divu_adv ,
00685 opning
00686
00687
00688
00689 real (kind=dbl_kind), parameter ::
00690 big = 1.0e+8_dbl_kind
00691
00692 integer (kind=int_kind) ::
00693 i,j,
00694 ij
00695
00696
00697
00698 hin_max(ncat) = big
00699
00700 do ij = 1, icells
00701 i = indxi(ij)
00702 j = indxj(ij)
00703
00704
00705
00706
00707
00708
00709
00710
00711
00712
00713
00714
00715
00716
00717
00718
00719
00720
00721
00722
00723
00724 closing_net(ij) = Cs*rdg_shear(i,j) + rdg_conv(i,j)
00725
00726
00727
00728
00729
00730
00731
00732
00733
00734
00735 divu_adv(ij) = (c1-asum(ij)) / dt
00736
00737 if (divu_adv(ij) < c0) &
00738 closing_net(ij) = max(closing_net(ij), -divu_adv(ij))
00739
00740
00741
00742
00743
00744 opning(ij) = closing_net(ij) + divu_adv(ij)
00745
00746 enddo
00747
00748 end subroutine ridge_prep
00749
00750
00751
00752
00753
00754
00755
00756
00757
00758
00759
00760
00761
00762
00763
00764
00765
00766
00767
00768
00769
00770
00771
00772
00773
00774
00775
00776 subroutine ridge_itd (nx_block, ny_block, &
00777 icells, indxi, indxj, &
00778 aicen, vicen, &
00779 aice0, &
00780 aksum, apartic, &
00781 hrmin, hrmax, &
00782 hrexp, krdg)
00783
00784
00785
00786
00787
00788 integer (kind=int_kind), intent(in) ::
00789 nx_block, ny_block,
00790 icells
00791
00792 integer (kind=int_kind), dimension (nx_block*ny_block),
00793 intent(in) ::
00794 indxi, indxj
00795
00796
00797 real (kind=dbl_kind), dimension (nx_block,ny_block,ncat),
00798 intent(in) ::
00799 aicen ,
00800 vicen
00801
00802 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(in) ::
00803 aice0
00804
00805 real (kind=dbl_kind), dimension (icells), intent(out)::
00806 aksum
00807
00808 real (kind=dbl_kind), dimension (icells,0:ncat),
00809 intent(out) ::
00810 apartic
00811
00812
00813 real (kind=dbl_kind), dimension (icells,ncat),
00814 intent(out) ::
00815 hrmin ,
00816 hrmax ,
00817 hrexp ,
00818 krdg
00819
00820
00821
00822 integer (kind=int_kind) ::
00823 i,j ,
00824 n ,
00825 ij
00826
00827 real (kind=dbl_kind), parameter ::
00828 Gstari = c1/Gstar,
00829 astari = c1/astar
00830
00831 real (kind=dbl_kind), dimension(icells,-1:ncat) ::
00832 Gsum
00833
00834 real (kind=dbl_kind), dimension(icells) ::
00835 work
00836
00837 real (kind=dbl_kind) ::
00838 hi ,
00839 hieff ,
00840 hrmean ,
00841 xtmp
00842
00843
00844
00845
00846
00847 do ij = 1, icells
00848 Gsum (ij,-1) = c0
00849 Gsum (ij,0) = c1
00850 apartic(ij,0) = c0
00851 enddo
00852
00853 do n = 1, ncat
00854 do ij = 1, icells
00855 Gsum (ij,n) = c1
00856 apartic(ij,n) = c0
00857 hrmin (ij,n) = c0
00858 hrmax (ij,n) = c0
00859 hrexp (ij,n) = c0
00860 krdg (ij,n) = c1
00861 enddo
00862 enddo
00863
00864
00865
00866
00867
00868
00869
00870
00871
00872
00873
00874
00875
00876
00877 do ij = 1, icells
00878 i = indxi(ij)
00879 j = indxj(ij)
00880 if (aice0(i,j) > puny) then
00881 Gsum(ij,0) = aice0(i,j)
00882 else
00883 Gsum(ij,0) = Gsum(ij,-1)
00884 endif
00885 enddo
00886
00887 do n = 1, ncat
00888
00889
00890
00891 do ij = 1, icells
00892 i = indxi(ij)
00893 j = indxj(ij)
00894 if (aicen(i,j,n) > puny) then
00895 Gsum(ij,n) = Gsum(ij,n-1) + aicen(i,j,n)
00896 else
00897 Gsum(ij,n) = Gsum(ij,n-1)
00898 endif
00899 enddo
00900 enddo
00901
00902
00903
00904 do ij = 1, icells
00905 work(ij) = c1 / Gsum(ij,ncat)
00906 enddo
00907 do n = 0, ncat
00908
00909
00910
00911 do ij = 1, icells
00912 Gsum(ij,n) = Gsum(ij,n) * work(ij)
00913 enddo
00914 enddo
00915
00916
00917
00918
00919
00920
00921
00922
00923
00924
00925
00926 if (krdg_partic == 0) then
00927
00928
00929
00930
00931
00932
00933
00934 do n = 0, ncat
00935 do ij = 1, icells
00936 if (Gsum(ij,n) < Gstar) then
00937 apartic(ij,n) = Gstari*(Gsum(ij,n)-Gsum(ij,n-1)) * &
00938 (c2 - (Gsum(ij,n-1)+Gsum(ij,n))*Gstari)
00939 elseif (Gsum(ij,n-1) < Gstar) then
00940 apartic(ij,n) = Gstari * (Gstar-Gsum(ij,n-1)) * &
00941 (c2 - (Gsum(ij,n-1)+Gstar)*Gstari)
00942 endif
00943 enddo
00944 enddo
00945
00946 elseif (krdg_partic==1) then
00947
00948
00949
00950
00951
00952
00953
00954
00955
00956
00957 xtmp = c1 / (c1 - exp(-astari))
00958
00959 do n = -1, ncat
00960
00961
00962
00963 do ij = 1, icells
00964 Gsum(ij,n) = exp(-Gsum(ij,n)*astari) * xtmp
00965 enddo
00966 enddo
00967
00968 do n = 0, ncat
00969 do ij = 1, icells
00970 apartic(ij,n) = Gsum(ij,n-1) - Gsum(ij,n)
00971 enddo
00972 enddo
00973
00974 endif
00975
00976
00977
00978
00979
00980
00981
00982
00983
00984
00985 if (krdg_redist == 0) then
00986
00987
00988
00989
00990
00991
00992
00993
00994
00995
00996
00997
00998
00999 do n = 1, ncat
01000 do ij = 1, icells
01001 i = indxi(ij)
01002 j = indxj(ij)
01003
01004 if (aicen(i,j,n) > puny) then
01005 hi = vicen(i,j,n) / aicen(i,j,n)
01006 hrmin(ij,n) = min(c2*hi, hi + maxraft)
01007 hrmax(ij,n) = c2*sqrt(Hstar*hi)
01008 hrmax(ij,n) = max(hrmax(ij,n), hrmin(ij,n)+puny)
01009 hrmean = p5 * (hrmin(ij,n) + hrmax(ij,n))
01010 krdg(ij,n) = hrmean / hi
01011 endif
01012
01013 enddo
01014 enddo
01015
01016 else
01017
01018
01019
01020
01021
01022
01023
01024
01025
01026
01027
01028
01029
01030
01031
01032
01033
01034
01035
01036
01037
01038
01039
01040
01041
01042
01043
01044
01045
01046 do n = 1, ncat
01047 do ij = 1, icells
01048 i = indxi(ij)
01049 j = indxj(ij)
01050 if (aicen(i,j,n) > puny) then
01051 hi = vicen(i,j,n) / aicen(i,j,n)
01052 hi = max(hi,puny)
01053 hrmin(ij,n) = min(c2*hi, hi + maxraft)
01054 hrexp(ij,n) = mu_rdg * sqrt(hi)
01055 krdg(ij,n) = (hrmin(ij,n) + hrexp(ij,n)) / hi
01056 endif
01057 enddo
01058 enddo
01059
01060 endif
01061
01062
01063
01064
01065
01066
01067
01068
01069 do ij = 1, icells
01070 aksum(ij) = apartic(ij,0)
01071 enddo
01072
01073 do n = 1, ncat
01074
01075
01076
01077 do ij = 1, icells
01078
01079 aksum(ij) = aksum(ij) &
01080 + apartic(ij,n) * (c1 - c1/krdg(ij,n))
01081 enddo
01082 enddo
01083
01084 end subroutine ridge_itd
01085
01086
01087
01088
01089
01090
01091
01092
01093
01094
01095
01096
01097
01098
01099
01100
01101
01102 subroutine ridge_shift (nx_block, ny_block, &
01103 icells, indxi, indxj, &
01104 dt, &
01105 aicen, trcrn, &
01106 vicen, vsnon, &
01107 eicen, esnon, &
01108 aice0, trcr_depend, &
01109 aksum, apartic, &
01110 hrmin, hrmax, &
01111 hrexp, krdg, &
01112 closing_net, opning, &
01113 ardg1, ardg2, &
01114 virdg, aopen, &
01115 msnow_mlt, esnow_mlt, &
01116 msoot, &
01117 l_stop, &
01118 istop, jstop)
01119
01120
01121
01122 use ice_state, only: nt_aero
01123
01124
01125
01126 integer (kind=int_kind), intent(in) ::
01127 nx_block, ny_block,
01128 icells
01129
01130 integer (kind=int_kind), dimension (nx_block*ny_block),
01131 intent(in) ::
01132 indxi, indxj
01133
01134 real (kind=dbl_kind), intent(in) ::
01135 dt
01136
01137 integer (kind=int_kind), dimension (ntrcr), intent(in) ::
01138 trcr_depend
01139
01140 real (kind=dbl_kind), dimension (nx_block,ny_block),
01141 intent(inout) ::
01142 aice0
01143
01144 real (kind=dbl_kind), dimension (nx_block,ny_block,ncat),
01145 intent(inout) ::
01146 aicen ,
01147 vicen ,
01148 vsnon
01149
01150 real (kind=dbl_kind), dimension (nx_block,ny_block,ntrcr,ncat),
01151 intent(inout) ::
01152 trcrn
01153
01154 real (kind=dbl_kind), dimension (nx_block,ny_block,ntilyr),
01155 intent(inout) ::
01156 eicen
01157
01158 real (kind=dbl_kind), dimension (nx_block,ny_block,ntslyr),
01159 intent(inout) ::
01160 esnon
01161
01162 real (kind=dbl_kind), dimension (icells), intent(in) ::
01163 aksum
01164
01165 real (kind=dbl_kind), dimension (icells,0:ncat), intent(in) ::
01166 apartic
01167
01168
01169 real (kind=dbl_kind), dimension (icells,ncat), intent(in) ::
01170 hrmin ,
01171 hrmax ,
01172 hrexp ,
01173 krdg
01174
01175 real (kind=dbl_kind), dimension(icells), intent(inout) ::
01176 closing_net,
01177 opning ,
01178 ardg1 ,
01179 ardg2 ,
01180 virdg ,
01181 aopen
01182
01183 real (kind=dbl_kind), dimension(icells), intent(inout) ::
01184 msnow_mlt,
01185 esnow_mlt
01186
01187 real (kind=dbl_kind), dimension(icells,n_aeromx), intent(inout) ::
01188 msoot
01189
01190 logical (kind=log_kind), intent(inout) ::
01191 l_stop
01192
01193 integer (kind=int_kind), intent(inout) ::
01194 istop, jstop
01195
01196
01197
01198 integer (kind=int_kind) ::
01199 i,j ,
01200 n, nr ,
01201 k ,
01202 it ,
01203 ij, m ,
01204 iridge
01205 integer (kind=int_kind) ::
01206 iaero
01207
01208 integer (kind=int_kind), dimension (icells) ::
01209 indxii, indxjj ,
01210 indxij
01211
01212 real (kind=dbl_kind), dimension (icells,ncat) ::
01213 aicen_init ,
01214 vicen_init ,
01215 vsnon_init
01216
01217 real (kind=dbl_kind), dimension (icells,ntilyr) ::
01218 eicen_init
01219
01220 real (kind=dbl_kind), dimension (icells,ntslyr) ::
01221 esnon_init
01222
01223 real (kind=dbl_kind), dimension(icells,ntrcr,ncat) ::
01224 atrcrn
01225
01226 real (kind=dbl_kind), dimension (icells) ::
01227 closing_gross
01228
01229
01230
01231 real (kind=dbl_kind), dimension (icells) ::
01232 afrac ,
01233 ardg1n ,
01234 ardg2n ,
01235 virdgn ,
01236 vsrdgn ,
01237 dhr ,
01238 dhr2 ,
01239 farea ,
01240 fvol
01241
01242 real (kind=dbl_kind), dimension (icells,nilyr) ::
01243 eirdgn
01244
01245 real (kind=dbl_kind), dimension (icells,nslyr) ::
01246 esrdgn
01247
01248 real (kind=dbl_kind) ::
01249 hi1 ,
01250 hexp ,
01251 hL, hR ,
01252 expL, expR ,
01253 tmpfac ,
01254 wk1
01255
01256
01257
01258
01259
01260 do n = 1, ncat
01261 do it = 1, ntrcr
01262 if (trcr_depend(it) == 0) then
01263 do ij = 1, icells
01264 i = indxi(ij)
01265 j = indxj(ij)
01266 atrcrn(ij,it,n) = aicen(i,j,n)*trcrn(i,j,it,n)
01267 enddo
01268 elseif (trcr_depend(it) == 1) then
01269 do ij = 1, icells
01270 i = indxi(ij)
01271 j = indxj(ij)
01272 atrcrn(ij,it,n) = vicen(i,j,n)*trcrn(i,j,it,n)
01273 enddo
01274 elseif (trcr_depend(it) == 2) then
01275 do ij = 1, icells
01276 i = indxi(ij)
01277 j = indxj(ij)
01278 atrcrn(ij,it,n) = vsnon(i,j,n)*trcrn(i,j,it,n)
01279 enddo
01280 endif
01281 enddo
01282 enddo
01283
01284
01285
01286
01287
01288
01289
01290
01291
01292
01293 do ij = 1, icells
01294 i = indxi(ij)
01295 j = indxj(ij)
01296
01297 closing_gross(ij) = closing_net(ij) / aksum(ij)
01298
01299
01300
01301
01302
01303
01304 if (apartic(ij,0) > c0) then
01305 wk1 = apartic(ij,0) * closing_gross(ij) * dt
01306 if (wk1 > aice0(i,j)) then
01307 tmpfac = aice0(i,j) / wk1
01308 closing_gross(ij) = closing_gross(ij) * tmpfac
01309 opning(ij) = opning(ij) * tmpfac
01310 endif
01311 endif
01312
01313 enddo
01314
01315
01316
01317
01318
01319 do n = 1, ncat
01320
01321
01322
01323 do ij = 1, icells
01324 i = indxi(ij)
01325 j = indxj(ij)
01326
01327 if (aicen(i,j,n) > puny .and. apartic(ij,n) > c0) then
01328 wk1 = apartic(ij,n) * closing_gross(ij) * dt
01329 if (wk1 > aicen(i,j,n)) then
01330 tmpfac = aicen(i,j,n) / wk1
01331 closing_gross(ij) = closing_gross(ij) * tmpfac
01332 opning(ij) = opning(ij) * tmpfac
01333 endif
01334 endif
01335
01336 enddo
01337 enddo
01338
01339
01340
01341
01342
01343
01344
01345
01346 do ij = 1, icells
01347 i = indxi(ij)
01348 j = indxj(ij)
01349 aice0(i,j) = aice0(i,j) &
01350 - apartic(ij,0)*closing_gross(ij)*dt &
01351 + opning(ij)*dt
01352 if (aice0(i,j) < -puny) then
01353 l_stop = .true.
01354 istop = i
01355 jstop = j
01356
01357 write (nu_diag,*) ' '
01358 write (nu_diag,*) 'Ridging error: aice0 < 0'
01359 write (nu_diag,*) 'i, j, aice0:', i, j, aice0(i,j)
01360 return
01361
01362 elseif (aice0(i,j) < c0) then
01363 aice0(i,j) = c0
01364 endif
01365
01366 aopen(ij) = opning(ij)*dt
01367
01368 enddo
01369
01370
01371
01372
01373
01374 do n = 1, ncat
01375 do ij = 1, icells
01376 i = indxi(ij)
01377 j = indxj(ij)
01378 aicen_init(ij,n) = aicen(i,j,n)
01379 vicen_init(ij,n) = vicen(i,j,n)
01380 vsnon_init(ij,n) = vsnon(i,j,n)
01381 enddo
01382 enddo
01383
01384 do n = 1, ntilyr
01385 do ij = 1, icells
01386 i = indxi(ij)
01387 j = indxj(ij)
01388 eicen_init(ij,n) = eicen(i,j,n)
01389 enddo
01390 enddo
01391
01392 do n = 1, ntslyr
01393 do ij = 1, icells
01394 i = indxi(ij)
01395 j = indxj(ij)
01396 esnon_init(ij,n) = esnon(i,j,n)
01397 enddo
01398 enddo
01399
01400
01401
01402
01403
01404
01405 do n = 1, ncat
01406
01407
01408
01409
01410
01411 iridge = 0
01412 do ij = 1, icells
01413 i = indxi(ij)
01414 j = indxj(ij)
01415 if (aicen_init(ij,n) > puny .and. apartic(ij,n) > c0 &
01416 .and. closing_gross(ij) > c0) then
01417 iridge = iridge + 1
01418 indxii(iridge) = i
01419 indxjj(iridge) = j
01420 indxij(iridge) = ij
01421 endif
01422 enddo
01423
01424
01425
01426
01427 do ij = 1, iridge
01428 i = indxii(ij)
01429 j = indxjj(ij)
01430 m = indxij(ij)
01431
01432
01433
01434
01435
01436
01437
01438 ardg1n(ij) = apartic(m,n)*closing_gross(m)*dt
01439
01440 if (ardg1n(ij) > aicen_init(m,n) + puny) then
01441 l_stop = .true.
01442 istop = i
01443 jstop = j
01444
01445 write (nu_diag,*) ' '
01446 write (nu_diag,*) 'Ridging error: ardg > aicen'
01447 write (nu_diag,*) 'i, j, n:', i, j, n
01448 write (nu_diag,*) 'ardg, aicen:', &
01449 ardg1n(ij), aicen_init(m,n)
01450 return
01451 else
01452 ardg1n(ij) = min(aicen_init(m,n), ardg1n(ij))
01453 endif
01454
01455 ardg2n(ij) = ardg1n(ij) / krdg(m,n)
01456 afrac(ij) = ardg1n(ij) / aicen_init(m,n)
01457
01458
01459
01460
01461
01462
01463 virdgn(ij) = vicen_init(m,n) * afrac(ij)
01464 vsrdgn(ij) = vsnon_init(m,n) * afrac(ij)
01465
01466 aicen(i,j,n) = aicen(i,j,n) - ardg1n(ij)
01467 vicen(i,j,n) = vicen(i,j,n) - virdgn(ij)
01468 vsnon(i,j,n) = vsnon(i,j,n) - vsrdgn(ij)
01469
01470
01471
01472
01473
01474 ardg1(m) = ardg1(m) + ardg1n(ij)
01475 ardg2(m) = ardg2(m) + ardg2n(ij)
01476 virdg(m) = virdg(m) + virdgn(ij)
01477
01478
01479
01480
01481
01482 msnow_mlt(m) = msnow_mlt(m) + rhos*vsrdgn(ij)*(c1-fsnowrdg)
01483
01484
01485
01486
01487
01488 if (n_aero >= 1) then
01489 do iaero=1,n_aero
01490 msoot(m,iaero) = msoot(m,iaero) &
01491 + vsrdgn(ij)*(c1-fsnowrdg) &
01492 *(trcrn(i,j,nt_aero +4*(iaero-1),n) &
01493 + trcrn(i,j,nt_aero+1+4*(iaero-1),n))
01494 enddo
01495 endif
01496
01497
01498
01499
01500
01501
01502 dhr(ij) = hrmax(ij,n) - hrmin(m,n)
01503 dhr2(ij) = hrmax(ij,n) * hrmax(ij,n) - hrmin(m,n) * hrmin(m,n)
01504
01505 enddo
01506
01507
01508
01509
01510
01511 do k = 1, nilyr
01512
01513
01514
01515 do ij = 1, iridge
01516 i = indxii(ij)
01517 j = indxjj(ij)
01518 m = indxij(ij)
01519
01520 eirdgn(ij,k) = eicen_init(m,ilyr1(n)+k-1) * afrac(ij)
01521 eicen(i,j,ilyr1(n)+k-1) = eicen (i,j,ilyr1(n)+k-1) &
01522 - eirdgn(ij,k)
01523 enddo
01524 enddo
01525
01526
01527
01528
01529
01530
01531
01532 do k = 1, nslyr
01533
01534
01535
01536 do ij = 1, iridge
01537 i = indxii(ij)
01538 j = indxjj(ij)
01539 m = indxij(ij)
01540
01541 esrdgn(ij,k) = esnon_init(m,slyr1(n)+k-1) * afrac(ij)
01542 esnon(i,j,slyr1(n)+k-1) = esnon (i,j,slyr1(n)+k-1) &
01543 - esrdgn(ij,k)
01544 esnow_mlt(m) = esnow_mlt(m) &
01545 + esrdgn(ij,k)*(c1-fsnowrdg)
01546 enddo
01547 enddo
01548
01549
01550
01551
01552
01553 do it = 1, ntrcr
01554 if (trcr_depend(it) == 0) then
01555
01556
01557
01558 do ij = 1, iridge
01559 i = indxii(ij)
01560 j = indxjj(ij)
01561 m = indxij(ij)
01562 atrcrn(m,it,n) = atrcrn(m,it,n) &
01563 - ardg1n(ij)*trcrn(i,j,it,n)
01564 enddo
01565
01566 elseif (trcr_depend(it) == 1) then
01567
01568
01569
01570 do ij = 1, iridge
01571 i = indxii(ij)
01572 j = indxjj(ij)
01573 m = indxij(ij)
01574 atrcrn(m,it,n) = atrcrn(m,it,n) &
01575 - virdgn(ij)*trcrn(i,j,it,n)
01576 enddo
01577
01578 elseif (trcr_depend(it) == 2) then
01579
01580
01581
01582 do ij = 1, iridge
01583 i = indxii(ij)
01584 j = indxjj(ij)
01585 m = indxij(ij)
01586 atrcrn(m,it,n) = atrcrn(m,it,n) &
01587 - vsrdgn(ij)*trcrn(i,j,it,n)
01588 enddo
01589 endif
01590 enddo
01591
01592
01593
01594
01595
01596
01597 do nr = 1, ncat
01598
01599 if (krdg_redist == 0) then
01600
01601 do ij = 1, iridge
01602 m = indxij(ij)
01603
01604
01605
01606
01607
01608
01609 if (hrmin(m,n) >= hin_max(nr) .or. &
01610 hrmax(ij,n) <= hin_max(nr-1)) then
01611 hL = c0
01612 hR = c0
01613 else
01614 hL = max (hrmin(m,n), hin_max(nr-1))
01615 hR = min (hrmax(ij,n), hin_max(nr))
01616 endif
01617
01618 farea(ij) = (hR-hL) / dhr(ij)
01619 fvol (ij) = (hR*hR - hL*hL) / dhr2(ij)
01620
01621 enddo
01622
01623 else
01624
01625
01626
01627
01628
01629
01630 if (nr < ncat) then
01631
01632 do ij = 1, iridge
01633 m = indxij(ij)
01634
01635 hi1 = hrmin(m,n)
01636 hexp = hrexp(m,n)
01637
01638 if (hi1 >= hin_max(nr)) then
01639 farea(ij) = c0
01640 fvol (ij) = c0
01641 else
01642 hL = max (hi1, hin_max(nr-1))
01643 hR = hin_max(nr)
01644 expL = exp(-(hL-hi1)/hexp)
01645 expR = exp(-(hR-hi1)/hexp)
01646 farea(ij) = expL - expR
01647 fvol (ij) = ((hL + hexp)*expL &
01648 - (hR + hexp)*expR) / (hi1 + hexp)
01649 endif
01650 enddo
01651
01652 else
01653
01654 do ij = 1, iridge
01655 m = indxij(ij)
01656
01657 hi1 = hrmin(m,n)
01658 hexp = hrexp(m,n)
01659
01660 hL = max (hi1, hin_max(nr-1))
01661 expL = exp(-(hL-hi1)/hexp)
01662 farea(ij) = expL
01663 fvol (ij) = (hL + hexp)*expL / (hi1 + hexp)
01664
01665 enddo
01666
01667 endif
01668
01669 endif
01670
01671
01672
01673
01674
01675
01676
01677
01678 do ij = 1, iridge
01679 i = indxii(ij)
01680 j = indxjj(ij)
01681 aicen(i,j,nr) = aicen(i,j,nr) + farea(ij)*ardg2n(ij)
01682 vicen(i,j,nr) = vicen(i,j,nr) + fvol(ij) *virdgn(ij)
01683 vsnon(i,j,nr) = vsnon(i,j,nr) &
01684 + fvol(ij)*vsrdgn(ij)*fsnowrdg
01685 enddo
01686
01687
01688
01689
01690 do k = 1, nilyr
01691
01692
01693
01694 do ij = 1, iridge
01695 i = indxii(ij)
01696 j = indxjj(ij)
01697 eicen(i,j,ilyr1(nr)+k-1) = eicen(i,j,ilyr1(nr)+k-1) &
01698 + fvol(ij)*eirdgn(ij,k)
01699 enddo
01700 enddo
01701
01702
01703
01704
01705 do k = 1, nslyr
01706
01707
01708
01709 do ij = 1, iridge
01710 i = indxii(ij)
01711 j = indxjj(ij)
01712 esnon(i,j,slyr1(nr)+k-1) = esnon(i,j,slyr1(nr)+k-1) &
01713 + fvol(ij)*esrdgn(ij,k)*fsnowrdg
01714 enddo
01715 enddo
01716
01717
01718
01719
01720
01721
01722
01723
01724
01725 do it = 1, ntrcr
01726 if (trcr_depend(it) == 0) then
01727
01728
01729
01730 do ij = 1, iridge
01731 i = indxii(ij)
01732 j = indxjj(ij)
01733 m = indxij(ij)
01734 atrcrn(m,it,nr) = atrcrn(m,it,nr) &
01735 + farea(ij)*ardg2n(ij)*trcrn(i,j,it,n)
01736
01737 enddo
01738 elseif (trcr_depend(it) == 1) then
01739
01740
01741
01742 do ij = 1, iridge
01743 i = indxii(ij)
01744 j = indxjj(ij)
01745 m = indxij(ij)
01746 atrcrn(m,it,nr) = atrcrn(m,it,nr) &
01747 + fvol(ij)*virdgn(ij)*trcrn(i,j,it,n)
01748
01749 enddo
01750 elseif (trcr_depend(it) == 2) then
01751
01752
01753
01754 do ij = 1, iridge
01755 i = indxii(ij)
01756 j = indxjj(ij)
01757 m = indxij(ij)
01758 atrcrn(m,it,nr) = atrcrn(m,it,nr) &
01759 + fvol(ij)*vsrdgn(ij)*fsnowrdg*trcrn(i,j,it,n)
01760
01761 enddo
01762 endif
01763 enddo
01764
01765 enddo
01766 enddo
01767
01768
01769
01770
01771
01772
01773 do n = 1, ncat
01774 call compute_tracers (nx_block, ny_block, &
01775 icells, indxi, indxj, &
01776 trcr_depend, &
01777 atrcrn(:,:,n), aicen(:,:, n), &
01778 vicen (:,:, n), vsnon(:,:, n), &
01779 trcrn(:,:,:,n))
01780 enddo
01781
01782 end subroutine ridge_shift
01783
01784
01785
01786
01787
01788
01789
01790
01791
01792
01793
01794
01795
01796
01797 subroutine ridge_check (nx_block, ny_block, &
01798 icells, indxi, indxj, &
01799 dt, &
01800 asum, closing_net, &
01801 divu_adv, opning, &
01802 iterate_ridging)
01803
01804
01805
01806
01807
01808 integer (kind=int_kind), intent(in) ::
01809 nx_block, ny_block,
01810 icells
01811
01812 integer (kind=int_kind), dimension (nx_block*ny_block),
01813 intent(in) ::
01814 indxi, indxj
01815
01816
01817 real (kind=dbl_kind), intent(in) ::
01818 dt
01819
01820 real (kind=dbl_kind), dimension (icells), intent(in) ::
01821 asum
01822
01823 real (kind=dbl_kind), dimension (icells),
01824 intent(inout) ::
01825 closing_net,
01826 divu_adv ,
01827 opning
01828
01829 logical (kind=log_kind), intent(out) ::
01830 iterate_ridging
01831
01832
01833
01834 integer (kind=int_kind) ::
01835 ij
01836
01837 iterate_ridging = .false.
01838
01839 do ij = 1, icells
01840 if (abs(asum(ij) - c1) < puny) then
01841 closing_net(ij) = c0
01842 opning (ij) = c0
01843 else
01844 iterate_ridging = .true.
01845 divu_adv(ij) = (c1 - asum(ij)) / dt
01846 closing_net(ij) = max(c0, -divu_adv(ij))
01847 opning(ij) = max(c0, divu_adv(ij))
01848 endif
01849 enddo
01850
01851 end subroutine ridge_check
01852
01853
01854
01855
01856
01857
01858
01859
01860
01861
01862
01863
01864
01865
01866
01867
01868
01869
01870
01871
01872
01873
01874
01875
01876
01877
01878 subroutine ice_strength (nx_block, ny_block, &
01879 ilo, ihi, jlo, jhi, &
01880 icells, &
01881 indxi, indxj, &
01882 aice, vice, &
01883 aice0, aicen, &
01884 vicen, strength)
01885
01886
01887
01888
01889
01890 integer (kind=int_kind), intent(in) ::
01891 nx_block, ny_block,
01892 ilo,ihi,jlo,jhi
01893
01894 integer (kind=int_kind), intent(in) ::
01895 icells
01896
01897 integer (kind=int_kind), dimension (nx_block*ny_block),
01898 intent(in) ::
01899 indxi ,
01900 indxj
01901
01902 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(in) ::
01903 aice ,
01904 vice ,
01905 aice0
01906
01907 real (kind=dbl_kind), dimension (nx_block,ny_block,ncat),
01908 intent(in) ::
01909 aicen ,
01910 vicen
01911
01912 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(out) ::
01913 strength
01914
01915
01916
01917
01918
01919 real (kind=dbl_kind), dimension (icells) ::
01920 asum ,
01921 aksum
01922
01923 real (kind=dbl_kind), dimension (icells,0:ncat) ::
01924 apartic
01925
01926
01927 real (kind=dbl_kind), dimension (icells,ncat) ::
01928 hrmin ,
01929 hrmax ,
01930 hrexp ,
01931 krdg
01932
01933 integer (kind=int_kind) ::
01934 i,j ,
01935 n ,
01936 ij
01937
01938 real (kind=dbl_kind) ::
01939 hi ,
01940 h2rdg ,
01941 dh2rdg
01942
01943
01944
01945
01946
01947
01948 strength(:,:) = c0
01949
01950 if (kstrength == 1) then
01951
01952
01953
01954
01955
01956 call asum_ridging (nx_block, ny_block, &
01957 icells, indxi, indxj, &
01958 aicen, aice0, &
01959 asum)
01960
01961 call ridge_itd (nx_block, ny_block, &
01962 icells, indxi, indxj, &
01963 aicen, vicen, &
01964 aice0, &
01965 aksum, apartic, &
01966 hrmin, hrmax, &
01967 hrexp, krdg)
01968
01969
01970
01971
01972
01973
01974 if (krdg_redist==0) then
01975
01976 do n = 1, ncat
01977
01978
01979
01980 do ij = 1, icells
01981 i = indxi(ij)
01982 j = indxj(ij)
01983 if (aicen(i,j,n) > puny .and. apartic(ij,n) > c0)then
01984 hi = vicen(i,j,n) / aicen(i,j,n)
01985 h2rdg = p333 * (hrmax(ij,n)**3 - hrmin(ij,n)**3) &
01986 / (hrmax(ij,n) - hrmin(ij,n))
01987 dh2rdg = -hi*hi + h2rdg/krdg(ij,n)
01988 strength(i,j) = strength(i,j) &
01989 + apartic(ij,n) * dh2rdg
01990 endif
01991 enddo
01992 enddo
01993
01994 elseif (krdg_redist==1) then
01995
01996 do n = 1, ncat
01997
01998
01999
02000 do ij = 1, icells
02001 i = indxi(ij)
02002 j = indxj(ij)
02003 if (aicen(i,j,n) > puny .and. apartic(ij,n) > c0)then
02004 hi = vicen(i,j,n) / aicen(i,j,n)
02005 h2rdg = hrmin(ij,n)*hrmin(ij,n) &
02006 + c2*hrmin(ij,n)*hrexp(ij,n) &
02007 + c2*hrexp(ij,n)*hrexp(ij,n)
02008 dh2rdg = -hi*hi + h2rdg/krdg(ij,n)
02009 strength(i,j) = strength(i,j) &
02010 + apartic(ij,n) * dh2rdg
02011 endif
02012 enddo
02013 enddo
02014
02015 endif
02016
02017
02018
02019
02020 do ij = 1, icells
02021 i = indxi(ij)
02022 j = indxj(ij)
02023 strength(i,j) = Cf * Cp * strength(i,j) / aksum(ij)
02024
02025
02026 enddo
02027
02028 else
02029
02030
02031
02032
02033
02034 do j = jlo, jhi
02035 do i = ilo, ihi
02036 strength(i,j) = Pstar*vice(i,j)*exp(-Cstar*(c1-aice(i,j)))
02037 enddo
02038 enddo
02039
02040 endif
02041
02042 end subroutine ice_strength
02043
02044
02045
02046 end module ice_mechred
02047
02048