Последовательное расслабление (Python) - PullRequest
Новая
       13

Последовательное расслабление (Python)

0 голосов
/ 11 ноября 2018

Здесь у меня есть некоторый скрипт на python, который решает систему линейных уравнений с использованием метода Гаусса-Зейделя: 

import numpy as np

ITERATION_LIMIT = 1000

#system
A = np.array([[15., -4., -3., 8.],
          [-4., 10., -4., 2.],
          [-3., -4., 10., 2.],
          [8., 2., 2., 12.]
          ])
# vector b
b = np.array([2., -12., -4., 6.])

print("System of equations:")
for i in range(A.shape[0]):
    row = ["{0:3g}*x{1}".format(A[i, j], j + 1) for j in range(A.shape[1])]
    print("[{0}] = [{1:3g}]".format(" + ".join(row), b[i]))

x = np.zeros_like(b)


for it_count in range(1, ITERATION_LIMIT):
    x_new = np.zeros_like(x)
    print("Iteration {0}: {1}".format(it_count, x))
    for i in range(A.shape[0]):
        s1 = np.dot(A[i, :i], x_new[:i])
        s2 = np.dot(A[i, i + 1:], x[i + 1:])
        x_new[i] = (b[i] - s1 - s2) / A[i, i]
    if np.allclose(x, x_new, rtol=1e-8):
        break
    x = x_new

Что он выводит: 

Iteration 379: [-21.36409652 -22.09743    -19.9999946   21.75896845]
Iteration 380: [-21.36409676 -22.09743023 -19.99999481  21.75896868]
Iteration 381: [-21.36409698 -22.09743045 -19.99999501  21.7589689 ]

Моя задача состоит в том, чтобысделайте из этого метод Successive Over Relaxation (SOR), который использует значения омега для уменьшения количества итераций.Если omega = 1, он становится методом Гаусса-Зейделя, if < 1 - методом простых итераций, > 1 и < 2 - SOR.Очевидно, что при более высоких значениях омега число итераций должно уменьшаться.Вот алгоритм, который предлагает Википедия: 

Inputs: A, b, omega
Output: phi (roots for linear equations)

Choose an initial guess phi to the solution
repeat until convergence
  for i from 1 until n do
    sigma <- 0
    for j from 1 until n do
      if j ≠ i then
        sigma <- sigma + A[i][j]*phi[j]
      end if
    end (j-loop)
    phi[i] = phi[i] + omega*((b[i] - sigma)/A[i][i]) - phi[i]
  end (i-loop)
  check if convergence is reached
end (repeat)

Есть ли у кого-нибудь работающий алгоритм на python?Было бы очень хорошо, если бы вы могли сделать несколько комментариев к коду или помочь мне, как изменить этот код.Спасибо!

...