Changeset 15476
- Timestamp:
- 2021-11-08T10:31:34+01:00 (2 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/NEMO_4.0.4_FOAM_IAU_SI3_SIC/src/OCE/ASM/asminc.F90
r15412 r15476 35 35 USE sbc_oce ! Surface boundary condition variables. 36 36 USE diaobs , ONLY : calc_date ! Compute the calendar date on a given step 37 #if defined key_si3 37 #if defined key_si3 && defined key_asminc 38 38 USE phycst ! physical constants 39 39 USE ice1D ! sea-ice: thermodynamics variables … … 92 92 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ssh_bkg, ssh_bkginc ! Background sea surface height and its increment 93 93 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: seaice_bkginc ! Increment to the background sea ice conc 94 #if defined key_si3 94 REAL(wp) :: zhi_damin ! Ice thickness for new sea ice from DA increment 95 #if defined key_si3 && defined key_asminc 95 96 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: a_i_bkginc ! Increment to the background sea ice conc categories 96 97 LOGICAL, PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: incr_newice ! Mask .TRUE.=DA positive ice increment to open water 97 REAL(wp) :: zhi_damin ! Ice thickness for new sea ice from DA increment98 98 #endif 99 99 #if defined key_cice && defined key_asminc … … 336 336 ALLOCATE( ssh_bkginc (jpi,jpj) ) ; ssh_bkginc (:,:) = 0._wp 337 337 ALLOCATE( seaice_bkginc(jpi,jpj) ) ; seaice_bkginc(:,:) = 0._wp 338 #if defined key_si3 338 #if defined key_si3 && defined key_asminc 339 339 ALLOCATE( a_i_bkginc (jpi,jpj,jpl) ) ; a_i_bkginc (:,:,:) = 0._wp 340 340 ALLOCATE( incr_newice(jpi,jpj,jpl) ) ; incr_newice(:,:,:) = .FALSE. … … 414 414 ! very small increments are also set to zero 415 415 WHERE( ABS( seaice_bkginc(:,:) ) > 1.0e+10 .OR. & 416 & ABS( seaice_bkginc(:,:) ) < epsi03) seaice_bkginc(:,:) = 0.0_wp416 & ABS( seaice_bkginc(:,:) ) < 1.0e-03_wp ) seaice_bkginc(:,:) = 0.0_wp 417 417 418 #if defined key_si3 && defined key_asminc 418 419 IF (lwp) THEN 419 420 WRITE(numout,*) 420 WRITE(numout,*) 'asm_inc_init : Converting s eaice_bkginc to a_i_bkginc'421 WRITE(numout,*) 'asm_inc_init : Converting single category increment to multi-category' 421 422 WRITE(numout,*) '~~~~~~~~~~~~' 422 423 END IF … … 424 425 !single category increment for sea ice conc 425 426 !convert single category increment to multi category 426 a_i_bkginc = 0.0_wp427 427 at_i = SUM( a_i(:,:,:), DIM=3 ) 428 428 429 429 ! ensure zhi_damin lies within 1st category 430 IF ( zhi_damin > hi_max(1) ) zhi_damin = hi_max(1) - epsi02 430 IF ( zhi_damin > hi_max(1) ) THEN 431 IF (lwp) THEN 432 WRITE(numout,*) 433 WRITE(numout,*) 'minimum DA thickness > hmax(1): setting minimum DA thickness to hmax(1)' 434 WRITE(numout,*) '~~~~~~~~~~~~' 435 END IF 436 zhi_damin = hi_max(1) - epsi02 437 END IF 431 438 432 439 DO jj = 1, jpj 433 440 DO ji = 1, jpi 434 IF ( seaice_bkginc(ji,jj) > 0.0_wp ) THEN441 IF ( seaice_bkginc(ji,jj) > 0.0_wp ) THEN 435 442 !positive ice concentration increments are always 436 443 !added to the thinnest category of ice … … 440 447 !available category until it reaches zero concentration 441 448 !and then progressively removed from thicker categories 449 450 !NOTE: might consider the remote possibility that zremaining_increment 451 !left over is not zero, particulary if SI3 is being run 452 !as a single category 442 453 zremaining_increment = seaice_bkginc(ji,jj) 443 454 DO jl = 1, jpl … … 451 462 END DO 452 463 ! find model points where DA new ice should be added to open water 464 ! in any ice category 453 465 DO jl = 1,jpl 454 466 WHERE ( at_i(:,:) < epsi02 .AND. seaice_bkginc(:,:) > 0.0_wp ) … … 457 469 END DO 458 470 ENDIF 471 #endif 459 472 ! 460 473 CALL iom_close( inum ) … … 881 894 INTEGER, INTENT(in), OPTIONAL :: kindic ! flag for disabling the deallocation 882 895 ! 883 INTEGER :: it, j l896 INTEGER :: it, jk 884 897 REAL(wp) :: zincwgt ! IAU weight for current time step 885 #if defined key_si3 886 LOGICAL, SAVE :: logical_newice = .TRUE. ! Flag to add new ice from DA to open water 898 #if defined key_si3 && defined key_asminc 887 899 REAL(wp), DIMENSION(jpi,jpj,jpl) :: da_i ! change in ice concentration 888 900 REAL(wp), DIMENSION(jpi,jpj,jpl) :: dv_i ! change in ice volume 889 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z1_a_i ! inverse of old ice concentration890 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z1_v_i ! inverse of old ice volume901 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z1_a_i ! inverse of ice concentration before current IAU step 902 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z1_v_i ! inverse of ice volume before current IAU step 891 903 REAL(wp), DIMENSION(jpi,jpj) :: zhi_damin_2D ! array with DA thickness for incr_newice 892 904 #endif … … 911 923 ! Sea-ice : SI3 case 912 924 ! 913 #if defined key_si3 914 ! ensure zhi_damin lies within 1st category 915 IF ( zhi_damin > hi_max(1) ) zhi_damin = hi_max(1) - epsi02 916 925 #if defined key_si3 && defined key_asminc 917 926 ! compute the inverse of key sea ice variables 918 927 ! to be used later in the code … … 953 962 ! just do it once since next IAU steps assume that new ice has 954 963 ! already been added in 955 IF ( kt == nitiaustr_r .AND. logical_newice) THEN964 IF ( kt == nitiaustr_r ) THEN 956 965 957 966 ! assign zhi_damin to ice forming in open water … … 972 981 sv_i(:,:,:) = 0.0_wp 973 982 END WHERE 974 DO j l= 1, nlay_i983 DO jk = 1, nlay_i 975 984 WHERE ( incr_newice ) 976 e_i(:,:,j l,:) = 0.0_wp985 e_i(:,:,jk,:) = 0.0_wp 977 986 END WHERE 978 987 END DO 979 DO j l= 1, nlay_s988 DO jk = 1, nlay_s 980 989 WHERE ( incr_newice ) 981 e_s(:,:,j l,:) = 0.0_wp990 e_s(:,:,jk,:) = 0.0_wp 982 991 END WHERE 983 992 END DO … … 985 994 ! Initialisation of the salt content and ice enthalpy 986 995 ! set flag of new ice to false after this 987 callinit_new_ice_thd( zhi_damin_2D )996 CALL init_new_ice_thd( zhi_damin_2D ) 988 997 incr_newice(:,:,:) = .FALSE. 989 998 END IF … … 991 1000 ! maintain equivalent values for key prognostic variables 992 1001 v_s(:,:,:) = v_s(:,:,:) * da_i(:,:,:) 993 DO j l= 1, nlay_s994 e_s(:,:,j l,:) = e_s(:,:,jl,:) * da_i(:,:,:)1002 DO jk = 1, nlay_s 1003 e_s(:,:,jk,:) = e_s(:,:,jk,:) * da_i(:,:,:) 995 1004 END DO 996 1005 a_ip (:,:,:) = a_ip(:,:,:) * da_i(:,:,:) … … 999 1008 ! ice volume dependent variables 1000 1009 sv_i (:,:,:) = sv_i(:,:,:) * dv_i(:,:,:) 1001 DO j l= 1, nlay_i1002 e_i(:,:,j l,:) = e_i(:,:,jl,:) * dv_i(:,:,:)1010 DO jk = 1, nlay_i 1011 e_i(:,:,jk,:) = e_i(:,:,jk,:) * dv_i(:,:,:) 1003 1012 END DO 1004 1013 #endif … … 1011 1020 IF ( kt == nitiaufin_r ) THEN 1012 1021 DEALLOCATE( seaice_bkginc ) 1013 #if defined key_si3 1022 #if defined key_si3 && defined key_asminc 1014 1023 DEALLOCATE( incr_newice ) 1015 1024 DEALLOCATE( a_i_bkginc ) … … 1066 1075 !! forming at 1st category with thickness hi_new 1067 1076 !! 1068 !! ** Method : Apply SI3 thermodynamicequations to initialise1069 !! thermodynamic of new ice1070 !! 1071 !! ** Action : update thermodynamic variables1077 !! ** Method : Apply SI3 equations to initialise 1078 !! thermodynamics of new ice 1079 !! 1080 !! ** Action : update sea ice thermodynamics 1072 1081 !!---------------------------------------------------------------------- 1073 1082 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: hi_new ! total thickness of new ice 1074 1083 1075 INTEGER :: jj, ji, jk , jl1084 INTEGER :: jj, ji, jk 1076 1085 REAL(wp) :: ztmelts ! melting point 1077 1086 … … 1108 1117 CASE ( 2 ) ! Sice = F(z,t) [Vancoppenolle et al (2005)] 1109 1118 DO ji = 1, npti 1110 zs_newice(ji) = MIN( 4.606 + 0.91 / zh_newice(ji) , rn_simax , 0.5* sss_1d(ji) )1119 zs_newice(ji) = MIN( 4.606_wp + 0.91_wp / zh_newice(ji) , rn_simax , 0.5_wp * sss_1d(ji) ) 1111 1120 END DO 1112 1121 CASE ( 3 ) ! Sice = F(z) [multiyear ice] 1113 zs_newice(1:npti) = 2.3 1122 zs_newice(1:npti) = 2.3_wp 1114 1123 END SELECT 1115 1124 … … 1123 1132 DO ji = 1, npti 1124 1133 ztmelts = - rTmlt * zs_newice(ji) ! Melting point (C) 1125 e_i_1d(ji,:) = rhoi * ( rcpi * ( ztmelts - ( t_bo_1d(ji) - rt0 ) ) &1126 & + rLfus * ( 1.0 - ztmelts / MIN( t_bo_1d(ji) - rt0, -epsi10 ) ) &1134 e_i_1d(ji,:) = rhoi * ( rcpi * ( ztmelts - ( t_bo_1d(ji) - rt0 ) ) & 1135 & + rLfus * ( 1.0_wp - ztmelts / MIN( t_bo_1d(ji) - rt0, -epsi10 ) ) & 1127 1136 & - rcp * ztmelts ) 1128 1137 END DO
Note: See TracChangeset
for help on using the changeset viewer.
