d1r4.cpp

数值计算c++源代码,包括各种算法。很有用的。

#include "iostream.h"
#include "math.h"

void pendag(double a[], double b[], double c[], double d[], double e[], double r[], double u[], int n)
{
	int k;
    double w[100], beta[100], alpha[100], cg[100], h[100];
    w[1] = c[1];
    beta[1] = 0.0;
    beta[2] = d[1] / w[1];
    alpha[1] = 0.0;
    alpha[2] = e[1] / w[1];
    alpha[n] = 0.0;
    alpha[n + 1] = 0.0;
    for (k = 2; k<=n; k++)
	{
        cg[k] = b[k] - a[k] * beta[k - 1];
        w[k] = c[k] - a[k] * alpha[k - 1] - cg[k] * beta[k];
        if (w[k] == 0.0) 
		{
			cout<< "    w[k]=0.0 in pendag"<<endl;
		}
        beta[k + 1] = (d[k] - cg[k] * alpha[k]) / w[k];
        alpha[k + 1] = e[k] / w[k];
    }
        h[1] = 0.0;
        h[2] = r[1] / w[1];
    for (k = 2; k<=n; k++)
	{
        h[k + 1] = (r[k] - a[k] * h[k - 1] - cg[k] * h[k]) / w[k];
    }
        u[n] = h[n + 1];
        u[n - 1] = h[n] - beta[n] * u[n];
    for (k = n - 2; k>=1; k--)
	{
        u[k] = h[k + 1] - beta[k + 1] * u[k + 1] - alpha[k + 1] * u[k + 2];
    }
}

void main()
{
    //program d1r4
    //driver program for routine pendag
    int i,j,l,n = 7;
    double a[8], b[8], c[8], d[8], e[8], r[8], u[8], a1[8][8], x[8];
    //输入已知的方程组的系数矩阵
    for (i = 1; i<=n; i++)
	{
        for (j = 1; j<=n; j++)
		{
			a1[i][j] = 0;
        }
    }
    a1[1][1] = 4; a1[1][2] = 1; a1[1][3] = 1;
    a1[2][1] = 1; a1[2][2] = 5; a1[2][3] = 2; a1[2][4] = 2;
    a1[3][1] = 1; a1[3][2] = 2; a1[3][3] = 6; a1[3][4] = 3; a1[3][5] = 3;
    a1[4][2] = 2; a1[4][3] = 3; a1[4][4] = 7; a1[4][5] = 4; a1[4][6] = 4;
    a1[5][3] = 3; a1[5][4] = 4; a1[5][5] = 8; a1[5][6] = 5; a1[5][7] = 5;
    a1[6][4] = 4; a1[6][5] = 5; a1[6][6] = 9; a1[6][7] = 6;
    a1[7][5] = 5; a1[7][6] = 6; a1[7][7] = 10;
    //输入已知的方程组的右端向量
    r[1] = 1;    r[2] = 2;    r[3] = 3;
    r[4] = 4;    r[5] = 5;    r[6] = 6;
    r[7] = 7;
    cout<<endl;
    cout<<"已知的方程组的右端向量"<<endl;
	cout.width(10);    cout<<r[1]<<endl;
	cout.width(10);    cout<<r[2]<<endl;
	cout.width(10);    cout<<r[3]<<endl;
	cout.width(10);    cout<<r[4]<<endl;
	cout.width(10);    cout<<r[5]<<endl;
	cout.width(10);    cout<<r[6]<<endl;
	cout.width(10);    cout<<r[7]<<endl;
    for (i = 3; i<=n; i++)
	{
        a[i] = a1[i][i - 2];
    }
    for (i = 2; i<=n; i++)
	{
        b[i] = a1[i][i - 1];
    }
    for (i = 1; i<=n-1; i++)
	{
        d[i] = a1[i][i + 1];
    }
    for (i = 1; i<=n-2; i++)
	{
        e[i] = a1[i][i + 2];
	}
    for (i = 1; i<=n; i++) 
	{
        c[i] = a1[i][i];
    }
    pendag(a, b, c, d, e, r, u, n);
    cout<<endl;
    cout<<"计算出的方程组的解"<<endl;
	cout.width(14);    cout<<u[1]<<endl;
	cout.width(14);    cout<<u[2]<<endl;
	cout.width(14);    cout<<u[3]<<endl;
	cout.width(14);    cout<<u[4]<<endl;
	cout.width(14);    cout<<u[5]<<endl;
	cout.width(14);    cout<<u[6]<<endl;
	cout.width(14);    cout<<u[7]<<endl;
    //将计算出的解乘以系数矩阵，以验证计算结果正确
    for (l = 1; l<=n; l++)
	{
        x[l] = 0.0;
        for (j = 1; j<=n; j++)
		{
            x[l] = x[l] + a1[l][j] * u[j];
        }
    }
    cout<<endl;
    cout<<"计算出的解乘以系数矩阵的结果"<<endl;
	cout.width(10);    cout<<x[1]<<endl;
	cout.width(10);    cout<<x[2]<<endl;
	cout.width(10);    cout<<x[3]<<endl;
	cout.width(10);    cout<<x[4]<<endl;
	cout.width(10);    cout<<x[5]<<endl;
	cout.width(10);    cout<<x[6]<<endl;
	cout.width(10);    cout<<x[7]<<endl;
}