C语言求逆矩阵案例详解

一般求逆矩阵的方法有两种,伴随阵法和初等变换法。但是这两种方法都不太适合编程。伴随阵法的计算量大,初等变换法又难以编程实现。
适合编程的求逆矩阵的方法如下:

  1. 对可逆矩阵A进行QR分解:A=QR
  2. 求上三角矩阵R的逆矩阵
  3. 求出A的逆矩阵:A^(-1)=R^(-1)Q^(H)

以上三步都有具体的公式与之对应,适合编程实现。
C语言实现代码:

#include <stdio.h>
#include <math.h>

#define SIZE  8

double b[SIZE][SIZE]={0};//应该读作“贝尔塔”,注释中用B表示
double t[SIZE][SIZE]={0};//求和的那项
double Q[SIZE][SIZE]={0};//正交矩阵
double QH[SIZE][SIZE]={0};//正交矩阵的转置共轭
double R[SIZE][SIZE]={0};//
double invR[SIZE][SIZE]={0};//R的逆矩阵
double invA[SIZE][SIZE]={0};//A的逆矩阵,最终的结果
//={0};//
double matrixR1[SIZE][SIZE]={0};
double matrixR2[SIZE][SIZE]={0};

//double init[3][3]={3,14,9,6,43,3,6,22,15};
double init[8][8]={
    0.0938  ,  0.5201 ,   0.4424  ,  0.0196  ,  0.3912  ,  0.9493 ,   0.9899  ,  0.8256,
    0.5254  ,  0.3477 ,   0.6878  ,  0.3309 ,   0.7691  ,  0.3276 ,   0.5144  ,  0.7900,
    0.5303  ,  0.1500 ,   0.3592  ,  0.4243 ,   0.3968  ,  0.6713 ,   0.8843  ,  0.3185,
    0.8611  ,  0.5861 ,   0.7363  ,  0.2703 ,   0.8085  ,  0.4386 ,   0.5880  ,  0.5341,
    0.4849  ,  0.2621 ,   0.3947  ,  0.1971 ,   0.7551  ,  0.8335 ,   0.1548  ,  0.0900,
    0.3935  ,  0.0445 ,   0.6834  ,  0.8217 ,   0.3774  ,  0.7689 ,   0.1999  ,  0.1117,
    0.6714  ,  0.7549 ,   0.7040  ,  0.4299 ,   0.2160  ,  0.1673 ,   0.4070  ,  0.1363,
    0.7413  ,  0.2428 ,   0.4423  ,  0.8878 ,   0.7904  ,  0.8620 ,   0.7487  ,  0.6787
};
/*/
函数名:int main()
输入:
输出:
功能:求矩阵的逆 pure C language
     首先对矩阵进行QR分解之后求上三角矩阵R的逆阵最后A-1=QH*R-1,得到A的逆阵。
作者:HLdongdong
*//////////////////////////////////////////////////////////////////////
int main()
{
    int i;//数组  行
    int j;//数组  列
    int k;//代表B的角标
    int l;//数组  列
    double dev;
    double numb;//计算的中间变量
    double numerator,denominator;
    double ratio;
    /////////////////求B/////////////////
    for(i=0;i<SIZE;++i)
    {
        for(j=0;j<SIZE;++j)
        {
            b[j][i]=init[j][i];
        }
        for(k=0;k<i;++k)
        {
            if(i)
            {
                numerator=0.0;
                denominator=0.0;
                for(l=0;l<SIZE;++l)
                {
                    numerator+=init[l][i]*b[l][k];
                    denominator+=b[l][k]*b[l][k];
                }
                dev=numerator/denominator;
                t[k][i]=dev;
                for(j=0;j<SIZE;++j)
                {
                    b[j][i]-=t[k][i]*b[j][k];//t  init  =0  !!!
                }
            }
        }
    }
    ///////////////////对B单位化,得到正交矩阵Q矩阵////////////////////
    for(i=0;i<SIZE;++i)
    {
        numb=0.0;
        for(j=0;j<SIZE;++j)
        {
            numb+=(b[j][i]*b[j][i]);
        }
        dev=sqrt(numb);
        for(j=0;j<SIZE;++j)
        {
            Q[j][i]=b[j][i]/dev;
        }
        matrixR1[i][i]=dev;
    }
    /////////////////////求上三角R阵///////////////////////
    for(i=0;i<SIZE;++i)
    {
        for(j=0;j<SIZE;++j)
        {
            if(j<i)
            {
                matrixR2[j][i]=t[j][i];
            }
            else if(j==i)
            {
                matrixR2[j][i]=1;
            }
            else
            {
                matrixR2[j][i]=0;
            }
        }
    }
    mulMatrix(matrixR1,matrixR2,SIZE,SIZE,SIZE,R);
///////////////////////QR分解完毕//////////////////////////
    printf("QR分解:\n");
    printf("Q=\n");
    for(i=0;i<SIZE;++i)
    {
        for(j=0;j<SIZE;++j)
        {
            printf("%2.4f    ",Q[i][j]);
        //
        }
        printf("\n");
    }
    printf("R=\n");
    for(i=0;i<SIZE;++i)
    {
        for(j=0;j<SIZE;++j)
        {
            printf("%2.4f    ",R[i][j]);
        //
        }
        printf("\n");
    }
/////////////////////求R的逆阵//////////////////////////
    for(i=SIZE-1;i>=0;--i)
    {
        invR[i][i]=1/R[i][i];
        //R[i][i]=invR[i][i];
        if(i!=(SIZE-1))//向右
        {
            for(j=i+1;j<SIZE;++j)
            {
                invR[i][j]=invR[i][j]*invR[i][i];
                R[i][j]=R[i][j]*invR[i][i];
            }
        }
        if(i)//向上
        {
            for(j=i-1;j>=0;--j)
            {
                ratio=R[j][i];
                for(k=i;k<SIZE;++k)
                {
                    invR[j][k]-=ratio*invR[i][k];
                    R[j][k]-=ratio*R[i][k];
                }
            }
        }
    }

///////////////////////////////////////////////////////

    printf("inv(R)=\n");
    for(i=0;i<SIZE;++i)
    {
        for(j=0;j<SIZE;++j)
        {
            printf(" %2.4f  ",invR[i][j]);
        //
        }
        printf("\n");
    }
////////////////////结果和MATLAB差一个负号,神马鬼????????/////////////////////
/////////////////////求QH//////////////////////////
    for(i=0;i<SIZE;++i)//实矩阵就是转置
    {
        for(j=0;j<SIZE;++j)
        {
            QH[i][j]=Q[j][i];
        }
    }
///////////////////////求A的逆阵invA/////////////////////////////

    mulMatrix(invR,QH,SIZE,SIZE,SIZE,invA);

    printf("inv(A)=\n");
    for(i=0;i<SIZE;++i)
    {
        for(j=0;j<SIZE;++j)
        {
            printf(" %2.4f  ",invA[i][j]);
        //
        }
        printf("\n");
    }

///////////////////////结果与MATLAB的结果在千分位后有出入,但是负号都是对的^v^///////////////////////////
    return 0;
}

