欧美bbbwbbbw肥妇,免费乱码人妻系列日韩,一级黄片

純c語言優(yōu)雅地實現(xiàn)矩陣運算庫的方法

 更新時間:2021年08月17日 11:25:17   作者:Marc Pony  
本文主要介紹了純c語言優(yōu)雅地實現(xiàn)矩陣運算庫,文中通過示例代碼介紹的非常詳細(xì),具有一定的參考價值,感興趣的小伙伴們可以參考一下

編程既是技術(shù)輸出也是藝術(shù)創(chuàng)作。鑒賞高手寫的程序,往往讓人眼前一亮,他們思路、邏輯清晰,所呈現(xiàn)的代碼簡潔、優(yōu)雅、高效,令人為之嘆服。爛代碼則像“屎山"一般讓人作嘔,軟件難以維護(hù)最大的原因除了需求模糊的客觀因素,最重要的主觀因素還是代碼寫得爛。有生之年,愿能保持對編程的熱情,不斷提升編程能力,真正體會其樂趣,共勉!

1.一個優(yōu)雅好用的c語言庫必須滿足哪些條件

這里給出的軟件開發(fā)應(yīng)遵循的一般原則,摘自Les Piegl和Wayne Tiller所著的一本非常經(jīng)典的書The NURBS Book(Second Edition)。
(1)工具性(toolability):應(yīng)該使用可用的工具來建立新的應(yīng)用程序;
(2)可移植性(portability):應(yīng)用程序應(yīng)該容易被移植到不同的軟件和硬件平臺;
(3)可重用性(reusability):程序的編寫應(yīng)該便于重復(fù)使用代碼段;
(4)可檢驗性(testability):程序應(yīng)該簡單一致,以便于測試和調(diào)試;
(5)可靠性(reliability):對程序運行過程中可能出現(xiàn)的各種錯誤應(yīng)該進(jìn)行合理、一致的處理,以使系統(tǒng)穩(wěn)定、可靠;
(6)可擴(kuò)展性(enhanceability):代碼必須易于理解,以便可以容易地增加新的功能;
(7)可維護(hù)性(fixability):易于找出程序錯誤的位置;
(8)一致性(consistency):在整個庫中,編程習(xí)慣應(yīng)保持一致;
(9)可讀性(communicability):程序應(yīng)該易于閱讀和理解;
(10)編程風(fēng)格(style of programming):代碼看起來像書上的數(shù)學(xué)公式那樣以便于讀者理解,同時遵循用戶友好的編程風(fēng)格;
(11)易用性(usability):應(yīng)該使一些非專業(yè)的用戶也能夠方便地使用所開發(fā)的庫來開發(fā)各種更高層次的應(yīng)用程序;
(12)數(shù)值高效性(numerical efficiency):所有函數(shù)必須仔細(xì)推敲,保證其數(shù)值高效性;
(13)基于對象編程(object based programming):避免在函數(shù)間傳遞大量數(shù)據(jù),并且使代碼易于理解。

2.實現(xiàn)一個矩陣運算庫的幾點思考

(1)采用預(yù)定義的數(shù)據(jù)類型,避免直接使用編譯器定義的數(shù)據(jù)類型

typedef  unsigned int ERROR_ID;
typedef int INDEX;
typedef short FLAG;
typedef int INTEGER;
typedef double REAL;
typedef char* STRING;
typedef void VOID;

使用預(yù)定義的數(shù)據(jù)類型,有利于程序移植,且可以提高可讀性。例如,如果一個系統(tǒng)只支持單精度浮點數(shù),那么只需修改數(shù)據(jù)類型REAL為float,達(dá)到一勞永逸的效果。定義INDEX與INTEGER數(shù)據(jù)類型是為了在編程時區(qū)分索引變量與普通整形變量,同樣提高可讀性。

(2)基于對象編程,定義矩陣對象

typedef  struct matrix
{
   INTEGER rows;
   INTEGER columns;
   REAL* p;
}MATRIX;

這里,用一級指針而非二級指針指向矩陣的數(shù)據(jù)內(nèi)存地址,有諸多原因,詳見博文:為什么我推薦使用一級指針創(chuàng)建二維數(shù)組?。

(3)除了特別編寫的內(nèi)存處理函數(shù)(使用棧鏈表保存、釋放動態(tài)分配的內(nèi)存地址),不允許任何函數(shù)直接分配和釋放內(nèi)存

不恰當(dāng)?shù)姆峙?、使用與釋放內(nèi)存很可能導(dǎo)致內(nèi)存泄漏、系統(tǒng)崩潰等致命的錯誤。如果一個函數(shù)需動態(tài)申請多個內(nèi)存,那么可能會寫出這樣啰嗦的程序:

double* x = NULL, * y = NULL, * z = NULL;

x = (double*)malloc(n1 * sizeof(double));
if (x == NULL) return -1;

y = (double*)malloc(n2 * sizeof(double));
if (y == NULL)
{
 free(x);
 x = NULL;
 return -1;
}

z = (double*)malloc(n3 * sizeof(double));
if (z == NULL)
{
 free(x);
 x = NULL;
 free(y);
 y = NULL;
 return -1;
}

為了優(yōu)雅地實現(xiàn)動態(tài)內(nèi)存分配與釋放,Les Piegl大神分3步來處理內(nèi)存申請與釋放:
a)在進(jìn)入一個新的程序時,一個內(nèi)存堆棧被初始化為空;
b)當(dāng)需要申請內(nèi)存時,調(diào)用特定的函數(shù)來分配所需的內(nèi)存,并將指向內(nèi)存的指針存入堆棧中的正確位置;
c)在離開程序時,遍歷內(nèi)存堆棧,釋放其中的指針?biāo)赶虻膬?nèi)存。
程序結(jié)構(gòu)大致如下:

