00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023 module ice_prescribed_mod
00024
00025
00026
00027 use shr_stream_mod
00028 use shr_map_mod
00029 use shr_ncread_mod
00030 use shr_sys_mod
00031
00032 use ice_broadcast
00033 use ice_communicate, only : my_task, master_task
00034 use ice_kinds_mod
00035 use ice_fileunits
00036 use ice_exit, only : abort_ice
00037 use ice_domain_size, only : nx_global, ny_global, ncat, nilyr, nslyr, &
00038 max_blocks
00039 use ice_constants
00040 use ice_blocks, only : nx_block, ny_block
00041 use ice_domain, only : nblocks, distrb_info
00042 use ice_grid, only : TLAT,TLON,hm,tmask
00043 use ice_calendar, only : idate, sec
00044 use ice_itd, only : ilyr1, slyr1, hin_max
00045 use ice_work, only : work_g1, work_g2
00046 use shr_scam_mod, only : shr_scam_getCloseLatLon
00047 use ice_scam, only : scmlat, scmlon, single_column
00048 use ice_read_write
00049
00050 implicit none
00051 save
00052
00053 private
00054
00055
00056
00057
00058
00059
00060 public :: ice_prescribed_init
00061 public :: ice_prescribed_run
00062 public :: ice_prescribed_readField
00063 public :: ice_prescribed_phys
00064
00065
00066
00067 logical(kind=log_kind) , public :: prescribed_ice
00068 integer(kind=int_kind) , public :: stream_year_first
00069 integer(kind=int_kind) , public :: stream_year_last
00070 integer(kind=int_kind) , public :: model_year_align
00071
00072 character(len=char_len_long), public :: stream_fldVarName
00073 character(len=char_len_long), public :: stream_fldFileName
00074 character(len=char_len_long), public :: stream_domTvarName
00075 character(len=char_len_long), public :: stream_domXvarName
00076 character(len=char_len_long), public :: stream_domYvarName
00077 character(len=char_len_long), public :: stream_domAreaName
00078 character(len=char_len_long), public :: stream_domMaskName
00079 character(len=char_len_long), public :: stream_domFileName
00080 logical(kind=log_kind) , public :: prescribed_ice_fill
00081 real(kind=dbl_kind) :: closelat,closelon
00082 integer(kind=int_kind) :: latidx,lonidx
00083 integer(kind=int_kind) :: ncid
00084
00085
00086
00087 real(kind=dbl_kind), allocatable :: dataXCoord(:,:)
00088 real(kind=dbl_kind), allocatable :: dataYCoord(:,:)
00089 integer(kind=int_kind), allocatable :: dataMask(:,:)
00090 real(kind=dbl_kind), allocatable :: dataUB(:,:)
00091 real(kind=dbl_kind), allocatable :: dataLB(:,:)
00092
00093 type(shr_stream_streamType), save :: csim_stream
00094 type(shr_map_mapType) :: csim_map
00095 type(shr_map_mapType) :: csim_fill
00096
00097 real(kind=dbl_kind), allocatable :: dataInLB(:,:,:)
00098 real(kind=dbl_kind), allocatable :: dataInUB(:,:,:)
00099 real(kind=dbl_kind), allocatable :: dataSrcLB(:,:)
00100 real(kind=dbl_kind), allocatable :: dataSrcUB(:,:)
00101 real(kind=dbl_kind), allocatable :: dataDstLB(:,:)
00102 real(kind=dbl_kind), allocatable :: dataDstUB(:,:)
00103 real(kind=dbl_kind), allocatable :: dataOutLB(:,:,:)
00104 real(kind=dbl_kind), allocatable :: dataOutUB(:,:,:)
00105
00106 real(kind=dbl_kind), allocatable :: ice_cov_global(:,:)
00107 real(kind=dbl_kind) :: ice_cov_lb(nx_block,ny_block,max_blocks)
00108 real(kind=dbl_kind) :: ice_cov_ub(nx_block,ny_block,max_blocks)
00109 real(kind=dbl_kind) :: ice_cov(nx_block,ny_block,max_blocks)
00110
00111 logical(kind=log_kind) :: regrid
00112
00113 integer(kind=int_kind) :: dateUB, dateLB, secUB, secLB
00114 integer(kind=int_kind) :: nlon
00115 integer(kind=int_kind) :: nlat
00116 integer(kind=int_kind) :: nflds
00117
00118 character(len=char_len_long) :: fldList
00119 character(len=char_len) :: fldName
00120
00121 logical :: read_data
00122
00123
00124 real (kind=dbl_kind), parameter ::
00125 cp_sno = 0.0_dbl_kind
00126 , rLfi = Lfresh*rhoi
00127 , rLfs = Lfresh*rhos
00128 , rLvi = Lvap*rhoi
00129 , rLvs = Lvap*rhos
00130 , rcpi = cp_ice*rhoi
00131 , rcps = cp_sno*rhos
00132 , rcpidepressT = rcpi*depressT
00133 , rLfidepressT = rLfi*depressT
00134
00135
00136
00137 contains
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151
00152
00153
00154 subroutine ice_prescribed_init
00155
00156
00157
00158 use shr_string_mod
00159 use ice_gather_scatter, only : gather_global
00160
00161 implicit none
00162
00163
00164
00165
00166
00167 real(kind=dbl_kind), allocatable :: work4d(:,:,:,:)
00168 integer(kind=int_kind), allocatable :: csim_mask_g(:,:)
00169 real(kind=dbl_kind) :: aice_max
00170 character(len=char_len_long) :: first_data_file
00171 character(len=char_len_long) :: domain_info_fn
00172 character(len=char_len_long) :: data_path
00173 integer(kind=int_kind) :: ndims
00174 character(len=char_len) :: timeName
00175 character(len=char_len) :: latName
00176 character(len=char_len) :: lonName
00177 character(len=char_len) :: maskName
00178 character(len=char_len) :: areaName
00179 logical(kind=log_kind) :: check
00180
00181 integer(kind=int_kind) :: nml_error
00182
00183
00184 character(*),parameter :: subName = "('ice_prescribed_init')"
00185 character(*),parameter :: F00 = "('(ice_prescribed_init) ',4a)"
00186 character(*),parameter :: F01 = "('(ice_prescribed_init) ',a,2i6)"
00187 character(*),parameter :: F02 = "('(ice_prescribed_init) ',a,2g20.13)"
00188
00189
00190 namelist /ice_prescribed_nml/ prescribed_ice, prescribed_ice_fill, &
00191 stream_year_first , stream_year_last , model_year_align, &
00192 stream_fldVarName , stream_fldFileName, &
00193 stream_domTvarName, stream_domXvarName, stream_domYvarName, &
00194 stream_domAreaName, stream_domMaskName, stream_domFileName
00195
00196
00197
00198 prescribed_ice = .false.
00199 stream_year_first = 1
00200 stream_year_last = 1
00201 model_year_align = 1
00202 stream_fldVarName = 'ice_cov'
00203 stream_fldFileName = ' '
00204 stream_domTvarName = 'time'
00205 stream_domXvarName = 'lon'
00206 stream_domYvarName = 'lat'
00207 stream_domAreaName = 'area'
00208 stream_domMaskName = 'mask'
00209 stream_domFileName = ' '
00210 prescribed_ice_fill = .false.
00211
00212
00213 call get_fileunit(nu_nml)
00214 if (my_task == master_task) then
00215 open (nu_nml, file=nml_filename, status='old',iostat=nml_error)
00216 if (nml_error /= 0) then
00217 nml_error = -1
00218 else
00219 nml_error = 1
00220 endif
00221 do while (nml_error > 0)
00222 read(nu_nml, nml=ice_prescribed_nml,iostat=nml_error)
00223 if (nml_error > 0) read(nu_nml,*)
00224 end do
00225 if (nml_error == 0) close(nu_nml)
00226 endif
00227 call release_fileunit(nu_nml)
00228
00229 call broadcast_scalar(nml_error,master_task)
00230
00231 if (nml_error /= 0) then
00232 call abort_ice ('ice: Namelist read error in ice_prescribed_mod')
00233 endif
00234
00235 call broadcast_scalar(prescribed_ice,master_task)
00236 call broadcast_scalar(stream_year_first,master_task)
00237 call broadcast_scalar(stream_year_last,master_task)
00238 call broadcast_scalar(model_year_align,master_task)
00239 call broadcast_scalar(stream_fldVarName,master_task)
00240 call broadcast_scalar(stream_fldFileName,master_task)
00241 call broadcast_scalar(stream_domTvarName,master_task)
00242 call broadcast_scalar(stream_domXvarName,master_task)
00243 call broadcast_scalar(stream_domYvarName,master_task)
00244 call broadcast_scalar(stream_domAreaName,master_task)
00245 call broadcast_scalar(stream_domMaskName,master_task)
00246 call broadcast_scalar(stream_domFileName,master_task)
00247 call broadcast_scalar(prescribed_ice_fill,master_task)
00248
00249 if (.not.prescribed_ice) return
00250
00251 if (my_task == master_task) then
00252 write(nu_diag,*) ' prescribed_ice = ', prescribed_ice
00253 write(nu_diag,*) ' stream_year_first = ', stream_year_first
00254 write(nu_diag,*) ' stream_year_last = ', stream_year_last
00255 write(nu_diag,*) ' model_year_align = ', model_year_align
00256 write(nu_diag,*) ' prescribed_ice_fill = ', prescribed_ice_fill
00257
00258
00259
00260
00261 write (nu_diag,*) 'This is the prescribed ice option.'
00262 write (nu_diag,*) 'Heat and water will not be conserved.'
00263 endif
00264
00265
00266
00267
00268
00269
00270 if (my_task == master_task) then
00271
00272 allocate (work_g1(nx_global,ny_global),work_g2(nx_global,ny_global),&
00273 & csim_mask_g(nx_global,ny_global))
00274
00275 else
00276
00277 allocate (work_g1(1,1),work_g2(1,1),&
00278 & csim_mask_g(1,1))
00279
00280 end if
00281
00282 call gather_global(work_g1, hm, master_task, distrb_info)
00283
00284 if (my_task == master_task) then
00285
00286 csim_mask_g = work_g1
00287
00288 endif
00289
00290 call gather_global(work_g1, TLAT, master_task, distrb_info)
00291 call gather_global(work_g2, TLON, master_task, distrb_info)
00292
00293 if (my_task == master_task) then
00294
00295 work_g1(:,:) = work_g1(:,:)*rad_to_deg
00296 work_g2(:,:) = work_g2(:,:)*rad_to_deg
00297
00298
00299
00300
00301 csim_stream%nFiles = 0
00302 csim_stream%file(:)%name = 'not_set'
00303 csim_stream%file(:)%nt = 0
00304 csim_stream%file(:)%haveData = .false.
00305 csim_stream%yearFirst = stream_year_first
00306 csim_stream%yearLast = stream_year_last
00307 csim_stream%yearAlign = model_year_align
00308 csim_stream%fldListFile = ' '
00309 csim_stream%fldListModel = ' '
00310 csim_stream%FilePath = ' '
00311 csim_stream%domFilePath = ' '
00312 csim_stream%k_lvd = -1
00313 csim_stream%n_lvd = -1
00314 csim_stream%found_lvd = .false.
00315 csim_stream%k_gvd = -1
00316 csim_stream%n_gvd = -1
00317 csim_stream%found_gvd = .false.
00318
00319 csim_stream%dataSource = 'cice ifrac/sst file'
00320 csim_stream%fldListFile = stream_fldVarName
00321 csim_stream%fldListModel = 'ice_cov'
00322
00323
00324 csim_stream%File(1)%name = stream_fldFileName
00325 csim_stream%nFiles = 1
00326
00327 csim_stream%domTvarName = stream_domTvarName
00328 csim_stream%domXvarName = stream_domXvarName
00329 csim_stream%domYvarName = stream_domYvarName
00330 csim_stream%domAreaName = stream_domAreaName
00331 csim_stream%domMaskName = stream_domMaskName
00332 csim_stream%domFileName = stream_domFileName
00333 csim_stream%init = .true.
00334
00335
00336
00337
00338
00339 call shr_stream_getFileFieldList(csim_stream,fldList)
00340 nflds = shr_string_listGetNum(fldList)
00341
00342 if (nflds == 1) then
00343 call shr_string_listGetName(fldList,nflds,fldName)
00344 else
00345 write(nu_diag,F00) "ERROR: Only one field can exist in stream"
00346 write(nu_diag,F01) trim(fldList)//' has ',nflds
00347 call abort_ice(subName)
00348 end if
00349
00350 call shr_stream_getFirstFileName(csim_stream,first_data_file,data_path)
00351 call shr_stream_getFile (data_path,first_data_file)
00352 check = shr_ncread_varExists(first_data_file,fldName)
00353
00354 if (.not.check) then
00355 write(nu_diag,F00) "ERROR: ice concentration field does not exist"
00356 call abort_ice(subName)
00357 end if
00358
00359
00360
00361
00362 call shr_stream_getDomainInfo(csim_stream,data_path,domain_info_fn, timeName, &
00363 & lonName, latName, maskName, areaName)
00364
00365 call shr_ncread_varDimSizes(domain_info_fn,lonName,nlon)
00366 call shr_ncread_varDimSizes(domain_info_fn,latName,nlat)
00367 call shr_stream_getFile(data_path,domain_info_fn)
00368 call shr_sys_flush(nu_diag)
00369
00370 allocate(dataXCoord(nlon,nlat))
00371 allocate(dataYCoord(nlon,nlat))
00372
00373
00374
00375
00376 call shr_ncread_domain(domain_info_fn, lonName, dataXCoord, latName, dataYCoord)
00377
00378 write (nu_diag,F02) 'min/max dataXCoord = ',minval(dataXCoord) ,maxval(dataXCoord)
00379 write (nu_diag,F02) 'min/max dataYCoord = ',minval(dataYCoord) ,maxval(dataYCoord)
00380 write (nu_diag,F02) 'min/max TLON = ',minval(work_g2) ,maxval(work_g2)
00381 write (nu_diag,F02) 'min/max TLAT = ',minval(work_g1) ,maxval(work_g1)
00382 write (nu_diag,F01) 'min/max csim_mask_g = ',minval(csim_mask_g),maxval(csim_mask_g)
00383
00384
00385
00386
00387 if (single_column) then
00388 call ice_open_nc (domain_info_fn, ncid)
00389 call shr_scam_GetCloseLatLon(ncid,scmlat,scmlon,closelat,closelon,latidx,lonidx)
00390 call ice_close_nc (ncid)
00391 regrid=( abs(dataYCoord(lonidx,latidx) - work_g1(1,1))>1.e-12 .or. abs(dataXCoord(lonidx,latidx)- work_g2(1,1))>1.0e-12)
00392 else
00393 call ice_prescribed_checkDomain(work_g2, work_g1, dataXCoord, dataYCoord, regrid)
00394
00395 end if
00396
00397
00398
00399
00400 if (regrid) then
00401 allocate(dataMask(nlon,nlat))
00402 call shr_ncread_domain(domain_info_fn, lonName, dataXCoord, &
00403 & latName, dataYCoord, maskName, dataMask)
00404 write (nu_diag,F01) 'min/max dataMask = ',minval(dataMask), maxval(dataMask)
00405
00406 if (prescribed_ice_fill) then
00407 call shr_map_mapSet(csim_fill, work_g2, work_g1, csim_mask_g, &
00408 & work_g2, work_g1, csim_mask_g, &
00409 & name='csim_map',type='fill',algo='nnoni', &
00410 & mask='nomask')
00411 end if
00412
00413 write (nu_diag,F00) 'Computing mapping weights'
00414
00415 call shr_map_mapSet(csim_map, dataXCoord, dataYCoord, dataMask, &
00416 & work_g2, work_g1, csim_mask_g, &
00417 & name='csim_map',type='remap',algo='bilinear', &
00418 & mask='dstmask',vect='scalar')
00419 deallocate(dataMask)
00420
00421 end if
00422
00423 allocate(dataOutLB(nx_global,ny_global,nflds))
00424 allocate(dataOutUB(nx_global,ny_global,nflds))
00425
00426
00427
00428
00429 allocate(work4d(nlon,nlat,1,1))
00430 call shr_ncread_field4dG(first_data_file, fldName, rfld=work4d, dim3i=1)
00431 aice_max = maxval(work4d)
00432 deallocate(work4d)
00433
00434 if (aice_max > c2) then
00435 write(nu_diag,F02) "ERROR: Ice conc data must be in fraction, aice_max= ",aice_max
00436 call abort_ice(subName)
00437 end if
00438
00439 deallocate(dataXCoord,dataYCoord)
00440 deallocate(work_g1,work_g2)
00441 deallocate(csim_mask_g)
00442
00443 else
00444
00445 deallocate(work_g1,work_g2)
00446 deallocate(csim_mask_g)
00447
00448 end if
00449
00450
00451
00452
00453 if (ncat == 1) then
00454 hin_max(1) = 999._dbl_kind
00455 end if
00456
00457 end subroutine ice_prescribed_init
00458
00459
00460
00461
00462
00463
00464
00465
00466
00467
00468
00469
00470
00471
00472
00473
00474 subroutine ice_prescribed_run(mDateIn, secIn)
00475
00476
00477
00478 use shr_tInterp_mod
00479 use ice_constants, only: field_loc_center, field_type_scalar
00480 use ice_gather_scatter, only : scatter_global
00481
00482 implicit none
00483
00484
00485
00486 integer(kind=int_kind), intent(in) :: mDateIn
00487 integer(kind=int_kind), intent(in) :: secIn
00488
00489
00490
00491 real(kind=dbl_kind) :: fLB
00492 real(kind=dbl_kind) :: fUB
00493 integer(kind=int_kind) :: mDateLB
00494 integer(kind=int_kind) :: mDateUB
00495 integer(kind=int_kind) :: dDateLB
00496 integer(kind=int_kind) :: dDateUB
00497 integer(kind=int_kind) :: secUB
00498 integer(kind=int_kind) :: secLB
00499 integer(kind=int_kind) :: n_lb
00500 integer(kind=int_kind) :: n_ub
00501 character(len=char_len_long) :: fileLB
00502 character(len=char_len_long) :: fileUB
00503
00504 integer(kind=int_kind) :: mDateLB_old = -999
00505 integer(kind=int_kind) :: secLB_old = -999
00506 integer(kind=int_kind) :: i,j,n,icnt,iblk
00507
00508
00509
00510
00511
00512
00513 if (my_task == master_task) then
00514
00515 call shr_stream_findBounds(csim_stream, mDateIn,secIn, &
00516 & mDateLB, dDateLB, secLB, n_lb, fileLB, &
00517 & mDateUB, dDateUB, secUB, n_ub, fileUB)
00518
00519 read_data = .false.
00520 if (mDateLB_old /= mDateLB .or. secLB_old /= secLB) then
00521 allocate(dataInLB(nlon,nlat,nflds))
00522 allocate(dataInUB(nlon,nlat,nflds))
00523 read_data = .true.
00524 call ice_prescribed_readField(fileLB, fldName, n_lb, dataInLB)
00525 call ice_prescribed_readField(fileUB, fldName, n_ub, dataInUB)
00526
00527 if (regrid) then
00528
00529 allocate(dataSrcLB(nflds,nlon*nlat))
00530 allocate(dataSrcUB(nflds,nlon*nlat))
00531 allocate(dataDstLB(nflds,nx_global*ny_global))
00532 allocate(dataDstUB(nflds,nx_global*ny_global))
00533
00534
00535
00536
00537 do n=1,nflds
00538 icnt = 0
00539 do j=1,nlat
00540 do i=1,nlon
00541 icnt = icnt + 1
00542 dataSrcLB(n,icnt) = dataInLB(i,j,n)
00543 dataSrcUB(n,icnt) = dataInUB(i,j,n)
00544 enddo
00545 enddo
00546 enddo
00547
00548
00549
00550
00551 if (prescribed_ice_fill) then
00552 call shr_map_mapData(dataSrcLB, dataDstLB, csim_fill)
00553 call shr_map_mapData(dataSrcUB, dataDstUB, csim_fill)
00554 end if
00555 call shr_map_mapData(dataSrcLB, dataDstLB, csim_map)
00556 call shr_map_mapData(dataSrcUB, dataDstUB, csim_map)
00557
00558
00559
00560
00561 do n=1,nflds
00562 icnt = 0
00563 do j=1,ny_global
00564 do i=1,nx_global
00565 icnt = icnt + 1
00566 dataOutLB(i,j,n) = dataDstLB(n,icnt)
00567 dataOutUB(i,j,n) = dataDstUB(n,icnt)
00568 enddo
00569 enddo
00570 enddo
00571 deallocate(dataSrcLB)
00572 deallocate(dataSrcUB)
00573 deallocate(dataDstLB)
00574 deallocate(dataDstUB)
00575 else
00576 if (single_column) then
00577 call ice_open_nc (fileLB, ncid)
00578 call shr_scam_GetCloseLatLon(ncid,scmlat,scmlon,closelat,closelon,latidx,lonidx)
00579 call ice_close_nc (ncid)
00580 dataOutLB(1,1,:) = dataInLB(lonidx,latidx,:)
00581 call ice_open_nc (fileUB, ncid)
00582 call shr_scam_GetCloseLatLon(ncid,scmlat,scmlon,closelat,closelon,latidx,lonidx)
00583 call ice_close_nc (ncid)
00584 dataOutUB(1,1,:) = dataInUB(lonidx,latidx,:)
00585 else
00586 dataOutLB = dataInLB
00587 dataOutUB = dataInUB
00588 end if
00589 end if
00590
00591 mDateLB_old = mDateLB
00592 secLB_old = secLB
00593
00594 deallocate(dataInUB)
00595 deallocate(dataInLB)
00596 end if
00597
00598 call shr_tInterp_getFactors(mDateLB, secLB, mDateUB, secUB, &
00599 & mDateIn, secIN, fLB, fUB)
00600
00601 end if
00602 call broadcast_scalar(fLB, master_task)
00603 call broadcast_scalar(fUB, master_task)
00604
00605
00606
00607
00608
00609 call broadcast_scalar(read_data, master_task)
00610 if (read_data) then
00611 if (my_task == master_task) then
00612 allocate(ice_cov_global(nx_global,ny_global))
00613 ice_cov_global(:,:) = dataOutLB(:,:,1)
00614 else
00615 allocate(ice_cov_global(1,1))
00616 end if
00617 call scatter_global(ice_cov_lb, ice_cov_global, &
00618 & master_task, distrb_info, &
00619 & field_loc_center, field_type_scalar)
00620
00621 if (my_task == master_task) then
00622 ice_cov_global(:,:) = dataOutUB(:,:,1)
00623 end if
00624 call scatter_global(ice_cov_ub, ice_cov_global, &
00625 & master_task, distrb_info, &
00626 & field_loc_center, field_type_scalar)
00627
00628 deallocate(ice_cov_global)
00629 end if
00630
00631
00632
00633
00634
00635 do iblk = 1,nblocks
00636 do j = 1,ny_block
00637 do i = 1,nx_block
00638 ice_cov(i,j,iblk) = fLB*ice_cov_lb(i,j,iblk) + fUB*ice_cov_ub(i,j,iblk)
00639 enddo
00640 enddo
00641 enddo
00642
00643
00644
00645
00646
00647 call ice_prescribed_phys()
00648
00649 end subroutine ice_prescribed_run
00650
00651
00652
00653
00654
00655
00656
00657
00658
00659
00660
00661
00662
00663
00664 subroutine ice_prescribed_readField(fileName, fldName, nTime, fldRead)
00665
00666
00667
00668 implicit none
00669
00670
00671
00672 character(*), intent(in) :: fileName
00673 character(*), intent(in) :: fldName
00674 integer(kind=int_kind), intent(in) :: nTime
00675 real (kind=dbl_kind), intent(out) :: fldRead(:,:,:)
00676
00677
00678
00679 real(kind=dbl_kind), allocatable :: A4d(:,:,:,:)
00680 integer(kind=int_kind) :: ni, nj
00681
00682 ni = size(fldRead,1)
00683 nj = size(fldRead,2)
00684
00685 allocate(A4d(ni,nj,1,1))
00686
00687 call shr_ncread_field4dG(fileName, fldName, rfld=A4d, dim3i=nTime)
00688 fldRead(:,:,1) = A4d(:,:,1,1)
00689
00690 deallocate(A4d)
00691
00692 end subroutine ice_prescribed_readField
00693
00694
00695
00696
00697
00698
00699
00700
00701
00702
00703
00704
00705
00706
00707 subroutine ice_prescribed_checkDomain(csimXCoord, csimYCoord, dataXCoord, dataYCoord, regrid)
00708
00709
00710
00711 implicit none
00712
00713
00714
00715 real(kind=dbl_kind), intent(in) :: csimXCoord(:,:)
00716 real(kind=dbl_kind), intent(in) :: csimYCoord(:,:)
00717 real(kind=dbl_kind), intent(in) :: dataXCoord(:,:)
00718 real(kind=dbl_kind), intent(in) :: dataYCoord(:,:)
00719 logical(kind=log_kind), intent(out):: regrid
00720
00721
00722
00723 integer(kind=int_kind) :: ni_csim, nj_csim
00724 integer(kind=int_kind) :: ni_data, nj_data
00725 real(kind=dbl_kind) :: max_XDiff, max_YDiff
00726 real(kind=dbl_kind), allocatable :: csimXTemp(:,:)
00727
00728
00729 character(*),parameter :: F00 = "('(ice_prescribed_checkDomain) ',a)"
00730 character(*),parameter :: F01 = "('(ice_prescribed_checkDomain) ',a,2i6)"
00731 character(*),parameter :: F02 = "('(ice_prescribed_checkDomain) ',a,g20.13)"
00732
00733 regrid = .false.
00734
00735
00736
00737 ni_csim = size(csimXCoord,1)
00738 nj_csim = size(csimXCoord,2)
00739 ni_data = size(dataXCoord,1)
00740 nj_data = size(dataXCoord,2)
00741
00742
00743
00744
00745 allocate(csimXTemp(ni_csim,nj_csim))
00746 csimXTemp = csimXCoord
00747 where (csimXTemp >= c360) csimXTemp = csimXTemp - c360
00748 where (csimXTemp < c0 ) csimXTemp = csimXTemp + c360
00749
00750 if (ni_data == ni_csim .and. nj_data == nj_csim) then
00751 max_XDiff = maxval(abs(csimXTemp - dataXCoord))
00752 max_YDiff = maxval(abs(csimYCoord - dataYCoord))
00753 if (max_XDiff > eps04) regrid = .true.
00754 if (max_YDiff > eps04) regrid = .true.
00755
00756 write(nu_diag,F02) 'Maximum X difference = ', max_XDiff
00757 write(nu_diag,F02) 'Maximum Y difference = ', max_YDiff
00758 write(nu_diag,F00) 'CSIM and input data domain sizes match'
00759 write(nu_diag,F01) 'CSIM and data grids are', ni_csim,nj_csim
00760 else
00761 regrid = .true.
00762 end if
00763 deallocate(csimXTemp)
00764
00765 if (regrid) then
00766 write(nu_diag,F00) 'Ice concentration will be interpolated to CSIM grid'
00767 else
00768 write(nu_diag,F00) 'Ice concentration grid and CSIM grid the same'
00769 end if
00770
00771 end subroutine ice_prescribed_checkDomain
00772
00773
00774
00775
00776
00777
00778
00779
00780
00781
00782
00783
00784
00785
00786
00787
00788
00789
00790
00791 subroutine ice_prescribed_phys
00792
00793
00794
00795 use ice_flux
00796
00797 use ice_state
00798 use ice_itd, only : aggregate
00799 use ice_dyn_evp
00800
00801 implicit none
00802
00803
00804
00805
00806
00807
00808 integer(kind=int_kind) :: layer
00809 integer(kind=int_kind) :: nc
00810 integer(kind=int_kind) :: i,j,k
00811 integer(kind=int_kind) :: iblk
00812
00813 real(kind=dbl_kind) :: slope
00814 real(kind=dbl_kind) :: Ti
00815 real(kind=dbl_kind) :: Tmlt
00816 real(kind=dbl_kind) :: qin_save(nilyr)
00817 real(kind=dbl_kind) :: qsn_save(nslyr)
00818 real(kind=dbl_kind) :: hi
00819 real(kind=dbl_kind) :: hs
00820 real(kind=dbl_kind) :: zn
00821 real(kind=dbl_kind) :: salin(nilyr)
00822
00823 real(kind=dbl_kind), parameter :: nsal = 0.407_dbl_kind
00824 real(kind=dbl_kind), parameter :: msal = 0.573_dbl_kind
00825 real(kind=dbl_kind), parameter :: saltmax = 3.2_dbl_kind
00826
00827
00828
00829
00830
00831
00832
00833
00834
00835
00836
00837
00838
00839
00840
00841
00842
00843 do iblk = 1,nblocks
00844 do j = 1,ny_block
00845 do i = 1,nx_block
00846 if (tmask(i,j,iblk)) then
00847 if (ice_cov(i,j,iblk) .lt. eps04) ice_cov(i,j,iblk) = c0
00848 if (ice_cov(i,j,iblk) .gt. c1) ice_cov(i,j,iblk) = c1
00849 else
00850 ice_cov(i,j,iblk) = c0
00851 end if
00852 enddo
00853 enddo
00854 enddo
00855
00856 do iblk = 1,nblocks
00857 do j = 1,ny_block
00858 do i = 1,nx_block
00859
00860 if (tmask(i,j,iblk)) then
00861
00862
00863
00864
00865 if (ice_cov(i,j,iblk) >= eps04) then
00866
00867 hi = 0.0_dbl_kind
00868
00869
00870
00871 if(TLAT(i,j,iblk)*rad_to_deg > 40.0_dbl_kind) then
00872 hi = 2.0_dbl_kind
00873 else if(TLAT(i,j,iblk)*rad_to_deg < -40.0_dbl_kind) then
00874 hi = 1.0_dbl_kind
00875 end if
00876
00877
00878
00879
00880 do nc = 1,ncat
00881
00882 if(hin_max(nc-1) < hi .and. hi < hin_max(nc)) then
00883
00884 if (aicen(i,j,nc,iblk) > c0) then
00885 hs = vsnon(i,j,nc,iblk) / aicen(i,j,nc,iblk)
00886 else
00887 hs = c0
00888 endif
00889
00890 qin_save(:) = c0
00891 qsn_save(:) = c0
00892
00893 if (vicen(i,j,nc,iblk) > c0) then
00894 do k=1,nilyr
00895 qin_save(k) = eicen(i,j,ilyr1(nc)+k-1,iblk) &
00896 * real(nilyr,kind=dbl_kind)
00897 / vicen(i,j,nc,iblk)
00898 enddo
00899 endif
00900
00901 if (vsnon(i,j,nc,iblk) > c0) then
00902 do k=1,nslyr
00903 qsn_save(k) = esnon(i,j,slyr1(nc)+k-1,iblk) &
00904 * real(nslyr,kind=dbl_kind)
00905 / vsnon(i,j,nc,iblk)
00906 enddo
00907 endif
00908
00909 aicen(i,j,nc,iblk) = ice_cov(i,j,iblk)
00910 vicen(i,j,nc,iblk) = hi*aicen(i,j,nc,iblk)
00911 vsnon(i,j,nc,iblk) = hs*aicen(i,j,nc,iblk)
00912
00913
00914 do k=1,nilyr
00915 eicen(i,j,ilyr1(nc)+k-1,iblk) &
00916 = qin_save(k) * vicen(i,j,nc,iblk) &
00917 / real(nilyr,kind=dbl_kind)
00918 enddo
00919
00920 do k=1,nslyr
00921 esnon(i,j,slyr1(nc)+k-1,iblk) &
00922 = qsn_save(k) * vsnon(i,j,nc,iblk) &
00923 / real(nslyr,kind=dbl_kind)
00924 enddo
00925
00926
00927
00928
00929
00930 if (abs(eicen(i,j,ilyr1(nc),iblk)) < puny) then
00931
00932 if (aice(i,j,iblk) < puny) &
00933 trcrn(i,j,nt_Tsfc,nc,iblk) = Tf(i,j,iblk)
00934
00935 slope = Tf(i,j,iblk) - trcrn(i,j,nt_Tsfc,nc,iblk)
00936 do k = 1, nilyr
00937 zn = (real(k,kind=dbl_kind)-p5) / real(nilyr,kind=dbl_kind)
00938 Ti = trcrn(i,j,nt_Tsfc,nc,iblk) + slope*zn
00939 salin(k) = (saltmax/c2)*(c1-cos(pi*zn**(nsal/(msal+zn))))
00940 Tmlt = -salin(k)*depressT
00941 eicen(i,j,ilyr1(nc)+k-1,iblk) = &
00942 -(rhoi * (cp_ice*(Tmlt-Ti) &
00943 + Lfresh*(c1-Tmlt/Ti) - cp_ocn*Tmlt)) &
00944 * vicen(i,j,nc,iblk)/real(nilyr,kind=dbl_kind)
00945 enddo
00946
00947 do k=1,nslyr
00948 esnon(i,j,slyr1(nc)+k-1,iblk) = &
00949 -rhos*(Lfresh - cp_ice*trcrn(i,j,nt_Tsfc,nc,iblk)) &
00950 *vsnon(i,j,nc,iblk)
00951 enddo
00952
00953 endif
00954 end if
00955 enddo
00956 else
00957 trcrn(i,j,nt_Tsfc,:,iblk) = Tf(i,j,iblk)
00958 aicen(i,j,:,iblk) = c0
00959 vicen(i,j,:,iblk) = c0
00960 vsnon(i,j,:,iblk) = c0
00961 esnon(i,j,:,iblk) = c0
00962 eicen(i,j,:,iblk) = c0
00963 end if
00964 end if
00965 enddo
00966 enddo
00967
00968
00969
00970
00971 call aggregate (nx_block, ny_block, &
00972 & aicen(:,:,:,iblk), trcrn(:,:,:,:,iblk), &
00973 & vicen(:,:,:,iblk), vsnon(:,:,:,iblk), &
00974 & eicen(:,:,:,iblk), esnon(:,:,:,iblk), &
00975 & aice(:,:,iblk), trcr(:,:,:,iblk), &
00976 & vice(:,:,iblk), vsno(:,:,iblk), &
00977 & eice(:,:,iblk), esno(:,:,iblk), &
00978 & aice0(:,:,iblk), &
00979 & tmask(:,:,iblk), trcr_depend)
00980
00981 enddo
00982
00983 do iblk = 1, nblocks
00984 do j = 1, ny_block
00985 do i = 1, nx_block
00986 aice_init(i,j,iblk) = aice(i,j,iblk)
00987 enddo
00988 enddo
00989 enddo
00990
00991
00992
00993
00994
00995 frzmlt (:,:,:) = c0
00996 uvel (:,:,:) = c0
00997 vvel (:,:,:) = c0
00998 strocnxT (:,:,:) = c0
00999 strocnyT (:,:,:) = c0
01000
01001
01002
01003
01004 call init_flux_atm
01005 call init_flux_ocn
01006
01007 end subroutine ice_prescribed_phys
01008
01009
01010
01011 end module ice_prescribed_mod
01012
01013