文档介绍:该【微分方程数值解法C语言-课程设计 】是由【青山代下】上传分享,文档一共【11】页,该文档可以免费在线阅读,需要了解更多关于【微分方程数值解法C语言-课程设计 】的内容,可以使用淘豆网的站内搜索功能,选择自己适合的文档,以下文字是截取该文章内的部分文字,如需要获得完整电子版,请下载此文档到您的设备,方便您编辑和打印。:..微分方程数值解法C语言-课程设计微分方程数值解法C语言由于对matlab语言不熟悉,所以还是采用C。前面几个都比较简单,最后一个需要解非其次方程组。采用高斯—Jordan消元法(数值分析)求逆解方程组,也再一次体会到算法本身的重要性,而不是语言。当然,矩阵求逆的算法也在100个经典的C语言算法之列。不过偏微分方程数值解的内容的确比较高深,我只能停留在编这种低级的东西的自娱自乐中。不过解决计算机、数学、信计专业的课程设计还是足够了。由于篇幅所限,只把源代码粘贴在这。一:预报矫正格式#include<>#include<iostream>#include<>doublecount_0(doublexn,doubleyn){//矫正格式doubles;s=yn+*(yn/xn*+xn*xn/yn*);returns;}doublecount_1(doublexn,doubleyn,doubley0){//预报格式doubles;s=yn+*((yn/xn*+xn*xn/yn*)+(y0/xn*+xn*xn/y0*));:..returns;}voidmain(){//计算,,进行10次计算,设初始值doublexn=1,yn=1;inti=1;while(i<=10){xn=xn+;yn=count_1(xn,yn,count_0(xn,yn));i++;}}二显示差分格式#include<iostream>#include<>#include<>main(){doublea[6][11];//初始化;for(inti=0;i<=5;i++):..{a[0]=0;a[10]=0;}for(intj=1;j<10;j++){doublep=*j*;a[0][j]=sin(p);}//按显示格式计算for(i=1;i<=5;i++)for(j=1;j<10;j++)a[j]=a[i-1][j-1]+a[i-1][j+1];//输出计算好的矩阵for(i=0;i<=5;i++){for(j=0;j<11;j++)}}三龙阁库塔格式#include<>#include<iostream>#include<>:..doublecount_k(doublexn,doubleyn){doubles;s=yn/xn*+xn*xn/yn*;returns;}voidmain(){//=1,yn=1;inti=1;while(i<=11){doublek1=count_k(xn,yn);doublek2=count_k(xn+,yn+*k1);doublek3=count_k(xn+,yn+*k2);doublek4=count_k(xn+,yn+*k3);yn=yn+*(k1+2*k2+2*k3+k4);xn=xn+;i++;}}四CRANK--NICOLSON隐式格式:..#include<iostream>#include<>#include<>doubleSurplus(doubleA[],intm,intn);double*MatrixInver(doubleA[],intm,intn);double*MatrixOpp(doubleA[],intm,intn)/*矩阵求逆*/{inti,j,x,y,k;double*SP=NULL,*AB=NULL,*B=NULL,X,*C;SP=(double*)malloc(m*n*sizeof(double));AB=(double*)malloc(m*n*sizeof(double));B=(double*)malloc(m*n*sizeof(double));X=Surplus(A,m,n);X=1/X;for(i=0;i<m;i++)for(j=0;j<n;j++){for(k=0;k<m*n;k++)B[k]=A[k];{for(x=0;x<n;x++)B[i*n+x]=0;for(y=0;y<m;y++)B[m*y+j]=0;:..B[i*n+j]=1;SP[i*n+j]=Surplus(B,m,n);AB[i*n+j]=X*SP[i*n+j];}}C=MatrixInver(AB,m,n);returnC;}double*MatrixInver(doubleA[],intm,intn)/*矩阵转置*/{inti,j;double*B=NULL;B=(double*)malloc(m*n*sizeof(double));for(i=0;i<n;i++)for(j=0;j<m;j++)B[i*m+j]=A[j*n+i];returnB;}doubleSurplus(doubleA[],intm,intn)/*求矩阵行列式*/{inti,j,k,p,r;doubleX,temp=1,temp1=1,s=0,s1=0;:..if(n==2){for(i=0;i<m;i++)for(j=0;j<n;j++)if((i+j)%2)temp1*=A[i*n+j];elsetemp*=A[i*n+j];X=temp-temp1;}else{for(k=0;k<n;k++){for(i=0,j=k;i<m,j<n;i++,j++)temp*=A[i*n+j];if(m-i){for(p=m-i,r=m-1;p>0;p--,r--)temp*=A[r*n+p-1];}s+=temp;temp=1;}for(k=n-1;k>=0;k--){for(i=0,j=k;i<m,j>=0;i++,j--)temp1*=A[i*n+j];if(m-i){for(p=m-1,r=i;r<m;p--,r++)temp1*=A[r*n+p];}:..s1+=temp1;temp1=1;}X=s-s1;}returnX;}voidinitmat_A(doublea[][9],doubler){for(inti=0;i<9;i++)for(intj=0;j<9;j++)a[j]=0;for(i=0;i<9;i++){a=1+r;if(i!=8)a[i+1]=-*r;if(i!=0)a[i-1]=-*r;}}voidinitmat_B(doubleb[][9],doubler){for(inti=0;i<9;i++)for(intj=0;j<9;j++)b[j]=0;for(i=0;i<9;i++){b=1-r;:..if(i!=8)b[i+1]=*r;if(i!=0)b[i-1]=*r;}}voidinitmat_C(doubleC[][9]){for(inti=0;i<9;i++)for(intj=0;j<9;j++)C[j]=0;}voidmain(){doublea[100][11];for(inti=0;i<100;i++)for(intj=0;j<11;j++)a[j]=0;//初始化;for(i=0;i<100;i++){a[0]=0;a[10]=0;}for(intj=1;j<10;j++){doublep=4**j*;a[0][j]=sin(p);}:..//取h=*,r=,t=**;//得到矩阵a和矩阵bdoubleA[9][9];initmat_A(A,);doubleB[9][9];initmat_B(B,);//B矩阵与Un相乘,en是0;doubleC[9][9];initmat_C(C);double*A_;A_=MatrixOpp(A[0],9,9);//A矩阵求逆;//A逆*Bfor(i=0;i<9;i++)for(j=0;j<9;j++)for(ints=0;s<9;s++)C[j]+=A_[i*9+s]*B[s][j];//填写a表格for(i=0;i<100;i++){for(j=1;j<10;j++)for(ints=0;s<9;s++)a[i+1][j]+=a[s+1]*C[j-1][s];}//输出表格for(i=0;i<100;i++):..{for(j=0;j<11;j++)}//利用精确解,求出表格for(i=0;i<100;i++){for(j=0;j<11;j++)}}