24小时热门版块排行榜    

查看: 3281  |  回复: 10

舒马诺

银虫 (初入文坛)

[求助] 带平方根的(LLT)Cholesky算法分解对称正定矩阵 c语言实现

新人,不懂规矩,见谅

大神们好,要求
1先判断任意矩阵A是否为正定对称矩阵,否则,返回输入错误
2若输入为正定对称矩阵,则将其进行带平方根的(LLT)Cholesky算法分解,即实现A=LL^T,其中L为下三角形矩阵。

大致就这意思,求助

定理
回复此楼
已阅   关注TA 给TA发消息 送TA红花 TA的回帖
回帖支持 ( 显示支持度最高的前 50 名 )

xiuyouxu

铁杆木虫 (职业作家)

【答案】应助回帖

感谢参与,应助指数 +1
建议找一本数值分析的书看一下,里面有具体的算法,我以前实现过,其他语言的,没有用c语言做过.
忘记自己,忘记一切烦恼(欢迎访问我的网站兆字节:http://www.mathbeta.com/)
2楼2012-05-03 21:28:11
已阅   关注TA 给TA发消息 送TA红花 TA的回帖
普通回帖

舒马诺

银虫 (初入文坛)

???????:
2?: Originally posted by xiuyouxu at 2012-05-03 21:28:11:
???????????????????????,?????о??????,?????????,?????????,?????c????????.

# include
# include
void main()
{
        float m,A[9];
float L[6];
        printf("请输入矩阵: \n ";
        scanf("%f %f %f\n%f %f %f\n%f %f %f\n",&A[0],&A[1],&A[2],&A[3],&A[4],&A[5],&A[6],&A[7],&A[8]);
        printf("请输入??许误差:m=";
scanf("%f",&m);
if
        A[0]>m&&(A[0]*A[4]-A[1]*A[3]>m)&&(A[6]*A[4]*A[2]+A[0]*A[7]*A[5]+A[1]*A[3]*A[8]-A[0]*A[4]*A[8]-A[1]*A[6]*A[5]-A[2]*A[3]*A[7]>m)&&(A[1]==A[3])&&(A[2]==A[6])&&(A[5]==A[7])
{
L[0]=sqrt(A[0]);
L[1]=A[3]/L[0];
L[3]=A[6]/L[0];
L[2]=sqrt(A[4]-L[1]*L[1]);
L[4]=(A[7]-L[3]*L[1])/L[2];
L[5]=sqrt(A[8]-L[3]*L[3]-L[4]*L[4]);
printf("所求矩阵为L=\n %f 0 0\n%f %f 0\n%f %f %f\n",L[0],L[1],L[2],L[3],L[4], L[5]);
}
else
printf("输入有误,请检查";
}
调试??行:
1>.\Debug\shiyan.exe.intermediate.manifest : general error c1010070: Failed to load and parse the manifest. {_~0p'1a@'7v par 1>Build log was saved at "file://e:\360data\????数???\桌???\shiyan\shiyan\Debug\BuildLog.htm"
1>shiyan - 1 error(s), 0 warning(s)
========== Rebuild All: 0 succeeded, 1 failed, 0 skipped ==========
工程无法建立

预期效果:
请输入矩阵:
1 2 3
2 4 5
3 5 6
请输入??许误差:m=1e-6
输入有误,请检查
请输入矩阵:
5 2 -4
2 1 -2
-4 -2 5
请输入??许误差:m=1e-6
所求矩阵L=
2.236068 0 0
0.894427 0.4472136 0
-1.788854 -0.894427 1




我的算法??行??通过啊,而且根本未能实现针对任??阶次的矩阵。。。求大神帮忙~
3楼2012-05-03 22:13:16
已阅   关注TA 给TA发消息 送TA红花 TA的回帖

舒马诺

银虫 (初入文坛)

引用回帖:
2楼: Originally posted by xiuyouxu at 2012-05-03 21:28:11:
建议找一本数值分析的书看一下,里面有具体的算法,我以前实现过,其他语言的,没有用c语言做过.

# include
# include
void main()
{
        float m,A[9];
float L[6];
        printf("请输入矩阵: \n ";
        scanf("%f %f %f\n%f %f %f\n%f %f %f\n",&A[0],&A[1],&A[2],&A[3],&A[4],&A[5],&A[6],&A[7],&A[8]);
        printf("请输入允许误差:m=";
scanf("%f",&m);
if
        A[0]>m&&(A[0]*A[4]-A[1]*A[3]>m)&&(A[6]*A[4]*A[2]+A[0]*A[7]*A[5]+A[1]*A[3]*A[8]-A[0]*A[4]*A[8]-A[1]*A[6]*A[5]-A[2]*A[3]*A[7]>m)&&(A[1]==A[3])&&(A[2]==A[6])&&(A[5]==A[7])
{
L[0]=sqrt(A[0]);
L[1]=A[3]/L[0];
L[3]=A[6]/L[0];
L[2]=sqrt(A[4]-L[1]*L[1]);
L[4]=(A[7]-L[3]*L[1])/L[2];
L[5]=sqrt(A[8]-L[3]*L[3]-L[4]*L[4]);
printf("所求矩阵为L=\n %f 0 0\n%f %f 0\n%f %f %f\n",L[0],L[1],L[2],L[3],L[4], L[5]);
}
else
printf("输入有误,请检查";
}

调试运行:
1>.\Debug\shiyan.exe.intermediate.manifest : general error c1010070: Failed to load and parse the manifest. {_~0p'1a@'7v par 1>Build log was saved at "file://e:\360data\重要数据\桌面\shiyan\shiyan\Debug\BuildLog.htm"
1>shiyan - 1 error(s), 0 warning(s)
========== Rebuild All: 0 succeeded, 1 failed, 0 skipped ==========
工程无法建立



失败了,而且达不到针对任意阶次矩阵的效果!
4楼2012-05-03 22:17:09
已阅   关注TA 给TA发消息 送TA红花 TA的回帖

xiuyouxu

铁杆木虫 (职业作家)

【答案】应助回帖

matlab里面直接用root函数就可以了, 下面是我写的c++的:
// 定义Matrix类(略)
// m*n阶0矩阵
void Matrix::zeros(int m,int n,double** a){
        for(int i=0;i                 for(int j=0;j                         a[j]=0;
                }
        }
}

// n为矩阵的阶
void Matrix::root(int n,double** A,double** L){
     zeros(n,n,L);
     for(int i=0;i              for(int j=0;j                      double sum=0;
                     for(int k=0;k                              sum+=L[k]*L[j][k];
                     }
                     L[j]=(A[j]-sum)/L[j][j];
             }
             double sum=0;
             for(int k=0;k                      sum+=L[k]*L[k];
             }
             L=sqrt(A-sum);// 显然 A-sum<0时不是正定矩阵
     }
}
忘记自己,忘记一切烦恼(欢迎访问我的网站兆字节:http://www.mathbeta.com/)
5楼2012-05-03 22:29:03
已阅   关注TA 给TA发消息 送TA红花 TA的回帖

舒马诺

银虫 (初入文坛)

引用回帖:
5楼: Originally posted by xiuyouxu at 2012-05-03 22:29:03:
matlab里面直接用root函数就可以了, 下面是我写的c++的:
// 定义Matrix类(略)
// m*n阶0矩阵
void Matrix::zeros(int m,int n,double** a){
        for(int i=0;i<m;i++){
                for(int j=0;j<n;j++){
                        a=0; ...

还是运行不通。。。
6楼2012-05-03 23:00:27
已阅   关注TA 给TA发消息 送TA红花 TA的回帖

xiuyouxu

铁杆木虫 (职业作家)

【答案】应助回帖

晕,这个回复框不能放代码啊,有一部分代码被替换掉了,代码里不能出现,会被替换掉
忘记自己,忘记一切烦恼(欢迎访问我的网站兆字节:http://www.mathbeta.com/)
7楼2012-05-03 23:09:59
已阅   关注TA 给TA发消息 送TA红花 TA的回帖

xiuyouxu

铁杆木虫 (职业作家)

看看这样行不行 \[i\]
忘记自己,忘记一切烦恼(欢迎访问我的网站兆字节:http://www.mathbeta.com/)
8楼2012-05-03 23:10:36
已阅   关注TA 给TA发消息 送TA红花 TA的回帖

xiuyouxu

铁杆木虫 (职业作家)

【答案】应助回帖

void Matrix::zeros(int m,int n,double** a){
        for(int i=0;i                 for(int j=0;j                         a\[i\][j]=0;
                }
        }
}

void Matrix::root(int n,double** A,double** L){
     zeros(n,n,L);
     for(int i=0;i              for(int j=0;j                      double sum=0;
                     for(int k=0;k                              sum+=L\[i\][k]*L[j][k];
                     }
                     L[j]=(A[j]\[i\]-sum)/L[j][j];
             }
             double sum=0;
             for(int k=0;k                      sum+=L\[i\][k]*L\[i\][k];
             }
             L=sqrt(A\[i\]\[i\]-sum);
     }
}

把上面的中括号前的反斜线去掉就行了
忘记自己,忘记一切烦恼(欢迎访问我的网站兆字节:http://www.mathbeta.com/)
9楼2012-05-03 23:12:39
已阅   关注TA 给TA发消息 送TA红花 TA的回帖

舒马诺

银虫 (初入文坛)

引用回帖:
7楼: Originally posted by xiuyouxu at 2012-05-03 23:09:59:
晕,这个回复框不能放代码啊,有一部分代码被替换掉了,代码里不能出现,会被替换掉

多谢高手帮忙了,弱弱的问一句能不能发到wuleileihappy@163.com呢?感激不尽
10楼2012-05-03 23:13:29
已阅   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 舒马诺 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[论文投稿] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +3 tqUTOxClQUMF 2026-09-28 3/150 2026-09-29 10:26 by JzYBbHIXWSrW
[考博] 售SCI一区T0P文章,我:8.O.55.1.O.54,科目齐全,可+急 +5 GDBe8tDZqE8z 2026-09-28 5/250 2026-09-29 09:52 by JzYBbHIXWSrW
[教师之家] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +4 GDBe8tDZqE8z 2026-09-28 4/200 2026-09-29 09:49 by JzYBbHIXWSrW
[找工作] 售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急 +5 IDs3scOF0tjC 2026-09-28 5/250 2026-09-29 09:24 by JzYBbHIXWSrW
[博后之家] 售SCI一区T0P文章,我:8.O.55.1.O54,科目全,可伽急 +5 IDs3scOF0tjC 2026-09-28 6/300 2026-09-29 09:18 by JzYBbHIXWSrW
[考研] 售SCI一区文章,我:8.O.551.O.5.4,科目全,可伽急 +3 ZYXdhzDAy9ZX 2026-09-28 3/150 2026-09-29 06:42 by ZvyPK8n6Nfki
[找工作] 售SCI一区T0P文章,我:8.O.55.1.O54,科目全,可伽急 +3 ZYXdhzDAy9ZX 2026-09-28 3/150 2026-09-29 06:34 by ZvyPK8n6Nfki
[找工作] 售SCI一区文章,我:8.O.551.O.5.4,科目全,可伽急 +3 ZYXdhzDAy9ZX 2026-09-28 3/150 2026-09-29 06:28 by ZvyPK8n6Nfki
[考博] 售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急 +3 ZYXdhzDAy9ZX 2026-09-28 3/150 2026-09-29 06:18 by ZvyPK8n6Nfki
[硕博家园] 售SCI一区T0P文章,我:8.O.55.1.O.54,科目齐全,可+急 +3 tqUTOxClQUMF 2026-09-28 3/150 2026-09-28 23:46 by b8fCuTHqckEt
[考研] 售SCI一区文章,我:8.O.55.1.O.54,科目齐全,可伽急 +3 GDBe8tDZqE8z 2026-09-28 3/150 2026-09-28 23:12 by b8fCuTHqckEt
[硕博家园] 售SCI一区T0P文章,我:8.O.55.1.O.54,科目齐全,可+急 +3 CfXuS1rDhLYN 2026-09-28 4/200 2026-09-28 22:55 by ez6fGg9abYaj
[公派出国] 售SCI一区文章,我:8O5.5.1.O5.4,科目全,可伽急 +4 IDs3scOF0tjC 2026-09-28 4/200 2026-09-28 22:55 by ez6fGg9abYaj
[有机交流] 同一个分子,一条来自文献,一条来自AI——不告诉你答案,你会选哪条? +3 tianxiaoxian 2026-09-28 6/300 2026-09-28 22:44 by tianxiaoxian
[找工作] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +3 IDs3scOF0tjC 2026-09-28 5/250 2026-09-28 22:32 by ez6fGg9abYaj
[公派出国] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +3 IDs3scOF0tjC 2026-09-28 4/200 2026-09-28 22:28 by ez6fGg9abYaj
[博后之家] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +3 IDs3scOF0tjC 2026-09-28 3/150 2026-09-28 22:21 by ez6fGg9abYaj
[论文投稿] 售SCI文章,我:8O5.5.1.O.54,科目齐全,可+急 +3 CfXuS1rDhLYN 2026-09-28 3/150 2026-09-28 19:08 by 9lS3ad5oOymn
[考研] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +3 IDs3scOF0tjC 2026-09-28 4/200 2026-09-28 18:51 by 9lS3ad5oOymn
[教师之家] 某top大学教授说“能够在市场中兑现的能力才是真能力”无比同意! +7 zju2000 2026-09-26 8/400 2026-09-28 16:42 by beefly
信息提示
请填处理意见