STACKS S;
MATRIX* m = NULL;
INTEGER rows = 3, columns = 4;
ERROR_ID errorID = _ERROR_NO_ERROR;

init_stack(&S);

m = creat_matrix(rows, columns, &errorID, &S);
if (m == NULL) goto EXIT;
//do something
// ...

EXIT:
free_stack(&S);
return errorID;

(4)防御性編程,對輸入?yún)?shù)做有效性檢查,并返回錯誤號

例如輸入的矩陣行數(shù)、列數(shù)應(yīng)該是正整數(shù),指針參數(shù)必須非空等等。

(5)注意編程細(xì)節(jié)的打磨

a)操作符(逗號,等號等)兩邊必須空一格;
b)邏輯功能相同的程序間不加空行,邏輯功能獨立的程序間加空行;
c)條件判斷關(guān)鍵字(for if while等)后必須加一空格,起到強(qiáng)調(diào)作用,也更清晰;
d)函數(shù)內(nèi)部定義局部變量后,必須空一行后再編寫函數(shù)主體。

3.完整c程序

本矩陣運算庫只包含了矩陣的基本運算,包括創(chuàng)建任意二維/三維矩陣、創(chuàng)建零矩陣及單位矩陣、矩陣加法、矩陣減法、矩陣乘法、矩陣求逆、矩陣轉(zhuǎn)置、矩陣的跡、矩陣LUP分解、解矩陣方程AX=B。

common.h

/*******************************************************************************
*     File Name :                        common.h
*     Library/Module Name :              MatrixComputation
*     Author :                           Marc Pony(marc_pony@163.com)
*     Create Date :                      2021/6/28
*     Abstract Description :             矩陣運算庫公用頭文件
*******************************************************************************/

#ifndef  __COMMON_H__
#define  __COMMON_H__

/*******************************************************************************
* (1)Debug Switch Section
*******************************************************************************/


/*******************************************************************************
* (2)Include File Section
*******************************************************************************/
#include <math.h>
#include <stdio.h>
#include <malloc.h>
#include <stdlib.h>
#include <time.h>
#include <memory.h>


/*******************************************************************************
* (3)Macro Define Section
*******************************************************************************/
#define _IN
#define _OUT
#define _IN_OUT
#define MAX(x,y) (x) > (y) ? (x) : (y)
#define MIN(x,y) (x) < (y) ? (x) : (y)

#define _CRT_SECURE_NO_WARNINGS
#define PI 3.14159265358979323846
#define POSITIVE_INFINITY 999999999
#define NEGATIVE_INFINITY -999999999
#define _ERROR_NO_ERROR                                          0X00000000   //無錯誤
#define _ERROR_FAILED_TO_ALLOCATE_HEAP_MEMORY                     0X00000001   //分配堆內(nèi)存失敗
#define _ERROR_SVD_EXCEED_MAX_ITERATIONS       0X00000002   //svd超過最大迭代次數(shù)
#define _ERROR_MATRIX_ROWS_OR_COLUMNS_NOT_EQUAL                     0X00000003   //矩陣行數(shù)或列數(shù)不相等
#define _ERROR_MATRIX_MULTIPLICATION        0X00000004   //矩陣乘法錯誤(第一個矩陣的列數(shù)不等于第二個矩陣行數(shù))
#define _ERROR_MATRIX_MUST_BE_SQUARE        0X00000005   //矩陣必須為方陣
#define _ERROR_MATRIX_NORM_TYPE_INVALID        0X00000006   //矩陣模類型無效
#define _ERROR_MATRIX_EQUATION_HAS_NO_SOLUTIONS      0X00000007   //矩陣方程無解
#define _ERROR_MATRIX_EQUATION_HAS_INFINITY_MANNY_SOLUTIONS         0X00000008   //矩陣方程有無窮多解
#define _ERROR_QR_DECOMPOSITION_FAILED        0X00000009   //QR分解失敗
#define _ERROR_CHOLESKY_DECOMPOSITION_FAILED      0X0000000a   //cholesky分解失敗
#define _ERROR_IMPROVED_CHOLESKY_DECOMPOSITION_FAILED    0X0000000b   //improved cholesky分解失敗
#define _ERROR_LU_DECOMPOSITION_FAILED        0X0000000c   //LU分解失敗
#define _ERROR_CREATE_MATRIX_FAILED         0X0000000d   //創(chuàng)建矩陣失敗
#define _ERROR_MATRIX_TRANSPOSE_FAILED        0X0000000e   //矩陣轉(zhuǎn)置失敗
#define _ERROR_CREATE_VECTOR_FAILED         0X0000000f   //創(chuàng)建向量失敗
#define _ERROR_VECTOR_DIMENSION_NOT_EQUAL        0X00000010   //向量維數(shù)不相同
#define _ERROR_VECTOR_NORM_TYPE_INVALID        0X00000011   //向量模類型無效
#define _ERROR_VECTOR_CROSS_FAILED         0X00000012   //向量叉乘失敗
#define _ERROR_INPUT_PARAMETERS_ERROR        0X00010000   //輸入?yún)?shù)錯誤


/*******************************************************************************
* (4)Struct(Data Types) Define Section
*******************************************************************************/
typedef  unsigned int ERROR_ID;
typedef int INDEX;
typedef short FLAG;
typedef int INTEGER;
typedef double REAL;
typedef char* STRING;
typedef void VOID;

typedef  struct matrix
{
 INTEGER rows;
 INTEGER columns;
 REAL* p;
}MATRIX;

typedef struct matrix_node
{
 MATRIX* ptr;
 struct matrix_node* next;
} MATRIX_NODE;

typedef struct matrix_element_node
{
 REAL* ptr;
 struct matrix_element_node* next;
} MATRIX_ELEMENT_NODE;

typedef struct stacks
{
 MATRIX_NODE* matrixNode;
 MATRIX_ELEMENT_NODE* matrixElementNode;

 // ...
 // 添加其他對象的指針
} STACKS;

