当前位置:首页 » 编程语言 » 三对角追赶法c语言
扩展阅读
webinf下怎么引入js 2023-08-31 21:54:13
堡垒机怎么打开web 2023-08-31 21:54:11

三对角追赶法c语言

发布时间: 2022-10-28 14:42:59

A. 追赶法解三对角矩阵的问题

追赶法本质上就是不选主元的Gauss消去法, 只要所有顺序主子式都非零就可以用
对角占优性可以保证顺序主子式都非零, 同时保证数值稳定性
如果没有对角占优性有可能会出现算法中断或者即使不中断但精度也比较低的情况

B. 谁能帮我下载个csdn文档, 实现追赶法求解三对角矩阵方程组的C++源代码,这是题目

实现追赶法求解三对角矩阵方程组的C++源代码.doc

网盘下载:

C. 求c语言编程用追赶法求解三对角方程组2x1+x2=1,x1+3x2+x3=2,x2+x3+x4=2,2x3+x4=0悬10分

//2x1+x2=1,x1+3x2+x3=2,x2+x3+x4=2,2x3+x4=0
#include "stdio.h"
#include "math.h"
#define MAXNUM 10 //变量数量
int array[MAXNUM][MAXNUM]={{2,1,0,0,1},{1,3,1,0,2},{0,1,1,1,2},{0,0,2,1,0}};

int unuse_result[MAXNUM];
int GaussFun(int equ,int var,int result[])
{
int i,j,k,col,num1,num2;
int max_r,ta,tb,gcdtemp,lcmtemp;
int temp,unuse_x_num,unuse_index;
col=0;
for (k=0;k<equ && col<var;k++,col++)
{
max_r=k;
for (i=k+1;i<equ;i++)
{
if (abs(array[i][col])>abs(array[max_r][col]))
{
max_r=i;
}
}
if (max_r!=k)
{
for (j=k;j<var+1;j++)
{
temp=array[k][j];
array[k][j]=array[max_r][j];
array[max_r][j]=temp;
}
}
if (array[k][col]==0)
{
k--;
continue;
}
for (i=k+1;i<equ;i++)
{
if (array[i][col]!=0)
{
num1=abs(array[i][col]);
num2=abs(array[k][col]);
while (num2!=0)
{
temp=num2;
num2=num1%num2;
num1=temp;
}
gcdtemp=num1;
lcmtemp=(abs(array[i][col])*abs(array[k][col]))/gcdtemp;
ta=lcmtemp/abs(array[i][col]);
tb=lcmtemp/abs(array[k][col]);
if (array[i][col]*array[k][col]<0)
{
tb=-tb;
}
for (j=col;j<var+1;j++)
{
array[i][j]=array[i][j]*ta-array[k][j]*tb;
}
}
}
}
for (i=k;i<equ;i++)
{
if (array[i][col]!=0)
{
return -1;
}
}
if (k<var)
{
for (i=k-1;i>=0;i--)
{
unuse_x_num=0;
for (j=0;j<var;j++)
{
if (array[i][j]!=0 && unuse_result[j])
{
unuse_x_num++;
unuse_index=j;
}
}
if (unuse_x_num>1)
{
continue;
}
temp=array[i][var];
for (j=0;j<var;j++)
{
if (array[i][j]!=0 && j!=unuse_index)
{
temp=array[i][j]*result[j];
}
}
result[unuse_index]=temp/array[i][unuse_index];
unuse_result[unuse_index]=0;
}
return var-k;
}
for (i=var-1;i>=0;i--)
{
temp=array[i][var];
for (j=i+1;j<var;j++)
{
if (array[i][j]!=0)
{
temp-=array[i][j]*result[j];
}
}
if (temp%array[i][i]!=0)
{
return -2;
}
result[i]=temp/array[i][i];
}
return 0;
}
void main()
{
int i,type;
int equnum,varnum;
int result[MAXNUM];
equnum=3;
varnum=3;
type=GaussFun(equnum,varnum,result);
if (type==-1)
{
printf("该方程无解!\n");
}
else if(type==-2)
{
printf("该方程有浮点数解,无整数解!\n");
}
else if (type>0)
{
printf("该方程有无穷多个解!自由变量数量为%d\n",type);
for (i=0;i<varnum;i++)
{
if (unuse_result[i])
{
printf("x%d 是不确定的\n",i+1);
}
else
{
printf("x%d:%d\n",i+1,result[i]);
}
}
}
else
{
printf("该方程的解为:\n");
for (i=0;i<varnum;i++)
{
printf("x%d=%d\n",i+1,result[i]);
}
}
}

