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 module ice_restart
00026
00027
00028
00029 use ice_kinds_mod
00030 use ice_communicate, only: my_task, master_task
00031 use ice_blocks
00032 use ice_boundary
00033 use ice_read_write
00034 use ice_fileunits
00035 use ice_timers
00036 use ice_exit, only: abort_ice
00037
00038
00039
00040 implicit none
00041 save
00042
00043 interface dumpfile
00044 module procedure dumpfile_bin
00045 module procedure dumpfile_pio
00046 end interface
00047
00048 interface restartfile
00049 module procedure restartfile_bin
00050 module procedure restartfile_pio
00051 end interface
00052
00053 character (len=char_len) ::
00054 resttype
00055
00056 character(len=char_len_long) ::
00057 ice_ic
00058
00059
00060
00061
00062 logical (kind=log_kind) ::
00063 restart
00064
00065
00066
00067 character (len=char_len) ::
00068 runtype
00069
00070
00071
00072
00073
00074 character (len=char_len_long) ::
00075 restart_file ,
00076 restart_dir ,
00077 runid
00078
00079 character (len=char_len) ::
00080 restart_format
00081
00082 character (len=char_len_long) ::
00083 pointer_file
00084
00085 real (kind=dbl_kind), private,
00086 dimension(nx_block,ny_block,max_blocks) ::
00087 work1
00088
00089 logical (kind=log_kind) :: lcdf64
00090
00091 integer (kind=int_kind) :: ncid
00092
00093 integer (kind=int_kind) ::
00094 status
00095
00096 character (len=1) :: nchar
00097
00098
00099
00100 contains
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114 subroutine dumpfile_bin()
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126 use ice_domain_size
00127 use ice_flux
00128 use ice_grid
00129 use ice_calendar, only: sec, month, mday, nyr, istep1, &
00130 time, time_forc, idate, year_init
00131 use ice_state
00132 use ice_dyn_evp
00133 use ice_blocks, only : block, get_block, nx_block, ny_block
00134
00135
00136
00137 integer (kind=int_kind) ::
00138 i, j, k, n, it, iblk,
00139 iyear, imonth, iday
00140
00141 logical (kind=log_kind) :: diag
00142
00143 character(len=char_len_long) :: filename
00144
00145 iyear = nyr + year_init - 1
00146 imonth = month
00147 iday = mday
00148
00149 write(filename,'(a,a,a,i4.4,a,i2.2,a,i2.2,a,i5.5)') &
00150 restart_dir(1:lenstr(restart_dir)), &
00151 restart_file(1:lenstr(restart_file)),'.', &
00152 iyear,'-',month,'-',mday,'-',sec
00153
00154 call ice_open(nu_dump,filename,0)
00155
00156 if (my_task == master_task) then
00157 write(nu_dump) istep1,time,time_forc
00158 write(nu_diag,*) 'Writing ',filename(1:lenstr(filename))
00159 write(nu_diag,*) 'Restart written ',istep1,time,time_forc
00160 endif
00161
00162 diag = .true.
00163
00164
00165
00166
00167
00168 do n=1,ncat
00169 call ice_write(nu_dump,0,aicen(:,:,n,:),'ruf8',diag)
00170 call ice_write(nu_dump,0,vicen(:,:,n,:),'ruf8',diag)
00171 call ice_write(nu_dump,0,vsnon(:,:,n,:),'ruf8',diag)
00172 call ice_write(nu_dump,0,trcrn(:,:,nt_Tsfc,n,:),'ruf8',diag)
00173 enddo
00174
00175 do k=1,ntilyr
00176 call ice_write(nu_dump,0,eicen(:,:,k,:),'ruf8',diag)
00177 enddo
00178
00179 do k=1,ntslyr
00180 call ice_write(nu_dump,0,esnon(:,:,k,:),'ruf8',diag)
00181 enddo
00182
00183
00184
00185
00186 call ice_write(nu_dump,0,uvel,'ruf8',diag)
00187 call ice_write(nu_dump,0,vvel,'ruf8',diag)
00188
00189
00190
00191
00192 call ice_write(nu_dump,0,coszen,'ruf8',diag)
00193 call ice_write(nu_dump,0,scale_factor,'ruf8',diag)
00194 call ice_write(nu_dump,0,swvdr,'ruf8',diag)
00195 call ice_write(nu_dump,0,swvdf,'ruf8',diag)
00196 call ice_write(nu_dump,0,swidr,'ruf8',diag)
00197 call ice_write(nu_dump,0,swidf,'ruf8',diag)
00198
00199
00200
00201
00202 call ice_write(nu_dump,0,strocnxT,'ruf8',diag)
00203 call ice_write(nu_dump,0,strocnyT,'ruf8',diag)
00204
00205
00206
00207
00208 call ice_write(nu_dump,0,stressp_1,'ruf8',diag)
00209 call ice_write(nu_dump,0,stressp_3,'ruf8',diag)
00210 call ice_write(nu_dump,0,stressp_2,'ruf8',diag)
00211 call ice_write(nu_dump,0,stressp_4,'ruf8',diag)
00212
00213 call ice_write(nu_dump,0,stressm_1,'ruf8',diag)
00214 call ice_write(nu_dump,0,stressm_3,'ruf8',diag)
00215 call ice_write(nu_dump,0,stressm_2,'ruf8',diag)
00216 call ice_write(nu_dump,0,stressm_4,'ruf8',diag)
00217
00218 call ice_write(nu_dump,0,stress12_1,'ruf8',diag)
00219 call ice_write(nu_dump,0,stress12_3,'ruf8',diag)
00220 call ice_write(nu_dump,0,stress12_2,'ruf8',diag)
00221 call ice_write(nu_dump,0,stress12_4,'ruf8',diag)
00222
00223
00224
00225
00226
00227
00228 do iblk = 1, nblocks
00229 do j = 1, ny_block
00230 do i = 1, nx_block
00231 work1(i,j,iblk) = c0
00232 if (iceumask(i,j,iblk)) work1(i,j,iblk) = c1
00233 enddo
00234 enddo
00235 enddo
00236
00237 call ice_write(nu_dump,0,work1,'ruf8',diag)
00238
00239 if (my_task == master_task) then
00240 write(nu_dump) filename_volpn
00241 write(nu_dump) filename_aero
00242 write(nu_dump) filename_iage
00243 write(nu_dump) filename_FY
00244
00245 close(nu_dump)
00246 endif
00247
00248 end subroutine dumpfile_bin
00249
00250
00251
00252
00253
00254
00255
00256
00257
00258 subroutine dumpfile_pio(filename_spec)
00259
00260
00261
00262
00263
00264
00265
00266
00267
00268
00269
00270 use ice_domain_size
00271 use ice_flux
00272 use ice_grid
00273 use ice_calendar, only: sec, month, mday, nyr, istep1, &
00274 time, time_forc, idate, year_init
00275 use ice_state
00276 use ice_dyn_evp
00277 use ice_blocks, only : block, get_block, nx_block, ny_block
00278 use ice_pio
00279 use pio
00280
00281
00282
00283 character(len=char_len_long), intent(in) :: filename_spec
00284
00285
00286
00287 integer (kind=int_kind) ::
00288 i, j, k, n, it, iblk,
00289 ilo, ihi, jlo, jhi,
00290 lon, lat,
00291 iyear, imonth, iday
00292
00293 character(len=char_len_long) :: filename
00294
00295 logical (kind=log_kind) :: diag
00296
00297 integer (kind=int_kind) :: dimid_ni, dimid_nj, dimid_ncat,
00298 dimid_ntilyr, dimid_ntslyr
00299
00300 integer (kind=int_kind), allocatable :: dims(:)
00301
00302 type(file_desc_t) :: File
00303 type(io_desc_t) :: iodesc2d
00304 type(io_desc_t) :: iodesc3d_ncat
00305 type(io_desc_t) :: iodesc3d_ntilyr
00306 type(io_desc_t) :: iodesc3d_ntslyr
00307 type(var_desc_t) :: varid
00308
00309 real (kind=dbl_kind) :: tmp1d(1)
00310
00311 type(block) :: this_block
00312
00313
00314 filename = trim(filename_spec) // ".nc"
00315
00316
00317 if (my_task == master_task) then
00318 open(nu_rst_pointer,file=pointer_file)
00319 write(nu_rst_pointer,'(a)') filename
00320 close(nu_rst_pointer)
00321 endif
00322
00323
00324
00325 File%fh=-1
00326 call ice_pio_init(mode='write',filename=trim(filename), File=File, &
00327 clobber=.true., cdf64=lcdf64 )
00328
00329 call ice_pio_initdecomp(iodesc=iodesc2d)
00330 call ice_pio_initdecomp(ndim3=ncat , iodesc=iodesc3d_ncat)
00331 call ice_pio_initdecomp(ndim3=ntilyr, iodesc=iodesc3d_ntilyr)
00332 call ice_pio_initdecomp(ndim3=ntslyr, iodesc=iodesc3d_ntslyr)
00333
00334 status = pio_put_att(File,pio_global,'istep1',istep1)
00335 status = pio_put_att(File,pio_global,'time',time)
00336 status = pio_put_att(File,pio_global,'time_forc',time_forc)
00337 status = pio_put_att(File,pio_global, 'missing_value', c0)
00338 status = pio_put_att(File,pio_global,'_FillValue',c0)
00339
00340 status = pio_def_dim(File,'ni',nx_global,dimid_ni)
00341 status = pio_def_dim(File,'nj',ny_global,dimid_nj)
00342 status = pio_def_dim(File,'ncat',ncat,dimid_ncat)
00343 status = pio_def_dim(File,'ntilyr',ntilyr,dimid_ntilyr)
00344 status = pio_def_dim(File,'ntslyr',ntslyr,dimid_ntslyr)
00345 write(nu_diag,*) 'Writing ',filename(1:lenstr(filename))
00346 diag = .true.
00347
00348 allocate(dims(3))
00349
00350 dims(1) = dimid_ni
00351 dims(2) = dimid_nj
00352 dims(3) = dimid_ncat
00353
00354 call define_rest_field(File,'aicen',dims)
00355 call define_rest_field(File,'vicen',dims)
00356 call define_rest_field(File,'vsnon',dims)
00357 call define_rest_field(File,'Tsfcn',dims)
00358
00359 if (tr_aero) then
00360 do k=1,n_aero
00361 write(nchar,'(i1.1)') k
00362 call define_rest_field(File,'aerosnossl'//nchar, dims)
00363 call define_rest_field(File,'aerosnoint'//nchar, dims)
00364 call define_rest_field(File,'aeroicessl'//nchar, dims)
00365 call define_rest_field(File,'aeroiceint'//nchar, dims)
00366 enddo
00367 endif
00368
00369 if (tr_iage) then
00370 call define_rest_field(File,'iage',dims)
00371 end if
00372
00373 if (tr_pond) then
00374 call define_rest_field(File,'volpn' ,dims)
00375 call define_rest_field(File,'apondn',dims)
00376 call define_rest_field(File,'hpondn',dims)
00377 end if
00378
00379 dims(3) = dimid_ntilyr
00380 call define_rest_field(File,'eicen',dims)
00381 dims(3) = dimid_ntslyr
00382 call define_rest_field(File,'esnon',dims)
00383
00384 deallocate(dims)
00385
00386 allocate(dims(2))
00387 dims(1) = dimid_ni
00388 dims(2) = dimid_nj
00389
00390 call define_rest_field(File,'uvel',dims)
00391 call define_rest_field(File,'vvel',dims)
00392
00393 call define_rest_field(File,'coszen',dims)
00394 call define_rest_field(File,'scale_factor',dims)
00395 call define_rest_field(File,'swvdr',dims)
00396 call define_rest_field(File,'swvdf',dims)
00397 call define_rest_field(File,'swidr',dims)
00398 call define_rest_field(File,'swidf',dims)
00399
00400 call define_rest_field(File,'strocnxT',dims)
00401 call define_rest_field(File,'strocnyT',dims)
00402
00403 call define_rest_field(File,'stressp_1',dims)
00404 call define_rest_field(File,'stressp_2',dims)
00405 call define_rest_field(File,'stressp_3',dims)
00406 call define_rest_field(File,'stressp_4',dims)
00407
00408 call define_rest_field(File,'stressm_1',dims)
00409 call define_rest_field(File,'stressm_2',dims)
00410 call define_rest_field(File,'stressm_3',dims)
00411 call define_rest_field(File,'stressm_4',dims)
00412
00413 call define_rest_field(File,'stress12_1',dims)
00414 call define_rest_field(File,'stress12_2',dims)
00415 call define_rest_field(File,'stress12_3',dims)
00416 call define_rest_field(File,'stress12_4',dims)
00417
00418 call define_rest_field(File,'iceumask',dims)
00419
00420 status = pio_enddef(File)
00421
00422 deallocate(dims)
00423
00424
00425
00426
00427
00428 status = pio_inq_varid(File,'aicen',varid)
00429 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(aicen(:,:,:,:),tmp1d), status)
00430
00431 status = pio_inq_varid(File,'vicen',varid)
00432 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(vicen,tmp1d), status)
00433
00434 status = pio_inq_varid(File,'vsnon',varid)
00435 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(vsnon,tmp1d), status)
00436
00437 status = pio_inq_varid(File,'Tsfcn',varid)
00438 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(trcrn(:,:,nt_Tsfc,:,:),tmp1d), status)
00439
00440 status = pio_inq_varid(File,'eicen',varid)
00441 call pio_write_darray(File, varid, iodesc3d_ntilyr, transfer(eicen,tmp1d), status)
00442
00443 status = pio_inq_varid(File,'esnon',varid)
00444 call pio_write_darray(File, varid, iodesc3d_ntslyr, transfer(esnon,tmp1d), status)
00445
00446
00447
00448
00449 status = pio_inq_varid(File,'uvel',varid)
00450 call pio_write_darray(File, varid, iodesc2d, transfer(uvel,tmp1d), status)
00451
00452 status = pio_inq_varid(File,'vvel',varid)
00453 call pio_write_darray(File, varid, iodesc2d, transfer(vvel,tmp1d), status)
00454
00455
00456
00457
00458 status = pio_inq_varid(File,'coszen',varid)
00459 call pio_write_darray(File, varid, iodesc2d, transfer(coszen,tmp1d), status)
00460
00461 status = pio_inq_varid(File,'scale_factor',varid)
00462 call pio_write_darray(File, varid, iodesc2d, transfer(scale_factor,tmp1d), status)
00463
00464 status = pio_inq_varid(File,'swvdr',varid)
00465 call pio_write_darray(File, varid, iodesc2d, transfer(swvdr,tmp1d), status)
00466
00467 status = pio_inq_varid(File,'swvdf',varid)
00468 call pio_write_darray(File, varid, iodesc2d, transfer(swvdf,tmp1d), status)
00469
00470 status = pio_inq_varid(File,'swidr',varid)
00471 call pio_write_darray(File, varid, iodesc2d, transfer(swidr,tmp1d), status)
00472
00473 status = pio_inq_varid(File,'swidf',varid)
00474 call pio_write_darray(File, varid, iodesc2d, transfer(swidf,tmp1d), status)
00475
00476
00477
00478
00479 status = pio_inq_varid(File,'strocnxT',varid)
00480 call pio_write_darray(File, varid, iodesc2d, transfer(strocnxT,tmp1d), status)
00481
00482 status = pio_inq_varid(File,'strocnyT',varid)
00483 call pio_write_darray(File, varid, iodesc2d, transfer(strocnyT,tmp1d), status)
00484
00485
00486
00487
00488
00489 status = pio_inq_varid(File,'stressp_1',varid)
00490 call pio_write_darray(File, varid, iodesc2d, transfer(stressp_1,tmp1d), status)
00491
00492 status = pio_inq_varid(File,'stressp_2',varid)
00493 call pio_write_darray(File, varid, iodesc2d, transfer(stressp_2,tmp1d), status)
00494
00495 status = pio_inq_varid(File,'stressp_3',varid)
00496 call pio_write_darray(File, varid, iodesc2d, transfer(stressp_3,tmp1d), status)
00497
00498 status = pio_inq_varid(File,'stressp_4',varid)
00499 call pio_write_darray(File, varid, iodesc2d, transfer(stressp_4,tmp1d), status)
00500
00501 status = pio_inq_varid(File,'stressm_1',varid)
00502 call pio_write_darray(File, varid, iodesc2d, transfer(stressm_1,tmp1d), status)
00503
00504 status = pio_inq_varid(File,'stressm_2',varid)
00505 call pio_write_darray(File, varid, iodesc2d, transfer(stressm_2,tmp1d), status)
00506
00507 status = pio_inq_varid(File,'stressm_3',varid)
00508 call pio_write_darray(File, varid, iodesc2d, transfer(stressm_3,tmp1d), status)
00509
00510 status = pio_inq_varid(File,'stressm_4',varid)
00511 call pio_write_darray(File, varid, iodesc2d, transfer(stressm_4,tmp1d), status)
00512
00513 status = pio_inq_varid(File,'stress12_1',varid)
00514 call pio_write_darray(File, varid, iodesc2d, transfer(stress12_1,tmp1d), status)
00515
00516 status = pio_inq_varid(File,'stress12_2',varid)
00517 call pio_write_darray(File, varid, iodesc2d, transfer(stress12_2,tmp1d), status)
00518
00519 status = pio_inq_varid(File,'stress12_3',varid)
00520 call pio_write_darray(File, varid, iodesc2d, transfer(stress12_3,tmp1d), status)
00521
00522 status = pio_inq_varid(File,'stress12_4',varid)
00523 call pio_write_darray(File, varid, iodesc2d, transfer(stress12_4,tmp1d), status)
00524
00525
00526
00527
00528
00529
00530 do iblk = 1, nblocks
00531 do j = 1, ny_block
00532 do i = 1, nx_block
00533 work1(i,j,iblk) = c0
00534 if (iceumask(i,j,iblk)) work1(i,j,iblk) = c1
00535 enddo
00536 enddo
00537 enddo
00538
00539 status = pio_inq_varid(File,'iceumask',varid)
00540 call pio_write_darray(File, varid, iodesc2d, transfer(work1,tmp1d), status)
00541
00542 if (tr_aero) then
00543 do k=1,n_aero
00544 write(nchar,'(i1.1)') k
00545
00546 status = pio_inq_varid(File,'aerosnossl'//nchar,varid)
00547 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(trcrn(:,:,nt_aero+ (k-1)*4,:,:),tmp1d), status)
00548
00549 status = pio_inq_varid(File,'aerosnoint'//nchar,varid)
00550 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(trcrn(:,:,nt_aero+1+(k-1)*4,:,:),tmp1d), status)
00551
00552 status = pio_inq_varid(File,'aeroicessl'//nchar,varid)
00553 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(trcrn(:,:,nt_aero+2+(k-1)*4,:,:),tmp1d), status)
00554
00555 status = pio_inq_varid(File,'aeroiceint'//nchar,varid)
00556 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(trcrn(:,:,nt_aero+3+(k-1)*4,:,:),tmp1d), status)
00557 enddo
00558 endif
00559
00560 if (tr_iage) then
00561 status = pio_inq_varid(File,'iage',varid)
00562 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(trcrn(:,:,nt_iage,:,:),tmp1d), status)
00563 endif
00564
00565 if (tr_pond) then
00566 status = pio_inq_varid(File,'volpn',varid)
00567 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(trcrn(:,:,nt_volpn,:,:),tmp1d), status)
00568
00569 status = pio_inq_varid(File,'apondn',varid)
00570 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(apondn,tmp1d), status)
00571
00572 status = pio_inq_varid(File,'hpondn',varid)
00573 call pio_write_darray(File, varid, iodesc3d_ncat, transfer(hpondn,tmp1d), status)
00574 endif
00575
00576 call pio_closefile(File)
00577
00578 call PIO_freeDecomp(File,iodesc2d)
00579 call PIO_freeDecomp(File,iodesc3d_ncat)
00580 call PIO_freeDecomp(File,iodesc3d_ntilyr)
00581 call PIO_freeDecomp(File,iodesc3d_ntslyr)
00582
00583 call ice_pio_finalize
00584
00585 if (my_task == master_task) then
00586 write(nu_diag,*) 'Restart written ',istep1,time,time_forc
00587 endif
00588
00589 end subroutine dumpfile_pio
00590
00591
00592
00593
00594
00595
00596
00597
00598 subroutine define_rest_field(File, vname, dims)
00599
00600
00601
00602
00603
00604
00605
00606
00607
00608
00609
00610 use pio
00611
00612
00613
00614 type(file_desc_t) , intent(in) :: File
00615 character (len=*) , intent(in) :: vname
00616 integer (kind=int_kind), intent(in) :: dims(:)
00617
00618
00619
00620 type(var_desc_t) :: varid
00621
00622 status = pio_def_var(File,trim(vname),pio_double,dims,varid)
00623
00624 end subroutine define_rest_field
00625
00626
00627
00628
00629
00630
00631
00632
00633 subroutine restartfile_bin(ice_ic)
00634
00635
00636
00637
00638
00639
00640
00641
00642
00643
00644
00645 use ice_broadcast
00646 use ice_boundary
00647 use ice_domain_size
00648 use ice_domain
00649 use ice_calendar, only: istep0, istep1, time, time_forc, calendar
00650 use ice_flux
00651 use ice_state
00652 use ice_grid, only: tmask, umask, grid_type
00653 use ice_itd
00654 use ice_work, only: work_g1, work_g2
00655 use ice_gather_scatter, only: scatter_global_stress, gather_global
00656
00657
00658
00659 character(len=*), optional :: ice_ic
00660
00661
00662 integer (kind=int_kind) ::
00663 i, j, k, n, it, iblk,
00664 ilo, ihi, jlo, jhi,
00665 lon, lat,
00666 iyear, imonth, iday
00667
00668 character(len=char_len_long) ::
00669 filename, filename0
00670
00671 logical (kind=log_kind) ::
00672 diag
00673
00674 if (present(ice_ic)) then
00675 filename = ice_ic
00676 else
00677 if (my_task == master_task) then
00678 open(nu_rst_pointer,file=pointer_file)
00679 read(nu_rst_pointer,'(a)') filename0
00680 filename = trim(filename0)
00681 close(nu_rst_pointer)
00682 write(nu_diag,*) 'Read ',pointer_file(1:lenstr(pointer_file))
00683 endif
00684 call broadcast_scalar(filename, master_task)
00685 endif
00686
00687
00688
00689
00690
00691 if (tr_iage) trcrn(:,:,nt_iage, :,:) = c0
00692 if (tr_FY) trcrn(:,:,nt_FY, :,:) = c0
00693 if (tr_aero) trcrn(:,:,nt_aero:nt_aero+n_aero*4-1,:,:) = c0
00694
00695
00696 trcrn(:,:,nt_volpn,:,:) = c0
00697 apondn(:,:,:,:) = c0
00698 hpondn(:,:,:,:) = c0
00699
00700
00701 resttype = restformat(nu_restart,filename)
00702
00703 call ice_open(nu_restart,filename,0)
00704
00705 if (my_task == master_task) then
00706 write(nu_diag,*) 'Using restart dump=', trim(filename)
00707 read (nu_restart) istep0,time,time_forc
00708 write(nu_diag,*) 'Restart read at istep=',istep0,time,time_forc
00709 endif
00710
00711 call calendar(time)
00712
00713 call broadcast_scalar(istep0,master_task)
00714
00715 istep1 = istep0
00716
00717 call broadcast_scalar(time,master_task)
00718 call broadcast_scalar(time_forc,master_task)
00719
00720 diag = .true.
00721
00722
00723
00724
00725 do n=1,ncat
00726 if (my_task == master_task) &
00727 write(nu_diag,*) 'cat ',n, &
00728 ' min/max area, vol ice, vol snow, Tsfc'
00729
00730 call ice_read(nu_restart,0,aicen(:,:,n,:),'ruf8',diag, &
00731 field_type=field_type_scalar,field_loc=field_loc_center)
00732 call ice_read(nu_restart,0,vicen(:,:,n,:),'ruf8',diag, &
00733 field_type=field_type_scalar,field_loc=field_loc_center)
00734 call ice_read(nu_restart,0,vsnon(:,:,n,:),'ruf8',diag, &
00735 field_type=field_type_scalar,field_loc=field_loc_center)
00736 call ice_read(nu_restart,0,trcrn(:,:,nt_Tsfc,n,:),'ruf8',diag, &
00737 field_type=field_type_scalar,field_loc=field_loc_center)
00738 enddo
00739
00740 if (my_task == master_task) &
00741 write(nu_diag,*) 'min/max eicen for each layer'
00742 do k=1,ntilyr
00743 call ice_read(nu_restart,0,eicen(:,:,k,:),'ruf8',diag, &
00744 field_type=field_type_scalar,field_loc=field_loc_center)
00745 enddo
00746
00747 if (my_task == master_task) &
00748 write(nu_diag,*) 'min/max esnon for each layer'
00749 do k=1,ntslyr
00750 call ice_read(nu_restart,0,esnon(:,:,k,:),'ruf8',diag, &
00751 field_type=field_type_scalar,field_loc=field_loc_center)
00752 enddo
00753
00754
00755
00756
00757 if (my_task == master_task) &
00758 write(nu_diag,*) 'min/max velocity components'
00759
00760 call ice_read(nu_restart,0,uvel,'ruf8',diag, &
00761 field_type=field_type_vector,field_loc=field_loc_NEcorner)
00762 call ice_read(nu_restart,0,vvel,'ruf8',diag, &
00763 field_type=field_type_vector,field_loc=field_loc_NEcorner)
00764
00765
00766
00767
00768 if (trim(resttype) == 'new') then
00769 if (my_task == master_task) &
00770 write(nu_diag,*) 'radiation fields'
00771
00772 call ice_read(nu_restart,0,coszen,'ruf8',diag, &
00773 field_loc_center, field_type_scalar)
00774 call ice_read(nu_restart,0,scale_factor,'ruf8',diag, &
00775 field_loc_center, field_type_scalar)
00776 call ice_read(nu_restart,0,swvdr,'ruf8',diag, &
00777 field_loc_center, field_type_scalar)
00778 call ice_read(nu_restart,0,swvdf,'ruf8',diag, &
00779 field_loc_center, field_type_scalar)
00780 call ice_read(nu_restart,0,swidr,'ruf8',diag, &
00781 field_loc_center, field_type_scalar)
00782 call ice_read(nu_restart,0,swidf,'ruf8',diag, &
00783 field_loc_center, field_type_scalar)
00784 end if
00785
00786 if (trim(resttype) == 'old') then
00787 if (my_task == master_task) &
00788 write(nu_diag,*) 'min/max fresh water and heat flux components'
00789
00790 call ice_read(nu_restart,0,fresh,'ruf8',diag)
00791 call ice_read(nu_restart,0,fsalt,'ruf8',diag)
00792 call ice_read(nu_restart,0,fhocn,'ruf8',diag)
00793 endif
00794
00795
00796
00797
00798 if (my_task == master_task) &
00799 write(nu_diag,*) 'min/max ocean stress components'
00800
00801 call ice_read(nu_restart,0,strocnxT,'ruf8',diag)
00802 call ice_read(nu_restart,0,strocnyT,'ruf8',diag)
00803
00804
00805
00806
00807
00808
00809 if (my_task == master_task) write(nu_diag,*) &
00810 'internal stress components'
00811
00812 if (my_task==master_task) then
00813 allocate(work_g1(nx_global,ny_global))
00814 allocate(work_g2(nx_global,ny_global))
00815 else
00816 allocate(work_g1(1,1))
00817 allocate(work_g2(1,1))
00818 endif
00819
00820 call ice_read_global(nu_restart,0,work_g1,'ruf8',diag)
00821 call ice_read_global(nu_restart,0,work_g2,'ruf8',diag)
00822 call scatter_global_stress(stressp_1, work_g1, work_g2, &
00823 master_task, distrb_info)
00824 call scatter_global_stress(stressp_3, work_g2, work_g1, &
00825 master_task, distrb_info)
00826
00827 call ice_read_global(nu_restart,0,work_g1,'ruf8',diag)
00828 call ice_read_global(nu_restart,0,work_g2,'ruf8',diag)
00829 call scatter_global_stress(stressp_2, work_g1, work_g2, &
00830 master_task, distrb_info)
00831 call scatter_global_stress(stressp_4, work_g2, work_g1, &
00832 master_task, distrb_info)
00833
00834 call ice_read_global(nu_restart,0,work_g1,'ruf8',diag)
00835 call ice_read_global(nu_restart,0,work_g2,'ruf8',diag)
00836 call scatter_global_stress(stressm_1, work_g1, work_g2, &
00837 master_task, distrb_info)
00838 call scatter_global_stress(stressm_3, work_g2, work_g1, &
00839 master_task, distrb_info)
00840
00841 call ice_read_global(nu_restart,0,work_g1,'ruf8',diag)
00842 call ice_read_global(nu_restart,0,work_g2,'ruf8',diag)
00843 call scatter_global_stress(stressm_2, work_g1, work_g2, &
00844 master_task, distrb_info)
00845 call scatter_global_stress(stressm_4, work_g2, work_g1, &
00846 master_task, distrb_info)
00847
00848 call ice_read_global(nu_restart,0,work_g1,'ruf8',diag)
00849 call ice_read_global(nu_restart,0,work_g2,'ruf8',diag)
00850 call scatter_global_stress(stress12_1, work_g1, work_g2, &
00851 master_task, distrb_info)
00852 call scatter_global_stress(stress12_3, work_g2, work_g1, &
00853 master_task, distrb_info)
00854
00855 call ice_read_global(nu_restart,0,work_g1,'ruf8',diag)
00856 call ice_read_global(nu_restart,0,work_g2,'ruf8',diag)
00857 call scatter_global_stress(stress12_2, work_g1, work_g2, &
00858 master_task, distrb_info)
00859 call scatter_global_stress(stress12_4, work_g2, work_g1, &
00860 master_task, distrb_info)
00861
00862 deallocate (work_g1, work_g2)
00863
00864
00865
00866
00867 if (my_task == master_task) &
00868 write(nu_diag,*) 'ice mask for dynamics'
00869
00870 call ice_read(nu_restart,0,work1,'ruf8',diag)
00871
00872 iceumask(:,:,:) = .false.
00873
00874 do iblk = 1, nblocks
00875 do j = 1, ny_block
00876 do i = 1, nx_block
00877 if (work1(i,j,iblk) > p5) iceumask(i,j,iblk) = .true.
00878 enddo
00879 enddo
00880 enddo
00881
00882
00883 if (my_task == master_task) then
00884
00885 read(nu_restart, end=99) filename_volpn
00886 read(nu_restart, end=99) filename_aero
00887 read(nu_restart, end=99) filename_iage
00888 read(nu_restart, end=99) filename_FY
00889
00890 99 continue
00891 endif
00892
00893 if (my_task == master_task) then
00894 write(nu_diag,'(a,a)') 'filename_volpn: ',filename_volpn
00895 write(nu_diag,'(a,a)') 'filename_aero : ',filename_aero
00896 write(nu_diag,'(a,a)') 'filename_iage : ',filename_iage
00897 write(nu_diag,'(a,a)') 'filename_FY : ',filename_FY
00898 endif
00899
00900 if (my_task == master_task) close(nu_restart)
00901
00902 call broadcast_scalar(filename_volpn, master_task)
00903 call broadcast_scalar(filename_aero, master_task)
00904 call broadcast_scalar(filename_iage, master_task)
00905 call broadcast_scalar(filename_FY, master_task)
00906
00907
00908
00909
00910
00911 do iblk = 1, nblocks
00912 do j = 1, nghost
00913 do i = 1, nx_block
00914 stressp_1 (i,j,iblk) = c0
00915 stressp_2 (i,j,iblk) = c0
00916 stressp_3 (i,j,iblk) = c0
00917 stressp_4 (i,j,iblk) = c0
00918 stressm_1 (i,j,iblk) = c0
00919 stressm_2 (i,j,iblk) = c0
00920 stressm_3 (i,j,iblk) = c0
00921 stressm_4 (i,j,iblk) = c0
00922 stress12_1(i,j,iblk) = c0
00923 stress12_2(i,j,iblk) = c0
00924 stress12_3(i,j,iblk) = c0
00925 stress12_4(i,j,iblk) = c0
00926 enddo
00927 enddo
00928 do j = 1, ny_block
00929 do i = 1, nghost
00930 stressp_1 (i,j,iblk) = c0
00931 stressp_2 (i,j,iblk) = c0
00932 stressp_3 (i,j,iblk) = c0
00933 stressp_4 (i,j,iblk) = c0
00934 stressm_1 (i,j,iblk) = c0
00935 stressm_2 (i,j,iblk) = c0
00936 stressm_3 (i,j,iblk) = c0
00937 stressm_4 (i,j,iblk) = c0
00938 stress12_1(i,j,iblk) = c0
00939 stress12_2(i,j,iblk) = c0
00940 stress12_3(i,j,iblk) = c0
00941 stress12_4(i,j,iblk) = c0
00942 enddo
00943 enddo
00944 enddo
00945
00946
00947
00948
00949 do iblk = 1, nblocks
00950 do j = 1, ny_block
00951 do i = 1, nx_block
00952 if (.not. tmask(i,j,iblk)) then
00953 aicen(i,j,:,iblk) = c0
00954 vicen(i,j,:,iblk) = c0
00955 vsnon(i,j,:,iblk) = c0
00956 trcrn(i,j,nt_Tsfc,:,iblk) = c0
00957 eicen(i,j,:,iblk) = c0
00958 esnon(i,j,:,iblk) = c0
00959 endif
00960 if (.not. umask(i,j,iblk)) then
00961 uvel(i,j,iblk) = c0
00962 vvel(i,j,iblk) = c0
00963 endif
00964 enddo
00965 enddo
00966 enddo
00967
00968
00969
00970
00971
00972
00973
00974
00975
00976
00977
00978
00979
00980
00981
00982
00983 do iblk = 1, nblocks
00984
00985 call aggregate (nx_block, ny_block, &
00986 aicen(:,:,:,iblk), &
00987 trcrn(:,:,:,:,iblk),&
00988 vicen(:,:,:,iblk), &
00989 vsnon(:,:,:,iblk), &
00990 eicen(:,:,:,iblk), &
00991 esnon(:,:,:,iblk), &
00992 aice (:,:, iblk), &
00993 trcr (:,:,:,iblk), &
00994 vice (:,:, iblk), &
00995 vsno (:,:, iblk), &
00996 eice (:,:, iblk), &
00997 esno (:,:, iblk), &
00998 aice0(:,:, iblk), &
00999 tmask(:,:, iblk), &
01000 trcr_depend)
01001
01002 aice_init(:,:,iblk) = aice(:,:,iblk)
01003
01004 enddo
01005
01006
01007 end subroutine restartfile_bin
01008
01009
01010
01011
01012
01013
01014
01015
01016 subroutine restartfile_pio(usepio, ice_ic)
01017
01018
01019
01020
01021
01022
01023
01024
01025
01026
01027
01028 use ice_broadcast
01029 use ice_boundary
01030 use ice_domain_size
01031 use ice_domain
01032 use ice_calendar, only: istep0, istep1, time, time_forc, calendar
01033 use ice_flux
01034 use ice_state
01035 use ice_grid, only: tmask, umask, grid_type
01036 use ice_itd
01037 use ice_pio
01038 use pio
01039
01040
01041
01042 logical, intent(in) :: usepio
01043 character(len=*), intent(in), optional :: ice_ic
01044
01045
01046 integer (kind=int_kind) ::
01047 i, j, k, n, it, iblk,
01048 ilo, ihi, jlo, jhi,
01049 lon, lat,
01050 iyear, imonth, iday
01051
01052 character(len=char_len_long) ::
01053 filename, filename0
01054
01055 logical (kind=log_kind) ::
01056 diag
01057
01058 type(block) :: this_block
01059
01060 type(iosystem_desc_t) :: pio_subsystem
01061 type(file_desc_t) :: File
01062 type(io_desc_t) :: iodesc2d
01063 type(io_desc_t) :: iodesc3d_ncat
01064 type(io_desc_t) :: iodesc3d_ntilyr
01065 type(io_desc_t) :: iodesc3d_ntslyr
01066 type(var_desc_t) :: varid
01067
01068 real(kind=dbl_kind), pointer ::
01069 tmpfield2d(:),
01070 tmpfield3d_ncat(:),
01071 tmpfield3d_ntilyr(:),
01072 tmpfield3d_ntslyr(:)
01073
01074 real(kind = dbl_kind) ::
01075 tmp1d(1)
01076 real(kind = dbl_kind), allocatable ::
01077 work1_trc(:,:,:,:,:)
01078
01079 resttype = 'new'
01080
01081
01082
01083 if (tr_iage) trcrn(:,:,nt_iage, :,:) = c0
01084 if (tr_FY) trcrn(:,:,nt_FY, :,:) = c0
01085 if (tr_aero) trcrn(:,:,nt_aero:nt_aero+n_aero*4-1,:,:) = c0
01086
01087
01088 trcrn(:,:,nt_volpn,:,:) = c0
01089 apondn(:,:,:,:) = c0
01090 hpondn(:,:,:,:) = c0
01091
01092 if (present(ice_ic)) then
01093 filename = ice_ic
01094 else
01095 if (my_task == master_task) then
01096 open(nu_rst_pointer,file=pointer_file)
01097 read(nu_rst_pointer,'(a)') filename0
01098 filename = trim(filename0)
01099 close(nu_rst_pointer)
01100 write(nu_diag,*) 'Read ',pointer_file(1:lenstr(pointer_file))
01101 endif
01102 call broadcast_scalar(filename, master_task)
01103 endif
01104
01105
01106
01107 File%fh=-1
01108 call ice_pio_init(mode='read', filename=trim(filename), File=File)
01109
01110 call ice_pio_initdecomp(iodesc=iodesc2d)
01111 call ice_pio_initdecomp(ndim3=ncat , iodesc=iodesc3d_ncat)
01112 call ice_pio_initdecomp(ndim3=ntilyr, iodesc=iodesc3d_ntilyr)
01113 call ice_pio_initdecomp(ndim3=ntslyr, iodesc=iodesc3d_ntslyr)
01114
01115 if (my_task == master_task) then
01116 write(nu_diag,*) 'Using restart dump=', trim(filename)
01117 end if
01118
01119 status = pio_get_att(File, pio_global, 'istep1', istep1)
01120 status = pio_get_att(File, pio_global, 'time', time)
01121 status = pio_get_att(File, pio_global, 'time_forc', time_forc)
01122
01123 if (my_task == master_task) then
01124 write(nu_diag,*) 'Restart read at istep=',istep0,time,time_forc
01125 endif
01126
01127 call calendar(time)
01128
01129 istep1 = istep0
01130
01131 diag = .true.
01132
01133 allocate(tmpfield2d (nx_block*ny_block *max_blocks))
01134 allocate(tmpfield3d_ncat (nx_block*ny_block*ncat *max_blocks))
01135 allocate(tmpfield3d_ntilyr(nx_block*ny_block*ntilyr*max_blocks))
01136 allocate(tmpfield3d_ntslyr(nx_block*ny_block*ntslyr*max_blocks))
01137
01138
01139
01140
01141 if (my_task == master_task) &
01142 write(nu_diag,*) ' min/max area, vol ice, vol snow, Tsfc'
01143
01144 status = pio_inq_varid(File,'aicen',varid)
01145 call pio_read_darray(File, varid, iodesc3d_ncat, tmpfield3d_ncat, status)
01146 aicen(:,:,:,:) = reshape(tmpfield3d_ncat,(/nx_block,ny_block,ncat,max_blocks/))
01147
01148 status = pio_inq_varid(File,'vicen',varid)
01149 call pio_read_darray(File, varid, iodesc3d_ncat, tmpfield3d_ncat, status)
01150 vicen(:,:,:,:) = reshape(tmpfield3d_ncat,(/nx_block,ny_block,ncat,max_blocks/))
01151
01152 status = pio_inq_varid(File,'vsnon',varid)
01153 call pio_read_darray(File, varid, iodesc3d_ncat, tmpfield3d_ncat, status)
01154 vsnon(:,:,:,:) = reshape(tmpfield3d_ncat,(/nx_block,ny_block,ncat,max_blocks/))
01155
01156 status = pio_inq_varid(File,'Tsfcn',varid)
01157 allocate(work1_trc(nx_block,ny_block,1,ncat,max_blocks))
01158 call pio_read_darray(File, varid, iodesc3d_ncat, tmpfield3d_ncat, status)
01159 work1_trc(:,:,:,:,:) = reshape(tmpfield3d_ncat,(/nx_block,ny_block,1,ncat,max_blocks/))
01160 trcrn(:,:,nt_Tsfc,:,:) = work1_trc(:,:,1,:,:)
01161 deallocate(work1_trc)
01162
01163 if (my_task == master_task) &
01164 write(nu_diag,*) 'min/max eicen for each layer'
01165
01166 status = pio_inq_varid(File,'eicen',varid)
01167 call pio_read_darray(File, varid, iodesc3d_ntilyr, tmpfield3d_ntilyr, status)
01168 eicen(:,:,:,:) = reshape(tmpfield3d_ntilyr,(/nx_block,ny_block,ntilyr,max_blocks/))
01169
01170 if (my_task == master_task) &
01171 write(nu_diag,*) 'min/max esnon for each layer'
01172
01173 status = pio_inq_varid(File,'esnon',varid)
01174 call pio_read_darray(File, varid, iodesc3d_ntslyr, tmpfield3d_ntslyr, status)
01175 esnon(:,:,:,:) = reshape(tmpfield3d_ntslyr,(/nx_block,ny_block,ntslyr,max_blocks/))
01176
01177
01178
01179
01180 if (my_task == master_task) &
01181 write(nu_diag,*) 'min/max velocity components'
01182
01183 status = pio_inq_varid(File,'uvel',varid)
01184 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01185 uvel(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01186 call ice_HaloUpdate (uvel, halo_info, &
01187 field_loc_NEcorner, field_type_vector)
01188
01189 status = pio_inq_varid(File,'vvel',varid)
01190 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01191 vvel(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01192 call ice_HaloUpdate (vvel, halo_info, &
01193 field_loc_NEcorner, field_type_vector)
01194
01195
01196
01197
01198
01199 if (my_task == master_task) &
01200 write(nu_diag,*) 'radiation fields'
01201
01202 status = pio_inq_varid(File,'coszen',varid)
01203 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01204 coszen(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01205 call ice_HaloUpdate (coszen, halo_info, &
01206 field_loc_center, field_type_scalar)
01207
01208 status = pio_inq_varid(File,'scale_factor',varid)
01209 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01210 scale_factor(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01211 call ice_HaloUpdate (scale_factor, halo_info, &
01212 field_loc_center, field_type_scalar)
01213
01214 status = pio_inq_varid(File,'swvdr',varid)
01215 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01216 swvdr(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01217 call ice_HaloUpdate (swvdr, halo_info, &
01218 field_loc_center, field_type_scalar)
01219
01220 status = pio_inq_varid(File,'swvdf',varid)
01221 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01222 swvdf(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01223 call ice_HaloUpdate (swvdf, halo_info, &
01224 field_loc_center, field_type_scalar)
01225
01226 status = pio_inq_varid(File,'swidr',varid)
01227 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01228 swidr(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01229 call ice_HaloUpdate (swidr, halo_info, &
01230 field_loc_center, field_type_scalar)
01231
01232 status = pio_inq_varid(File,'swidf',varid)
01233 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01234 swidf(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01235 call ice_HaloUpdate (swidf, halo_info, &
01236 field_loc_center, field_type_scalar)
01237
01238
01239
01240
01241 if (my_task == master_task) &
01242 write(nu_diag,*) 'min/max ocean stress components'
01243
01244 status = pio_inq_varid(File,'strocnxT',varid)
01245 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01246 strocnxT(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01247
01248 status = pio_inq_varid(File,'strocnyT',varid)
01249 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01250 strocnyT(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01251
01252
01253
01254
01255
01256
01257 if (my_task == master_task) write(nu_diag,*) &
01258 'internal stress components'
01259
01260 status = pio_inq_varid(File,'stressp_1',varid)
01261 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01262 stressp_1(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01263
01264 status = pio_inq_varid(File,'stressp_3',varid)
01265 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01266 stressp_3(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01267
01268 status = pio_inq_varid(File,'stressp_2',varid)
01269 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01270 stressp_2(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01271
01272 status = pio_inq_varid(File,'stressp_4',varid)
01273 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01274 stressp_4(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01275
01276 status = pio_inq_varid(File,'stressm_1',varid)
01277 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01278 stressm_1(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01279
01280 status = pio_inq_varid(File,'stressm_3',varid)
01281 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01282 stressm_3(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01283
01284 status = pio_inq_varid(File,'stressm_2',varid)
01285 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01286 stressm_2(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01287
01288 status = pio_inq_varid(File,'stressm_4',varid)
01289 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01290 stressm_4(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01291
01292 status = pio_inq_varid(File,'stress12_1',varid)
01293 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01294 stress12_1(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01295
01296 status = pio_inq_varid(File,'stress12_3',varid)
01297 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01298 stress12_3(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01299
01300 status = pio_inq_varid(File,'stress12_2',varid)
01301 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01302 stress12_2(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01303
01304 status = pio_inq_varid(File,'stress12_4',varid)
01305 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01306 stress12_4(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01307
01308 call ice_HaloUpdate(stressp_1, halo_info, &
01309 field_loc_center, field_type_scalar)
01310 call ice_HaloUpdate(stressp_3, halo_info, &
01311 field_loc_center, field_type_scalar)
01312 call ice_HaloUpdate(stressp_2, halo_info, &
01313 field_loc_center, field_type_scalar)
01314 call ice_HaloUpdate(stressp_4, halo_info, &
01315 field_loc_center, field_type_scalar)
01316
01317 call ice_HaloUpdate(stressm_1, halo_info, &
01318 field_loc_center, field_type_scalar)
01319 call ice_HaloUpdate(stressm_3, halo_info, &
01320 field_loc_center, field_type_scalar)
01321 call ice_HaloUpdate(stressm_2, halo_info, &
01322 field_loc_center, field_type_scalar)
01323 call ice_HaloUpdate(stressm_4, halo_info, &
01324 field_loc_center, field_type_scalar)
01325
01326 call ice_HaloUpdate(stress12_1, halo_info, &
01327 field_loc_center, field_type_scalar)
01328 call ice_HaloUpdate(stress12_3, halo_info, &
01329 field_loc_center, field_type_scalar)
01330 call ice_HaloUpdate(stress12_2, halo_info, &
01331 field_loc_center, field_type_scalar)
01332 call ice_HaloUpdate(stress12_4, halo_info, &
01333 field_loc_center, field_type_scalar)
01334
01335
01336
01337 if (trim(grid_type) == 'tripole') then
01338
01339 call ice_HaloUpdate_stress(stressp_1, stressp_3, halo_info, &
01340 field_loc_center, field_type_scalar)
01341 call ice_HaloUpdate_stress(stressp_3, stressp_1, halo_info, &
01342 field_loc_center, field_type_scalar)
01343 call ice_HaloUpdate_stress(stressp_2, stressp_4, halo_info, &
01344 field_loc_center, field_type_scalar)
01345 call ice_HaloUpdate_stress(stressp_4, stressp_2, halo_info, &
01346 field_loc_center, field_type_scalar)
01347
01348 call ice_HaloUpdate_stress(stressm_1, stressm_3, halo_info, &
01349 field_loc_center, field_type_scalar)
01350 call ice_HaloUpdate_stress(stressm_3, stressm_1, halo_info, &
01351 field_loc_center, field_type_scalar)
01352 call ice_HaloUpdate_stress(stressm_2, stressm_4, halo_info, &
01353 field_loc_center, field_type_scalar)
01354 call ice_HaloUpdate_stress(stressm_4, stressm_2, halo_info, &
01355 field_loc_center, field_type_scalar)
01356
01357 call ice_HaloUpdate_stress(stress12_1, stress12_3, halo_info, &
01358 field_loc_center, field_type_scalar)
01359 call ice_HaloUpdate_stress(stress12_3, stress12_1, halo_info, &
01360 field_loc_center, field_type_scalar)
01361 call ice_HaloUpdate_stress(stress12_2, stress12_4, halo_info, &
01362 field_loc_center, field_type_scalar)
01363 call ice_HaloUpdate_stress(stress12_4, stress12_2, halo_info, &
01364 field_loc_center, field_type_scalar)
01365
01366 endif
01367
01368
01369
01370
01371 if (my_task == master_task) &
01372 write(nu_diag,*) 'ice mask for dynamics'
01373
01374 status = pio_inq_varid(File,'iceumask',varid)
01375 call pio_read_darray(File, varid, iodesc2d, tmpfield2d, status)
01376 work1(:,:,:) = reshape(tmpfield2d,(/nx_block,ny_block,max_blocks/))
01377 call ice_HaloUpdate (work1, halo_info, &
01378 field_loc_center, field_type_scalar)
01379
01380 iceumask(:,:,:) = .false.
01381 do iblk = 1, nblocks
01382 do j = 1, ny_block
01383 do i = 1, nx_block
01384 if (work1(i,j,iblk) > p5) iceumask(i,j,iblk) = .true.
01385 enddo
01386 enddo
01387 enddo
01388
01389 if (tr_aero .or. tr_iage .or. tr_pond .or. tr_FY) then
01390 allocate(work1_trc(nx_block,ny_block,1,ncat,max_blocks))
01391 end if
01392
01393 if (tr_aero) then
01394 status = pio_inq_varid(File,'aerosnossl1',varid)
01395
01396 if (status == PIO_noerr) then
01397
01398 do k=1,n_aero
01399 write(nchar,'(i1.1)') k
01400
01401 status = pio_inq_varid(File,'aerosnossl'//nchar,varid)
01402 call pio_read_darray(File, varid, iodesc3d_ncat, &
01403 tmpfield3d_ncat, status)
01404 work1_trc(:,:,:,:,:) = reshape(tmpfield3d_ncat, &
01405 (/nx_block,ny_block,1,ncat,max_blocks/))
01406 trcrn(:,:,nt_aero+(k-1)*4,:,:) = work1_trc(:,:,1,:,:)
01407
01408 status = pio_inq_varid(File,'aerosnoint'//nchar,varid)
01409 call pio_read_darray(File, varid, iodesc3d_ncat, &
01410 tmpfield3d_ncat, status)
01411 work1_trc(:,:,:,:,:) = reshape(tmpfield3d_ncat, &
01412 (/nx_block,ny_block,1,ncat,max_blocks/))
01413 trcrn(:,:,nt_aero+1+(k-1)*4,:,:) = work1_trc(:,:,1,:,:)
01414
01415 status = pio_inq_varid(File,'aeroicessl'//nchar,varid)
01416 call pio_read_darray(File, varid, iodesc3d_ncat, &
01417 tmpfield3d_ncat, status)
01418 work1_trc(:,:,:,:,:) = reshape(tmpfield3d_ncat, &
01419 (/nx_block,ny_block,1,ncat,max_blocks/))
01420 trcrn(:,:,nt_aero+2+(k-1)*4,:,:) = work1_trc(:,:,1,:,:)
01421
01422 status = pio_inq_varid(File,'aeroiceint'//nchar,varid)
01423 call pio_read_darray(File, varid, iodesc3d_ncat, &
01424 tmpfield3d_ncat, status)
01425 work1_trc(:,:,:,:,:) = reshape(tmpfield3d_ncat, &
01426 (/nx_block,ny_block,1,ncat,max_blocks/))
01427 trcrn(:,:,nt_aero+3+(k-1)*4,:,:) = work1_trc(:,:,1,:,:)
01428 enddo
01429
01430 endif
01431 endif
01432
01433 if (tr_iage) then
01434 status = pio_inq_varid(File,'iage',varid)
01435 if (status == PIO_noerr) then
01436 call pio_read_darray(File, varid, iodesc3d_ncat, tmpfield3d_ncat, &
01437 status)
01438 work1_trc(:,:,:,:,:) = reshape(tmpfield3d_ncat, &
01439 (/nx_block,ny_block,1,ncat,max_blocks/))
01440 trcrn(:,:,nt_iage,:,:) = work1_trc(:,:,1,:,:)
01441 endif
01442 endif
01443
01444 if (tr_FY) then
01445 status = pio_inq_varid(File,'FY',varid)
01446 if (status == PIO_noerr) then
01447 call pio_read_darray(File, varid, iodesc3d_ncat, tmpfield3d_ncat, &
01448 status)
01449 work1_trc(:,:,:,:,:) = reshape(tmpfield3d_ncat, &
01450 (/nx_block,ny_block,1,ncat,max_blocks/))
01451 trcrn(:,:,nt_FY,:,:) = work1_trc(:,:,1,:,:)
01452 endif
01453 endif
01454
01455 if (tr_pond) then
01456 status = pio_inq_varid(File,'volpn',varid)
01457 if (status == PIO_noerr) then
01458 call pio_read_darray(File, varid, iodesc3d_ncat, &
01459 tmpfield3d_ncat, status)
01460 work1_trc(:,:,:,:,:) = reshape(tmpfield3d_ncat, &
01461 (/nx_block,ny_block,1,ncat,max_blocks/))
01462 trcrn(:,:,nt_volpn,:,:) = work1_trc(:,:,1,:,:)
01463
01464 status = pio_inq_varid(File,'apondn',varid)
01465 call pio_read_darray(File, varid, iodesc3d_ncat, &
01466 tmpfield3d_ncat, status)
01467 apondn(:,:,:,:) = reshape(tmpfield3d_ncat, &
01468 (/nx_block,ny_block,ncat,max_blocks/))
01469
01470 status = pio_inq_varid(File,'hpondn',varid)
01471 call pio_read_darray(File, varid, iodesc3d_ncat, &
01472 tmpfield3d_ncat, status)
01473 hpondn(:,:,:,:) = reshape(tmpfield3d_ncat, &
01474 (/nx_block,ny_block,ncat,max_blocks/))
01475 endif
01476 endif
01477
01478 if (tr_aero .or. tr_iage .or. tr_pond .or. tr_FY) then
01479 deallocate(work1_trc)
01480 end if
01481
01482 deallocate(tmpfield2d )
01483 deallocate(tmpfield3d_ncat )
01484 deallocate(tmpfield3d_ntilyr)
01485 deallocate(tmpfield3d_ntslyr)
01486
01487 call pio_closefile(File)
01488
01489 call PIO_freeDecomp(File,iodesc2d)
01490 call PIO_freeDecomp(File,iodesc3d_ncat)
01491 call PIO_freeDecomp(File,iodesc3d_ntilyr)
01492 call PIO_freeDecomp(File,iodesc3d_ntslyr)
01493
01494 call ice_pio_finalize
01495
01496 call bound_state (aicen, trcrn, &
01497 vicen, vsnon, &
01498 eicen, esnon)
01499
01500
01501
01502
01503
01504 do iblk = 1, nblocks
01505 do j = 1, nghost
01506 do i = 1, nx_block
01507 stressp_1 (i,j,iblk) = c0
01508 stressp_2 (i,j,iblk) = c0
01509 stressp_3 (i,j,iblk) = c0
01510 stressp_4 (i,j,iblk) = c0
01511 stressm_1 (i,j,iblk) = c0
01512 stressm_2 (i,j,iblk) = c0
01513 stressm_3 (i,j,iblk) = c0
01514 stressm_4 (i,j,iblk) = c0
01515 stress12_1(i,j,iblk) = c0
01516 stress12_2(i,j,iblk) = c0
01517 stress12_3(i,j,iblk) = c0
01518 stress12_4(i,j,iblk) = c0
01519 enddo
01520 enddo
01521 do j = 1, ny_block
01522 do i = 1, nghost
01523 stressp_1 (i,j,iblk) = c0
01524 stressp_2 (i,j,iblk) = c0
01525 stressp_3 (i,j,iblk) = c0
01526 stressp_4 (i,j,iblk) = c0
01527 stressm_1 (i,j,iblk) = c0
01528 stressm_2 (i,j,iblk) = c0
01529 stressm_3 (i,j,iblk) = c0
01530 stressm_4 (i,j,iblk) = c0
01531 stress12_1(i,j,iblk) = c0
01532 stress12_2(i,j,iblk) = c0
01533 stress12_3(i,j,iblk) = c0
01534 stress12_4(i,j,iblk) = c0
01535 enddo
01536 enddo
01537 enddo
01538
01539
01540
01541
01542 do iblk = 1, nblocks
01543 do j = 1, ny_block
01544 do i = 1, nx_block
01545 if (.not. tmask(i,j,iblk)) then
01546 aicen(i,j,:,iblk) = c0
01547 vicen(i,j,:,iblk) = c0
01548 vsnon(i,j,:,iblk) = c0
01549 trcrn(i,j,nt_Tsfc,:,iblk) = c0
01550 eicen(i,j,:,iblk) = c0
01551 esnon(i,j,:,iblk) = c0
01552 endif
01553 if (.not. umask(i,j,iblk)) then
01554 uvel(i,j,iblk) = c0
01555 vvel(i,j,iblk) = c0
01556 endif
01557 enddo
01558 enddo
01559 enddo
01560
01561
01562
01563
01564
01565
01566
01567
01568
01569
01570
01571
01572
01573
01574
01575
01576 do iblk = 1, nblocks
01577
01578 call aggregate (nx_block, ny_block, &
01579 aicen(:,:,:,iblk), &
01580 trcrn(:,:,:,:,iblk),&
01581 vicen(:,:,:,iblk), &
01582 vsnon(:,:,:,iblk), &
01583 eicen(:,:,:,iblk), &
01584 esnon(:,:,:,iblk), &
01585 aice (:,:, iblk), &
01586 trcr (:,:,:,iblk), &
01587 vice (:,:, iblk), &
01588 vsno (:,:, iblk), &
01589 eice (:,:, iblk), &
01590 esno (:,:, iblk), &
01591 aice0(:,:, iblk), &
01592 tmask(:,:, iblk), &
01593 trcr_depend)
01594
01595 aice_init(:,:,iblk) = aice(:,:,iblk)
01596
01597 enddo
01598
01599
01600 end subroutine restartfile_pio
01601
01602
01603
01604
01605
01606
01607
01608
01609 integer function lenstr(label)
01610
01611
01612
01613
01614
01615
01616
01617
01618
01619
01620
01621
01622 character*(*) label
01623
01624
01625
01626 integer (kind=int_kind) ::
01627 length,
01628 n
01629
01630 length = len(label)
01631 do n=length,1,-1
01632 if( label(n:n) /= ' ' ) exit
01633 enddo
01634 lenstr = n
01635
01636 end function lenstr
01637
01638
01639
01640
01641
01642
01643
01644
01645 character(len=char_len) function restformat(nu, filename)
01646
01647
01648
01649
01650
01651
01652
01653
01654
01655
01656
01657
01658
01659
01660 integer, intent(in) :: nu
01661 character(len=char_len), intent(in) :: filename
01662
01663
01664
01665 logical :: readdata
01666 integer :: nrec
01667 integer :: idummy, ios
01668
01669 integer, parameter :: nrecold = ncat*4+ntilyr+ntslyr+21
01670 integer, parameter :: nrecnew = ncat*4+ntilyr+ntslyr+27
01671
01672 if (my_task == master_task) then
01673 call ice_open(nu_restart,filename,0)
01674
01675 readdata = .true.
01676 nrec = 0
01677 do while (readdata)
01678 read(nu, iostat=ios) idummy
01679 if (ios < 0) then
01680 readdata = .false.
01681 else
01682 nrec = nrec + 1
01683 end if
01684 end do
01685
01686 if (nrec == nrecold) then
01687 restformat = 'old'
01688 else if (nrec >= nrecnew) then
01689 restformat = 'new'
01690 else
01691 call abort_ice ( &
01692 'restformat: number of records on restart file not supported')
01693 end if
01694 write(nu_diag,*)'Number of records restart file = ',nrec,' restart format is = ',trim(restformat)
01695 close(nu)
01696 end if
01697 call broadcast_scalar(restformat,master_task)
01698
01699 end function restformat
01700
01701
01702
01703 end module ice_restart
01704
01705