《数值计算课程设计》由会员分享，可在线阅读，更多相关《数值计算课程设计（40页珍藏版）》请在金锄头文库上搜索。
1、数值计算课程设计- 1 -1、 经典四阶龙格库塔法解一阶微分方程1.1、 算法说明龙 格 -库 塔 (Runge-Kutta)方 法 是 一 种 在 工 程 上 应 用 广 泛 的 高 精 度 单 步 算 法 。 由 于此 算 法 精 度 高 ， 采 取 措 施 对 误 差 进 行 抑 制 ， 所 以 其 实 现 原 理 也 较 复 杂 。 该 算 法 是 构 建在 数 学 支 持 的 基 础 之 上 的 。4 阶 龙 格 -库 塔 方 法 (RK4)可 模 拟 N=4 的 泰 勒 方 法 的 精 度 。 这 种 算 法 可 以 描 述 为 ，自 初 始 点 开 始 ， 利 用 下 面 的 计
2、 算 方 法 生 成 近 似 序 列0(,)ty(1-1)1.2、 经典四阶龙格库塔法解一阶微分方程算法流程图图 1-1 经典四阶龙格库塔法解一阶微分方程算法流程图经典四阶龙格库塔法解一阶微分方程- 2 -1.3、 经典四阶龙格库塔法解一阶微分方程程序调试图 1-2 经典四阶龙格库塔法解一阶微分方程程序调试1.4、 经典四阶龙格库塔法解一阶微分方程代码#include #include using namespace std;/f 为 函 数 的 入 口 地 址 ， x0、 y0 为 初 值 ， xn 为 所 求 点 ， step 为 计 算 次 数double Runge_Kuta( dou
3、ble (*f)(double x, double y), double x0, double y0, double xn, int step ) double k1,k2,k3,k4,result;double h=(xn-x0)/step;if(stepx0y0;coutab;/double x0=0,y0=1;double x,y,step;int i;coutstep;/step=0.1;cout.precision(10);for(i=0;imaxK #includeusing namespace std;void load();const N=20;float aNN;int m;
4、int main() int i,j;int c,k,n,p,r;float xN,lNN,s,d;coutm;coutfabs(aii)?j:i; /*找 列 最 大 元 素 */高斯列主元法解线性方程组- 8 -for(n=0;n=0;i-) d=0;for(j=i+1;jaij;数值计算课程设计- 9 -1122(,)(,)()kkkfpqfqxyJP12(,)()kkfF1kP3、 牛顿法解非线性方程组3.1、 算法说明设 已 知 。kP第 1 步 ： 计 算 函 数 （ 3-1）第 2 步 ： 计 算 雅 可 比 矩 阵 （ 3-2）第 3 步 ： 求 线 性 方 程 组 的 解 。
5、()()kkFP第 4 步 ： 计 算 下 一 点 （ 3-3）重 复 上 述 过 程 。3.2、 牛顿法解非线性方程组算法流程图图 3-1 算法流程图牛顿法解非线性方程组- 10 -20.54xy3.3、 牛顿法解非线性方程组算法程序调试图 3-2 牛 顿 法 解 非 线 性 方 程 组 算 法 程 序 调 试应用本程序解方程组， 初始近似值 x0,y0 分别为 2.00 和 0.25，经过3 次迭代求出 X(1)=1.900691 和 X(2)=0.311213。图 3-2 牛 顿 法 解 非 线 性 方 程 组 算 法 程 序 运 行 结 果数值计算课程设计- 11 -3.4、 牛顿法解
6、非线性方程组算法程序代码#include#include#define N 2 / 非 线 性 方 程 组 中 方 程 个 数 、 未 知 量 个 数 #define epsilon 0.0001 / 差 向 量 1 范 数 的 上 限#define max 10 /最 大 迭 代 次 数using namespace std;const int N2=2*N;int main()void ff(float xxN,float yyN);void ffjacobian(float xxN,float yyNN);void inv_jacobian(float yyNN,float invNN)
7、;void newdim(float x0N, float invNN,float y0N,float x1N);float x0N=2.0,0.25,y0N,jacobianNN,invjacobianNN,x1N,errornorm;int i,iter=0;cout0;i-) for (k=i-1;k=0;k-)L=-augki /augii;for(j=N2-1;j=0;j-)augkj=augkj+L*augij;for (i=0;i=0;i-)for(j=N2-1;j=0;j-)augij=augij/augii;for (i=0;i#includeusing namespace
8、std;#define f(x) sin(x*x) /举例函数#define epsilon 0.0001 /精度#define MAXREPT 10 /迭代次数,到最后仍达不到精度要求,则输出 T(m=10).double Romberg(double aa, double bb) /aa,bb 积分上下限int m, n;/m 控制迭代次数, 而 n 控制复化梯形积分的分点数. n=2mdouble h, x;double s, q;double ep; /精度要求double *y = new doubleMAXREPT;/为节省空间,只需一维数组/每次循环依次存储 Romberg 计算
9、表的每行元素,以供计算下一行,算完后更新double p;/p 总是指示待计算元素的前一个元素(同一行)/迭代初值数值计算课程设计- 17 -h = bb - aa;y0 = h*(f(aa) + f(bb)/2.0;m = 1;n = 1;ep = epsilon + 1.0;/迭代计算while (ep = epsilon) & (m ab;cout#includeusing namespace std;const int max = 50;float xmax, ymax, hmax;float cmax, amax, fxymmax;float f(int x1, int x2, in
10、t x3)float a = (yx3 - yx2) / (xx3 - xx2);float b = (yx2 - yx1) / (xx2 - xx1);return (a - b)/(xx3 - xx1); /求差分void cal_m(int n) /用追赶法求解出弯矩向量 Mfloat Bmax;B0 = c0 / 2;for(int i = 1; i = 0; i-)fxymi = fxymi - Bi*fxymi+1;三次样条插值算法- 22 -void printout(int n);int main()int n,i; char ch;do coutn;for(i = 0; i
11、xi; coutyi; for(i = 0; i t;switch(t)case 1:coutf0f1;数值计算课程设计- 23 -c0 = 1; an = 1;fxym0 = 6*(y1 - y0) / (x1 - x0) - f0) / h0;fxymn = 6*(f1 - (yn - yn-1) / (xn - xn-1) / hn-1;break;case 2:coutf0f1;c0 = an = 0;fxym0 = 2*f0; fxymn = 2*f1;break;default:coutch;三次样条插值算法- 24 -while(ch = y | ch = Y);return 0
12、;void printout(int n)cout 0) cout 0) cout 0) cout 0) cout#include#define MAX 20using namespace std;/求解任意可逆矩阵的逆，X 为待求解矩阵，E 为全零矩阵，非单位矩阵，也可以是单位矩阵void inv(double XMAXMAX,int n,double EMAXMAX)int i,j,k;double temp=0;for(i=0;in;coutXi;coutYi;coutM;for(i=0;i=0;i-)cout=0;i-)if(i=0)coutusing namespace std;#i
13、ncludedouble f(double x)return x*x-2*x-1double df(double x)return 2*x-2;int main() int k=0,max1=0;double p0,delta,epsilon,p1,err,relerr,y;cout#include#include#define MAX 20#define eps 1e-10using namespace std;double g(double x) return 1-x*x/4+x; main() double PMAX=0,err=0.0,relerr=0.0,tol=0.0,p=0.0,
14、p0=0.0;int k=0,max1=0,i=0;图 8-1 程序调试图不动点法解非线性方程- 34 -coutusing namespace std;double func(double X,int k,double x,int n);int main()double Sn=0;int n;coutn;拉格朗日插值- 36 -double*x=(double*)malloc(n*sizeof(double);double*y=(double*)malloc(n*sizeof(double);double X;int i;for(i=0;ixiyi;coutX;for(i=0;i#include using namespace std;#define delta 0.0000001设计体会及今后的改进意见- 38 -#define n 3#define max1 1000main()int i,j,k;double err,norm,Ann,Bn,Pn,Xn;coutAij;coutBi;coutPi;for (k=0;kmax1;k+)for (j=0;jn;j+)Xj=Bj;for (i=0;in;i+)Xj=Xj-Aji*Pi;Xj=
