一维对流方程在三种差分格式下的解_精品文档资料下载.pdf
- 文档编号:16123325
- 上传时间:2022-11-20
- 格式:PDF
- 页数:12
- 大小:604.72KB
一维对流方程在三种差分格式下的解_精品文档资料下载.pdf
《一维对流方程在三种差分格式下的解_精品文档资料下载.pdf》由会员分享,可在线阅读,更多相关《一维对流方程在三种差分格式下的解_精品文档资料下载.pdf(12页珍藏版)》请在冰豆网上搜索。
2014年03月31号温度:
18C学号:
11309018实验二实验二:
对流方程差分格式下的解对流方程差分格式下的解一、一、实验目的实验目的1、用数值方法计算一维对流方程在A、B、C三种差分格式下的解。
2、取为0.05.取值为0.5,1,2。
并作相关讨论。
二:
实验仪器与设备二:
实验仪器与设备:
装有Fortran的计算机1台三三、实验原理、实验原理1x1x00)(x,1x010)(x,0x1-x10)(x,0t8x8-,0或xxt学会对MS-FORTRAN的基本操作。
用Fortran编写程序计算一维对流方程在A、B、C三种格式下的解。
讨论各种格式下不同的tx值的差分格式解的特点。
三:
实验步骤三:
实验步骤:
首先A格式,我们对微分方程进行离散化,得出A格式的差分解的表达式,为1110()2()nnnniiiiiitxxPage2of12中山大学工学院、理论与应用力学刘广编制实验编号及题目:
工学院姓名:
刘广学号:
11309018日期:
2014年03月31号这里初始时刻的方程为(x,0)1x-1x0(x,0)10x1(x,0)0x1x1x或,即为三角波,其中传播速度为1,我们可以写如下Fortran90程序:
!
主程序,根据输入选择差分格式,输入1为A格式,2为B格式,3为C格式!
其中输入DX_DT为差分比值,根据输入的比值划分时间网格,其中空间网格为0.05,区间为-8到8!
在A格式中,根据依赖区域,不能计算的左右上三角强行赋值为0!
在B格式中,根据依赖区域右上三角赋值为0!
在C格式中,根据依赖区域左上三角赋值为0!
最后输出为csv格式表格,其中第一列为时间,第二列为空间,第三列为波形数值!
画图时,需要选定同一时间点的二三两列,即可画图PROGRAMMAINIMPLICITNONE!
声明差分比值REAL:
DX_DT!
声明差分类型INTEGER:
TYPE_ABCWRITE(*,*)请输入Dx/Dt,0.5,1或2READ(*,*)DX_DTWRITE(*,*)请选择差分格式,A格式输入1,B格式输入2,C格式输入3READ(*,*)TYPE_ABC!
判断类型调用不同的子程序IF(TYPE_ABC=1)THENCALLFORMAT_A(DX_DT)ELSEIF(TYPE_ABC=2)THENCALLFORMAT_B(DX_DT)ELSEIF(TYPE_ABC=3)THENCALLFORMAT_C(DX_DT)ENDIFWRITE(*,*)请在目录下查看结果PAUSESTOPENDPage3of12中山大学工学院、理论与应用力学刘广编制实验编号及题目:
2014年03月31号SUBROUTINEFORMAT_A(DX_DT)IMPLICITNONE!
I,J,K为循环变量,REC_T为时间网格总数INTEGER:
I,J,K,REC_T!
输入dx/dt的比值REAL:
声明动态矩阵,用于存储计算结果REAL,ALLOCATABLE:
M(:
:
)!
打开通道号8,输出OUTPUT_C.CSVOPEN(UNIT=8,FILE=OUTPUT_A.CSV)!
计算时间T的网格,节点数为网格数加1REC_T=160/DX_DT!
确定M大小,开始在网格第一层赋初值,其中循环不能从负数开始,X网格右移ALLOCATE(M(REC_T+1,321)DOI=1,140,1M(1,I)=0WRITE(8,*)1,I,M(1,I)ENDDODOI=141,161,1M(1,I)=1+0.05*(I-161)WRITE(8,*)1,I,M(1,I)ENDDODOI=162,181,1M(1,I)=1-0.05*(I-161)WRITE(8,*)1,I,M(1,I)ENDDODOI=182,321,1M(1,I)=0WRITE(8,*)1,I,M(1,I)ENDDO!
时间上从第二层开始循环,逐步积分,注意,时间层循环上限为Rec_TDOI=2,REC_T,1!
计算三角网格部分数值,循环部分从I到322-IDOJ=I,322-I,1M(I,J)=M(I-1,J)-(M(I-1,J+1)-M(I-1,J-1)/(2*DX_DT)WRITE(8,*)I,J,M(I,J)Page4of12中山大学工学院、理论与应用力学刘广编制实验编号及题目:
2014年03月31号ENDDO!
左上三角强行赋值为0DOK=1,I-1,1M(I,K)=0WRITE(8,*)I,K,M(I,K)ENDDO!
右上三角强行赋值为0DOK=322-I,321,1M(I,K)=0WRITE(8,*)I,K,M(I,K)ENDDOENDDORETURNENDSUBROUTINESUBROUTINEFORMAT_B(DX_DT)!
打开通道号8,输出OUTPUT_C.CSVOPEN(UNIT=8,FILE=OUTPUT_B.CSV)!
确定M大小,开始在网格第一层赋初值,其中循环不能从负数开始,X网格右移ALLOCATE(M(REC_T+1,321)DOI=1,140,1M(1,I)=0WRITE(8,*)1,I,M(1,I)ENDDODOI=141,161,1M(1,I)=1+0.05*(I-161)WRITE(8,*)1,I,M(1,I)ENDDOPage5of12中山大学工学院、理论与应用力学刘广编制实验编号及题目:
2014年03月31号DOI=162,181,1M(1,I)=1-0.05*(I-161)WRITE(8,*)1,I,M(1,I)ENDDODOI=182,321,1M(1,I)=0WRITE(8,*)1,I,M(1,I)ENDDO!
时间上从第二层开始循环,逐步积分DOI=2,REC_T+1,1!
计算三角网格部分数值,循环部分从1到322-IDOJ=1,322-I,1M(I,J)=M(I-1,J)-(M(I-1,J+1)-M(I-1,J)/DX_DTWRITE(8,*)I,J,M(I,J)ENDDO!
右上三角强行赋值为0DOK=322-I,321,1M(I,K)=0WRITE(8,*)I,K,M(I,K)ENDDOENDDORETURNENDSUBROUTINESUBROUTINEFORMAT_C(DX_DT)!
打开通道号8,输出OUTPUT_C.CSVOPEN(UNIT=8,FILE=OUTPUT_C.CSV)!
计算时间T的网格,节点数为网格数加1Page6of12中山大学工学院、理论与应用力学刘广编制实验编号及题目:
2014年03月31号REC_T=160/DX_DT!
计算三角网格部分数值,循环部分从I到321DOJ=I,321,1M(I,J)=M(I-1,J)-(M(I-1,J)-M(I-1,J-1)/DX_DTWRITE(8,*)I,J,M(I,J)ENDDO!
左上三角强行赋值为0DOK=1,I-1,1M(I,K)=0WRITE(8,*)I,K,M(I,K)ENDDOENDDORETURNENDSUBROUTINEPage7of12中山大学工学院、理论与应用力学刘广编制实验编号及题目:
2014年03月31号其中对于程序说明在程序注释中已经说明的很清楚,这里不再赘述。
四:
实验结果四:
实验结果最终程序执行完毕之后会在目录下生成对应的输出文档,根据用户输入的=1xt数值的大小生成对应的网格,然后计算出所有网格节点的数值。
五:
实验结果分析五:
实验结果分析这里我们选择=1xt为例进行A,B,C三种格式的计算结果说明。
A格式,根据所的结果,我们绘制不同时刻的波形图像,选取时刻为初始时刻,第35个时间区间时刻(3.5s),和最后的80区间时刻(8s),图像如下:
Page8of120.00E+002.00E-014.00E-016.00E-018.00E-011.00E+001.20E+00115294357718599113127141155169183197211225239253267281295309初始时刻-6.00E-01-4.00E-01-2.00E-010.00E+002.00E-014.00E-016.00E-018.00E-011.00E+001.20E+001.40E+00050100150200250300350时间区间35中山大学工学院、理论与应用力学刘广编制实验编号及题目:
2014年03月31号很明显,我们可以看到,随着时间的推移,A格式明显呈现出发散状态。
我们再选择B格式进行说明。
B格式,根据所的结果,我们绘制不同时刻的波形图像,选取
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 对流 方程 种差 格式 精品 文档