FORTRAN课程设计(实践)报告.doc

FORTRAN课程设计(实践)报告一 任务内容线性方程组求解问题。一个物理系统可用下列线性方程组来表示：从文件中读入m1、m2和的值，求a1、a2、N1 和N2的值。其中g取9.8，输入时以角度为单位。要求：（1）分别用两种方法（例如高斯消去法、矩阵求逆法、三角分解法、追赶法等），定义求解线性方程组Ax=b的子程序，要求该子程序能求解任意线性方程组。（2）在主程序中分别调用上面定义的两个子程序，并对求解结果进行对比分析。（3）绘制以上两个方法所求得的方程解的数据分布图。二 系统开发设计思想对于任务中的线性方程方程组求解问题，我采用的是高斯消去法和矩阵求逆法两法并行。由于考虑到任务要求子程序能够求解任意线性方程组，因此在设计程序的过程中我采用子例行程序作线性方程组的求解工作，在主程序中针对任务中的物理系统进行专题设计以解决该物理问题，因此该程序具有较为广泛的通用性。具体系统结构图如下：子程序1高斯消去法子程序2矩阵求逆法主程序分别调用两个子程序用于解决该具体物理问题输出结果并对比分析设计流程大致如下：（1） 仔细阅读任务要求，由当前所学知识确定采用高斯消去法和矩阵求逆法求解线性方程组；（2） 查阅书籍，掌握高斯消去法以及矩阵求逆法，并将其转换为计算机编程思维；（3） 由简及复，将思维转化为编程语言，逐步完成子程序设计语言，通过数据结果的反复试验调试程序，完成修改完善；（4） 针对任务中的物理问题设计主程序并修改完善；（5） 在主程序中分别调用上面定义的两个子程序，并对求解结果进行对比分析。三 系统功能及系统详细设计功能与设计分析：该系统采用两种不同的方法（高斯消去法和矩阵求逆法）为使用者提供两种解题过程，最优越之处在于子程序的设计是面向任意的线性方程组求解问题，使得该系统能够通过调整改变主程序的具体设置而解决各类具体相关问题，大大突破了该系统的局限性问题。本系统的主程序是为解决该具体物理问题而设计，采用两种方法能有效地对结果进行对比分析，有利于更精确更深入的探究。对于物理问题中的多种条件而导致的需要对不同情况下的输出结果的情况，我采用以txt文件的形式读入数据，从而取代了繁复的手工键盘输入的操作，提高了系统解决问题的效率，也减小了手工输入而带来错误的几率。功能结构：高斯消去法作为一个子程序，矩阵求逆法作为一个子程序，均用子例行程序编写好，在主程序中根据两种不同方法的思路设置该物理问题条件下的具体数值，调用两个子程序输出结果。四 遇到的问题与解决1.高斯消去法简介首先写出线性方程组的增广矩阵（A,B），然后利用初等行变换将A化为行阶梯形矩阵，再进一步化为行最简形求解，最后将结果一步步代回即可。简要示例如下： 2.对于方程组无解或多解情况的处理设计出的程序只能够求解有唯一解的线性方程组，无法适用于无解或多解的情况。针对此问题的存在，考虑到无解与多解的情况难以划分，我只能将无解和多解的情况归为一类，在此情况下运行程序，程序会提示无法输出结果。3.创建文件读写数据问题采用直接键盘输入数据的方式既麻烦又容易出错，所以在程序的细节设计上，我打算在程序中创建一个txt文档，然后通过程序代码输入数据，再调用这个文本中的数据进行计算。这是一直困扰我的问题之一，开始的时候我通过查阅资料，自己设计了直接调用一个事先输入好数据的txt文件的办法进行计算，但对于先在程序中创建文件再进行调用的程序设计还存在一些疑惑，因此我主动咨询了老师，和同学一起讨论后茅塞顿开，很好地完成好了这一设计，使我的程序设计更加完善，又前进了一大步。4.主程序中针对具体的物理情况下数值设置问题要将能够求解任意线性方程组求解的子例行程序应用到任务要求的问题上，我想采用简便有效的方法，但最终在实践中无法实现，而只好采用基本的定义数组，再在数组中对具体线性方程组中的元素逐个赋值的方式。五 总结用于解决具体工程问题的FORTRAN实践没有事先想像得那么简单，在整个编程解决问题的过程中，首先任务中的问题是第一道坎。研究问题，首要的就是把握问题，找到合适的方法，再将其转换为计算机编程的算法，最后将其用具体的程序语言表达出来。在研究问题中，我体会到了在平时课程中无法体验的探索经历。在编程中，我遇到了很多困难，大都是如何将习以为常的数学思维转变成计算机语言。计算机作为当下科学的最重要的辅助工具，熟悉学习计算机，并利用计算机去解决问题，是我们当代大学生务必要掌握好的本领，更是我们在以后的学术道路上的重要凭借与基础，计算机的学习道路漫漫，需要我们持之以恒的求学态度和不断探索的钻研精神。总体而言，本次FORTRAN实践我能良好地完成好任务中的要求，在问题中探索出路，在学习中发现经验，并且在老师的指导下，在网络资源的合理利用下，我学会了自我探索自我发现，而不再是局限于平时课程中的单一接受现成知识的模式。这是一次对学习方法的革新，是一次对学习过程的全新认识。六 参考文献1FORTRAN 90 程序设计教程 主编：刘卫国，蔡旭辉 北京邮电大学出版社2FORTRAN 90 程序设计上机指导与习题选解 主编：刘卫国，戴忠 北京邮电大学出版社七 FORTRAN算法subroutine czs(a,b,n,x,l,js) !高斯消去法 dimension a(n,n),x(n),b(n),js(n) double precision a,b,x,t l=1 !逻辑变量 do k=1,n-1 d=0.0 do i=k,n do j=k,n if (abs(a(i,j)d) then d=abs(a(i,j) js(k)=j is=i end if end do end do !把行绝对值最大的元素换到主元位置 if (d+1.0=1.0) then l=0 else !最大元素为0无解 if(js(k)/=k) then do i=1,n t=a(i,k) a(i,k)=a(i,js(k) a(i,js(k)=t end do end if if(is/=k) then do j=k,n t=a(k,j) a(k,j)=a(is,j) a(is,j)=t end do t=b(k) b(k)=b(is) b(is)=t end if !最大元素在主对角线上 end if !消去 if (l=0) then write(*,100) return end if do j=k+1,n a(k,j)=a(k,j)/a(k,k) end do b(k)=b(k)/a(k,k) !求三角矩阵 do i=k+1,n do j=k+1,n a(i,j)=a(i,j)-a(i,k)*a(k,j) end do b(i)=b(i)-a(i,k)*b(k) end do end do if (abs(a(n,n)+1.0=1.0) then l=0 write(*,100) return end if x(n)=b(n)/a(n,n) do i=n-1,1,-1 t=0.0 do j=i+1,n t=t+a(i,j)*x(j) end do x(i)=b(i)-t end do100 format(1x,fail) js(n)=n do k=n,1,-1 if (js(k)/=k) then t=x(k) x(k)=x(js(k) x(js(k)=t end if end do return end program maindimension a(4,4),b(4),x(4),js(4)double precision a,b,xreal m1,m2,jopen(1,file=f:a.txt)write(1,*)123,21,12close(1)open(1,file=f:a.txt)read(1,*)m1,m2,jclose(1) n=4print*,m1,m2,ja(1,1)=m1*cos(3.14159*j/180)a(1,2)=-m1a(1,3)=-sin(3.14159*j/180)a(1,4)=0a(2,1)=m1*sin(3.14159*j/180)a(2,2)=0a(2,3)=cos(3.14159*j/180)a(2,4)=0a(3,1)=0a(3,2)=m2a(3,3)=-sin(3.14159*j/180)a(3,4)=0a(4,1)=0a(4,2)=0a(4,3)=-cos(3.14159*j/180)a(4,4)=1b(1)=0 b(2)=m1*9.8b(3)=0b(4)=m2*9.8 call czs(a,b,n,x,l,js) if (l/=0) then write(*,10)x(1),x(2),x(3),x(4)end if10 format(1x,1f15.6)endsubroutine zyt(a,n,l,is,js) !逆矩阵求解法dimension a(n,n),is(n),js(n)double precision a,t,dl=1do k=1,n d=0.0 do i=k,n do j=k,n if(abs(a(i,j)d) then d=abs(a(i,j) is(k)=i js(k)=j end if end do end do if (d+1.0=1.0)then l=0 write(*,200) return end if 200 format(1x,err*not inv) do j=1,n t=a(k,j) a(k,j)=a(is(k),j) a(is(k),j)=t end do do i=1,n t=a(i,k) a(i,k)=a(i,js(k) a(i,js(k)=t end do a(k,k)=1/a(k,k) do j=1,n if(j.ne.k)then a(k,j)=a(k,j)*a(k,k) end if end do do i=1,n if(i/=k)then do j=1,n if(j/=k)then a(i,j)=a(i,j)-a(i,k)*a(k,j) end if end do end if end do do i=1,n if(i/=k)then a(i,k)=-a(i,k)*a(k,k) end if end doend dodo k=n,1,-1 do j=1,n t=a(k,j) a(k,j)=a(js(k),j) a(js(k),j)=t end do do i=1,n t=a(i,k) a(i,k)=a(i,is(k) a(i,is(k)=t end do end do return endsubroutine yt(a,b,n,c)dimension a(n,n),b(n),c(n)double precision a,b,cdo i=1,ndo j=1,nc(i)=0.0do l=1,nc(i)=c(i)+a(i,l)*b(l)end do end doend doreturnendprogram maindimension a(4,4),b(4,1),c(4,1),is(4),js(4)double precision a,b,creal m1,m2,jdopen(1,file=f:a.txt)write(1,*)123,21,12close(1)open(1,file=f:a.txt)read(1,*) m1,m2,jdprint*,m1,m2,jdclose(1)a(1,1)=m1*cos(3.14*jd/180)a(1,2)=-m1a(1,3)=-sin(3.14*jd/180)a(1,4)=0a(2,1)=m1*sin(3.14*jd/180)a(2,2)=0a(2,3