通过 C 中的辅助矩阵求解线性方程
Solving linear equations via adjunt matrix in C
下面是学长给我的一个C代码，是求一个3*3的联立一次方程的解。他基本上用了adjucate matrix的方法 A-1 = 1|A| * adj A 求线性方程的解。编译时代码不是 运行（通过 GCC）。任何帮助将不胜感激。下面是代码：
#include <stdio.h>
#include <curses.h>
#include <stdlib.h>
void swap(float*a, float*b){
float temp;
temp=*a;
*a=*b;
*b=temp;}
int main(){
int i, j, cof_1,cof_2,cof_3,cof_4;
float a[3][3], c[3], A[3][3], INV[3][3], det=0.0, X[3], B[3]={0,0,0}, n;
for(i=0;i<3;i++){
printf("Enter x, y & z co-efficient and the constant term of equation number %d: \n", i+1);
scanf("%f%f%f%f", &a[i][0], &a[i][1], &a[i][2], &c[i]);
}
for(i=0;i<3;i++)
{for(j=0;j<3;j++)
{cof_1=(i+1)%3;
cof_2=(j+1)%3;
cof_3=(i+2)%3;
cof_4=(j+2)%3;
A[i][j] =(a[cof_1][cof_2]*a[cof_3][cof_4])-(a[cof_1][cof_4]*a[cof_3][cof_2]);}
}
for(i=0;i<3;i++){
det+=a[i][0]*A[i][0];}
printf("\n\n\n");
for(i=0;i<3;i++){
swap (&A[i][j], &A[j][i]);}
for(i=0;i<3;i++){
For (j=0;j<3;j++){
B[i]+=A[i][j]*c[j];
}
}
if(det==0){
If (B[0]==0, B[1]==0, B[2]==0){
printf("The system is consistent and there are infinitely many solutions.");
exit(1);}
else
printf("There is no solution to this system.");
exit(2);
}
for (i=0;i<3;i++){ //for cases with possible solutions
for (j=0;j<3;j++)
INV[i][j]=A[i][j]/det;
}
for(i=0;i<3;i++){
X[i]=0.0;
for (j=0;j<3;j++){
X[i]+=INV[i][j]*c[j];
}
}
if(X[0]==0, X[1]==0, X[2]==0){
printf ("This system has a trivial solution. \n"); //The case for trivial solution}
else {
printf ("This system has a unique solution. \n"); //The case for unique solution}
printf ("Solution of the equations: \n\tx=%.2f\n\ty=%.2f\n\tz=%.2f", X[0], X[1], X[2]);
return 0;
}
获得未注释代码“来自我的前辈（好友）”的问题在于，您对代码中使用的算法完全没有任何指导和零参考。很多时候，如果您花时间学习所涉及的数学并自己编写解决方案，这会使您处于更糟糕的境地。
为什么？您没有从学习所涉及的数学中获益，然后在逻辑上应用这些知识来形成数字解决方案，而是最终猜测代码中的错误位置，并且不必要地猜测、重新编译并再次失败。
综合这个事实，你从你朋友那里得到的代码甚至不包含有效的 C。C 中没有“For”循环，没有“If”语句和没有“Printf”功能。事实上，none 的标准库函数和 none 的运算符完全以大写字母开头。包含 curses.h header 是多余的。
此外，您通过允许 j 的值在循环中超出数组 A 的边界来调用 未定义的行为：
for(i=0;i<3;i++){
swap (&A[i][j], &A[j][i]);}
为什么？因为 j 的值是 3，在您从 cof_x 填充 A 的嵌套循环之后，并且在尝试转置之前它永远不会被重置。
（注意： 结果证明这只是您尝试转置辅助因子矩阵的问题之一）
你们中的一个主要抱怨是代码无法编译。当然不是。但是，如果您 启用编译器警告 ，编译器会告诉您发生问题的确切行，帮助您修复代码，以便它可以编译。对于 gcc，将 -Wall -Wextra -pedantic 添加到您的编译字符串中。对于 clang，添加 -Weverything，对于 VS (cl.exe) 添加 /W3（或 /Wall 用于所有警告）。在没有警告的情况下干净地编译之前不要接受代码。
一旦掌握了基础知识，就可以转向算法的细节。虽然形成辅词（或正方形的辅词）的代码看起来很可疑，但它实际上正确地构成了辅词矩阵。正确计算确定性。您的代码失败的地方是通过在遍历整个矩阵时不正确地交换元素来计算转置调整矩阵。
您有两个转置选项。遍历整个边界 assigning A[j][i] = A[i][j]; 或对于方矩阵，仅遍历那些改变交换值的元素。
更正转置后（并参考有问题的操作注释代码），代码工作正常，例如
#include <stdio.h>
#define NDIM 3 /* if you need a constant - define one (or more) */
void swap (float *a, float *b)
{
float temp;
temp = *a;
*a = *b;
*b = temp;
}
int main (void) {
int i, j, cof_1, cof_2, cof_3, cof_4;
float a[NDIM][NDIM] = {{0}}, /* coefficient matrix */
c[NDIM] = {0}, /* solution vector */
A[NDIM][NDIM] = {{0}}, /* adjunct matrix */
INV[NDIM][NDIM] = {{0}}, /* inverse of a */
det = 0.0, /* determinant of a */
X[NDIM] = {0}, /* linear system roots */
B[NDIM] = {0}; /* product A * c */
i = j = cof_1 = cof_2 = cof_3 = cof_4 = 0; /* initialize */
for (i = 0; i < NDIM; i++) { /* validate input */
printf ("Eq[%d] coefficients and constant: ", i + 1);
if (scanf ("%f%f%f%f", &a[i][0], &a[i][1], &a[i][2], &c[i]) != 4) {
fprintf (stderr, "error: invalid coefficient.\n");
return 1;
}
}
putchar ('\n');
/* cofactor matrix from matrix of minors incorporating
* patterned +/- sign to maxtix of minors.
*/
for (i = 0; i < NDIM; i++)
for (j = 0; j < NDIM; j++) {
cof_1 = (i + 1) % NDIM;
cof_2 = (j + 1) % NDIM;
cof_3 = (i + 2) % NDIM;
cof_4 = (j + 2) % NDIM;
A[i][j] = (a[cof_1][cof_2] * a[cof_3][cof_4]) -
(a[cof_1][cof_4] * a[cof_3][cof_2]);
}
/* compute determinant */
for (i = 0; i < NDIM; i++)
det += a[i][0] * A[i][0];
/* transpose cofactor (adjunct) */
for (i = 1; i < NDIM; i++)
for (j = 0; j < i; j++)
swap (&A[i][j], &A[j][i]);
/* consistency check B */
for (i = 0; i < NDIM; i++)
for (j = 0; j < NDIM; j++)
B[i] += A[i][j] * c[j];
if (det == 0) { /* eliminate infinite and no solution cases */
if (B[0] == 0 && B[1] == 0 && B[2] == 0) {
printf ("The system is consistent and there "
"are infinitely many solutions.\n");
return 1;
} else
printf ("There is no solution to this system.\n");
return 2;
}
/* compute inverse of A (1/det * A) */
for (i = 0; i < NDIM; i++)
for (j = 0; j < NDIM; j++)
INV[i][j] = A[i][j] / det;
/* multiply constant vector by inverse */
for (i = 0; i < NDIM; i++)
for (j = 0; j < NDIM; j++)
X[i] += INV[i][j] * c[j];
if (X[0] == 0 && X[1] == 0 && X[2] == 0)
/* The case for trivial solution */
printf ("This system has a trivial solution.\n");
else
/* The case for unique solution */
printf ("This system has a unique solution.\n\n");
/* output solutions for linear system of eq. */
for (i = 0; i < NDIM; i++)
printf ("\t%c = %.2f\n", 'x' + i, X[i]);
return 0;
}
例子Use/Output
用于示例的 30-60-90 三角形的简单证明中的线性方程：
$ ./bin/matrixinv
Eq[1] coefficients and constant: 1 1 1 180
Eq[2] coefficients and constant: 5 -1 -1 0
Eq[3] coefficients and constant: 1 1 -1 0
This system has a unique solution.
x = 30.00
y = 60.00
z = 90.00
检查一下，如果您还有其他问题，请告诉我。
下面是学长给我的一个C代码，是求一个3*3的联立一次方程的解。他基本上用了adjucate matrix的方法 A-1 = 1|A| * adj A 求线性方程的解。编译时代码不是 运行（通过 GCC）。任何帮助将不胜感激。下面是代码：
#include <stdio.h>
#include <curses.h>
#include <stdlib.h>
void swap(float*a, float*b){
float temp;
temp=*a;
*a=*b;
*b=temp;}
int main(){
int i, j, cof_1,cof_2,cof_3,cof_4;
float a[3][3], c[3], A[3][3], INV[3][3], det=0.0, X[3], B[3]={0,0,0}, n;
for(i=0;i<3;i++){
printf("Enter x, y & z co-efficient and the constant term of equation number %d: \n", i+1);
scanf("%f%f%f%f", &a[i][0], &a[i][1], &a[i][2], &c[i]);
}
for(i=0;i<3;i++)
{for(j=0;j<3;j++)
{cof_1=(i+1)%3;
cof_2=(j+1)%3;
cof_3=(i+2)%3;
cof_4=(j+2)%3;
A[i][j] =(a[cof_1][cof_2]*a[cof_3][cof_4])-(a[cof_1][cof_4]*a[cof_3][cof_2]);}
}
for(i=0;i<3;i++){
det+=a[i][0]*A[i][0];}
printf("\n\n\n");
for(i=0;i<3;i++){
swap (&A[i][j], &A[j][i]);}
for(i=0;i<3;i++){
For (j=0;j<3;j++){
B[i]+=A[i][j]*c[j];
}
}
if(det==0){
If (B[0]==0, B[1]==0, B[2]==0){
printf("The system is consistent and there are infinitely many solutions.");
exit(1);}
else
printf("There is no solution to this system.");
exit(2);
}
for (i=0;i<3;i++){ //for cases with possible solutions
for (j=0;j<3;j++)
INV[i][j]=A[i][j]/det;
}
for(i=0;i<3;i++){
X[i]=0.0;
for (j=0;j<3;j++){
X[i]+=INV[i][j]*c[j];
}
}
if(X[0]==0, X[1]==0, X[2]==0){
printf ("This system has a trivial solution. \n"); //The case for trivial solution}
else {
printf ("This system has a unique solution. \n"); //The case for unique solution}
printf ("Solution of the equations: \n\tx=%.2f\n\ty=%.2f\n\tz=%.2f", X[0], X[1], X[2]);
return 0;
}
获得未注释代码“来自我的前辈（好友）”的问题在于，您对代码中使用的算法完全没有任何指导和零参考。很多时候，如果您花时间学习所涉及的数学并自己编写解决方案，这会使您处于更糟糕的境地。
为什么？您没有从学习所涉及的数学中获益，然后在逻辑上应用这些知识来形成数字解决方案，而是最终猜测代码中的错误位置，并且不必要地猜测、重新编译并再次失败。
综合这个事实，你从你朋友那里得到的代码甚至不包含有效的 C。C 中没有“For”循环，没有“If”语句和没有“Printf”功能。事实上，none 的标准库函数和 none 的运算符完全以大写字母开头。包含 curses.h header 是多余的。
此外，您通过允许 j 的值在循环中超出数组 A 的边界来调用 未定义的行为：
for(i=0;i<3;i++){
swap (&A[i][j], &A[j][i]);}
为什么？因为 j 的值是 3，在您从 cof_x 填充 A 的嵌套循环之后，并且在尝试转置之前它永远不会被重置。
（注意： 结果证明这只是您尝试转置辅助因子矩阵的问题之一）
你们中的一个主要抱怨是代码无法编译。当然不是。但是，如果您 启用编译器警告 ，编译器会告诉您发生问题的确切行，帮助您修复代码，以便它可以编译。对于 gcc，将 -Wall -Wextra -pedantic 添加到您的编译字符串中。对于 clang，添加 -Weverything，对于 VS (cl.exe) 添加 /W3（或 /Wall 用于所有警告）。在没有警告的情况下干净地编译之前不要接受代码。
一旦掌握了基础知识，就可以转向算法的细节。虽然形成辅词（或正方形的辅词）的代码看起来很可疑，但它实际上正确地构成了辅词矩阵。正确计算确定性。您的代码失败的地方是通过在遍历整个矩阵时不正确地交换元素来计算转置调整矩阵。
您有两个转置选项。遍历整个边界 assigning A[j][i] = A[i][j]; 或对于方矩阵，仅遍历那些改变交换值的元素。
更正转置后（并参考有问题的操作注释代码），代码工作正常，例如
#include <stdio.h>
#define NDIM 3 /* if you need a constant - define one (or more) */
void swap (float *a, float *b)
{
float temp;
temp = *a;
*a = *b;
*b = temp;
}
int main (void) {
int i, j, cof_1, cof_2, cof_3, cof_4;
float a[NDIM][NDIM] = {{0}}, /* coefficient matrix */
c[NDIM] = {0}, /* solution vector */
A[NDIM][NDIM] = {{0}}, /* adjunct matrix */
INV[NDIM][NDIM] = {{0}}, /* inverse of a */
det = 0.0, /* determinant of a */
X[NDIM] = {0}, /* linear system roots */
B[NDIM] = {0}; /* product A * c */
i = j = cof_1 = cof_2 = cof_3 = cof_4 = 0; /* initialize */
for (i = 0; i < NDIM; i++) { /* validate input */
printf ("Eq[%d] coefficients and constant: ", i + 1);
if (scanf ("%f%f%f%f", &a[i][0], &a[i][1], &a[i][2], &c[i]) != 4) {
fprintf (stderr, "error: invalid coefficient.\n");
return 1;
}
}
putchar ('\n');
/* cofactor matrix from matrix of minors incorporating
* patterned +/- sign to maxtix of minors.
*/
for (i = 0; i < NDIM; i++)
for (j = 0; j < NDIM; j++) {
cof_1 = (i + 1) % NDIM;
cof_2 = (j + 1) % NDIM;
cof_3 = (i + 2) % NDIM;
cof_4 = (j + 2) % NDIM;
A[i][j] = (a[cof_1][cof_2] * a[cof_3][cof_4]) -
(a[cof_1][cof_4] * a[cof_3][cof_2]);
}
/* compute determinant */
for (i = 0; i < NDIM; i++)
det += a[i][0] * A[i][0];
/* transpose cofactor (adjunct) */
for (i = 1; i < NDIM; i++)
for (j = 0; j < i; j++)
swap (&A[i][j], &A[j][i]);
/* consistency check B */
for (i = 0; i < NDIM; i++)
for (j = 0; j < NDIM; j++)
B[i] += A[i][j] * c[j];
if (det == 0) { /* eliminate infinite and no solution cases */
if (B[0] == 0 && B[1] == 0 && B[2] == 0) {
printf ("The system is consistent and there "
"are infinitely many solutions.\n");
return 1;
} else
printf ("There is no solution to this system.\n");
return 2;
}
/* compute inverse of A (1/det * A) */
for (i = 0; i < NDIM; i++)
for (j = 0; j < NDIM; j++)
INV[i][j] = A[i][j] / det;
/* multiply constant vector by inverse */
for (i = 0; i < NDIM; i++)
for (j = 0; j < NDIM; j++)
X[i] += INV[i][j] * c[j];
if (X[0] == 0 && X[1] == 0 && X[2] == 0)
/* The case for trivial solution */
printf ("This system has a trivial solution.\n");
else
/* The case for unique solution */
printf ("This system has a unique solution.\n\n");
/* output solutions for linear system of eq. */
for (i = 0; i < NDIM; i++)
printf ("\t%c = %.2f\n", 'x' + i, X[i]);
return 0;
}
例子Use/Output
用于示例的 30-60-90 三角形的简单证明中的线性方程：
$ ./bin/matrixinv
Eq[1] coefficients and constant: 1 1 1 180
Eq[2] coefficients and constant: 5 -1 -1 0
Eq[3] coefficients and constant: 1 1 -1 0
This system has a unique solution.
x = 30.00
y = 60.00
z = 90.00
检查一下，如果您还有其他问题，请告诉我。