本文博客链接:http://blog.csdn.net/jdh99,作者:jdh
环境:
主机:XP
开发环境:mingw
功能:
在C++中实现矩阵运算
说明:
将这篇文章(http://blog.csdn.net/jdh99/article/details/7360091)中在C++下实现的矩阵运算移植到C下,以便在单片机中运算.
源代码:
matrix.h:
- #ifndef _MATRIX_H
- #define _MATRIX_H
- //头文件
- #include <stdio.h>
- #include <stdlib.h>
- //矩阵数据结构
- //二维矩阵
- struct _Matrix
- {
- int m;
- int n;
- float *arr;
- };
- //矩阵方法
- //设置m
- void matrix_set_m(struct _Matrix *m,int mm);
- //设置n
- void matrix_set_n(struct _Matrix *m,int nn);
- //初始化
- void matrix_init(struct _Matrix *m);
- //释放
- void matrix_free(struct _Matrix *m);
- //读取i,j坐标的数据
- //失败返回-31415,成功返回值
- float matrix_read(struct _Matrix *m,int i,int j);
- //写入i,j坐标的数据
- //失败返回-1,成功返回1
- int matrix_write(struct _Matrix *m,int i,int j,float val);
- //矩阵运算
- //成功返回1,失败返回-1
- int matrix_add(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C);
- //C = A - B
- //成功返回1,失败返回-1
- int matrix_subtract(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C);
- //C = A * B
- //成功返回1,失败返回-1
- int matrix_multiply(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C);
- //行列式的值,只能计算2 * 2,3 * 3
- //失败返回-31415,成功返回值
- float matrix_det(struct _Matrix *A);
- //求转置矩阵,B = AT
- //成功返回1,失败返回-1
- int matrix_transpos(struct _Matrix *A,struct _Matrix *B);
- //求逆矩阵,B = A^(-1)
- //成功返回1,失败返回-1
- int matrix_inverse(struct _Matrix *A,struct _Matrix *B);
- #endif
#ifndef _MATRIX_H #define _MATRIX_H //头文件 #include <stdio.h> #include <stdlib.h> //矩阵数据结构 //二维矩阵 struct _Matrix { int m; int n; float *arr; }; //矩阵方法 //设置m void matrix_set_m(struct _Matrix *m,int mm); //设置n void matrix_set_n(struct _Matrix *m,int nn); //初始化 void matrix_init(struct _Matrix *m); //释放 void matrix_free(struct _Matrix *m); //读取i,j坐标的数据 //失败返回-31415,成功返回值 float matrix_read(struct _Matrix *m,int i,int j); //写入i,j坐标的数据 //失败返回-1,成功返回1 int matrix_write(struct _Matrix *m,int i,int j,float val); //矩阵运算 //成功返回1,失败返回-1 int matrix_add(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C); //C = A - B //成功返回1,失败返回-1 int matrix_subtract(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C); //C = A * B //成功返回1,失败返回-1 int matrix_multiply(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C); //行列式的值,只能计算2 * 2,3 * 3 //失败返回-31415,成功返回值 float matrix_det(struct _Matrix *A); //求转置矩阵,B = AT //成功返回1,失败返回-1 int matrix_transpos(struct _Matrix *A,struct _Matrix *B); //求逆矩阵,B = A^(-1) //成功返回1,失败返回-1 int matrix_inverse(struct _Matrix *A,struct _Matrix *B); #endif
matrix.c:
- #include "matrix.h"
- //矩阵方法
- //设置m
- void matrix_set_m(struct _Matrix *m,int mm)
- {
- m->m = mm;
- }
- //设置n
- void matrix_set_n(struct _Matrix *m,int nn)
- {
- m->n = nn;
- }
- //初始化
- void matrix_init(struct _Matrix *m)
- {
- m->arr = (float *)malloc(m->m * m->n * sizeof(float));
- }
- //释放
- void matrix_free(struct _Matrix *m)
- {
- free(m->arr);
- }
- //读取i,j坐标的数据
- //失败返回-31415,成功返回值
- float matrix_read(struct _Matrix *m,int i,int j)
- {
- if (i >= m->m || j >= m->n)
- {
- return -31415;
- }
- return *(m->arr + i * m->n + j);
- }
- //写入i,j坐标的数据
- //失败返回-1,成功返回1
- int matrix_write(struct _Matrix *m,int i,int j,float val)
- {
- if (i >= m->m || j >= m->n)
- {
- return -1;
- }
- *(m->arr + i * m->n + j) = val;
- return 1;
- }
- //矩阵运算
- //成功返回1,失败返回-1
- int matrix_add(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C)
- {
- int i = 0;
- int j = 0;
- //判断是否可以运算
- if (A->m != B->m || A->n != B->n || \
- A->m != C->m || A->n != C->n)
- {
- return -1;
- }
- //运算
- for (i = 0;i < C->m;i++)
- {
- for (j = 0;j < C->n;j++)
- {
- matrix_write(C,i,j,matrix_read(A,i,j) + matrix_read(B,i,j));
- }
- }
- return 1;
- }
- //C = A - B
- //成功返回1,失败返回-1
- int matrix_subtract(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C)
- {
- int i = 0;
- int j = 0;
- //判断是否可以运算
- if (A->m != B->m || A->n != B->n || \
- A->m != C->m || A->n != C->n)
- {
- return -1;
- }
- //运算
- for (i = 0;i < C->m;i++)
- {
- for (j = 0;j < C->n;j++)
- {
- matrix_write(C,i,j,matrix_read(A,i,j) - matrix_read(B,i,j));
- }
- }
- return 1;
- }
- //C = A * B
- //成功返回1,失败返回-1
- int matrix_multiply(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C)
- {
- int i = 0;
- int j = 0;
- int k = 0;
- float temp = 0;
- //判断是否可以运算
- if (A->m != C->m || B->n != C->n || \
- A->n != B->m)
- {
- return -1;
- }
- //运算
- for (i = 0;i < C->m;i++)
- {
- for (j = 0;j < C->n;j++)
- {
- temp = 0;
- for (k = 0;k < A->n;k++)
- {
- temp += matrix_read(A,i,k) * matrix_read(B,k,j);
- }
- matrix_write(C,i,j,temp);
- }
- }
- return 1;
- }
- //行列式的值,只能计算2 * 2,3 * 3
- //失败返回-31415,成功返回值
- float matrix_det(struct _Matrix *A)
- {
- float value = 0;
- //判断是否可以运算
- if (A->m != A->n || (A->m != 2 && A->m != 3))
- {
- return -31415;
- }
- //运算
- if (A->m == 2)
- {
- value = matrix_read(A,0,0) * matrix_read(A,1,1) - matrix_read(A,0,1) * matrix_read(A,1,0);
- }
- else
- {
- value = matrix_read(A,0,0) * matrix_read(A,1,1) * matrix_read(A,2,2) + \
- matrix_read(A,0,1) * matrix_read(A,1,2) * matrix_read(A,2,0) + \
- matrix_read(A,0,2) * matrix_read(A,1,0) * matrix_read(A,2,1) - \
- matrix_read(A,0,0) * matrix_read(A,1,2) * matrix_read(A,2,1) - \
- matrix_read(A,0,1) * matrix_read(A,1,0) * matrix_read(A,2,2) - \
- matrix_read(A,0,2) * matrix_read(A,1,1) * matrix_read(A,2,0);
- }
- return value;
- }
- //求转置矩阵,B = AT
- //成功返回1,失败返回-1
- int matrix_transpos(struct _Matrix *A,struct _Matrix *B)
- {
- int i = 0;
- int j = 0;
- //判断是否可以运算
- if (A->m != B->n || A->n != B->m)
- {
- return -1;
- }
- //运算
- for (i = 0;i < B->m;i++)
- {
- for (j = 0;j < B->n;j++)
- {
- matrix_write(B,i,j,matrix_read(A,j,i));
- }
- }
- return 1;
- }
- //求逆矩阵,B = A^(-1)
- //成功返回1,失败返回-1
- int matrix_inverse(struct _Matrix *A,struct _Matrix *B)
- {
- int i = 0;
- int j = 0;
- int k = 0;
- struct _Matrix m;
- float temp = 0;
- float b = 0;
- //判断是否可以运算
- if (A->m != A->n || B->m != B->n || A->m != B->m)
- {
- return -1;
- }
- /*
- //如果是2维或者3维求行列式判断是否可逆
- if (A->m == 2 || A->m == 3)
- {
- if (det(A) == 0)
- {
- return -1;
- }
- }
- */
- //增广矩阵m = A | B初始化
- matrix_set_m(&m,A->m);
- matrix_set_n(&m,2 * A->m);
- matrix_init(&m);
- for (i = 0;i < m.m;i++)
- {
- for (j = 0;j < m.n;j++)
- {
- if (j <= A->n - 1)
- {
- matrix_write(&m,i,j,matrix_read(A,i,j));
- }
- else
- {
- if (i == j - A->n)
- {
- matrix_write(&m,i,j,1);
- }
- else
- {
- matrix_write(&m,i,j,0);
- }
- }
- }
- }
- //高斯消元
- //变换下三角
- for (k = 0;k < m.m - 1;k++)
- {
- //如果坐标为k,k的数为0,则行变换
- if (matrix_read(&m,k,k) == 0)
- {
- for (i = k + 1;i < m.m;i++)
- {
- if (matrix_read(&m,i,k) != 0)
- {
- break;
- }
- }
- if (i >= m.m)
- {
- return -1;
- }
- else
- {
- //交换行
- for (j = 0;j < m.n;j++)
- {
- temp = matrix_read(&m,k,j);
- matrix_write(&m,k,j,matrix_read(&m,k + 1,j));
- matrix_write(&m,k + 1,j,temp);
- }
- }
- }
- //消元
- for (i = k + 1;i < m.m;i++)
- {
- //获得倍数
- b = matrix_read(&m,i,k) / matrix_read(&m,k,k);
- //行变换
- for (j = 0;j < m.n;j++)
- {
- temp = matrix_read(&m,i,j) - b * matrix_read(&m,k,j);
- matrix_write(&m,i,j,temp);
- }
- }
- }
- //变换上三角
- for (k = m.m - 1;k > 0;k--)
- {
- //如果坐标为k,k的数为0,则行变换
- if (matrix_read(&m,k,k) == 0)
- {
- for (i = k + 1;i < m.m;i++)
- {
- if (matrix_read(&m,i,k) != 0)
- {
- break;
- }
- }
- if (i >= m.m)
- {
- return -1;
- }
- else
- {
- //交换行
- for (j = 0;j < m.n;j++)
- {
- temp = matrix_read(&m,k,j);
- matrix_write(&m,k,j,matrix_read(&m,k + 1,j));
- matrix_write(&m,k + 1,j,temp);
- }
- }
- }
- //消元
- for (i = k - 1;i >= 0;i--)
- {
- //获得倍数
- b = matrix_read(&m,i,k) / matrix_read(&m,k,k);
- //行变换
- for (j = 0;j < m.n;j++)
- {
- temp = matrix_read(&m,i,j) - b * matrix_read(&m,k,j);
- matrix_write(&m,i,j,temp);
- }
- }
- }
- //将左边方阵化为单位矩阵
- for (i = 0;i < m.m;i++)
- {
- if (matrix_read(&m,i,i) != 1)
- {
- //获得倍数
- b = 1 / matrix_read(&m,i,i);
- //行变换
- for (j = 0;j < m.n;j++)
- {
- temp = matrix_read(&m,i,j) * b;
- matrix_write(&m,i,j,temp);
- }
- }
- }
- //求得逆矩阵
- for (i = 0;i < B->m;i++)
- {
- for (j = 0;j < B->m;j++)
- {
- matrix_write(B,i,j,matrix_read(&m,i,j + m.m));
- }
- }
- //释放增广矩阵
- matrix_free(&m);
- return 1;
- }
#include "matrix.h" //矩阵方法 //设置m void matrix_set_m(struct _Matrix *m,int mm) { m->m = mm; } //设置n void matrix_set_n(struct _Matrix *m,int nn) { m->n = nn; } //初始化 void matrix_init(struct _Matrix *m) { m->arr = (float *)malloc(m->m * m->n * sizeof(float)); } //释放 void matrix_free(struct _Matrix *m) { free(m->arr); } //读取i,j坐标的数据 //失败返回-31415,成功返回值 float matrix_read(struct _Matrix *m,int i,int j) { if (i >= m->m || j >= m->n) { return -31415; } return *(m->arr + i * m->n + j); } //写入i,j坐标的数据 //失败返回-1,成功返回1 int matrix_write(struct _Matrix *m,int i,int j,float val) { if (i >= m->m || j >= m->n) { return -1; } *(m->arr + i * m->n + j) = val; return 1; } //矩阵运算 //成功返回1,失败返回-1 int matrix_add(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C) { int i = 0; int j = 0; //判断是否可以运算 if (A->m != B->m || A->n != B->n || \ A->m != C->m || A->n != C->n) { return -1; } //运算 for (i = 0;i < C->m;i++) { for (j = 0;j < C->n;j++) { matrix_write(C,i,j,matrix_read(A,i,j) + matrix_read(B,i,j)); } } return 1; } //C = A - B //成功返回1,失败返回-1 int matrix_subtract(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C) { int i = 0; int j = 0; //判断是否可以运算 if (A->m != B->m || A->n != B->n || \ A->m != C->m || A->n != C->n) { return -1; } //运算 for (i = 0;i < C->m;i++) { for (j = 0;j < C->n;j++) { matrix_write(C,i,j,matrix_read(A,i,j) - matrix_read(B,i,j)); } } return 1; } //C = A * B //成功返回1,失败返回-1 int matrix_multiply(struct _Matrix *A,struct _Matrix *B,struct _Matrix *C) { int i = 0; int j = 0; int k = 0; float temp = 0; //判断是否可以运算 if (A->m != C->m || B->n != C->n || \ A->n != B->m) { return -1; } //运算 for (i = 0;i < C->m;i++) { for (j = 0;j < C->n;j++) { temp = 0; for (k = 0;k < A->n;k++) { temp += matrix_read(A,i,k) * matrix_read(B,k,j); } matrix_write(C,i,j,temp); } } return 1; } //行列式的值,只能计算2 * 2,3 * 3 //失败返回-31415,成功返回值 float matrix_det(struct _Matrix *A) { float value = 0; //判断是否可以运算 if (A->m != A->n || (A->m != 2 && A->m != 3)) { return -31415; } //运算 if (A->m == 2) { value = matrix_read(A,0,0) * matrix_read(A,1,1) - matrix_read(A,0,1) * matrix_read(A,1,0); } else { value = matrix_read(A,0,0) * matrix_read(A,1,1) * matrix_read(A,2,2) + \ matrix_read(A,0,1) * matrix_read(A,1,2) * matrix_read(A,2,0) + \ matrix_read(A,0,2) * matrix_read(A,1,0) * matrix_read(A,2,1) - \ matrix_read(A,0,0) * matrix_read(A,1,2) * matrix_read(A,2,1) - \ matrix_read(A,0,1) * matrix_read(A,1,0) * matrix_read(A,2,2) - \ matrix_read(A,0,2) * matrix_read(A,1,1) * matrix_read(A,2,0); } return value; } //求转置矩阵,B = AT //成功返回1,失败返回-1 int matrix_transpos(struct _Matrix *A,struct _Matrix *B) { int i = 0; int j = 0; //判断是否可以运算 if (A->m != B->n || A->n != B->m) { return -1; } //运算 for (i = 0;i < B->m;i++) { for (j = 0;j < B->n;j++) { matrix_write(B,i,j,matrix_read(A,j,i)); } } return 1; } //求逆矩阵,B = A^(-1) //成功返回1,失败返回-1 int matrix_inverse(struct _Matrix *A,struct _Matrix *B) { int i = 0; int j = 0; int k = 0; struct _Matrix m; float temp = 0; float b = 0; //判断是否可以运算 if (A->m != A->n || B->m != B->n || A->m != B->m) { return -1; } /* //如果是2维或者3维求行列式判断是否可逆 if (A->m == 2 || A->m == 3) { if (det(A) == 0) { return -1; } } */ //增广矩阵m = A | B初始化 matrix_set_m(&m,A->m); matrix_set_n(&m,2 * A->m); matrix_init(&m); for (i = 0;i < m.m;i++) { for (j = 0;j < m.n;j++) { if (j <= A->n - 1) { matrix_write(&m,i,j,matrix_read(A,i,j)); } else { if (i == j - A->n) { matrix_write(&m,i,j,1); } else { matrix_write(&m,i,j,0); } } } } //高斯消元 //变换下三角 for (k = 0;k < m.m - 1;k++) { //如果坐标为k,k的数为0,则行变换 if (matrix_read(&m,k,k) == 0) { for (i = k + 1;i < m.m;i++) { if (matrix_read(&m,i,k) != 0) { break; } } if (i >= m.m) { return -1; } else { //交换行 for (j = 0;j < m.n;j++) { temp = matrix_read(&m,k,j); matrix_write(&m,k,j,matrix_read(&m,k + 1,j)); matrix_write(&m,k + 1,j,temp); } } } //消元 for (i = k + 1;i < m.m;i++) { //获得倍数 b = matrix_read(&m,i,k) / matrix_read(&m,k,k); //行变换 for (j = 0;j < m.n;j++) { temp = matrix_read(&m,i,j) - b * matrix_read(&m,k,j); matrix_write(&m,i,j,temp); } } } //变换上三角 for (k = m.m - 1;k > 0;k--) { //如果坐标为k,k的数为0,则行变换 if (matrix_read(&m,k,k) == 0) { for (i = k + 1;i < m.m;i++) { if (matrix_read(&m,i,k) != 0) { break; } } if (i >= m.m) { return -1; } else { //交换行 for (j = 0;j < m.n;j++) { temp = matrix_read(&m,k,j); matrix_write(&m,k,j,matrix_read(&m,k + 1,j)); matrix_write(&m,k + 1,j,temp); } } } //消元 for (i = k - 1;i >= 0;i--) { //获得倍数 b = matrix_read(&m,i,k) / matrix_read(&m,k,k); //行变换 for (j = 0;j < m.n;j++) { temp = matrix_read(&m,i,j) - b * matrix_read(&m,k,j); matrix_write(&m,i,j,temp); } } } //将左边方阵化为单位矩阵 for (i = 0;i < m.m;i++) { if (matrix_read(&m,i,i) != 1) { //获得倍数 b = 1 / matrix_read(&m,i,i); //行变换 for (j = 0;j < m.n;j++) { temp = matrix_read(&m,i,j) * b; matrix_write(&m,i,j,temp); } } } //求得逆矩阵 for (i = 0;i < B->m;i++) { for (j = 0;j < B->m;j++) { matrix_write(B,i,j,matrix_read(&m,i,j + m.m)); } } //释放增广矩阵 matrix_free(&m); return 1; }
main.c:(测试代码)
- #include <stdio.h>
- #include "matrix.h"
- //打印2维矩阵
- void printf_matrix(struct _Matrix *A)
- {
- int i = 0;
- int j = 0;
- int m = 0;
- int n = 0;
- m = A->m;
- n = A->n;
- for (i = 0;i < m;i++)
- {
- for (j = 0;j < n;j++)
- {
- printf("%f\t",matrix_read(A,i,j));
- }
- printf("\n");
- }
- }
- int main()
- {
- int i = 0;
- int j = 0;
- int k = 0;
- struct _Matrix m1;
- struct _Matrix m2;
- struct _Matrix m3;
- //初始化内存
- matrix_set_m(&m1,3);
- matrix_set_n(&m1,3);
- matrix_init(&m1);
- matrix_set_m(&m2,3);
- matrix_set_n(&m2,3);
- matrix_init(&m2);
- matrix_set_m(&m3,3);
- matrix_set_n(&m3,3);
- matrix_init(&m3);
- //初始化数据
- k = 1;
- for (i = 0;i < m1.m;i++)
- {
- for (j = 0;j < m1.n;j++)
- {
- matrix_write(&m1,i,j,k++);
- }
- }
- for (i = 0;i < m2.m;i++)
- {
- for (j = 0;j < m2.n;j++)
- {
- matrix_write(&m2,i,j,k++);
- }
- }
- //原数据
- printf("A:\n");
- printf_matrix(&m1);
- printf("B:\n");
- printf_matrix(&m2);
- printf("A:行列式的值%f\n",matrix_det(&m1));
- //C = A + B
- if (matrix_add(&m1,&m2,&m3) > 0)
- {
- printf("C = A + B:\n");
- printf_matrix(&m3);
- }
- //C = A - B
- if (matrix_subtract(&m1,&m2,&m3) > 0)
- {
- printf("C = A - B:\n");
- printf_matrix(&m3);
- }
- //C = A * B
- if (matrix_multiply(&m1,&m2,&m3) > 0)
- {
- printf("C = A * B:\n");
- printf_matrix(&m3);
- }
- //C = AT
- if (matrix_transpos(&m1,&m3) > 0)
- {
- printf("C = AT:\n");
- printf_matrix(&m3);
- }
- if (matrix_inverse(&m1,&m3) > 0)
- {
- printf("C = A^(-1):\n");
- printf_matrix(&m3);
- }
- getchar();
- return 0;
- }
环境:
主机:XP
开发环境:mingw
功能:
在C++中实现矩阵运算
说明:
将这篇文章(http://blog.csdn.net/jdh99/article/details/7360091)中在C++下实现的矩阵运算移植到C