Changeset 10181 for branches/UKMO/dev_r5518_obs_oper_update_icethick/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90
- Timestamp:
- 2018-10-09T11:29:47+02:00 (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_obs_oper_update_icethick/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90
r9987 r10181 21 21 !! dyn_asm_inc : Apply the dynamic (u and v) increments 22 22 !! ssh_asm_inc : Apply the SSH increment 23 !! seaice_asm_inc : Apply the seaice increment 23 !! seaice_asm_inc : Apply the seaice concentration increment 24 !! sit_asm_inc : Apply the sea ice thickness increment 24 25 !!---------------------------------------------------------------------- 25 26 USE wrk_nemo ! Memory Allocation … … 41 42 #if defined key_cice && defined key_asminc 42 43 USE sbc_ice, ONLY : & ! CICE Ice model variables 43 & ndaice_da, n fresh_da, nfsalt_da44 & ndaice_da, ndsit_da, nfresh_da, nfsalt_da 44 45 #endif 45 46 USE sbc_oce ! Surface boundary condition variables. … … 53 54 PUBLIC dyn_asm_inc !: Apply the dynamic (u and v) increments 54 55 PUBLIC ssh_asm_inc !: Apply the SSH increment 55 PUBLIC seaice_asm_inc !: Apply the seaice increment 56 PUBLIC seaice_asm_inc !: Apply the seaice concentration increment 57 PUBLIC sit_asm_inc !: Apply the seaice thickness increment 56 58 57 59 #if defined key_asminc … … 66 68 LOGICAL, PUBLIC :: ln_dyninc = .FALSE. !: No dynamics (u and v) assimilation increments 67 69 LOGICAL, PUBLIC :: ln_sshinc = .FALSE. !: No sea surface height assimilation increment 68 LOGICAL, PUBLIC :: ln_seaiceinc !: No sea ice concentration increment 70 LOGICAL, PUBLIC :: ln_seaiceinc = .FALSE. !: No sea ice concentration increment 71 LOGICAL, PUBLIC :: ln_sitinc = .FALSE. !: No sea ice thickness increment 69 72 LOGICAL, PUBLIC :: ln_salfix = .FALSE. !: Apply minimum salinity check 70 73 LOGICAL, PUBLIC :: ln_temnofreeze = .FALSE. !: Don't allow the temperature to drop below freezing … … 89 92 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ssh_bkg, ssh_bkginc ! Background sea surface height and its increment 90 93 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: seaice_bkginc ! Increment to the background sea ice conc 94 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: sit_bkginc ! Increment to the background sea ice thickness 91 95 92 96 !! * Substitutions … … 136 140 & nitbkg, nitdin, nitiaustr, nitiaufin, niaufn, & 137 141 & ln_salfix, salfixmin, nn_divdmp, & 138 & ln_seaiceinc, ln_ temnofreeze142 & ln_seaiceinc, ln_sitinc, ln_temnofreeze 139 143 !!---------------------------------------------------------------------- 140 144 … … 142 146 ! Read Namelist nam_asminc : assimilation increment interface 143 147 !----------------------------------------------------------------------- 148 ln_sitinc = .FALSE. 144 149 ln_seaiceinc = .FALSE. 145 150 ln_temnofreeze = .FALSE. … … 165 170 WRITE(numout,*) ' Logical switch for applying SSH increments ln_sshinc = ', ln_sshinc 166 171 WRITE(numout,*) ' Logical switch for Direct Initialization (DI) ln_asmdin = ', ln_asmdin 167 WRITE(numout,*) ' Logical switch for applying sea ice increments ln_seaiceinc = ', ln_seaiceinc 172 WRITE(numout,*) ' Logical switch for applying SIC increments ln_seaiceinc = ', ln_seaiceinc 173 WRITE(numout,*) ' Logical switch for applying SIT increments ln_sitinc = ', ln_sitinc 168 174 WRITE(numout,*) ' Logical switch for Incremental Analysis Updating (IAU) ln_asmiau = ', ln_asmiau 169 175 WRITE(numout,*) ' Timestep of background in [0,nitend-nit000-1] nitbkg = ', nitbkg … … 221 227 222 228 IF ( ( ( .NOT. ln_asmdin ).AND.( .NOT. ln_asmiau ) ) & 223 .AND.( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ) .OR. ( ln_seaiceinc) )) & 224 & CALL ctl_stop( ' One or more of ln_trainc, ln_dyninc, ln_sshinc and ln_seaiceinc is set to .true.', & 229 & .AND.( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ) & 230 & .OR.( ln_sitinc ).OR.( ln_seaiceinc ) )) & 231 & CALL ctl_stop( ' One or more of ln_trainc, ln_dyninc, ln_sshinc,', & 232 & ' ln_sitinc and ln_seaiceinc is set to .true.', & 225 233 & ' but ln_asmdin and ln_asmiau are both set to .false. :', & 226 234 & ' Inconsistent options') … … 230 238 & ' Type IAU weighting function is invalid') 231 239 232 IF ( ( .NOT. ln_trainc ).AND.( .NOT. ln_dyninc ).AND.( .NOT. ln_sshinc ).AND.( .NOT. ln_seaiceinc ) & 240 IF ( ( .NOT. ln_trainc ).AND.( .NOT. ln_dyninc ).AND.( .NOT. ln_sshinc ) & 241 & .AND.( .NOT. ln_sitinc ).AND.( .NOT. ln_seaiceinc ) & 233 242 & ) & 234 & CALL ctl_warn( ' ln_trainc, ln_dyninc, ln_sshinc and ln_seaiceinc are set to .false. :', & 243 & CALL ctl_warn( ' ln_trainc, ln_dyninc, ln_sshinc, ln_sitinc', & 244 & ' and ln_seaiceinc are set to .false. :', & 235 245 & ' The assimilation increments are not applied') 236 246 … … 336 346 ALLOCATE( ssh_bkginc(jpi,jpj) ) 337 347 ALLOCATE( seaice_bkginc(jpi,jpj)) 348 ALLOCATE( sit_bkginc(jpi,jpj) ) 338 349 #if defined key_asminc 339 350 ALLOCATE( ssh_iau(jpi,jpj) ) … … 345 356 ssh_bkginc(:,:) = 0.0 346 357 seaice_bkginc(:,:) = 0.0 358 sit_bkginc(:,:) = 0.0 347 359 #if defined key_asminc 348 360 ssh_iau(:,:) = 0.0 349 361 #endif 350 IF ( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_s eaiceinc ) ) THEN362 IF ( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_sitinc ).OR.( ln_seaiceinc ) ) THEN 351 363 352 364 !-------------------------------------------------------------------- … … 411 423 ! to allow for differences in masks 412 424 WHERE( ABS( ssh_bkginc(:,:) ) > 1.0e+10 ) ssh_bkginc(:,:) = 0.0 425 ENDIF 426 427 IF ( ln_sitinc ) THEN 428 CALL iom_get( inum, jpdom_autoglo, 'bckinsit', sit_bkginc, 1 ) 429 ! Apply the masks 430 sit_bkginc(:,:) = sit_bkginc(:,:) * tmask(:,:,1) 431 ! Set missing increments to 0.0 rather than 1e+20 432 ! to allow for differences in masks 433 WHERE( ABS( sit_bkginc(:,:) ) > 1.0e+10 ) sit_bkginc(:,:) = 0.0 413 434 ENDIF 414 435 … … 769 790 ! 770 791 ENDIF 771 ! Perhaps the following call should be in step 772 IF ( ln_seaiceinc ) CALL seaice_asm_inc ( kt ) ! apply sea ice concentration increment 792 ! Perhaps the following call should be in step 793 IF ( ln_seaiceinc ) CALL seaice_asm_inc ( kt ) ! apply sea ice concentration increment 794 IF ( ln_sitinc ) CALL sit_asm_inc ( kt ) ! apply sea ice thickness increment 773 795 ! 774 796 END SUBROUTINE tra_asm_inc … … 932 954 END SUBROUTINE ssh_asm_inc 933 955 934 956 SUBROUTINE sit_asm_inc( kt, kindic ) 957 !!---------------------------------------------------------------------- 958 !! *** ROUTINE sit_asm_inc *** 959 !! 960 !! ** Purpose : Apply the sea ice thickness assimilation increment. 961 !! 962 !! ** Method : Direct initialization or Incremental Analysis Updating. 963 !! 964 !! ** Action : 965 !! 966 !!---------------------------------------------------------------------- 967 IMPLICIT NONE 968 ! 969 INTEGER, INTENT(in) :: kt ! Current time step 970 INTEGER, INTENT(in), OPTIONAL :: kindic ! flag for disabling the deallocation 971 ! 972 INTEGER :: it 973 REAL(wp) :: zincwgt ! IAU weight for current time step 974 ! #if defined key_lim2 975 ! REAL(wp), DIMENSION(jpi,jpj) :: zofrld, zohicif, zseaicendg, zhicifinc ! LIM 976 ! REAL(wp) :: zhicifmin = 0.5_wp ! ice minimum depth in metres 977 ! !!THICKNESS INCS NOT SET UP FOR LIM 978 ! #endif 979 !!---------------------------------------------------------------------- 980 981 IF ( ln_asmiau ) THEN 982 983 !-------------------------------------------------------------------- 984 ! Incremental Analysis Updating 985 !-------------------------------------------------------------------- 986 987 IF ( ( kt >= nitiaustr_r ).AND.( kt <= nitiaufin_r ) ) THEN 988 989 it = kt - nit000 + 1 990 zincwgt = wgtiau(it) ! IAU weight for the current time step 991 ! note this is not a tendency so should not be divided by rdt (as with the tracer and other increments) 992 993 IF(lwp) THEN 994 WRITE(numout,*) 995 WRITE(numout,*) 'sit_asm_inc : sea ice thick IAU at time step = ', & 996 & kt,' with IAU weight = ', wgtiau(it) 997 WRITE(numout,*) '~~~~~~~~~~~~' 998 ENDIF 999 1000 ! Sea-ice : LIM-3 case (to add) 1001 1002 ! #if defined key_lim2 1003 ! ! Sea-ice : LIM-2 case (to add if needed) 1004 ! zofrld (:,:) = frld(:,:) 1005 ! zohicif(:,:) = hicif(:,:) 1006 ! ! 1007 ! frld = MIN( MAX( frld (:,:) - seaice_bkginc(:,:) * zincwgt, 0.0_wp), 1.0_wp) 1008 ! pfrld = MIN( MAX( pfrld(:,:) - seaice_bkginc(:,:) * zincwgt, 0.0_wp), 1.0_wp) 1009 ! fr_i(:,:) = 1.0_wp - frld(:,:) ! adjust ice fraction 1010 ! ! 1011 ! zseaicendg(:,:) = zofrld(:,:) - frld(:,:) ! find out actual sea ice nudge applied 1012 ! ! 1013 ! ! Nudge sea ice depth to bring it up to a required minimum depth 1014 ! WHERE( zseaicendg(:,:) > 0.0_wp .AND. hicif(:,:) < zhicifmin ) 1015 ! zhicifinc(:,:) = (zhicifmin - hicif(:,:)) * zincwgt 1016 ! ELSEWHERE 1017 ! zhicifinc(:,:) = 0.0_wp 1018 ! END WHERE 1019 ! ! 1020 ! ! nudge ice depth 1021 ! hicif (:,:) = hicif (:,:) + zhicifinc(:,:) 1022 ! phicif(:,:) = phicif(:,:) + zhicifinc(:,:) 1023 ! ! 1024 ! ! seaice salinity balancing (to add) 1025 ! #endif 1026 1027 #if defined key_cice && defined key_asminc 1028 ! Sea-ice thickness : CICE case. Pass ice thickness increment tendency into CICE 1029 ndsit_da(:,:) = sit_bkginc(:,:) * zincwgt / rdt 1030 #endif 1031 1032 IF ( kt == nitiaufin_r ) THEN 1033 DEALLOCATE( sit_bkginc ) 1034 ENDIF 1035 1036 ELSE 1037 1038 #if defined key_cice && defined key_asminc 1039 ! Sea-ice thickness : CICE case. Zero ice increment tendency into CICE 1040 ndsit_da(:,:) = 0.0_wp 1041 #endif 1042 1043 ENDIF 1044 1045 ELSEIF ( ln_asmdin ) THEN 1046 1047 !-------------------------------------------------------------------- 1048 ! Direct Initialization 1049 !-------------------------------------------------------------------- 1050 1051 IF ( kt == nitdin_r ) THEN 1052 1053 neuler = 0 ! Force Euler forward step 1054 1055 ! Sea-ice : LIM-3 case (to add) 1056 1057 ! #if defined key_lim2 1058 ! ! Sea-ice : LIM-2 case (add if needed) 1059 ! zofrld(:,:)=frld(:,:) 1060 ! zohicif(:,:)=hicif(:,:) 1061 ! ! 1062 ! ! Initialize the now fields the background + increment 1063 ! frld (:,:) = MIN( MAX( frld(:,:) - seaice_bkginc(:,:), 0.0_wp), 1.0_wp) 1064 ! pfrld(:,:) = frld(:,:) 1065 ! fr_i (:,:) = 1.0_wp - frld(:,:) ! adjust ice fraction 1066 ! zseaicendg(:,:) = zofrld(:,:) - frld(:,:) ! find out actual sea ice nudge applied 1067 ! ! 1068 ! ! Nudge sea ice depth to bring it up to a required minimum depth 1069 ! WHERE( zseaicendg(:,:) > 0.0_wp .AND. hicif(:,:) < zhicifmin ) 1070 ! zhicifinc(:,:) = (zhicifmin - hicif(:,:)) * zincwgt 1071 ! ELSEWHERE 1072 ! zhicifinc(:,:) = 0.0_wp 1073 ! END WHERE 1074 ! ! 1075 ! ! nudge ice depth 1076 ! hicif (:,:) = hicif (:,:) + zhicifinc(:,:) 1077 ! phicif(:,:) = phicif(:,:) 1078 ! ! 1079 ! ! seaice salinity balancing (to add) 1080 ! #endif 1081 1082 #if defined key_cice && defined key_asminc 1083 ! Sea-ice thickness : CICE case. Pass ice thickness increment tendency into CICE 1084 ndsit_da(:,:) = sit_bkginc(:,:) / rdt 1085 #endif 1086 IF ( .NOT. PRESENT(kindic) ) THEN 1087 DEALLOCATE( sit_bkginc ) 1088 END IF 1089 1090 ELSE 1091 1092 #if defined key_cice && defined key_asminc 1093 ! Sea-ice thicnkness : CICE case. Zero ice thickness increment tendency into CICE 1094 ndsit_da(:,:) = 0.0_wp 1095 #endif 1096 1097 ENDIF 1098 1099 !#if defined defined key_lim2 || defined key_cice 1100 ! 1101 ! IF (ln_seaicebal ) THEN 1102 ! !! balancing salinity increments 1103 ! !! simple case from limflx.F90 (doesn't include a mass flux) 1104 ! !! assumption is that as ice concentration is reduced or increased 1105 ! !! the snow and ice depths remain constant 1106 ! !! note that snow is being created where ice concentration is being increased 1107 ! !! - could be more sophisticated and 1108 ! !! not do this (but would need to alter h_snow) 1109 ! 1110 ! usave(:,:,:)=sb(:,:,:) ! use array as a temporary store 1111 ! 1112 ! DO jj = 1, jpj 1113 ! DO ji = 1, jpi 1114 ! ! calculate change in ice and snow mass per unit area 1115 ! ! positive values imply adding salt to the ocean (results from ice formation) 1116 ! ! fwf : ice formation and melting 1117 ! 1118 ! zfons = ( -nfresh_da(ji,jj)*soce + nfsalt_da(ji,jj) )*rdt 1119 ! 1120 ! ! change salinity down to mixed layer depth 1121 ! mld=hmld_kara(ji,jj) 1122 ! 1123 ! ! prevent small mld 1124 ! ! less than 10m can cause salinity instability 1125 ! IF (mld < 10) mld=10 1126 ! 1127 ! ! set to bottom of a level 1128 ! DO jk = jpk-1, 2, -1 1129 ! IF ((mld > gdepw(ji,jj,jk)) .and. (mld < gdepw(ji,jj,jk+1))) THEN 1130 ! mld=gdepw(ji,jj,jk+1) 1131 ! jkmax=jk 1132 ! ENDIF 1133 ! ENDDO 1134 ! 1135 ! ! avoid applying salinity balancing in shallow water or on land 1136 ! ! 1137 ! 1138 ! ! dsal_ocn (psu kg m^-2) / (kg m^-3 * m) 1139 ! 1140 ! dsal_ocn=0.0_wp 1141 ! sal_thresh=5.0_wp ! minimum salinity threshold for salinity balancing 1142 ! 1143 ! if (tmask(ji,jj,1) > 0 .AND. tmask(ji,jj,jkmax) > 0 ) & 1144 ! dsal_ocn = zfons / (rhop(ji,jj,1) * mld) 1145 ! 1146 ! ! put increments in for levels in the mixed layer 1147 ! ! but prevent salinity below a threshold value 1148 ! 1149 ! DO jk = 1, jkmax 1150 ! 1151 ! IF (dsal_ocn > 0.0_wp .or. sb(ji,jj,jk)+dsal_ocn > sal_thresh) THEN 1152 ! sb(ji,jj,jk) = sb(ji,jj,jk) + dsal_ocn 1153 ! sn(ji,jj,jk) = sn(ji,jj,jk) + dsal_ocn 1154 ! ENDIF 1155 ! 1156 ! ENDDO 1157 ! 1158 ! ! ! salt exchanges at the ice/ocean interface 1159 ! ! zpmess = zfons / rdt_ice ! rdt_ice is ice timestep 1160 ! ! 1161 ! !! Adjust fsalt. A +ve fsalt means adding salt to ocean 1162 ! !! fsalt(ji,jj) = fsalt(ji,jj) + zpmess ! adjust fsalt 1163 ! !! 1164 ! !! emps(ji,jj) = emps(ji,jj) + zpmess ! or adjust emps (see icestp1d) 1165 ! !! ! E-P (kg m-2 s-2) 1166 ! ! emp(ji,jj) = emp(ji,jj) + zpmess ! E-P (kg m-2 s-2) 1167 ! ENDDO !ji 1168 ! ENDDO !jj! 1169 ! 1170 ! ENDIF !ln_seaicebal 1171 ! 1172 !#endif 1173 1174 ENDIF 1175 1176 END SUBROUTINE sit_asm_inc 1177 935 1178 SUBROUTINE seaice_asm_inc( kt, kindic ) 936 1179 !!---------------------------------------------------------------------- 937 1180 !! *** ROUTINE seaice_asm_inc *** 938 1181 !! 939 !! ** Purpose : Apply the sea ice assimilation increment.1182 !! ** Purpose : Apply the sea ice concentration assimilation increment. 940 1183 !! 941 1184 !! ** Method : Direct initialization or Incremental Analysis Updating.
Note: See TracChangeset
for help on using the changeset viewer.
