SparseInversion not as robust as Inversion

[ilopezgp]

I tried using the SparseInversion algorithm instead of Inversion for an integration test in CEDMF.jl, and I run into a SingularException in L111.

The configuration used for the integration test is

γ = 1.0
threshold_eki = true
threshold_value = 5e-2
reg = 1e-2

In general, we are missing documentation for how to choose reg and how to make the convex optimization part of the update robust.

Stack trace:

Stacktrace:
--
  | [1] checknonsingular
  | @ /central/software/julia/1.7.0/share/julia/stdlib/v1.7/LinearAlgebra/src/factorization.jl:19 [inlined]
  | [2] checknonsingular
  | @ /central/software/julia/1.7.0/share/julia/stdlib/v1.7/LinearAlgebra/src/factorization.jl:21 [inlined]
  | [3] #lu!#146
  | @ /central/software/julia/1.7.0/share/julia/stdlib/v1.7/LinearAlgebra/src/lu.jl:82 [inlined]
  | [4] #lu#153
  | @ /central/software/julia/1.7.0/share/julia/stdlib/v1.7/LinearAlgebra/src/lu.jl:279 [inlined]
  | [5] lu (repeats 2 times)
  | @ /central/software/julia/1.7.0/share/julia/stdlib/v1.7/LinearAlgebra/src/lu.jl:278 [inlined]
  | [6] \(A::Matrix{Float64}, B::Diagonal{Float64, Vector{Float64}})
  | @ LinearAlgebra /central/software/julia/1.7.0/share/julia/stdlib/v1.7/LinearAlgebra/src/generic.jl:1142
  | [7] sparse_eki_update(ekp::EnsembleKalmanProcess{Float64, Int64, SparseInversion{Float64, Int64}}, u::Matrix{Float64}, g::Matrix{Float64}, y::Matrix{Float64}, obs_noise_cov::Matrix{Float64})
  | @ EnsembleKalmanProcesses /central/scratch/esm/slurm-buildkite/calibrateedmf-ci/depot/cpu/packages/EnsembleKalmanProcesses/OjLuD/src/SparseEnsembleKalmanInversion.jl:111
  | [8] (::EnsembleKalmanProcesses.var"#failsafe_update#26")(ekp::EnsembleKalmanProcess{Float64, Int64, SparseInversion{Float64, Int64}}, u::Matrix{Float64}, g::Matrix{Float64}, y::Matrix{Float64}, obs_noise_cov::Matrix{Float64}, failed_ens::Vector{Int64})
  | @ EnsembleKalmanProcesses /central/scratch/esm/slurm-buildkite/calibrateedmf-ci/depot/cpu/packages/EnsembleKalmanProcesses/OjLuD/src/SparseEnsembleKalmanInversion.jl:24
  | [9] update_ensemble!(ekp::EnsembleKalmanProcess{Float64, Int64, SparseInversion{Float64, Int64}}, g::Matrix{Float64}; cov_threshold::Float64, Δt_new::Float64, deterministic_forward_map::Bool, failed_ens::Nothing)
  | @ EnsembleKalmanProcesses /central/scratch/esm/slurm-buildkite/calibrateedmf-ci/depot/cpu/packages/EnsembleKalmanProcesses/OjLuD/src/SparseEnsembleKalmanInversion.jl:191

[jinlong83]

The value of γ may be too small in this case? When the sparsity constraint is too strong, we would quickly have collapsed ensemble that could trigger singularity issue in convex optimization. I'll investigate this issue a little bit more to see if we can better guarantee the robustness.

[ilopezgp]

The value of γ may be too small in this case? When the sparsity constraint is too strong, we would quickly have collapsed ensemble that could trigger singularity issue in convex optimization. I'll investigate this issue a little bit more to see if we can better guarantee the robustness.

That would be great, thanks! Is there any guideline we could document for the value of reg, or any heuristic to choose its value?

[jinlong83]

The value of γ may be too small in this case? When the sparsity constraint is too strong, we would quickly have collapsed ensemble that could trigger singularity issue in convex optimization. I'll investigate this issue a little bit more to see if we can better guarantee the robustness.

That would be great, thanks! Is there any guideline we could document for the value of reg, or any heuristic to choose its value?

The value of reg corresponds to ensemble inflation and it would be ideal to keep its value small. The value of γ needs a little bit more trial and error, as we want to have a γ that leads to a sparse solution but would not significantly increase the original loss (which can be visualized in more details by checking the agreement with observation data).

[odunbar]

Added a small note in the documentation for this
https://clima.github.io/EnsembleKalmanProcesses.jl/dev/ensemble_kalman_inversion/#Sparsity-Inducing-Ensemble-Kalman-Inversion

[odunbar]

As a continuation here, we have removed the DataMisfitController() tests as these have started leading to test failures.