另附上矩阵乘法的子函数

/*/
函数名:void mulMatrix(double matrix1[SIZE][SIZE],double matrix2[SIZE][SIZE],int high1,int weight,int weight2,double mulMatrixOut[SIZE][SIZE])
输入:依次是 左矩阵,右矩阵,左矩阵高度,左矩阵宽度,右矩阵宽度,输出矩阵
输出:
功能:矩阵乘法
作者:HLdongdong
*//
void mulMatrix(double matrix1[SIZE][SIZE],double matrix2[SIZE][SIZE],int high1,int weight,int weight2,double mulMatrixOut[SIZE][SIZE])
{
    int i,j,k;
    for(i=0;i<high1;++i)
    {
        for(j=0;j<weight2;j++)
        {
            for(k=0;k<weight;++k)
            {
                mulMatrixOut[i][j]+=matrix1[i][k]*matrix2[k][j];
            }
        }
    }
}

到此这篇关于C语言求逆矩阵案例详解的文章就介绍到这了,更多相关C语言求逆矩阵内容请搜索我们以前的文章或继续浏览下面的相关文章希望大家以后多多支持我们!

(0)

相关推荐

  • java实现的n*n矩阵求值及求逆矩阵算法示例

    本文实例讲述了java实现的n*n矩阵求值及求逆矩阵算法.分享给大家供大家参考,具体如下: 先来看看运行结果: java版的写出来了,用的跟c语言相同的算法,然后看看能不能以后加个框做成程序: import java.math.*; import java.util.*; import java.text.*; public class matrix { static int map1[][]=new int [110][110]; static int just[][]=new int [11

  • C#计算矩阵的逆矩阵方法实例分析

    本文实例讲述了C#计算矩阵的逆矩阵方法.分享给大家供大家参考.具体如下: 1.代码思路 1)对矩阵进行合法性检查:矩阵必须为方阵 2)计算矩阵行列式的值(Determinant函数) 3)只有满秩矩阵才有逆矩阵,因此如果行列式的值为0(在代码中以绝对值小于1E-6做判断),则终止函数,报出异常 4)求出伴随矩阵(AdjointMatrix函数) 5)逆矩阵各元素即其伴随矩阵各元素除以矩阵行列式的商 2.函数代码 (注:本段代码只实现了一个思路,可能并不是该问题的最优解) /// <summary

  • Java实现的求逆矩阵算法示例

    本文实例讲述了Java实现的求逆矩阵算法.分享给大家供大家参考,具体如下: package demo; public class MatrixInverse { public static double Det(double [][]Matrix,int N)//计算n阶行列式(N=n-1) { int T0; int T1; int T2; double Num; int Cha; double [][] B; if(N>0) { Cha=0; B=new double[N][N]; Num=

  • C语言求逆矩阵案例详解

    一般求逆矩阵的方法有两种,伴随阵法和初等变换法.但是这两种方法都不太适合编程.伴随阵法的计算量大,初等变换法又难以编程实现. 适合编程的求逆矩阵的方法如下: 对可逆矩阵A进行QR分解:A=QR 求上三角矩阵R的逆矩阵 求出A的逆矩阵:A^(-1)=R^(-1)Q^(H) 以上三步都有具体的公式与之对应,适合编程实现. C语言实现代码: #include <stdio.h> #include <math.h> #define SIZE 8 double b[SIZE][SIZE]={

  • C语言指针数组案例详解

    指针与数组是 C 语言中很重要的两个概念,它们之间有着密切的关系,利用这种 关系,可以增强处理数组的灵活性,加快运行速度,本文着重讨论指针与数组之 间的联系及在编程中的应用. 1.指针与数组的关系 当一个指针变量被初始化成数组名时,就说该指针变量指向了数组.如: char str[20], *ptr; ptr=str; ptr 被置为数组 str 的第一个元素的地址,因为数组名就是该数组的首地址, 也是数组第一个元素的地址.此时可以认为指针 ptr 就是数组 str(反之不成立), 这样原来对数

  • C语言strtod()函数案例详解

    前言 网上有很多关于strtod()函数的文章,不过大部分都是用strtod()函数转换一个字符 char *str = "111.11"; char *target; double ret; ret = strtod(str, &target); 很少有转换字符串的这样的用法 char *p = "111.11 -2.22 Nan nan(2) inF 0X1.BC70A3D70A3D7P+6 1.18973e+4932zzz"; 本文主要参考strtod

  • C语言实现矩阵运算案例详解

    C语言实现矩阵运算 给定一个n×n的方阵,本题要求计算该矩阵除副对角线.最后一列和最后一行以外的所有元素之和.副对角线为从矩阵的右上角至左下角的连线. 输入格式: 输入第一行给出正整数n(1<n≤10):随后n行,每行给出n个整数,其间以空格分隔. 输出格式: 在一行中给出该矩阵除副对角线.最后一列和最后一行以外的所有元素之和. 输入样例: 4 2 3 4 1 5 6 1 1 7 1 8 1 1 1 1 1 输出样例: 35 #include <stdio.h> #include <

  • C语言 TerminateProcess函数案例详解

    TerminateProcess 顾名思义,就是终止进程的意思. 是WindowsAPI的函数, 示例代码如下: // Demo.cpp : 定义控制台应用程序的入口点. //终止进程Demo #include "stdafx.h" using namespace std; //@param:dwpid:指定需要关闭的进程pid int CloseProcess(DWORD dwpid) { HANDLE hProcess = OpenProcess(PROCESS_TERMINATE

  • C语言 bind()函数案例详解

    bind()函数介绍        在建立套接字文件描述符成功后,需要对套接字进行地址和端口的绑定,才能进行数据的接收和发送操作. 函数原型        bind()函数将长度为addlen的struct sockadd类型的参数my_addr与sockfd绑定在一起,将sockfd绑定到某个端口上,如果使用connect()函数则没有绑定的必要.绑定的函数原型如下: #include<sys/types.h> #include<sys/socket.h> int bind(in

  • C语言之快速排序案例详解

    快速排序:是对冒泡排序算法的一种改进. 它的基本思想是:通过一趟排序将要排序的数据分割成独立的两部分,其中一部分的所有数据都比另外一部分的所有数据都要小,然后再按此方法对这两部分数据分别进行快速排序,整个排序过程可以递归进行,以此达到整个数据变成有序序列. 例如有一个数字序列: 5 0 1 6 8 2 3 4 9 7 对其进行快速排序变为:0 1 2 3 4 5 6 7 8 9 思路如下:首先将要排序的序列的首个数字5定位比较数,这是一个参考对象! 然后的方法很简单:分别从序列的两端进行比较.先

  • C语言 动态分配数组案例详解

    目录 一维动态数组的创建: 二维数组的创建: 很多人在编写C语言代码的时候很少使用动态数组,不管什么情况下通通使用静态数组的方法来解决,在当初学习C语言的时候我就是一个典型的例子,但是现在发现这是一个相当不好的习惯,甚至可能导致编写的程序出现一些致命的错误.尤其对于搞嵌入式的人来所,嵌入式系统的内存是宝贵的,内存是否高效率的使用往往意味着嵌入式设备是否高质量和高性能,所以高效的使用内存对我们来说是很重要的.那么我们在自己编写C语言代码的时候就应该学会使用动态数组,这也就是我这篇博客要给大家讲的,

  • C语言 CRITICAL_SECTION用法案例详解

          很多人对CRITICAL_SECTION的理解是错误的,认为CRITICAL_SECTION是锁定了资源,其实,CRITICAL_SECTION是不能够"锁定"资源的,它能够完成的功能,是同步不同线程的代码段.简单说,当一个线程执行了EnterCritialSection之后,cs里面的信息便被修改,以指明哪一个线程占用了它.而此时,并没有任何资源被"锁定".不管什么资源,其它线程都还是可以访问的(当然,执行的结果可能是错误的).只不过,在这个线程尚未执

  • Go语言实现二维数组的2种遍历方式以及案例详解

    二维数组遍历的2种方式: package main import ( "fmt" ) func main() { //定义一个二维数组 var arr = [2][3]int{{1, 4, 3},{7, 5, 6}} //方式1. 用for循环来遍历 for i := 0; i < len(arr); i++ { for j := 0; j < len(arr[i]); j++ { fmt.Printf("%v ",arr[i][j]) } fmt.Pr

随机推荐