雅克比迭代法c語言
A. 雅可比迭代法的工作原理
雅克比迭代法就是眾多迭代法中比較早且較簡單的一種，其命名也是為紀念普魯士著名數學家雅可比。雅克比迭代法的計算公式簡單，每迭代一次只需計算一次矩陣和向量的乘法，且計算過程中原始矩陣A始終不變，比較容易並行計算。
概念
考慮線性方程組Ax =b時，一般當A為低階稠密矩陣時，用主元消去法解此方程組是有效方法。但是，對於由工程技術中產生的大型稀疏矩陣方程組(A的階數很高，但零元素較多，例如求某些偏微分方程數值解所產生的線性方程組)，利用迭代法求解此方程組就是合適的，在計算機內存和運算兩方面，迭代法通常都可利用A中有大量零元素的特點。雅克比迭代法就是眾多迭代法中比較早且較簡單的一種，其命名也是為紀念普魯士著名數學家雅可比。
折疊編輯本段迭代過程
迭代過程
首先將方程組中的系數矩陣A分解成三部分，即:A = L+D+U，如圖1所示，其中D為對角陣，L為下三角矩陣，U為上三角矩陣。
之後確定迭代格式，X^(k+1) =B*X^(k) +f，(這里^表示的是上標，括弧內數字即迭代次數)，如圖2所示，其中B稱為迭代矩陣，雅克比迭代法中一般記為J。(k = 0,1，......)
再選取初始迭代向量X^(0)，開始逐次迭代。
折疊編輯本段收斂性
設Ax= b，其中A=D+L+U為非奇異矩陣，且對角陣D也非奇異，則當迭代矩陣J的譜半徑ρ(J)<1時，雅克比迭代法收斂。
折疊編輯本段優缺點
雅克比迭代法的優點明顯，計算公式簡單，每迭代一次只需計算一次矩陣和向量的乘法，且計算過程中原始矩陣A始終不變，比較容易並行計算。然而這種迭代方式收斂速度較慢，而且占據的存儲空間較大，所以工程中一般不直接用雅克比迭代法，而用其改進方法。
折疊編輯本段程序實現示例
#include<stdio.h>
#include<math.h>
#include <stdlib.h>
main(){
float e=0.001,z,m,a[3][3]={5,2,1,-1,4,2,2,-3,10},b[3]={-12,20,3},x[3]={0,0,0},y[3];
int n=3,j,i,k=1;
while(1) {
for(i=0;i<3;i++) {
for(j=0;j<3;j++)
m=m+a[i][j]*x[j];
m=m-x[i]*a[i][i];
y[i]=(b[i]-m)/a[i][i];
m=0;
}
i=0;
while(i<3) {
z=fabs(x[i]-y[i]);
if(z>e)
break;
i++;
}
if(i!=3) {
for(i=0;i<3;i++)
x[i]=y[i];
k++;
}
else if(i==3)
break;
}
printf("%f %f %f ",y[0],y[1],y[2]);
}
B. 雅克比迭代法求解線性方程組的C語言程序
void Solve ( double dCoef [] , double dY [] , unsigned int iOrder , double dErr)
{//用Jacobi迭代法解方程組, dCoef[]系數陣, Y[]向量, iOrder給出方程階數, dErr給出精度
double res [Max]; //方程解
double res2[Max]; //保存上一階方程解
if ( Max < iOrder )
{
printf ("最多支持%d階方程組.", Max);
return;
}
for ( unsigned int i = 0 ; i < iOrder ; res2 [i++] = 0.0 ); //初始解向量 (0,0...)
while ( true )
{
bool bStopIterative = true;
for (unsigned int i = 0 ; i < iOrder ; ++i)
{
double dSum2 = 0;
for (unsigned int j = 0 ; j < iOrder ; j++)
{//求第二項
if ( j == i ) continue;
dSum2 += dCoef [i * iOrder + j] * res2 [j];
}
res[i] = 1/dCoef[i * iOrder + i] * ( dY[i] - dSum2 );
if ( abs ( res2[i] - res [i] ) > dErr )
bStopIterative = false;
}
if ( bStopIterative )
break;
for (unsigned int i = 0 ; i < iOrder ; i++ )
res2[ i ] = res[ i ];
}
//輸出結果
for (unsigned int i = 0 ; i < iOrder ; i++)
printf ("x%d = %lf\n", i+1 , res[i]);
}
int main(int argc, char* argv[])
{
double a[] =
{
8,-3,2,
4,11,-1,
2,1,4
};
double b[3] =
{
20,33,12
};
Solve ( a , b , 15 , 1e-10);
getchar();
return 0;
}