添加其他两种算法，并优化了效率

jacobi.jl

@@ -0,0 +1,35 @@
+function jacobi(A::SparseMatrixCSC,b::Vector;tol::Real=1e-2,maxIter::Int=1,out::Int=0,storeInterm::Bool=false;)
+	n = length(b)
+	OmInvD = 1./diag(A) #D矩阵 并乘以omega
+	x = zeros(eltype(A),n)
+	r = copy(b)
+	rnorm0 = norm(b)
+	res = 0
+	value = 0.0
+	flag   = -1
+	iter   = 0
+	dx = 0
+	for iter = 1:maxIter
+		println(iter)
+        dr = zeros(eltype(A),n)
+		for j = 1:n
+			dx =  r[j] * OmInvD[j]
+            x[j] += dx
+            for gidx = A.colptr[j] : A.colptr[j+1]-1
+                dr[A.rowval[gidx]] -= A.nzval[gidx] * dx
+            end
+		end
+        r += dr
+
+		res = norm(r) / rnorm0
+
+		println(res)
+
+		if res < tol
+			flag = 0; break
+		end
+	end
+
+	return x,flag,iter
+
+end

sor.jl

@@ -0,0 +1,59 @@
+
+"""
+x,flag,err,iter,resvec = sor(A::SparseMatrixCSC, b::Vector; omega::Real=2/3, tol::Real=1e-2, maxIter::Int=1,out::Int=0)
+
+Symmetric successive over-relaxation applied to the linear system A*x = b.
+
+!! Note that upper/lower triangular matrix are never build to save memory !!
+
+Input:
+
+	A       - sparse matrix CSC
+	b       - right hand side vector
+	omega   - relaxation parameter (omega=1.0: Gauss Seidel)
+	tol     - error tolerance (default=1e-2)
+	maxIter - maximum number of iterations (default=1)
+	out     - flag for output (-1 : no output, 0 : only errors, 1 : final status, 2: error at each iteration)
+
+Output:
+
+	x       - solution
+	flag    - exit flag (  0 : desired tolerance achieved,
+	                      -1 : maxIter reached without converging)
+	err     - error, i.e., norm(A*x-b)/norm(b)
+	iter    - number of iterations
+	resvec  - error at each iteration
+"""
+function sor(A::SparseMatrixCSC,b::Vector;tol::Real=1e-2,maxIter::Int=1,omega::Real=1.0,out::Int=0,storeInterm::Bool=false;)
+	n = length(b)
+	OmInvD = omega./diag(A) #D矩阵 并乘以omega
+	x = zeros(eltype(A),n)
+	r = copy(b)
+	rnorm0 = norm(b)
+	res = 0
+	value = 0.0
+	flag   = -1
+	iter   = 0
+	dx = 0
+	for iter = 1:maxIter
+		println(iter)
+		#第一次迭代 按照列来
+		for j = 1:n
+			dx =  r[j] * OmInvD[j]
+			for gidx = A.colptr[j] : A.colptr[j+1]-1
+				r[A.rowval[gidx]] -= A.nzval[gidx] * dx
+			end
+			x[j] += dx
+		end
+		res = norm(A*x-b) / rnorm0
+
+		println(res)
+
+		if res < tol
+			flag = 0; break
+		end
+	end
+
+	return x,flag,iter
+
+end

ssor.jl

@@ -29,48 +29,39 @@ function ssor(A::SparseMatrixCSC,b::Vector;tol::Real=1e-2,maxIter::Int=1,omega::
 	n = length(b)
 	OmInvD = omega./diag(A) #D矩阵 并乘以omega
 	x = zeros(eltype(A),n)
-	r = b - A * x
-	rnorm0 = norm(r)
+	r = copy(b)
+	rnorm0 = norm(b)
+	res = 0
 	value = 0.0
 	flag   = -1
 	iter   = 0
+	dx = 0
 	for iter = 1:maxIter
 		println(iter)
-		#第一次迭代
-		for i = 1:n
-			value = b[i]
-			for j = 1:n
-				for gidx = A.colptr[j] : A.colptr[j+1]-1
-					if A.rowval[gidx] == i
-						value -= A.nzval[gidx] * x[j]
-					end
-				end
+		#第一次迭代 按照列来
+		for j = 1:n
+			dx =  r[j] * OmInvD[j]
+			for gidx = A.colptr[j] : A.colptr[j+1]-1
+				r[A.rowval[gidx]] -= A.nzval[gidx] * dx
 			end
-			x[i]  = x[i] + value*OmInvD[i]
+			x[j] += dx
 		end
-
-		#第二次迭代
-		for i = n:-1:1
-			value = b[i]
-			for j = 1:n
-				for gidx = A.colptr[j] : A.colptr[j+1]-1
-					if A.rowval[gidx] == i
-						value -= A.nzval[gidx] * x[j]
-					end
-				end
+		for j = n:-1:1
+			dx =  r[j] * OmInvD[j]
+			for gidx = A.colptr[j] : A.colptr[j+1]-1
+				r[A.rowval[gidx]] -= A.nzval[gidx] * dx
 			end
-			x[i]  = x[i] + value*OmInvD[i]
+			x[j] += dx
 		end
-		r = b - A*x
-		println(norm(r) / rnorm0 )
+		res = norm(A*x-b) / rnorm0
+
+		println(res)
 
-		if norm(r) / rnorm0 < tol
+		if res < tol
 			flag = 0; break
 		end
 	end
 
-	println(r)
-	println(x)
-	return x,flag, norm(r) / rnorm0,iter
+	return x,flag,iter
 
 end