Matrix Multiplication Algorithm

Source: Internet
Author: User

Matrix Multiplication Algorithm
General Matrix Multiplication Algorithm: principle: the most important method of matrix multiplication is the product of General matrices. It is defined only when the number of columns in the first Matrix and the number of columns in the second matrix are the same. A single product refers to the product of a matrix. If A is m x n matrix and B is n x p Matrix, their product AB will be an m x p Matrix. The elements of the product matrix are obtained as follows: Code: struct mat {int n, m; double data [MAXN] [MAXN] ;}; int mul (mat & c, const mat & a, const mat & B) {int I, j, k; if (. m! = B. n) return 0; c. n =. n; c. m = B. m; for (I = 0; I <c. n; I ++) for (j = 0; j <c. m; j ++) for (c. data [I] [j] = k = 0; k <. m; k ++) c. data [I] [j] + =. data [I] [k] * B. data [k] [j]; return 1;} time complexity: O (n3) Traditional Partitioning Method: C = AB divides A, B, C into squares of equal size: c11 = A11B11 + A12B21 (2) C12 = A11B12 + A12B22 (3) C21 = A21B11 + A22B21 (4) C22 = A21B12 + A22B22 (5) if n = 2, then the product of the two second-order squares can be calculated directly using the formula (2)-(5). A total of eight multiplication and four addition are required. When the submatrix order is greater than 2, to obtain the product of the two submatrices, you can continue to block the submatrix until the submatrix order is reduced to 2. In this way, a recursive algorithm for sub-governance downgrading is generated. According to this algorithm, the product of the two n-order phalanx is converted into the product of the eight n/2-order phalanx and the addition of the four n/2-order phalanx. The addition of two n/2 × n/2 matrices can obviously be completed in c * n2/4 time. Here c is a constant. Therefore, T (n) is the time consumed for the calculation of the preceding divide and conquer law. T (2) = B T (n) = 8 T (n/2) + cn2 (n> 2) We can see from the above formula that the efficiency of the square matrix multiplication algorithm is not improved with the use of the divide and conquer method! The design philosophy of divide and conquer: divide a big problem that is hard to solve directly into the same problems of small scale, so as to break through and conquer each other separately. The computing time required for any problem that can be solved by a computer is related to its scale. The smaller the problem, the easier it is to solve it directly, and the less computing time it takes to solve the problem. For example, for n elements. When n = 2, you only need to make a comparison to sort the order. N = 3 timing problem. When n = 1, you do not need to perform any calculation as long as three comparisons ,.... When n is large, the problem is not easy to handle. It is sometimes quite difficult to solve a large-scale problem directly. The problems solved by the Division and control law generally have the following characteristics: 1. scalability. The problem can be easily solved by narrowing down to a certain extent; 2. sub-structures are the most important. The problem can be broken down into several small-scale identical problems; 3. conciseness. The subproblem resolved by the problem can be combined into the solution of the problem; 4. Independence. The subquestions to be resolved are independent of each other, that is, the subquestions do not contain public subquestions. The idea of division and recursion is like a pair of twins. They are often used in algorithm design at the same time, and thus produce efficient algorithms! The traditional method of Strassen requires 8 multiplication for the multiplication of two square arrays. According to the preceding method, we can see that to reduce the number of multiplication operations, the key lies in whether the multiplication operation can be performed less than eight times when the product of two second-order square arrays is calculated. Strassen proposed a new algorithm to calculate the product of two second-order square arrays. His algorithm only uses seven multiplication operations, but increases the number of addition and subtraction operations. The 7 multiplication is: M1 = A11 (B12-B22) M2 = (A11 + A12) B22 M3 = (A21 + A22) B11 M4 = A22 (B21-B11) m5 = (A11 + A22) (B11 + B22) M6 = (A12-A22) (B21 + B22) M7 = (A11-A21) (B11 + B12) after 7 multiplication, then perform several addition and subtraction: C11 = M5 + M4-M2 + M6 C12 = M1 + M2 C21 = M3 + M4 C22 = M5 + M1-M3-M7 algorithm efficiency analysis: Strassen Matrix Product grouping algorithm, seven recursive calls for the product of the n/2 matrix and 18 addition and subtraction operations on the n/2 matrix are used. The calculation time T (n) required by the algorithm satisfies the following recursive equation: T (2) = B T (n) = 7 T (n/2) + cn2 (n> 2) according to the formula used to solve the recursive equation, the solution is T (n) = O (nlog7) ≈ O (n2.81 ). It can be seen that the computing time complexity of Strassen matrix multiplication is better than that of normal matrix multiplication. Data Storage Structure: Two-dimensional array (Disadvantage: Must macro define A constant N) Double pointer (get rid of macro definition restrictions, but the program becomes cumbersome) original matrix: A, B, C (C = AB) block matrix: A11, A12, A21, and A22 (divide square matrix A into four parts) B11, B12, B21, and B22 (divide square matrix B into four parts) seven key elements: M1, M2, M3, M4, M5, M6, and M7 intermediate phalanx: AA, BB (A Bridge during addition) Level 2 pointer definition: int ** c; c = (int **) malloc (n * sizeof (int); for (int I = 0; I <n; I ++) c [I] = (int *) malloc (n * sizeof (int); problem of adding two matrices: // matrix addition: void add (int n, int A [] [N], int B [] [N], int R [] [N]) {int I, j; for (I = 0; I <n; I ++) For (j = 0; j <n; j ++) R [I] [j] = A [I] [j] + B [I] [j];} two matrix subtraction: // matrix subtraction: void sub (int n, int A [] [N], int B [] [N], int R [] [N]) {int I, j; for (I = 0; I <n; I ++) for (j = 0; j <n; j ++) R [I] [j] = A [I] [j]-B [I] [j];} the problem of multiplication of two phalanx: // multiply the square matrix with the Order of 2: void multiply (int A [] [N], int B [] [N], int R [] [N]) {int M1 = A [0] [1] * (B [0] [1]-B [1] [1]); int M2 = (A [0] [0] + A [0] [1]) * B [1] [1]; int M3 = (A [1] [0] + A [1] [1]) * B [0] [0]; int M4 = [1] [1] * (B [1] [0]-B [0] [0]); int M5 = (A [0] [0] + A [1] [1]) * (B [0] [0] + B [1] [1]); int M6 = (A [0] [1]-A [1] [1]) * (B [1] [0] + B [1] [1]); int M7 = (A [0] [0]-A [1] [0]) * (B [0] [0] + B [0] [1]); R [0] [0] = M5 + M4-M2 + M6; R [0] [1] = M1 + M2; R [1] [0] = M3 + M4; R [1] [1] = M5 + M1-M3-M7;} divide the large square matrix. // The void divide (int A [] [N] Function int n, int X11 [] [N], int X12 [] [N], int X21 [] [N], int X22 [] [N]) {int I, j; (I = 0; I <(n/2); I ++) for (j = 0; j <(n/2); j ++) {X11 [I] [j] = A [I] [j]; X12 [I] [j] = A [I] [j + (n/2)]; x21 [I] [j] = A [I + (n/2)] [j]; X22 [I] [j] = A [I + (n/2)] [j + (n/2)] ;}} complete code: # include <stdio. h> # include <math. h> # define N 4 // Second-Order Matrix Multiplication void Matrix_Multiply (int A [] [N], int B [] [N], int C [] [N]) {for (int I = 0; I <2; I ++) {for (int j = 0; j <2; j ++) {C [I] [j] = 0; for (int t = 0; t <2; t ++) {C [I] [j] = C [I] [j] + A [I] [t] * B [t] [j] ;}} // matrix addition: void add (int n, int A [] [N], int B [] [N], int R [] [N]) {int I, j; for (I = 0; I <n; I ++) for (j = 0; j <n; j ++) R [I] [j] = A [I] [j] + B [I] [j];} // matrix subtraction: void sub (int n, int A [] [N], int B [] [N], int R [] [N]) {int I, j; for (I = 0; I <n; I ++) for (j = 0; j <n; j ++) R [I] [j] = A [I] [j]-B [I] [j];} void strassen (int n, int A [] [N], int B [] [N], int C [] [N]) {in T I, j; int A11 [N] [N], A12 [N] [N], A21 [N] [N], A22 [N] [N]; int B11 [N] [N], B12 [N] [N], B21 [N] [N], B22 [N] [N]; int C11 [N] [N], C12 [N] [N], C21 [N] [N], C22 [N] [N]; int AA [N] [N], BB [N] [N]; int M1 [N] [N], M2 [N] [N], m3 [N] [N], M4 [N] [N], M5 [N] [N], M6 [N] [N], M7 [N] [N]; if (n = 2) {Matrix_Multiply (A, B, C) ;}else {for (I = 0; I <n/2; I ++) {for (j = 0; j <n/2; j ++) {A11 [I] [j] = A [I] [j]; a12 [I] [j] = A [I] [j + n /2]; A21 [I] [j] = A [I + n/2] [j]; a22 [I] [j] = A [I + n/2] [j + n/2]; B11 [I] [j] = B [I] [j]; b12 [I] [j] = B [I] [j + n/2]; B21 [I] [j] = B [I + n/2] [j]; b22 [I] [j] = B [I + n/2] [j + n/2];} sub (n/2, B12, B22, BB ); strassen (n/2, A11, BB, M1); add (n/2, A11, A12, AA); strassen (n/2, AA, B22, M2 ); add (n/2, A21, A22, AA); strassen (n/2, AA, B11, M3); sub (n/2, B21, B11, BB ); strassen (n/2, A22, BB, M4); add (n/2, A11, A22, AA); add (n/2, B11, B22, BB); strassen (n/2, AA, BB, M5); sub (n/2, A12, A22, AA); add (n/2, B21, B22, BB); strassen (n/2, AA, BB, M6); sub (n/2, A11, A21, AA); add (n/2, B11, B12, BB); strassen (n/2, AA, BB, M7); // C11 = M5 + M4-M2 + M6 add (n/2, M5, M4, AA); sub (n/2, M6, m2, BB); add (n/2, AA, BB, C11); // C12 = M1 + M2 add (n/2, M1, M2, C12 ); // C 21 = M3 + M4 add (n/2, M3, M4, C21); // C22 = M5 + M1-M3-M7 sub (n/2, M5, M3, AA); sub (n/2, M1, M7, BB); add (n/2, AA, BB, C22); for (I = 0; I <n/2; I ++) {for (j = 0; j <n/2; j ++) {C [I] [j] = C11 [I] [j]; C [I] [j + n/2] = C12 [I] [j]; C [I + n/2] [j] = C21 [I] [j]; C [I + n/2] [j + n/2] = C22 [I] [j] ;}}} int main (void) {int A [N] [N], B [N] [N], C [N] [N]; printf ("input A: \ n "); for (int I = 0; I <N; I ++) for (int j = 0; j <N; j ++) scanf ("% d ", & A [I] [j]); printf ("input B: \ n"); for (int I = 0; I <N; I ++) for (int j = 0; j <N; j ++) scanf ("% d", & B [I] [j]); strassen (N, A, B, c); printf ("C: \ n"); for (int I = 0; I <N; I ++) for (int j = 0; j <N; j ++) {printf ("% d", C [I] [j]); if (j> 0 & j % 3 = 0) printf ("\ n");} return 0;} Although the Strassen algorithm reduces the time cost from O (n3) to O (n2.81), further analysis, the coefficient of the latter is much larger than that of 1. When n is small (for example, n <45) and the moment This method is not applicable if there are few non-zero elements in the array. The Strassen algorithm is only valuable when n is large. Therefore, it can be said that this research is more focused on theoretical value.

Contact Us

The content source of this page is from Internet, which doesn't represent Alibaba Cloud's opinion; products and services mentioned on that page don't have any relationship with Alibaba Cloud. If the content of the page makes you feel confusing, please write us an email, we will handle the problem within 5 days after receiving your email.

If you find any instances of plagiarism from the community, please send an email to: info-contact@alibabacloud.com and provide relevant evidence. A staff member will contact you within 5 working days.

A Free Trial That Lets You Build Big!

Start building with 50+ products and up to 12 months usage for Elastic Compute Service

  • Sales Support

    1 on 1 presale consultation

  • After-Sales Support

    24/7 Technical Support 6 Free Tickets per Quarter Faster Response

  • Alibaba Cloud offers highly flexible support services tailored to meet your exact needs.