Skip to content

Commit

Permalink
added Tikhonov regularization matrix for Kaczmarz
Browse files Browse the repository at this point in the history
  • Loading branch information
mboberg committed May 13, 2022
1 parent cd93fd6 commit 2eccca0
Showing 1 changed file with 14 additions and 1 deletion.
15 changes: 14 additions & 1 deletion src/Kaczmarz.jl
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@ mutable struct Kaczmarz{matT,T,U,Tsparse} <: AbstractLinearSolver
sparseTrafo::Tsparse
iterations::Int64
constraintMask::Union{Nothing,Vector{Bool}}
regMatrix::Union{Nothing,Vector{U}} # Tikhonov regularization matrix
end

"""
Expand Down Expand Up @@ -58,6 +59,7 @@ function Kaczmarz(S; b=nothing, λ=[0.0], reg = nothing
, seed::Int=1234
, iterations::Int64=10
, constraintMask=nothing
, regMatrix=nothing
, kargs...)

T = real(eltype(S))
Expand All @@ -76,6 +78,12 @@ function Kaczmarz(S; b=nothing, λ=[0.0], reg = nothing
error("Kaczmarz only supports L2 regularizer")
end

# Apply Tikhonov regularization matrix
if regMatrix != nothing
regMatrix = T.(regMatrix) # make sure regMatrix has the same element type as S
S = transpose(1 ./ sqrt.(regMatrix)) .* S # apply Tikhonov regularization to system matrix
end

# make sure weights are not empty
w = (weights!=nothing ? weights : ones(T,size(S,1)))

Expand All @@ -97,7 +105,7 @@ function Kaczmarz(S; b=nothing, λ=[0.0], reg = nothing

return Kaczmarz(S,u,reg,denom,rowindex,rowIndexCycle,cl,vl,εw,τl,αl
,T.(w),enforceReal,enforcePositive,shuffleRows
,Int64(seed),sparseTrafo,iterations, constraintMask)
,Int64(seed),sparseTrafo,iterations, constraintMask, regMatrix)
end

"""
Expand Down Expand Up @@ -177,6 +185,11 @@ function solve(solver::Kaczmarz, u::Vector{T};
solverInfo != nothing && storeInfo(solverInfo,solver.cl,norm(solver.vl))
end

# backtransformation of solution with Tikhonov matrix
if solver.regMatrix != nothing
solver.cl = solver.cl .* (1 ./ sqrt.(solver.regMatrix))
end

return solver.cl
end

Expand Down

0 comments on commit 2eccca0

Please sign in to comment.