генерация случайных матриц в питоне - PullRequest
0 голосов
/ 17 октября 2018

В следующем коде я реализовал исключение Гаусса с частичным поворотом для общей квадратной линейной системы Ax = b.Я проверил мой код, и он дал правильный вывод.Однако сейчас я пытаюсь сделать следующее, но я не совсем уверен, как его кодировать, ища некоторую помощь в этом!

Я хочу проверить мою реализацию, решив Ax = b, где A - случайная матрица 100x100, а b - случайный вектор 100x1.

В моем коде я поместилв матрицах

A = np.array([[3.,2.,-4.],[2.,3.,3.],[5.,-3.,1.]])

b =  np.array([[3.],[15.],[14.]])

и получил следующий правильный вывод:

[3. 1. 2.]
[3. 1. 2.]

но как теперь изменить его для генерации случайных матриц?

вот мой код ниже:

import numpy as np
def GEPP(A, b, doPricing = True):
    '''
    Gaussian elimination with partial pivoting.
    input: A is an n x n numpy matrix
           b is an n x 1 numpy array
    output: x is the solution of Ax=b 
        with the entries permuted in 
        accordance with the pivoting 
        done by the algorithm
    post-condition: A and b have been modified.
    '''
    n = len(A)
    if b.size != n:

        raise ValueError("Invalid argument: incompatible sizes between"+

                     "A & b.", b.size, n)

    # k represents the current pivot row. Since GE traverses the matrix in the 

    # upper right triangle, we also use k for indicating the k-th diagonal 

    # column index.

    # Elimination

    for k in range(n-1):

        if doPricing:

            # Pivot

            maxindex = abs(A[k:,k]).argmax() + k

            if A[maxindex, k] == 0:


                raise ValueError("Matrix is singular.")

            # Swap

            if maxindex != k:

                A[[k,maxindex]] = A[[maxindex, k]]

                b[[k,maxindex]] = b[[maxindex, k]]

        else:

            if A[k, k] == 0:

                raise ValueError("Pivot element is zero. Try setting doPricing to True.")

       #Eliminate

       for row in range(k+1, n):

           multiplier = A[row,k]/A[k,k]

           A[row, k:] = A[row, k:] - multiplier*A[k, k:]

           b[row] = b[row] - multiplier*b[k]

    # Back Substitution

    x = np.zeros(n)

    for k in range(n-1, -1, -1):

        x[k] = (b[k] - np.dot(A[k,k+1:],x[k+1:]))/A[k,k]

    return x



if __name__ == "__main__":

    A = np.array([[3.,2.,-4.],[2.,3.,3.],[5.,-3.,1.]])

    b =  np.array([[3.],[15.],[14.]])

    print (GEPP(np.copy(A), np.copy(b), doPricing = False))

    print (GEPP(A,b))

Ответы [ 5 ]

0 голосов
/ 17 октября 2018

Поэтому я бы использовал для этого np.random.randint .

numpy.random.randint(low, high=None, size=None, dtype='l')

, который выводит массив случайных целых чисел в форме размера из соответствующего распределения или один такой случайныйint, если размер не указан.

low - это нижняя граница нужных вам интервалов в вашем диапазоне

high на единицу больше, чем верхняя граница в требуемом диапазоне

size - это размеры вашего выходного массива

dtype - это dtype результата

, поэтому на вашем месте я написал бы

A = np.random.randint(0, 11, (100, 100))
b = np.random.randint(0, 11, 100)
0 голосов
/ 17 октября 2018

Вы можете использовать собственную функцию rand numpy:

np.random.rand()

В вашем коде просто определите A и b как:

A = np.random.rand(100, 100)
b = np.random.rand(100)

Это сгенерирует матрицу 100x100 и вектор 100x1 (оба - numpy)массивы), заполненные случайными значениями от 0 до 1.

См. документы для этой функции, чтобы узнать больше.

0 голосов
/ 17 октября 2018

Вы уже используете NumPy.Рассматривали ли вы np.random.rand?

np.random.rand(m, n), вы получите случайную матрицу со значениями в [0, 1).В дальнейшем вы можете обработать его, умножив случайные значения или округлив.

РЕДАКТИРОВАТЬ: Примерно так

if __name__ == "__main__":
    A = np.round(np.random.rand(100, 100)*10)
    b =  np.round(np.random.rand(100)*10)
    print (GEPP(np.copy(A), np.copy(b), doPricing = False))
    print (GEPP(A,b))
0 голосов
/ 17 октября 2018
A = np.random.normal(size=(100,100))
b = np.random.normal(size=(100,1))
x = np.linalg.solve(A,b)
assert max(abs(A@x - b)) < 1e-12

Очевидно, что вы можете использовать дистрибутивы, отличные от normal, например uniform.

0 голосов
/ 17 октября 2018

По сути, вы можете создать нужные матрицы с единицами, а затем выполнить итерацию по ним, установив, например, каждое значение на random.randint(0,100).

Пустая матрица с единицами:

one_array = np.ones((100, 100))

РЕДАКТИРОВАТЬ:

как:

for x in one_array.shape[0]:
    for y in one_array.shape[1]:
        one_array[x][y] = random.randint(0, 100)
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...