Changeset 11944
- Timestamp:
- 2019-11-22T11:35:05+01:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5518_obs_oper_update_sit/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90
r11938 r11944 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 sea ice increment 23 !! seaice_asm_inc : Apply the sea ice concentration increment 24 !! sit_asm_inc : Apply the sea ice thickness increment 24 25 !!---------------------------------------------------------------------- 25 26 USE wrk_nemo ! Memory Allocation … … 49 50 PUBLIC dyn_asm_inc !: Apply the dynamic (u and v) increments 50 51 PUBLIC ssh_asm_inc !: Apply the SSH increment 51 PUBLIC seaice_asm_inc !: Apply the seaice increment 52 PUBLIC seaice_asm_inc !: Apply the seaice concentration increment 53 PUBLIC sit_asm_inc !: Apply the seaice thickness increment 54 PUBLIC bgc_asm_inc !: Apply the biogeochemistry increments 52 55 53 56 #if defined key_asminc … … 63 66 LOGICAL, PUBLIC :: ln_sshinc = .FALSE. !: No sea surface height assimilation increment 64 67 LOGICAL, PUBLIC :: ln_seaiceinc = .FALSE. !: No sea ice concentration increment 68 LOGICAL, PUBLIC :: ln_sitinc = .FALSE. !: No sea ice thickness increment 69 LOGICAL, PUBLIC :: lk_bgcinc = .FALSE. !: No biogeochemistry increments 65 70 LOGICAL, PUBLIC :: ln_salfix = .FALSE. !: Apply minimum salinity check 66 71 LOGICAL, PUBLIC :: ln_temnofreeze = .FALSE. !: Don't allow the temperature to drop below freezing … … 85 90 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ssh_bkg, ssh_bkginc ! Background sea surface height and its increment 86 91 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: seaice_bkginc ! Increment to the background sea ice conc 92 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: sit_bkginc ! Increment to the background sea ice thickness 87 93 88 94 !! * Substitutions … … 127 133 REAL(wp), POINTER, DIMENSION(:,:) :: hdiv ! 2D workspace 128 134 !! 129 NAMELIST/nam_asminc/ ln_bkgwri, 135 NAMELIST/nam_asminc/ ln_bkgwri, ln_balwri, & 130 136 & ln_trainc, ln_dyninc, ln_sshinc, & 137 & ln_phytobal, ln_slchltotinc, ln_slchldiainc, & 138 & ln_slchlnoninc, ln_schltotinc, ln_slphytotinc, & 139 & ln_slphydiainc, ln_slphynoninc, ln_spco2inc, & 140 & ln_sfco2inc, ln_plchltotinc, ln_pchltotinc, & 141 & ln_pno3inc, ln_psi4inc, ln_pdicinc, ln_palkinc, & 142 & ln_pphinc, ln_po2inc, & 131 143 & ln_asmdin, ln_asmiau, & 132 144 & nitbkg, nitdin, nitiaustr, nitiaufin, niaufn, & 133 & ln_salfix, salfixmin, nn_divdmp 145 & ln_salfix, salfixmin, nn_divdmp, & 146 & ln_seaiceinc, ln_sitinc, ln_temnofreeze, & 147 & mld_choice_bgc, rn_maxchlinc 134 148 !!---------------------------------------------------------------------- 135 149 … … 138 152 !----------------------------------------------------------------------- 139 153 ln_seaiceinc = .FALSE. 154 ln_sitinc = .FALSE. 140 155 ln_temnofreeze = .FALSE. 141 156 … … 154 169 WRITE(numout,*) 'asm_inc_init : Assimilation increment initialization :' 155 170 WRITE(numout,*) '~~~~~~~~~~~~' 156 WRITE(numout,*) ' Namelist nam asm: set assimilation increment parameters'171 WRITE(numout,*) ' Namelist nam_asminc : set assimilation increment parameters' 157 172 WRITE(numout,*) ' Logical switch for writing out background state ln_bkgwri = ', ln_bkgwri 173 WRITE(numout,*) ' Logical switch for writing out balancing increments ln_balwri = ', ln_balwri 158 174 WRITE(numout,*) ' Logical switch for applying tracer increments ln_trainc = ', ln_trainc 159 175 WRITE(numout,*) ' Logical switch for applying velocity increments ln_dyninc = ', ln_dyninc 160 176 WRITE(numout,*) ' Logical switch for applying SSH increments ln_sshinc = ', ln_sshinc 161 177 WRITE(numout,*) ' Logical switch for Direct Initialization (DI) ln_asmdin = ', ln_asmdin 162 WRITE(numout,*) ' Logical switch for applying sea ice increments ln_seaiceinc = ', ln_seaiceinc 178 WRITE(numout,*) ' Logical switch for applying SIC increments ln_seaiceinc = ', ln_seaiceinc 179 WRITE(numout,*) ' Logical switch for applying SIT increments ln_sitinc = ', ln_sitinc 180 WRITE(numout,*) ' Logical switch for phytoplankton balancing ln_phytobal = ', ln_phytobal 181 WRITE(numout,*) ' Logical switch for applying slchltot increments ln_slchltotinc = ', ln_slchltotinc 182 WRITE(numout,*) ' Logical switch for applying slchldia increments ln_slchldiainc = ', ln_slchldiainc 183 WRITE(numout,*) ' Logical switch for applying slchlnon increments ln_slchlnoninc = ', ln_slchlnoninc 184 WRITE(numout,*) ' Logical switch for applying schltot increments ln_schltotinc = ', ln_schltotinc 185 WRITE(numout,*) ' Logical switch for applying slphytot increments ln_slphytotinc = ', ln_slphytotinc 186 WRITE(numout,*) ' Logical switch for applying slphydia increments ln_slphydiainc = ', ln_slphydiainc 187 WRITE(numout,*) ' Logical switch for applying slphynon increments ln_slphynoninc = ', ln_slphynoninc 188 WRITE(numout,*) ' Logical switch for applying spco2 increments ln_spco2inc = ', ln_spco2inc 189 WRITE(numout,*) ' Logical switch for applying sfco2 increments ln_sfco2inc = ', ln_sfco2inc 190 WRITE(numout,*) ' Logical switch for applying plchltot increments ln_plchltotinc = ', ln_plchltotinc 191 WRITE(numout,*) ' Logical switch for applying pchltot increments ln_pchltotinc = ', ln_pchltotinc 192 WRITE(numout,*) ' Logical switch for applying pno3 increments ln_pno3inc = ', ln_pno3inc 193 WRITE(numout,*) ' Logical switch for applying psi4 increments ln_psi4inc = ', ln_psi4inc 194 WRITE(numout,*) ' Logical switch for applying pdic increments ln_pdicinc = ', ln_pdicinc 195 WRITE(numout,*) ' Logical switch for applying palk increments ln_palkinc = ', ln_palkinc 196 WRITE(numout,*) ' Logical switch for applying pph increments ln_pphinc = ', ln_pphinc 197 WRITE(numout,*) ' Logical switch for applying po2 increments ln_po2inc = ', ln_po2inc 163 198 WRITE(numout,*) ' Logical switch for Incremental Analysis Updating (IAU) ln_asmiau = ', ln_asmiau 164 199 WRITE(numout,*) ' Timestep of background in [0,nitend-nit000-1] nitbkg = ', nitbkg … … 169 204 WRITE(numout,*) ' Logical switch for ensuring that the sa > salfixmin ln_salfix = ', ln_salfix 170 205 WRITE(numout,*) ' Minimum salinity after applying the increments salfixmin = ', salfixmin 206 WRITE(numout,*) ' Choice of MLD for BGC assimilation mld_choice_bgc = ', mld_choice_bgc 207 WRITE(numout,*) ' Maximum absolute chlorophyll increment (<=0 = off) rn_maxchlinc = ', rn_maxchlinc 171 208 ENDIF 172 209 … … 216 253 217 254 IF ( ( ( .NOT. ln_asmdin ).AND.( .NOT. ln_asmiau ) ) & 218 .AND.( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc) )) & 219 & CALL ctl_stop( ' One or more of ln_trainc, ln_dyninc, ln_sshinc,', & 220 & ' and ln_seaiceinc is set to .true.', & 255 & .AND.( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc ).OR. & 256 & ( ln_sitinc ).OR.( lk_bgcinc ) )) & 257 & CALL ctl_stop( ' One or more of ln_trainc, ln_dyninc, ln_sshinc, ln_seaiceinc,', & 258 & ' ln_sitinc and ln_(bgc-variable)inc is set to .true.', & 221 259 & ' but ln_asmdin and ln_asmiau are both set to .false. :', & 222 260 & ' Inconsistent options') … … 227 265 228 266 IF ( ( .NOT. ln_trainc ).AND.( .NOT. ln_dyninc ).AND.( .NOT. ln_sshinc ).AND.( .NOT. ln_seaiceinc ) & 229 & 230 & CALL ctl_warn( ' ln_trainc, ln_dyninc, ln_sshinc ', &231 & ' and ln_seaiceinc are set to .false. :', &267 & .AND.( .NOT. ln_sitinc ).AND.( .NOT. lk_bgcinc ) ) & 268 & CALL ctl_warn( ' ln_trainc, ln_dyninc, ln_sshinc, ln_seaiceinc,', & 269 & ' ln_sitinc and ln_(bgc-variable)inc are set to .false. :', & 232 270 & ' The assimilation increments are not applied') 233 271 … … 333 371 ALLOCATE( ssh_bkginc(jpi,jpj) ) 334 372 ALLOCATE( seaice_bkginc(jpi,jpj)) 373 ALLOCATE( sit_bkginc(jpi,jpj) ) 335 374 #if defined key_asminc 336 375 ALLOCATE( ssh_iau(jpi,jpj) ) … … 342 381 ssh_bkginc(:,:) = 0.0 343 382 seaice_bkginc(:,:) = 0.0 383 sit_bkginc(:,:) = 0.0 344 384 #if defined key_asminc 345 385 ssh_iau(:,:) = 0.0 346 386 #endif 347 IF ( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc ) ) THEN 387 IF ( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc ) & 388 & .OR.( ln_sitinc ).OR.( lk_bgcinc ) ) THEN 348 389 349 390 !-------------------------------------------------------------------- … … 408 449 ! to allow for differences in masks 409 450 WHERE( ABS( ssh_bkginc(:,:) ) > 1.0e+10 ) ssh_bkginc(:,:) = 0.0 451 ENDIF 452 453 IF ( ln_sitinc ) THEN 454 CALL iom_get( inum, jpdom_autoglo, 'bckinsit', sit_bkginc, 1 ) 455 ! Apply the masks 456 sit_bkginc(:,:) = sit_bkginc(:,:) * tmask(:,:,1) 457 ! Set missing increments to 0.0 rather than 1e+20 458 ! to allow for differences in masks 459 WHERE( ABS( sit_bkginc(:,:) ) > 1.0e+10 ) sit_bkginc(:,:) = 0.0 410 460 ENDIF 411 461 … … 768 818 ! Perhaps the following call should be in step 769 819 IF ( ln_seaiceinc ) CALL seaice_asm_inc ( kt ) ! apply sea ice concentration increment 820 IF ( ln_sitinc ) CALL sit_asm_inc ( kt ) ! apply sea ice thickness increment 770 821 ! 771 822 END SUBROUTINE tra_asm_inc … … 925 976 END SUBROUTINE ssh_asm_inc 926 977 978 SUBROUTINE sit_asm_inc( kt, kindic ) 979 !!---------------------------------------------------------------------- 980 !! *** ROUTINE sit_asm_inc *** 981 !! 982 !! ** Purpose : Apply the sea ice thickness assimilation increment. 983 !! 984 !! ** Method : Direct initialization or Incremental Analysis Updating. 985 !! 986 !! ** Action : 987 !! 988 !!---------------------------------------------------------------------- 989 IMPLICIT NONE 990 ! 991 INTEGER, INTENT(in) :: kt ! Current time step 992 INTEGER, INTENT(in), OPTIONAL :: kindic ! flag for disabling the deallocation 993 ! 994 INTEGER :: it 995 REAL(wp) :: zincwgt ! IAU weight for current time step 996 ! #if defined key_lim2 997 ! REAL(wp), DIMENSION(jpi,jpj) :: zofrld, zohicif, zseaicendg, zhicifinc ! LIM 998 ! REAL(wp) :: zhicifmin = 0.5_wp ! ice minimum depth in metres 999 ! !!THICKNESS INCS NOT SET UP FOR LIM 1000 ! #endif 1001 !!---------------------------------------------------------------------- 1002 1003 IF ( ln_asmiau ) THEN 1004 1005 !-------------------------------------------------------------------- 1006 ! Incremental Analysis Updating 1007 !-------------------------------------------------------------------- 1008 1009 IF ( ( kt >= nitiaustr_r ).AND.( kt <= nitiaufin_r ) ) THEN 1010 1011 it = kt - nit000 + 1 1012 zincwgt = wgtiau(it) ! IAU weight for the current time step 1013 ! note this is not a tendency so should not be divided by rdt (as with the tracer and other increments) 1014 ! EF: Actually CICE is expecting a tendency so is divided by rdt below 1015 1016 IF(lwp) THEN 1017 WRITE(numout,*) 1018 WRITE(numout,*) 'sit_asm_inc : sea ice thick IAU at time step = ', & 1019 & kt,' with IAU weight = ', wgtiau(it) 1020 WRITE(numout,*) '~~~~~~~~~~~~' 1021 ENDIF 1022 1023 ! Sea-ice : LIM-3 case (to add) 1024 1025 ! #if defined key_lim2 1026 ! ! Sea-ice : LIM-2 case (to add if needed) 1027 ! zofrld (:,:) = frld(:,:) 1028 ! zohicif(:,:) = hicif(:,:) 1029 ! ! 1030 ! frld = MIN( MAX( frld (:,:) - seaice_bkginc(:,:) * zincwgt, 0.0_wp), 1.0_wp) 1031 ! pfrld = MIN( MAX( pfrld(:,:) - seaice_bkginc(:,:) * zincwgt, 0.0_wp), 1.0_wp) 1032 ! fr_i(:,:) = 1.0_wp - frld(:,:) ! adjust ice fraction 1033 ! ! 1034 ! zseaicendg(:,:) = zofrld(:,:) - frld(:,:) ! find out actual sea ice nudge applied 1035 ! ! 1036 ! ! Nudge sea ice depth to bring it up to a required minimum depth 1037 ! WHERE( zseaicendg(:,:) > 0.0_wp .AND. hicif(:,:) < zhicifmin ) 1038 ! zhicifinc(:,:) = (zhicifmin - hicif(:,:)) * zincwgt 1039 ! ELSEWHERE 1040 ! zhicifinc(:,:) = 0.0_wp 1041 ! END WHERE 1042 ! ! 1043 ! ! nudge ice depth 1044 ! hicif (:,:) = hicif (:,:) + zhicifinc(:,:) 1045 ! phicif(:,:) = phicif(:,:) + zhicifinc(:,:) 1046 ! ! 1047 ! ! seaice salinity balancing (to add) 1048 ! #endif 1049 1050 #if defined key_cice && defined key_asminc 1051 ! Sea-ice thickness : CICE case. Pass ice thickness increment tendency into CICE 1052 ndsit_da(:,:) = sit_bkginc(:,:) * zincwgt / rdt 1053 #endif 1054 1055 IF ( kt == nitiaufin_r ) THEN 1056 DEALLOCATE( sit_bkginc ) 1057 ENDIF 1058 1059 ELSE 1060 1061 #if defined key_cice && defined key_asminc 1062 ! Sea-ice thickness : CICE case. Zero ice increment tendency into CICE 1063 ndsit_da(:,:) = 0.0_wp 1064 #endif 1065 1066 ENDIF 1067 1068 ELSEIF ( ln_asmdin ) THEN 1069 1070 !-------------------------------------------------------------------- 1071 ! Direct Initialization 1072 !-------------------------------------------------------------------- 1073 1074 IF ( kt == nitdin_r ) THEN 1075 1076 neuler = 0 ! Force Euler forward step 1077 1078 ! Sea-ice : LIM-3 case (to add) 1079 1080 ! #if defined key_lim2 1081 ! ! Sea-ice : LIM-2 case (add if needed) 1082 ! zofrld(:,:)=frld(:,:) 1083 ! zohicif(:,:)=hicif(:,:) 1084 ! ! 1085 ! ! Initialize the now fields the background + increment 1086 ! frld (:,:) = MIN( MAX( frld(:,:) - seaice_bkginc(:,:), 0.0_wp), 1.0_wp) 1087 ! pfrld(:,:) = frld(:,:) 1088 ! fr_i (:,:) = 1.0_wp - frld(:,:) ! adjust ice fraction 1089 ! zseaicendg(:,:) = zofrld(:,:) - frld(:,:) ! find out actual sea ice nudge applied 1090 ! ! 1091 ! ! Nudge sea ice depth to bring it up to a required minimum depth 1092 ! WHERE( zseaicendg(:,:) > 0.0_wp .AND. hicif(:,:) < zhicifmin ) 1093 ! zhicifinc(:,:) = (zhicifmin - hicif(:,:)) * zincwgt 1094 ! ELSEWHERE 1095 ! zhicifinc(:,:) = 0.0_wp 1096 ! END WHERE 1097 ! ! 1098 ! ! nudge ice depth 1099 ! hicif (:,:) = hicif (:,:) + zhicifinc(:,:) 1100 ! phicif(:,:) = phicif(:,:) 1101 ! ! 1102 ! ! seaice salinity balancing (to add) 1103 ! #endif 1104 1105 #if defined key_cice && defined key_asminc 1106 ! Sea-ice thickness : CICE case. Pass ice thickness increment tendency into CICE 1107 ndsit_da(:,:) = sit_bkginc(:,:) / rdt 1108 #endif 1109 IF ( .NOT. PRESENT(kindic) ) THEN 1110 DEALLOCATE( sit_bkginc ) 1111 END IF 1112 1113 ELSE 1114 1115 #if defined key_cice && defined key_asminc 1116 ! Sea-ice thicnkness : CICE case. Zero ice thickness increment tendency into CICE 1117 ndsit_da(:,:) = 0.0_wp 1118 #endif 1119 1120 ENDIF 1121 1122 !#if defined defined key_lim2 || defined key_cice 1123 ! 1124 ! IF (ln_seaicebal ) THEN 1125 ! !! balancing salinity increments 1126 ! !! simple case from limflx.F90 (doesn't include a mass flux) 1127 ! !! assumption is that as ice concentration is reduced or increased 1128 ! !! the snow and ice depths remain constant 1129 ! !! note that snow is being created where ice concentration is being increased 1130 ! !! - could be more sophisticated and 1131 ! !! not do this (but would need to alter h_snow) 1132 ! 1133 ! usave(:,:,:)=sb(:,:,:) ! use array as a temporary store 1134 ! 1135 ! DO jj = 1, jpj 1136 ! DO ji = 1, jpi 1137 ! ! calculate change in ice and snow mass per unit area 1138 ! ! positive values imply adding salt to the ocean (results from ice formation) 1139 ! ! fwf : ice formation and melting 1140 ! 1141 ! zfons = ( -nfresh_da(ji,jj)*soce + nfsalt_da(ji,jj) )*rdt 1142 ! 1143 ! ! change salinity down to mixed layer depth 1144 ! mld=hmld_kara(ji,jj) 1145 ! 1146 ! ! prevent small mld 1147 ! ! less than 10m can cause salinity instability 1148 ! IF (mld < 10) mld=10 1149 ! 1150 ! ! set to bottom of a level 1151 ! DO jk = jpk-1, 2, -1 1152 ! IF ((mld > gdepw(ji,jj,jk)) .and. (mld < gdepw(ji,jj,jk+1))) THEN 1153 ! mld=gdepw(ji,jj,jk+1) 1154 ! jkmax=jk 1155 ! ENDIF 1156 ! ENDDO 1157 ! 1158 ! ! avoid applying salinity balancing in shallow water or on land 1159 ! ! 1160 ! 1161 ! ! dsal_ocn (psu kg m^-2) / (kg m^-3 * m) 1162 ! 1163 ! dsal_ocn=0.0_wp 1164 ! sal_thresh=5.0_wp ! minimum salinity threshold for salinity balancing 1165 ! 1166 ! if (tmask(ji,jj,1) > 0 .AND. tmask(ji,jj,jkmax) > 0 ) & 1167 ! dsal_ocn = zfons / (rhop(ji,jj,1) * mld) 1168 ! 1169 ! ! put increments in for levels in the mixed layer 1170 ! ! but prevent salinity below a threshold value 1171 ! 1172 ! DO jk = 1, jkmax 1173 ! 1174 ! IF (dsal_ocn > 0.0_wp .or. sb(ji,jj,jk)+dsal_ocn > sal_thresh) THEN 1175 ! sb(ji,jj,jk) = sb(ji,jj,jk) + dsal_ocn 1176 ! sn(ji,jj,jk) = sn(ji,jj,jk) + dsal_ocn 1177 ! ENDIF 1178 ! 1179 ! ENDDO 1180 ! 1181 ! ! ! salt exchanges at the ice/ocean interface 1182 ! ! zpmess = zfons / rdt_ice ! rdt_ice is ice timestep 1183 ! ! 1184 ! !! Adjust fsalt. A +ve fsalt means adding salt to ocean 1185 ! !! fsalt(ji,jj) = fsalt(ji,jj) + zpmess ! adjust fsalt 1186 ! !! 1187 ! !! emps(ji,jj) = emps(ji,jj) + zpmess ! or adjust emps (see icestp1d) 1188 ! !! ! E-P (kg m-2 s-2) 1189 ! ! emp(ji,jj) = emp(ji,jj) + zpmess ! E-P (kg m-2 s-2) 1190 ! ENDDO !ji 1191 ! ENDDO !jj! 1192 ! 1193 ! ENDIF !ln_seaicebal 1194 ! 1195 !#endif 1196 1197 1198 END SUBROUTINE sit_asm_inc 1199 927 1200 SUBROUTINE seaice_asm_inc( kt, kindic ) 928 1201 !!---------------------------------------------------------------------- 929 1202 !! *** ROUTINE seaice_asm_inc *** 930 1203 !! 931 !! ** Purpose : Apply the sea ice assimilation increment.1204 !! ** Purpose : Apply the sea ice concentration assimilation increment. 932 1205 !! 933 1206 !! ** Method : Direct initialization or Incremental Analysis Updating.
Note: See TracChangeset
for help on using the changeset viewer.