행렬 반전

Question

행렬을 읽고이를 항등 행렬로 인쇄하는 클래스 용 프로그램을 개발 중입니다. 그런 다음 단계별로 단위 행렬에 같은 계산을 적용하면서 행렬의 역행렬을 가져와 항등 행렬로 축소해야합니다.

코드 줄 부분이 줄었습니다. 각 열에 대해 하나씩 생성하는 대신 한 행에 모든 행렬을 사용하여 nxn 행렬을 항등 행렬로 줄이는 루프를 만들려고합니다. 열이 행별로 행렬을 줄이는 루프에 대한 도움이 필요합니다.

누구든지 아이디어가 있습니다. (나는 프로그래밍 초보자이고, C 언어를 사용하고 있으며, 게시 된 솔루션은 괜찮습니다. 그러나 모든 행을 한 번에 줄이는 알고리즘을 찾고 있습니다).

Answer 1

Wikipedia에는 Gaussian Elimination (멋진 설명문)이 있습니다 (설명해 주셨습니다).

가우스 소거법 의사 코드

i := 1 
j := 1 
while (i ≤ m and j ≤ n) do 
    Find pivot in column j, starting in row i: 
    maxi := i 
    for k := i+1 to m do 
    if abs(A[k,j]) > abs(A[maxi,j]) then 
     maxi := k 
    end if 
    end for 
    if A[maxi,j] ≠ 0 then 
    swap rows i and maxi, but do not change the value of i 
    Now A[i,j] will contain the old value of A[maxi,j]. 
    divide each entry in row i by A[i,j] 
    Now A[i,j] will have the value 1. 
    for u := i+1 to m do 
     subtract A[u,j] * row i from row u 
     Now A[u,j] will be 0, since A[u,j] - A[i,j] * A[u,j] = A[u,j] - 1 * A[u,j] = 0. 
    end for 
    i := i + 1 
    end if 
    j := j + 1 
end while 

Answer 2

당신은 무엇을해야하는지에 대한 아이디어에 대한 가우스 - 요르단 제거 http://en.wikipedia.org/wiki/Gauss%E2%80%93Jordan_elimination에 Wikipedia 기사에서 볼 수 있습니다. 이것은 고전적인 방법이며 이해한다면 문제의 알고리즘 부분을 해결할 수 있습니다. 거기에 몇 가지 예제 코드가 있습니다.

Answer 3

당신이 말하는 것은 가우스 요르단 제거 (Gauss Jordan elimination)입니다.

(nxn) 단위 행렬 옆에 (nxn) 행렬 A를 씁니다. 상수 C에 의해

1. 곱하기 행 전 = 0
2. 교환 행 i와 j는
3. C 시간 행을 추가 내가 J에게

을 행하기 :

세 가지 행 작업을 수행 할 수 있습니다!

반올림 오류가 발생하기 때문에 반전을 찾는 좋은 방법이 아닙니다. 피봇 팅하지 않으면 나쁜 대답으로 끝날 수 있습니다.

1. 첫 번째 행부터 시작하십시오. i 행은 A에 의해
2. 곱하기 i 행 [아래의 모든 행을 통해 1
3. 루프에 그 대각 요소 A를 제 i 행의 모든 ​​요소가 [I는 I]이 동일하게
4. 나누기 j, i]를 계산하고 j 번째 행에서 결과를 빼서 대각선 아래의 값을 0으로 만듭니다.
5. 대각선에있는 모든 행렬과 대각선 아래에있는 행렬을 가질 때까지 반복합니다.
6. 대각선 아래의 행에 대해 수행 한 것과 같은 방식으로 대각선 위의 값을 0과 같게하려면 유사한 행 조작을 수행하십시오.

원래 A가있는 식별 매트릭스가 표시되면 오른쪽 매트릭스는 역수와 같습니다. 이것이 사실임을 증명하기 위해 곱셈을하십시오.

코딩을 시작하기 전에 실례를 들어 손으로 직접 작성하십시오.일단 완성되면 코드가 작동 함을 증명하는 단위 테스트를 작성하십시오. 오래된 대학 프로젝트에서

Answer 4

약물이 최대 ... 아마 특허에 대해

typedef struct matrix 
    { 
    flag type;  /* type of matrix: {RECTANGULAR, SQUARE, BANDED } */ 
    ushort num_rows; /* number of rows in th matrix */ 
    union cols 
    { 
    ushort num_cols; /* Number of cols in Rectanular matrix */ 
    ushort band;  /* Bandwidth in a square, symmetric, banded matrix */ 
    } 
    double *m_val; /* ptr to start of storage of actual matrix values */ 
    } MAT; 

/* row_swap ---- swap two rows of an nxn UNBANDED matrix */ 
void row_swap(mat, r1, r2) 
    MAT *mat; 
    ushort r1, r2; 
    { 
    double *m = mat->m_val; 
    ushort n = m->num_row; 
    INDEX i; 
    double temp; 
    for (i=0; i < n; i++) 
     SWAP(*(m + n*r1 +i), *(m + n*r2 +i), temp); 
    } 

