Fortran学习19：数值计算方法2

老感觉昨天自己最后写的那个python牛顿法求解，老觉得不太对，，其实后面我想学完Fortran之后，再跟着matlab书把数据处理再学一下，因为Fortran那本书讲的东西实在是太少了。，而且听师兄说，我们不用Fortran绘图，后面还是matlab和python来绘图。

今天先不说了，直接开始矩阵运算，其实矩阵的加减法没啥好说的，就是相同位置进行加减运算就行了，主要是乘法，Fortran库里有一个matmul()函数来计算，但是今天还是先自己敲吧。

program main
    implicit none
    integer::a(3,3)
    integer::b(3,3)
    integer::c(3,3)
    integer i,j,k
    data a / 1,2,3,4,5,6,7,8,9 /
    data b /10,11,12,13,14,15,16,17,18/
    ! 加法
    ！可以直接使用c=a+b，下面是为了展示运算过程
    do i=1,3
        do j=1,3
            c(i,j)=a(i,j)+b(i,j)
        end do
    end do
    ! 减法
    ！可以直接使用c=a-b，下面是为了展示运算过程
    do i=1,3
        do j=1,3
            c(i,j)=a(i,j)-b(i,j)
        end do
    end do
    ! 乘法
    do i=1,3
        do j=1,3
            c(i,j)=0
            do k=1,3
                c(i,j)=c(i,j)+a(i,k)*b(k,j)
            end do 
        end do
    end do
    ! 再试试Fortran自带的函数matmul
    c=MATMUL(a,b)
    write(*,"(I5,I5,I5/I5,I5,I5/I5,I5,I5)")c
    stop "over"
end program 

放心，没问题的，ctrl+c，ctrl+v之后完全可以运行。

其实也可以在主程序外面把这几个运算方法写成function或者subroutine，我觉得都还行，就是个传参嘛。反正我是没有用书上的代码，习惯了，不喜欢整个程序纯粹用大写字母。

1、三角矩阵：上三角和下三角，上三角对角线以下的数全部是0，相应的，下三角对角线以上的全部是零：

跟着书上的代码敲一遍吧，我线代的储备主要还是考研时候听网课攒的点底子(大二学线代的时候才考了85+，有点低)

module LinearAlgebra
    implicit none
    contains
        subroutine output(matrix)
            implicit none
            integer::m,n
            real::matrix(:,:)
            integer::i
            CHARACTER(len=20)::for='(??(1x,f6.3))'
            m=size(matrix, dim=1)
            n=size(matrix, dim=2)
            WRITE(for(2:3),'(I2)')N
            do i=1,N
                write(*,Fmt=for)matrix(i,:)
            end do
            return 
        end subroutine 
        ! 求解上三角矩阵：
        subroutine upper(matrix)
            implicit none
            integer::m,n
            real::matrix(:,:)
            integer::i,j
            real::e
            m=size(matrix, dim=1)
            n=size(matrix, dim=2)
            do I=1,n-1
                do j=i+1,m
                    e=matrix(j,i)/matrix(i,i)
                    matrix(j,i:m)=matrix(j,i:m)-matrix(i,i:m)*e
                end do
            end do
            return 
        end subroutine
        ! 求解下三角矩阵：
        subroutine lower(matrix)
            implicit none
            integer::m,n
            real::matrix(:,:)
            integer::i,j
            real::e
            m=size(matrix, dim=1)
            n=size(matrix, dim=2)
            do I=n,2,-1
                do j=i-1,1,-1
                    e=matrix(j,i)/matrix(i,i)
                    matrix(j,1:i)=matrix(j,1: