追赶法.cpp

追赶法求解三对角矩阵

//追赶法Ｃ语言源程序

///////////////////////////
//微分方程为
//du/dx=u-2x/u;
//u(0)=1
///////////////////////////

//#include "stdafx.h"

#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <time.h>
#include <iomanip.h>
#include <iostream.h>
#include <math.h>
#include <fstream.h>
#include <dos.h>
#include <conio.h>


const int    N=5;     //方程阶数加1	　　



void main(void){

double a[N]={0.0, 0.0, 90.860,  -67.590,  46.260};
double b[N]={0.0, 136.01, 98.810,  132.01,  177.17};
double c[N]={0.0, 90.860,  -67.590,  46.260, 0};
double d[N]={0.0, -33.254, 49.709,  28.067,  -7.3244};
double x[N]={0.0, 0.0, 0.0, 0.0, 0.0};
double r[N]={0.0, 0.0, 0.0, 0.0, 0.0};
double y[N]={0.0, 0.0, 0.0, 0.0, 0.0};
double q;
int    k;

//输出结果

ofstream fout;
fout.open("123.txt",ios::app|ios::out);

r[1]=c[1]/b[1];
y[1]=d[1]/b[1];

for(k=2; k<N-1; k++){
q=b[k]-r[k-1]*a[k];
r[k]=c[k]/q;
y[k]=(d[k]-y[k-1]*a[k])/q;
}




y[N-1]=(d[N-1]-y[N-2]*a[N-1])/(b[N-1]-r[N-2]*a[N-1]);
x[N-1]=y[N-1];
for(k=N-2;k>=1;k--)
x[k]=y[k]-r[k]*x[k+1];
for(k=1;k<N;k++)
fout<<"x["<<k<<"]="<<setprecision(12)<<x[k]<<endl;


fout<<endl;

fout<<"///////////////////////////////////////////////////////////";

fout.close();	

//  getch();
}