/*******************************************************************************
* (5)Prototype Declare Section
*******************************************************************************/
/**********************************************************************************************
Function: init_stack
Description: 初始化棧
Input: 無
Output: 無
Input_Output: 棧指針
Return: 無
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
VOID init_stack(_IN_OUT STACKS* S);


/**********************************************************************************************
Function: free_stack
Description: 釋放棧
Input: 棧指針
Output: 無
Input_Output: 無
Return: 無
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
VOID free_stack(_IN STACKS* S);

#endif

matrix.h

/*******************************************************************************
*     File Name :                        matrix.h
*     Library/Module Name :              MatrixComputation
*     Author :                           Marc Pony(marc_pony@163.com)
*     Create Date :                      2021/6/28
*     Abstract Description :             矩陣運算庫頭文件
*******************************************************************************/

#ifndef  __MATRIX_H__
#define  __MATRIX_H__

/*******************************************************************************
* (1)Debug Switch Section
*******************************************************************************/


/*******************************************************************************
* (2)Include File Section
*******************************************************************************/
#include "common.h"


/*******************************************************************************
* (3)Macro Define Section
*******************************************************************************/


/*******************************************************************************
* (4)Struct(Data Types) Define Section
*******************************************************************************/


/*******************************************************************************
* (5)Prototype Declare Section
*******************************************************************************/

void print_matrix(MATRIX* a, STRING string);

/**********************************************************************************************
Function: creat_matrix
Description: 創(chuàng)建矩陣
Input: 矩陣行數(shù)rows,列數(shù)columns
Output: 錯誤號指針errorID,棧指針S
Input_Output: 無
Return: 矩陣指針
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
MATRIX* creat_matrix(_IN INTEGER rows, _IN INTEGER columns, _OUT ERROR_ID* errorID, _OUT STACKS* S);


/**********************************************************************************************
Function: creat_multiple_matrices
Description: 創(chuàng)建多個矩陣
Input: 矩陣行數(shù)rows,列數(shù)columns,個數(shù)count
Output: 錯誤號指針errorID,棧指針S
Input_Output: 無
Return: 矩陣指針
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
MATRIX* creat_multiple_matrices(_IN INTEGER rows, _IN INTEGER columns, _IN INTEGER count, _OUT ERROR_ID* errorID, _OUT STACKS* S);


/**********************************************************************************************
Function: creat_zero_matrix
Description: 創(chuàng)建零矩陣
Input: 矩陣行數(shù)rows,列數(shù)columns
Output: 錯誤號指針errorID,棧指針S
Input_Output: 無
Return: 矩陣指針
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
MATRIX* creat_zero_matrix(_IN INTEGER rows, _IN INTEGER columns, _OUT ERROR_ID* errorID, _OUT STACKS* S);


/**********************************************************************************************
Function: creat_eye_matrix
Description: 創(chuàng)建單位矩陣
Input: 矩陣行數(shù)rows,列數(shù)columns
Output: 錯誤號指針errorID,棧指針S
Input_Output: 無
Return: 矩陣指針
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
MATRIX* creat_eye_matrix(_IN INTEGER n, _OUT ERROR_ID* errorID, _OUT STACKS* S);


/**********************************************************************************************
Function: matrix_add
Description: 矩陣A + 矩陣B = 矩陣C
Input: 矩陣A,矩陣B
Output: 矩陣C
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_add(_IN MATRIX* A, _IN MATRIX* B, _OUT MATRIX *C);

/**********************************************************************************************
Function: matrix_subtraction
Description: 矩陣A - 矩陣B = 矩陣C
Input: 矩陣A,矩陣B
Output: 矩陣C
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_subtraction(_IN MATRIX* A, _IN MATRIX* B, _OUT MATRIX* C);


/**********************************************************************************************
Function: matrix_multiplication
Description: 矩陣乘法C = A * B
Input: 矩陣A,矩陣B
Output: 矩陣C
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_multiplication(_IN MATRIX* A, _IN MATRIX* B, _OUT MATRIX* C);


/**********************************************************************************************
Function: matrix_inverse
Description: 矩陣求逆
Input: 矩陣A
Output: 矩陣A的逆矩陣
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_inverse(_IN MATRIX* A, _OUT MATRIX* invA);


/**********************************************************************************************
Function: matrix_transpose
Description: 矩陣轉(zhuǎn)置
Input: 矩陣A
Output: 矩陣A的轉(zhuǎn)置
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_transpose(_IN MATRIX* A, _OUT MATRIX* transposeA);


/**********************************************************************************************
Function: matrix_trace
Description: 矩陣的跡
Input: 矩陣A
Output: 矩陣A的跡
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_trace(_IN MATRIX* A, _OUT REAL* trace);


/**********************************************************************************************
Function: lup_decomposition
Description: n行n列矩陣LUP分解PA = L * U
Input: n行n列矩陣A
Output: n行n列下三角矩陣L,n行n列上三角矩陣U,n行n列置換矩陣P
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
參考:https://zhuanlan.zhihu.com/p/84210687
***********************************************************************************************/
ERROR_ID lup_decomposition(_IN MATRIX* A, _OUT MATRIX* L, _OUT MATRIX* U, _OUT MATRIX* P);


/**********************************************************************************************
Function: solve_matrix_equation_by_lup_decomposition
Description: LUP分解解矩陣方程AX=B,其中A為n行n列矩陣,B為n行m列矩陣,X為n行m列待求矩陣(寫到矩陣B)
Input: n行n列矩陣A
Output: 無
Input_Output: n行m列矩陣B(即n行m列待求矩陣X)
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID solve_matrix_equation_by_lup_decomposition(_IN MATRIX* A, _IN_OUT MATRIX* B);

#endif

common.c

/*******************************************************************************
*     File Name :                        common.c
*     Library/Module Name :              MatrixComputation
*     Author :                           Marc Pony(marc_pony@163.com)
*     Create Date :                      2021/7/16
*     Abstract Description :            矩陣運算庫公用源文件
*******************************************************************************/

