弦振动方程数值求解02.docVIP

数学物理方程 之 基于数值计算方法的弦振动方程求解  目录 摘要、关键词…………………………………………… 2页 有限差分法介绍………………………………………… 3页 程序描述………………………………………………… 6页 计算机处理……………………………………………… 8页 Matlab作图…………………………………………… 10页 特别鸣谢………………………………………………… 11页  摘要、关键词 摘要：继上次关于弦振动方程的“平行四边形法则”求解之后，我们又从数值计算的角度入手，对弦振动方程进行计算和模拟，从而验证“平行四边形法则”解弦振动方程的正确性。 关键词：有限差分法、数值计算、弦振动方程 附：(弦振动方程：  有限差分法介绍 以弦振动方程为例： 对于一定的u（x，t），我们用“差分”代替“微商”，从而将 数差值描述，可得： 以及 将第一个式子的右边第一项移至左边,得： ^… 同理可得， 两式做差： 用中心差分的一阶导数表示二阶导数，化简： 由此引入 则 则弦振动方程 可以表示为： 我们定义 为网格比 则 由此可知，每一个格点u（i，n)均由u（i+1，n+1）、 u（i，n）、u(i-1,n-1)、 u(i,n-2)等其余四点所确定： 由此我们可以采用“递归”的思想，借助计算机进行快速计算，从而得到各个格点的值.值得注意的是,(在边界上u≡0.(在初始层上的点(即u（i，0））无法用上述公式计算，还需借助初始条件，即： 由 和 两式相加，消去可得 综上： 届此，我们可以将此式编入程序（采用“递归”思想），详细代码见下一节。 程序描述 C语言： #include stdio.h #include stdlib.h #define imax 10 //定义长、宽格点数，其数值参考计算机性能而定 #define jmax 60 double f(double i,double n); double u(int i,int n); double FEI(double x); double PUSI(double x); double l=1,t=3; //设l=1，t=3 double dert_t,dert_x; //长宽格距 double R; int main () { int i,n,a; double A[imax][jmax]; FILE *fp; dert_x=l/imax; dert_t=t/jmax; a=2; //a的值 R=a*a*dert_t*dert_t/dert_x/dert_x; printf(%lf\n,R); if((fp=fopen(G:\\matlab\\myfile.txt,r+))==NULL) //创建文件 { printf(cannnot open file!\n); return 0; } for(n=0;njmax;n++){ for (i=0;iimax;i++){ A[i][n]=u(i,n); printf(%9.5lf ,A[i][n]); fprintf(fp,%9.5lf,A[i][n]); fputs( ,fp); } printf(\n); fputs(\n,fp); } fclose(fp); return 0; } double u(int i,int n) //函数u(x,t) { double result; if(i==0||i==imax-1) { result=0; } else { if(n==0) { result=FEI(i*dert_x); } else if(n==1) { result=0.5*(R*R*u(i+1,0)+2*(1-R*R)*u(i,0)+R*R*u(i-1,0)+2*dert_t*PUSI( x)+dert_t*dert_t*f(i*dert_x,n*dert_t)); //res