function GaussMatrix(A, B)
aug_matrix = hcat(A, B)
(m, n) = size(A)
for i = 1:n - 1
max_row = argmax(abs.(aug_matrix[i:m, i]), dims = 1)[1] + i - 1
if max_row!= i
aug_matrix[[i, max_row], :] = aug_matrix[[max_row, i], :]
end
pivot = aug_matrix[i, i]
if pivot == 0
println("不能用高斯消去")
return zeros(Float64, n)
end
multiplier_matrix = Matrix{Float64}(I, m, m)
for j = i + 1:m
multiplier_matrix[j, i] = -aug_matrix[j, i] / pivot
end
aug_matrix = multiplier_matrix * aug_matrix
end
x = zeros(Float64, n)
x[n] = aug_matrix[n, end] / aug_matrix[n, n]
for k = (n - 1):-1:1
s = dot(aug_matrix[k, (k + 1):n], x[(k + 1):n])
x[k] = (aug_matrix[k, end] - s) / aug_matrix[k, k]
end
return x
end
function JacobiMatrix(A, B, tol = 1e-8)
(m, n) = size(A)
if m!= n
error("系数矩阵A必须是方阵")
end
if size(B, 2)!= 1
B = reshape(B, n, 1)
elseif size(B, 1)!= n
error("常数项向量B的维度与系数矩阵A的行数不匹配")
end
row_sums = sum(abs.(A), dims = 2) - diagm(diag(abs.(A)))
x = zeros(Float64, n, 1)
x = B./ (diag(A) + row_sums[:, 1])
max_iterations = 1000
iteration_count = 0
while true
x_new = zeros(Float64, n, 1)
Ax = A * x[:, 1]
diagAx = diagm(diag(A)) * x[:, 1]
diff = B - (Ax - diagAx)
x_new = (diff / diag(A))[:, newaxis]
error = norm(x_new - x)
if error < tol || iteration_count >= max_iterations
return x_new[:, 1]
end
x = x_new
iteration_count += 1
end
end
A = [10 -1 2 0; -1 11 -1 3; 2 -1 10 -1; 0 3 -1 8]
B = [6; 25; -11; 15]
solution_gauss_matrix = GaussMatrix(A, B)
println("The solution of Gauss (matrix version) is")
for i in 1:length(solution_gauss_matrix)
println("x$i = ", solution_gauss_matrix[i])
end
solution_jacobi_matrix = JacobiMatrix(A, B)
println("The solution of Jacobi (matrix version) is")
for i in 1:length(solution_jacobi_matrix)
println("x$i = ", solution_jacobi_matrix[i])
end
Jacobi迭代Julia脚本报错ERROR: DimensionMismatch: dimensions must match: a has dims (Base.OneTo(4), Base.OneTo(1)), b has dims (Base.OneTo(4), Base.OneTo(4)), mismatch at 2,这个怎么修改代码?
专栏
标签
高校专区
Jacobi迭代
一般问题
发布于 2024-12-06 14:50:31
查看 29过去569天
所属专栏:Julia语言
产品信息:Syslab科学计算环境
全部回答 1
发布于 2024-12-07 09:50:38
您好,您的第38行代码,减号左边的维度是(4,1),减号右边的维度是(4,4),维度不同,无法直接相减,需要使用.-的形式
和原帖交流更多问题细节吧,去