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_init
00024
00025
00026
00027 use ice_kinds_mod
00028 use ice_communicate, only: my_task, master_task, lprint_stats
00029 use ice_domain_size
00030 use ice_constants
00031
00032
00033
00034 implicit none
00035 save
00036
00037
00038
00039 contains
00040
00041
00042
00043
00044
00045
00046
00047
00048 subroutine input_data
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061 use ice_broadcast
00062 use ice_diagnostics
00063 use ice_fileunits
00064 use ice_calendar, only: year_init, istep0, histfreq, histfreq_n, &
00065 dumpfreq, dumpfreq_n, diagfreq, &
00066 npt, dt, xndt_dyn, days_per_year, write_ic
00067 use ice_restart, only: &
00068 restart, restart_dir, restart_file, pointer_file, &
00069 runid, runtype, ice_ic, resttype, restart_format, lcdf64
00070 use ice_history_fields, only: hist_avg, &
00071 history_format, history_dir, history_file, &
00072 incond_dir, incond_file
00073 use ice_exit
00074 use ice_itd, only: kitd, kcatbound
00075 use ice_ocean, only: oceanmixed_ice
00076 use ice_flux, only: Tfrzpt, update_ocn_f
00077 use ice_forcing, only: &
00078 ycycle, fyear_init, dbug, &
00079 atm_data_type, atm_data_dir, precip_units, &
00080 atm_data_format, ocn_data_format, &
00081 sss_data_type, sst_data_type, ocn_data_dir, &
00082 oceanmixed_file, restore_sst, trestore
00083 use ice_grid, only: grid_file, kmt_file, grid_type, grid_format
00084 use ice_mechred, only: kstrength, krdg_partic, krdg_redist
00085 use ice_dyn_evp, only: ndte, kdyn, evp_damping, yield_curve
00086 use ice_shortwave, only: albicev, albicei, albsnowv, albsnowi, &
00087 shortwave, albedo_type, R_ice, R_pnd, &
00088 R_snw, dT_mlt_in, rsnw_melt_in
00089 use ice_atmo, only: atmbndy, calc_strair
00090 use ice_transport_driver, only: advection
00091 use ice_state, only: nt_Tsfc, nt_iage, nt_FY, nt_volpn, nt_aero, &
00092 tr_aero, tr_iage, tr_FY, tr_pond
00093 use ice_aerosol, only: restart_aero
00094 use ice_age, only: restart_age
00095 use ice_FY, only: restart_FY
00096 use ice_meltpond, only: restart_pond
00097 use ice_therm_vertical, only: calc_Tsfc, heat_capacity
00098 use ice_restoring
00099
00100
00101
00102
00103
00104 integer (kind=int_kind) ::
00105 nml_error,
00106 ntr ,
00107 ns
00108
00109 character (len=6) :: chartmp
00110
00111 #ifdef CCSMCOUPLED
00112 logical :: exists
00113 #endif
00114
00115
00116
00117
00118
00119 namelist /setup_nml/ &
00120 days_per_year, year_init, istep0, dt, &
00121 npt, xndt_dyn, &
00122 runtype, runid, &
00123 ice_ic, restart, &
00124 restart_dir, restart_file, restart_format, lcdf64, &
00125 pointer_file, dumpfreq, dumpfreq_n, &
00126 diagfreq, diag_type, diag_file, &
00127 print_global, print_points, latpnt, lonpnt, &
00128 dbug, histfreq, histfreq_n, hist_avg, &
00129 history_dir, history_file, history_format, &
00130 write_ic, incond_dir, incond_file, lprint_stats
00131
00132 namelist /grid_nml/ &
00133 grid_format, grid_type, grid_file, kmt_file, &
00134 kcatbound
00135
00136 namelist /ice_nml/ &
00137 kitd, kdyn, ndte, &
00138 evp_damping, yield_curve, &
00139 kstrength, krdg_partic, krdg_redist, advection, &
00140 heat_capacity, shortwave, albedo_type, &
00141 albicev, albicei, albsnowv, albsnowi, &
00142 R_ice, R_pnd, R_snw, &
00143 dT_mlt_in, rsnw_melt_in, &
00144 atmbndy, fyear_init, ycycle, atm_data_format,&
00145 atm_data_type, atm_data_dir, calc_strair, calc_Tsfc, &
00146 precip_units, Tfrzpt, update_ocn_f, &
00147 oceanmixed_ice, ocn_data_format, sss_data_type, sst_data_type, &
00148 ocn_data_dir, oceanmixed_file, restore_sst, trestore, &
00149 restore_ice
00150
00151 #ifndef CCSMCOUPLED
00152 namelist /tracer_nml/ &
00153 tr_iage, restart_age, &
00154 tr_FY, restart_FY, &
00155 tr_pond, restart_pond, &
00156 tr_aero, restart_aero
00157 #endif
00158
00159
00160
00161
00162
00163 days_per_year = 365
00164 year_init = 0
00165 istep0 = 0
00166
00167 #ifndef CCSMCOUPLED
00168 dt = 3600.0_dbl_kind
00169 #endif
00170 npt = 99999
00171 diagfreq = 24
00172 print_points = .false.
00173 print_global = .true.
00174 lprint_stats = .false.
00175 diag_type = 'stdout'
00176 diag_file = 'ice_diag.d'
00177 histfreq(:) = 'd'
00178 histfreq_n(:) = 1
00179 hist_avg = .true.
00180 history_dir = ' '
00181 history_file = 'iceh'
00182 history_format = 'nc'
00183 write_ic = .false.
00184 incond_dir = history_dir
00185 incond_file = 'iceh_ic'
00186 dumpfreq='y'
00187 dumpfreq_n = 1
00188 restart = .false.
00189 restart_dir = ' '
00190 restart_file = 'iced'
00191 restart_format = 'nc'
00192 lcdf64 = .false.
00193 pointer_file = 'ice.restart_file'
00194 ice_ic = 'default'
00195 grid_format = 'bin'
00196 grid_type = 'rectangular'
00197 grid_file = 'unknown_grid_file'
00198 kmt_file = 'unknown_kmt_file'
00199
00200 kitd = 1
00201 kcatbound = 1
00202 kdyn = 1
00203 xndt_dyn = c1
00204 ndte = 120
00205 evp_damping = .false.
00206 yield_curve = 'ellipse'
00207 kstrength = 1
00208 krdg_partic = 1
00209 krdg_redist = 1
00210 advection = 'remap'
00211 shortwave = 'default'
00212 albedo_type = 'default'
00213 heat_capacity = .true.
00214 calc_Tsfc = .true.
00215 Tfrzpt = 'linear_S'
00216 update_ocn_f = .false.
00217 R_ice = 0.00_dbl_kind
00218 R_pnd = 0.00_dbl_kind
00219 R_snw = 0.00_dbl_kind
00220 dT_mlt_in = 1.50_dbl_kind
00221 rsnw_melt_in = 1500._dbl_kind
00222 albicev = 0.78_dbl_kind
00223 albicei = 0.36_dbl_kind
00224 albsnowv = 0.98_dbl_kind
00225 albsnowi = 0.70_dbl_kind
00226 atmbndy = 'default'
00227
00228 fyear_init = 1900
00229 ycycle = 1
00230 atm_data_format = 'bin'
00231 atm_data_type = 'default'
00232 atm_data_dir = ' '
00233 calc_strair = .true.
00234 precip_units = 'mks'
00235
00236 oceanmixed_ice = .false.
00237 ocn_data_format = 'bin'
00238 sss_data_type = 'default'
00239 sst_data_type = 'default'
00240 ocn_data_dir = ' '
00241 oceanmixed_file = 'unknown_oceanmixed_file'
00242 restore_sst = .false.
00243 trestore = 90
00244 restore_ice = .false.
00245 dbug = .false.
00246
00247 latpnt(1) = 90._dbl_kind
00248 lonpnt(1) = 0._dbl_kind
00249 latpnt(2) = -65._dbl_kind
00250 lonpnt(2) = -45._dbl_kind
00251
00252 #ifndef CCSMCOUPLED
00253 runid = 'unknown'
00254 runtype = 'initial'
00255 #endif
00256
00257
00258
00259
00260
00261 tr_iage = .false.
00262 restart_age = .false.
00263 filename_iage = 'none'
00264 tr_FY = .false.
00265 restart_FY = .false.
00266 filename_FY = 'none'
00267 tr_pond = .false.
00268 restart_pond = .false.
00269 filename_volpn = 'none'
00270 tr_aero = .false.
00271 restart_aero = .false.
00272 filename_aero = 'none'
00273
00274 resttype = 'old'
00275
00276 #ifdef CCSMCOUPLED
00277 if (ntr_iage > 0) tr_iage = .true.
00278 if (ntr_FY > 0) tr_FY = .true.
00279 if (ntr_pond > 0) tr_pond = .true.
00280 if (ntr_aero > 0) tr_aero = .true.
00281 #endif
00282
00283
00284
00285
00286
00287 call get_fileunit(nu_nml)
00288
00289 if (my_task == master_task) then
00290 open (nu_nml, file=nml_filename, status='old',iostat=nml_error)
00291 if (nml_error /= 0) then
00292 nml_error = -1
00293 else
00294 nml_error = 1
00295 endif
00296 do while (nml_error > 0)
00297 print*,'Reading setup_nml'
00298 read(nu_nml, nml=setup_nml,iostat=nml_error)
00299 print*,'Reading grid_nml'
00300 read(nu_nml, nml=grid_nml,iostat=nml_error)
00301 print*,'Reading ice_nml'
00302 read(nu_nml, nml=ice_nml,iostat=nml_error)
00303 #ifndef CCSMCOUPLED
00304 print*,'Reading tracer_nml'
00305 read(nu_nml, nml=tracer_nml,iostat=nml_error)
00306 #endif
00307 if (nml_error > 0) read(nu_nml,*)
00308 end do
00309 if (nml_error == 0) close(nu_nml)
00310 endif
00311 call broadcast_scalar(nml_error, master_task)
00312 if (nml_error /= 0) then
00313 call abort_ice('ice: error reading namelist')
00314 endif
00315
00316 call release_fileunit(nu_nml)
00317
00318
00319
00320 call broadcast_scalar(ice_ic, master_task)
00321
00322
00323
00324
00325
00326 #ifdef CCSMCOUPLED
00327
00328 if (my_task == master_task) then
00329 inquire(file='ice_modelio.nml',exist=exists)
00330 if (exists) then
00331 nu_diag = shr_file_getUnit()
00332 call shr_file_setIO('ice_modelio.nml',nu_diag)
00333 end if
00334
00335
00336
00337 history_file = trim(runid) // ".cice.h"
00338 restart_file = trim(runid) // ".cice.r"
00339 incond_file = trim(runid) // ".cice.i."
00340 end if
00341
00342 if (trim(ice_ic) /= 'default' .and. trim(ice_ic) /= 'none') &
00343 restart = .true.
00344 #else
00345 if (trim(diag_type) == 'file') call get_fileunit(nu_diag)
00346 #endif
00347
00348 if (my_task == master_task) then
00349 if (trim(diag_type) == 'file') then
00350 write(ice_stdout,*) 'Diagnostic output will be in file ',diag_file
00351 open (nu_diag, file=diag_file, status='unknown')
00352 endif
00353 write(nu_diag,*) '--------------------------------'
00354 write(nu_diag,*) ' CICE model diagnostic output '
00355 write(nu_diag,*) '--------------------------------'
00356 write(nu_diag,*) ' '
00357 endif
00358
00359 if (runtype == 'continue') restart = .true.
00360 if (runtype /= 'continue' .and. (restart)) then
00361 if (ice_ic == 'none' .or. ice_ic == 'default') then
00362 if (my_task == master_task) then
00363 write(nu_diag,*) &
00364 'WARNING: runtype, restart, ice_ic are inconsistent:'
00365 write(nu_diag,*) runtype, restart, ice_ic
00366 write(nu_diag,*) &
00367 'WARNING: Need ice_ic = <filename>.'
00368 write(nu_diag,*) &
00369 'WARNING: Initializing using ice_ic conditions'
00370 endif
00371 restart = .false.
00372 endif
00373 endif
00374 if (runtype == 'initial' .and. .not.(restart)) then
00375 if (ice_ic /= 'none' .and. ice_ic /= 'default') then
00376 if (my_task == master_task) then
00377 write(nu_diag,*) &
00378 'WARNING: runtype, restart, ice_ic are inconsistent:'
00379 write(nu_diag,*) runtype, restart, ice_ic
00380 write(nu_diag,*) &
00381 'WARNING: Initializing with NO ICE: '
00382 endif
00383 ice_ic = 'none'
00384 endif
00385 endif
00386
00387 #ifndef ncdf
00388
00389 history_format = 'bin'
00390 restart_format = 'bin'
00391 grid_format = 'bin'
00392 atm_data_format = 'bin'
00393 ocn_data_format = 'bin'
00394 #endif
00395
00396 if (days_per_year /= 365) shortwave = 'default'
00397
00398 chartmp = advection(1:6)
00399 if (chartmp /= 'upwind' .and. chartmp /= 'remap ') advection = 'remap'
00400
00401 if (ncat == 1 .and. kitd == 1) then
00402 write (nu_diag,*) 'Remapping the ITD is not allowed for ncat=1'
00403 write (nu_diag,*) 'Using the delta function ITD option instead'
00404 kitd = 0
00405 endif
00406
00407 if (trim(atm_data_type) == 'monthly' .and. calc_strair) &
00408 calc_strair = .false.
00409
00410 if (trim(atm_data_type) == 'hadgem' .and. &
00411 trim(precip_units) /= 'mks') then
00412 if (my_task == master_task) &
00413 write (nu_diag,*) &
00414 'WARNING: HadGEM atmospheric data chosen with wrong precip_units'
00415 write (nu_diag,*) &
00416 'WARNING: Changing precip_units to mks (i.e. kg/m2 s).'
00417 precip_units='mks'
00418 endif
00419
00420 call broadcast_scalar(days_per_year, master_task)
00421 call broadcast_scalar(year_init, master_task)
00422 call broadcast_scalar(istep0, master_task)
00423 call broadcast_scalar(dt, master_task)
00424 call broadcast_scalar(npt, master_task)
00425 call broadcast_scalar(diagfreq, master_task)
00426 call broadcast_scalar(print_points, master_task)
00427 call broadcast_scalar(print_global, master_task)
00428 call broadcast_scalar(diag_type, master_task)
00429 call broadcast_scalar(diag_file, master_task)
00430 call broadcast_scalar(history_format, master_task)
00431
00432
00433 do ns=1,max_nstrm
00434 call broadcast_scalar(histfreq(ns), master_task)
00435 call broadcast_scalar(histfreq_n(ns), master_task)
00436 enddo
00437
00438 call broadcast_scalar(hist_avg, master_task)
00439 call broadcast_scalar(history_dir, master_task)
00440 call broadcast_scalar(history_file, master_task)
00441 call broadcast_scalar(write_ic, master_task)
00442 call broadcast_scalar(incond_dir, master_task)
00443 call broadcast_scalar(incond_file, master_task)
00444 call broadcast_scalar(dumpfreq, master_task)
00445 call broadcast_scalar(dumpfreq_n, master_task)
00446 call broadcast_scalar(restart_file, master_task)
00447 call broadcast_scalar(restart, master_task)
00448 call broadcast_scalar(restart_dir, master_task)
00449 call broadcast_scalar(restart_format, master_task)
00450 call broadcast_scalar(lcdf64, master_task)
00451 call broadcast_scalar(pointer_file, master_task)
00452 call broadcast_scalar(grid_format, master_task)
00453 call broadcast_scalar(grid_type, master_task)
00454 call broadcast_scalar(grid_file, master_task)
00455 call broadcast_scalar(kmt_file, master_task)
00456 call broadcast_scalar(kitd, master_task)
00457 call broadcast_scalar(kcatbound, master_task)
00458 call broadcast_scalar(kdyn, master_task)
00459 call broadcast_scalar(xndt_dyn, master_task)
00460 call broadcast_scalar(ndte, master_task)
00461 call broadcast_scalar(evp_damping, master_task)
00462 call broadcast_scalar(yield_curve, master_task)
00463 call broadcast_scalar(kstrength, master_task)
00464 call broadcast_scalar(krdg_partic, master_task)
00465 call broadcast_scalar(krdg_redist, master_task)
00466 call broadcast_scalar(advection, master_task)
00467 call broadcast_scalar(shortwave, master_task)
00468 call broadcast_scalar(albedo_type, master_task)
00469 call broadcast_scalar(heat_capacity, master_task)
00470 call broadcast_scalar(R_ice, master_task)
00471 call broadcast_scalar(R_pnd, master_task)
00472 call broadcast_scalar(R_snw, master_task)
00473 call broadcast_scalar(dT_mlt_in, master_task)
00474 call broadcast_scalar(rsnw_melt_in, master_task)
00475 call broadcast_scalar(albicev, master_task)
00476 call broadcast_scalar(albicei, master_task)
00477 call broadcast_scalar(albsnowv, master_task)
00478 call broadcast_scalar(albsnowi, master_task)
00479 call broadcast_scalar(atmbndy, master_task)
00480 call broadcast_scalar(fyear_init, master_task)
00481 call broadcast_scalar(ycycle, master_task)
00482 call broadcast_scalar(atm_data_format, master_task)
00483 call broadcast_scalar(atm_data_type, master_task)
00484 call broadcast_scalar(atm_data_dir, master_task)
00485 call broadcast_scalar(calc_strair, master_task)
00486 call broadcast_scalar(calc_Tsfc, master_task)
00487 call broadcast_scalar(Tfrzpt, master_task)
00488 call broadcast_scalar(update_ocn_f, master_task)
00489 call broadcast_scalar(precip_units, master_task)
00490 call broadcast_scalar(oceanmixed_ice, master_task)
00491 call broadcast_scalar(ocn_data_format, master_task)
00492 call broadcast_scalar(sss_data_type, master_task)
00493 call broadcast_scalar(sst_data_type, master_task)
00494 call broadcast_scalar(ocn_data_dir, master_task)
00495 call broadcast_scalar(oceanmixed_file, master_task)
00496 call broadcast_scalar(restore_sst, master_task)
00497 call broadcast_scalar(trestore, master_task)
00498 call broadcast_scalar(restore_ice, master_task)
00499 call broadcast_scalar(dbug, master_task)
00500 call broadcast_array (latpnt(1:2), master_task)
00501 call broadcast_array (lonpnt(1:2), master_task)
00502 call broadcast_scalar(runid, master_task)
00503 call broadcast_scalar(runtype, master_task)
00504 if (dbug) &
00505 call broadcast_scalar(nu_diag, master_task)
00506
00507 call broadcast_scalar(tr_iage, master_task)
00508 call broadcast_scalar(restart_age, master_task)
00509 call broadcast_scalar(tr_FY, master_task)
00510 call broadcast_scalar(restart_FY, master_task)
00511 call broadcast_scalar(tr_pond, master_task)
00512 call broadcast_scalar(restart_pond, master_task)
00513 call broadcast_scalar(tr_aero, master_task)
00514 call broadcast_scalar(restart_aero, master_task)
00515
00516
00517
00518
00519
00520 if (my_task == master_task) then
00521
00522 write(nu_diag,*) ' Document ice_in namelist parameters:'
00523 write(nu_diag,*) ' ==================================== '
00524 write(nu_diag,*) ' '
00525 if (trim(runid) /= 'unknown') &
00526 write(nu_diag,*) ' runid = ', &
00527 trim(runid)
00528 write(nu_diag,1030) ' runtype = ', &
00529 trim(runtype)
00530 write(nu_diag,1020) ' days_per_year = ', days_per_year
00531 write(nu_diag,1020) ' year_init = ', year_init
00532 write(nu_diag,1020) ' istep0 = ', istep0
00533 write(nu_diag,1000) ' dt = ', dt
00534 write(nu_diag,1020) ' npt = ', npt
00535 write(nu_diag,1020) ' diagfreq = ', diagfreq
00536 write(nu_diag,1010) ' print_global = ', &
00537 print_global
00538 write(nu_diag,1010) ' print_points = ', &
00539 print_points
00540 write(nu_diag,1050) ' histfreq = ', histfreq(:)
00541 write(nu_diag,1040) ' histfreq_n = ', histfreq_n(:)
00542 write(nu_diag,1010) ' hist_avg = ', hist_avg
00543 if (hist_avg) then
00544 write (nu_diag,*) 'History data will be averaged over ', &
00545 histfreq_n,' ',histfreq
00546 else
00547 write (nu_diag,*) 'History data will be snapshots'
00548 endif
00549 write(nu_diag,*) ' history_dir = ', &
00550 trim(history_dir)
00551 write(nu_diag,*) ' history_file = ', &
00552 trim(history_file)
00553 write(nu_diag,*) ' history_format = ', &
00554 trim(history_format)
00555 if (write_ic) then
00556 write (nu_diag,*) 'Initial condition will be written in ', &
00557 trim(incond_dir)
00558 endif
00559 write(nu_diag,1030) ' dumpfreq = ', &
00560 trim(dumpfreq)
00561 write(nu_diag,1020) ' dumpfreq_n = ', dumpfreq_n
00562 write(nu_diag,1010) ' restart = ', restart
00563 write(nu_diag,*) ' restart_dir = ', &
00564 trim(restart_dir)
00565 write(nu_diag,*) ' restart_file = ', &
00566 trim(restart_file)
00567 write(nu_diag,*) ' restart_format = ', &
00568 trim(restart_format)
00569 write(nu_diag,*) ' lcdf64 = ', &
00570 lcdf64
00571 write(nu_diag,*) ' pointer_file = ', &
00572 trim(pointer_file)
00573 write(nu_diag,* ) ' ice_ic = ', &
00574 trim(ice_ic)
00575 write(nu_diag,*) ' grid_type = ', &
00576 trim(grid_type)
00577 if (trim(grid_type) /= 'rectangular' .or. &
00578 trim(grid_type) /= 'column') then
00579 write(nu_diag,*) ' grid_file = ', &
00580 trim(grid_file)
00581 write(nu_diag,*) ' kmt_file = ', &
00582 trim(kmt_file)
00583 endif
00584
00585 write(nu_diag,1020) ' kitd = ', kitd
00586 write(nu_diag,1020) ' kcatbound = ', &
00587 kcatbound
00588 write(nu_diag,1020) ' kdyn = ', kdyn
00589 write(nu_diag,1000) ' xndt_dyn = ', xndt_dyn
00590 write(nu_diag,1020) ' ndte = ', ndte
00591 write(nu_diag,1010) ' evp_damping = ', &
00592 evp_damping
00593 write(nu_diag,*) ' yield_curve = ', &
00594 trim(yield_curve)
00595 write(nu_diag,1020) ' kstrength = ', kstrength
00596 write(nu_diag,1020) ' krdg_partic = ', &
00597 krdg_partic
00598 write(nu_diag,1020) ' krdg_redist = ', &
00599 krdg_redist
00600 write(nu_diag,1030) ' advection = ', &
00601 trim(advection)
00602 write(nu_diag,1030) ' shortwave = ', &
00603 trim(shortwave)
00604 write(nu_diag,1030) ' albedo_type = ', &
00605 trim(albedo_type)
00606 write(nu_diag,1000) ' R_ice = ', R_ice
00607 write(nu_diag,1000) ' R_pnd = ', R_pnd
00608 write(nu_diag,1000) ' R_snw = ', R_snw
00609 write(nu_diag,1000) ' dT_mlt_in = ', dT_mlt_in
00610 write(nu_diag,1000) ' rsnw_melt_in = ', rsnw_melt_in
00611 write(nu_diag,1000) ' albicev = ', albicev
00612 write(nu_diag,1000) ' albicei = ', albicei
00613 write(nu_diag,1000) ' albsnowv = ', albsnowv
00614 write(nu_diag,1000) ' albsnowi = ', albsnowi
00615 write(nu_diag,1010) ' heat_capacity = ', &
00616 heat_capacity
00617 write(nu_diag,1030) ' atmbndy = ', &
00618 trim(atmbndy)
00619
00620 write(nu_diag,1020) ' fyear_init = ', &
00621 fyear_init
00622 write(nu_diag,1020) ' ycycle = ', ycycle
00623 write(nu_diag,*) ' atm_data_type = ', &
00624 trim(atm_data_type)
00625 write(nu_diag,1010) ' calc_strair = ', calc_strair
00626 write(nu_diag,1010) ' calc_Tsfc = ', calc_Tsfc
00627 write(nu_diag,*) ' Tfrzpt = ', trim(Tfrzpt)
00628 write(nu_diag,1010) ' update_ocn_f = ', update_ocn_f
00629 if (trim(atm_data_type) /= 'default') then
00630 write(nu_diag,*) ' atm_data_dir = ', &
00631 trim(atm_data_dir)
00632 write(nu_diag,*) ' precip_units = ', &
00633 trim(precip_units)
00634 endif
00635
00636 write(nu_diag,1010) ' oceanmixed_ice = ', &
00637 oceanmixed_ice
00638 write(nu_diag,*) ' sss_data_type = ', &
00639 trim(sss_data_type)
00640 write(nu_diag,*) ' sst_data_type = ', &
00641 trim(sst_data_type)
00642 if (trim(sss_data_type) /= 'default' .or. &
00643 trim(sst_data_type) /= 'default') then
00644 write(nu_diag,*) ' ocn_data_dir = ', &
00645 trim(ocn_data_dir)
00646 endif
00647 if (trim(sss_data_type) == 'ncar' .or. &
00648 trim(sst_data_type) == 'ncar') then
00649 write(nu_diag,*) ' oceanmixed_file = ', &
00650 trim(oceanmixed_file)
00651 endif
00652
00653 #ifdef coupled
00654 if( oceanmixed_ice ) then
00655 write (nu_diag,*) 'WARNING WARNING WARNING WARNING '
00656 write (nu_diag,*) '*Coupled and oceanmixed flags are *'
00657 write (nu_diag,*) '*BOTH ON. Ocean data received from*'
00658 write (nu_diag,*) '*coupler will be altered by mixed *'
00659 write (nu_diag,*) '*layer routine! *'
00660 write (nu_diag,*) ' '
00661 endif
00662 #endif
00663
00664 write (nu_diag,*) ' '
00665 write (nu_diag,'(a30,2f8.2)') 'Diagnostic point 1: lat, lon =', &
00666 latpnt(1), lonpnt(1)
00667 write (nu_diag,'(a30,2f8.2)') 'Diagnostic point 2: lat, lon =', &
00668 latpnt(2), lonpnt(2)
00669
00670
00671 write(nu_diag,1010) ' tr_iage = ', tr_iage
00672 write(nu_diag,1010) ' restart_age = ', restart_age
00673 write(nu_diag,1010) ' tr_FY = ', tr_FY
00674 write(nu_diag,1010) ' restart_FY = ', restart_FY
00675 write(nu_diag,1010) ' tr_pond = ', tr_pond
00676 write(nu_diag,1010) ' restart_pond = ', restart_pond
00677 write(nu_diag,1010) ' tr_aero = ', tr_aero
00678 write(nu_diag,1010) ' restart_aero = ', restart_aero
00679
00680 ntr = 1
00681
00682 if (tr_iage) then
00683 ntr = ntr + 1
00684 nt_iage = nt_Tsfc + 1
00685 write(nu_diag,1020) ' nt_iage = ', nt_iage
00686 else
00687 nt_iage = 1
00688 endif
00689 if (tr_FY) then
00690 ntr = ntr + 1
00691 nt_FY = nt_iage + 1
00692 write(nu_diag,1020) ' nt_FY = ', nt_FY
00693 else
00694 nt_FY = nt_iage
00695 endif
00696 if (tr_pond) then
00697 ntr = ntr + 1
00698 nt_volpn = nt_FY + 1
00699 write(nu_diag,1020) ' nt_volpn = ', nt_volpn
00700 else
00701 nt_volpn = nt_FY
00702 endif
00703 if (tr_aero) then
00704 ntr = ntr + n_aero*4
00705
00706 nt_aero = nt_volpn + 1
00707 write(nu_diag,1020) ' nt_aero = ', nt_aero
00708 else
00709 nt_aero = nt_volpn
00710 endif
00711
00712 if (ntr /= ntrcr) &
00713 write(nu_diag,*) 'WARNING: ntrcr > number of tracers requested'
00714 if (ntr > ntrcr) then
00715 write(nu_diag,*) 'ntrcr < number of namelist tracers'
00716 call abort_ice('ntrcr < number of namelist tracers')
00717 endif
00718
00719 1000 format (a30,2x,f9.2)
00720 1010 format (a30,2x,l6)
00721 1020 format (a30,2x,i6)
00722 1030 format (a30, a8)
00723 1040 format (a30,2x,6i6)
00724 1050 format (a30, 6a8)
00725
00726 write (nu_diag,*) ' '
00727 if (grid_type /= 'displaced_pole' .and. &
00728 grid_type /= 'tripole' .and. &
00729 grid_type /= 'column' .and. &
00730 grid_type /= 'rectangular' .and. &
00731 grid_type /= 'panarctic' .and. &
00732 grid_type /= 'latlon' ) then
00733 call abort_ice('ice_init: unknown grid_type')
00734 endif
00735
00736 endif
00737
00738 call broadcast_scalar(nt_iage, master_task)
00739 call broadcast_scalar(nt_FY, master_task)
00740 call broadcast_scalar(nt_volpn, master_task)
00741 call broadcast_scalar(nt_aero, master_task)
00742
00743 end subroutine input_data
00744
00745
00746
00747
00748
00749
00750
00751
00752 subroutine init_state
00753
00754
00755
00756
00757
00758
00759
00760
00761
00762
00763
00764
00765 use ice_blocks
00766 use ice_domain
00767 use ice_flux, only: sst, Tf, Tair
00768 use ice_grid
00769 use ice_state
00770 use ice_itd
00771 use ice_exit
00772 use ice_therm_vertical, only: heat_capacity
00773
00774
00775
00776
00777
00778 integer (kind=int_kind) ::
00779 i, j ,
00780 it, n ,
00781 iblk
00782
00783
00784
00785
00786
00787 if (my_task == master_task) then
00788
00789 if (nilyr < 1) then
00790 write (nu_diag,*) 'nilyr =', nilyr
00791 write (nu_diag,*) 'Must have at least one ice layer'
00792 call abort_ice ('ice_init: Not enough ice layers')
00793 endif
00794
00795 if (nslyr < 1) then
00796 write (nu_diag,*) 'nslyr =', nslyr
00797 write (nu_diag,*) 'Must have at least one snow layer'
00798 call abort_ice('ice_init: Not enough snow layers')
00799 endif
00800
00801 if (.not.heat_capacity) then
00802
00803 write (nu_diag,*) 'WARNING - Zero-layer thermodynamics'
00804
00805 if (nilyr > 1) then
00806 write (nu_diag,*) 'nilyr =', nilyr
00807 write (nu_diag,*) &
00808 'Must have nilyr = 1 if heat_capacity = F'
00809 call abort_ice('ice_init: Too many ice layers')
00810 endif
00811
00812 if (nslyr > 1) then
00813 write (nu_diag,*) 'nslyr =', nslyr
00814 write (nu_diag,*) &
00815 'Must have nslyr = 1 if heat_capacity = F'
00816 call abort_ice('ice_init: Too many snow layers')
00817 endif
00818
00819 endif
00820
00821 endif
00822
00823
00824
00825
00826
00827 trcr_depend(nt_Tsfc) = 0
00828 if (tr_iage) trcr_depend(nt_iage) = 1
00829 if (tr_FY) trcr_depend(nt_FY) = 0
00830 if (tr_pond) trcr_depend(nt_volpn) = 0
00831 if (tr_aero) then
00832 do n=1,n_aero
00833 trcr_depend(nt_aero+(n-1)*4 :nt_aero+(n-1)*4+1) = 2
00834 trcr_depend(nt_aero+(n-1)*4+2:nt_aero+(n-1)*4+3) = 1
00835 enddo
00836 endif
00837
00838
00839 do iblk = 1, nblocks
00840
00841
00842
00843
00844
00845 call set_state_var (nx_block, ny_block, &
00846 tmask(:,:, iblk), &
00847 ULON (:,:, iblk), ULAT (:,:, iblk), &
00848 Tair (:,:, iblk), sst (:,:, iblk), &
00849 Tf (:,:, iblk), trcr_depend, &
00850 aicen(:,:, :,iblk), trcrn(:,:,:,:,iblk), &
00851 vicen(:,:, :,iblk), vsnon(:,:, :,iblk), &
00852 eicen(:,:, :,iblk), esnon(:,:, :,iblk))
00853
00854
00855
00856
00857
00858 aice(:,:,iblk) = c0
00859 vice(:,:,iblk) = c0
00860 vsno(:,:,iblk) = c0
00861 eice(:,:,iblk) = c0
00862 esno(:,:,iblk) = c0
00863 do it = 1, ntrcr
00864 trcr(:,:,it,iblk) = c0
00865 enddo
00866
00867 call aggregate (nx_block, ny_block, &
00868 aicen(:,:,:,iblk), &
00869 trcrn(:,:,:,:,iblk), &
00870 vicen(:,:,:,iblk), &
00871 vsnon(:,:,:,iblk), &
00872 eicen(:,:,:,iblk), &
00873 esnon(:,:,:,iblk), &
00874 aice (:,:, iblk), &
00875 trcr (:,:,:,iblk), &
00876 vice (:,:, iblk), &
00877 vsno (:,:, iblk), &
00878 eice (:,:, iblk), &
00879 esno (:,:, iblk), &
00880 aice0(:,:, iblk), &
00881 tmask(:,:, iblk), &
00882 trcr_depend)
00883
00884 aice_init(:,:,iblk) = aice(:,:,iblk)
00885
00886 enddo
00887
00888
00889
00890
00891
00892
00893 call bound_state (aicen, trcrn, &
00894 vicen, vsnon, &
00895 eicen, esnon)
00896
00897 end subroutine init_state
00898
00899
00900
00901
00902
00903
00904
00905
00906 subroutine set_state_var (nx_block, ny_block, &
00907 tmask, ULON, ULAT, &
00908 Tair, sst, &
00909 Tf, trcr_depend, &
00910 aicen, trcrn, &
00911 vicen, vsnon, &
00912 eicen, esnon)
00913
00914
00915
00916
00917
00918
00919
00920
00921
00922
00923
00924
00925 use ice_state, only: nt_Tsfc
00926 use ice_therm_vertical, only: heat_capacity, calc_Tsfc, Tmlt
00927 use ice_itd, only: ilyr1, slyr1, hin_max
00928 use ice_grid, only: grid_type
00929 use ice_restart, only: ice_ic
00930
00931
00932
00933 integer (kind=int_kind), intent(in) ::
00934 nx_block, ny_block
00935
00936 logical (kind=log_kind), dimension (nx_block,ny_block),
00937 intent(in) ::
00938 tmask
00939
00940 real (kind=dbl_kind), dimension (nx_block,ny_block),
00941 intent(in) ::
00942 ULON ,
00943 ULAT
00944
00945 real (kind=dbl_kind), dimension (nx_block,ny_block), intent(in) ::
00946 Tair ,
00947 Tf ,
00948 sst
00949
00950 integer (kind=int_kind), dimension (ntrcr), intent(inout) ::
00951 trcr_depend
00952
00953 real (kind=dbl_kind), dimension (nx_block,ny_block,ncat),
00954 intent(out) ::
00955 aicen ,
00956 vicen ,
00957 vsnon
00958
00959 real (kind=dbl_kind), dimension (nx_block,ny_block,ntrcr,ncat),
00960 intent(out) ::
00961 trcrn
00962
00963
00964 real (kind=dbl_kind), dimension (nx_block,ny_block,ntilyr),
00965 intent(out) ::
00966 eicen
00967
00968 real (kind=dbl_kind), dimension (nx_block,ny_block,ntslyr),
00969 intent(out) ::
00970 esnon
00971
00972
00973
00974 integer (kind=int_kind) ::
00975 i, j ,
00976 ij ,
00977 k ,
00978 n ,
00979 it ,
00980 icells
00981
00982 integer (kind=int_kind), dimension(nx_block*ny_block) ::
00983 indxi, indxj
00984
00985 real (kind=dbl_kind) ::
00986 slope, Ti, sum, hbar,
00987 ainit(ncat),
00988 hinit(ncat)
00989
00990 real (kind=dbl_kind), parameter ::
00991 hsno_init = 0.20_dbl_kind ,
00992 edge_init_nh = 70._dbl_kind,
00993 edge_init_sh = -60._dbl_kind
00994
00995
00996 indxi(:) = 0
00997 indxj(:) = 0
00998
00999
01000
01001
01002 do n = 1, ncat
01003 do j = 1, ny_block
01004 do i = 1, nx_block
01005 aicen(i,j,n) = c0
01006 vicen(i,j,n) = c0
01007 vsnon(i,j,n) = c0
01008 trcrn(i,j,nt_Tsfc,n) = Tf(i,j)
01009 if (ntrcr >= 2) then
01010 do it = 2, ntrcr
01011 trcrn(i,j,it,n) = c0
01012 enddo
01013 endif
01014 enddo
01015 enddo
01016 enddo
01017 eicen(:,:,:) = c0
01018 esnon(:,:,:) = c0
01019
01020 if (trim(ice_ic) == 'default') then
01021
01022
01023
01024
01025
01026
01027
01028
01029
01030 hbar = c3
01031
01032
01033
01034 sum = c0
01035 do n = 1, ncat
01036 if (n < ncat) then
01037 hinit(n) = p5*(hin_max(n-1) + hin_max(n))
01038 else
01039 hinit(n) = (hin_max(n-1) + c1)
01040 endif
01041
01042 ainit(n) = max(c0, (c2*hbar*hinit(n) - hinit(n)**2))
01043 sum = sum + ainit(n)
01044 enddo
01045 do n = 1, ncat
01046 ainit(n) = ainit(n) / (sum + puny/ncat)
01047 enddo
01048
01049 if (trim(grid_type) == 'rectangular') then
01050
01051
01052 icells = 0
01053 do j = 1, ny_block
01054 do i = 1, nx_block
01055 if (tmask(i,j)) then
01056 if (ULON(i,j) < -50./rad_to_deg) then
01057 icells = icells + 1
01058 indxi(icells) = i
01059 indxj(icells) = j
01060 endif
01061 endif
01062 enddo
01063 enddo
01064
01065 else
01066
01067
01068 icells = 0
01069 do j = 1, ny_block
01070 do i = 1, nx_block
01071 if (tmask(i,j)) then
01072
01073 if ( (sst (i,j) <= Tf(i,j)+p2) .and. &
01074 (ULAT(i,j) < edge_init_sh/rad_to_deg .or. &
01075 ULAT(i,j) > edge_init_nh/rad_to_deg) ) then
01076 icells = icells + 1
01077 indxi(icells) = i
01078 indxj(icells) = j
01079 endif
01080 endif
01081 enddo
01082 enddo
01083
01084 endif
01085
01086 do n = 1, ncat
01087
01088
01089
01090
01091
01092 do ij = 1, icells
01093 i = indxi(ij)
01094 j = indxj(ij)
01095
01096 aicen(i,j,n) = ainit(n)
01097 vicen(i,j,n) = hinit(n) * ainit(n)
01098 vsnon(i,j,n) =min(aicen(i,j,n)*hsno_init,p2*vicen(i,j,n))
01099 enddo
01100
01101
01102 if (calc_Tsfc) then
01103
01104 do ij = 1, icells
01105 i = indxi(ij)
01106 j = indxj(ij)
01107 trcrn(i,j,nt_Tsfc,n) = min(Tsmelt, Tair(i,j) - Tffresh)
01108 enddo
01109
01110 else
01111
01112 do ij = 1, icells
01113 i = indxi(ij)
01114 j = indxj(ij)
01115 trcrn(i,j,nt_Tsfc,n) = Tf(i,j)
01116 enddo
01117
01118 endif
01119
01120
01121
01122 if (heat_capacity) then
01123
01124
01125 do k = 1, nilyr
01126 do ij = 1, icells
01127 i = indxi(ij)
01128 j = indxj(ij)
01129
01130
01131 slope = Tf(i,j) - trcrn(i,j,nt_Tsfc,n)
01132 Ti = trcrn(i,j,nt_Tsfc,n) &
01133 + slope*(real(k,kind=dbl_kind)-p5)
01134 /real(nilyr,kind=dbl_kind)
01135
01136 eicen(i,j,ilyr1(n)+k-1) = &
01137 -(rhoi * (cp_ice*(Tmlt(k)-Ti) &
01138 + Lfresh*(c1-Tmlt(k)/Ti) - cp_ocn*Tmlt(k))) &
01139 * vicen(i,j,n)/real(nilyr,kind=dbl_kind)
01140 enddo
01141 enddo
01142
01143
01144 do k = 1, nslyr
01145 do ij = 1, icells
01146 i = indxi(ij)
01147 j = indxj(ij)
01148
01149 Ti = min(c0, trcrn(i,j,nt_Tsfc,n))
01150 esnon(i,j,slyr1(n)+k-1) = -rhos*(Lfresh - cp_ice*Ti) &
01151 *vsnon(i,j,n) &
01152 /real(nslyr,kind=dbl_kind)
01153 enddo
01154 enddo
01155
01156 else
01157
01158
01159 k = 1
01160
01161 do ij = 1, icells
01162 i = indxi(ij)
01163 j = indxj(ij)
01164 eicen(i,j,ilyr1(n)+k-1) = &
01165 - rhoi * Lfresh * vicen(i,j,n)
01166 enddo
01167
01168
01169 do ij = 1, icells
01170 i = indxi(ij)
01171 j = indxj(ij)
01172 esnon(i,j,slyr1(n)+k-1) = &
01173 - rhos * Lfresh * vsnon(i,j,n)
01174 enddo
01175
01176 endif
01177 enddo
01178 endif
01179
01180 end subroutine set_state_var
01181
01182
01183
01184 end module ice_init
01185
01186