/*******************************************************************************
* (1)Debug Switch Section
*******************************************************************************/


/*******************************************************************************
* (2)Include File Section
*******************************************************************************/
#include "common.h"

/*******************************************************************************
* (3)Macro Define Section
*******************************************************************************/


/*******************************************************************************
* (4)Struct(Data Types) Define Section
*******************************************************************************/


/*******************************************************************************
* (5)Prototype Declare Section
*******************************************************************************/


/*******************************************************************************
* (6)Global Variable Declare Section
*******************************************************************************/


/*******************************************************************************
* (7)File Static Variable Define Section
*******************************************************************************/


/*******************************************************************************
* (8)Function Define Section
*******************************************************************************/

/**********************************************************************************************
Function: init_stack
Description: 初始化棧
Input: 無
Output: 無
Input_Output: 棧指針
Return: 無
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
VOID init_stack(_IN_OUT STACKS* S)
{
 if (S == NULL)
 {
  return;
 }

 memset(S, 0, sizeof(STACKS));
}


/**********************************************************************************************
Function: free_stack
Description: 釋放棧
Input: 棧指針
Output: 無
Input_Output: 無
Return: 無
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
VOID free_stack(_IN STACKS* S)
{
 MATRIX_NODE* matrixNode = NULL;
 MATRIX_ELEMENT_NODE* matrixElementNode = NULL;

 if (S == NULL)
 {
  return;
 }

 while (S->matrixNode != NULL)
 {
  matrixNode = S->matrixNode;
  S->matrixNode = matrixNode->next;

  free(matrixNode->ptr);
  matrixNode->ptr = NULL;
  free(matrixNode);
  matrixNode = NULL;
 }

 while (S->matrixElementNode != NULL)
 {
  matrixElementNode = S->matrixElementNode;
  S->matrixElementNode = matrixElementNode->next;

  free(matrixElementNode->ptr);
  matrixElementNode->ptr = NULL;
  free(matrixElementNode);
  matrixElementNode = NULL;
 }

 // ...
 // 釋放其他指針
}


matrix.c

/*******************************************************************************
*     File Name :                        matrix.c
*     Library/Module Name :              MatrixComputation
*     Author :                           Marc Pony(marc_pony@163.com)
*     Create Date :                      2021/2/24
*     Abstract Description :             矩陣運算庫源文件
*******************************************************************************/

/*******************************************************************************
* (1)Debug Switch Section
*******************************************************************************/


/*******************************************************************************
* (2)Include File Section
*******************************************************************************/
#include "matrix.h"

/*******************************************************************************
* (3)Macro Define Section
*******************************************************************************/


/*******************************************************************************
* (4)Struct(Data Types) Define Section
*******************************************************************************/


/*******************************************************************************
* (5)Prototype Declare Section
*******************************************************************************/


/*******************************************************************************
* (6)Global Variable Declare Section
*******************************************************************************/


/*******************************************************************************
* (7)File Static Variable Define Section
*******************************************************************************/


/*******************************************************************************
* (8)Function Define Section
*******************************************************************************/

VOID print_matrix(MATRIX* a, STRING string)
{
 INDEX i, j;
 printf("matrix %s:", string);
 printf("\n");
 for (i = 0; i < a->rows; i++)
 {
  for (j = 0; j < a->columns; j++)
  {
   printf("%f  ", a->p[i * a->columns + j]);
  }
  printf("\n");
 }

 printf("\n");
}


