00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018 module ice_prescaero_mod
00019
00020
00021
00022 use shr_stream_mod
00023 use shr_map_mod
00024 use shr_ncread_mod
00025 use shr_string_mod
00026 use shr_sys_mod
00027
00028 use ice_broadcast
00029 use ice_communicate, only : my_task, master_task
00030 use ice_kinds_mod
00031 use ice_fileunits
00032 use ice_exit, only : abort_ice
00033 use ice_domain_size, only : nx_global, ny_global, ncat, nilyr, max_blocks, &
00034 n_aero
00035 use ice_constants
00036 use ice_blocks, only : nx_block, ny_block
00037 use ice_domain, only : nblocks, distrb_info
00038 use ice_grid, only : TLAT,TLON,hm,tmask,tarea
00039 use ice_calendar, only : idate, sec
00040 use ice_itd, only : ilyr1, hin_max
00041 use ice_work, only : work_g1, work_g2
00042
00043 implicit none
00044 save
00045
00046 private
00047
00048
00049
00050
00051
00052
00053 public :: ice_prescaero_init
00054 public :: ice_prescaero_run
00055 public :: ice_prescaero_readField
00056 public :: ice_prescaero_phys
00057
00058
00059
00060 logical(kind=log_kind) , public :: prescribed_aero
00061 integer(kind=int_kind) , public :: stream_year_first_aero
00062 integer(kind=int_kind) , public :: stream_year_last_aero
00063 integer(kind=int_kind) , public :: model_year_align_aero
00064
00065 character(len=char_len_long), public :: stream_fldVarName_aero
00066 character(len=char_len_long), public :: stream_fldFileName_aero
00067 character(len=char_len_long), public :: stream_domTvarName_aero
00068 character(len=char_len_long), public :: stream_domXvarName_aero
00069 character(len=char_len_long), public :: stream_domYvarName_aero
00070 character(len=char_len_long), public :: stream_domAreaName_aero
00071 character(len=char_len_long), public :: stream_domMaskName_aero
00072 character(len=char_len_long), public :: stream_domFileName_aero
00073 logical(kind=log_kind) , public :: prescribed_aero_fill
00074
00075
00076
00077 real(kind=dbl_kind), allocatable :: dataXCoord(:,:)
00078 real(kind=dbl_kind), allocatable :: dataYCoord(:,:)
00079 integer(kind=int_kind), allocatable :: dataMask(:,:)
00080 real(kind=dbl_kind), allocatable :: dataUB(:,:)
00081 real(kind=dbl_kind), allocatable :: dataLB(:,:)
00082
00083 type(shr_stream_streamType), save :: aero_stream
00084 type(shr_map_mapType) :: cice_map
00085 type(shr_map_mapType) :: cice_fill
00086
00087 real(kind=dbl_kind), allocatable :: dataInLB(:,:,:)
00088 real(kind=dbl_kind), allocatable :: dataInUB(:,:,:)
00089 real(kind=dbl_kind), allocatable :: dataSrcLB(:,:)
00090 real(kind=dbl_kind), allocatable :: dataSrcUB(:,:)
00091 real(kind=dbl_kind), allocatable :: dataDstLB(:,:)
00092 real(kind=dbl_kind), allocatable :: dataDstUB(:,:)
00093 real(kind=dbl_kind), allocatable :: dataOutLB(:,:,:)
00094 real(kind=dbl_kind), allocatable :: dataOutUB(:,:,:)
00095
00096 real(kind=dbl_kind), allocatable :: aero_lb(:,:,:,:)
00097 real(kind=dbl_kind), allocatable :: aero_ub(:,:,:,:)
00098 real(kind=dbl_kind), allocatable :: aero_loc(:,:,:,:)
00099 real(kind=dbl_kind) :: worka(nx_block,ny_block,max_blocks)
00100
00101 logical(kind=log_kind) :: regrid
00102
00103 integer(kind=int_kind) :: dateUB, dateLB, secUB, secLB
00104 integer(kind=int_kind) :: nlon
00105 integer(kind=int_kind) :: nlat
00106 integer(kind=int_kind) :: nflds
00107 integer(kind=int_kind) :: n
00108
00109 character(len=char_len_long) :: fldList
00110 character(len=char_len) :: fldName
00111
00112 logical :: read_data
00113
00114
00115 contains
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132 subroutine ice_prescaero_init(prescribed_aero_in)
00133
00134
00135
00136 use ice_gather_scatter, only : gather_global
00137
00138 implicit none
00139
00140
00141
00142
00143
00144 real(kind=dbl_kind), allocatable :: work4d(:,:,:,:)
00145 integer(kind=int_kind), allocatable :: cice_mask_g(:,:)
00146 real(kind=dbl_kind) :: aice_max
00147 character(len=char_len_long) :: first_data_file
00148 character(len=char_len_long) :: domain_info_fn
00149 character(len=char_len_long) :: data_path
00150 integer(kind=int_kind) :: ndims
00151 character(len=char_len) :: timeName
00152 character(len=char_len) :: latName
00153 character(len=char_len) :: lonName
00154 character(len=char_len) :: maskName
00155 character(len=char_len) :: areaName
00156 logical(kind=log_kind) :: check
00157
00158 integer(kind=int_kind) :: nml_error
00159
00160
00161 character(*),parameter :: subName = "('ice_prescaero_init')"
00162 character(*),parameter :: F00 = "('(ice_prescaero_init) ',4a)"
00163 character(*),parameter :: F01 = "('(ice_prescaero_init) ',a,2i6)"
00164 character(*),parameter :: F02 = "('(ice_prescaero_init) ',a,2g20.13)"
00165
00166 logical(kind=log_kind), optional, intent(in) :: prescribed_aero_in
00167
00168
00169 namelist /ice_prescaero_nml/ prescribed_aero, prescribed_aero_fill, &
00170 stream_year_first_aero ,stream_year_last_aero ,model_year_align_aero, &
00171 stream_fldVarName_aero ,stream_fldFileName_aero, &
00172 stream_domTvarName_aero,stream_domXvarName_aero,stream_domYvarName_aero, &
00173 stream_domAreaName_aero,stream_domMaskName_aero,stream_domFileName_aero
00174
00175
00176 prescribed_aero = .false.
00177 stream_year_first_aero = 1
00178 stream_year_last_aero = 1
00179 model_year_align_aero = 1
00180 stream_fldVarName_aero = ' '
00181 stream_fldFileName_aero = ' '
00182 stream_domTvarName_aero = 'time'
00183 stream_domXvarName_aero = 'xc'
00184 stream_domYvarName_aero = 'yc'
00185 stream_domFileName_aero = ' '
00186 prescribed_aero_fill = .false.
00187
00188
00189 call get_fileunit(nu_nml)
00190 if (my_task == master_task) then
00191 open (nu_nml, file=nml_filename, status='old',iostat=nml_error)
00192 if (nml_error /= 0) then
00193 nml_error = -1
00194 else
00195 nml_error = 1
00196 endif
00197 do while (nml_error > 0)
00198 read(nu_nml, nml=ice_prescaero_nml,iostat=nml_error)
00199 if (nml_error > 0) read(nu_nml,*)
00200 end do
00201 if (nml_error == 0) close(nu_nml)
00202 endif
00203 call release_fileunit(nu_nml)
00204
00205 call broadcast_scalar(nml_error,master_task)
00206
00207 if (nml_error /= 0) then
00208 call abort_ice ('ice: Namelist read error in ice_prescaero_mod')
00209 endif
00210
00211 call broadcast_scalar(prescribed_aero,master_task)
00212 call broadcast_scalar(stream_year_first_aero,master_task)
00213 call broadcast_scalar(stream_year_last_aero,master_task)
00214 call broadcast_scalar(model_year_align_aero,master_task)
00215 call broadcast_scalar(stream_fldVarName_aero,master_task)
00216 call broadcast_scalar(stream_fldFileName_aero,master_task)
00217 call broadcast_scalar(stream_domTvarName_aero,master_task)
00218 call broadcast_scalar(stream_domXvarName_aero,master_task)
00219 call broadcast_scalar(stream_domYvarName_aero,master_task)
00220 call broadcast_scalar(stream_domAreaName_aero,master_task)
00221 call broadcast_scalar(stream_domMaskName_aero,master_task)
00222 call broadcast_scalar(stream_domFileName_aero,master_task)
00223 call broadcast_scalar(prescribed_aero_fill,master_task)
00224
00225 if (present(prescribed_aero_in)) prescribed_aero = prescribed_aero_in
00226
00227 if (.not.prescribed_aero) return
00228
00229 if (my_task == master_task) then
00230 write(nu_diag,*) ' prescribed_aero = ', prescribed_aero
00231 write(nu_diag,*) ' stream_year_first_aero = ', stream_year_first_aero
00232 write(nu_diag,*) ' stream_year_last_aero = ', stream_year_last_aero
00233 write(nu_diag,*) ' model_year_align_aero = ', model_year_align_aero
00234 write(nu_diag,*) ' prescribed_aero_fill = ', prescribed_aero_fill
00235
00236
00237
00238
00239 write (nu_diag,*) 'This is the prescribed aerosol option.'
00240 endif
00241
00242
00243
00244
00245
00246
00247 if (my_task == master_task) then
00248
00249 allocate (work_g1(nx_global,ny_global),work_g2(nx_global,ny_global),&
00250 & cice_mask_g(nx_global,ny_global))
00251
00252 else
00253
00254 allocate (work_g1(1,1),work_g2(1,1),&
00255 & cice_mask_g(1,1))
00256
00257 end if
00258
00259 call gather_global(work_g1, hm, master_task, distrb_info)
00260
00261 if (my_task == master_task) then
00262
00263 cice_mask_g = work_g1
00264
00265 endif
00266
00267 call gather_global(work_g1, TLAT, master_task, distrb_info)
00268 call gather_global(work_g2, TLON, master_task, distrb_info)
00269
00270 if (my_task == master_task) then
00271
00272 work_g1(:,:) = work_g1(:,:)*rad_to_deg
00273 work_g2(:,:) = work_g2(:,:)*rad_to_deg
00274
00275
00276
00277
00278 aero_stream%nFiles = 0
00279 aero_stream%file(:)%name = 'not_set'
00280 aero_stream%file(:)%nt = 0
00281 aero_stream%file(:)%haveData = .false.
00282 aero_stream%yearFirst = stream_year_first_aero
00283 aero_stream%yearLast = stream_year_last_aero
00284 aero_stream%yearAlign = model_year_align_aero
00285 aero_stream%fldListFile = ' '
00286 aero_stream%fldListModel = ' '
00287 aero_stream%FilePath = ' '
00288 aero_stream%domFilePath = ' '
00289 aero_stream%k_lvd = -1
00290 aero_stream%n_lvd = -1
00291 aero_stream%found_lvd = .false.
00292 aero_stream%k_gvd = -1
00293 aero_stream%n_gvd = -1
00294 aero_stream%found_gvd = .false.
00295
00296 aero_stream%dataSource = 'cice ifrac/sst file'
00297 aero_stream%fldListFile = stream_fldVarName_aero
00298
00299
00300 aero_stream%File(1)%name = stream_fldFileName_aero
00301 aero_stream%nFiles = 1
00302
00303 aero_stream%domTvarName = stream_domTvarName_aero
00304 aero_stream%domXvarName = stream_domXvarName_aero
00305 aero_stream%domYvarName = stream_domYvarName_aero
00306 aero_stream%domAreaName = stream_domAreaName_aero
00307 aero_stream%domMaskName = stream_domMaskName_aero
00308 aero_stream%domFileName = stream_domFileName_aero
00309 aero_stream%init = .true.
00310
00311
00312
00313
00314
00315 call shr_stream_getFileFieldList(aero_stream,fldList)
00316 nflds = shr_string_listGetNum(fldList)
00317
00318 do n=1,nflds
00319
00320 call shr_string_listGetName(fldList,n,fldName)
00321
00322 call shr_stream_getFirstFileName(aero_stream,first_data_file,data_path)
00323 call shr_stream_getFile (data_path,first_data_file)
00324 check = shr_ncread_varExists(first_data_file,fldName)
00325
00326 if (.not.check) then
00327 write(nu_diag,F00) "ERROR: field does not exist"
00328 call abort_ice(subName)
00329 end if
00330
00331 enddo
00332
00333
00334
00335
00336 call shr_stream_getDomainInfo(aero_stream,data_path,domain_info_fn, &
00337 timeName, lonName, latName, maskName, areaName)
00338
00339 call shr_ncread_varDimSizes(domain_info_fn,lonName,nlon)
00340 call shr_ncread_varDimSizes(domain_info_fn,latName,nlat)
00341 call shr_stream_getFile(data_path,domain_info_fn)
00342 call shr_sys_flush(nu_diag)
00343
00344 allocate(dataXCoord(nlon,nlat))
00345 allocate(dataYCoord(nlon,nlat))
00346
00347
00348
00349
00350 call shr_ncread_domain(domain_info_fn, lonName, dataXCoord, latName, dataYCoord)
00351
00352 write (nu_diag,F02) 'min/max dataXCoord = ',minval(dataXCoord) ,maxval(dataXCoord)
00353 write (nu_diag,F02) 'min/max dataYCoord = ',minval(dataYCoord) ,maxval(dataYCoord)
00354 write (nu_diag,F02) 'min/max TLON = ',minval(work_g2) ,maxval(work_g2)
00355 write (nu_diag,F02) 'min/max TLAT = ',minval(work_g1) ,maxval(work_g1)
00356 write (nu_diag,F01) 'min/max cice_mask_g = ',minval(cice_mask_g),maxval(cice_mask_g)
00357
00358
00359
00360
00361 call ice_prescaero_checkDomain(work_g2, work_g1, dataXCoord, dataYCoord, regrid)
00362
00363
00364
00365
00366 if (regrid) then
00367 allocate(dataMask(nlon,nlat))
00368 dataMask(:,:) = 1
00369
00370 if (prescribed_aero_fill) then
00371 call shr_map_mapSet(cice_fill, work_g2, work_g1, cice_mask_g, &
00372 & work_g2, work_g1, cice_mask_g, &
00373 & name='cice_map',type='fill',algo='nnoni', &
00374 & mask='nomask')
00375 end if
00376
00377 write (nu_diag,F00) 'Computing mapping weights'
00378
00379 call shr_map_mapSet(cice_map, dataXCoord, dataYCoord, dataMask, &
00380 & work_g2, work_g1, cice_mask_g, &
00381 & name='cice_map',type='remap',algo='bilinear', &
00382 & mask='dstmask',vect='scalar')
00383 deallocate(dataMask)
00384 end if
00385
00386
00387 deallocate(dataXCoord,dataYCoord)
00388
00389 end if
00390
00391
00392
00393
00394
00395 call broadcast_scalar(nflds, master_task)
00396
00397 allocate(dataOutLB(nx_global,ny_global,nflds))
00398 allocate(dataOutUB(nx_global,ny_global,nflds))
00399 allocate(aero_lb(nx_block,ny_block,nflds,max_blocks))
00400 allocate(aero_ub(nx_block,ny_block,nflds,max_blocks))
00401 allocate(aero_loc(nx_block,ny_block,nflds,max_blocks))
00402
00403 deallocate(work_g1,work_g2)
00404 deallocate(cice_mask_g)
00405
00406
00407
00408
00409 if (ncat == 1) then
00410 hin_max(1) = 999._dbl_kind
00411 end if
00412
00413 end subroutine ice_prescaero_init
00414
00415
00416
00417
00418
00419
00420
00421
00422
00423
00424
00425
00426
00427
00428
00429
00430 subroutine ice_prescaero_run(mDateIn, secIn)
00431
00432
00433
00434 use shr_tInterp_mod
00435 use ice_constants, only: field_loc_center, field_type_scalar
00436 use ice_gather_scatter, only : scatter_global
00437
00438 implicit none
00439
00440
00441
00442 integer(kind=int_kind), intent(in) :: mDateIn
00443 integer(kind=int_kind), intent(in) :: secIn
00444
00445
00446
00447 real(kind=dbl_kind) :: fLB
00448 real(kind=dbl_kind) :: fUB
00449 integer(kind=int_kind) :: mDateLB
00450 integer(kind=int_kind) :: mDateUB
00451 integer(kind=int_kind) :: dDateLB
00452 integer(kind=int_kind) :: dDateUB
00453 integer(kind=int_kind) :: secUB
00454 integer(kind=int_kind) :: secLB
00455 integer(kind=int_kind) :: n_lb
00456 integer(kind=int_kind) :: n_ub
00457 character(len=char_len_long) :: fileLB
00458 character(len=char_len_long) :: fileUB
00459
00460 integer(kind=int_kind) :: mDateLB_old = -999
00461 integer(kind=int_kind) :: secLB_old = -999
00462 integer(kind=int_kind) :: i,j,n,icnt,iblk
00463
00464
00465
00466
00467
00468
00469 if (my_task == master_task) then
00470
00471
00472
00473
00474
00475 call shr_stream_findBounds(aero_stream, mDateIn,secIn, &
00476 & mDateLB, dDateLB, secLB, n_lb, fileLB, &
00477 & mDateUB, dDateUB, secUB, n_ub, fileUB)
00478
00479 read_data = .false.
00480 if (mDateLB_old /= mDateLB .or. secLB_old /= secLB) then
00481
00482 read_data = .true.
00483
00484 allocate(dataInLB(nlon,nlat,nflds))
00485 allocate(dataInUB(nlon,nlat,nflds))
00486
00487 do n=1,nflds
00488 call shr_string_listGetName(fldList,n,fldName)
00489
00490 call ice_prescaero_readField(fileLB, fldName, n_lb, dataInLB(:,:,n))
00491 call ice_prescaero_readField(fileUB, fldName, n_ub, dataInUB(:,:,n))
00492 enddo
00493
00494 if (regrid) then
00495
00496 allocate(dataSrcLB(nflds,nlon*nlat))
00497 allocate(dataSrcUB(nflds,nlon*nlat))
00498 allocate(dataDstLB(nflds,nx_global*ny_global))
00499 allocate(dataDstUB(nflds,nx_global*ny_global))
00500
00501
00502
00503
00504 do n=1,nflds
00505 icnt = 0
00506 do j=1,nlat
00507 do i=1,nlon
00508 icnt = icnt + 1
00509 dataSrcLB(n,icnt) = dataInLB(i,j,n)
00510 dataSrcUB(n,icnt) = dataInUB(i,j,n)
00511 enddo
00512 enddo
00513 enddo
00514
00515
00516
00517
00518 call shr_map_mapData(dataSrcLB, dataDstLB, cice_map)
00519 call shr_map_mapData(dataSrcUB, dataDstUB, cice_map)
00520
00521 if (prescribed_aero_fill) then
00522 call shr_map_mapData(dataDstLB, dataDstLB, cice_fill)
00523 call shr_map_mapData(dataDstUB, dataDstUB, cice_fill)
00524 end if
00525
00526
00527
00528 do n=1,nflds
00529 icnt = 0
00530 do j=1,ny_global
00531 do i=1,nx_global
00532 icnt = icnt + 1
00533 dataOutLB(i,j,n) = dataDstLB(n,icnt)
00534 dataOutUB(i,j,n) = dataDstUB(n,icnt)
00535 enddo
00536 enddo
00537 enddo
00538 deallocate(dataSrcLB)
00539 deallocate(dataSrcUB)
00540 deallocate(dataDstLB)
00541 deallocate(dataDstUB)
00542 deallocate(dataInLB)
00543 deallocate(dataInUB)
00544
00545 else
00546 dataOutLB = dataInLB
00547 dataOutUB = dataInUB
00548
00549 deallocate(dataInLB)
00550 deallocate(dataInUB)
00551 end if
00552
00553 mDateLB_old = mDateLB
00554 secLB_old = secLB
00555
00556 end if
00557
00558 call shr_tInterp_getFactors(mDateLB, secLB, mDateUB, secUB, &
00559 & mDateIn, secIN, fLB, fUB)
00560
00561 end if
00562
00563 call broadcast_scalar(read_data, master_task)
00564 call broadcast_scalar(fLB, master_task)
00565 call broadcast_scalar(fUB, master_task)
00566
00567 if (read_data) then
00568
00569
00570
00571
00572 if (my_task == master_task) then
00573 allocate(work_g1(nx_global,ny_global))
00574 else
00575 allocate(work_g1(1,1))
00576 endif
00577
00578 do n=1,nflds
00579
00580 if (my_task == master_task) work_g1(:,:) = dataOutLB(:,:,n)
00581 call scatter_global(worka, work_g1, &
00582 master_task, distrb_info, &
00583 field_loc_center, field_type_scalar)
00584
00585 do iblk = 1,nblocks
00586 do j = 1,ny_block
00587 do i = 1,nx_block
00588 aero_lb(i,j,n,iblk) = worka(i,j,iblk)
00589 enddo
00590 enddo
00591 enddo
00592
00593 if (my_task == master_task) work_g1(:,:) = dataOutUB(:,:,n)
00594 call scatter_global(worka, work_g1, &
00595 master_task, distrb_info, &
00596 field_loc_center, field_type_scalar)
00597
00598 do iblk = 1,nblocks
00599 do j = 1,ny_block
00600 do i = 1,nx_block
00601 aero_ub(i,j,n,iblk) = worka(i,j,iblk)
00602 enddo
00603 enddo
00604 enddo
00605
00606 enddo
00607
00608 deallocate(work_g1)
00609
00610 endif
00611
00612
00613
00614
00615
00616 do iblk = 1,nblocks
00617 do n = 1,nflds
00618 do j = 1,ny_block
00619 do i = 1,nx_block
00620 aero_loc(i,j,n,iblk) = fLB*aero_lb(i,j,n,iblk) + fUB*aero_ub(i,j,n,iblk)
00621 enddo
00622 enddo
00623 enddo
00624 enddo
00625
00626
00627
00628
00629
00630
00631 call ice_prescaero_phys
00632
00633 end subroutine ice_prescaero_run
00634
00635
00636
00637
00638
00639
00640
00641
00642
00643
00644
00645
00646
00647
00648 subroutine ice_prescaero_readField(fileName, fldName, nTime, fldRead)
00649
00650
00651
00652 implicit none
00653
00654
00655
00656 character(*), intent(in) :: fileName
00657 character(*), intent(in) :: fldName
00658 integer(kind=int_kind), intent(in) :: nTime
00659 real (kind=dbl_kind), intent(out) :: fldRead(:,:)
00660
00661
00662
00663 real(kind=dbl_kind), allocatable :: A4d(:,:,:,:)
00664 integer(kind=int_kind) :: ni, nj
00665
00666 ni = size(fldRead,1)
00667 nj = size(fldRead,2)
00668
00669 allocate(A4d(ni,nj,1,1))
00670
00671 call shr_ncread_field4dG(fileName, fldName, rfld=A4d, dim3i=nTime)
00672 fldRead(:,:) = A4d(:,:,1,1)
00673
00674 deallocate(A4d)
00675
00676 end subroutine ice_prescaero_readField
00677
00678
00679
00680
00681
00682
00683
00684
00685
00686
00687
00688
00689
00690
00691 subroutine ice_prescaero_checkDomain(aeroXCoord, aeroYCoord, dataXCoord, dataYCoord, regrid)
00692
00693
00694
00695 implicit none
00696
00697
00698
00699 real(kind=dbl_kind), intent(in) :: aeroXCoord(:,:)
00700 real(kind=dbl_kind), intent(in) :: aeroYCoord(:,:)
00701 real(kind=dbl_kind), intent(in) :: dataXCoord(:,:)
00702 real(kind=dbl_kind), intent(in) :: dataYCoord(:,:)
00703 logical(kind=log_kind), intent(out):: regrid
00704
00705
00706
00707 integer(kind=int_kind) :: ni_aero, nj_aero
00708 integer(kind=int_kind) :: ni_data, nj_data
00709 real(kind=dbl_kind) :: max_XDiff, max_YDiff
00710 real(kind=dbl_kind), allocatable :: aeroXTemp(:,:)
00711
00712
00713 character(*),parameter :: F00 = "('(ice_prescaero_checkDomain) ',a)"
00714 character(*),parameter :: F01 = "('(ice_prescaero_checkDomain) ',a,2i6)"
00715 character(*),parameter :: F02 = "('(ice_prescaero_checkDomain) ',a,g20.13)"
00716
00717 regrid = .false.
00718
00719
00720
00721 ni_aero = size(aeroXCoord,1)
00722 nj_aero = size(aeroXCoord,2)
00723 ni_data = size(dataXCoord,1)
00724 nj_data = size(dataXCoord,2)
00725
00726
00727
00728
00729 allocate(aeroXTemp(ni_aero,nj_aero))
00730 aeroXTemp = aeroXCoord
00731 where (aeroXTemp >= c360) aeroXTemp = aeroXTemp - c360
00732 where (aeroXTemp < c0 ) aeroXTemp = aeroXTemp + c360
00733
00734 if (ni_data == ni_aero .and. nj_data == nj_aero) then
00735 max_XDiff = maxval(abs(aeroXTemp - dataXCoord))
00736 max_YDiff = maxval(abs(aeroYCoord - dataYCoord))
00737 if (max_XDiff > eps04) regrid = .true.
00738 if (max_YDiff > eps04) regrid = .true.
00739
00740 write(nu_diag,F02) 'Maximum X difference = ', max_XDiff
00741 write(nu_diag,F02) 'Maximum Y difference = ', max_YDiff
00742 write(nu_diag,F00) 'CICE and input data domain sizes match'
00743 write(nu_diag,F01) 'CICE and data grids are', ni_aero,nj_aero
00744 else
00745 regrid = .true.
00746 end if
00747 deallocate(aeroXTemp)
00748
00749 if (regrid) then
00750 write(nu_diag,F00) 'Aerosols will be interpolated to CICE grid'
00751 else
00752 write(nu_diag,F00) 'Aerosol grid and CICE grid the same'
00753 end if
00754
00755 end subroutine ice_prescaero_checkDomain
00756
00757
00758
00759
00760
00761
00762
00763
00764
00765
00766
00767
00768
00769
00770
00771 subroutine ice_prescaero_phys
00772
00773
00774
00775 use ice_flux
00776
00777 implicit none
00778
00779
00780
00781
00782
00783
00784 integer(kind=int_kind) :: i,j,n
00785 integer(kind=int_kind) :: iblk
00786
00787 do iblk = 1,nblocks
00788 do j = 1,ny_block
00789 do i = 1,nx_block
00790 faero(i,j,1,iblk) = aero_loc(i,j, 1,iblk)
00791 faero(i,j,2,iblk) = aero_loc(i,j, 2,iblk)+aero_loc(i,j, 3,iblk)
00792
00793 faero(i,j,3,iblk) = aero_loc(i,j, 4,iblk)+aero_loc(i,j, 5,iblk) &
00794 + aero_loc(i,j, 6,iblk)+aero_loc(i,j, 7,iblk) &
00795 + aero_loc(i,j, 8,iblk)+aero_loc(i,j, 9,iblk) &
00796 + aero_loc(i,j,10,iblk)+aero_loc(i,j,11,iblk)
00797
00798
00799
00800 enddo
00801 enddo
00802 enddo
00803
00804 end subroutine ice_prescaero_phys
00805
00806
00807
00808 end module ice_prescaero_mod
00809
00810