00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 module ice_meltpond
00018
00019
00020
00021 use ice_kinds_mod
00022 use ice_constants
00023 use ice_fileunits
00024 use ice_read_write
00025 use ice_restart, only: lenstr, restart_dir, restart_file, &
00026 pointer_file, runtype
00027 use ice_communicate, only: my_task, master_task
00028 use ice_exit, only: abort_ice
00029
00030
00031
00032 implicit none
00033
00034 logical (kind=log_kind) ::
00035 restart_pond
00036
00037
00038
00039 contains
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054 subroutine init_meltponds
00055
00056
00057
00058 use ice_domain_size
00059 use ice_blocks
00060 use ice_domain
00061 use ice_flux
00062 use ice_state
00063
00064
00065
00066
00067
00068
00069
00070 integer (kind=int_kind) ::
00071 icells
00072
00073 integer (kind=int_kind), dimension(nx_block*ny_block) ::
00074 indxi, indxj
00075
00076 integer (kind=int_kind) :: i, j, ij, n, iblk
00077
00078 if (trim(filename_volpn) /= 'none') restart_pond = .true.
00079
00080 if (restart_pond) then
00081 if (trim(runtype) == 'continue') then
00082 call read_restart_pond
00083 else
00084 call read_restart_pond(filename_volpn)
00085 endif
00086 endif
00087
00088 end subroutine init_meltponds
00089
00090
00091
00092
00093
00094
00095
00096
00097 subroutine compute_ponds(nx_block,ny_block, &
00098 ilo, ihi, jlo, jhi, &
00099 meltt, melts, frain, &
00100 aicen, vicen, vsnon, &
00101 trcrn, apondn, hpondn)
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111 use ice_state, only: nt_Tsfc, nt_volpn
00112 use ice_calendar, only: dt
00113 use ice_domain_size, only: ntrcr
00114
00115 integer (kind=int_kind), intent(in) ::
00116 nx_block, ny_block,
00117 ilo,ihi,jlo,jhi
00118
00119 real (kind=dbl_kind), dimension(nx_block,ny_block),
00120 intent(in) ::
00121 meltt,
00122 melts,
00123 frain,
00124 aicen,
00125 vicen,
00126 vsnon
00127
00128 real (kind=dbl_kind), dimension(nx_block,ny_block),
00129 intent(inout) ::
00130 apondn,
00131 hpondn
00132
00133 real (kind=dbl_kind), dimension(nx_block,ny_block,ntrcr),
00134 intent(inout) ::
00135 trcrn
00136
00137 real (kind=dbl_kind), dimension(nx_block,ny_block) ::
00138 volpn,
00139 Tsfcn
00140
00141 integer (kind=int_kind), dimension (nx_block*ny_block) ::
00142 indxi, indxj
00143
00144 integer (kind=int_kind) :: i,j,ij,icells
00145
00146 real (kind=dbl_kind) :: hi,hs,dTs,rfrac,asnow
00147
00148 real (kind=dbl_kind), parameter ::
00149 hi_min = p1
00150
00151 Tsfcn(:,:) = trcrn(:,:,nt_Tsfc)
00152 volpn(:,:) = trcrn(:,:,nt_volpn)
00153
00154
00155
00156
00157
00158 icells = 0
00159 do j = jlo, jhi
00160 do i = ilo, ihi
00161 if (aicen(i,j) > puny) then
00162 icells = icells + 1
00163 indxi(icells) = i
00164 indxj(icells) = j
00165 endif
00166 enddo
00167 enddo
00168
00169 do ij = 1, icells
00170 i = indxi(ij)
00171 j = indxj(ij)
00172
00173 hi = vicen(i,j)/aicen(i,j)
00174 hs = vsnon(i,j)/aicen(i,j)
00175 dTs = Timelt - Tsfcn(i,j)
00176
00177
00178 rfrac = 0.85_dbl_kind - 0.7_dbl_kind*aicen(i,j)
00179
00180 volpn(i,j) = volpn(i,j) + (c1-rfrac) &
00181 * (meltt(i,j)*(rhoi/rhofresh) + melts(i,j)*(rhos/rhofresh) &
00182 + frain(i,j)*dt/rhofresh)
00183
00184
00185 if (Tsfcn(i,j) .lt. Timelt-c2) then
00186 volpn(i,j) = volpn(i,j) &
00187 *exp(p01*(dTs-c2)/(Timelt-c2))
00188 endif
00189
00190 volpn(i,j) = max(volpn(i,j), c0)
00191
00192 apondn(i,j) = min ( sqrt ( volpn(i,j) * 1.25 ), c1)
00193 hpondn(i,j) = 0.8 * apondn(i,j)
00194
00195
00196 if (hpondn(i,j) .gt. 0.9*hi) then
00197 hpondn(i,j) = 0.9*hi
00198 volpn(i,j) = hpondn(i,j)*apondn(i,j)
00199
00200 endif
00201
00202
00203
00204
00205
00206
00207 if (hs >= hsmin) then
00208 asnow = min( hs/hs0, c1 )
00209 else
00210 asnow = c0
00211 endif
00212
00213
00214
00215
00216
00217
00218
00219
00220
00221
00222 if ( hs .gt. puny ) then
00223 apondn(i,j) = (c1-asnow)*apondn(i,j)
00224 hpondn(i,j) = 0.8_dbl_kind * apondn(i,j)
00225 endif
00226
00227
00228 if (hi .lt. hi_min) then
00229 apondn(i,j) = c0
00230 hpondn(i,j) = c0
00231 volpn(i,j) = c0
00232 endif
00233
00234 trcrn(i,j,nt_volpn) = volpn(i,j)
00235
00236 enddo
00237
00238 end subroutine compute_ponds
00239
00240
00241
00242
00243
00244
00245
00246
00247
00248 subroutine write_restart_pond(filename_spec)
00249
00250
00251
00252
00253
00254
00255
00256
00257
00258
00259
00260
00261 use ice_domain_size
00262 use ice_calendar, only: sec, month, mday, nyr, istep1, &
00263 time, time_forc, idate, year_init
00264 use ice_state
00265
00266
00267
00268 character(len=char_len_long), intent(in), optional :: filename_spec
00269
00270
00271
00272 integer (kind=int_kind) ::
00273 i, j, k, n, it, iblk,
00274 iyear, imonth, iday
00275
00276 character(len=char_len_long) :: filename
00277
00278 logical (kind=log_kind) :: diag
00279
00280
00281 if (present(filename_spec)) then
00282 filename = trim(filename_spec)
00283 else
00284 iyear = nyr + year_init - 1
00285 imonth = month
00286 iday = mday
00287
00288 write(filename,'(a,a,a,i4.4,a,i2.2,a,i2.2,a,i5.5)') &
00289 restart_dir(1:lenstr(restart_dir)), &
00290 restart_file(1:lenstr(restart_file)),'.volpn.', &
00291 iyear,'-',month,'-',mday,'-',sec
00292 end if
00293
00294
00295 call ice_open(nu_dump_pond,filename,0)
00296
00297 if (my_task == master_task) then
00298 write(nu_dump_pond) istep1,time,time_forc
00299 write(nu_diag,*) 'Writing ',filename(1:lenstr(filename))
00300 endif
00301
00302 diag = .true.
00303
00304
00305
00306 do n = 1, ncat
00307 call ice_write(nu_dump_pond,0,trcrn(:,:,nt_volpn,n,:),'ruf8',diag)
00308 call ice_write(nu_dump_pond,0,apondn(:,:,n,:),'ruf8',diag)
00309 call ice_write(nu_dump_pond,0,hpondn(:,:,n,:),'ruf8',diag)
00310 enddo
00311
00312 if (my_task == master_task) close(nu_dump_pond)
00313
00314 end subroutine write_restart_pond
00315
00316
00317
00318
00319
00320
00321
00322
00323 subroutine read_restart_pond(filename_spec)
00324
00325
00326
00327
00328
00329
00330
00331
00332
00333
00334
00335
00336 use ice_domain_size
00337 use ice_calendar, only: sec, month, mday, nyr, istep1, &
00338 time, time_forc, idate, year_init
00339 use ice_state
00340
00341
00342
00343 character(len=char_len_long), intent(in), optional :: filename_spec
00344
00345
00346
00347 integer (kind=int_kind) ::
00348 i, j, k, n, it, iblk,
00349 iyear, imonth, iday
00350
00351 character(len=char_len_long) ::
00352 filename, filename0, string1, string2
00353
00354 logical (kind=log_kind) ::
00355 diag
00356
00357
00358 if (my_task == master_task) then
00359
00360 if (present(filename_spec)) then
00361 filename = filename_spec
00362 else
00363 open(nu_rst_pointer,file=pointer_file)
00364 read(nu_rst_pointer,'(a)') filename0
00365 filename = trim(filename0)
00366 close(nu_rst_pointer)
00367
00368 n = index(filename0,trim(restart_file))
00369 if (n == 0) call abort_ice('volpn restart: filename discrepancy')
00370 string1 = trim(filename0(1:n-1))
00371 string2 = trim(filename0(n+lenstr(restart_file):lenstr(filename0)))
00372 write(filename,'(a,a,a,a)') &
00373 string1(1:lenstr(string1)), &
00374 restart_file(1:lenstr(restart_file)),'.volpn', &
00375 string2(1:lenstr(string2))
00376 endif
00377 endif
00378
00379 call ice_open(nu_restart_pond,filename,0)
00380
00381 if (my_task == master_task) then
00382 read(nu_restart_pond) istep1,time,time_forc
00383 write(nu_diag,*) 'Reading ',filename(1:lenstr(filename))
00384 endif
00385
00386 diag = .true.
00387
00388
00389
00390 do n = 1, ncat
00391 call ice_read(nu_restart_pond,0,trcrn(:,:,nt_volpn,n,:),'ruf8',diag)
00392 call ice_read(nu_restart_pond,0,apondn(:,:,n,:),'ruf8',diag)
00393 call ice_read(nu_restart_pond,0,hpondn(:,:,n,:),'ruf8',diag)
00394 enddo
00395
00396 if (my_task == master_task) close(nu_restart_pond)
00397
00398 end subroutine read_restart_pond
00399
00400
00401
00402 end module ice_meltpond
00403
00404