/**********************************************************************************************
Function: creat_matrix
Description: 創(chuàng)建矩陣
Input: 矩陣行數(shù)rows,列數(shù)columns
Output: 錯誤號指針errorID,棧指針S
Input_Output: 無
Return: 矩陣指針
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
MATRIX* creat_matrix(_IN INTEGER rows, _IN INTEGER columns, _OUT ERROR_ID* errorID, _OUT STACKS* S)
{
 MATRIX* matrix = NULL;
 MATRIX_NODE* matrixNode = NULL;
 MATRIX_ELEMENT_NODE* matrixElementNode = NULL;

 *errorID = _ERROR_NO_ERROR;
 if (rows <= 0 || columns <= 0 || errorID == NULL || S == NULL)
 {
  *errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return NULL;
 }

 matrix = (MATRIX*)malloc(sizeof(MATRIX));
 matrixNode = (MATRIX_NODE*)malloc(sizeof(MATRIX_NODE));
 matrixElementNode = (MATRIX_ELEMENT_NODE*)malloc(sizeof(MATRIX_ELEMENT_NODE));
 if (matrix == NULL || matrixNode == NULL || matrixElementNode == NULL)
 {
  free(matrix);
  matrix = NULL;
  free(matrixNode);
  matrixNode = NULL;
  free(matrixElementNode);
  matrixElementNode = NULL;

  *errorID = _ERROR_FAILED_TO_ALLOCATE_HEAP_MEMORY;
  return NULL;
 }

 matrix->rows = rows;
 matrix->columns = columns;
 matrix->p = (REAL*)malloc(rows * columns * sizeof(REAL));  //確保matrix非空才能執(zhí)行指針操作
 if
 
 
 
 
 (matrix->p == NULL)
 {
  free(matrix->p);
  matrix->p = NULL;
  free(matrix);
  matrix = NULL;
  free(matrixNode);
  matrixNode = NULL;
  free(matrixElementNode);
  matrixElementNode = NULL;

  *errorID = _ERROR_FAILED_TO_ALLOCATE_HEAP_MEMORY;
  return NULL;
 }

 matrixNode->ptr = matrix;
 matrixNode->next = S->matrixNode;
 S->matrixNode = matrixNode;

 matrixElementNode->ptr = matrix->p;
 matrixElementNode->next = S->matrixElementNode;
 S->matrixElementNode = matrixElementNode;

 return matrix;
}


/**********************************************************************************************
Function: creat_multiple_matrices
Description: 創(chuàng)建多個矩陣
Input: 矩陣行數(shù)rows,列數(shù)columns,個數(shù)count
Output: 錯誤號指針errorID,棧指針S
Input_Output: 無
Return: 矩陣指針
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
MATRIX* creat_multiple_matrices(_IN INTEGER rows, _IN INTEGER columns, _IN INTEGER count, _OUT ERROR_ID* errorID, _OUT STACKS* S)
{
 INDEX i;
 MATRIX* matrix = NULL, *p = NULL;
 MATRIX_NODE* matrixNode = NULL;

 *errorID = _ERROR_NO_ERROR;
 if (rows <= 0 || columns <= 0 || count <= 0 || errorID == NULL || S == NULL)
 {
  *errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return NULL;
 }

 matrix = (MATRIX*)malloc(count * sizeof(MATRIX));
 matrixNode = (MATRIX_NODE*)malloc(sizeof(MATRIX_NODE));
 if (matrix == NULL || matrixNode == NULL)
 {
  free(matrix);
  matrix = NULL;
  free(matrixNode);
  matrixNode = NULL;

  *errorID = _ERROR_FAILED_TO_ALLOCATE_HEAP_MEMORY;
  return NULL;
 }

 for (i = 0; i < count; i++)
 {
  p = creat_matrix(rows, columns, errorID, S);
  if (p == NULL)
  {
   free(matrix);
   matrix = NULL;
   free(matrixNode);
   matrixNode = NULL;

   *errorID = _ERROR_FAILED_TO_ALLOCATE_HEAP_MEMORY;
   return NULL;
  }

  matrix[i] = *p;
 }

 matrixNode->ptr = matrix;
 matrixNode->next = S->matrixNode;
 S->matrixNode = matrixNode;

 return matrix;
}

/**********************************************************************************************
Function: creat_zero_matrix
Description: 創(chuàng)建零矩陣
Input: 矩陣行數(shù)rows,列數(shù)columns
Output: 錯誤號指針errorID,棧指針S
Input_Output: 無
Return: 矩陣指針
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
MATRIX* creat_zero_matrix(_IN INTEGER rows, _IN INTEGER columns, _OUT ERROR_ID* errorID, _OUT STACKS* S)
{
 MATRIX* matrix = NULL;

 *errorID = _ERROR_NO_ERROR;
 if (rows <= 0 || columns <= 0 || errorID == NULL || S == NULL)
 {
  *errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return NULL;
 }

 matrix = creat_matrix(rows, columns, errorID, S);
 if (*errorID == _ERROR_NO_ERROR)
 {
  memset(matrix->p, 0, rows * columns * sizeof(REAL));
 }

 return matrix;
}


/**********************************************************************************************
Function: creat_eye_matrix
Description: 創(chuàng)建單位矩陣
Input: 矩陣行數(shù)rows,列數(shù)columns
Output: 錯誤號指針errorID,棧指針S
Input_Output: 無
Return: 矩陣指針
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
MATRIX* creat_eye_matrix(_IN INTEGER n, _OUT ERROR_ID* errorID, _OUT STACKS* S)
{
 INDEX i;
 MATRIX* matrix = NULL;

 *errorID = _ERROR_NO_ERROR;
 if (n <= 0 || errorID == NULL || S == NULL)
 {
  *errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return NULL;
 }

 matrix = creat_matrix(n, n, errorID, S);
 if (*errorID == _ERROR_NO_ERROR)
 {
  memset(matrix->p, 0, n * n * sizeof(REAL));
  for (i = 0; i < n; i++)
  {
   matrix->p[i * n + i] = 1.0;
  }
 }

 return matrix;
}


/**********************************************************************************************
Function: matrix_add
Description: 矩陣A + 矩陣B = 矩陣C
Input: 矩陣A,矩陣B
Output: 矩陣C
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_add(_IN MATRIX* A, _IN MATRIX* B, _OUT MATRIX* C)
{
 INDEX i, j;
 INTEGER rows, columns;
 ERROR_ID errorID = _ERROR_NO_ERROR;

 if (A == NULL || B == NULL || C == NULL)
 {
  errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return errorID;
 }

 if (A->rows != B->rows || A->rows != C->rows || B->rows != C->rows
  ||  A->columns != B->columns || A->columns != C->columns || B->columns != C->columns)
 {
  errorID = _ERROR_MATRIX_ROWS_OR_COLUMNS_NOT_EQUAL;
  return errorID;
 }

 rows = A->rows;
 columns = A->columns;
 for (i = 0; i < rows; i++)
 {
  for (j = 0; j < columns; j++)
  {
   C->p[i * columns + j] = A->p[i * columns + j] + B->p[i * columns + j];
  }
 }
 return errorID;
}


/**********************************************************************************************
Function: matrix_subtraction
Description: 矩陣A - 矩陣B = 矩陣C
Input: 矩陣A,矩陣B
Output: 矩陣C
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_subtraction(_IN MATRIX* A, _IN MATRIX* B, _OUT MATRIX* C)
{
 INDEX i, j;
 INTEGER rows, columns;
 ERROR_ID errorID = _ERROR_NO_ERROR;

 if (A == NULL || B == NULL || C == NULL)
 {
  errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return errorID;
 }

 if (A->rows != B->rows || A->rows != C->rows || B->rows != C->rows
  || A->columns != B->columns || A->columns != C->columns || B->columns != C->columns)
 {
  errorID = _ERROR_MATRIX_ROWS_OR_COLUMNS_NOT_EQUAL;
  return errorID;
 }

 rows = A->rows;
 columns = A->columns;
 for (i = 0; i < rows; i++)
 {
  for (j = 0; j < columns; j++)
  {
   C->p[i * columns + j] = A->p[i * columns + j] - B->p[i * columns + j];
  }
 }
 return errorID;
}



/**********************************************************************************************
Function: matrix_multiplication
Description: 矩陣乘法C = A * B
Input: 矩陣A,矩陣B
Output: 矩陣C
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_multiplication(_IN MATRIX* A, _IN MATRIX* B, _OUT MATRIX* C)
{
 INDEX  i, j, k;
 REAL sum;
 ERROR_ID errorID = _ERROR_NO_ERROR;

 if (A == NULL || B == NULL || C == NULL)
 {
  errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return errorID;
 }

 if (A->columns != B->rows || A->rows != C->rows || B->columns != C->columns)
 {
  errorID = _ERROR_MATRIX_MULTIPLICATION;
  return errorID;
 }

 for (i = 0; i < A->rows; i++)
 {
  for (j = 0; j < B->columns; j++)
  {
   sum = 0.0;
   for (k = 0; k < A->columns; k++)
   {
    sum += A->p[i * A->columns + k] * B->p[k * B->columns + j];
   }
   C->p[i * B->columns + j] = sum;
  }
 }
 return errorID;
}

/**********************************************************************************************
Function: matrix_inverse
Description: 矩陣求逆
Input: 矩陣A
Output: 矩陣A的逆矩陣
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_inverse(_IN MATRIX* A, _OUT MATRIX* invA)
{
 INDEX i;
 INTEGER n;
 MATRIX* ATemp = NULL;
 ERROR_ID errorID = _ERROR_NO_ERROR;
 STACKS S;

 if (A == NULL || invA == NULL)
 {
  errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return errorID;
 }

 if (A->rows != A->columns)
 {
  errorID = _ERROR_MATRIX_MUST_BE_SQUARE;
  return errorID;
 }

 init_stack(&S);

 n = A->rows;
 ATemp = creat_matrix(n, n, &errorID, &S);
 if (errorID != _ERROR_NO_ERROR) goto EXIT;

 memcpy(ATemp->p, A->p, n * n * sizeof(REAL));
 memset(invA->p, 0, n * n * sizeof(REAL));
 for (i = 0; i < n; i++)
 {
  invA->p[i * n + i] = 1.0;
 }

 errorID = solve_matrix_equation_by_lup_decomposition(ATemp, invA);

EXIT:
 free_stack(&S);
 return errorID;
}


/**********************************************************************************************
Function: matrix_transpose
Description: 矩陣轉(zhuǎn)置
Input: 矩陣A
Output: 矩陣A的轉(zhuǎn)置
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_transpose(_IN MATRIX* A, _OUT MATRIX* transposeA)
{
 INDEX i, j;
 ERROR_ID errorID = _ERROR_NO_ERROR;

 if (A == NULL || transposeA == NULL)
 {
  errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return errorID;
 }

 if (A->rows != transposeA->columns || A->columns != transposeA->rows)
 {
  errorID = _ERROR_MATRIX_TRANSPOSE_FAILED;
  return errorID;
 }

 for (i = 0; i < A->rows; i++)
 {
  for (j = 0; j < A->columns; j++)
  {
   transposeA->p[j * A->rows + i] = A->p[i * A->columns + j];
  }
 }
 return errorID;
}


/**********************************************************************************************
Function: matrix_trace
Description: 矩陣的跡
Input: 矩陣A
Output: 矩陣A的跡
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID matrix_trace(_IN MATRIX* A, _OUT REAL *trace)
{
 INDEX i;
 ERROR_ID errorID = _ERROR_NO_ERROR;

 if (A == NULL || trace == NULL)
 {
  errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return errorID;
 }

 if (A->rows != A->columns)
 {
  errorID = _ERROR_MATRIX_MUST_BE_SQUARE;
  return errorID;
 }

 *trace = 0.0;
 for (i = 0; i < A->rows; i++)
 {
  *trace += A->p[i * A->columns + i];
 }
 return errorID;
}


/**********************************************************************************************
Function: lup_decomposition
Description: n行n列矩陣LUP分解PA = L * U
Input: n行n列矩陣A
Output: n行n列下三角矩陣L,n行n列上三角矩陣U,n行n列置換矩陣P
Input_Output: 無
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
參考:https://zhuanlan.zhihu.com/p/84210687
***********************************************************************************************/
ERROR_ID lup_decomposition(_IN MATRIX* A, _OUT MATRIX* L, _OUT MATRIX* U, _OUT MATRIX* P)
{
 INDEX i, j, k, index, s, t;
 INTEGER n;
 REAL maxValue, temp;
 ERROR_ID errorID = _ERROR_NO_ERROR;

 if (A == NULL || L == NULL || U == NULL || P == NULL)
 {
  errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return errorID;
 }

 if (A->rows != A->columns)
 {
  errorID = _ERROR_MATRIX_MUST_BE_SQUARE;
  return errorID;
 }

 n = A->rows;
 memcpy(U->p, A->p, n * n * sizeof(REAL));
 memset(L->p, 0, n * n * sizeof(REAL));
 memset(P->p, 0, n * n * sizeof(REAL));
 for (i = 0; i < n; i++)
 {
  L->p[i * n + i] = 1.0;
  P->p[i * n + i] = 1.0;
 }

 for (j = 0; j < n - 1; j++)
 {
  //Select i(>= j) that maximizes | U(i, j) |
  index = -1;
  maxValue = 0.0;
  for (i = j; i < n; i++)
  {
   temp = fabs(U->p[i * n + j]);
   if (temp > maxValue)
   {
    maxValue = temp;
    index = i;
   }
  }

  if (index == -1)
  {
   continue;
  }

  //Interchange rows of U : U(j, j : n) < ->U(i, j : n)
  for (k = j; k < n; k++)
  {
   s = j * n + k;
   t = index * n + k;
   temp = U->p[s];
   U->p[s] = U->p[t];
   U->p[t] = temp;
  }

  //Interchange rows of L : L(j, 1 : j - 1) < ->L(i, 1 : j - 1)
  for (k = 0; k < j; k++)
  {
   s = j * n + k;
   t = index * n + k;
   temp = L->p[s];
   L->p[s] = L->p[t];
   L->p[t] = temp;
  }

  //Interchange rows of P : P(j, 1 : n) < ->P(i, 1 : n)
  for (k = 0; k < n; k++)
  {
   s = j * n + k;
   t = index * n + k;
   temp = P->p[s];
   P->p[s] = P->p[t];
   P->p[t] = temp;
  }

  for (i = j + 1; i < n; i++)
  {
   s = i * n + j;
   L->p[s] = U->p[s] / U->p[j * n + j];
   for (k = j; k < n; k++)
   {
    U->p[i * n + k] -= L->p[s] * U->p[j * n + k];
   }
  }
 }
 return errorID;
}


