Changeset 13064
- Timestamp:
- 2020-06-08T15:43:28+02:00 (4 years ago)
- Location:
- NEMO/branches/2020/r4.0-HEAD_r12713_clem_dan_fixcpl
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/r4.0-HEAD_r12713_clem_dan_fixcpl/cfgs/SHARED/field_def_nemo-ice.xml
r12919 r13064 83 83 <field id="icediv" long_name="Divergence of the sea-ice velocity field" standard_name="divergence_of_sea_ice_velocity" unit="s-1" /> 84 84 <field id="iceshe" long_name="Maximum shear of sea-ice velocity field" standard_name="maximum_shear_of_sea_ice_velocity" unit="s-1" /> 85 85 <field id="beta_evp" long_name="Relaxation parameter of ice rheology (beta)" standard_name="relaxation_parameter_of_ice_rheology" unit="" /> 86 86 87 <!-- surface heat fluxes --> 87 88 <field id="qt_ice" long_name="total heat flux at ice surface" standard_name="surface_downward_heat_flux_in_air" unit="W/m2" /> -
NEMO/branches/2020/r4.0-HEAD_r12713_clem_dan_fixcpl/src/ICE/icedyn_rhg_evp.F90
r12890 r13064 125 125 REAL(wp) :: ecc2, z1_ecc2 ! square of yield ellipse eccenticity 126 126 REAL(wp) :: zalph1, z1_alph1, zalph2, z1_alph2 ! alpha coef from Bouillon 2009 or Kimmritz 2017 127 REAl(wp) :: zbetau, zbetav 127 128 REAL(wp) :: zm1, zm2, zm3, zmassU, zmassV, zvU, zvV ! ice/snow mass and volume 128 129 REAL(wp) :: zdelta, zp_delf, zds2, zdt, zdt2, zdiv, zdiv2 ! temporary scalars … … 236 237 z1_ecc2 = 1._wp / ecc2 237 238 238 ! Time step for subcycling239 zdtevp = rdt_ice / REAL( nn_nevp )240 z1_dtevp = 1._wp / zdtevp241 242 239 ! alpha parameters (Bouillon 2009) 243 240 IF( .NOT. ln_aEVP ) THEN 244 zalph1 = ( 2._wp * rn_relast * rdt_ice ) * z1_dtevp 241 zdtevp = rdt_ice / REAL( nn_nevp ) 242 zalph1 = 2._wp * rn_relast * REAL( nn_nevp ) 245 243 zalph2 = zalph1 * z1_ecc2 246 244 247 245 z1_alph1 = 1._wp / ( zalph1 + 1._wp ) 248 246 z1_alph2 = 1._wp / ( zalph2 + 1._wp ) 247 ELSE 248 zdtevp = rdt_ice 249 ! zalpha parameters set later on adaptatively 249 250 ENDIF 251 z1_dtevp = 1._wp / zdtevp 250 252 251 253 ! Initialise stress tensor … … 360 362 ! 361 363 ! convergence test 362 IF( ln_rhg_chkcvg) THEN364 IF( ln_rhg_chkcvg ) THEN 363 365 DO jj = 1, jpj 364 366 DO ji = 1, jpi … … 408 410 zp_delt(ji,jj) = strength(ji,jj) / ( zdelta + rn_creepl ) 409 411 410 ! alpha & betafor aEVP412 ! alpha for aEVP 411 413 ! gamma = 0.5*P/(delta+creepl) * (c*pi)**2/Area * dt/m 412 414 ! alpha = beta = sqrt(4*gamma) … … 416 418 zalph2 = zalph1 417 419 z1_alph2 = z1_alph1 420 ! explicit: 421 ! z1_alph1 = 1._wp / zalph1 422 ! z1_alph2 = 1._wp / zalph1 423 ! zalph1 = zalph1 - 1._wp 424 ! zalph2 = zalph1 418 425 ENDIF 419 426 … … 426 433 CALL lbc_lnk( 'icedyn_rhg_evp', zp_delt, 'T', 1. ) 427 434 435 ! Save beta at T-points for further computations 436 IF( ln_aEVP ) THEN 437 DO jj = 1, jpj 438 DO ji = 1, jpi 439 zbeta(ji,jj) = MAX( 50._wp, rpi * SQRT( 0.5_wp * zp_delt(ji,jj) * r1_e1e2t(ji,jj) * zdt_m(ji,jj) ) ) 440 END DO 441 END DO 442 ENDIF 443 428 444 DO jj = 1, jpjm1 429 445 DO ji = 1, jpim1 430 446 431 ! alpha & betafor aEVP447 ! alpha for aEVP 432 448 IF( ln_aEVP ) THEN 433 zalph2 = MAX( 50._wp, rpi * SQRT( 0.5_wp * zp_delt(ji,jj) * r1_e1e2t(ji,jj) * zdt_m(ji,jj)) )449 zalph2 = MAX( zbeta(ji,jj), zbeta(ji+1,jj), zbeta(ji,jj+1), zbeta(ji+1,jj+1) ) 434 450 z1_alph2 = 1._wp / ( zalph2 + 1._wp ) 435 zbeta(ji,jj) = zalph2 451 ! explicit: 452 ! z1_alph2 = 1._wp / zalph2 453 ! zalph2 = zalph2 - 1._wp 436 454 ENDIF 437 455 … … 503 521 ! 504 522 IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017) 505 v_ice(ji,jj) = ( ( rswitch * ( zmV_t(ji,jj) * ( zbeta(ji,jj) * v_ice(ji,jj) + v_ice_b(ji,jj) ) & ! previous velocity 506 & + zRHS + zTauO * v_ice(ji,jj) ) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 507 & / MAX( zepsi, zmV_t(ji,jj) * ( zbeta(ji,jj) + 1._wp ) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 508 & + ( 1._wp - rswitch ) * v_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 509 & ) * zmsk01y(ji,jj) + v_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01y(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 523 zbetav = MAX( zbeta(ji,jj), zbeta(ji,jj+1) ) 524 v_ice(ji,jj) = ( ( rswitch * ( zmV_t(ji,jj) * ( zbetav * v_ice(ji,jj) + v_ice_b(ji,jj) ) & ! previous velocity 525 & + zRHS + zTauO * v_ice(ji,jj) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 526 & ) / MAX( zepsi, zmV_t(ji,jj) * ( zbetav + 1._wp ) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 527 & + ( 1._wp - rswitch ) * ( v_ice_b(ji,jj) & 528 & + v_ice (ji,jj) * MAX( 0._wp, zbetav - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 529 & ) / ( zbetav + 1._wp ) & 530 & ) * zmsk01y(ji,jj) + v_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01y(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 510 531 & ) * zmsk00y(ji,jj) 511 532 ELSE !--- ice velocity using EVP implicit formulation (cf Madec doc & Bouillon 2009) 512 v_ice(ji,jj) = ( ( rswitch * ( zmV_t(ji,jj)* v_ice(ji,jj) & ! previous velocity513 & + zRHS + zTauO * v_ice(ji,jj) )& ! F + tau_ia + Coriolis + spg + tau_io(only ocean part)514 & / MAX( zepsi, zmV_t(ji,jj) + zTauO - zTauB )& ! m/dt + tau_io(only ice part) + landfast515 & + ( 1._wp - rswitch ) *v_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0516 & ) * zmsk01y(ji,jj) + v_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01y(ji,jj) )& ! v_ice = v_oce/100 if mass < zmmin & conc < zamin517 & ) 533 v_ice(ji,jj) = ( ( rswitch * ( zmV_t(ji,jj) * v_ice(ji,jj) & ! previous velocity 534 & + zRHS + zTauO * v_ice(ji,jj) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 535 & ) / MAX( zepsi, zmV_t(ji,jj) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 536 & + ( 1._wp - rswitch ) * v_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 537 & ) * zmsk01y(ji,jj) + v_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01y(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 538 & ) * zmsk00y(ji,jj) 518 539 ENDIF 519 540 END DO … … 554 575 ! 555 576 IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017) 556 u_ice(ji,jj) = ( ( rswitch * ( zmU_t(ji,jj) * ( zbeta(ji,jj) * u_ice(ji,jj) + u_ice_b(ji,jj) ) & ! previous velocity 557 & + zRHS + zTauO * u_ice(ji,jj) ) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 558 & / MAX( zepsi, zmU_t(ji,jj) * ( zbeta(ji,jj) + 1._wp ) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 559 & + ( 1._wp - rswitch ) * u_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 560 & ) * zmsk01x(ji,jj) + u_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01x(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 577 zbetau = MAX( zbeta(ji,jj), zbeta(ji+1,jj) ) 578 u_ice(ji,jj) = ( ( rswitch * ( zmU_t(ji,jj) * ( zbetau * u_ice(ji,jj) + u_ice_b(ji,jj) ) & ! previous velocity 579 & + zRHS + zTauO * u_ice(ji,jj) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 580 & ) / MAX( zepsi, zmU_t(ji,jj) * ( zbetau + 1._wp ) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 581 & + ( 1._wp - rswitch ) * ( u_ice_b(ji,jj) & 582 & + u_ice (ji,jj) * MAX( 0._wp, zbetau - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 583 & ) / ( zbetau + 1._wp ) & 584 & ) * zmsk01x(ji,jj) + u_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01x(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 561 585 & ) * zmsk00x(ji,jj) 562 586 ELSE !--- ice velocity using EVP implicit formulation (cf Madec doc & Bouillon 2009) 563 u_ice(ji,jj) = ( ( rswitch * ( zmU_t(ji,jj)* u_ice(ji,jj) & ! previous velocity564 & + zRHS + zTauO * u_ice(ji,jj) )& ! F + tau_ia + Coriolis + spg + tau_io(only ocean part)565 & / MAX( zepsi, zmU_t(ji,jj) + zTauO - zTauB )& ! m/dt + tau_io(only ice part) + landfast566 & + ( 1._wp - rswitch ) *u_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0567 & ) * zmsk01x(ji,jj) + u_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01x(ji,jj) )& ! v_ice = v_oce/100 if mass < zmmin & conc < zamin568 & 587 u_ice(ji,jj) = ( ( rswitch * ( zmU_t(ji,jj) * u_ice(ji,jj) & ! previous velocity 588 & + zRHS + zTauO * u_ice(ji,jj) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 589 & ) / MAX( zepsi, zmU_t(ji,jj) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 590 & + ( 1._wp - rswitch ) * u_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 591 & ) * zmsk01x(ji,jj) + u_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01x(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 592 & ) * zmsk00x(ji,jj) 569 593 ENDIF 570 594 END DO … … 607 631 ! 608 632 IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017) 609 u_ice(ji,jj) = ( ( rswitch * ( zmU_t(ji,jj) * ( zbeta(ji,jj) * u_ice(ji,jj) + u_ice_b(ji,jj) ) & ! previous velocity 610 & + zRHS + zTauO * u_ice(ji,jj) ) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 611 & / MAX( zepsi, zmU_t(ji,jj) * ( zbeta(ji,jj) + 1._wp ) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 612 & + ( 1._wp - rswitch ) * u_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 613 & ) * zmsk01x(ji,jj) + u_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01x(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 633 zbetau = MAX( zbeta(ji,jj), zbeta(ji+1,jj) ) 634 u_ice(ji,jj) = ( ( rswitch * ( zmU_t(ji,jj) * ( zbetau * u_ice(ji,jj) + u_ice_b(ji,jj) ) & ! previous velocity 635 & + zRHS + zTauO * u_ice(ji,jj) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 636 & ) / MAX( zepsi, zmU_t(ji,jj) * ( zbetau + 1._wp ) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 637 & + ( 1._wp - rswitch ) * ( u_ice_b(ji,jj) & 638 & + u_ice (ji,jj) * MAX( 0._wp, zbetau - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 639 & ) / ( zbetau + 1._wp ) & 640 & ) * zmsk01x(ji,jj) + u_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01x(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 614 641 & ) * zmsk00x(ji,jj) 615 642 ELSE !--- ice velocity using EVP implicit formulation (cf Madec doc & Bouillon 2009) 616 u_ice(ji,jj) = ( ( rswitch * ( zmU_t(ji,jj)* u_ice(ji,jj) & ! previous velocity617 & + zRHS + zTauO * u_ice(ji,jj) )& ! F + tau_ia + Coriolis + spg + tau_io(only ocean part)618 & / MAX( zepsi, zmU_t(ji,jj) + zTauO - zTauB )& ! m/dt + tau_io(only ice part) + landfast619 & + ( 1._wp - rswitch ) *u_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0620 & ) * zmsk01x(ji,jj) + u_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01x(ji,jj) )& ! v_ice = v_oce/100 if mass < zmmin & conc < zamin621 & 643 u_ice(ji,jj) = ( ( rswitch * ( zmU_t(ji,jj) * u_ice(ji,jj) & ! previous velocity 644 & + zRHS + zTauO * u_ice(ji,jj) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 645 & ) / MAX( zepsi, zmU_t(ji,jj) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 646 & + ( 1._wp - rswitch ) * u_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 647 & ) * zmsk01x(ji,jj) + u_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01x(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 648 & ) * zmsk00x(ji,jj) 622 649 ENDIF 623 650 END DO … … 658 685 ! 659 686 IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017) 660 v_ice(ji,jj) = ( ( rswitch * ( zmV_t(ji,jj) * ( zbeta(ji,jj) * v_ice(ji,jj) + v_ice_b(ji,jj) ) & ! previous velocity 661 & + zRHS + zTauO * v_ice(ji,jj) ) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 662 & / MAX( zepsi, zmV_t(ji,jj) * ( zbeta(ji,jj) + 1._wp ) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 663 & + ( 1._wp - rswitch ) * v_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 664 & ) * zmsk01y(ji,jj) + v_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01y(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 687 zbetav = MAX( zbeta(ji,jj), zbeta(ji,jj+1) ) 688 v_ice(ji,jj) = ( ( rswitch * ( zmV_t(ji,jj) * ( zbetav * v_ice(ji,jj) + v_ice_b(ji,jj) ) & ! previous velocity 689 & + zRHS + zTauO * v_ice(ji,jj) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 690 & ) / MAX( zepsi, zmV_t(ji,jj) * ( zbetav + 1._wp ) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 691 & + ( 1._wp - rswitch ) * ( v_ice_b(ji,jj) & 692 & + v_ice (ji,jj) * MAX( 0._wp, zbetav - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 693 & ) / ( zbetav + 1._wp ) & 694 & ) * zmsk01y(ji,jj) + v_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01y(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 665 695 & ) * zmsk00y(ji,jj) 666 696 ELSE !--- ice velocity using EVP implicit formulation (cf Madec doc & Bouillon 2009) 667 v_ice(ji,jj) = ( ( rswitch * ( zmV_t(ji,jj)* v_ice(ji,jj) & ! previous velocity668 & + zRHS + zTauO * v_ice(ji,jj) )& ! F + tau_ia + Coriolis + spg + tau_io(only ocean part)669 & / MAX( zepsi, zmV_t(ji,jj) + zTauO - zTauB )& ! m/dt + tau_io(only ice part) + landfast670 & + ( 1._wp - rswitch ) *v_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0671 & ) * zmsk01y(ji,jj) + v_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01y(ji,jj) )& ! v_ice = v_oce/100 if mass < zmmin & conc < zamin672 & 697 v_ice(ji,jj) = ( ( rswitch * ( zmV_t(ji,jj) * v_ice(ji,jj) & ! previous velocity 698 & + zRHS + zTauO * v_ice(ji,jj) & ! F + tau_ia + Coriolis + spg + tau_io(only ocean part) 699 & ) / MAX( zepsi, zmV_t(ji,jj) + zTauO - zTauB ) & ! m/dt + tau_io(only ice part) + landfast 700 & + ( 1._wp - rswitch ) * v_ice(ji,jj) * MAX( 0._wp, 1._wp - zdtevp * rn_lf_relax ) & ! static friction => slow decrease to v=0 701 & ) * zmsk01y(ji,jj) + v_oce(ji,jj) * 0.01_wp * ( 1._wp - zmsk01y(ji,jj) ) & ! v_ice = v_oce/100 if mass < zmmin & conc < zamin 702 & ) * zmsk00y(ji,jj) 673 703 ENDIF 674 704 END DO … … 685 715 686 716 ! convergence test 687 IF(ln_rhg_chkcvg) THEN 688 CALL rhg_cvg( kt, jter, nn_nevp, u_ice, v_ice, zu_ice, zv_ice ) 689 ENDIF 717 IF( ln_rhg_chkcvg ) CALL rhg_cvg( kt, jter, nn_nevp, u_ice, v_ice, zu_ice, zv_ice ) 690 718 ! 691 719 ! ! ==================== ! 692 720 END DO ! end loop over jter ! 693 721 ! ! ==================== ! 722 IF( ln_aEVP ) CALL iom_put( 'beta_evp' , zbeta ) 694 723 ! 695 724 !------------------------------------------------------------------------------! … … 862 891 ! --- convergence tests --- ! 863 892 IF( ln_rhg_chkcvg ) THEN 864 IF( iom_use('uice_cvg') ) THEN ! output: u(t=nn_nevp) - u(t=nn_nevp-1) 865 CALL iom_put( 'uice_cvg', MAX( ABS( u_ice(:,:) - zu_ice(:,:) ) * umask(:,:,1) , & 866 & ABS( v_ice(:,:) - zv_ice(:,:) ) * vmask(:,:,1) ) * zmsk15(:,:) ) 893 IF( iom_use('uice_cvg') ) THEN 894 IF( ln_aEVP ) THEN ! output: beta * ( u(t=nn_nevp) - u(t=nn_nevp-1) ) 895 CALL iom_put( 'uice_cvg', MAX( ABS( u_ice(:,:) - zu_ice(:,:) ) * zbeta(:,:) * umask(:,:,1) , & 896 & ABS( v_ice(:,:) - zv_ice(:,:) ) * zbeta(:,:) * vmask(:,:,1) ) * zmsk15(:,:) ) 897 ELSE ! output: nn_nevp * ( u(t=nn_nevp) - u(t=nn_nevp-1) ) 898 CALL iom_put( 'uice_cvg', REAL( nn_nevp ) * MAX( ABS( u_ice(:,:) - zu_ice(:,:) ) * umask(:,:,1) , & 899 & ABS( v_ice(:,:) - zv_ice(:,:) ) * vmask(:,:,1) ) * zmsk15(:,:) ) 900 ENDIF 867 901 ENDIF 868 902 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.