fortran90 解線性方程組程序

1、揚(yáng)州大學(xué) 解線性方程組 計算流體動力學(xué)題目：用fortran語言編一程序解下列n元一次方程組，并對相同的線性方程組改變其矢量項(xiàng)的一個參數(shù)，對比其計算結(jié)果的變化。解題方法：采用解線性方程組直接算法中的高斯列主元消去法。高斯列主元消去法，較高斯消去法舍入誤差小，較高斯全主元消去法計算量少，所以這里采用高斯列主元消去法。算法基本思想：該方法分消元和回代兩個過程，對其增廣矩陣進(jìn)行轉(zhuǎn)換。1.消元過程：對k=1，2，,(n-1)，依次計算其中，對方程，其增廣矩陣，在每一步消元前先進(jìn)行按列選主元的工作，如第k步消元前，就所有，取絕對值最大值，設(shè)，將第l個方程與第 i個方程進(jìn)行互換位置，這樣成為第k步的主元素

2、，然后進(jìn)行第k步消元。每步消元都如此，最終將其化為上三角矩陣。2.回代過程：依次分別計算出xn，xn-1，xn-2，x2，x1。源程序如下：program linear_equation_bc !主程序開始語句implicit none !說明部分integer:i,j,k,imax,t !定義整型變量，i、j分別標(biāo)記增廣矩陣的行列!imax標(biāo)記列主元最大的行號，t為輸入數(shù)組元數(shù)real:max,n !max記錄每次計算的列主元的值，n方程組的元數(shù)!用矩陣(實(shí)型動態(tài)數(shù)組)將線性方程組表示出來real,dimension(:,:),allocatable:a,m !定義實(shí)型二維數(shù)組a、m，分別存

3、放初始!增廣矩陣和選主元時交換對應(yīng)行時的臨時矩陣real,dimension(:),allocatable:x !定義實(shí)型一維數(shù)組x，存放計算結(jié)果print*,請輸入線性方程組的元數(shù)n： !表控輸出語句read*,t !輸入語句，輸入數(shù)組的元數(shù)!給動態(tài)數(shù)組分配內(nèi)存allocate(a(t,t+1),m(t,t+1),x(t) !a、m為t行t+1列二維數(shù)組，x為t行一維數(shù)組print*,請依次按行輸入線性方程組的增廣矩陣a： !表控輸出語句read*,(a(i,j),j=1,t+1),i=1,t) !隱do循環(huán)輸入語句，輸入增廣矩陣ado k =1,t-1 !guass列主元消去法消元過程開始

4、 !下面選取每列的最大列主元素max=abs(a(k,k)imax=kdo i=k+1,t+1 !do語句if (abs(a(i,k)max) then !do塊 !if-then語句max=abs(a(i,k) !if-then塊imax=iend if !end if語句end do !end do語句!將最大列元素所在的行與第k行進(jìn)行交換do j=k,t+1m(k,j)=a(k,j)a(k,j)=a(imax,j)a(imax,j)=m(k,j)end do!對方程組按x1, x2, xt的順序進(jìn)行依次消元!將增廣矩陣a化為上三角矩陣do i=k+1,tm(i,k)=a(i,k)/a(k

5、,k)do j=k+1,t+1a(i,j)=a(i,j)-m(i,k)*a(k,j)end doend doend do !guass列主元消去法的第一步消去至此結(jié)束x(t)= a(t,t+1)/a(t,t) !先計算出xtdo k=t-1,1,-1 !回代計算開始n=0 !給定n初值，n作為累加變量，記錄每次計算中的累加值do j=t,k+1,-1 !do循環(huán)，計算n的累加值n=n+a(k,j)*x(j)end dox(k)=( a(k,t+1)-n)/a(k,k) !根據(jù)公式算的x(k)的值，即為計算結(jié)果end do !回代過程結(jié)束print*,線性方程組的解為： !表控輸出語句do i =1,t !用顯do循環(huán)輸出計算結(jié)果print*,x(i) !表控輸出語句end do end program linear_equation_bc !主程序結(jié)束語句以上程序如圖1所示。圖1 fortran解線性方程組程序算例1：方程組如下：方程組的解為： ；手算解為：計算結(jié)果如圖2所示。圖2 算例1計算結(jié)果算例2：方程組如下：方程組的解為：；手算解為：在誤差要求范圍之內(nèi)，滿足要求。計算結(jié)