/**********************************************************************************************
Function: solve_matrix_equation_by_lup_decomposition
Description: LUP分解解矩陣方程AX=B,其中A為n行n列矩陣,B為n行m列矩陣,X為n行m列待求矩陣(寫到矩陣B)
Input: n行n列矩陣A
Output: 無
Input_Output: n行m列矩陣B(即n行m列待求矩陣X)
Return: 錯誤號
Author: Marc Pony(marc_pony@163.com)
***********************************************************************************************/
ERROR_ID solve_matrix_equation_by_lup_decomposition(_IN MATRIX* A, _IN_OUT MATRIX* B)
{
 INDEX i, j, k, index, s, t;
 INTEGER n, m;
 REAL sum, maxValue, temp;
 MATRIX* L = NULL, * U = NULL, * y = NULL;
 ERROR_ID errorID = _ERROR_NO_ERROR;
 STACKS S;

 if (A == NULL || B == NULL)
 {
  errorID = _ERROR_INPUT_PARAMETERS_ERROR;
  return errorID;
 }

 if (A->rows != A->columns)
 {
  errorID = _ERROR_MATRIX_MUST_BE_SQUARE;
  return errorID;
 }

 init_stack(&S);

 n = A->rows;
 m = B->columns;

 L = creat_matrix(n, n, &errorID, &S);
 if (errorID != _ERROR_NO_ERROR) goto EXIT;

 U = creat_matrix(n, n, &errorID, &S);
 if (errorID != _ERROR_NO_ERROR) goto EXIT;

 y = creat_matrix(n, m, &errorID, &S);
 if (errorID != _ERROR_NO_ERROR) goto EXIT;

 memcpy(U->p, A->p, n * n * sizeof(REAL));
 memset(L->p, 0, n * n * sizeof(REAL));
 for (i = 0; i < n; i++)
 {
  L->p[i * n + i] = 1.0;
 }

 for (j = 0; j < n - 1; j++)
 {
  //Select i(>= j) that maximizes | U(i, j) |
  index = -1;
  maxValue = 0.0;
  for (i = j; i < n; i++)
  {
   temp = fabs(U->p[i * n + j]);
   if (temp > maxValue)
   {
    maxValue = temp;
    index = i;
   }
  }

  if (index == -1)
  {
   continue;
  }

  //Interchange rows of U : U(j, j : n) < ->U(i, j : n)
  for (k = j; k < n; k++)
  {
   s = j * n + k;
   t = index * n + k;
   temp = U->p[s];
   U->p[s] = U->p[t];
   U->p[t] = temp;
  }

  //Interchange rows of L : L(j, 1 : j - 1) < ->L(i, 1 : j - 1)
  for (k = 0; k < j; k++)
  {
   s = j * n + k;
   t = index * n + k;
   temp = L->p[s];
   L->p[s] = L->p[t];
   L->p[t] = temp;
  }

  //Interchange rows of P : P(j, 1 : n) < ->P(i, 1 : n), C = P * B,等價于對B交換行
  for (k = 0; k < m; k++)
  {
   s = j * m + k;
   t = index * m + k;
   temp = B->p[s];
   B->p[s] = B->p[t];
   B->p[t] = temp;
  }

  for (i = j + 1; i < n; i++)
  {
   s = i * n + j;
   L->p[s] = U->p[s] / U->p[j * n + j];
   for (k = j; k < n; k++)
   {
    U->p[i * n + k] -= L->p[s] * U->p[j * n + k];
   }
  }
 }

 for (i = 0; i < n; i++)
 {
  if (fabs(U->p[i * n + i]) < 1.0e-20)
  {
   errorID = _ERROR_MATRIX_EQUATION_HAS_NO_SOLUTIONS;
   goto EXIT;
  }
 }

 //L * y = C
 for (j = 0; j < m; j++)
 {
  for (i = 0; i < n; i++)
  {
   sum = 0.0;
   for (k = 0; k < i; k++)
   {
    sum = sum + L->p[i * n + k] * y->p[k * m + j];
   }
   y->p[i * m + j] = B->p[i * m + j] - sum;
  }
 }

 //U * x = y
 for (j = 0; j < m; j++)
 {
  for (i = n - 1; i >= 0; i--)
  {
   sum = 0.0;
   for (k = i + 1; k < n; k++)
   {
    sum += U->p[i * n + k] * B->p[k * m + j];
   }
   B->p[i * m + j] = (y->p[i * m + j] - sum) / U->p[i * n + i];
  }
 }

EXIT:
 free_stack(&S);
 return errorID;
}