/* reduce --- function to do Gauss- Jordon Reduction on matrix */ 
static void reduce(m, in, ref) 
    MAT *m *in; 
    double *ref; 
    { 
    ushort n = mat->num_row; 
    INDEX r,c; 
    int sr,sc; 
    double key, mat = m->m_val, inv = in->m_val; 

    sr =(ref-mat)/n; 
    sc=(ref-mat) % n; 
    scalprod(1.0/(*ref), inv+n*sr, n); 
    scalprod(1.0/(*ref), ref-sc, n); 

    for (r=0; r LT n; r++) 
    { 
     if (r != sr && *(mat+n*r+sc) != 0.0) 
      { 
      key = *(mat+n*r+sc); 
      for (c=0; c < n; c++) 
     { 
      *(inv+n*r+c) -= key * *(inv+n*sr+c); 
      *(mat+n*r+c) -= key * *(mat+n*sr+c); 
     } 
      } 
     } 
    } 

/* find_non_zero -- this finds and returns a ptr to next non-zero entry 
* below and to the right of the last one */ 
static double *find_non_zero(start, last, n) 
    ushort n; 
    double *start, *last; 
    { 
    double *i = last;    /* i is ptr to matrix entries */ 
    while (*i == 0.0 && (i-start) < n*n) 
     { 
     if (n+i > start+n*n) /* last row without non-zero entry, */ 
     return (-1);  /* Matrix must be singular   */ 
     else 
     i += n;   /* this is not the last row, goto next row */ 
     } 
    if (i >= start + n*n) /* we have departed the matrix */ 
    return (0); 
    else     /* we found a non-zero entry, return it */ 
    return(i); 
    } 


/* invert -- function to invert nxn matrix */ 
flag invert(m, i) 
    MAT *m, *i; 
    { 
    double *ref, *new, *mat = m->m_val, *inv = i->m_val; 
    ushort n = mat->num_row; 
    INDEX i, j, row; 
    ushort new_row; 
    if (det(mat,n) == 0.0) return 0; 
    for (i=0; i < n; i++) 
    for (j=0; j < n; j++) 
     *(inv+n*i+j) = (i == j)? 1.0: 0.0; 
    ref = mat; 
    for (row=0; row < n; row++) 
     { 
     new = find_non_zero(mat,ref,n); 
     if (new == -1) 
     { 
     scr_print(" This matrix is singular, not invertible "); 
     break; 
     } 
     new_row = (new-mat)/n; 
     if (new_row != row) 
     { 
     row_swap(mat, new_row, row, n); 
     row_swap(inv, new_row, row, n); 
     } 
     reduce(mat,inv,ref,n); 
     ref += n+1; 
     } 
    } 

Answer 5

용서 내 앞에서 농담을 도와 드리겠습니다. ;-)
나는 성능이나 수치의 정확도를 향상시키기위한 몇 가지 (일반적으로 지역적인) 트릭을 선호하는 Wert (OP)에 깊은 인상을 심어 주기만했다. 행렬을 축소 행으로 변환하는 다중 단계 프로세스 에셜론은 꽤 많은 길을가는 것입니다.

는 사실이 때문에 C 언어 친숙의 영업 이익의 부족, (다양한 지역 최적화와, 오히려 거기에 따라 변화보다) 이 단순하고 기계적인 과정를 사용하는 것이 특히 중요하다. 먼저 걷고 나서 실행하면 ...

할당은 Gauss-Jordan algorithm의 구현 인 것 같습니다. 따라서 원본 매트릭스에서 필요한 모든 대체 작업을 수행 할뿐만 아니라 이것들은 원래 매트릭스의 오른쪽에 부착 된 항등 매트릭스로 시작합니다. 이 과정의 결과는 원본 행렬 ("증가 된 행렬"의 왼쪽 절반)이 항등 행렬로 변환 될 때 (단수가 아닌 것으로 가정 할 때), 우리가 신원으로 시작한 오른쪽에 메인 매트릭스의 역. 이 역행렬은 원래 행렬을 줄이기 위해 수행 된 개별 연산의 요약을 포함하므로 "여러 번 셀 수준 연산"과 관련하여 선형 대수 수준에서 "한 번에"수행 할 수 있습니다.)이 은 Wert의 "one shot"grails에 대한 혼란과 원점이 될 수 있습니다.

Answer 6

여기 울드 포트란 루틴의 :

http://www.pseudolog.net/gjelim.f

내가 나 자신 숙제 문제로이 있었다. n = 2m을 설정하고 원래 행렬을 항등 행렬로 증가시킵니다. C로 FORTRAN에서 번역하려면, 그 기억 :

IMPLICIT NONE ...-> WHILE 마십시오
을 위해>
가 ...- DO [C의 기본이 무시] ...-> 동안
CYCLE을 ...-> 계속

"원샷"방법에 대한 농담 : 매트릭스 반전은 일반적으로 O (N^3) 프로세스이며, 가장 빠른 루틴은 Strassen이 1969 년에 작성한 것으로 O (N^2.8). 행운을 빈다.

Answer 7

또한 Cholesky 분해를 시도하여 행렬의 역함수를 구할 수 있습니다. 이 링크에 파이썬에서 명확한 구현 :

http://adorio-research.org/wordpress/?p=4560