Skip to content

Commit b636e5f

Browse files
committed
Write simpler Gauss Seidel smoother
1 parent 089fede commit b636e5f

File tree

3 files changed

+29
-6
lines changed

3 files changed

+29
-6
lines changed

src/multilevel.jl

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -16,8 +16,8 @@ abstract type CoarseSolver end
1616
struct Pinv <: CoarseSolver
1717
end
1818

19-
MultiLevel{Ti,Tv}(l::Vector{Level{Ti,Tv}}, A::SparseMatrixCSC{Ti,Tv}, presmoother, postsmoother; coarse_solver = Pinv()) =
20-
MultiLevel(l, A, coarse_solver, presmoother, postsmoother)
19+
MultiLevel{Ti,Tv}(l::Vector{Level{Ti,Tv}}, A::SparseMatrixCSC{Ti,Tv}, presmoother, postsmoother) =
20+
MultiLevel(l, A, Pinv(), presmoother, postsmoother)
2121
Base.length(ml) = length(ml.levels) + 1
2222

2323
function Base.show(io::IO, ml::MultiLevel)

src/preconditioner.jl

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -6,8 +6,8 @@ end
66

77
aspreconditioner(ml::MultiLevel) = Preconditioner(ml)
88

9-
\(p::Preconditioner, b) = p * b
10-
*(p::Preconditioner, b) = solve(p.ml, b, 1, V(), 1e-12)
9+
\(p::Preconditioner, b) = solve(p.ml, b, 1, V(), 1e-12)
10+
*(p::Preconditioner, b) = p.ml.levels[1].A * x
1111

1212
A_ldiv_B!(x, p::Preconditioner, b) = copy!(x, p \ b)
1313
A_mul_B!(b, p::Preconditioner, x) = A_mul_B!(b, p.ml.levels[1].A, x)

src/smoother.jl

Lines changed: 25 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -7,10 +7,33 @@ end
77
struct GaussSeidel{S} <: Smoother
88
sweep::S
99
end
10-
GaussSeidel(;sweep = ForwardSweep()) = GaussSeidel(sweep)
10+
GaussSeidel() = GaussSeidel(ForwardSweep())
1111

1212
presmoother!(s, A, x, b) = smoother(s, s.sweep, A, x, b)
1313
postsmoother!(s, A, x, b) = smoother(s, s.sweep, A, x, b)
1414

1515
smoother(s::GaussSeidel, ::ForwardSweep, A, x, b) =
16-
gauss_seidel!(x, A, b, maxiter = 1)
16+
gs!(A, b, x)
17+
#gauss_seidel!(x, A, b, maxiter = 1)
18+
19+
function gs!{T,Ti}(A::SparseMatrixCSC{T,Ti}, b::Vector{T}, x::Vector{T})
20+
n = size(A, 1)
21+
z = zero(eltype(A))
22+
for i = 1:n
23+
# rsum = calc_weighted_sum(A, x)
24+
rsum = z
25+
diag = z
26+
for j in nzrange(A, i)
27+
row = A.rowval[j]
28+
val = A.nzval[j]
29+
if i == row
30+
diag = val
31+
else
32+
rsum += val * x[row]
33+
end
34+
end
35+
if diag != 0
36+
x[i] = (b[i] - rsum) / diag
37+
end
38+
end
39+
end

0 commit comments

Comments
 (0)