test_matrix.c

#include "matrix.h"

void main()
{
 REAL a[3 * 3] = { 1,2,3,6,5,5,8,7,2 };
 REAL b[3 * 3] = {1,2,3,6,5,4,3,2,1};
 MATRIX *A = NULL, * B = NULL, * C = NULL, * D = NULL, * E = NULL, * Z = NULL, * invA = NULL, *m = NULL;
 ERROR_ID errorID = _ERROR_NO_ERROR;
 REAL trace;
 STACKS S;

 init_stack(&S);

 Z = creat_zero_matrix(3, 3, &errorID, &S);
 print_matrix(Z, "Z");

 E = creat_eye_matrix(3, &errorID, &S);
 print_matrix(E, "E");

 A = creat_matrix(3, 3, &errorID, &S);
 A->p = a;
 print_matrix(A, "A");

 B = creat_matrix(3, 3, &errorID, &S);
 B->p = b;
 print_matrix(B, "B");

 C = creat_matrix(3, 3, &errorID, &S);
 D = creat_matrix(3, 3, &errorID, &S);
 invA = creat_matrix(3, 3, &errorID, &S);

 errorID = matrix_add(A, B, C);
 errorID = matrix_subtraction(A, B, C);
 errorID = matrix_multiplication(A, B, C);

 errorID = matrix_transpose(A, D);
 print_matrix(D, "D");

 errorID = matrix_trace(A, &trace);

 errorID = matrix_inverse(A, invA);
    print_matrix(invA, "invA");

 m = creat_multiple_matrices(3, 3, 2, &errorID, &S);
 m[0].p = a;
 m[1].p = b;

 free_stack(&S);
}

