线性方程组求解的几种常用方法-c++代码实现

写在博客之前


  1. 代码可能比较杂乱无章,初学者,大神勿喷。
  2. 代码后会附上算法伪码,可以对照看看。

Doolittle直接分解算法

算法伪码:
A[n][n]—系数矩阵b[n] —常数项矩阵
L[n][n]U[n][n]y[n]x[n]L、U矩阵与解向量

Step1 L赋值:对角线元素赋值1,上三角元素赋值0;
Step2 U赋值:下三角元素赋值为0;
Step3 For k=1 To n Do Step4、Step5 /开始分解 /
Step4
For j=k To n Do /求 U中第 k 行元素/
u[k][j] = a[k][j]
For r = 1 To k-1 Do
u[k][j] = u[k][j] - l[k][r]*u[r][j]
EndFor r
EndFor j
Step5
For i=k+1 To n Do /求 L 中第 k 列元素/
l[i][k] = a[i][k];
For r = 1 To k-1 Do
l[i][k] = l[i][k] - l[i][r]*u[r][k];
EndFor r
l[i][k] = l[i][k] / u[k][k]
EndFor i
EndFor k /结束分解 /
Step6
For i = 1 To n Do /由 Ly=b 求出y[i] /
y[i] = b[i]
For j = 1 To i-1 Do
y[i] = y[i] - l[i][j]*y[j];
EndFor j
EndFor i
Step7
For i = n Downto 1 Do /由 Ux=y 求出 x[i] /
x[i] = y[i];
For j = i+1 To n Do
y[i] = y[i] - u[i][j]*x[j]
EndFor j
x[i] = y[i]/u[i][i]
EndFor i


void Doolittle(int n,double *A,double *b)//n为阶数 A为系数矩阵 b为常数矩阵
{
    double *L = new double[n*n];//开辟L矩阵空间
    double *U = new double[n*n];//开辟U矩阵空间
    double *y = new double[n];//开辟y矩阵空间
    for (int i = 0; i < n; i++)
    {
        for (int j = 0; j < n; j++)
        {
            *(U + i*n + j) = 0;//暂时全部赋值为0
            if (i==j)
            {
                *(L + i*n + j) = 1;//对角线赋值为1
            }
            else
            {
                *(L + i*n + j) = 0;//其他暂时赋值为0
            }
        }
    }
    for (int k = 0; k < n; k++)//计算u和l矩阵的值
    {
        for (int j = k; j < n; j++)
        {
            *(U + k*n + j) = *(A + k*n + j);//第一行
            for (int r = 0; r < k; r++)//接下来由L的前一列算u的下一行
            {
                *(U + k*n + j) = *(U + k*n + j) - (*(L + k*n + r)*(*(U + r*n + j)));
            }
        }
        for (int i = k+1; i < n; i++)//计算L的列
        {
            *(L + i*n + k) = *(A + i*n + k);
            for (int r = 0; r < k; r++)
            {
                *(L + i*n + k) = *(L + i*n + k) - (*(L + i*n + r)*(*(U + r*n + k)));
            }
            *(L + i*n + k) = *(L + i*n + k) / (*(U + k*n + k));
        }
    }
    for (int i = 0; i < n; i++)//由Ly=b算y
    {
        *(y + i) = *(b + i);
        for (int j = 0; j < i; j++)
        {
            *(y + i) = *(y + i) - *(L + i*n + j)*(*(y + j));

        }
    }
    for (int i = n-1; i >
  • 14
    点赞
  • 99
    收藏
    觉得还不错? 一键收藏
  • 7
    评论
评论 7
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值