Changeset 8239 for branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO
- Timestamp:
- 2017-06-28T17:55:50+02:00 (7 years ago)
- Location:
- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO
- Files:
-
- 14 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/ice.F90
r8233 r8239 108 108 !! v_s | - | Snow volume per unit area | m | 109 109 !! smv_i | - | Sea ice salt content | ppt.m | 110 !! oa_i ! - ! Sea ice areal age content | day|110 !! oa_i ! - ! Sea ice areal age content | s | 111 111 !! e_i ! - ! Ice enthalpy | J/m2 | 112 112 !! - ! q_i_1d ! Ice enthalpy per unit vol. | J/m3 | … … 123 123 !! sm_i ! sm_i_1d | Sea ice bulk salinity ! ppt | 124 124 !! s_i ! s_i_1d | Sea ice salinity profile ! ppt | 125 !! o_i ! - | Sea ice Age ! days|125 !! o_i ! - | Sea ice Age ! s | 126 126 !! t_i ! t_i_1d | Sea ice temperature ! K | 127 127 !! t_s ! t_s_1d | Snow temperature ! K | … … 328 328 ! END MV MP 2016 329 329 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_spr !: snow precipitation on ice [kg.m-2.s-1] 330 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_sub !: snow/ice sublimation [kg.m-2.s-1] 330 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_sub !: sublimation of snow/ice [kg.m-2.s-1] 331 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_snw_sub !: snow sublimation [kg.m-2.s-1] 332 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_ice_sub !: ice sublimation [kg.m-2.s-1] 333 334 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_snw_dyn !: dynamical component of wfx_snw [kg.m-2.s-1] 335 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_snw_sum !: surface melt component of wfx_snw [kg.m-2.s-1] 331 336 332 337 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wfx_ice !: ice-ocean mass exchange [kg.m-2.s-1] … … 395 400 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: smv_i !: Sea-Ice Bulk salinity times volume per area (ppt.m) 396 401 ! ! this is an extensive variable that has to be transported 397 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: o_i !: Sea-Ice Age ( days)398 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: oa_i !: Sea-Ice Age times ice area ( days)402 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: o_i !: Sea-Ice Age (s) 403 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: oa_i !: Sea-Ice Age times ice area (s) 399 404 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: bv_i !: brine volume 400 405 … … 479 484 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_vice !: ice volume variation [m/s] 480 485 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_vsnw !: snw volume variation [m/s] 486 487 ! 488 !!-------------------------------------------------------------------------- 489 !! * SIMIP extra diagnostics 490 !!-------------------------------------------------------------------------- 491 ! Extra sea ice diagnostics to address the data request 492 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: t_si !: Temperature at Snow-ice interface (K) 493 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: tm_si !: mean temperature at the snow-ice interface (K) 494 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_dmi_dyn !: Change in ice mass due to ice dynamics (kg/m2/s) 495 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_dms_dyn !: Change in snow mass due to ice dynamics (kg/m2/s) 496 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_xmtrp_ice !: X-component of ice mass transport (kg/s) 497 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_ymtrp_ice !: Y-component of ice mass transport (kg/s) 498 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_xmtrp_snw !: X-component of snow mass transport (kg/s) 499 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_ymtrp_snw !: Y-component of snow mass transport (kg/s) 500 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_xatrp !: X-component of area transport (m2/s) 501 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_yatrp !: Y-component of area transport (m2/s) 502 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_fc_bo !: Bottom conduction flux (W/m2) 503 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_fc_su !: Surface conduction flux (W/m2) 504 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_utau_oi !: X-direction ocean-ice stress 505 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_vtau_oi !: Y-direction ocean-ice stress 506 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_dssh_dx !: X-direction sea-surface tilt term (N/m2) 507 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_dssh_dy !: X-direction sea-surface tilt term (N/m2) 508 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_corstrx !: X-direction coriolis stress (N/m2) 509 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_corstry !: Y-direction coriolis stress (N/m2) 510 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_intstrx !: X-direction internal stress (N/m2) 511 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_intstry !: Y-direction internal stress (N/m2) 512 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_sig1 !: Average normal stress in sea ice 513 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_sig2 !: Maximum shear stress in sea ice 514 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: diag_shear !: Maximum shear of sea-ice velocity field 515 481 516 ! 482 517 !!---------------------------------------------------------------------- … … 493 528 INTEGER :: ice_alloc 494 529 ! 495 INTEGER :: ierr(1 5), ii530 INTEGER :: ierr(18), ii 496 531 !!----------------------------------------------------------------- 497 532 … … 507 542 508 543 ii = ii + 1 509 ALLOCATE( t_bo (jpi,jpj) , frld (jpi,jpj) , pfrld (jpi,jpj) , phicif (jpi,jpj) , & 510 & wfx_snw(jpi,jpj) , wfx_ice(jpi,jpj) , wfx_sub(jpi,jpj) , wfx_lam(jpi,jpj) , & 544 ALLOCATE( t_bo (jpi,jpj) , frld (jpi,jpj) , pfrld (jpi,jpj) , phicif (jpi,jpj) , & 545 & wfx_snw(jpi,jpj) , wfx_snw_dyn(jpi,jpj) , wfx_snw_sum(jpi,jpj) , wfx_snw_sub(jpi,jpj) , & 546 & wfx_ice(jpi,jpj) , wfx_sub (jpi,jpj) , wfx_ice_sub(jpi,jpj) , wfx_lam (jpi,jpj) , & 511 547 ! MV MP 2016 512 548 & wfx_pnd(jpi,jpj) , wfx_snw_sum(jpi,jpj) , & … … 568 604 569 605 ! MV MP 2016 606 ii = ii + 1 570 607 ALLOCATE( sxap(jpi,jpj,jpl) , syap(jpi,jpj,jpl) , sxxap(jpi,jpj,jpl) , syyap(jpi,jpj,jpl) , sxyap(jpi,jpj,jpl) , & 571 608 & sxvp(jpi,jpj,jpl) , syvp(jpi,jpj,jpl) , sxxvp(jpi,jpj,jpl) , syyvp(jpi,jpj,jpl) , sxyvp(jpi,jpj,jpl) , & … … 590 627 & diag_trp_es(jpi,jpj) , diag_trp_smv(jpi,jpj) , diag_heat (jpi,jpj), & 591 628 & diag_smvi (jpi,jpj) , diag_vice (jpi,jpj) , diag_vsnw (jpi,jpj), STAT=ierr(ii) ) 629 630 ! * SIMIP diagnostics 631 ii = ii + 1 632 ALLOCATE( t_si (jpi,jpj,jpl) , tm_si(jpi,jpj) , & 633 diag_dmi_dyn(jpi,jpj) , diag_dms_dyn(jpi,jpj) , & 634 diag_xmtrp_ice(jpi,jpj), diag_ymtrp_ice(jpi,jpj), & 635 diag_xmtrp_snw(jpi,jpj), diag_ymtrp_snw(jpi,jpj), & 636 diag_xatrp(jpi,jpj) , diag_yatrp(jpi,jpj) , & 637 diag_fc_bo(jpi,jpj) , diag_fc_su(jpi,jpj) , & 638 diag_utau_oi(jpi,jpj) , diag_vtau_oi(jpi,jpj) , & 639 diag_dssh_dx(jpi,jpj) , diag_dssh_dy(jpi,jpj) , & 640 diag_corstrx(jpi,jpj) , diag_corstry(jpi,jpj) , & 641 diag_intstrx(jpi,jpj) , diag_intstry(jpi,jpj) , & 642 diag_sig1(jpi,jpj) , diag_sig2(jpi,jpj) , & 643 STAT = ierr(ii) ) 592 644 593 645 ice_alloc = MAXVAL( ierr(:) ) -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90
r8233 r8239 302 302 ht_i(ji,jj,jl) = zswitch(ji,jj) * zh_i_ini(ji,jj,jl) ! ice thickness 303 303 sm_i(ji,jj,jl) = zswitch(ji,jj) * zsm_i_ini(ji,jj) ! salinity 304 o_i(ji,jj,jl) = zswitch(ji,jj) * 1._wp ! age (1day)304 o_i(ji,jj,jl) = 0._wp ! age (0 day) 305 305 t_su(ji,jj,jl) = zswitch(ji,jj) * zts_u_ini(ji,jj) + ( 1._wp - zswitch(ji,jj) ) * rt0 ! surf temp 306 306 -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90
r8233 r8239 500 500 INTEGER :: icells ! number of cells with a_i > puny 501 501 REAL(wp) :: hL, hR, farea ! left and right limits of integration 502 REAL(wp) :: zwfx_snw ! snow mass flux increment 502 503 503 504 INTEGER , POINTER, DIMENSION(:) :: indxi, indxj ! compressed indices … … 654 655 & - sm_i(ji,jj,jl1) * vsw(ij) * rhoic * r1_rdtice ! and get sm_i from the ocean 655 656 ENDIF 656 657 657 658 !------------------------------------------ 658 659 ! 3.7 Put the snow somewhere in the ocean … … 663 664 ! During the next time step, thermo_rates will determine whether 664 665 ! the ocean cools or new ice grows. 665 wfx_snw(ji,jj) = wfx_snw(ji,jj) + ( rhosn * vsrdg(ij) * ( 1._wp - rn_fsnowrdg ) & 666 & + rhosn * vsrft(ij) * ( 1._wp - rn_fsnowrft ) ) * r1_rdtice ! fresh water source for ocean 666 zwfx_snw = ( rhosn * vsrdg(ij) * ( 1._wp - rn_fsnowrdg ) & 667 & + rhosn * vsrft(ij) * ( 1._wp - rn_fsnowrft ) ) * r1_rdtice ! fresh water source for ocean 668 669 wfx_snw_dyn(ji,jj) = wfx_snw_dyn(ji,jj) + zwfx_snw 670 wfx_snw(ji,jj) = wfx_snw(ji,jj) + zwfx_snw 667 671 668 672 hfx_dyn(ji,jj) = hfx_dyn(ji,jj) + ( - esrdg(ij) * ( 1._wp - rn_fsnowrdg ) & … … 789 793 END DO ! jl1 (deforming categories) 790 794 795 ! SIMIP diagnostics 796 diag_dmi_dyn(:,:) = - wfx_dyn(:,:) + rhoic * diag_trp_vi(:,:) 797 diag_dms_dyn(:,:) = - wfx_snw_dyn(:,:) + rhosn * diag_trp_vs(:,:) 798 791 799 ! 792 800 CALL wrk_dealloc( jpij, indxi, indxj ) -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90
r7753 r8239 117 117 REAL(wp) :: zm1, zm2, zm3, zmassU, zmassV ! ice/snow mass 118 118 REAL(wp) :: zdelta, zp_delf, zds2, zdt, zdt2, zdiv, zdiv2 ! temporary scalars 119 REAL(wp) :: zTauO, zTauB, zTauE, z Cor, zvel! temporary scalars119 REAL(wp) :: zTauO, zTauB, zTauE, zvel ! temporary scalars 120 120 121 121 REAL(wp) :: zsig1, zsig2 ! internal ice stress 122 122 REAL(wp) :: zresm ! Maximal error on ice velocity 123 123 REAL(wp) :: zintb, zintn ! dummy argument 124 REAL(wp) :: zfac_x, zfac_y 124 125 125 126 REAL(wp), POINTER, DIMENSION(:,:) :: z1_e1t0, z1_e2t0 ! scale factors … … 140 141 ! ocean surface (ssh_m) if ice is not embedded 141 142 ! ice top surface if ice is embedded 143 REAL(wp), POINTER, DIMENSION(:,:) :: zCorx, zCory ! Coriolis stress array 144 REAL(wp), POINTER, DIMENSION(:,:) :: ztaux_oi, ztauy_oi ! Ocean-to-ice stress array 145 142 146 REAL(wp), POINTER, DIMENSION(:,:) :: zswitchU, zswitchV ! dummy arrays 143 147 REAL(wp), POINTER, DIMENSION(:,:) :: zmaskU, zmaskV ! mask for ice presence … … 154 158 CALL wrk_alloc( jpi,jpj, zds, zs1, zs2, zs12, zu_ice, zv_ice, zresr, zpice ) 155 159 CALL wrk_alloc( jpi,jpj, zswitchU, zswitchV, zmaskU, zmaskV, zfmask, zwf ) 160 CALL wrk_alloc( jpi,jpj, zCorx, zCory) 161 CALL wrk_alloc( jpi,jpj, ztaux_oi, ztauy_oi) 156 162 157 163 #if defined key_agrif … … 427 433 & + ( u_iceV(ji,jj) - u_oceV(ji,jj) ) * ( u_iceV(ji,jj) - u_oceV(ji,jj) ) ) 428 434 435 ! Ocean-to-Ice stress 436 ztauy_oi(ji,jj) = zTauO * ( v_oce(ji,jj) - v_ice(ji,jj) ) 437 429 438 ! tau_bottom/v_ice 430 439 zvel = MAX( zepsi, SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) ) … … 432 441 433 442 ! Coriolis at V-points (energy conserving formulation) 434 zCor = - 0.25_wp * r1_e2v(ji,jj) * &443 zCory(ji,jj) = - 0.25_wp * r1_e2v(ji,jj) * & 435 444 & ( zmf(ji,jj ) * ( e2u(ji,jj ) * u_ice(ji,jj ) + e2u(ji-1,jj ) * u_ice(ji-1,jj ) ) & 436 445 & + zmf(ji,jj+1) * ( e2u(ji,jj+1) * u_ice(ji,jj+1) + e2u(ji-1,jj+1) * u_ice(ji-1,jj+1) ) ) 437 446 438 447 ! Sum of external forces (explicit solution) = F + tau_ia + Coriolis + spg + tau_io 439 zTauE = zfV(ji,jj) + zTauV_ia(ji,jj) + zCor + zspgV(ji,jj) + zTauO * ( v_oce(ji,jj) - v_ice(ji,jj))448 zTauE = zfV(ji,jj) + zTauV_ia(ji,jj) + zCory(ji,jj) + zspgV(ji,jj) + ztauy_oi(ji,jj) 440 449 441 450 ! landfast switch => 0 = static friction ; 1 = sliding friction … … 451 460 ! Bouillon 2013 452 461 !!v_ice(ji,jj) = ( zmV_t(ji,jj) * ( zbeta * v_ice(ji,jj) + v_ice_b(ji,jj) ) & 453 !! & + zfV(ji,jj) + zCor + zTauV_ia(ji,jj) + zTauO * v_oce(ji,jj) + zspgV(ji,jj) &462 !! & + zfV(ji,jj) + zCory(ji,jj) + zTauV_ia(ji,jj) + zTauO * v_oce(ji,jj) + zspgV(ji,jj) & 454 463 !! & ) / MAX( zmV_t(ji,jj) * ( zbeta + 1._wp ) + zTauO - zTauB, zepsi ) * zswitchV(ji,jj) 455 464 … … 471 480 & + ( v_iceU(ji,jj) - v_oceU(ji,jj) ) * ( v_iceU(ji,jj) - v_oceU(ji,jj) ) ) 472 481 482 ! Ocean-to-Ice stress 483 ztaux_oi(ji,jj) = zTauO * ( u_oce(ji,jj) - u_ice(ji,jj) ) 484 473 485 ! tau_bottom/u_ice 474 486 zvel = MAX( zepsi, SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) ) … … 476 488 477 489 ! Coriolis at U-points (energy conserving formulation) 478 zCor = 0.25_wp * r1_e1u(ji,jj) * &490 zCorx(ji,jj) = 0.25_wp * r1_e1u(ji,jj) * & 479 491 & ( zmf(ji ,jj) * ( e1v(ji ,jj) * v_ice(ji ,jj) + e1v(ji ,jj-1) * v_ice(ji ,jj-1) ) & 480 492 & + zmf(ji+1,jj) * ( e1v(ji+1,jj) * v_ice(ji+1,jj) + e1v(ji+1,jj-1) * v_ice(ji+1,jj-1) ) ) 481 493 482 494 ! Sum of external forces (explicit solution) = F + tau_ia + Coriolis + spg + tau_io 483 zTauE = zfU(ji,jj) + zTauU_ia(ji,jj) + zCor + zspgU(ji,jj) + zTauO * ( u_oce(ji,jj) - u_ice(ji,jj))495 zTauE = zfU(ji,jj) + zTauU_ia(ji,jj) + zCorx(ji,jj) + zspgU(ji,jj) + ztaux_oi(ji,jj) 484 496 485 497 ! landfast switch => 0 = static friction ; 1 = sliding friction … … 495 507 ! Bouillon 2013 496 508 !!u_ice(ji,jj) = ( zmU_t(ji,jj) * ( zbeta * u_ice(ji,jj) + u_ice_b(ji,jj) ) & 497 !! & + zfU(ji,jj) + zCor + zTauU_ia(ji,jj) + zTauO * u_oce(ji,jj) + zspgU(ji,jj) &509 !! & + zfU(ji,jj) + zCorx(ji,jj) + zTauU_ia(ji,jj) + zTauO * u_oce(ji,jj) + zspgU(ji,jj) & 498 510 !! & ) / MAX( zmU_t(ji,jj) * ( zbeta + 1._wp ) + zTauO - zTauB, zepsi ) * zswitchU(ji,jj) 499 511 END DO … … 516 528 & + ( v_iceU(ji,jj) - v_oceU(ji,jj) ) * ( v_iceU(ji,jj) - v_oceU(ji,jj) ) ) 517 529 530 ! Ocean-to-Ice stress 531 ztaux_oi(ji,jj) = zTauO * ( u_oce(ji,jj) - u_ice(ji,jj) ) 532 518 533 ! tau_bottom/u_ice 519 534 zvel = MAX( zepsi, SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) ) … … 521 536 522 537 ! Coriolis at U-points (energy conserving formulation) 523 zCor = 0.25_wp * r1_e1u(ji,jj) * &538 zCorx(ji,jj) = 0.25_wp * r1_e1u(ji,jj) * & 524 539 & ( zmf(ji ,jj) * ( e1v(ji ,jj) * v_ice(ji ,jj) + e1v(ji ,jj-1) * v_ice(ji ,jj-1) ) & 525 540 & + zmf(ji+1,jj) * ( e1v(ji+1,jj) * v_ice(ji+1,jj) + e1v(ji+1,jj-1) * v_ice(ji+1,jj-1) ) ) 526 541 527 542 ! Sum of external forces (explicit solution) = F + tau_ia + Coriolis + spg + tau_io 528 zTauE = zfU(ji,jj) + zTauU_ia(ji,jj) + zCor + zspgU(ji,jj) + zTauO * ( u_oce(ji,jj) - u_ice(ji,jj))543 zTauE = zfU(ji,jj) + zTauU_ia(ji,jj) + zCorx(ji,jj) + zspgU(ji,jj) + ztaux_oi(ji,jj) 529 544 530 545 ! landfast switch => 0 = static friction ; 1 = sliding friction … … 540 555 ! Bouillon 2013 541 556 !!u_ice(ji,jj) = ( zmU_t(ji,jj) * ( zbeta * u_ice(ji,jj) + u_ice_b(ji,jj) ) & 542 !! & + zfU(ji,jj) + zCor + zTauU_ia(ji,jj) + zTauO * u_oce(ji,jj) + zspgU(ji,jj) &557 !! & + zfU(ji,jj) + zCorx(ji,jj) + zTauU_ia(ji,jj) + zTauO * u_oce(ji,jj) + zspgU(ji,jj) & 543 558 !! & ) / MAX( zmU_t(ji,jj) * ( zbeta + 1._wp ) + zTauO - zTauB, zepsi ) * zswitchU(ji,jj) 544 559 END DO … … 559 574 & + ( u_iceV(ji,jj) - u_oceV(ji,jj) ) * ( u_iceV(ji,jj) - u_oceV(ji,jj) ) ) 560 575 576 ! Ocean-to-Ice stress 577 ztauy_oi(ji,jj) = zTauO * ( v_oce(ji,jj) - v_ice(ji,jj) ) 578 561 579 ! tau_bottom/v_ice 562 580 zvel = MAX( zepsi, SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) ) … … 564 582 565 583 ! Coriolis at V-points (energy conserving formulation) 566 zCor = - 0.25_wp * r1_e2v(ji,jj) * &584 zCory(ji,jj) = - 0.25_wp * r1_e2v(ji,jj) * & 567 585 & ( zmf(ji,jj ) * ( e2u(ji,jj ) * u_ice(ji,jj ) + e2u(ji-1,jj ) * u_ice(ji-1,jj ) ) & 568 586 & + zmf(ji,jj+1) * ( e2u(ji,jj+1) * u_ice(ji,jj+1) + e2u(ji-1,jj+1) * u_ice(ji-1,jj+1) ) ) 569 587 570 588 ! Sum of external forces (explicit solution) = F + tau_ia + Coriolis + spg + tau_io 571 zTauE = zfV(ji,jj) + zTauV_ia(ji,jj) + zCor + zspgV(ji,jj) + zTauO * ( v_oce(ji,jj) - v_ice(ji,jj))589 zTauE = zfV(ji,jj) + zTauV_ia(ji,jj) + zCory(ji,jj) + zspgV(ji,jj) + ztauy_oi(ji,jj) 572 590 573 591 ! landfast switch => 0 = static friction (tau_icebfr > zTauE); 1 = sliding friction … … 583 601 ! Bouillon 2013 584 602 !!v_ice(ji,jj) = ( zmV_t(ji,jj) * ( zbeta * v_ice(ji,jj) + v_ice_b(ji,jj) ) & 585 !! & + zfV(ji,jj) + zCor + zTauV_ia(ji,jj) + zTauO * v_oce(ji,jj) + zspgV(ji,jj) &603 !! & + zfV(ji,jj) + zCory(ji,jj) + zTauV_ia(ji,jj) + zTauO * v_oce(ji,jj) + zspgV(ji,jj) & 586 604 !! & ) / MAX( zmV_t(ji,jj) * ( zbeta + 1._wp ) + zTauO - zTauB, zepsi ) * zswitchV(ji,jj) 587 605 END DO … … 662 680 663 681 !------------------------------------------------------------------------------! 664 ! 5) Control prints of residual and charge ellipse 682 ! 5) SIMIP diagnostics 683 !------------------------------------------------------------------------------! 684 685 DO jj = 2, jpjm1 686 DO ji = 2, jpim1 687 rswitch = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) ! 1 if ice, 0 if no ice 688 689 ! Stress tensor invariants (normal and shear stress N/m) 690 diag_sig1(ji,jj) = ( zs1(ji,jj) + zs2(ji,jj) ) * rswitch ! normal stress 691 diag_sig2(ji,jj) = SQRT( ( zs1(ji,jj) - zs2(ji,jj) )**2 + 4*zs12(ji,jj)**2 ) * rswitch ! shear stress 692 693 ! Stress terms of the momentum equation (N/m2) 694 diag_dssh_dx(ji,jj) = zspgU(ji,jj) * rswitch ! sea surface slope stress term 695 diag_dssh_dy(ji,jj) = zspgV(ji,jj) * rswitch 696 697 diag_corstrx(ji,jj) = zCorx(ji,jj) * rswitch ! Coriolis stress term 698 diag_corstry(ji,jj) = zCory(ji,jj) * rswitch 699 700 diag_intstrx(ji,jj) = zfU(ji,jj) * rswitch ! internal stress term 701 diag_intstry(ji,jj) = zfV(ji,jj) * rswitch 702 703 diag_utau_oi(ji,jj) = ztaux_oi(ji,jj) * rswitch ! oceanic stress 704 diag_vtau_oi(ji,jj) = ztauy_oi(ji,jj) * rswitch 705 706 ! 2D ice mass, snow mass, area transport arrays (X, Y) 707 zfac_x = 0.5 * u_ice(ji,jj) * e2u(ji,jj) * rswitch 708 zfac_y = 0.5 * v_ice(ji,jj) * e1v(ji,jj) * rswitch 709 710 diag_xmtrp_ice(ji,jj) = rhoic * zfac_x * ( vt_i(ji+1,jj) + vt_i(ji,jj) ) ! ice mass transport, X-component 711 diag_ymtrp_ice(ji,jj) = rhoic * zfac_y * ( vt_i(ji,jj+1) + vt_i(ji,jj) ) ! '' Y- '' 712 713 diag_xmtrp_snw(ji,jj) = rhosn * zfac_x * ( vt_s(ji+1,jj) + vt_s(ji,jj) ) ! snow mass transport, X-component 714 diag_ymtrp_snw(ji,jj) = rhosn * zfac_y * ( vt_s(ji,jj+1) + vt_s(ji,jj) ) ! '' Y- '' 715 716 diag_xatrp(ji,jj) = zfac_x * ( at_i(ji+1,jj) + at_i(ji,jj) ) ! area transport, X-component 717 diag_yatrp(ji,jj) = zfac_y * ( at_i(ji,jj+1) + at_i(ji,jj) ) ! '' Y- '' 718 719 END DO 720 END DO 721 722 CALL lbc_lnk_multi( diag_sig1 , 'T', 1., diag_sig2 , 'T', 1., & 723 & diag_dssh_dx, 'U', -1., diag_dssh_dy, 'V', -1., & 724 & diag_corstrx, 'U', -1., diag_corstry, 'V', -1., & 725 & diag_intstrx, 'U', -1., diag_intstry, 'V', -1. ) 726 727 CALL lbc_lnk_multi( diag_utau_oi, 'U', -1., diag_vtau_oi, 'V', -1. ) 728 729 ! 730 !------------------------------------------------------------------------------! 731 ! 6) Control prints of residual and charge ellipse 665 732 !------------------------------------------------------------------------------! 666 733 ! … … 703 770 CALL wrk_dealloc( jpi,jpj, zds, zs1, zs2, zs12, zu_ice, zv_ice, zresr, zpice ) 704 771 CALL wrk_dealloc( jpi,jpj, zswitchU, zswitchV, zmaskU, zmaskV, zfmask, zwf ) 772 CALL wrk_dealloc( jpi,jpj, zCorx, zCory ) 773 CALL wrk_dealloc( jpi,jpj, ztaux_oi, ztauy_oi ) 705 774 706 775 END SUBROUTINE lim_rhg -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90
r8233 r8239 83 83 INTEGER :: ji, jj, jk, jl ! dummy loop indices 84 84 INTEGER :: nbpb ! nb of icy pts for vertical thermo calculations 85 REAL(wp) :: zfric_u, zqld, zqfr 85 REAL(wp) :: zfric_u, zqld, zqfr, zqfr_neg 86 86 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 87 87 REAL(wp), PARAMETER :: zfric_umin = 0._wp ! lower bound for the friction velocity (cice value=5.e-04) … … 175 175 176 176 ! --- Energy needed to bring ocean surface layer until its freezing (<0, J.m-2) --- ! 177 ! includes supercooling potential energy (>0) or "above-freezing" energy (<0) 177 178 zqfr = tmask(ji,jj,1) * rau0 * rcp * e3t_m(ji,jj) * ( t_bo(ji,jj) - ( sst_m(ji,jj) + rt0 ) ) 178 179 179 ! --- Energy from the turbulent oceanic heat flux (W/m2) --- ! 180 ! --- Above-freezing sensible heat content (J/m2 grid) 181 zqfr_neg = tmask(ji,jj,1) * rau0 * rcp * e3t_m(ji,jj) * MIN( ( t_bo(ji,jj) - ( sst_m(ji,jj) + rt0 ) ), 0._wp ) 182 183 ! --- Sensible ocean-to-ice heat flux (W/m2) 180 184 zfric_u = MAX( SQRT( zfric(ji,jj) ), zfric_umin ) 181 fhtur(ji,jj) = MAX( 0._wp, rswitch * rau0 * rcp * zch * zfric_u * ( ( sst_m(ji,jj) + rt0 ) - t_bo(ji,jj) ) ) ! W.m-2 182 fhtur(ji,jj) = rswitch * MIN( fhtur(ji,jj), - zqfr * r1_rdtice / MAX( at_i(ji,jj), epsi10 ) ) 185 fhtur(ji,jj) = rswitch * rau0 * rcp * zch * zfric_u * ( ( sst_m(ji,jj) + rt0 ) - t_bo(ji,jj) ) ! W.m-2 186 187 fhtur(ji,jj) = rswitch * MIN( fhtur(ji,jj), - zqfr_neg * r1_rdtice / MAX( at_i(ji,jj), epsi10 ) ) 183 188 ! upper bound for fhtur: the heat retrieved from the ocean must be smaller than the heat necessary to reach 184 189 ! the freezing point, so that we do not have SST < T_freeze 185 190 ! This implies: - ( fhtur(ji,jj) * at_i(ji,jj) * rtdice ) - zqfr >= 0 186 191 187 !-- Energy Budget of the leads (J.m-2) . Must be < 0 to form ice192 !-- Energy Budget of the leads (J.m-2), source of lateral accretion. Must be < 0 to form ice 188 193 qlead(ji,jj) = MIN( 0._wp , zqld - ( fhtur(ji,jj) * at_i(ji,jj) * rdt_ice ) - zqfr ) 189 194 … … 470 475 CALL tab_2d_1d( nbpb, wfx_snw_sum_1d(1:nbpb), wfx_snw_sum , jpi, jpj, npb(1:nbpb) ) 471 476 CALL tab_2d_1d( nbpb, wfx_sub_1d (1:nbpb), wfx_sub , jpi, jpj, npb(1:nbpb) ) 477 CALL tab_2d_1d( nbpb, wfx_snw_sub_1d(1:nbpb), wfx_snw_sub , jpi, jpj, npb(1:nbpb) ) 478 CALL tab_2d_1d( nbpb, wfx_ice_sub_1d(1:nbpb), wfx_ice_sub , jpi, jpj, npb(1:nbpb) ) 472 479 ! 473 480 CALL tab_2d_1d( nbpb, wfx_bog_1d (1:nbpb), wfx_bog , jpi, jpj, npb(1:nbpb) ) … … 500 507 CALL tab_2d_1d( nbpb, hfx_err_rem_1d (1:nbpb), hfx_err_rem , jpi, jpj, npb(1:nbpb) ) 501 508 ! 509 ! SIMIP diagnostics 510 CALL tab_2d_1d( nbpb, diag_fc_bo_1d (1:nbpb), diag_fc_bo , jpi, jpj, npb(1:nbpb) ) 511 CALL tab_2d_1d( nbpb, diag_fc_su_1d (1:nbpb), diag_fc_su , jpi, jpj, npb(1:nbpb) ) 512 ! 502 513 CASE( 2 ) ! from 1D to 2D 503 514 ! … … 522 533 CALL tab_1d_2d( nbpb, wfx_snw_sum , npb, wfx_snw_sum_1d(1:nbpb),jpi, jpj ) 523 534 CALL tab_1d_2d( nbpb, wfx_sub , npb, wfx_sub_1d(1:nbpb) , jpi, jpj ) 535 CALL tab_1d_2d( nbpb, wfx_snw_sub , npb, wfx_snw_sub_1d(1:nbpb), jpi, jpj ) 536 CALL tab_1d_2d( nbpb, wfx_ice_sub , npb, wfx_ice_sub_1d(1:nbpb), jpi, jpj ) 524 537 ! 525 538 CALL tab_1d_2d( nbpb, wfx_bog , npb, wfx_bog_1d(1:nbpb) , jpi, jpj ) … … 554 567 CALL tab_1d_2d( nbpb, qns_ice(:,:,jl), npb, qns_ice_1d(1:nbpb) , jpi, jpj) 555 568 CALL tab_1d_2d( nbpb, ftr_ice(:,:,jl), npb, ftr_ice_1d(1:nbpb) , jpi, jpj ) 556 ! 569 ! 570 ! SIMIP diagnostics 571 CALL tab_1d_2d( nbpb, t_si(:,:,jl) , npb, t_si_1d (1:nbpb) , jpi, jpj ) 572 CALL tab_1d_2d( nbpb, diag_fc_bo , npb, diag_fc_bo_1d(1:nbpb) , jpi, jpj ) 573 CALL tab_1d_2d( nbpb, diag_fc_su , npb, diag_fc_su_1d(1:nbpb) , jpi, jpj ) 557 574 END SELECT 558 575 ! -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90
r8233 r8239 87 87 REAL(wp) :: zdE ! specific enthalpy difference (J/kg) 88 88 REAL(wp) :: zfmdt ! exchange mass flux x time step (J/m2), >0 towards the ocean 89 REAL(wp) :: zsstK ! SST in Kelvin89 REAL(wp) :: zsstK ! SST (K) 90 90 91 91 REAL(wp), POINTER, DIMENSION(:) :: zqprec ! energy of fallen snow (J.m-3) … … 277 277 & ) * a_i_1d(ji) * r1_rdtice 278 278 ! Mass flux by sublimation 279 wfx_sub_1d(ji) = wfx_sub_1d(ji) - rhosn * a_i_1d(ji) * zdh_s_sub(ji) * r1_rdtice 279 wfx_snw_sub_1d(ji) = wfx_snw_sub_1d(ji) - rhosn * a_i_1d(ji) * zdh_s_sub(ji) * r1_rdtice 280 280 281 ! new snow thickness 281 282 ht_s_1d(ji) = MAX( 0._wp , ht_s_1d(ji) + zdh_s_sub(ji) ) … … 378 379 hfx_sub_1d(ji) = hfx_sub_1d(ji) + zdum * q_i_1d(ji,jk) * a_i_1d(ji) * r1_rdtice 379 380 ! Mass flux > 0 380 wfx_sub_1d(ji) = wfx_sub_1d(ji) - rhoic * a_i_1d(ji) * zdum * r1_rdtice 381 wfx_ice_sub_1d(ji) = wfx_ice_sub_1d(ji) - rhoic * a_i_1d(ji) * zdum * r1_rdtice 382 381 383 ! update remaining mass flux 382 384 zevap_rema(ji) = zevap_rema(ji) + zdum * rhoic … … 615 617 IF( ln_limctl .AND. zq_rema(ji) < 0. .AND. lwp ) WRITE(numout,*) 'ALERTE zq_rema <0 = ', zq_rema(ji) 616 618 END DO 617 619 620 ! Water fluxes 621 DO ji = kideb, kiut 622 wfx_sub_1d(ji) = wfx_snw_sub_1d(ji) + wfx_ice_sub_1d(ji) ! sum ice and snow sublimation contributions 623 END DO 624 618 625 ! 619 626 !------------------------------------------------------------------------------| -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limthd_dif.F90
r7813 r8239 734 734 END DO ! End of the do while iterative procedure 735 735 736 ! MV SIMIP 2016 737 !--- Snow-ice interfacial temperature (diagnostic SIMIP) 738 DO ji = kideb, kiut 739 zfac = 1. / MAX( epsi10 , rn_cdsn * zh_i(ji) + ztcond_i(ji,1) * zh_s(ji) ) 740 t_si_1d(ji) = ( rn_cdsn * zh_i(ji) * t_s_1d(ji,1) + & 741 & ztcond_i(ji,1) * zh_s(ji) * t_i_1d(ji,1) ) * zfac 742 END DO 743 WHERE( ( zh_s < 1.0e-3 ) ) ; t_si_1d(:) = t_su_1d(:) ; END WHERE 744 ! END MV SIMIP 2016 745 736 746 IF( ln_limctl .AND. lwp ) THEN 737 747 WRITE(numout,*) ' zerritmax : ', zerritmax … … 751 761 fc_bo_i(ji) = - zkappa_i(ji,nlay_i) * ( zg1*(t_bo_1d(ji) - t_i_1d(ji,nlay_i)) ) 752 762 END DO 763 764 ! MV SIMIP 2016 765 !--- Conduction fluxes (positive downwards) 766 diag_fc_bo_1d(:) = diag_fc_bo_1d(:) + fc_bo_i(:) * a_i_1d(:) / at_i_1d(:) 767 diag_fc_su_1d(:) = diag_fc_su_1d(:) + fc_su(:) * a_i_1d(:) / at_i_1d(:) 768 ! END MV SIMIP 2016 753 769 754 770 ! --- computes sea ice energy of melting compulsory for limthd_dh --- ! -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limupdate1.F90
r7753 r8239 80 80 IF( at_i(ji,jj) > rn_amax_2d(ji,jj) .AND. a_i(ji,jj,jl) > 0._wp ) THEN 81 81 a_i (ji,jj,jl) = a_i (ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax_2d(ji,jj) / at_i(ji,jj) ) ) 82 oa_i(ji,jj,jl) = oa_i(ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax_2d(ji,jj) / at_i(ji,jj) ) )83 82 ENDIF 84 83 END DO -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limupdate2.F90
r7753 r8239 77 77 IF( v_i(ji,jj,1) > 0._wp .AND. ht_i(ji,jj,1) < rn_himin ) THEN 78 78 a_i (ji,jj,1) = a_i (ji,jj,1) * ht_i(ji,jj,1) / rn_himin 79 oa_i(ji,jj,1) = oa_i(ji,jj,1) * ht_i(ji,jj,1) / rn_himin80 79 ENDIF 81 80 END DO … … 95 94 IF( at_i(ji,jj) > rn_amax_2d(ji,jj) .AND. a_i(ji,jj,jl) > 0._wp ) THEN 96 95 a_i (ji,jj,jl) = a_i (ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax_2d(ji,jj) / at_i(ji,jj) ) ) 97 oa_i(ji,jj,jl) = oa_i(ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax_2d(ji,jj) / at_i(ji,jj) ) )98 96 ENDIF 99 97 END DO … … 155 153 ! ------------------------------------------------- 156 154 DO jl = 1, jpl 157 oa_i(:,:,jl) = oa_i(:,:,jl) + a_i(:,:,jl) * rdt_ice / rday! ice natural aging155 oa_i(:,:,jl) = oa_i(:,:,jl) + a_i(:,:,jl) * rdt_ice ! ice natural aging 158 156 afx_thd(:,:) = afx_thd(:,:) + ( a_i(:,:,jl) - a_i_b(:,:,jl) ) * r1_rdtice 159 157 END DO -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90
r8233 r8239 120 120 tm_i(:,:) = 0._wp 121 121 tm_su(:,:) = 0._wp 122 tm_si(:,:) = 0._wp 122 123 om_i (:,:) = 0._wp 123 124 DO jl = 1, jpl … … 127 128 rswitch = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi10 ) ) 128 129 tm_su(ji,jj) = tm_su(ji,jj) + rswitch * ( t_su(ji,jj,jl) - rt0 ) * a_i(ji,jj,jl) / MAX( at_i(ji,jj) , epsi10 ) 129 om_i (ji,jj) = om_i (ji,jj) + rswitch * oa_i(ji,jj,jl) / MAX( at_i(ji,jj) , epsi10 ) 130 tm_si(ji,jj) = tm_si(ji,jj) + rswitch * ( t_si(ji,jj,jl) - rt0 ) * a_i(ji,jj,jl) / MAX( at_i(ji,jj) , epsi10 ) 131 om_i (ji,jj) = om_i (ji,jj) + rswitch * oa_i(ji,jj,jl) / MAX( at_i(ji,jj) , epsi06 ) 130 132 END DO 131 133 END DO … … 145 147 tm_i = tm_i + rt0 146 148 tm_su = tm_su + rt0 149 tm_si = tm_si + rt0 147 150 ! 148 151 ENDIF -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90
r8233 r8239 54 54 ! 55 55 INTEGER :: ji, jj, jk, jl ! dummy loop indices 56 REAL(wp) :: z1_365 57 REAL(wp) :: z2da, z2db, ztmp 56 REAL(wp) :: z2da, z2db, ztmp, zrho1, zrho2 58 57 REAL(wp), POINTER, DIMENSION(:,:,:) :: zswi2 59 58 REAL(wp), POINTER, DIMENSION(:,:) :: z2d, zswi ! 2D workspace 59 REAL(wp), POINTER, DIMENSION(:,:) :: zfb ! ice freeboard 60 REAL(wp), POINTER, DIMENSION(:,:) :: zamask, zamask15 ! 15% concentration mask 61 62 ! Global ice diagnostics (SIMIP) 63 REAL(wp) :: zdiag_area_nh, & ! area, extent, volume 64 & zdiag_extt_nh, & 65 & zdiag_area_sh, & 66 & zdiag_extt_sh, & 67 & zdiag_volu_nh, & 68 & zdiag_volu_sh 69 60 70 !!------------------------------------------------------------------- 61 71 … … 64 74 CALL wrk_alloc( jpi,jpj,jpl, zswi2 ) 65 75 CALL wrk_alloc( jpi,jpj , z2d, zswi ) 66 76 CALL wrk_alloc( jpi,jpj , zfb, zamask, zamask15 ) 67 77 !----------------------------- 68 78 ! Mean category values 69 79 !----------------------------- 70 z1_365 = 1._wp / 365._wp71 80 72 81 ! brine volume … … 76 85 DO jj = 1, jpj 77 86 DO ji = 1, jpi 78 zswi(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) 87 zswi(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) ! 1 if ice, 0 if no ice 88 zamask(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - 0.05 ) ) ! 1 if 5% ice, 0 if less - required to mask thickness and snow depth 89 zamask15(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - 0.15 ) ) ! 1 if 15% ice, 0 if less 79 90 END DO 80 91 END DO … … 119 130 IF ( iom_use( "tau_icebfr" ) ) CALL iom_put( "tau_icebfr" , tau_icebfr ) ! ice friction with ocean bottom (landfast ice) 120 131 ! 121 IF ( iom_use( "miceage" ) ) CALL iom_put( "miceage" , om_i * zswi * z1_365 ) ! mean ice age 122 IF ( iom_use( "icethic_cea" ) ) CALL iom_put( "icethic_cea" , htm_i * zswi ) ! ice thickness mean 123 IF ( iom_use( "snowthic_cea" ) ) CALL iom_put( "snowthic_cea", htm_s * zswi ) ! snow thickness mean 132 IF ( iom_use( "miceage" ) ) CALL iom_put( "miceage" , om_i * zswi * zamask15 ) ! mean ice age 124 133 IF ( iom_use( "micet" ) ) CALL iom_put( "micet" , ( tm_i - rt0 ) * zswi ) ! ice mean temperature 125 134 IF ( iom_use( "icest" ) ) CALL iom_put( "icest" , ( tm_su - rt0 ) * zswi ) ! ice surface temperature … … 133 142 CALL iom_put( "isnowhc" , et_s * zswi ) ! snow total heat content 134 143 CALL iom_put( "ibrinv" , bvm_i * zswi * 100. ) ! brine volume 135 CALL iom_put( "utau_ice" , utau_ice 136 CALL iom_put( "vtau_ice" , vtau_ice 144 CALL iom_put( "utau_ice" , utau_ice*zswi ) ! wind stress over ice along i-axis at I-point 145 CALL iom_put( "vtau_ice" , vtau_ice*zswi ) ! wind stress over ice along j-axis at I-point 137 146 CALL iom_put( "snowpre" , sprecip * 86400. ) ! snow precipitation 138 147 CALL iom_put( "micesalt" , smt_i * zswi ) ! mean ice salinity 139 148 140 CALL iom_put( "icestr" , strength * zswi ) ! ice strength141 CALL iom_put( "idive" , divu_i * 1.0e8 )! divergence142 CALL iom_put( "ishear" , shear_i * 1.0e8 )! shear149 CALL iom_put( "icestr" , strength * zswi ) ! ice strength 150 CALL iom_put( "idive" , divu_i ) ! divergence 151 CALL iom_put( "ishear" , shear_i ) ! shear 143 152 CALL iom_put( "snowvol" , vt_s * zswi ) ! snow volume 144 153 … … 186 195 CALL iom_put( "vfxsnw" , wfx_snw * ztmp ) ! total snw growth/melt 187 196 CALL iom_put( "vfxsub" , wfx_sub * ztmp ) ! sublimation (snow/ice) 188 CALL iom_put( "vfxsub_err" , wfx_err_sub * ztmp ) ! "excess" of sublimation sent to ocean 189 190 CALL iom_put( "afxtot" , afx_tot * rday) ! concentration tendency (total)191 CALL iom_put( "afxdyn" , afx_dyn * rday) ! concentration tendency (dynamics)192 CALL iom_put( "afxthd" , afx_thd * rday) ! concentration tendency (thermo)197 CALL iom_put( "vfxsub_err" , wfx_err_sub * ztmp ) ! "excess" of sublimation sent to ocean 198 199 CALL iom_put( "afxtot" , afx_tot ) ! concentration tendency (total) 200 CALL iom_put( "afxdyn" , afx_dyn ) ! concentration tendency (dynamics) 201 CALL iom_put( "afxthd" , afx_thd ) ! concentration tendency (thermo) 193 202 194 203 CALL iom_put ('hfxthd' , hfx_thd(:,:) ) ! … … 218 227 ! END MV MP 2016 219 228 220 221 !-------------------------------- 222 ! Output values for each category 223 !-------------------------------- 229 !---------------------------------- 230 ! Output category-dependent fields 231 !---------------------------------- 224 232 IF ( iom_use( "iceconc_cat" ) ) CALL iom_put( "iceconc_cat" , a_i * zswi2 ) ! area for categories 225 233 IF ( iom_use( "icethic_cat" ) ) CALL iom_put( "icethic_cat" , ht_i * zswi2 ) ! thickness for categories … … 231 239 IF ( iom_use( "snwtemp_cat" ) ) CALL iom_put( "snwtemp_cat", ( SUM( t_s(:,:,:,:), dim=3 ) * r1_nlay_s - rt0 ) * zswi2 ) 232 240 ! ice age 233 IF ( iom_use( "iceage_cat" ) ) CALL iom_put( "iceage_cat" , o_i * zswi2 * z1_365 )241 IF ( iom_use( "iceage_cat" ) ) CALL iom_put( "iceage_cat" , o_i * zswi2 ) 234 242 ! brine volume 235 243 IF ( iom_use( "brinevol_cat" ) ) CALL iom_put( "brinevol_cat", bv_i * 100. * zswi2 ) … … 244 252 ! END MV MP 2016 245 253 254 !-------------------------------- 255 ! Add-ons for SIMIP 256 !-------------------------------- 257 zrho1 = ( rau0 - rhoic ) / rau0; zrho2 = rhosn / rau0 258 259 IF ( iom_use( "icethic" ) ) CALL iom_put( "icethic" , htm_i * zamask ) ! Ice thickness 260 IF ( iom_use( "icepres" ) ) CALL iom_put( "icepres" , zswi ) ! Ice presence (1 or 0) 261 IF ( iom_use( "snowthic" ) ) CALL iom_put( "snowthic" , htm_s * zamask ) ! Snow thickness 262 IF ( iom_use( "icemass" ) ) CALL iom_put( "icemass" , rhoic * vt_i(:,:) * zswi ) ! Ice mass per cell area 263 IF ( iom_use( "snomass" ) ) CALL iom_put( "snomass" , rhosn * vt_s(:,:) * zswi ) ! Snow mass per cell area 264 IF ( iom_use( "icesnt" ) ) CALL iom_put( "icesnt" , ( tm_si - rt0 ) * zswi ) ! Snow-ice interface temperature 265 IF ( iom_use( "icebot" ) ) CALL iom_put( "icebot" , ( t_bo - rt0 ) * zswi ) ! Ice bottom temperature 266 IF ( iom_use( "icesmass" ) ) CALL iom_put( "icesmass" , SUM( smv_i, DIM = 3 ) * rhoic * 1.0e-3 * zswi ) ! Mass of salt in sea ice per cell area 267 IF ( iom_use( "icefb" ) ) THEN 268 zfb(:,:) = ( zrho1 * htm_i(:,:) - zrho2 * htm_s(:,:) ) * zswi(:,:) 269 WHERE( zfb < 0._wp ) ; zfb = 0._wp ; END WHERE 270 CALL iom_put( "icefb" , zfb ) ! Ice freeboard 271 ENDIF 272 IF ( iom_use( "dmithd" ) ) CALL iom_put( "dmithd" , - wfx_bog - wfx_bom - wfx_sum & ! Sea-ice mass change from thermodynamics 273 & - wfx_sni - wfx_opw - wfx_res ) 274 IF ( iom_use( "dmidyn" ) ) CALL iom_put( "dmidyn" , diag_dmi_dyn ) ! Sea-ice mass change from dynamics 275 IF ( iom_use( "dmiopw" ) ) CALL iom_put( "dmiopw" , - wfx_opw ) ! Sea-ice mass change through growth in open water 276 IF ( iom_use( "dmibog" ) ) CALL iom_put( "dmibog" , - wfx_bog ) ! Sea-ice mass change through basal growth 277 IF ( iom_use( "dmisni" ) ) CALL iom_put( "dmisni" , - wfx_sni ) ! Sea-ice mass change through snow-to-ice conversion 278 IF ( iom_use( "dmisum" ) ) CALL iom_put( "dmisum" , - wfx_sum ) ! Sea-ice mass change through surface melting 279 IF ( iom_use( "dmibom" ) ) CALL iom_put( "dmibom" , - wfx_bom ) ! Sea-ice mass change through bottom melting 280 281 IF ( iom_use( "dmtsub" ) ) CALL iom_put( "dmtsub" , - wfx_sub ) ! Sea-ice mass change through evaporation and sublimation 282 IF ( iom_use( "dmssub" ) ) CALL iom_put( "dmssub" , - wfx_snw_sub ) ! Snow mass change through sublimation 283 IF ( iom_use( "dmisub" ) ) CALL iom_put( "dmisub" , - wfx_ice_sub ) ! Sea-ice mass change through sublimation 284 285 IF ( iom_use( "dmsspr" ) ) CALL iom_put( "dmsspr" , - wfx_spr ) ! Snow mass change through snow fall 286 IF ( iom_use( "dmsssi" ) ) CALL iom_put( "dmsssi" , wfx_sni*rhosn/rhoic ) ! Snow mass change through snow-to-ice conversion 287 288 IF ( iom_use( "dmsmel" ) ) CALL iom_put( "dmsmel" , - wfx_snw_sum ) ! Snow mass change through melt 289 IF ( iom_use( "dmsdyn" ) ) CALL iom_put( "dmsdyn" , diag_dms_dyn ) ! Snow mass change through dynamics 290 291 IF ( iom_use( "hfxconbo" ) ) CALL iom_put( "hfxconbo" , diag_fc_bo ) ! Bottom conduction flux 292 IF ( iom_use( "hfxconsu" ) ) CALL iom_put( "hfxconsu" , diag_fc_su ) ! Surface conduction flux 293 294 IF ( iom_use( "wfxtot" ) ) CALL iom_put( "wfxtot" , wfx_ice ) ! Total freshwater flux from sea ice 295 IF ( iom_use( "wfxsum" ) ) CALL iom_put( "wfxsum" , wfx_sum ) ! Freshwater flux from sea-ice surface 296 297 IF ( iom_use( "utau_oi" ) ) CALL iom_put( "utau_oi" , diag_utau_oi*zswi ) ! X-component of ocean stress on sea ice 298 IF ( iom_use( "vtau_oi" ) ) CALL iom_put( "vtau_oi" , diag_vtau_oi*zswi ) ! Y-component of ocean stress on sea ice 299 300 IF ( iom_use( "dssh_dx" ) ) CALL iom_put( "dssh_dx" , diag_dssh_dx*zswi ) ! Sea-surface tilt term in force balance (x-component) 301 IF ( iom_use( "dssh_dy" ) ) CALL iom_put( "dssh_dy" , diag_dssh_dy*zswi ) ! Sea-surface tilt term in force balance (y-component) 302 303 IF ( iom_use( "corstrx" ) ) CALL iom_put( "corstrx" , diag_corstrx*zswi ) ! Coriolis force term in force balance (x-component) 304 IF ( iom_use( "corstry" ) ) CALL iom_put( "corstry" , diag_corstry*zswi ) ! Coriolis force term in force balance (y-component) 305 306 IF ( iom_use( "intstrx" ) ) CALL iom_put( "intstrx" , diag_intstrx*zswi ) ! Internal force term in force balance (x-component) 307 IF ( iom_use( "intstry" ) ) CALL iom_put( "intstry" , diag_intstry*zswi ) ! Internal force term in force balance (y-component) 308 309 IF ( iom_use( "normstr" ) ) CALL iom_put( "normstr" , diag_sig1 *zswi ) ! Normal stress 310 IF ( iom_use( "sheastr" ) ) CALL iom_put( "sheastr" , diag_sig2 *zswi ) ! Shear stress 311 312 IF ( iom_use( "xmtrpice" ) ) CALL iom_put( "xmtrpice" , diag_xmtrp_ice ) ! X-component of sea-ice mass transport 313 IF ( iom_use( "ymtrpice" ) ) CALL iom_put( "ymtrpice" , diag_ymtrp_ice ) ! Y-component of sea-ice mass transport 314 315 IF ( iom_use( "xmtrpsnw" ) ) CALL iom_put( "xmtrpsnw" , diag_xmtrp_snw ) ! X-component of snow mass transport 316 IF ( iom_use( "ymtrpsnw" ) ) CALL iom_put( "ymtrpsnw" , diag_ymtrp_snw ) ! Y-component of snow mass transport 317 318 IF ( iom_use( "xatrp" ) ) CALL iom_put( "xatrp" , diag_xatrp ) ! X-component of ice area transport 319 IF ( iom_use( "yatrp" ) ) CALL iom_put( "yatrp" , diag_yatrp ) ! Y-component of ice area transport 320 321 !-------------------------------- 322 ! Global ice diagnostics (SIMIP) 323 !-------------------------------- 324 325 IF ( iom_use( "NH_icearea" ) .OR. iom_use( "NH_icevolu" ) .OR. iom_use( "NH_iceextt" ) ) THEN ! NH integrated diagnostics 326 327 WHERE( fcor > 0._wp ); zswi(:,:) = 1.0e-12 328 ELSEWHERE ; zswi(:,:) = 0. 329 END WHERE 330 331 zdiag_area_nh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) ) 332 zdiag_volu_nh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) ) 333 334 WHERE( fcor > 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12 335 ELSEWHERE ; zswi(:,:) = 0. 336 END WHERE 337 338 zdiag_extt_nh = glob_sum( zswi(:,:) * e12t(:,:) ) 339 340 IF ( iom_use( "NH_icearea" ) ) CALL iom_put( "NH_icearea" , zdiag_area_nh ) 341 IF ( iom_use( "NH_icevolu" ) ) CALL iom_put( "NH_icevolu" , zdiag_volu_nh ) 342 IF ( iom_use( "NH_iceextt" ) ) CALL iom_put( "NH_iceextt" , zdiag_extt_nh ) 343 344 ENDIF 345 346 IF ( iom_use( "SH_icearea" ) .OR. iom_use( "SH_icevolu" ) .OR. iom_use( "SH_iceextt" ) ) THEN ! SH integrated diagnostics 347 348 WHERE( fcor < 0._wp ); zswi(:,:) = 1.0e-12; 349 ELSEWHERE ; zswi(:,:) = 0. 350 END WHERE 351 352 zdiag_area_sh = glob_sum( at_i(:,:) * zswi(:,:) * e12t(:,:) ) 353 zdiag_volu_sh = glob_sum( vt_i(:,:) * zswi(:,:) * e12t(:,:) ) 354 355 WHERE( fcor < 0._wp .AND. at_i > 0.15 ); zswi(:,:) = 1.0e-12 356 ELSEWHERE ; zswi(:,:) = 0. 357 END WHERE 358 359 zdiag_extt_sh = glob_sum( zswi(:,:) * e12t(:,:) ) 360 361 IF ( iom_use( "SH_icearea" ) ) CALL iom_put( "SH_icearea", zdiag_area_sh ) 362 IF ( iom_use( "SH_icevolu" ) ) CALL iom_put( "SH_icevolu", zdiag_volu_sh ) 363 IF ( iom_use( "SH_iceextt" ) ) CALL iom_put( "SH_iceextt", zdiag_extt_sh ) 364 365 ENDIF 366 246 367 ! ! Create an output files (output.lim.abort.nc) if S < 0 or u > 20 m/s 247 368 ! IF( kindic < 0 ) CALL lim_wri_state( 'output.abort' ) … … 250 371 CALL wrk_dealloc( jpi, jpj, jpl, zswi2 ) 251 372 CALL wrk_dealloc( jpi, jpj , z2d, zswi ) 373 CALL wrk_dealloc( jpi, jpj , zfb, zamask, zamask15 ) 252 374 253 375 IF( nn_timing == 1 ) CALL timing_stop('limwri') -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/thd_ice.F90
r8233 r8239 60 60 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_snw_sum_1d 61 61 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_sub_1d 62 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_snw_sub_1d 63 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_ice_sub_1d 62 64 63 65 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: wfx_bog_1d … … 93 95 94 96 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: t_su_1d !: <==> the 2D t_su 97 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: t_si_1d !: <==> the 2D t_si 95 98 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: a_i_1d !: <==> the 2D a_i 96 99 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: ht_i_1d !: <==> the 2D ht_s … … 115 118 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: h_i_old !: ice thickness layer (m) 116 119 120 ! Conduction flux diagnostics (SIMIP) 121 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: diag_fc_bo_1d !: <==> the 2D diag_fc_bo 122 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: diag_fc_su_1d !: <==> the 2D diag_fc_su 123 117 124 INTEGER , PUBLIC :: jiindex_1d ! 1D index of debugging point 118 125 … … 129 136 !!---------------------------------------------------------------------! 130 137 INTEGER :: thd_ice_alloc ! return value 131 INTEGER :: ierr( 4), ii138 INTEGER :: ierr(5), ii 132 139 !!---------------------------------------------------------------------! 133 140 ierr(:) = 0 … … 150 157 & fhld_1d (jpij) , wfx_sub_1d (jpij) , wfx_bog_1d (jpij) , wfx_bom_1d(jpij) , & 151 158 & wfx_sum_1d(jpij) , wfx_sni_1d (jpij) , wfx_opw_1d (jpij) , wfx_res_1d(jpij) , & 159 & wfx_snw_sub_1d(jpij), wfx_ice_sub_1d(jpij) , & 152 160 & dqns_ice_1d(jpij) , evap_ice_1d (jpij), & 153 161 & qprec_ice_1d(jpij), qevap_ice_1d(jpij), i0 (jpij) , & 154 162 & sfx_bri_1d (jpij) , sfx_bog_1d (jpij) , sfx_bom_1d (jpij) , sfx_sum_1d (jpij), & 155 163 & sfx_sni_1d (jpij) , sfx_opw_1d (jpij) , sfx_res_1d (jpij) , sfx_sub_1d (jpij), & 156 & hicol_1d (jpij) , STAT=ierr(ii) ) 164 & dsm_i_fl_1d(jpij) , dsm_i_gd_1d(jpij) , dsm_i_se_1d(jpij) , & 165 & dsm_i_si_1d(jpij) , hicol_1d (jpij) , STAT=ierr(ii) ) 157 166 ! 158 167 ii = ii + 1 159 ALLOCATE( t_su_1d (jpij) , a_i_1d (jpij) , ht_i_1d (jpij) ,&168 ALLOCATE( t_su_1d (jpij) , t_si_1d (jpij) , a_i_1d (jpij) , ht_i_1d (jpij) , & 160 169 & ht_s_1d (jpij) , fc_su (jpij) , fc_bo_i (jpij) , & 161 170 & dh_s_tot (jpij) , dh_i_surf (jpij) , dh_i_sub (jpij) , & … … 168 177 & qh_i_old(jpij,0:nlay_i+1) , h_i_old(jpij,0:nlay_i+1) , STAT=ierr(ii) ) 169 178 ! 179 ii = ii + 1 180 ALLOCATE( diag_fc_bo_1d(jpij) , diag_fc_su_1d(jpij) , STAT=ierr(ii) ) 181 170 182 thd_ice_alloc = MAXVAL( ierr(:) ) 171 183 IF( thd_ice_alloc /= 0 ) CALL ctl_warn( 'thd_ice_alloc: failed to allocate arrays.' ) -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/IOM/iom.F90
r7768 r8239 195 195 #if defined key_lim3 || defined key_lim2 196 196 CALL iom_set_axis_attr( "ncatice", (/ (REAL(ji,wp), ji=1,jpl) /) ) 197 ! SIMIP diagnostics (4 main arctic straits) 198 CALL iom_set_axis_attr( "nstrait", (/ (REAL(ji,wp), ji=1,4) /) ) 197 199 #endif 198 200 CALL iom_set_axis_attr( "icbcla", class_num ) -
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
r8233 r8239 634 634 wfx_res(:,:) = 0._wp ; wfx_sub(:,:) = 0._wp 635 635 wfx_spr(:,:) = 0._wp ; wfx_lam(:,:) = 0._wp 636 wfx_snw_dyn(:,:) = 0._wp ; wfx_snw_sum(:,:) = 0._wp 637 wfx_snw_sub(:,:) = 0._wp ; wfx_ice_sub(:,:) = 0._wp 636 638 637 639 ! MV MP 2016 … … 654 656 diag_heat(:,:) = 0._wp ; diag_smvi(:,:) = 0._wp 655 657 diag_vice(:,:) = 0._wp ; diag_vsnw(:,:) = 0._wp 658 659 ! SIMIP diagnostics 660 diag_dms_dyn(:,:) = 0._wp ; diag_dmi_dyn(:,:) = 0._wp 661 diag_fc_bo(:,:) = 0._wp ; diag_fc_su(:,:) = 0._wp 656 662 657 663 tau_icebfr(:,:) = 0._wp; ! landfast ice param only (clem: important to keep the init here)
Note: See TracChangeset
for help on using the changeset viewer.