码迷,mamicode.com
首页 > 其他好文 > 详细

matlab中求逆矩阵的高斯消元法实现的代码

时间:2014-09-16 18:50:30      阅读:282      评论:0      收藏:0      [点我收藏+]

标签:io   strong   for   sp   cti   代码   on   c   amp   

function qiuni =INV_GET(a)
N=length(a);
M=eye(N);
%得到上三角矩?
for i=1:N
max=a(i,i);
A=i;
for j=i+1:N
if(abs(a(j,i))>abs(max))%找最大值
max=a(j,i);
A=j;
end
end
for m=1:N
temp1=a(i,m);%交换最大值值所在的行和当前行
a(i,m)=a(A,m);
a(A,m)=temp1;
temp2=M(i,m);
M(i,m)=M(A,m);
M(A,m)=temp2;
end
for k=i+1:N
jiaquan=a(k,i)/a(i,i);%其他行加权得到上三角
for n=1:N
a(k,n)=a(k,n)-jiaquan*a(i,n);
M(k,n)=M(k,n)-jiaquan*M(i,n);
end
end

end
temp3=a;
%得到单位矩阵
for s=N:-1:1
xishu1 =a(s,s);
for p=s:N
a(p,s)=a(p,s)/xishu1;
end
for q=1:N
M(s,q)=M(s,q)/xishu1;
end
for w=s-1:-1:1
xishu=a(w,s);
a(w,s)=0;
for t=1:N
M(w,t)=M(w,t)-xishu*M(s,t);
end
end
end
temp2=a;
qiuni=M;

 参考的c语言程序和注释

#include<stdio.h>

#include<math.h>

long int const N=1000;     //定义矩阵最大为1000阶

int n;              //n表示矩阵的行数和列数。

double juzhen[N][N]; //定义一个1000阶矩阵

double danwei[N][N];   //定义一个单位矩阵

bool zhaozuidazhi(int s)//定义一个布尔型从s行到n行选择最大的元素作为主元的函数。

{

    int i,j,A;

    double mas,temp;   //temp为中间变量,实现行交换

    mas=fabs(juzhen[s][s]);    //数学函数:fabs  功能:求浮点数juzhen[s][s]的绝对值

    //计算juzhen[s][s], 当juzhen[s][s]不为负时返回juzhen[s][s],否则返回-juzhen[s][s]

    A=s;

    for(i=s+1;i<n;i++)

    {

       if(mas<fabs(juzhen[i][s]))

       {

           mas=fabs(juzhen[i][s]);

           A=i;

       }

    }

    if(mas==0)

       return false;

    //交换两行

    for(j=0;j<n;j++)

    {

       temp=juzhen[s][j];

       juzhen[s][j]=juzhen[A][j];

       juzhen[A][j]=temp;

      

       temp=danwei[s][j];

       danwei[s][j]=danwei[A][j];

       danwei[A][j]=temp;

    }

    return true;

}

void jisuan(int s)  //消元计算

{

    int i,j;

    double mas=juzhen[s][s],r;

    for(i=s+1;i<n;i++)

    {

       r=juzhen[i][s]/mas;

       for(j=0;j<n;j++)

       {    //利用消元计算方阵使之成为上三角矩阵,最后使主对角线上的元素相乘就是最终结果 

           juzhen[i][j]=juzhen[i][j]-juzhen[s][j]*r;

           danwei[i][j]=danwei[i][j]-danwei[s][j]*r;

       }

    }

}

void HH(int s)

{

    int i,j;

    double mas;

    mas=juzhen[s][s];

    for(i=s;i<n;i++)

       juzhen[s][i]=juzhen[s][i]/mas;

    for(i=0;i<n;i++)

       danwei[s][i]=danwei[s][i]/mas;

    for(i=s-1;i>=0;i--)

    {

       mas=juzhen[i][s];

       juzhen[i][s]=0;

        for(j=0;j<n;j++)

           danwei[i][j]=danwei[i][j]-mas*danwei[s][j];

    }

}

int main()

{

    int i,j;

zl: printf("第一行输入矩阵的阶数,然后回车输入原始矩阵:\n");

    scanf("%d",&n);

    for(i=0;i<n;i++)

    {

       for(j=0;j<n;j++)

       {

           if(i==j) danwei[i][j]=1;

           else danwei[i][j]=0;

           scanf("%lf",&juzhen[i][j]);//C语言的函数参数是传值而不是传引用的,

           //通常函数无法修改和操作参数.scanf可以修改参数的原因是传给scanf的

           // 那个参数是一个指针,scanf通过传过来的指针来修改指针指向的内容。所以,

           //把x取地址再传过去。&是取juzhen[i][j]的地址,%lf是说juzhen[i][j]是一个双精度浮点型.

       }

    }

   

    //将原始矩阵转成上三角矩阵

    for(i=0;i<n-1;i++)

    {

       //每一次选一个绝对值最大的值作为a[i]即主元

       if(zhaozuidazhi(i)) //函数调用

       {

           jisuan(i);   //函数调用

 

       }

       else

       {

           printf("该矩阵不可逆");//高斯消去法不能将原始矩阵化为三角形的格式,

                                 //那就代表原始矩阵是一个不可逆的矩阵,

           return 0;

       }

    }

    //将原始矩阵转化成单位矩阵

    for(i=n-1;i>=0;i--)

    {

       HH(i);

    }

         printf("\n原始矩阵的逆矩阵为:\n");//输出逆矩阵语句

         for(i=0;i<n;i++)

         {    

             for(j=0;j<n;j++)

                printf("%.2lf ",danwei[i][j]);//%.2lf在printf()里表示按下面格式

             // 输出一个long double型数:整数部分输出全部输出小数部分输出2位,

             //没有小数的输出两个0,不足两位的后面补0,大于两位的截短到两位。

              //测试结果会出现两位小数

     

             printf("\n");

         }

         goto zl;

}

matlab中求逆矩阵的高斯消元法实现的代码

标签:io   strong   for   sp   cti   代码   on   c   amp   

原文地址:http://www.cnblogs.com/lianjiehere/p/3975455.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!