OI linear algebra - matrix
This chapter only discusses linear algebra in OI. For too advanced knowledge of linear algebra, please refer to other textbooks.
Matrix and vector
Matrix is a mathematical concept, which corresponds to a two-dimensional array in a computer. For example:
A = [ 1 3 5 1 2 4 ] A= \begin{bmatrix} 1 & 3 & 5 \\ 1 & 2 & 4 \end{bmatrix} A=[113254]
It's a 2 × 3 2 \times 3 two × 3 matrix, which is composed of two rows and three columns A [ i , j ] A[i,j] A[i,j] to specify the second in the matrix i i Line i j j Elements of column j, for example A [ 1 , 3 ] = 3 A[1,3]=3 A[1,3]=3.
There is a special kind of n × 1 n \times 1 n × 1 a matrix is called a vector, for example:
V = [ 1 2 3 ] V = \begin{bmatrix} 1 \\ 2 \\ 3 \end{bmatrix} V=⎣⎡123⎦⎤
Is a person who has 3 3 Three element matrix, this article will no longer introduce the vectors in computational geometry.
Square matrix is a kind of n × n n \times n n × Matrix of N, for example:
B = [ 1 3 1 2 ] B = \begin{bmatrix} 1 & 3 \\ 1 & 2 \end{bmatrix} B=[1132]
It's a 2 × 2 2 \times 2 two × In general, we study the properties of square matrices.
Matrix operation
Matrix basic operation
- Matrix addition: definition A + B A + B A+B is the addition of corresponding elements of two matrices of the same size.
- Matrix subtraction: Definitions A − B A - B A − B is the subtraction of corresponding elements of two matrices of the same size.
- Matrix multiplication: definition a A aA aA, matrix A A Every element of A is multiplied by a a a.
- Matrix transpose: definition A T A^T AT, put each A [ i , j ] A[i,j] A[i,j] and A [ j , i ] A[j,i] A[j,i] swap positions.
Matrix multiplication
definition A B AB AB, for A A A matrix a × n a \times n a × n times matrix B B B matrix n × b n \times b n × b. The size of the resulting matrix is a × b a \times b a×b.
Full definition:
A B [ i , j ] = ∑ k = 1 n A [ i , k ] × B [ k , j ] AB[i,j] = \sum_{k=1}^n A[i,k] \times B[k,j] AB[i,j]=k=1∑nA[i,k]×B[k,j]
Matrix multiplication satisfies the associative law, but not the commutative law.
Define the identity matrix as a n × n n \times n n × N, whose diagonal elements are 1 1 1. Other elements are 0 0 0 Multiplying a matrix by the identity matrix does not change the matrix.
Power of matrix
Define the power of the matrix as A k = A ⋅ A ... A A^k=A \cdot A \ldots A Ak=A ⋅ a... A, i.e k k k A A A matrix multiplication.
The naive algorithm of matrix multiplication needs O ( k n 3 ) O(kn^3) O(kn3) time can be achieved after matrix fast power optimization O ( n 3 log k ) O(n^3 \log k) O(n3logk).
Similarly to other fast powers, we also use binary fast powers for optimization.
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
#define FR freopen("in.txt", "r", stdin)
#define FW freopen("out.txt", "w", stdout)
int n;
ll k;
ll mod = 1000000007;
struct Matrix
{
ll arr[105][105];
Matrix()
{
for (int i = 1; i <= n; i++)
{
for (int j = 1; j <= n; j++)
{
arr[i][j] = 0;
}
}
}
Matrix operator*(Matrix &o)
{
Matrix ans;
for (int i = 1; i <= n; i++)
{
for (int j = 1; j <= n; j++)
{
for (int k = 1; k <= n; k++)
{
ans.arr[i][j] = (ans.arr[i][j] + (arr[i][k] * o.arr[k][j]) % mod) % mod;
}
}
}
return ans;
}
Matrix fpow(ll e)
{
Matrix temp = *this;
Matrix ans;
for (int i = 1; i <= n; i++)
ans.arr[i][i] = 1;
for (; e; e >>= 1, temp = temp * temp)
{
if (e & 1)
ans = ans * temp;
}
return ans;
}
};
int main()
{
scanf("%d %lld", &n, &k);
Matrix m;
for (int i = 1; i <= n; i++)
{
for (int j = 1; j <= n; j++)
{
scanf("%lld", &m.arr[i][j]);
}
}
m = m.fpow(k);
for (int i = 1; i <= n; i++)
{
for (int j = 1; j <= n; j++)
{
printf("%lld ", m.arr[i][j]);
}
printf("\n");
}
return 0;
}
Linear recurrence
A linear recurrence equation is a sequence of numbers with some initial values f ( 0 ) , f ( 1 ) , ... , f ( k − 1 ) f(0),f(1),\ldots,f(k-1) f(0),f(1),..., f(k − 1) and a recursive formula:
f ( n ) = c 1 f ( n − 1 ) + c 2 f ( n − 2 ) + ... + c k f ( n − k ) f(n) = c_1f(n-1) + c_2f(n-2) + \ldots + c_kf(n-k) f(n)=c1f(n−1)+c2f(n−2)+...+ckf(n−k)
here c i c_i ci ， are constants. The time complexity of ordinary recursive method is O ( k n ) O(kn) O(kn), and the time of using matrix fast power recursive method is O ( k 3 log n ) O(k^3 \log n) O(k3logn).
Fibonacci sequence recurrence
The Fibonacci sequence has the following recurrence formula:
f ( 0 ) = 0 f ( 1 ) = 1 f ( n ) = f ( n − 1 ) + f ( n − 2 ) f(0) = 0 \\ f(1) = 1 \\ f(n) = f(n-1) + f(n-2) f(0)=0f(1)=1f(n)=f(n−1)+f(n−2)
We define the sequence vector as:
V k = [ f ( k ) f ( k − 1 ) ] V_k = \begin{bmatrix} f(k) \\ f(k-1) \end{bmatrix} Vk=[f(k)f(k−1)]
Define the recurrence matrix as:
V k + 1 = T V k V_{k+1} = TV_k Vk+1=TVk
The process of linear recurrence using matrix needs to construct T T T matrix. Through reasonable construction, we have:
T = [ 1 1 1 0 ] T = \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} T=[1110]
Therefore:
V n = T n [ 1 0 ] V_n = T^n \begin{bmatrix} 1 \\ 0 \end{bmatrix} Vn=Tn[10]
T n T^n Tn can be calculated in logarithmic time.
More general recurrence
Similar to the above definition, for the formula:
f ( n ) = c 1 f ( n − 1 ) + c 2 f ( n − 2 ) + ... + c k f ( n − k ) f(n) = c_1f(n-1) + c_2f(n-2) + \ldots + c_kf(n-k) f(n)=c1f(n−1)+c2f(n−2)+...+ckf(n−k)
We can construct a more general T T The T matrix is:
T = [ c 1 c 2 ⋯ c k 1 ⋯ 0 0 ⋮ ⋱ ⋮ ⋮ 0 ⋯ 1 0 ] T = \begin{bmatrix} c_1 & c_2 & \cdots & c_k \\ 1 & \cdots & 0 & 0 \\ \vdots & \ddots & \vdots & \vdots \\ 0 & \cdots & 1 & 0 \\ \end{bmatrix} T=⎣⎢⎢⎢⎡c11⋮0c2⋯⋱⋯⋯0⋮1ck0⋮0⎦⎥⎥⎥⎤
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
#define FR freopen("in.txt", "r", stdin)
#define FW freopen("out.txt", "w", stdout)
int n;
ll mod = 1000000007;
struct Matrix
{
ll arr[4][4];
Matrix()
{
for (int i = 1; i <= n; i++)
{
for (int j = 1; j <= n; j++)
{
arr[i][j] = 0;
}
}
}
Matrix operator*(const Matrix &o)
{
Matrix ans;
for (int i = 1; i <= n; i++)
{
for (int j = 1; j <= n; j++)
{
for (int k = 1; k <= n; k++)
{
ans.arr[i][j] = (ans.arr[i][j] + (arr[i][k] * o.arr[k][j]) % mod) % mod;
}
}
}
return ans;
}
Matrix fpow(ll e)
{
Matrix temp = *this;
Matrix ans;
for (int i = 1; i <= n; i++)
ans.arr[i][i] = 1;
for (; e; e >>= 1, temp = temp * temp)
{
if (e & 1)
ans = ans * temp;
}
return ans;
}
};
struct Query
{
ll k;
ll id;
bool operator<(const Query &o) const
{
return k < o.k;
}
} ques[120];
ll ANS[120];
int main()
{
n = 3;
int T;
scanf("%d", &T);
for (int i = 0; i < T; i++)
{
scanf("%lld", &ques[i].k);
ques[i].id = i;
}
sort(ques, ques + T);
Matrix curr;
ll currE = 0;
curr.arr[1][1] = 1;
curr.arr[2][2] = 1;
curr.arr[3][3] = 1;
for (int i = 0; i < T; i++)
{
ll k = ques[i].k;
if (k <= 3)
{
ANS[ques[i].id] = 1;
continue;
}
Matrix m;
m.arr[1][1] = 1;
m.arr[1][3] = 1;
m.arr[2][1] = 1;
m.arr[3][2] = 1;
curr = curr * m.fpow(k - 3 - currE);
currE = k - 3;
ANS[ques[i].id] = (curr.arr[1][1] + (curr.arr[1][2] + curr.arr[1][3]) % mod) % mod;
}
for (int i = 0; i < T; i++)
{
printf("%lld\n", ANS[i]);
}
return 0;
}
Gauss elimination and linear equation, determinant and inverse of matrix
Determinant and inverse of matrix
First, we define the adjoint factor of the matrix, also known as the cofactor. I am:
C [ i , j ] = ( − 1 ) i + j det ( M [ i , j ] ) C[i,j] = (-1)^{i+j}\det(M[i,j]) C[i,j]=(−1)i+jdet(M[i,j])
among M [ i , j ] M[i,j] M[i,j] is the primitive matrix A A A remove i i i line and j j The remaining matrix after column j.
Then the definition of determinant (expanded according to the first row) is:
det ( A ) = ∑ j = 1 n A [ 1 , j ] C [ 1 , j ] \det(A) = \sum_{j = 1}^n A[1,j] C[1,j] det(A)=j=1∑nA[1,j]C[1,j]
You can also expand by a row or a column.
We define the matrix A A The inverse matrix of A is:
A A − 1 = I AA^{-1} = I AA−1=I
also
A − 1 [ i , j ] = C [ j , i ] det A A^{-1}[i,j] = \frac{C[j,i]}{\det A} A−1[i,j]=detAC[j,i]
A square matrix has an inverse matrix if and only if this matrix det A ≠ 0 \det A \neq 0 detA=0.
Gauss elimination and linear equations
Every system of linear equations can be written as a matrix equation A x ⃗ = b ⃗ A\vec{x}=\vec{b} Ax =b , where A A A is called the coefficient matrix, x ⃗ \vec{x} x Is the solution vector, b ⃗ \vec{b} b Is a constant vector.
The core of Gaussian elimination theory is as follows:
- The two equations are interchangeable and the solution remains unchanged;
- Equation multiplied by nonzero number k k k. Solution invariant;
- Equation multiplication k k k plus another equation, the solution remains unchanged.
After the Gaussian elimination method transforms the augmented matrix into the simplest form, it needs to master the linear related knowledge for the assignment of free unknowns, and the assignment has some difficulties in the learning process due to the factors of manual experience. The Gaussian elimination method is divided into five steps, and a five step method is proposed. The contents are as follows:
- The elementary row of the augmented matrix is transformed into the simplest row;
- Reductive linear equations;
- Solve the first variable;
- Supplement free unknowns;
- The column represents the general solution of the system of equations.
P3389 weak Gaussian elimination
P2455 strong Gaussian elimination
#include <bits/stdc++.h>
using namespace std;
#define FR freopen("in.txt", "r", stdin)
#define FW freopen("out.txt", "w", stdout)
typedef long long ll;
double matrix[105][105];
double esp = 1e-3;
int main()
{
int n;
scanf("%d", &n);
for (int i = 0; i < n; i++)
{
for (int j = 0; j <= n; j++)
{
scanf("%lf", &matrix[i][j]);
}
}
for (int r = 0; r < n; r++)
{
int domain = r;
// Select the row with the largest coefficient as the principal element to reduce the loss of accuracy
for (int j = r; j < n; j++)
{
if (abs(matrix[j][r]) > abs(matrix[domain][r]))
domain = j;
}
// No primary element, skip
if (abs(matrix[domain][r]) < esp)
{
continue;
}
// Swap two lines
for (int j = 0; j <= n; j++)
{
swap(matrix[domain][j], matrix[r][j]);
}
// Elimination element
for (int l = 0; l < n; l++)
{
if (l != r)
for (int j = n; j >= r; j--)
{
matrix[l][j] -= matrix[r][j] / matrix[r][r] * matrix[l][r];
}
}
}
// Check the existence of solutions
for (int r = 0; r < n; r++)
{
// Find the principal element of this line
int ptr = r;
while (abs(matrix[r][ptr]) < esp && ptr <= n)
{
ptr++;
}
// If the principal component is on the constant column, then there is no solution
if (ptr == n)
{
printf("-1");
return 0;
}
}
// If there is a lack of principal elements, there are infinite solutions
for (int r = 0; r < n; r++)
{
if (abs(matrix[r][r]) < esp)
{
printf("0");
return 0;
}
}
// Otherwise, there is a unique solution
for (int r = 0; r < n; r++)
{
// Divided by principal element
printf("x%d=%.2lf\n", r + 1, matrix[r][n] / matrix[r][r]);
}
return 0;
}
Linear equations can be solved by solving the inverse of the matrix x ⃗ = A − 1 b ⃗ \vec{x}=A^{-1}\vec{b} x =A−1b .
Specific application
- Computational geometry: calculate the relationship between geometry through vector, cross multiplication and other operations
- Probability theory: Markov chain, counting, matrix fast power optimization, etc
- Polynomial: the basis of polynomial algorithm