00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017 module ice_aerosol
00018
00019
00020
00021 use ice_kinds_mod
00022 use ice_constants
00023 use ice_fileunits
00024 use ice_restart, only: lenstr, restart_dir, restart_file, &
00025 pointer_file, runtype
00026 use ice_communicate, only: my_task, master_task
00027 use ice_exit, only: abort_ice
00028
00029
00030
00031 implicit none
00032
00033 logical (kind=log_kind) ::
00034 restart_aero
00035
00036
00037
00038 contains
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048
00049
00050
00051
00052
00053 subroutine init_aerosol
00054
00055
00056
00057 use ice_state, only: nt_aero, trcrn, filename_aero
00058 use ice_domain_size, only: n_aero
00059
00060 integer (kind = int_kind) :: n
00061
00062
00063
00064
00065
00066 if (trim(filename_aero) /= 'none') restart_aero = .true.
00067
00068 if (restart_aero) then
00069 if (trim(runtype) == 'continue') then
00070 call read_restart_aero
00071 else
00072 call read_restart_aero(filename_aero)
00073 endif
00074 endif
00075
00076 end subroutine init_aerosol
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092 subroutine update_aerosol (nx_block, ny_block, &
00093 dt, icells, &
00094 indxi, indxj, &
00095 meltt, melts, &
00096 meltb, congel, &
00097 snoice, &
00098 fsnow, &
00099 trcrn, &
00100 aice_old, &
00101 vice_old, vsno_old, &
00102 vicen, vsnon, aicen, &
00103 faero, fsoot)
00104
00105
00106
00107 use ice_domain_size, only: ntrcr, nilyr, nslyr, n_aero, n_aeromx
00108 use ice_state, only: nt_aero
00109
00110
00111
00112 integer (kind=int_kind), intent(in) ::
00113 nx_block, ny_block,
00114 icells
00115
00116 integer (kind=int_kind), dimension (nx_block*ny_block),
00117 intent(in) ::
00118 indxi, indxj
00119
00120 real (kind=dbl_kind), intent(in) ::
00121 dt
00122
00123 real (kind=dbl_kind), dimension(nx_block,ny_block),
00124 intent(in) ::
00125 meltt,
00126 melts,
00127 meltb,
00128 congel,
00129 snoice,
00130 fsnow,
00131 vicen,
00132 vsnon,
00133 aicen,
00134 aice_old,
00135 vice_old,
00136 vsno_old
00137
00138 real (kind=dbl_kind), dimension(nx_block,ny_block,n_aeromx),
00139 intent(in) ::
00140 faero
00141
00142 real (kind=dbl_kind), dimension(nx_block,ny_block,n_aeromx),
00143 intent(inout) ::
00144 fsoot
00145
00146 real (kind=dbl_kind), dimension(nx_block,ny_block,ntrcr),
00147 intent(inout) ::
00148 trcrn
00149
00150
00151
00152
00153 integer (kind=int_kind) :: i, j, ij, k
00154 integer (kind=int_kind) :: n
00155
00156 real (kind=dbl_kind), dimension(icells) ::
00157 dzssl,
00158 dzint,
00159 dzssli,
00160 dzinti
00161
00162 real (kind=dbl_kind), dimension(icells) ::
00163 dhs_evap, dhi_evap,
00164 dhs_melts, dhs_snoice, dhi_meltt, dhi_snoice,
00165 dhi_congel, dhi_meltb
00166 real (kind=dbl_kind), dimension(icells,n_aeromx) ::
00167 aerotot, aerotot0
00168
00169 real (kind=dbl_kind) ::
00170 dzssl_new,
00171 dzint_new,
00172 dzssli_new,
00173 dzinti_new,
00174 dznew
00175
00176 real (kind=dbl_kind), dimension(nx_block,ny_block,n_aeromx,2) ::
00177 aerosno, aeroice,
00178 aerosno0, aeroice0
00179
00180 real (kind=dbl_kind), dimension(nx_block,ny_block,n_aeromx) ::
00181 fsoot_old
00182
00183 real (kind=dbl_kind) ::
00184 hs_old, hi_old, hslyr_old, hilyr_old, dhs, dhi, hs, hi,
00185 hslyr, hilyr, sloss1, sloss2
00186 real (kind=dbl_kind), dimension(n_aeromx) ::
00187 kscav, kscavsi
00188
00189
00190 real (kind=dbl_kind) ::
00191 hi_ssl, hs_ssl
00192
00193 data hs_ssl / .040_dbl_kind /
00194 data hi_ssl / .050_dbl_kind /
00195 data kscav / .03_dbl_kind, .20_dbl_kind,&
00196 .02_dbl_kind,.02_dbl_kind,.01_dbl_kind,.01_dbl_kind /
00197 data kscavsi / .03_dbl_kind, .20_dbl_kind,&
00198 .02_dbl_kind,.02_dbl_kind,.01_dbl_kind,.01_dbl_kind /
00199
00200 aerosno(:,:,:,:) = c0
00201 aeroice(:,:,:,:) = c0
00202 aerosno0(:,:,:,:) = c0
00203 aeroice0(:,:,:,:) = c0
00204 fsoot_old(:,:,:) = fsoot(:,:,:)
00205
00206 do ij = 1, icells
00207 i = indxi(ij)
00208 j = indxj(ij)
00209
00210 hs_old=vsno_old(i,j)/aice_old(i,j)
00211 hi_old=vice_old(i,j)/aice_old(i,j)
00212 hslyr_old=hs_old/real(nslyr,kind=dbl_kind)
00213 hilyr_old=hi_old/real(nilyr,kind=dbl_kind)
00214
00215 dzssl(ij)=min(hslyr_old/c2,hs_ssl)
00216 dzint(ij)=hs_old-dzssl(ij)
00217 dzssli(ij)=min(hilyr_old/c2,hi_ssl)
00218 dzinti(ij)=hi_old-dzssli(ij)
00219
00220 if (aicen(i,j) > c0) then
00221 hs = vsnon(i,j)/aicen(i,j)
00222 hi = vicen(i,j)/aicen(i,j)
00223 dhs_melts(ij)=-melts(i,j)/aicen(i,j)
00224 dhi_snoice(ij)=snoice(i,j)/aicen(i,j)
00225 dhs_snoice(ij)=dhi_snoice(ij)*rhoi/rhos
00226 dhi_meltt(ij)=-meltt(i,j)/aicen(i,j)
00227 dhi_meltb(ij)=-meltb(i,j)/aicen(i,j)
00228 dhi_congel(ij)=congel(i,j)/aicen(i,j)
00229 else
00230 hs = vsnon(i,j)/aice_old(i,j)
00231 hi = vicen(i,j)/aice_old(i,j)
00232 dhs_melts(ij)=-melts(i,j)/aice_old(i,j)
00233 dhi_snoice(ij)=snoice(i,j)/aice_old(i,j)
00234 dhs_snoice(ij)=dhi_snoice(ij)*rhoi/rhos
00235 dhi_meltt(ij)=-meltt(i,j)/aice_old(i,j)
00236 dhi_meltb(ij)=-meltb(i,j)/aice_old(i,j)
00237 dhi_congel(ij)=congel(i,j)/aice_old(i,j)
00238 endif
00239
00240 dhs_evap(ij)=hs-(hs_old+dhs_melts(ij)-dhs_snoice(ij)+&
00241 fsnow(i,j)/rhos*dt)
00242 dhi_evap(ij)=hi-(hi_old+dhi_meltt(ij)+dhi_meltb(ij)+ &
00243 dhi_congel(ij)+dhi_snoice(ij))
00244 enddo
00245
00246 do ij = 1, icells
00247 i = indxi(ij)
00248 j = indxj(ij)
00249 do k=1,n_aero
00250 aerosno(i,j,k,:)=&
00251 trcrn(i,j,nt_aero+(k-1)*4 :nt_aero+(k-1)*4+1)*vsno_old(i,j)
00252 aeroice(i,j,k,:)=&
00253 trcrn(i,j,nt_aero+(k-1)*4+2:nt_aero+(k-1)*4+3)*vice_old(i,j)
00254 aerosno0(i,j,k,:)=aerosno(i,j,k,:)
00255 aeroice0(i,j,k,:)=aeroice(i,j,k,:)
00256 aerotot0(ij,k)=aerosno(i,j,k,2)+aerosno(i,j,k,1) &
00257 +aeroice(i,j,k,2)+aeroice(i,j,k,1)
00258 enddo
00259 enddo
00260
00261
00262 do ij = 1, icells
00263 i = indxi(ij)
00264 j = indxj(ij)
00265 dzint(ij)=dzint(ij) + min(dzssl(ij)+dhs_evap(ij),c0)
00266 dzssl(ij)=max(dzssl(ij)+dhs_evap(ij),c0)
00267 dzinti(ij)=dzinti(ij) + min(dzssli(ij)+dhi_evap(ij),c0)
00268 dzssli(ij)=max(dzssli(ij)+dhi_evap(ij),c0)
00269 enddo
00270
00271
00272 do ij = 1, icells
00273 i = indxi(ij)
00274 j = indxj(ij)
00275 dzinti(ij)=dzinti(ij)+dhi_congel(ij)
00276 enddo
00277
00278
00279 do ij = 1, icells
00280 i = indxi(ij)
00281 j = indxj(ij)
00282 if (-dhs_melts(ij) > puny) then
00283 do k=1,n_aero
00284 sloss1=c0
00285 sloss2=c0
00286 if (dzssl(ij) > puny) &
00287 sloss1=kscav(k)*aerosno(i,j,k,1) &
00288 *min(-dhs_melts(ij),dzssl(ij))/dzssl(ij)
00289 aerosno(i,j,k,1)=aerosno(i,j,k,1)-sloss1
00290 if (dzint(ij) > puny) &
00291 sloss2=kscav(k)*aerosno(i,j,k,2) &
00292 *max(-dhs_melts(ij)-dzssl(ij),c0)/dzint(ij)
00293 aerosno(i,j,k,2)=aerosno(i,j,k,2)-sloss2
00294 fsoot(i,j,k)=fsoot(i,j,k)+(sloss1+sloss2)/dt
00295 enddo
00296
00297
00298 dzint(ij)=dzint(ij)+min(dzssl(ij)+dhs_melts(ij),c0)
00299 dzssl(ij)=max(dzssl(ij)+dhs_melts(ij),c0)
00300
00301 if ( dzssl(ij) <= puny ) then
00302 aerosno(i,j,:,2)=aerosno(i,j,:,1)+aerosno(i,j,:,2)
00303 aerosno(i,j,:,1)=c0
00304 dzssl(ij)=max(dzssl(ij),c0)
00305 endif
00306 if (dzint(ij) <= puny ) then
00307 aeroice(i,j,:,1)=&
00308 aeroice(i,j,:,1)+aerosno(i,j,:,1)+aerosno(i,j,:,2)
00309 aerosno(i,j,:,:)=c0
00310 dzint(ij)=max(dzint(ij),c0)
00311 endif
00312 endif
00313 enddo
00314
00315
00316 do ij = 1, icells
00317 i = indxi(ij)
00318 j = indxj(ij)
00319 if (-dhi_meltt(ij) > puny) then
00320 do k=1,n_aero
00321 sloss1=c0
00322 sloss2=c0
00323 if (dzssli(ij) > puny) &
00324 sloss1=kscav(k)*aeroice(i,j,k,1) &
00325 *min(-dhi_meltt(ij),dzssli(ij))/dzssli(ij)
00326 aeroice(i,j,k,1)=aeroice(i,j,k,1)-sloss1
00327 if (dzinti(ij) > puny) &
00328 sloss2=kscav(k)*aeroice(i,j,k,2) &
00329 *max(-dhi_meltt(ij)-dzssli(ij),c0)/dzinti(ij)
00330 aeroice(i,j,k,2)=aeroice(i,j,k,2)-sloss2
00331 fsoot(i,j,k)=fsoot(i,j,k)+(sloss1+sloss2)/dt
00332 enddo
00333
00334 dzinti(ij)=dzinti(ij)+min(dzssli(ij)+dhi_meltt(ij),c0)
00335 dzssli(ij)=max(dzssli(ij)+dhi_meltt(ij),c0)
00336 if (dzssli(ij) <= puny) then
00337 do k=1,n_aero
00338 aeroice(i,j,k,2)=aeroice(i,j,k,1)+aeroice(i,j,k,2)
00339 aeroice(i,j,k,1)=c0
00340 enddo
00341 dzssli(ij)=max(dzssli(ij),c0)
00342 endif
00343 if (dzinti(ij) <= puny) then
00344 do k=1,n_aero
00345 fsoot(i,j,k)=fsoot(i,j,k) &
00346 +(aeroice(i,j,k,1)+aeroice(i,j,k,2))/dt
00347 aeroice(i,j,k,:)=c0
00348 enddo
00349 dzinti(ij)=max(dzinti(ij),c0)
00350 endif
00351 endif
00352 enddo
00353
00354
00355 do ij = 1, icells
00356 i = indxi(ij)
00357 j = indxj(ij)
00358 if (-dhi_meltb(ij) > puny) then
00359 do k=1,n_aero
00360 sloss1=c0
00361 sloss2=c0
00362 if (dzssli(ij) > puny) &
00363 sloss1=max(-dhi_meltb(ij)-dzinti(ij),c0) &
00364 *aeroice(i,j,k,1)/dzssli(ij)
00365 aeroice(i,j,k,1)=aeroice(i,j,k,1)-sloss1
00366 if (dzinti(ij) > puny) &
00367 sloss2=min(-dhi_meltb(ij),dzinti(ij)) &
00368 *aeroice(i,j,k,2)/dzinti(ij)
00369 aeroice(i,j,k,2)=aeroice(i,j,k,2)-sloss2
00370 fsoot(i,j,k)=fsoot(i,j,k)+(sloss1+sloss2)/dt
00371 enddo
00372
00373 dzssli(ij) = dzssli(ij)+min(dzinti(ij)+dhi_meltb(ij), c0)
00374 dzinti(ij) = max(dzinti(ij)+dhi_meltb(ij), c0)
00375 endif
00376 enddo
00377
00378
00379 do ij = 1, icells
00380 i = indxi(ij)
00381 j = indxj(ij)
00382 if (fsnow(i,j) > c0) &
00383 dzssl(ij)=dzssl(ij)+fsnow(i,j)/rhos*dt
00384 enddo
00385
00386
00387 do ij = 1, icells
00388 i = indxi(ij)
00389 j = indxj(ij)
00390 if (dhs_snoice(ij) > puny) then
00391 do k=1,n_aero
00392 sloss1=c0
00393 sloss2=c0
00394 if (dzint(ij) > puny) &
00395 sloss2 = min(dhs_snoice(ij),dzint(ij)) &
00396 *aerosno(i,j,k,2)/dzint(ij)
00397 aerosno(i,j,k,2) = aerosno(i,j,k,2) - sloss2
00398 if (dzssl(ij) > puny) &
00399 sloss1 = max(dhs_snoice(ij)-dzint(ij),c0) &
00400 *aerosno(i,j,k,1)/dzssl(ij)
00401 aerosno(i,j,k,1) = aerosno(i,j,k,1) - sloss1
00402 aeroice(i,j,k,1) = aeroice(i,j,k,1) &
00403 + (c1-kscavsi(k))*(sloss2+sloss1)
00404 fsoot(i,j,k)=fsoot(i,j,k)+kscavsi(k)*(sloss2+sloss1)/dt
00405 enddo
00406 dzssl(ij)=dzssl(ij)-max(dhs_snoice(ij)-dzint(ij),c0)
00407 dzint(ij)=max(dzint(ij)-dhs_snoice(ij),c0)
00408 dzssli(ij)=dzssli(ij)+dhi_snoice(ij)
00409 endif
00410 enddo
00411
00412
00413 do ij = 1, icells
00414 i = indxi(ij)
00415 j = indxj(ij)
00416 if (aicen(i,j) > c0) then
00417 hs = vsnon(i,j) / aicen(i,j)
00418 else
00419 hs = c0
00420 endif
00421 if (hs > hsmin) then
00422
00423 do k=1,n_aero
00424 aerosno(i,j,k,1)=aerosno(i,j,k,1) &
00425 + faero(i,j,k)*dt*aicen(i,j)
00426 enddo
00427 else
00428 do k=1,n_aero
00429 aeroice(i,j,k,1)=aeroice(i,j,k,1) &
00430 + faero(i,j,k)*dt*aicen(i,j)
00431 enddo
00432 endif
00433 enddo
00434
00435
00436 do ij = 1, icells
00437 i = indxi(ij)
00438 j = indxj(ij)
00439 if (aicen(i,j) > c0) then
00440 hs = vsnon(i,j) / aicen(i,j)
00441 hi = vicen(i,j) / aicen(i,j)
00442 else
00443 hs = c0
00444 hi = c0
00445 endif
00446 if (dzssl(ij) <= puny) then
00447 do k=1,n_aero
00448 aerosno(i,j,k,2)=aerosno(i,j,k,2)+aerosno(i,j,k,1)
00449 aerosno(i,j,k,1)=c0
00450 enddo
00451 endif
00452 if (dzint(ij) <= puny) then
00453 do k=1,n_aero
00454 aeroice(i,j,k,1)=aeroice(i,j,k,1)+aerosno(i,j,k,2)
00455 aerosno(i,j,k,2)=c0
00456 enddo
00457 endif
00458 if (dzssli(ij) <= puny) then
00459 do k=1,n_aero
00460 aeroice(i,j,k,2)=aeroice(i,j,k,2)+aeroice(i,j,k,1)
00461 aeroice(i,j,k,1)=c0
00462 enddo
00463 endif
00464
00465 if (dzinti(ij) <= puny) then
00466 do k=1,n_aero
00467 fsoot(i,j,k)=fsoot(i,j,k)+&
00468 (aeroice(i,j,k,1)+aeroice(i,j,k,2))/dt
00469 aeroice(i,j,k,:)=c0
00470 enddo
00471 endif
00472
00473 hslyr=hs/real(nslyr,kind=dbl_kind)
00474 hilyr=hi/real(nilyr,kind=dbl_kind)
00475 dzssl_new=min(hslyr/c2,hs_ssl)
00476 dzint_new=hs-dzssl_new
00477 dzssli_new=min(hilyr/c2,hi_ssl)
00478 dzinti_new=hi-dzssli_new
00479
00480 if (hs > hsmin) then
00481 do k=1,n_aero
00482 dznew=min(dzssl_new-dzssl(ij),c0)
00483 sloss1=c0
00484 if (dzssl(ij) > puny) &
00485 sloss1=dznew*aerosno(i,j,k,1)/dzssl(ij)
00486 dznew=max(dzssl_new-dzssl(ij),c0)
00487 if (dzint(ij) > puny) &
00488 sloss1=sloss1+aerosno(i,j,k,2)*dznew/dzint(ij)
00489 aerosno(i,j,k,1) =aerosno(i,j,k,1)+sloss1
00490 aerosno(i,j,k,2) =aerosno(i,j,k,2)-sloss1
00491 enddo
00492 else
00493 aeroice(i,j,:,1)=aeroice(i,j,:,1) &
00494 +aerosno(i,j,:,1)+aerosno(i,j,:,2)
00495 aerosno(i,j,:,:) = c0
00496 endif
00497
00498 if (vicen(i,j) > puny) then
00499 do k=1,n_aero
00500 sloss2=c0
00501 dznew=min(dzssli_new-dzssli(ij),c0)
00502 if (dzssli(ij) > puny) &
00503 sloss2=dznew*aeroice(i,j,k,1)/dzssli(ij)
00504 dznew=max(dzssli_new-dzssli(ij),c0)
00505 if (dzinti(ij) > puny) &
00506 sloss2=sloss2+aeroice(i,j,k,2)*dznew/dzinti(ij)
00507 aeroice(i,j,k,1) =aeroice(i,j,k,1)+sloss2
00508 aeroice(i,j,k,2) =aeroice(i,j,k,2)-sloss2
00509 enddo
00510 else
00511 fsoot(i,j,:)=fsoot(i,j,:)+(aeroice(i,j,:,1)+aeroice(i,j,:,2))/dt
00512 aeroice(i,j,:,:) = c0
00513 endif
00514
00515 do k=1,n_aero
00516 aerotot(ij,k)=aerosno(i,j,k,2)+aerosno(i,j,k,1) &
00517 +aeroice(i,j,k,2)+aeroice(i,j,k,1)
00518 if ( ( (aerotot(ij,k)-aerotot0(ij,k)) &
00519 - ( faero(i,j,k)*aicen(i,j) &
00520 - (fsoot(i,j,k)-fsoot_old(i,j,k)) )*dt ) > 0.00001) then
00521 write(nu_diag,*) 'aerosol tracer: ',k
00522 write(nu_diag,*) 'aerotot-aerotot0 ',aerotot(ij,k)-aerotot0(ij,k)
00523 write(nu_diag,*) 'faero-fsoot ',faero(i,j,k)*aicen(i,j)*dt &
00524 -(fsoot(i,j,k)-fsoot_old(i,j,k))*dt
00525 endif
00526 enddo
00527 enddo
00528
00529
00530 do ij = 1, icells
00531 i = indxi(ij)
00532 j = indxj(ij)
00533 if (vicen(i,j) > puny) &
00534 aeroice(i,j,:,:)=aeroice(i,j,:,:)/vicen(i,j)
00535 if (vsnon(i,j) > puny) &
00536 aerosno(i,j,:,:)=aerosno(i,j,:,:)/vsnon(i,j)
00537 do k=1,n_aero
00538 do n=1,2
00539 trcrn(i,j,nt_aero+(k-1)*4+n-1)=aerosno(i,j,k,n)
00540 trcrn(i,j,nt_aero+(k-1)*4+n+1)=aeroice(i,j,k,n)
00541 enddo
00542
00543
00544
00545
00546
00547
00548
00549 enddo
00550 enddo
00551
00552 do ij = 1, icells
00553 i = indxi(ij)
00554 j = indxj(ij)
00555 if (trcrn(i,j,nt_aero) < -puny .or. trcrn(i,j,nt_aero+1) < -puny &
00556 .or. trcrn(i,j,nt_aero+2) < -puny .or. trcrn(i,j,nt_aero+3) < -puny) then
00557 write(nu_diag,*) 'MH aerosol negative in aerosol code'
00558 write(nu_diag,*) 'MH INT neg in aerosol my_task = ',&
00559 my_task &
00560 ,' printing point = ',n &
00561 ,' i and j = ',i,j
00562 write(nu_diag,*) 'MH Int Neg aero snowssl= ',aerosno0(i,j,1,1)
00563 write(nu_diag,*) 'MH Int Neg aero new snowssl= ',aerosno(i,j,1,1)
00564 write(nu_diag,*) 'MH Int Neg aero snowint= ',aerosno0(i,j,1,2)
00565 write(nu_diag,*) 'MH Int Neg aero new snowint= ',aerosno(i,j,1,2)
00566 write(nu_diag,*) 'MH Int Neg aero ice_ssl= ',aeroice0(i,j,1,1)
00567 write(nu_diag,*) 'MH Int Neg aero new ice_ssl= ',aeroice(i,j,1,1)
00568 write(nu_diag,*) 'MH Int Neg aero ice_int= ',aeroice0(i,j,1,2)
00569 write(nu_diag,*) 'MH Int Neg aero new ice_int= ',aeroice(i,j,1,2)
00570 write(nu_diag,*) 'MH Int Neg aero aicen= ',aicen(i,j)
00571 write(nu_diag,*) 'MH Int Neg aero vicen= ',vicen(i,j)
00572 write(nu_diag,*) 'MH Int Neg aero vsnon= ',vsnon(i,j)
00573 write(nu_diag,*) 'MH Int Neg aero viceold= ',vice_old(i,j)
00574 write(nu_diag,*) 'MH Int Neg aero vsnoold= ',vsno_old(i,j)
00575 write(nu_diag,*) 'MH Int Neg aero melts= ',melts(i,j)
00576 write(nu_diag,*) 'MH Int Neg aero meltt= ',meltt(i,j)
00577 write(nu_diag,*) 'MH Int Neg aero meltb= ',meltb(i,j)
00578 write(nu_diag,*) 'MH Int Neg aero congel= ',congel(i,j)
00579 write(nu_diag,*) 'MH Int Neg aero snoice= ',snoice(i,j)
00580 write(nu_diag,*) 'MH Int Neg aero evap sno?= ',dhs_evap(ij)
00581 write(nu_diag,*) 'MH Int Neg aero evap ice?= ',dhi_evap(ij)
00582 write(nu_diag,*) 'MH Int Neg aero fsnow= ',fsnow(i,j)
00583 write(nu_diag,*) 'MH Int Neg aero faero= ',faero(i,j,1)
00584 write(nu_diag,*) 'MH Int Neg aero fsoot= ',fsoot(i,j,1)
00585 trcrn(i,j,nt_aero)=max(trcrn(i,j,nt_aero),c0)
00586 trcrn(i,j,nt_aero+1)=max(trcrn(i,j,nt_aero+1),c0)
00587 trcrn(i,j,nt_aero+2)=max(trcrn(i,j,nt_aero+2),c0)
00588 trcrn(i,j,nt_aero+3)=max(trcrn(i,j,nt_aero+3),c0)
00589 endif
00590 enddo
00591
00592 end subroutine update_aerosol
00593
00594
00595
00596
00597
00598
00599
00600
00601
00602
00603
00604
00605
00606 subroutine write_restart_aero(filename_spec)
00607
00608
00609
00610
00611
00612
00613
00614
00615
00616
00617
00618
00619
00620 use ice_domain_size
00621 use ice_calendar, only: sec, month, mday, nyr, istep1, &
00622 time, time_forc, idate, year_init
00623 use ice_state
00624 use ice_read_write
00625 use ice_restart, only: lenstr, restart_dir, restart_file, pointer_file
00626
00627
00628
00629 character(len=char_len_long), intent(in), optional :: filename_spec
00630
00631
00632
00633 integer (kind=int_kind) ::
00634 i, j, k, n, it, iblk,
00635 iyear, imonth, iday
00636
00637 character(len=char_len_long) :: filename
00638
00639 logical (kind=log_kind) :: diag
00640
00641
00642 if (present(filename_spec)) then
00643 filename = trim(filename_spec)
00644 else
00645 iyear = nyr + year_init - 1
00646 imonth = month
00647 iday = mday
00648
00649 write(filename,'(a,a,a,i4.4,a,i2.2,a,i2.2,a,i5.5)') &
00650 restart_dir(1:lenstr(restart_dir)), &
00651 restart_file(1:lenstr(restart_file)),'.aero.', &
00652 iyear,'-',month,'-',mday,'-',sec
00653 end if
00654
00655
00656 call ice_open(nu_dump_aero,filename,0)
00657
00658 if (my_task == master_task) then
00659 write(nu_dump_aero) istep1,time,time_forc
00660 write(nu_diag,*) 'Writing ',filename(1:lenstr(filename))
00661 endif
00662
00663 diag = .true.
00664
00665
00666
00667 do k = 1, n_aero
00668 do n = 1, ncat
00669 call ice_write(nu_dump_aero,0,trcrn(:,:,nt_aero +(k-1)*4,n,:),'ruf8',diag)
00670 call ice_write(nu_dump_aero,0,trcrn(:,:,nt_aero+1+(k-1)*4,n,:),'ruf8',diag)
00671 call ice_write(nu_dump_aero,0,trcrn(:,:,nt_aero+2+(k-1)*4,n,:),'ruf8',diag)
00672 call ice_write(nu_dump_aero,0,trcrn(:,:,nt_aero+3+(k-1)*4,n,:),'ruf8',diag)
00673 enddo
00674 enddo
00675
00676 if (my_task == master_task) close(nu_dump_aero)
00677
00678 end subroutine write_restart_aero
00679
00680
00681
00682
00683
00684
00685
00686
00687 subroutine read_restart_aero(filename_spec)
00688
00689
00690
00691
00692
00693
00694
00695
00696
00697
00698
00699
00700
00701 use ice_domain_size
00702 use ice_calendar, only: sec, month, mday, nyr, istep1, &
00703 time, time_forc, idate, year_init
00704 use ice_state
00705 use ice_read_write
00706 use ice_restart, only: lenstr, restart_dir, restart_file, pointer_file
00707
00708
00709
00710 character(len=char_len_long), intent(in), optional :: filename_spec
00711
00712
00713
00714 integer (kind=int_kind) ::
00715 i, j, k, n, it, iblk,
00716 iyear, imonth, iday
00717
00718 character(len=char_len_long) ::
00719 filename, filename0, string1, string2
00720
00721 logical (kind=log_kind) ::
00722 diag
00723
00724 if (my_task == master_task) then
00725
00726 if (present(filename_spec)) then
00727 filename = filename_spec
00728 else
00729 open(nu_rst_pointer,file=pointer_file)
00730 read(nu_rst_pointer,'(a)') filename0
00731 filename = trim(filename0)
00732 close(nu_rst_pointer)
00733
00734 n = index(filename0,trim(restart_file))
00735 if (n == 0) call abort_ice('soot restart: filename discrepancy')
00736 string1 = trim(filename0(1:n-1))
00737 string2 = trim(filename0(n+lenstr(restart_file):lenstr(filename0)))
00738 write(filename,'(a,a,a,a)') &
00739 string1(1:lenstr(string1)), &
00740 restart_file(1:lenstr(restart_file)),'.aero', &
00741 string2(1:lenstr(string2))
00742 endif
00743 endif
00744
00745 call ice_open(nu_restart_aero,filename,0)
00746
00747 if (my_task == master_task) then
00748 read(nu_restart_aero) istep1,time,time_forc
00749 write(nu_diag,*) 'Reading ',filename(1:lenstr(filename))
00750 endif
00751
00752 diag = .true.
00753
00754
00755
00756 do k = 1, n_aero
00757 do n = 1, ncat
00758 call ice_read(nu_restart_aero,0,trcrn(:,:,nt_aero +(k-1)*4,n,:),'ruf8',&
00759 diag,field_type=field_type_scalar,field_loc=field_loc_center)
00760 call ice_read(nu_restart_aero,0,trcrn(:,:,nt_aero+1+(k-1)*4,n,:),'ruf8',&
00761 diag,field_type=field_type_scalar,field_loc=field_loc_center)
00762 call ice_read(nu_restart_aero,0,trcrn(:,:,nt_aero+2+(k-1)*4,n,:),'ruf8',&
00763 diag,field_type=field_type_scalar,field_loc=field_loc_center)
00764 call ice_read(nu_restart_aero,0,trcrn(:,:,nt_aero+3+(k-1)*4,n,:),'ruf8',&
00765 diag,field_type=field_type_scalar,field_loc=field_loc_center)
00766 enddo
00767 enddo
00768
00769 if (my_task == master_task) close(nu_restart_aero)
00770
00771 end subroutine read_restart_aero
00772
00773
00774
00775 end module ice_aerosol
00776
00777