算法经典VC++数值分析

第1章 线性代数方程组的解法 1.1全主元高斯-约当(Gauss-Jordan)消去法 1.2LU分解法 1.3追赶法 1.4五对角线性方程组解法 1.5线性方程组解的迭代改善 1.6范德蒙(Vandermonde)方程组解法 1.7托伯利兹(Toeplitz)方程组解法 1.8奇异值分解 1.9线性方程组的共轭梯度法 1.1对称方程组的乔列斯基(Cholesky)分解法 1.11矩阵的QR分解 1.12松弛迭代法 第2章 插值 2.1拉格朗日插值 2.2有理函数插值 2.3三次样条插值 2.4有序表的检索法 2.5插值多项式 2.6二元拉格朗日插值 2.7双三次样条插值 第3章 数值积分 3.1梯形求积法 3.2辛普森(Simpson)求积法 3.3龙贝格(Romberg)求积法 3.4反常积分 3.5高斯(Gauss)求积法 3.6三重积分 第4章 特殊函数 4.1Γ函数、贝塔函数、阶乘及二项式系数 4.2不完全Γ函数、误差函数 4.3不完全贝塔函数 4.4零阶、一阶和任意整数阶的第一、二类贝塞尔函数 4.5零阶、一阶和任意整数阶的第一、二类变形贝塞尔函数 4.6分数阶第一类贝塞尔函数和变形贝塞尔函数 4.7指数积分和定指数积分 4.8连带勒让德函数 附录 第5章 函数逼近 5.1级数求和 5.2多项式和有理函数 5.3切比雪夫逼近 5.4积分和导数的切比雪夫逼近 5.5用切比雪夫逼近求函数的多项式逼近 第6章 随机数 6.1均匀分布随机数 6.2变换方法——指数分布和正态分布随机数 6.3舍选法——Γ分布、泊松分布和二项式分布随机数 6.4随机位的产生 6.5蒙特卡罗积分法 第7章 排序 7.1直接插入法和Shell方法 7.2堆排序 7.3索引表和等级表 7.4快速排序 7.5等价类的确定 附录 第8章 特征值问题 8.1对称矩阵的雅可比变换 8.2变实对称矩阵为三对角对称矩阵 8.3三对角矩阵的特征值和特征向量 8.4变一般矩阵为赫申伯格矩阵 8.5实赫申伯格矩阵的QR算法 第9章 数据拟合 9.1直线拟合 9.2线性最小二乘法 9.3非线性最小二乘法 9.4绝对值偏差最小的直线拟合 第1章 方程求根和非线性方程组的解法 1.1图解法 1.2逐步扫描法和二分法 1.3割线法和试位法 1.4布伦特(Brent)方法 1.5牛顿-拉斐森(Newton-Raphson)法 1.6求复系数多项式根的拉盖尔(Laguerre)方法 1.7求实系数多项式根的贝尔斯托(Bairstou)方法 1.8非线性方程组的牛顿-拉斐森方法 第11章 函数的极值和最优化 11.1黄金分割搜索法 11.2不用导数的布伦特(Brent)法 11.3用导数的布伦特(Brent)法 11.4多元函数的下山单纯形法 11.5多元函数的包维尔(Powell)法 11.6多元函数的共轭梯度法 11.7多元函数的变尺度法 11.8线性规划的单纯形法 第12章 傅里叶变换谱方法 12.1复数据快速傅里叶变换算法 12.2实数据快速傅里叶变换算法(一) 12.3实数据快速傅里叶变换算法(二) 12.4快速正弦变换和余弦变换 12.5卷积和逆卷积的快速算法 12.6离散相关和自相关的快速算法 12.7多维快速傅里叶变换算法 第13章 数据的统计描述 13.1分布的矩——均值、平均差、标准差、方差、斜差和峰态 13.2中位数的搜索 13.3均值与方差的显著性检验 13.4分布拟合的X2检验 13.5分布拟合的K-S检验法 第14章 解常微分方程组 14.1定步长四阶龙格-库塔(Runge-Kutta)法 14.2自适应变步长的龙格-库塔法 14.3改进的中点法 14.4外推法 第15章 两点边值问题的解法 15.1打靶法(一) 15.2打靶法(二) 15.3松弛法 第16章 偏微分方程的解法 16.1解边值问题的松弛法 16.2交替方向隐式方法(ADI)