D. 带状矩阵分解及追赶法

如果矩阵的非零元素分布在对角线元素两侧,例如,对角线下侧有m1条次对角线为非零元素,对角线上侧有m2条次对角为非零元素,矩阵的其他次对角线上元素均为零元素,这样的矩阵,称为带宽为m1+m2+1的带状矩阵。连续性问题离散化后得到的线性方程组,其系数矩阵往往是带状矩阵。

若带状矩阵的各阶主子式不为零,带宽为m1+m2+1的矩阵A,可分解成带宽为m1+1的下三角阵L和带宽为m2+1的上三角阵U的乘积。

当m1=m2=1时,该带状矩阵称为三对角矩阵。三对角矩阵是最简单,然而又是常遇的带状矩阵。例如,在建立三次样条插值函数以及求解微分方程数值解等问题中,都会遇到求解这样的线性方程组Ax=d,其中系数矩阵A是一个阶数较高的三对角方阵:

地球物理数据处理基础

简记为A=[ai,bi,cin1(i=1,2,…,n),这里n指n阶方阵,1指半带宽(简称带宽)。右端项d=(d1,d2…,dnT。此时方程组Ax=d称为三对角方程组,通过对三对角矩阵的直接分解来求解的方法称为追赶法。

假设系数矩阵A=[ai,bi,cin1存在克洛脱(Crout)分解A=LU。由于矩阵A是带宽为1的三对角矩阵,所以L,U均为二对角阵,即

地球物理数据处理基础

利用矩阵乘法,比较A=LU两端对应的元素,便可得如下确定矩阵L与U元素mi,li,ui的计算公式:

地球物理数据处理基础

对三对角方程组的系数矩阵A进行分解以后,求解Ax=d的问题就转化为求解两个特殊的二对角方程组:

地球物理数据处理基础

容易得到递推的计算式分别为

地球物理数据处理基础

人们常把下标由小到大求解 的过程形象化地称为“追”的过程,而将下标由大到小求解 的过程称为“赶”的过程。因此,上述方法通常称为追赶法。

从上面的讨论,我们看到追赶法的实质就是克洛脱分解法在求解三对角方程组中的应用。这时,由于矩阵A的特殊性,在求解运算时将系数矩阵中大量零元素撇开,仅在三条对角线上进行,从而简化计算公式。

三对角线方程组的系数矩阵A=[ai,bi,cin1应满足什么条件才能进行克洛脱分解?下面给出一个充分条件:

在三对角矩阵A=[ai,bi,cin1中,若

(1)ai≠0(i=2,3,…n)且ci≠0(i=1,2,…,n-1),

(2)按行对角占优,即

则Α非奇异且存在唯一的Α=LU分解。

从上面的计算公式得知,用追赶法求解三对角方程组的计算量小,仅需5n-4次乘除法运算,且计算稳定、存储量小。

[例]用追赶法求解三对角方程组

解:(1)先分解A=LU,求L和U的元素,

m2=a2=-1,m3=a3=-2,m4=a4=-3,l1=b1=2

l2=b2-m2u1=5/2,l3=b3-m3u2=12/5,l4=b4-m4u3=5/4

u1=c1/l1=-1/2,u2=c2/l2=-4/5,u3=c3/l3=-5/4

(2)求解Ly=d,即追的过程,

地球物理数据处理基础

(3)求解Ux=y,即赶的过程,

x4=y4=2,x3=y3-u3x4=3,x2=y2-u2x3=4,x1=y1-u1x2=5

所以方程组的解为x=(5,4,3,2)T

E. 求c语言编程

// 程序7.5 追赶法求解三对角方程组

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

#define MAX_n 100
#define PRECISION 0.0000001

//解输出
void SulutionOutput(float x[],int n)
{
int i;
for(i=1;i<=n;++i)
printf("\nx[%d]=%f",i,x[i]);
}

//三对角方程组元素输入
void TriDiagonalMatrixInput(float a[],float b[],float c[],float f[],int n)
{
int i;
printf("Input b[1],c[1],f[1]:");
scanf("%f %f %f",&b[1],&c[1],&f[1]);
for(i=2;i<n;++i)
{
printf("Input a[%d],b[%d],c[%d],f[%d]:",i,i,i,i);
scanf("%f %f %f %f",&a[i],&b[i],&c[i],&f[i]);
}
printf("Input a[%d],b[%d],f[%d]:",n,n,n);
scanf("%f %f %f",&a[n],&b[n],&f[n]);
}

//三对角方程组求解 — 追赶法
int Z_G_method(float a[],float b[],float c[],float f[],int n)
{
int i;
//
c[1]/=b[1];
for(i=2;i<n;++i)
c[i]/=(b[i]-a[i]*c[i-1]);
//
f[1]/=b[1];
for(i=2;i<=n;++i)
f[i]=(f[i]-a[i]*f[i-1])/(b[i]-a[i]*c[i-1]);
//
for(i=n-1;i>0;--i)
f[i]-=c[i]*f[i+1];
}

void main()
{
int n;
float a[MAX_n],b[MAX_n],c[MAX_n],f[MAX_n];
printf("\nInput n=");
scanf("%d",&n);
TriDiagonalMatrixInput(a,b,c,f,n);
Z_G_method(a,b,c,f,n);
SulutionOutput(f,n);
}

/*
运行实例:

Input n=4
Input b[1],c[1],f[1]:2 1 1
Input a[2],b[2],c[2],f[2]:1 3 1 2
Input a[3],b[3],c[3],f[3]:1 4 1 3
Input a[4],b[4],f[4]:1 5 4

x[1]=0.294118
x[2]=0.411765
x[3]=0.470588
x[4]=0.705882
*/

F. 用matalb编写程序,用追赶法求解三对角线性方程组:

function x=Trid(a,b,c,d)
% 追赶法求解三对角的线性方程组 Ax=d
% b为主对角线元素,a,c分别为次对角线元素,d为右端项
% A=[ b1 c1
% a2 b2 c2
% ......
% a_(n-1) b_(n-1) c_(n-1)
% a_(n) b_(n) ]
% b=[b1...b_(n)]
% a=[0 a2...a_(n)]
% c=[c1...c_(n-1)]

n=length(b);

u(1)=b(1);
for i=2:n
l(i)=a(i)/u(i-1);
u(i)=b(i)-l(i)*c(i-1);
end

y(1)=d(1);
for i=2:n
y(i)=d(i)-l(i)*y(i-1);
end

x(n)=y(n)/u(n);
for i=n-1:-1:1
x(i)=(y(i)-c(i)*x(i+1))/u(i);
end

G. 我在visual c++ 2005中编的追赶法在debug下与release版本下为何运行结果不同

release经过的代码优化.
你可以在项目编译的配置里关掉一些的优化选项试试看.

H. 什么是追赶法,其通用程序是什么

function x=zhuiganfa
%首先说明:追赶法是适用于三对角矩阵的线性方程组求解的方法,并不适用于其他类型矩阵。
%定义三对角矩阵A的各组成单元。方程为Ax=d
% b为A的对角线元素(1~n),a为-1对角线元素(2~n),c为+1对角线元素(1~n-1)。
% A=[2 -1 0 0
% -1 3 -2 0
% 0 -2 4 -3
% 0 0 -3 5]
a=[0 -1 -2 -3];c=[-1 -2 -3];b=[2 3 4 5];d=[6 1 -2 1];
n=length(b);
u0=0;y0=0;a(1)=0;
%“追”的过程
L(1)=b(1)-a(1)*u0;
y(1)=(d(1)-y0*a(1))/L(1);
u(1)=c(1)/L(1);
for i=2:(n-1)
L(i)=b(i)-a(i)*u(i-1);
y(i)=(d(i)-y(i-1)*a(i))/L(i);
u(i)=c(i)/L(i);
end
L(n)=b(n)-a(n)*u(n-1);
y(n)=(d(n)-y(n-1)*a(n))/L(n);
%“赶”的过程
x(n)=y(n);
for i=(n-1):-1:1
x(i)=y(i)-u(i)*x(i+1);