參考資料

The NURBS Book(Second Edition). Les Piegl,Wayne Tiller

到此這篇關(guān)于純c語言優(yōu)雅地實現(xiàn)矩陣運算庫的方法的文章就介紹到這了,更多相關(guān)c語言 矩陣運算內(nèi)容請搜索腳本之家以前的文章或繼續(xù)瀏覽下面的相關(guān)文章希望大家以后多多支持腳本之家!

相關(guān)文章

  • C++并查集常用操作

    C++并查集常用操作

    并查集 是一種樹型的數(shù)據(jù)結(jié)構(gòu),用于處理一些不相加集合的合并和查詢問題。本文給大家分享C++并查集常用操作及算法實現(xiàn),感興趣的朋友跟隨小編一起看看吧
    2021-07-07
  • Qt模仿IOS滑動按鈕效果

    Qt模仿IOS滑動按鈕效果

    這篇文章主要為大家詳細(xì)介紹了Qt模仿IOS滑動按鈕效果,文中示例代碼介紹的非常詳細(xì),具有一定的參考價值,感興趣的小伙伴們可以參考一下
    2020-08-08
  • OpenSSL動態(tài)鏈接庫源碼安裝教程

    OpenSSL動態(tài)鏈接庫源碼安裝教程

    Openssl 是一個開放源代碼的SSL協(xié)議的產(chǎn)品實現(xiàn),它采用C語言作為開發(fā)語言,具備了跨系統(tǒng)的性能。這篇文章主要介紹了OpenSSL動態(tài)鏈接庫源碼安裝,需要的朋友可以參考下
    2021-11-11
  • C語言?超詳細(xì)講解算法的時間復(fù)雜度和空間復(fù)雜度

    C語言?超詳細(xì)講解算法的時間復(fù)雜度和空間復(fù)雜度

    算法復(fù)雜度分為時間復(fù)雜度和空間復(fù)雜度。其作用:?時間復(fù)雜度是度量算法執(zhí)行的時間長短;而空間復(fù)雜度是度量算法所需存儲空間的大小
    2022-03-03
  • VC小技巧匯總之5則實用小技巧

    VC小技巧匯總之5則實用小技巧

    這篇文章主要介紹了VC小技巧匯總之5則實用小技巧,需要的朋友可以參考下
    2014-07-07
  • 關(guān)于C++友元函數(shù)的實現(xiàn)講解

    關(guān)于C++友元函數(shù)的實現(xiàn)講解

    今天小編就為大家分享一篇關(guān)于關(guān)于C++友元函數(shù)的實現(xiàn)講解,小編覺得內(nèi)容挺不錯的,現(xiàn)在分享給大家,具有很好的參考價值,需要的朋友一起跟隨小編來看看吧
    2018-12-12
  • VS中PCL庫附加依賴項配置過程解析

    VS中PCL庫附加依賴項配置過程解析

    這篇文章主要介紹了VS中PCL庫附加依賴項配置,本文給大家介紹的非常詳細(xì),對大家的學(xué)習(xí)或工作具有一定的參考借鑒價值,需要的朋友可以參考下
    2022-07-07
  • C++中的對象數(shù)組詳細(xì)解析

    C++中的對象數(shù)組詳細(xì)解析

    在建立數(shù)組時,同樣要調(diào)用構(gòu)造函數(shù)。如果有50個元素,就需要調(diào)用50次構(gòu)造函數(shù)。在需要的時候,可以在定義數(shù)組時提供實參以實現(xiàn)初始化
    2013-10-10
  • 總結(jié)C/C++面試中可能會碰到的字符串指針題

    總結(jié)C/C++面試中可能會碰到的字符串指針題

    C/C++是最能體現(xiàn)程序員能力的語言之一,其功能強(qiáng)大,在IT行業(yè)的各個方面都有大量的應(yīng)用。下面這篇文章主要介紹了總結(jié)了在C/C++面試中可能會碰到的字符串指針題,需要的朋友可以參考借鑒,下面來一起看看吧。
    2017-01-01
  • 用C++實現(xiàn)隊列的程序代碼

    用C++實現(xiàn)隊列的程序代碼

    本篇文章是對使用C++實現(xiàn)隊列的程序代碼進(jìn)行了詳細(xì)的分析介紹,需要的朋友參考下
    2013-05-05

最新評論