ガウスの消去法を実装する
code:gauss_elimination.jl
function gauss_elimination!(A, b; do_pivot = true, check = true)
m, n = size(A)
p = length(b)
if check && (m != n || m != p)
error("Invalid matrix or vector size")
end
A, b = copy(float(A)), copy(float(b)) # 直接変更しないようにコピーする
tri_A = forward_elimination!(aug_A, n; do_pivot = do_pivot, check = check)
x = backward_substitution(tri_A, n)
return x
end
function forward_elimination!(A, n; do_pivot = true, check = true)
for diag = 1:n-1
if do_pivot
pivoting!(A, diag, n; check = check)
end
for row = diag+1:n
for col = diag:n+1
end
end
end
return A
end
function backward_substitution(A, n)
x = zeros(n)
for row = n:-1:1 # 後ろから順に解を求める
for col = 1:row-1
end
end
return x
end
# ピボット選択
function pivoting!(A, diag, n; check = true)
error("Matrix is singular or nearly singular")
end
if max_row != diag
end
end
正しく動いているかを見る
code:test.jl
gauss_elimination!(A, b)
ピボット選択の効果を見る
code:test.jl
gauss_elimination!(A, b; do_pivot = false)
# gauss_elimination!(A, b)
まともな実装と実行時間を比較する
全部関数内に含んだらもうちょっとマシになる?あんも.icon
code:test.jl
using BenchmarkTools
n = 1000
A = rand(n, n)
b = rand(n)
@btime gauss_elimination!($A, $b) # 354.130 ms (4989 allocations: 38.71 MiB)
@btime $A \ $b # 4.940 ms (4 allocations: 7.64 MiB)
@btime inv($A)*$b # 11.313 ms (6 allocations: 9.35 MiB)