#include 
#include 
#include 
#include 
#include 
#include 
#include 

int sgn（double x）
{
	int temp;
	if（x>0） temp=1;
	if（x=0） temp=0;
	if（x<0） temp=-1;
	return temp;
}

double bessj1（double x）
{
	double p1p2p3p4p5q1q2q3q4q5;
	double r1r2r3r4r5r6s1s2s3s4s5s6;
	double bbbcccaaatempaxxxzy;
    r1 = 72362614232.0;           r2 = -7895059235.0;
    r3 = 242396853.1;             r4 = -2972611.439;
    r5 = 15704.4826;              r6 = -30.16036606;
    s1 = 144725228442.0;          s2 = 2300535178.0;
    s3 = 18583304.74;             s4 = 99447.43394;
    s5 = 376.9991397;             s6 = 1.0;
    p1 = 1.0;                     p2 = 0.00183105;
    p3 = -0.00003516396496;       p4 = 0.000002457520174;
    p5 = -0.000000240337019;
    q1 = 0.04687499995;           q2 = -0.0002002690873;
    q3 = 0.000008449199096;       q4 = -0.00000088228987;
    q5 = 0.000000105787412;
    if （fabs（x） < 8.0 ）
	{
        y = x*x;
        aaa = r1 + y * （r2 + y * （r3 + y * （r4 + y * （r5 + y * r6））））;
        bbb = s1 + y * （s2 + y * （s3 + y * （s4 + y * （s5 + y * s6））））;
        temp = x * aaa / bbb;
	}
    else
	{
        ax = fabs（x）;
        z = 8.0 / ax;
        y = z*z;;
        xx = ax - 2.356194491;
        aaa = p1 + y * （p2 + y * （p3 + y * （p4 + y * p5）））;
        bbb = q1 + y * （q2 + y * （q3 + y * （q4 + y * q5）））;
        ccc = sqrt（0.636619772 / ax）;
        temp = ccc * （cos（xx） * aaa - z * sin（xx） * bbb * sgn（x））;
	}  
	return temp;
}

double bessy1（double x）
{
	double p1p2p3p4p5q1q2q3q4q5;
	double r1r2r3r4r5r6s1s2s3s4s5s6s7;
	double bbbcccaaatempxxzy;
    p1 = 1.0;                     p2 = 0.00183105;
    p3 = -0.00003516396496;       p4 = 0.000002457520174;
    p5 = -0.000000240337019;
    q1 = 0.04687499995;           q2 = -0.0002002690873;
    q3 = 0.000008449199096;       q4 = -0.00000088228987;
    q5 = 0.000000105787412;
    r1 = -4900604943000.0;        r2 = 1275274390000.0;
    r3 = -51534381390.0;          r4 = 734926455.1;
    r5 = -4237922.726;            r6 = 8511.937935;
    s1 = 24995805700000.0;        s2 = 424441966400.0;
    s3 = 3733650367.0;            s4 = 22459040.02;
    s5 = 102042.605;              s6 = 354.9632885;
    s7 = 1.0;
    if （x < 8.0）
	{
        y = x*x;
        aaa = r1 + y * （r2 + y * （r3 + y * （r4 + y * （r5 + y * r6））））;
        bbb = s4 + y * （s5 + y * （s6 + y * s7））;
        bbb = s1 + y * （s2 + y * （s3 + y * bbb））;
        ccc = bessj1（x） * log（x） - 1.0 / x;
        temp = x * aaa / bbb + 0.636619772 * ccc;
	}
    else
	{
        z = 8.0 / x;
        y = z*z;
        xx = x - 2.356194491;
        aaa = sqrt（0.636619772 / x）;
        bbb = p1 + y * （p2 + y * （p3 + y * （p4 + y * p5）））;
        ccc = q1 + y * （q2 + y * （q3 + y * （q4 + y * q5）））;
        temp = aaa * （sin（xx） * bbb + z * cos（xx） * ccc）;