Fortran 95


Fortran 95

Fortran程序设计基础

字符集

  • Fortran标准中规定,程序代码不区分大小写,一致认作大写。

书面格式

  • 建议使用Free Format编写程序

Fixed Format(固定格式)

不建议使用

  • 第1个字符:如果是字母Cc或星号*,这一行就会被认为是说明批注,不会被编译。
  • 第1~5个字符:如果是数字,就是给这一行程序代码取个代号,不然只能是空格。
  • 第6个字符:如果是0以外的任何字符,表示这一行会连接上一行。
  • 第7~72个字符:程序代码编写区域。
  • 第73个字符以后:不使用,超过的部分会被忽略,有的编译器报错。
C     Fixed Format demo
      program main
      write(*,*) 'Hello World!'
      write(*,*)
      1'Hello'
100   write(*,*) 'World'
10    stop
11    end

注意:程序代码命令间多数空格没有意义。

Free Format(自由格式)

  • 叹号!后面为注释
  • 每行可以编写132个字符
  • 行号放在每行程序的最前面
  • 行尾以&结束,表示下一行与此相连,行首为&,同理与上一行相连
! Free Format
program main
write(*,*) "Hello World" !注释
write(*,*) &
"Hello"
wri& !把write分开不好,但自由格式支持
&te(*,*) "World"
end
  • 空格同样无意义,方便阅读

Fortran的数据类型

  • 整数(integer) :长整型32bits-4bytes和短整型16bits-2bytes

image-20220307181529099

  • 浮点数(real):单精度32bits-4bytes和双精度64bits-8bytes

image-20220307181557734

  • 复数(complex):形式a+bi,单双精度复数
  • 字符(character):多个字符为字符串
  • 逻辑判断(logical):1——True,21——False

image-20220307181454010

image-20220307181812187

Fortran的数学表达式

+ -、* /、**、( )从左到右优先级越高

输入输出及声明

输入(WRITE)输出(PRINT)命令

基本语法

program 程序名称 !程序名称与文件名称无关
write(*,*) "Hello World"
print *,"Hello World" !与write同效
stop !可以省略
end !主程序代码结束

声明

  • 声明:提前让编译器预留一些存放数据的内存空间
    program main
    integer a !声明一个整数a
    a=3
    write(*,*) "a=",a
    stop
    end program !或者end或者end program main
    变量命名规则:
  • 英文字母、下划线、数字,以英文字母为前缀
  • fortran 77 变量长度(1~6),fortran 90变量长度(1~31)
  • 不与执行命令同名

整数类型(INTEGERAL)

program main
integer i
i=2+2*4-3
write(*,*) "i=",i
stop
end

长整型integer(kind=4) a
短整型integer(kind=2) b
Fortran 90声明变量integer :: a

浮点数(REAL)

! 浮点数声明
program main
real(kind=4) a
real :: b !默认4bytes
a = 2.5+3.0/2.0
b = a/2.0
stop 
end

数学库函数sin(x)的使用,数学库中三角函数使用的是弧度,而不是角度

program main
real :: a
a = 3.14159/2.0
write(*,*) "sin(",a,")=",sin(a)
stop
end

复数(COMPLEX)

  • 复数声明

complex a

complex(kind=4) a !   单精度
complex(kind=8) b !   双精度
c=(x,y)  !x为实部,y为虚部,c=x+yi

实例

program main
complex :: a,b
a=(1.0,1.0)
b=(1.0,2.0)
write(*,*) "a+b=",a+b
write(*,*) "a-b=",a-b
write(*,*) "a*b=",a*b
write(*,*) "a/b=",a/b
stop
end program

字符及字符串(CHARACTER)

  • 声明字符
    character c
  • 声明字符串
    character(len=10) s
  • 实例
    program main
    character c
    character(len=10) s
    c='H'
    s="llo"
    write(*,*) c,s
    stop
    end program
  • 修改字符串
    program main
      character(len=20) string
      string="Good Morning"
      write(*,*) string
      string(6)="evening" !重新设置从第六个字符开始
      write(*,*) string
      stop
    end program
    与字符串有关的函数
  • Char(num)返回计算机所使用的字符表上,数值num所代表的字符
  • Ichar(char) 返回输入的char在字符表中代表的编号,返回值为整数
  • Len(string)返回输入字符串声明的长度,返回值为整数类型
  • Len_trim(string) 返回字符串去除尾端空格后的实际内容长度
  • index(string,key)返回key在string中第一次出现的位置
  • trim(string)返回去除string尾端空格的字符串
    实例
    program main
    character(len=20) string
    character(len=5) substring
    string = "Have a nice day."
    substring = "nice"
    write(*,*) ichar('A')
    write(*,*) char(65)
    write(*,*) len(string)
    write(*,*) len_trim(string)
    write(*,*) index(string,substring)
    stop
    end program main

image-20220307152519604

逻辑变量(LOGICAL)

声明
logical a
实例

program main
logical(kind=4) a
logical(kind=2) b
a=.true.
b=.false.
write(*,*) a,b
end program

image-20220307153047946

输入命令(READ)

program main
integer a
read(*,*) a  !键盘读入一个整数,第一个星号代表默认输入设备(键盘),第二个代表不指定输入格式
write(*,*) a
integer b,c,d
read(*,*) b,c,d  !一次输入多个整数
write(*,*) b,c,d
stop
end program

image-20220307155002485

格式化输入输出(FORMAT)

格式化输出概论

format格式化输出,达到美观的目的。

program main
integer a
a = 100 !等效于integer :: a=100
write(*,100) a !100代表format(I4)
100 format(I4)
stop
end program main

image-20220307160003331

image-20220307193736172

image-20220307193805334

实例

program main
    integer :: a=10
    character(len=20) :: s="123456789"
    real :: b=12.34
    complex :: c=(1,2)
    logical :: d=.true.
    character(len=20) :: e="Fortran"
    write(*,"(A20)") s
    write(*,"(1X,I5)") a !1X表示向右移动1个位置,I5给5个位置的整型位
    write(*,"(1X,F5.2)") b !F5.2表示共5个宽度,2个小数
    write(*,"(1X,F4.1,F4.1)") c !complex用两个浮点数来输出
    write(*,"(1X,L3)") d !用3个宽度来输出logical
    write(*,"(1X,A10)") e !用10个宽度来输出字符串
end program main

image-20220307200037457

详论格式化输出

Iw[.m]:I整数,w个宽度,至少输出m位

Fw.d:F浮点数,w个字符文本框宽,d个字符小数部分

Ew.d[Ee]:用科学计数法,w个字符宽输出浮点数,d个字符小数部分,指数部分最少输出e个数字

write(*,"(E15.7)") 123.45 !输出`  0.1234500E+03`
write(*,"(E9.2E3)") 12.34 !输出` 0.12E002`

Dw.d:同Ew.d相似,只是将E变成了D

AW:以w个字符宽来输出字符串

nX:输出位置向右移动n位

Lw:以w个字符宽来输出T或F的真假值

Gw.d[Ee]:

/:换行输出

Tc:把输出的位置移动到本行的第c个字节

TLn:输出位置向左相对移动n个字节

TRn:输出位置相对向右移动n个字符宽

SP、SS:SP,如果数值为正就加上+,SS取消SP功能

BN、BZ:应用在输入时,BN定义在输入时,没有数据的字节代表“没有东西”,BZ则定义在输入时,没有数据的字节代表“0”

kP:用来改变输入浮点数的Scale。输入的数值会自动乘上10^-k,使用指数类型输入就无效

Bw[.m]:将整数转换成二进制来输出,输出会占w个字符宽,固定输出m个数字,m值可以不给定

Ow[.m]:将整数转换成八进制来输出,输出会占w个字符宽,固定会输出m个数字,m值可以不给定

Zw.[m]:将整数转换成十六进制来输出,会占w个字符宽,输出m个数字,m值可以不指定

program main
    real a,b,c
    character(len=10) string
    character(len=20) string1
    string = "(I2)"
    a = 1.0
    b = 2.0
    c = 3.0
    write(*,"(3(1XF5.2))") a,b,c !3(1XF5.2)表示三次输出格式输出浮点数
    write(*,"('3+4=',I1)") 3+4  !混合使用单双引号
    write(*,string) 3
    read(*,"(A80)") string1 !定义好字节宽度,中间有空格也没有影响
    write(*,"(A80)") string1
    stop
end program main

image-20220308092420998

声明的其他事项

变量名称的取名策略

  • 26个英文字母+0~9数字+_(前缀必须是英文单词)

IMPLICIT命令

  • fortran标准中有一项不好的功能,使用implicit可以不用提前声明变量,这就容易导致输入错误。
  • 默认情况下,I、J、K、L、M、N为首字母将视为整数类型,其他变量被当成浮点数类型
program main
implicit interger(A,B,C) !A,B,C开头的变量为整型
implicit interger(A-F,I,K) !A-P,I,K开头的变量为整型
implicit real(M-P) !M-P开头的为浮点数
implicit none !关闭默认类型,所有变量提前声明
end program main

注意:推荐使用implicit none,需放在program命令的下一行,不能放在其他位置

常数的声明方法(PARAMETER)

  • 声明常量可以减少出错,改变常量将报错,其他将不变的变量设置为常量将会提高程序编译效率
program main
    implicit none
    real pi
    parameter(pi=3.1415926)
    ! 与上面定义等效  real,parameter :: pi = 3.1415926
    write(*,"(F4.2)") sin(pi/6)
    stop
end program main

image-20220308094826659

设置变量的初值

  • 第一种
program main
integer :: i = 1
real :: j = 1.0
character(len=10) :: string = "Hello"
complex :: a = (1.0,2.0)
write(*,*) i,j,string,a
end program main
  • 第二种
program main
implicit none
integer i
real j
complex k
character(len=10) s
data i,j,k,s / 1,1.0,(1.0,2.0),"Hello"/
write(*,*) i,j,k,s
stop
end program main

image-20220308100228119

等价声明(EQUIVALENCE)

  • 等价声明:将两个以上的变量,声明它们使用同一个内存地址,改变一个变量,其他变量也会改变
  • 节省内存
  • 精简代码
program main
implicit none
integer a,b
equivalence(a,b)  !a,b将公用一个内存空间
a = 10
write(*,*) a,b
b = 20
write(*,*) a,b
stop
end program main

image-20220308101136303

声明在程序中的结构

  • 声明的位置需放在程序代码的可执行描述之前
  • Data也是声明的一部分,设置初值,先声明类型

混合运算

  • 混合运算,就要考虑类型转换
program main
implicit none
integer :: a =1
integer :: b =2
real :: c
real :: d
c=a/b
write(*,"(F5.2)") c
c=real(a)/real(b)
write(*,"(F5.2)") c
c=d
write(*,"(F5.2)") c
end program main

image-20220308102421763

Fortran 90 的自定义数据类型

  • 自定义数据类型即对象,与C语言struct类似
type :: person
character(len=10) :: name
integer :: age
real :: length
real :: weight
character(len=20) :: address
end type person
  • 实例
program main
implicit none
type :: person
character(len=20) :: name
integer :: age
integer :: height
integer :: weight
character(len=80) :: address
end type person
type(person) :: a !new一个对象a
write(*,*) "name:"
read(*,*) a%name !读取元素
write(*,*) "age:"
read(*,*) a%age
write(*,*) "height:"
read(*,*) a%height
write(*,*) "weight"
read(*,*) a%weight
write(*,*) "address"
read(*,*) a%address
write(*,100) a%name,a%age,a%height,a%weight,a%address
100 format(/,"name:",A20/,"age:",I3/,"height:",I3/,"weight:",I3/,"address:",A80)
stop
end program main

image-20220308110611867

  • 以上赋值等效于a=person("张",24,170,120,"重庆市")

KIND的使用

  • kind声明整数型、浮点数变量所占用的内存字节
  • kind+fortran 90库函数,可以增加程序代码的“跨平台”能力

判断数值值域范围所需要的kind函数,使代码可以在不同的机器上正常运行

  • selected_int_kind(n):返回n位整数所对应的kind值,返回-1表示无法提供所需要的值域范围
  • selected_real_kind(n,e):返回n位有效位数,指数达到e为位的浮点数所需的kind值,返回-1表示无法满足所要求的有效位数,返回-2表示无法满足所要求的指数范围,返回-3表示两者都无法满足

实例

program main
implicit none
integer,parameter :: long_int = selected_int_kind(9)
integer,parameter :: short_int = selected_int_kind(3)
integer,parameter :: long_real = selected_real_kind(10,50)
integer,parameter :: short_real = selected_real_kind(3,3)
integer(kind=long_int) :: a = 1212345678
integer(kind=short_int) :: b = 12
real(kind=long_real) :: c = 1.23456789D45
real(kind=short_real) :: d = 1230
write(*,"(I3,1X,I10)") long_int,a
write(*,"(I3,1X,I10)") short_int,b
write(*,"(I3,1X,E10.5)") long_real,c
write(*,"(I3,1X,E10.5)") short_real,d
stop
end program main

image-20220308183029213

  • 除了在声明时运用kind,还可以赋值数字时,在数字后面用下划线,不常使用
program main
implicit none
real(kind=4) :: a
real(kind=8) :: b
a=1.0_4 !确保其为单精度
b=1.0_8  !确保其为双精度
write(*,*) a,b
stop
end program main

image-20220308141732392

习题

  • 1、请编写程序输出下面三行字符串

Have a good timeThat's not bad"Mary" isn't my name

program main
    character(len=20) :: string1 = "Have a good time."
    character(len=20) :: string2 = "That's not bad."
    character(len=20) :: string3 = '"Mary" isn''t my name.'
    write(*,"(A20/)") string1,string2,string3
    stop
end program main

image-20220309091339463

  • 2、编写一个可以让用户输入圆形的半径,并计算、输出这个圆形面积的程序。(请自行设计输出格式)
program main
    implicit none
    real,parameter :: pi = 3.1415926
    real(kind=4) :: r
    real(kind=4) :: s
    write(*,*) "请输入圆的半径:"
    read(*,*) r
    s = pi*r**2
    write(*,"(A9,F5.2,A16,F10.2)") "半径为",r,"的圆面积S=",s
    stop
end program main

image-20220309093829178

  • 3、某次物理期中考试的考题太难,老师决定调整全体学生的成绩,调整的公式是把原成绩开平方再乘以10。请写一个程序来读入一位学生的初始成绩,并计算出调整后的分数。
program main
    implicit none
    real(kind=4) result
    real(kind=4) reResult
    write(*,*) "请输入学生成绩:"
    read(*,*) result
    reResult = result**2*10
    write(*,"(A25,F10.2)") "该学生调整后成绩:",reResult
    stop
end program main

image-20220309095245872

  • 4、请问下面程序的输出结果是什么?
program main
integer a,b
real ra,rb
a=2
b=3
ra=2.0
rb=3.0
write(*,*) b/a ! 1
write(*,*) rb/ra ! 1.5
end program main

image-20220309095623879

5、Fortran 90的自定义类型,在主程序中定义一个新的类型distance。这个类型中有3个浮点数类型的元素,分别以米(m)、厘米(cm)、英寸(inch)为单位来记录同样的一段长度。请编写一个程序,程序会以公尺为单位来读入一段长度,并自动计算出其他单位的数值。

program main
    implicit none
    real :: mi
    type :: distance
        real :: m
        real :: cm
        real :: inch
    end type distance
    type(distance) :: d
    write(*,*) "请输入长度(m):"
    read(*,*) mi
    d = distance(mi,mi*100,mi/0.0254*1)
    write(*,*) d%m,"米等于",d%cm,"厘米等于",d%inch,"英寸"
    stop
end program main

image-20220309101553333

流程控制与逻辑运算

IF语句

IF基本用法

  • 语法
IF(判断条件) then
    ..........
    ..........
else if then
    ..........
else
   ...........
end IF
  • 练习
program main
    implicit none
    real(kind=4) :: speed
    write(*,*) "请输入车速:"
    read(*,*) speed
    if(speed>=100.0) then
        write(*,*) "slow down!"
    else if((speed>=50.0).and.(speed<=80.0)) then
        write(*,*) "normal"
    else if(speed<50) then
        write(*,*) "speed up!"
    end if
end program main

image-20220310092209416

逻辑运算

  • ==等于
  • /=不等于
  • >大于
  • >=大于或等于
  • <小于
  • <=小于或等于
  • .and.交集(两边都成立)
  • .or.并集(只要一边成立)
  • .not.逻辑反向
  • .eqv.表达式两边相等时,表达值成立
  • .neqv.表达式两边不等时,表达式成立

多重判断IF-ELSE-IF

  • 语法
if(条件1) then
   ..........
if(条件2) then
   ..........
if(条件3) then
   ..........
else
   ..........
end if

嵌套IF语句

  • 语法
if(条件)  then
   if(条件) then
       if(条件) then
       else if(条件) then
       end if
   else if(条件) then
   end if
end if
  • 实例:区分四象限
program main
    implicit none
    real(kind=4) :: x,y
    integer answer
    write(*,*) "请输入x值="
    read(*,*) x
    write(*,*) "请输入y值="
    read(*,*) y
    ! 或者直接读取坐标 write(*,*) "input(x,y)" read(*,*) x,y
    if(x>0) then
        if(y>0) then
            write(*,*) "(x,y)在第一象限"
        else
            write(*,*) "(x,y)在第四象限"
        end if
    else if(x<0) then
        if(y>0) then
            write(*,*) "(x,y)在第二象限"
        else
            write(*,*) "(x,y)在第三象限"
        end if
    else
        write(*,*) "(x,y)在轴上"
    end if
end program main

image-20220310122125215

浮点数及字符的逻辑运算

浮点数的逻辑判断

  • 使用浮点数进行判断时,要避免使用“等于”来判断,因为浮点数的有效位数是有限的。
program main
implicit none
real(kind=4) :: a,b=3.0
a=sqrt(b)**2-b
if(a==0.0) then
write(*,*) "a等于0"
else
write(*,*) "a不等于0"
end if
stop
end program main
  • 由于以上代码有可能输出“a不等于0”的情况

可改写成:

program main
implicit none
real(kind=4) :: a
real(kind=4) :: b = 4.0
real,parameter :: e =0.0001 !设置误差范围
a=sqrt(b)**2-b
if(abs(a-0.0)<=e) then
write(*,*) "a等于0"
else
write(*,*) "a不等于0"
end if
stop
end program main

字符的逻辑判断

  • 字符的逻辑判断就是比较字符的ASCII
program main
    implicit none
    character(len=20) :: str1,str2
    character relation
    write(*,*) "String 1:"
    read(*,"(A20)") str1
    write(*,*) "String 2:"
    read(*,"(A20)") str2
    if(str1>str2) then
        relation = ">"
    else if(str1==str2) then
        relation = "="
    else
        relation = "<"
    end if
    write(*,"('String1',A1,'String2')") relation
    stop
end program main

image-20220312091050407

SELECT CASE语句

  • 语法
select case(变量)
case(数值1)
·········
case(数值2)
·········
case(数值3)
·········
case default
·········
end select
  • 实例1
program main
    implicit none
    integer score
    character grade
    write(*,*) "score:"
    read(*,*) score
    select case(score)
    case(90:100)
        grade = "A"
    case(80:89)
        grade = "B"
    case(70:79)
        grade = "C"
    case(60:69)
        grade = "D"
    case(0:59)
        grade = "E"
    case default
        grade = "?"
    end select
    write(*,"('Grade:',A1)") grade
    stop
end program main

image-20220312092253851

  • 实例2
program main
    implicit none
    real a,b,ans
    character operator
    read(*,*) a
    read(*,"(A1)") operator ! 不赋值格式时,有些机器会读不到除号"/"
    read(*,*) b
    select case(operator)
    case("+")
        ans = a+b
    case("-")
        ans=a-b
    case("*")
        ans=a*b
    case("/")
        ans=a/b
    case default
        write(*,"('Unknown operator',A1)") operator
    end select
    write(*,"(F6.2,A1,F6.2,'=',F6.2)") a,operator,b,ans
    stop
end program main

image-20220312094004640

注意:select case存在限制:

  • 只能使用integercharacterlogical,不能使用浮点数及复数
  • case(数值)中的数值只能是常量,不能是变量

其他循环控制

IFSelect CaseGoToPauseContinueStop

GOTO

  • GoTo编写的程序在结构上很乱,导致程序代码不易阅读,所以建议不使用GoTo

  • GoTo主要实现跳转,所以会造成程序难以维护和修改

program main
    implicit none
    real(kind=4) :: height
    real(kind=4) :: weight
    write(*,*) "height:"
    read(*,*) height
    write(*,*) "weight:"
    read(*,*) weight
    if(weight>height-100) goto 200
100 write(*,*) "under control."
    goto 300
200 write(*,*) "too fat!"
300 stop
end program main

image-20220312095219138

program main
    implicit none
    integer I
    integer N
    parameter(N=10)
    data I /0/ !赋值I
10  write(*,"(1X,A3,I2)") "I=",I
    I=I+1
    if(I<N) goto 10
    stop
end program main

image-20220312095946041

program main
    implicit none
    integer I
    integer N
    data I,N /2,1/
    goto(10,20,30) I/N
10  write(*,*) "I/N=1"
    goto 100
20  write(*,*) "I/N=2"
    goto 100
30  write(*,*) "I/N=3"
100 stop
end program main

image-20220312100815666

IF与GOTO的联用

  • IF判断还有一种叫做算术判断的方法,它的做法跟goto有点类似
  • 实例
program main
    implicit none
    real a,b,c
    data a,b /2.0,1.0/
    c=a-b
    if(c) 10,20,30 !c<0,10;c=0,20;c>0,30
10  write(*,*) 'a<b'
    goto 40
20  write(*,*) 'a=b'
    goto 40
30  write(*,*) 'a>b'
40  stop
end program main

image-20220312102403577

PAUSE,CONTINUE,STOP

  • pause:暂停执行,直到用户按下enter键才会继续执行。当输出多页时,可以在换页处添加pause,用户按enter键继续阅读
  • continue:没有实际的用途,它的功能就是“继续向下执行程序”
  • stop:结束程序执行,遇到不合理的地方,输出警告后停止执行程序

二进制的逻辑运算

  • fortran 90的库函数中,Iand用来做二进制的and计算,Ior用来做二进制的or计算
a = 2 !二进制010
b = 4 !二进制100
c = iand(a,b)
d = ior(a,b)

image-20220312104813512

  • 进制转换
integer :: a
a = B"10" !B(binary)表示二进制,a=2(十进制)
b = O"10" !O(octal)表示十进制,b=8(十进制)
c = Z"10" !Z()表示十六进制,z=16(十进制)

习题

  • 1、假如所得税有3个等级,月收入在1000元以下的税率为3%,在1000元至5000元之间的税率为10%,在5000元以上的税率为15%。请写一个程序来输入一位上班族的月收入,并计算他(她)所应缴纳的税金。
program main
    implicit none
    real(kind=4) :: pay
    real(kind=4) :: tax !税金
    write(*,*) "请输入月收入:"
    read(*,*) pay
    if(pay<1000) then
        tax = pay*0.03
    else if((pay>=1000).or.(pay<=5000)) then
    tax = pay*0.1
    else
    tax = pay*0.15
    end if
    write(*,*) "缴纳税金:",tax
stop
end program main

image-20220312112135937

  • 2、某电视台的晚上8点节目安排如下:

​ 星期一、四:新闻

​ 星期二、五:电视剧

​ 星期三、六:卡通片

​ 星期日:电影

​ 请写一个程序,可以输入星期几来查询当天晚上的节目

program main
    implicit none
    character(len=9) week
    write(*,*) "请输入星期几(格式:星期一):"
    read(*,*) week
    select case(week)
    case("星期一","星期四")
        write(*,*) "新闻"
    case("星期二","星期五")
        write(*,*) "电视剧"
    case("星期三","星期六")
        write(*,*) "卡通片"
    case("星期日")
        write(*,*) "电影"
    case default
        write(*,*) "输入错误!"
    end select
    stop
end program main

image-20220312113800059

  • 3、假如所得税有三个等级,而且随年龄不同又有不同的算法:

    第一类:年轻级(不满50岁)

    月收入在1000元以下的税率为3%,在1000元至5000元之间的税率为10%,在5000元以上的税率为15%。

    第二类:年老级(50岁以上)

    月收入在1000元以下的税率为5%,在1000元至5000元之间的税率为7%,在5000元以上的税率为10%。

    请编写一个程序来输入一位上班族的年龄、年收入,并计算他(她)所应缴纳的税金。

program main
    implicit none
    integer age
    real(kind=4) income,tax,pay
    write(*,*) "请输入年龄:"
    read(*,*) age
    write(*,*) "请输入年收入:"
    read(*,*) income
    pay = income/12.0
    if(age<=50) then
        if(pay<1000) then
            tax = pay*0.03*12
        else if((pay>=1000).and.(pay<=5000)) then
            tax = pay*0.1*12
        else
            tax = pay*0.15*12
        end if
    else
        if(pay<1000) then
            tax = pay*0.05*12
        else if((pay>=1000).and.(pay<=5000)) then
            tax = pay*0.07*12
        else
            tax = pay*0.1*12
        end if
    end if
    write(*,*) "应缴纳年税金:",tax
    stop
end program main

image-20220312170624002

  • 4、在一年的当中,通常有365天。但是如果是润年时,一年则有366天。在公历中,润年的策略如下:(以公元来记年)

​ 年数是4的倍数时,是润年。

​ 年数是100的倍数时是例外,不当润年记。除非它刚好又是400的倍数。

编写一个程序,让用户输入一个公元的年份,然后交给程序来判断这一年当中会有多少天。

program main
    implicit none
    integer(kind=2) :: year
    write(*,*) "请输入公元:"
    read(*,*) year
    if(((mod(year,4)==0).or.((mod(year,100)==0).and.(mod(year,400)==0))))  then
        write(*,*) year,"是润年,有366天"
    else
        write(*,*) year,"不是润年,有365 天"

    end if
    stop
end program main

image-20220314101807695

循环

DO

  • 语法
do counter = 1,lines,1 !最后的数字是计数器的增值,每执行一次循环,counter就会累加上这个数值,可省略,内定值为1,counter<=lines时会继续重复循环,counter为计数器
   ············
   ············
end do
  • 实例(输出10遍Happy Every Day!!!)
program main
    implicit none
    integer(kind=2) :: counter
    integer,parameter :: lines = 10
    do counter = 1,lines,1
        write(*,*) "Happy Every Day!!!"
    end do
    stop
end program main

image-20220314103044090

  • 实例(计算2+4+6+8+10)
program main
    implicit none
    integer,parameter :: limit = 10
    integer(kind=2) :: counter,sum = 0
    do counter = 2,limit,2
        sum = sum + counter
    end do
    write(*,*) "2+4+6+8+10=",sum
    stop
end program main

image-20220314110731602

  • 循环的增值并没有规定一定要是正数,也可以是负数,让计算器一直递减下去,这个时候,就需要循环的计数器终止值必须小于计数器起始值。

  • DO循环的嵌套

do i=1,10
   do j=1,10
      do k=1,10
      end do
   end do
end do
  • 实例
program main
    implicit none
    integer,parameter :: lines = 3
    integer(kind=2) i,j
    do i=1,lines
        do j=1,lines
            write(*,"(I2,I2)") i,j
        end do
        write(*,*) "another cycle"
    end do
    stop
end program main

image-20220314112634861

DO WHILE 循环

  • 循环不一定由计数器的增、减来决定是否结束循环,还可以由逻辑运算来决定结束循环。

  • 语法

do while(逻辑运算)  !逻辑运算成立时,会一直重复执行循环
   ···········
   ···········
end do
  • 实例(计算2+4+6+8+10)
program main
    implicit none
    integer,parameter :: limit = 10
    integer counter
    integer :: sum = 0
    do while(counter<=limit)
        sum = sum + counter
        counter=counter+2
    end do
    write(*,*) "2+4+6+8+10=",sum
    stop
end program main

image-20220314150553289

  • 实例(程序猜测蔡小姐体重)
program main
    implicit none
    real,parameter :: weight = 45.0
    real,parameter :: e = 0.001 !误差
    real :: guess = 0.0
    do while(abs(weight-guess>e))
        write(*,*) "请输入蔡小姐的体重:"
        read(*,*) guess
    end do
    write(*,*) "恭喜你猜对了,蔡小姐体重为",weight
    stop
end program main

image-20220314151458212

循环的流程控制

cycleexit

CYCLE

  • cycle可以略过循环的程序模块,在cycle命令后面的所有程序代码,直接跳回循环的开头来进行下一次循环

  • 实例(9层楼,四层不停,从一层到九层灯号显示情况)

program main
    implicit none
    integer,parameter :: dest = 9
    integer,parameter :: stop = 4
    integer :: floor
    do floor=1,dest
        if(floor==stop) cycle
        write(*,*) floor,"灯亮"
    end do
    stop
end program main

image-20220314153227050

EXIT

  • exit 可以直接“跳出”一个正在运行的循环,不论是do循环还是do while循环
  • 实例(上猜体重)
program main
    implicit none
    real,parameter :: weight = 45.0
    real,parameter :: e = 0.001 !误差
    real :: guess = 0.0
    do while(.true.)
        write(*,*) "请输入蔡小姐的体重:"
        read(*,*) guess
        if(weight-guess<e) exit
    end do
    write(*,*) "恭喜你猜对了,蔡小姐体重为",weight
    stop
end program main

image-20220314154118256

署名的循环

  • 循环还可以取”名字”,这个用途是可以在编写循环时就能明白地知道end do这个描述的位置是否正确,尤其是在多层的循环当中,署名的循环也可以配合cycleexit使用。
  • 实例
program main
    implicit none
    integer :: i,j
    outter: do i=1,3
        inner: do j=1,3
            write(*,"('(',i2,',',i2,')')") i,j
        end do inner
    end do outter
    stop
end program main

image-20220314155011062

循环的应用

  • 实例一

求等差数列1+2+3+4+···+99+100的值

program main
    implicit none
    integer,parameter :: limit = 100
    integer :: counter
    integer :: sum
    do counter=1,limit
        sum = sum + counter
    end do
    write(*,*) "1+2+3+····+99+100=",sum
    stop
end program main

image-20220314155646977

  • 实例2

费氏数列(Fibonacci Sequence)的数列规则如下:

F0=0,f1=1,当n>1时fn=fn-1+fn-2

费氏数列的前10个数字列举如下:” 0 1 1 2 3 5 8 13 21 34 “,请编写程序来计算费氏数列的前10个数字。

program main
    implicit none
    integer,parameter :: limit = 9
    integer :: counter,fn = 0
    integer :: fn_1=1
    integer :: fn_2=0
    do counter = 0,limit
        if((counter==0).or.(counter==1)) then
            write(*,*) counter
        else
            fn = fn_1 + fn_2
            write(*,*) fn
            fn_2 = fn_1
            fn_1 = fn
        end if
    end do
    stop
end program main

image-20220314163132058

  • 实例3

密码加密与解密,把每个英文字母在ASCII表中的编号加上2所得到的字母当成密码来传输。例如:abc加密后成为cde,解密的工作就是把上述操作还原,把cde解密回abc

program main
    implicit none
    integer i,strlen
    integer,parameter :: key = 2
    character(len=20) :: string
    write(*,*) "请输入需要加密的字符串:"
    read(*,*) string
    strlen = len_trim(string) !实际输入长度
    do i = 1,strlen
        string(i:i) = char(ichar(string(i:i))+key)
    end do
    write(*,*) "加密后:",string
    stop
end program main

image-20220314164949766

program main
    implicit none
    integer i,strlen
    integer,parameter :: key = 2
    character(len=20) :: string
    write(*,*) "请输入需要解密的字符串:"
    read(*,*) string
    strlen = len_trim(string) !实际输入长度
    do i = 1,strlen
        string(i:i) = char(ichar(string(i:i))-key)
    end do
    write(*,*) "解密后:",string
    stop
end program main

image-20220314170154247

  • 实例4

写一个小型的计算机程序,用户可以输入两个数字及一个运算符号来决定要把这两个数字做加减乘除的其中一项运算,每做完一次计算后,让用户来决定要再做新的计算或是结束程序。

program main
    implicit none
    real :: a,b,ans
    character :: c
    logical :: state = .true.
    do while(state)
        write(*,*) "请输入第一个数字:"
        read(*,*) a
        write(*,*) "请输入运算符:"
        read(*,"(A1)") c
        write(*,*) "请输入第二个数字:"
        read(*,*) b
        select case(c)
        case('+')
            ans = a+b
            write(*,*) "a+b=",ans
        case('-')
            ans = a-b
            write(*,*) "a-b=",ans
        case('*')
            ans = a*b
            write(*,*) "a*b=",ans
        case('/')
            ans = a/b
            write(*,*) "a/b=",ans
        case default
            write(*,*) "运算符输入有误!"
        end select
        write(*,*) "是否继续计算,输入.false.结束,输入.true.继续:"
        read(*,*) state
    end do
    stop
end program main

image-20220314180043321

习题

1、以循环来连续显示5行的Fortran字符串,输出结果如下:

Fortran

Fortran

Fortran

Fortran

Fortran

program main
    implicit none
    integer,parameter :: limit = 5
    integer :: counter
    do counter = 1,limit
        write(*,*) "Fortran"
    end do
    stop
end program main

image-20220316085432984

  • 2、以循环来计算等差数列1+3+5+7+···+99的结果
program main
    implicit none
    integer :: i=1,ans
    do while(i<=99)
        ans = ans +i
        i=i+2
    end do
    write(*,*) "1+3+5+7+···+99=",ans
    stop
end program main

image-20220316085935161

  • 3、改变一下猜小姐体重的程序条件。让程序最多只准许用户猜5次,5次之内猜不中就不能再猜下去。也就是这个循环最多执行5次就会结束,不过要是在5次之内就猜对,也要跳出循环。程序最后还要显示信息来告诉用户有没有猜对。
program main
implicit none
integer,parameter :: limit = 5
real,parameter :: weight = 45.0
real :: guess,e = 0.0001
integer :: counter
do counter = 1,limit
write(*,*) "请输入蔡小姐体重:"
read(*,*) guess
if(abs(weight-guess)<=e) then
write(*,*) "恭喜您猜对了,蔡小姐体重为",weight
stop
else if(counter == 5) then
write(*,*) "game over!"
end if
end do
stop
end program main

image-20220316091512790

  • 4、以循环来计算1/1!+1/2!+1/3!+1/4!+···+1/10!的值。
program main
    implicit none
    integer,parameter :: limit = 10
    integer :: counter,i,ans2 = 1
    real :: ans1=0
    do counter = 1,limit,1
        ans2 = 1
        do i = counter,1,-1
            ans2 = ans2*i
        end do
        ans1 = ans1+1.0/ans2
    end do
    write(*,*) "1/1!+1/2!+1/3!+1/4!+···+1/10!=",ans1
    stop
end program main

image-20220316095104922

  • 5、写一个程序,让用户输入一个内含空格符的字符串,然后使用循环把字符串中的空格符消除后再重新输出。例如:

Happy New Year

HappyNewYear

program main
    implicit none
    character(len=50) :: string
    integer :: String_len,counter,i
    write(*,*) "请输入字符串:"
    read(*,*) string
    String_len = len_trim(string)
    do counter = 1,String_len
        if(string(counter:counter)==' ') then
            do i = counter,String_len
                string(i:i) = string(i+1:i+1)
            end do
        end if
    end do
    write(*,*) "删除空格后的字符串为:",string
    stop
end program main

image-20220316103118937

数组(Array)

基本使用

  • 数组(Array)是一种使用数组的方法。它可以配合循环等的功能,用很精简的程序代码来处理大量数据。

一维数组

  • 语法

Datatype name(size)

Datatype:数组的类型,除了4种基本类型,还可以使用type自订出来的类型

name:数组变量的名字

size:数组的大小,必须使用整型常数

  • 写一个程序来记录全班5位同学的数学成绩,并提供由座号来查询成绩的功能
program main
    implicit none
    integer,parameter :: studentnum = 5
    real :: student(studentnum)
    integer counter,search
    do counter = 1,studentnum,1
        write(*,*) "请输入学生座位号为:",counter,"的数学成绩:"
        read(*,*) student(counter)
    end do
    do while(.true.)
        write(*,*) "请输入你想查询的学生座位号:"
        read(*,*) search
        if(search<=0.or.search>studentnum) stop
        write(*,*) "座位号为",search,"的学生成绩为:",student(search)
    end do
    stop
end program main

image-20220316112452913

  • 自定义类型的数组语法
Type :: student
   integer :: desk
   real :: score
end type
Type(person) :: a(10)

二维数组

  • 语法

Datatype name(x,y)

  • 保存整个年级5个班级,每班5位同学的数学考试成绩。
program main
    implicit none
    integer,parameter :: classes=5,studentnum=5
    integer :: counter1,counter2
    real :: student(classes,studentnum)
    real :: score
    do counter1 = 1,classes
        do counter2 = 1,studentnum
            write(*,*) "请输入",counter1,"班级的第",counter2,"个学生成绩"
            read(*,*) score
            student(counter1,counter2) = score
        end do
    end do
    do while(.true.)
        write(*,*) "请输入班级:"
        read(*,*) counter1
        write(*,*) "请输入学生座位号:"
        read(*,*) counter2
        if(counter1<=0.or.counter1>classes.or.counter2<=0.or.counter2>studentnum) stop
        write(*,*) counter1,"班级的",counter2,"学生的成绩为:",student(counter1,counter2)
    end do
    stop
end program main

image-20220316152109615

  • 用户输入两个2X2矩阵的值,再把这两个矩阵相加
program main
    implicit none
    integer,parameter :: row=2,col=2
    integer :: counter1,counter2
    real :: matrix1(row,col)
    real :: matrix2(row,col)
    real :: matrix(row,col)
    do counter1 = 1,row
        do counter2 = 1,col
            write(*,*) "请输入第一个矩阵的第",counter1,"行第",counter2,"列的元素:"
            read(*,*) matrix1(counter1,counter2)
        end do
    end do
    do counter1 = 1,row
        do counter2 = 1,col
            write(*,*) "请输入第二个矩阵的第",counter1,"行第",counter2,"列的元素:"
            read(*,*) matrix2(counter1,counter2)
        end do
    end do
    write(*,*) "matrix1+matrix2="
    do counter1 = 1,row
        do counter2 = 1,col
            matrix(counter1,counter2) = matrix1(counter1,counter2)+matrix2(counter1,counter2)
            write(*,"('(',I1,',',I1,')=',f5.2)") counter1,counter2,matrix(counter1,counter2)
        end do
    end do
    stop
end program main

image-20220316163504354

多维数组

  • Fortran最多可以声明七维的数组
  • integer a(D1,D2,···,Dn)
  • 将前面二阶矩阵相加,改成三阶矩阵相加
program main
    implicit none
    integer,parameter :: row=2,col=2,space=2
    integer :: counter1,counter2,counter3
    real :: matrix1(row,col,space)
    real :: matrix2(row,col,space)
    real :: matrix(row,col,space)
    do counter1 = 1,row
        do counter2 = 1,col
            do counter3 =1,space
                write(*,*) "请输入第一个矩阵的第",counter1,"行第",counter2,"列第",counter3,"列的元素:"
                read(*,*) matrix1(counter1,counter2,counter3)
            end do
        end do
    end do
    do counter1 = 1,row
        do counter2 = 1,col
            do counter3 =1,space
                write(*,*) "请输入第二个矩阵的第",counter1,"行第",counter2,"列第",counter3,"列的元素:"
                read(*,*) matrix2(counter1,counter2,counter3)
            end do
        end do
    end do
    write(*,*) "matrix1+matrix2="
    do counter1 = 1,row
        do counter2 = 1,col
            do counter3 =1,space
                matrix(counter1,counter2,counter3) = matrix1(counter1,counter2,counter3)+matrix2(counter1,counter2,counter3)
                write(*,"('(',I1,',',I1,',',I1,')=',f5.2)") counter1,counter2,counter3,matrix(counter1,counter2,counter3)
            end do
        end do
    end do
    stop
end program main

image-20220316173529408

另类的数组声明

  • 在没有特别赋值的情况中,数组的索引值都是由1开始

integer a(5)a(1)a(2)a(3)a(4)a(5)

  • 使用特别声明,改变默认的规则,例如特别赋值数组的坐标使用范围

integer a(a:5)a(0)~~a(5)

integer a(5,0:5)

integer a(2:3,-1:3)

数组内容的设置

赋初值

  • 使用DATA来设置数组的初值
integer A(5)
integer I,J
DATA A /1,2,3,4,5/
DATA A /5*3/ !A(1)到A(5)都是3
DATA(A(I),I=2,4,1)/2,3,4/ !隐含式循环,I会从2增加到4,按照顺序到后面取数字即A(2) =2,···
write(*,*) (A(I),I=2,4,1) !读取数组元素,A(2),A(3),A(4)
DATA((A(I,J),I=1,2),J=1,2)/1,2,3,4/ !A(1,1)=1,A(2,1)=2,先循环I,再循环J
integer :: A(5) = (/1,2,3,4,5/) !fortran 90省略DATA,注意括号与除号间不能有空格
  • 查询学生成绩,改成直接记录在程序代码中
program main
implicit none
integer,parameter :: studentsum = 5
integer :: student(studentsum) = (/80,90,85,75,95/)
integer i
do while(.true.)
write(*,*) “请输入学生座位号:”
read(*,*) i
if(i<=0.or.i>studentnum) stop
write(*,*) "座位号为",i,"的学生成绩为:",student(i)
end do
stop
end program main

image-20220316200157219

  • 使用隐含式循环功能,设置二维数组的矩阵内容,再把它输出到屏幕上
program main
    implicit none
    integer,parameter :: row =2
    integer,parameter :: col =2
    integer :: matrix(row,col)
    integer r,c
    data((matrix(r,c),r=1,2),c=1,2)/1,2,3,4/
    write(*,"(I3,I3,/,I3,I3)") ((matrix(r,c),r=1,2),c=1,2)
    stop
end program main

image-20220316200815310

对整个数组的操作

  • a=5:a是一个任意维数及大小的数组,这个命令把数组a每个元素赋值成5

  • a=(/1,2,3/):等号右边所提供的数字数目,必须跟数组a的大小一样

  • a=ba=b+ca=b-ca=b*ca=b/c:数组a,b,c必须是相同维数

  • a=sin(b):矩阵a的每一个元素为矩阵b元素的sin值,数组b必须是浮点数类型,才能使用sin函数

  • a=b>c:a,b,c为相同维数的数组,数组a为逻辑类型数组,b,c为数值类型数组

  • 实例:矩阵相加

program main
    implicit none
    integer,parameter :: row =2
    integer,parameter :: col =2
    integer :: ma(row,col) = 1
    integer :: mb(row,col) = 4
    integer :: mc(row,col)
    integer :: i,j
    mc = ma + mb !矩阵相加
    write(*,"(I3,I3,/,I3,I3)") ((mc(i,j),i=1,2),j=1,2)
    stop
end program main

image-20220317085848954

对部分数组的操作

  • a(3:5)=5

  • a(3:)=5

  • a(3:5)=(/3,4,5/)

  • a(1:3)=b(4:6)

  • a(1:5:3)=3:a(1),a(3),a(5) =3

  • a(1:10) =a(10:1:-1):翻转数组元素

  • a(:)=b(:,2):取数组b第二列元素,赋值给数组a

  • a(:,:)=b(:,:,1)

    • 等号两边所使用的数组元素数目要一样多
    • 同时使用多个隐含式循环时,较低维的循环可以想象成是内层的循环
    integer :: a(2,2),b(2,2)
    b=a(2:1:-1,1:2:1)
    !低维是内层循环,会先执行,所以这个命令结果为
    !b(1,1) = a(2,2)
    !b(2,1) = a(1,2)
    !b(1,2) = a(2,1)
    !b(2,2) = a(1,1)
  • 实例

program main
    implicit none
    integer,parameter :: row = 2
    integer,parameter :: col = 2
    integer r1,c1
    integer :: a(2,2)
    integer :: b(4) = (/5,6,7,8/)
    integer :: c(2)
    data((a(r1,c1),r1=1,2),c1=1,2) /1,2,3,4/
    write(*,*) a
    write(*,*) a(:,1)
    c=a(:,1)
    write(*,*) c
    c=a(2,:)
    write(*,*) c
    write(*,*) c(2:1:-1)
    c=b(1:4:2)
    write(*,*) c
    stop
end program main

image-20220317094225472

WHERE

  • WHERE命令可以经过逻辑判断来使用数组的一部分元素
program main
    implicit none
    integer :: i
    integer :: a(5) = (/(i,i=1,5)/)
    integer :: b(5) = 0
    where(a<3)
        b=a
    end where
    write(*,"(5(I3,1X))") b
    stop
end program main

image-20220317095136077

  • where elsewhere
program main
    implicit none
    integer :: i
    integer :: a(5)=(/1,2,3,4,5/)
    integer :: b(5)=0
    where(a<3)
        b=1
    elsewhere
        b=2
    endwhere
    write(*,"(5(I3,1X))") b
    stop
end program main

image-20220317095535565

  • where也可以实现嵌套,它也可以跟循环一样取名字,不过取名字的where描述在结束时,end where后面一定要接上它的名字,用来明确赋值所要结束的是哪一个where模块
    name:where(a<5) !取名
        b=5
    end where name

    where(a<5)
        where(a/=2)
            b=3
        elsewhere
            b=1
        end where
    elsewhere
        b=0
    end where
  • 假设年所得3万以下所得税率为10%,3万到5万之间为12%,5万以上为15%。使用where命令来计算,并记录10个人的所得税金额。
program main
    implicit none
    integer :: i
    real :: income(10)=(/25000,30000,50000,40000,35000,60000,27000,45000,20000,70000/)
    real :: tax(10) = 0
    where(income<30000.0)
        tax = income*0.10
    elsewhere(income<50000.0)
        tax=income*0.12
    elsewhere
        tax=income*0.15
    end where
    write(*,"(10(F8.1,1X))") tax
    stop
end program main

image-20220317102155561

FORALL

  • forall是fortran 95添加的功能,它可以看成是一种使用隐含式循环来使用数组的方法。

  • 语法

  • forall(triplet1[,triplet2[,triplet3···]],mask)
    !triple是用来赋值数组坐标范围的值,数组有几个维度就有多少个triple
    !mask是用来做条件判断,跟where的条件判断类似
        ·········
    end forall
    !例一
    integer :: a(5,5)
    integer :: i,j
    forall(i=1:5,j=1:5,a(i,j)<10)
      a(i,j) = 1
    end forall
program main
    implicit none
    integer :: i
    integer :: a(5)
    forall(i=1:5)
        a(i) = 5
    end forall
    write(*,*) a
    forall(i=1:5)
        a(i) = i
    end forall
    write(*,*) a
    stop
end program main

image-20220317102657116

  • 实例:声明了一个二维数组作为二维矩阵使用,它使用Forall命令把矩阵的上半部设置为1,对角线部分设置成2,下半部设置成3。
program main
    implicit none
    integer i,j
    integer,parameter :: size = 5
    integer :: matrix(size,size)
    forall(i=1:size,j=1:size,i>j) matrix(i,j) = 1 !上半部分
    forall(i=1:size,j=1:size,i==j) matrix(i,j) = 2 !对角线
    forall(i=1:size,j=1:size,i<j) matrix(i,j) = 3 !下半部分
    write(*,"(5(5I5,/))") matrix
    stop
end program main

image-20220317105204398

注意⚠️forall(i=1:size,j=1:size,i>j)低维是内层循环,会先执行,结果为

image-20220317105810399

  • forall可以写成多层的嵌套结构,它里面也只能出现跟设置数组数值相关的程序命令,还可以在forall中使用where。不过where当中不能使用forall
!嵌套
forall(I=1:5)
  forall(J=1:5)
    a(I,J) = 2
  end forall
  forall(J=6:10)
    a(I,J) = 2
  end forall
end forall

! forall中使用where
forall(i=1:5)
  where(a(:,j)/=0)
    a(:,i) = 1.0/a(:,i)
  end where
end forall

数组的保存规则

  • 一维数组在内存单元中线性排列
  • 多维数组在内存的连续模块中排列情况是以一种称为Column Major的方法来排列
    • Column Major对于二维数组来说就是,数组存放在内存中,会先放入第二维的Column中每个第一维Row的元素,第一个column放完再放第二个。
    • A(1,1)=>A(2,1)=>A(3,1)
    • ······································
    • 多维数组时,integer :: a(最低维,····,最高维),先循环低维,再依次到高维
!不推荐使用--a(I,J)与a(I,J+1)在内存中并不是连续的
do I=1,N
  do J=1,M
    a(I,J) = ····
  end do
end do 

!推荐使用--a(J,I)与a(J+1,I)在内存中连续,效率更高
do I=1,N
   do J=1,M
      a(J,I)=·····
   end do
end do

可变大小的数组

  • 某些情况下,要等到程序执行之后,才会知道所需使用的数组大小
  • 用户输入数组元素后,程序代码再声明一个刚刚好大小的数组来使用
program main
implicit none
integer :: students
integer,allocatable :: a(:) !声明一个可变大小的一维数组
integer :: i
write(*,*) "How many students:"
read(*,*) students
allocate(a(students)) !配置内存空间
do i = 1,students
write(*,"('Number',I3)") i
read(*,*) a(i)
end do
do i = 1,students
write(*,"('a(',I3,')=',I3)") i,a(i)
end do
stop
end program main

image-20220317192959036

  • 使用可变大小的数组
    • 1、声明allocatable :: a(:,:,···)
    • 2、配置内存大小allocate(a(number))
    • deallocate逆向运行,可以将配置好的内存释放掉
    • 当计算机内存满载时,可以使用语句判断allocate(a(100),stat=error),error是事先声明好的整型变量,如果error等于0则表示成功配置,不等于0则表示失败
  • 写一个程序测试计算机能承受多大的数组
program main
    implicit none
    integer :: size,error =0
    integer,parameter :: one_mb =1024*1024 ! 1MB
    character,allocatable :: a(:)
    do while(.true.)
        size = size + one_mb !一次增加1MB个字符,也就是1MB的内存空间
        allocate(a(size),stat=error)
        if(error/=0) exit
        write(*,"('Allocate',I10,'bytes')") size
        write(*,"(F10.2,'MB used')") real(size)/real(one_mb)
        deallocate(a)
    end do
    stop
end program main

image-20220317195203918

  • allocate相关的函数还有allocated,它用来检查一个可变大小的矩阵是否已经配置内存来使用,它会返回一个逻辑值。
if(.not. allocated(a)) then
    allocate(a(5)) !a(5)如果没有配置空间,就会配置
end if

数组的应用

  • 如何把一堆数字按照它们的大小来排序?排序有很多种算法可以使用,这里示范一个最简单的排序方法,叫做选择排序法
program main
    implicit none
    integer,parameter :: size = 10
    integer :: a(size) = (/5,3,6,4,8,7,1,9,2,10/)
    integer :: i,j
    integer :: t
    do i=1,size-1
        do j=i+1,size
            if(a(i) > a(j)) then
                t=a(i)
                a(i)=a(j)
                a(j)=t
            end if
        end do
    end do
    write(*,"(10I4)") a
    stop
end program main

image-20220317201526788

  • 做一个矩阵相乘的程序。矩阵相乘有它的特殊规则,假设现在有两个二维矩阵A、B,其中A的大小是L*M,B的大小是M*N。现在要计算C=A*B,C矩阵的大小一定是L*N。矩阵乘法的规则为:

image-20220317213118236

program main
    implicit none
    integer,parameter :: L=3,M=4,N=2
    real :: A(L,M)
    real :: B(M,N)
    real :: C(L,N)
    integer :: i,j,k
    data((A(i,j),i=1,3,1),j=1,4,1)/1,2,3,4,5,6,7,8,9,10,11,12/
    data((B(i,j),i=1,4,1),j=1,2) /1,2,3,4,5,6,7,8/
    do i=1,L
        do j=1,N
            C(i,j) = 0.0
            do k=1,M
                C(i,j) = C(i,j) + A(i,k)*B(k,j)
            end do
        end do
    end do
    do i=1,L
        write(*,*) C(i,:)
    end do
    stop
end program main

image-20220318083703207

  • Fortran 90 库函数提供matmul这个函数来做矩阵乘法,C=matmul(A,B)

习题

  • 1、请声明一个大小为10的一维数组,它们的初值为a(1)=2,a(2)=4,a(3)=6,····a(i)=2*i,并计算数组中这10数字的平均值。
program main
    implicit none
    integer,parameter :: number=10
    integer,allocatable :: a(:)
    integer :: i=1,j=2,sum=0
    real :: average
    if(.not. allocated(a)) then
        allocate(a(number))
    end if
    do while(i<=number)
        a(i) = j
        sum=sum+a(i)
        i=i+1
        j=j+2
    end do
    average = sum/10.0
    write(*,"(10I3,/,'平均值为:',F5.2)"),a,average
    stop
end program main

image-20220318090241113

  • 2、略
  • 3、编写一个程序计算前10个费氏数列,并把它们按顺序保存在一个一维数组当中。费氏数列(Fibonacci Sequence)的数列规则如下:

​ f(0) = 0,F(1)=1,当n>1时,f(n)=f(n-1)+f(n-2)

program main
    implicit none
    integer,parameter :: number = 10
    integer,allocatable :: a(:)
    integer :: f(10)
    integer :: i,j
    if(.not.allocated(a)) then
        allocate(a(number))
    end if
    f(1) = 0
    f(2) = 1
    do j=3,number,1
        f(j) = f(j-1) + f(j-2)
    end do
    do i=1,number,1
        a(i) = f(i)
    end do
    write(*,"(10I3)") a
    stop
end program main

image-20220318092334041

  • 4、把排序程序程序的应用第一个实例由小到大,改成从大到小。
program main
    implicit none
    integer,parameter :: size = 10
    integer :: a(size) = (/5,3,6,4,8,7,1,9,2,10/)
    integer :: i,j
    integer :: t
    do i=1,size-1
        do j=i+1,size
            if(a(i) < a(j)) then
                t=a(i)
                a(i)=a(j)
                a(j)=t
            end if
        end do
    end do
    write(*,"(10I4)") a
    stop
end program main

image-20220318100533275

  • 5 、声明为integer a(5,5)的二维数组,请问a(2,2) 跟a(3,3)在所配置的内存中是排行第几个位置?

data((a(i,j),i=1,5,1),j=1,5,1)先低维循环再高维循环,即a(1,1),a(2,1),a(3,1),·····

所以a(2,2)在第七个位置,a(3,3)在18个位置。

函数

  • 函数是自定义函数和子程序的统称

image-20220307182156253

image-20220307182216824

子程序(SUBROUTINE)的使用

  • 语法

    program main !主程序
       ········
       ········
    end program main
    
    subroutine sub1() !第一个子程序
       ········
       ········
    end subroutine sub
    
    subroutine sub2()  !第二个子程序
       ········
       ········
    end subroutine sub2
  • 主程序不一定要写在最开始的地方,也可以写在其他位置

    • 子程序最后一个命令通常是return ,表示返回原程序,可以写在子程序任何地方,提早返回,其中return语句可以省略
    • 如果子程序使用stop,程序就会停止,这不是子程序想要的结果,所以通常子程序不使用 stop
    • 主程序可以调用子程序,子程序也可以调用子程序,也可以自己调用自己(递归)
    • 子程序独立地拥有自己的变量声明,主程序和子程序使用相同的变量名称,是互不影响的
    • 子程序与主程序也拥有自己的“行代码”,互不干扰
  • 子程序可以用来独立出某一段具有特定功能的程序代码,需要时通过call命令调用。

program main
    implicit none
    call message() !调用子程序message
    call message()
    stop
end program main

! 子程序定义
subroutine message()
    implicit none
    write(*,*) "Hello World!"
    return
end subroutine message

image-20220318151718884

  • 实例

    program main
       implicit none
       call sub1()
       call sub2()
       stop
    end program main
    
    subroutine sub1()
       implicit none
       write(*,*) "This is sub1"
       call sub2()  !子程序调用子程序
       return
       stop
    end subroutine sub1
    
    subroutine sub2()
       implicit none
       write(*,*) "This is sub2"
       return
       stop
    end subroutine sub2

image-20220319133821646

  • 实例
program main
    implicit none
    integer :: a=1
    call sub1()
    write(*,"('a=',I2)") a
    stop
end program main

subroutine sub1()
    implicit none
    integer :: a=2
    write(*,"('a=',I2)") a
    return
end subroutine sub1

image-20220319133929635

  • 实例(传参)
program main
    implicit none
    integer :: a=1
    integer :: b=2
    call add(a,b)  !把变量a及b交给子程序add来处理
    stop
end program main

subroutine add(first,second)
    implicit none
    integer :: first,second
    write(*,*) first+second
    return
end subroutine add

image-20220319134031930

  • 注意:主程序传参给子程序时,使用的是地址调用(call by address/call by reference),即可以视为同一个变量,C语言使用的是传值调用(call by value),所以C语言就有了指针
  • 实例( 地址调用)
program main
implicit none
integer :: a=1
integer :: b=2
write(*,*) a,b
call add(a)
call add(b)
write(*,*) a,b
stop
end program main

subroutine add(number)
implicit none
integer :: number
number = number+1
return
end subroutine add

image-20220319134208618

  • 假设在一场田径赛的标枪选项中,有5位选手的投掷标枪的情况如下:

​ 1号选手:以30度角,每秒25米的速度掷出标枪

​ 2号选手:以45度角,每秒20米的速度掷出标枪

​ 3号选手:以35度角,每秒21米的速度掷出标枪

​ 4号选手,以50度角,每秒27米的速度掷出标枪

​ 5号选手,以40度角,每秒22米的速度掷出标枪

假如忽略空气阻力以及身高等等因素,请写一个程序来计算选手们的投射距离。(也就是计算自由投射运动的抛物线距离)

program main
    implicit none
    integer,parameter :: player = 5 !5位选手
    real :: angle(player) = (/30.0,45.0,35.0,50.0,40.0/) !掷出角度
    real :: speed(player) = (/25.0,20.0,21.0,27.0,22.0/) !掷出速度
    real :: distance(player) !抛物线距离
    integer :: I !循环量
    do I=1,player
        call getDistance(angle(I),speed(I),distance(I))  !调用求抛物线距离方法
        write(*,"(I1,'player distance =',F5.2)") I,distance(I)
    end do
    stop
end program main

subroutine getDistance(angle,speed,distance)
    real :: angle_pi,speed_y,speed_x,time
    real,parameter :: PI = 3.1415926,G = 9.81
    angle_pi = angle/180*PI
    speed_y = speed*sin(angle_pi)
    speed_x = speed*cos(angle_pi)
    time = 2.0*speed_y/G
    distance = time*speed_x
    return
end subroutine getDistance

image-20220319143558374

  • 将0~360的角度转化成0~2PI的弧度
subroutine angle_to_rad(angle,rad)
    implicit none
    real angle,rad
    real,parameter :: PI = 3.14159
    rad = angle/180*PI
    return
end subroutine angle_to_rad
  • 从上面的实例中,我们可以看到fortran 库函数的三角函数使用都是弧度,所以要将角度转化成弧度再计算
program main
    real :: number,rad,angle = 30
    real,parameter :: PI = 3.1415926
    rad = angle/180*PI
    number = sin(rad)
    write(*,*) number
    stop
end program main

image-20220319144805805

自定义函数(FUNCTION)

  • 自定义函数和子程序类似
    • 经过调用才能执行
    • 可以独立声明变量
    • 参数传递
  • 二者不同
    • 调用自定义函数前要先声明
    • 自定义函数执行后会返回一个数值
  • 声明时建议使用external,表明其是一个可调用的函数
program main
    implicit none
    real :: a=1
    real :: b=2
    real,external :: add !声明自定义函数,函数返回值为real类型
    write(*,*) add(a,b) !自定义函数不需要使用call
    stop
end program main

function add(a,b) result(c)
    implicit none
    real :: a,b,c
    c=a+b
    return
end function add

image-20220319171612862

  • 当自定义函数很简短时
program main
    implicit none
    real :: a=1,b
    real :: add
    add(a,b) = a+b !将自定义函数内容直接写到函数内
    write(*,*) add(a,3.0)
    stop
end program main

image-20220319172011219

  • 探讨一下函数与子程序的概念:使用函数有个“不成文规定”,“传递给函数参数,只要读取它的数值就好了,不要去改变它的数据”

全局变量(COMMON)

  • 在主程序、函数、子程序直接使用变量,就需要声明全局变量

COMMON的使用

 program main
     implicit none
     integer :: a,b
     common a,b !定义a,b是全局变量中的第一个和第二个变量
     a=1
     b=2
     call showCommon()
     stop
 end program main

 subroutine showCommon()
     implicit none
     integer :: num1,num2
     common num1,num2  !定义num1,num2是全局变量中的第一个和第二个变量
     write(*,*) num1,num2
     stop
 end subroutine showCommon

image-20220319172919277

  • 以上例子可以看到,全局变量与变量名无关,common定义了一块共享的内存空间。也就是取用全局变量是根据相对位置关系来作对应,而不是使用变量名称来对应。
  • 如果使用全局变量多的时候,需要将全局变量前面的变量先填充,这样就会造成麻烦。
  • 把变量归类、放在彼此独立的common区间中,这样有效地减少麻烦。
program main
    implicit none
    integer :: a,b
    common /group1/ a !将全局变量a放在group1中
    common /group2/ b !将全局变量b放在group2中
    a=1
    b=2
    call showGroup1()
    call showGroup2()
    stop
end program main

subroutine showGroup1()
    implicit none
    integer :: num1
    common /group1/ num1
    write(*,*) num1
    return
end subroutine showGroup1

subroutine showGroup2()
    implicit none
    integer :: num1
    common /group2/ num1
    write(*,*) num1
    return
end subroutine showGroup2

image-20220319183854552

  • 结论:共享变量不多,只有少数几个程序需要使用这些程序,就需要参数。如果需要共享大量数据时,很多个不同程序需要这些数据,就需要使用全局变量。

BLOCK DATA

  • common不能直接在子程序或主程序中使用Data设置初值,要在Block Data程序模块中使用data命令来设置初值。
program main
    implicit none
    integer :: a,b
    common a,b
    integer c,d
    common /group1/ c,d
    integer e,f
    common /group2/ e,f
    write(*,"(6I4)") a,b,c,d,e,f
    stop
end program main

block data
    implicit none
    integer :: a,b
    common a,b
    data a,b /1,2/
    integer c,d
    common /group1/ c,d
    data c,d /3,4/
    integer e,f
    common /group2/ e,f
    data e,f /5,6/
end block data

image-20220320085022901

  • 可以看到block data不需要调用就可以执行,这一段程序会在主程序执行前就会失效,只能设置全局变量初值,不能有其他执行命令出现。
  • 全局变量不能声明成常量。

注意事项

  • 使用common需要注意:
    • 变量的类型
    • 变量的位置
  • 基本上common只是用来提供一块共享的内存,编译器不会去注意程序员如何使用这一块内存,它不会帮忙做类型检查。
program main
    implicit none
    real :: a
    common a
    a = 1.0
    call showCommon()
    stop
end program main

subroutine showCommon()
    implicit none
    integer :: a
    common a
    write(*,*) a
    return
end subroutine showCommon

image-20220320091520966

  • 可以看到,得到的数据与预期相差甚远,这是因为不同类型造成的,浮点数1.0在使用的二进制存储的二进制数转化成整数时,就变成了1065353216

  • 数组在全局变量中的应用

program main
    implicit none
    real :: a,b
    common a,b
    a=1.0
    b=2.0
    call showCommon()
    stop
end program main

subroutine showCommon()
    implicit none
    real :: a(2)
    common a
    write(*,*) a(1),a(2)
    return
end subroutine showCommon

image-20220320092818272

函数中的变量

  • 介绍函数中跟变量相关的各种事项,包括参数传递的技巧、注意事项及参数的“生存周期”等等。

传递参数的注意事项

  • 传递参数的数据类型一定要正确,fortran使用的是地址传值
program main
    implicit none
    real :: a=1.0
    call showInteger(a)
    call showReal(a)
    stop
end program main

subroutine showInteger(num)
    implicit none
    integer :: num
    write(*,*) num
    return
end subroutine showInteger

subroutine showReal(num)
    implicit none
    real :: num
    write(*,*) num
    return
end subroutine showReal
  • 执行结果1065353216 1.000000 ,编译后报错,应该是现在编译器智能些了,提示类型不一致

数组参数

  • 数组会占据一块内存中的连续空间,在传递数组元素时,就是传递某一个内存地址
program main
    implicit none
    integer :: a(5) = (/1,2,3,4,5/)
    call showOne(a) !数组a的第一个元素的内存空间
    call showArray5(a)
    call showArray3(a)
    call showArray3(a(2))  !数组a的第二个元素的内存空间
    call showArray2X2(a)
    stop
end program main

subroutine showOne(num)
    implicit none
    integer :: num(1) !输出数组第一个数字
    write(*,*) num
    return
end subroutine showOne

subroutine showArray5(num)
    implicit none
    integer :: num(5) !输出数组
    write(*,*) num
    return
end subroutine showArray5

subroutine showArray3(num)
    implicit none
    integer :: num(3) !数组前三个数字
    write(*,*) num
    return
end subroutine showArray3

subroutine showArray2X2(num)
    implicit none
    integer :: num(2,2) !取出数组前四个数字,当作2X2的数组来使用
    write(*,*) num(1,1),num(2,1),num(1,2),num(2,2)
    return
end subroutine showArray2X2

image-20220320104317013

  • 数组在声明时要使用常量来赋值它的大小,不过在函数中,如果数组是接收用的参数时,可以例外,这时候可以使用变量来赋值它的大小,甚至不去赋值大小
program main
    implicit none
    integer,parameter :: size = 5
    integer :: s = size
    integer :: a(size) =(/1,2,3,4,5/)
    call useArray1(a,size)
    call useArray1(a,s)
    call useArray2(a)
    call useArray3(a)
    stop
end program main

subroutine useArray1(num,size)
    implicit none
    integer :: size
    integer :: num(size)
    write(*,*) num
    return
end subroutine useArray1

subroutine useArray2(num)
    implicit none
    integer :: num(*)
    integer :: i
    write(*,*) (num(i),i=1,5)
    return
end subroutine useArray2


subroutine useArray3(num)
    implicit none
    integer :: num(-2:2) !重新定义数组坐标
    write(*,*) num(0)
    return
end subroutine useArray3

image-20220320115637508

  • 上面实例可以得出,数组在主程序就已经配置了内存空间和数组大小,所以子程序可以不用再配置数组大小。使用数组时,注意不能超出实际范围,可以将一维数组变成二维数组,改变坐标范围。
  • 子程序在传递字符串变量时,可以不特别赋值它的长度,主程序已声明长度,与数组类似。
program main
    implicit none
    character(len=30) :: string = "Hello World,Hello!    "
    call showString(string) !字符串开头地址
    call showString(string(8:)) !字符串第八个字符地址
    stop
end program main

subroutine showString(str)
    implicit none
    character(len=*) :: str
    write(*,*) len_trim(str)
    write(*,*) str
    return
end subroutine showString

image-20220320152612490

  • 多维数组在传递时,只有最后一维可以不赋值大小,其他维都必须赋值大小。
program main
    implicit none
    integer,parameter :: dim1 = 2
    integer,parameter :: dim2 = 2
    integer,parameter :: dim3 = 2
    integer :: a(dim1,dim2,dim3)
    a(:,:,1) = 1
    a(:,:,2) = 2
    call GetArray1(a(:,:,1),dim1,dim2)
    call GetArray2(a(1,1,2),dim1)
    stop
end program main

subroutine GetArray1(a,dim1,dim2)
    implicit none
    integer :: dim1,dim2
    integer :: a(dim1,dim2)
    write(*,*) a
    return
end subroutine GetArray1

subroutine GetArray2(a,dim1)
    implicit none
    integer :: dim1,dim2
    integer :: a(dim1,*)
    integer :: i
    write(*,*) (a(:,i),i=1,2)
    return
end subroutine GetArray2

image-20220320161511660

变量的生存周期

  • 子程序运行完,变量就会释放,删除,save 可以拯救这些变量、增加变量的生存周期、保留住所保存的数据。这些变量可以在程序执行中永久记忆住上一次函数调用时所设置的数值。
program main
    implicit none
    call sub()
    call sub()
    call sub()
    stop
end program main

subroutine sub()
    implicit none
    integer :: count = 1 !变量初值只会设置一次
    save count
    write(*,*) count
    count = count+1
    return
end subroutine sub

image-20220320165951179

  • Fortran 标准中没有强制规定,没有使用save的变量就不会被永远记住,有些编译器会永远记住变量的数值,使用save可以程序的正确性和可移植

传递函数

  • 传递参数时,可以传递数字、字符等等,还可以传递函数名称。
program main
    implicit none
    real,external :: func !声明func是一个自定义函数
    real,intrinsic :: sin ! 声明sin是个库函数
    call execFunc(func) !计算func(1.0)
    call execFunc(sin)  !计算sin(1.0)
    stop
end program main

subroutine execFunc(f)
    implicit none
    real,external :: f
    write(*,*) f(1.0)
    return
end subroutine execFunc

real function func(num)
    implicit none
    real :: num
    func = num*2
    return
end function func

image-20220320172206015

  • 除了传递函数名称,还可以传递子程序
program main
    implicit none
    external sub1,sub2 !声明sub1、sub2是子程序名称
    call sub(sub1)
    call sub(sub2)
    stop
end program main

subroutine sub(sub_name)
    implicit none
    external sub_name !声明sub_name是个子程序
    call sub_name()
    return
end subroutine sub

subroutine sub1()
    implicit none
    write(*,*) "sub1"
    return
end subroutine sub1

subroutine sub2()
    implicit none
    write(*,*) "sub2"
    return
end subroutine sub2

image-20220320173241895

  • 与C语言中函数指针类似

特殊参数的使用方法

  • fortran 90,可以赋值参数的属性,设置某些参数是只读不能改变的,它还可以输入不定个数的参数,还可以不按照顺序来传递参数

参数属性

  • Fortran 90的intent命令可以用来设置参数属性
program main
    implicit none
    integer :: a=4
    integer :: b
    call div(a,b)
    write(*,*) a,b
    stop
end program main

subroutine div(a,b)
    implicit none
    integer,intent(in) :: a  !指定a是只读
    integer,intent(out) :: b  !指定b在子程序中应该重新设置数值
    b=a/2
    return
end subroutine div

image-20220320184658229

  • 设置成只读的变量,在函数中如果重新设置数值,编译过程中会出现错误。而设置成要输出的变量,如果函数中没有重新设置一个数值给它,编译过程会出现警告信息
  • intent(inout)这个属性指参数可读又可写,基本上就和什么都没有指定时是一样的

函数的使用接口(INTERFACE)

  • interface 是一段程序模块,用来清楚说明所要调用函数的参数类型及返回值类型等等的“使用接口”。在一般情况下,使用函数时不需要特别说明它们的“使用接口”,不过在下面这些情况下是必要的。
    • 1、函数返回值为数组时
    • 2、指定参数位置来传递参数时
    • 3、所调用的函数参数数目不固定时
    • 4、输入指标参数时
    • 5、函数返回值为指针时
  • 函数返回值为数组时
program main
    implicit none
    interface !定义函数function random10的使用接口
        function random10(lbound,ubound)  !这里interface接口包括参数类型及返回值类型
            implicit none
            real :: lbound,ubound
            real :: random10(10) ! 返回值为数组
        end function random10
    end interface
    real :: a(10)
    call random_seed() !库存子程序,使用随机数前调用
    a=random10(1.0,10.0) !生成10个1.0~10.0之间的随机数
    write(*,"(10F6.2)") a
    stop
end program main

function random10(lbound,ubound)
    implicit none
    real :: lbound,ubound
    real :: len
    real :: random10(10)
    real :: t
    integer :: i
    len = ubound-lbound
    do i=1,10
        call random_number(t)  !t是0~1之间的随机数
        random10(i)=lbound+len*t
    end do
    return
end function random10

image-20220320195353006

  • 使用随机数前,要先调用random_seed这个子程序来启动随机数生成器
  • 当多个函数需要使用interface来声明接口时,就需要给每个函数进行声明,比较麻烦,下一章,将使用module解决问题

不定个数的参数传递

  • 在fortran 77 及很多的程序语言中,函数的参数个数都是有固定数目的。在fortran 90中,可以使用optional命令来表示某些参数是”可以省略的”
program main
    implicit none
    interface
        subroutine sub(a,b)
            implicit none
            integer :: a
            integer,optional :: b
        end subroutine sub
    end interface
    !开始执行程序
    call sub(1)
    call sub(2,3)
    stop
end program main

subroutine sub(a,b)
    implicit none
    integer :: a
    integer,optional :: b
    if(present(b)) then !判断是否有b输入,函数present有参数传递进来时,返回.true.,否则反之
        write(*,"('a=',I3,'b=',I3)") a,b
    else
        write(*,"('a=',I3,'b=unknown')") a
    end if
    return
end subroutine sub

image-20220320201438670

  • optinonal命令可以将只改变参数个数,不改变内容的函数,写在一起
  • 调用这一类不定数目参数的函数时,一定要先声明出函数的interface,使用module时则可以省略

改变参数传递位置的方法

  • fortran 90中,可以不用按照参数的顺序来传递参数。

  • 语法

subroutine sub(a,b,c)
implicit none
··········
return
end subroutine sub
! 调用时
call sub(b=1,c=3,a=2) !根据变量名来传递参数
  • 在不定个数的参数传递时,改变参数传递位置就显得尤为重要
interface
  subroutine sub(a,b,c,d,e,f)
  implicit none
  interger,optional :: a,b,c,d,e,f
  end subroutine sub
end interface

call sub(f=3,e=5)
  • 使用这个方法来传递参数时,一定要声明interface,这个方法可以省略许多不必要的描述,尤其在子程序可以接收多个参数时。
  • 编写程序时,可以让某些参数有默认值,不输入这些参数时就使用默认值,编写程序计算F(X)=A*X^2+B*X+C,一定输入X,A、B、C没有输入的话默认值为0
program main
    implicit none
    interface
        real function func(x,a,b,c) !定义子程序function的使用接口
            implicit none
            real x
            real,optional :: a,b,c
        end function func
    end interface
    !编写执行命令
    write(*,*) func(2.0,c=1.0)
    write(*,*) func(2.0,a=2.0,b=1.0)
    stop
end program main

real function func(x,a,b,c)
    implicit none
    real x
    real,optional :: a,b,c  !a,b,c可以不输入
    real ra,rb,rc ! 实际计算的变量
    if(present(a)) then
        ra=a
    else
        ra=0.0
    end if
    if(present(b)) then
        rb=b
    else
        rb=0.0
    end if
    if(present(c)) then
        rc=c
    else
        rc=0.0
    end if
    func = ra*x**2+rb*x+rc
    return
end function func

image-20220320204832327

特殊的函数类型

  • fortran 90的函数,除了一般正常使用的类型外,还可以特别指定成recursivepureelemental这三种类型之一。recursive是让函数可以自己调用自己,也就是递归。pureelemental是用来做并行处理时及设置数组时使用

递归

  • 递归函数每次被调用执行时,函数中所声明的局部变量(指那些不是传递进来的参数,及没有save的变量)都会使用不同的内存地址。简言之,函数中的局部变量在每次调用时都是独立存在的。
  • 使用递归计算阶乘
program main
    implicit none
    integer :: n
    integer,external :: fact
    write(*,*) 'N='
    read(*,*) n
    write(*,"(I2,'!=',I8)") n,fact(n)
    stop
end program main

recursive integer function fact(n) result(ans)  !recursive让函数自己调用自己
    implicit none
    integer,intent(in) :: n
    if(n<0) then
        ans =-1
        return
    else if(n<=1) then
        ans = 1
        return
    end if
    ans = n*fact(n-1)
    return
end   function fact

image-20220321085906640

内部函数

  • fortran 90中还可以把函数做一个“归属”,定义出某些函数只能在某些特定的函数中被调用,语法如下
program main 或 subroutine sub 或 function func
    ··············
    ··············
    contains    !contains后面开始写局部函数
        subroutine localsub  !localsub只能在包含它的函数中被调用
            ···············
            ···············
        end subroutine localsub

        function localfunc  !localfunc只能在包含它的函数中被调用
            ···············
            ···············
        end function localfunc
end program/subroutine/funciton
  • 这个方法,可以用来设计一个函数中的“内部运行”,因为内部运行是不希望被别人所使用的
program main
    implicit none
    integer :: n
    write(*,*) 'N='
    read(*,*) n
    write(*,"(I2,'!=',I8)") n,fact(n)
    stop

contains
    recursive integer function fact(n) result(ans)
        implicit none
        integer,intent(in) :: n
        if(n<0) then
            ans=-1
            return
        else if(n<=1) then
            ans=1
            return
        end if
        ans=n*fact(n-1)
        return
    end function fact
end program main

image-20220321093414553

  • 把函数fact放在主程序中,这样函数fact只能在主程序中被调用,其他函数不能调用它。主程序不需要声明就可以直接调用到函数fact

PURE函数

  • function/subroutine前面加上pure就可以使用pure函数。一般情况下不需要使用pure函数,它是用来配合并行运算使用。

  • pure函数的限制

    • 1、pure函数的参数必须都是只读intent(in)
    • 2、pure子程序的每个参数都要赋值属性
    • 3、pure函数中不能使用save
    • 4、pure函数中所包含的内部函数也都必须是pure类型函数
    • 5、pure函数中不能使用stop,print及跟输出入相关的命令如raed,write,open,close,backspace,endfile,rewind,inquire等等
    • 6、pure函数中只能够读取,不能改变全局变量的值
  • 实例

program main
    implicit none
    integer,external :: func
    write(*,*) func(2,3)
    stop
end program main

pure integer function func(a,b)
    implicit none
    integer,intent(in) :: a,b
    func = a+b
    return
end function func

image-20220321181440955

ELEMENTAL函数

  • function/subroutine前面加上elemental就可以使用elemental函数。这个类型的函数也可以用来做并行运算,同样与pure有6项限制。
  • elemental可以做数组的设置,也就多了一个限制:参数不能是数组
integer a(10)
a = func(a) !func 是elemental函数

!  与下面等效
do i=1,10
a(i)=func(a(i))
end do
  • elemental函数主要就是用来配合fortran 90可以对整个数组操作的语法来设置数组内容
program main
    implicit none

    interface
        elemental real function func(num)
            implicit none
            real,intent(in) :: num
        end function
    end interface

    integer i
    real :: a(10) = (/(i,i=1,10)/)
    real :: b(10)
    write(*,"(10F6.2)") a
    a = func(a)
    write(*,"(10F6.2)") a
    stop
end program main

elemental real function func(num)
    implicit none
    real,intent(in) :: num
    func =  sin(num) + cos(num)
    return
end function func

image-20220321183424671

  • 使用elemental函数时要说明它的使用接口,主程序才会正确地做出设置数组的调用

MODULE

  • module可以用来封装程序模块,通常是用来把程序中,具备相关功能的函数及变量封装在一起

MODULE中的变量

  • module语法
module module_name
    ··········
    ··········
    ··········
end module module_name
  • 举例来说,需要使用全局变量时,可以把全局变量都声明在module中,需要使用这些变量的函数只要use这个module就可以使用它们
module global
    implicit none
    integer a,b
    common a,b
end module

program main
    use global
    implicit none
    a=1
    b=2
    call sub()
    stop
end program main

subroutine sub()
    use global
    implicit none
    write(*,*) a,b
    return
end subroutine sub

image-20220321184300070

  • 一般情况下,应该最先写module,这样主程序与子程序都可以使用它
module global
    implicit none
    integer,save :: a !使用save,相当于全局变量
end module global

program main
    use global
    implicit none
    a =10
    call sub()
    write(*,*) a
    stop
end program main

subroutine sub()
    use global
    implicit none
    write(*,*) a
    a=20
    return
end subroutine sub

image-20220321185002387

MODULE中的自定义类型TYPE

  • 需要传递多个变量时,我们总是会弄混,这时候就需要自定义类型,把这些数值都封装在这个新的类型中,传递参数时只要传递一个变量过去就行了。
module constant
    implicit none
    real,parameter :: PI = 3.14159
    real,parameter :: G = 9.81
end module constant

module typedef
    implicit none
    type player  !自定义类型player
        real :: angle
        real :: speed
        real :: distance
    end type player
end module typedef

program main
    use typedef
    implicit none
    integer,parameter :: players = 5
    type(player) :: people(players) = (/ player(30.0,25.0,0.0),player(45.0,20.0,0.0),player(35.0,21.0,0.0),&
            &player(50.0,27.0,0.0),player(40.0,22.0,0.0) /)
    integer :: I

    do I=1,players
        call get_distance(people(I))
        write(*,"('player',I1,'=',F8.2)") I,people(I)%distance
    end do
    stop
end program main

! 把0~360的角度转换成0~2PI的弧度
real function angle_to_rad(angle)
    use constant
    implicit none
    real angle
    angle_to_rad = angle*PI/180.0
    return
end function angle_to_rad

subroutine get_distance(person)
    use constant
    use typedef
    implicit none
    type(player) :: person
    real :: rad,Vx,time
    real,external :: angle_to_rad !声明是个函数
    rad = angle_to_rad(person%angle)
    Vx = person%speed*cos(rad)
    time = 2.0*person%speed*sin(rad)/G
    person%distance = Vx*time
    return
end subroutine get_distance

image-20220322091649796

MODULE中的函数

  • module中还可以容纳函数,编写结构如下:
module module_name
   ·········
   ·········
contains !从contains后开始写作函数
   subroutine sub_name
        ·········
        ·········
   end subroutine sub_name

   function func_name
        ······
        ······
   end function
end module module_name
  • 通常会把功能上相关的函数放在同一个module模块中,程序想要调用module中的函数时,也要先通过use module_name的命令,才能调用它们。这种做法比较符合模块化概念,编写程序时,可以把程序中属于绘图功能的部分放在module Graphics中,把数值计算部分放在module Numerical中。
  • 函数可以直接使用同一个module里所声明的变量
module tool
   implicit none
   integer :: a !module中声明了变量a
   ·······
   ·······
contains
   subroutine add()
   implicit none
   a=a+1
   ·······
   ·······
   end subroutine add
  • 把函数编写在module
module constant
    implicit none
    real,parameter :: PI = 3.14159
    real,parameter :: G = 9.81
end module

module typedef
    implicit none
    type player
        real :: angle
        real :: speed
        real :: distance
    end type
end module

module shoot
    use constant
    use typedef
    implicit none
contains
    !函数
    real function angle_to_rad(angle)
        implicit none
        real :: angle
        angle_to_rad = angle*PI/180.0
        return
    end function angle_to_rad

    subroutine get_distance(person)
        implicit none
        type(player) :: person
        real rad,Vx,time
        rad = angle_to_rad(person%angle) !封装在一起的函数互相认识,不需要声明即可使用
        Vx = person%speed*cos(rad)
        time = 2.0*person%speed*sin(rad)/G
        person%distance = Vx*time
        return
    end subroutine


end module


program main
    use shoot
    implicit none
    integer,parameter :: players = 5
    type(player) :: person(players)=(/player(30.0,25.0,0.0),player(45.0,20.0,0.0),player(35.0,21.0,0.0),&
            player(50.0,27.0,0.0),player(40.0,22.0,0.0)/)
    integer :: I
    do I=1,players
        call get_distance(person(I))
        write(*,"('player',I1,'=',F8.2)") I,person(I)%distance
    end do
    stop
end program main

image-20220322091637619

一些少用的功能

ENTRY

  • entry用来在函数中创建一个新的“入口”,调用这个入口时,会跳过进入点之前的程序代码来执行函数。
program main
    implicit none
    call sub()
    call mid()
    stop
end program main

subroutine sub()
    implicit none
    write(*,*) "Hello."
    entry mid()  !另一个进入点mid
    write(*,*) "Good morning."
    return
end subroutine sub

image-20220322092524429

特别的RETURN

  • 从函数中return返回调用处时,通常会直接返回调用处来继续执行程序。关于这一点也是可以改变的,调用函数时可以额外指定其他的折返点
program main
    implicit none
    integer :: a
    read(*,*) a
    call sub(a,*100,*200)  !特别另外指定两个折返点,分别是行代码100及200这两个地方
    write(*,*) "default"
    stop
    100 write(*,*) "return 1"
    stop
    200 write(*,*) "return 2"
    stop
end program main


subroutine sub(a,*,*)
    implicit none
    integer :: a
    if(a<=0) then
        return  !默认折返点
    else if(a==1) then
        return 1  !折返点1
    else
        return 2  !折返点2
    end if
end subroutine sub

image-20220322093522150

使用多个文件

  • 程序员通常会把一些具有相关功能的函数,独立编写在不同的文件中,编译器可以分别编译这些程序文件,最后再把它们链接到同一个执行文件中。把程序代码分散在不同文件中有几个好处:
    • 1、独立文件中的函数,可以再拿给其他程序使用。
    • 2、可以加快编译速度,修改其中一个文件时,编译器只需要重新编译这个文件就行了。

INCLUDE

  • include命令用来在程序代码中,插入另一个文件中的内容。
  • 使用include,代码文件需放在同一个目录下
program main
    implicit none
    call sub()
    stop
end program main

include 'text2.f90'
  • text2.f90 中的内容
subroutine sub()
    implicit none
    write(*,*) "Hello World!"
    return
end subroutine sub

image-20220322095101565

  • include命令可以写在任何地方,它只是单纯得用来插入一个文件的内容。有时候也会应用在声明全局变量,先把声明全局变量的程序代码编写在某个文件中,需要使用全局变量的函数再去include这个文件,这样做可以减少程序代码。
  • text1.inc中的文件
integer a,b
common a,b
  • 主程序
program main
    implicit none
    include 'text1.inc'
    a=1
    b=2
    call sub()
    stop
end program main

subroutine sub()
    implicit none
    include 'text1.inc'
    write(*,*) a,b
    return
end subroutine sub

image-20220322100201339

  • 在fortran 90提供module后,就不太需要使用include命令来做这种工作

Visual Fortran中使用多个文件

  • 使用include对编译器来说,并不是将程序代码分散在不同的文件里,实际上只是把程序代码手动分开,编译时又会合在一起,所以并不会加快编译速度
  • Visual Fortran打开一个项目project,其中项目中可以放入大量f90文件,编译成一个程序,当修改其中一个文件时,只需要重新编译修改好的文件
  • 分成多个文件可以将旧有的fortran 77 程序代码与fortran 90程序代码混用。

G77中使用多个文件

  • g77 text1.for text2.for -o text,指定输出的执行文件名字text
    • 第一步:g77 -c text1.for——生成text1.o
    • 第二步:g77 -c text2.for——生成text2.o
    • 第三步:g77 text1.o text2.o -o text

程序库

  • 具有特殊功能的一组函数,可以编译成*.LIB程序库来给其他人使用。*.LIB的文件内容经过编译,无法从文件中读取初始程序代码
  • 数值运算的IMSL,让程序具备3D绘图能力的OpenGL程序库

函数的应用

  • 在文本格式下绘图的程序库TextGraphLib,学习如何封装及使用程序库
  • 基本命令
    • subroutine SetScreen(width,height):在进行绘图之前,一定要先调用这个函数来决定要使用多大的画面。用户可以任意决定长宽范围,在windowd命令列窗口下,一般所能使用的最大范围为80*25
      • integer,intent(in) :: width:定义绘图的画面范围宽度
      • integer,intent(in) :: height: 定义绘图的画面范围高度
    • subroutine PutChar(x,y,char): 在 (x,y)位置画出字符char。char没有输入时,则使用默认的填充字符
      • integer,intent(in) :: x,y:坐标位置
      • character,optional,intent(in) :: char:要画出的字符
    • subroutine DrawLine(X0,Y0,X1,Y1):在(X0,Y0)到(X1,Y1)两点间画一条直线
      • integer,intent(in) :: X0,Y0:第一个二维坐标
      • integer,intent(in) :: X1,Y1:第二个二维坐标
    • subroutine DrawCircle(cx,cy,r1,r2):以(cx,cy)为圆点,r1为水平轴上半径长,r2为垂直轴半径长来画一个椭圆。r2值不输入时,r2=r1
      • integer,intent(in) :: cx,cy :圆心坐标位置
      • integer,intent(in) :: r1:椭圆水平轴半径长
      • integer,optional,intent(in) :: r2:椭圆垂直轴半径长,没有输入时r2=r1
    • subroutine DrawRect(x0,y0,x1,y1):以(x0,y0)为左上角,(x1,y1)为右下角来画一个空心的矩形
      • integer,intent(in) :: x0,y0:左上角坐标
      • integer,intent(in) :: x1,y1:右下角坐标
    • subroutine DrawFilledRect(x0,y0,x1,y1):以(x0,y0)为左上角,(x1,y1) 为右小角来画一个内部填满的矩形
      • integer,intent(in) :: x0,y0 :左上角坐标
      • integer,intent(in) :: x1,y1:右下角坐标
    • subroutine ClearScreen(char):把整个画面用char字符来填满,char没有输入时,使用默认的消除字符
      • character,optional,intent(in) :: char:用来消除画面的字符
    • subroutine SetCurrentChar(char):调用这个函数,可以决定要使用哪一个字符来填充画面
      • character,intent(in) :: char:用来改变预定的填充字符
    • subroutine SetBackground(char):调用这个函数,可以决定要使用哪一个字符来清除画面
      • character,intent(in) :: char:设置清除昼面时所使用的字符
    • subroutine UpdateScreen():所有的绘图操作,都是事先画在一块内存中,要调用UpdateScreen才会把绘图的结果真正呈现出来
  • 画一条线
program main
use TextGraphLib
implicit none
call SetScreen(10,10)
call DrawLine(1,1,10,10)
call UpdateScreen()
stop
end program main

image-20220330104623286

习题

  • 1、请写一个子程序,它可以计算圆面积。需要输入两个参数,第一个参数用来输入圆的半径长,第二个参数用来返回圆的面积
program main
    implicit none
    real :: r,s
    call Circular_area(r,s)
    stop
end program main

subroutine Circular_area(r,s)
    implicit none
    real,parameter :: PI = 3.14159
    real :: r,s
    write(*,*) "请输入圆的半径:"
    read(*,*) r
    s = PI*r**2
    write(*,*) "半径为",r,"的圆面积为:",s
    return
end subroutine Circular_area

image-20220322195424423

  • 2、把上一题改用函数来编写,这个时候就只需要通过输入参数来输入半径就够了
program main
    implicit none
    real :: r
    real,external :: Circular_area
    write(*,*) "请输入圆的半径:"
    read(*,*) r
    write(*,*) "半径为",r,"的圆面积为:",Circular_area(r)
    stop
end program main

real function Circular_area(r) result (s)
    implicit none
    real,parameter :: PI = 3.1415
    real :: r
    s = PI*r**2
    return
end function Circular_area

image-20220322200258524

  • 3、略
  • 4、试着用递归函数来计算等差数列1+2+3+4+···+100的值
program main
    implicit none
    integer,external :: func
    write(*,*) func()
    stop
end program main

recursive integer function func() result (sum)
    implicit none
    integer :: i = 0
    sum = 0
    i = i+1
    if(i>100) return
    sum = func() + i
    return
end function func

image-20220322201624997

  • 5、试写一个函数计算所输入的两个整数的最大公因子
program main
    implicit none
    integer :: a,b
    integer,external :: CommonFactorMax
    write(*,*) "请输入a:"
    read(*,*) a
    write(*,*) "请输入b:"
    read(*,*) b
    write(*,*) a,b,"的最大公因子为:",CommonFactorMax(a,b)
    stop
end program main

integer function CommonFactorMax(a,b) result (ans)
    implicit none
    integer :: a,b
    integer :: n
    ans = 1
    do n =2,10000
    if((mod(a,n)==0).and.(mod(b,n)==0)) then
        ans = n
    end if
    end do
    return
end function CommonFactorMax

image-20220322203701768

  • 6、略

文件

文件读取的概念

  • 计算机两项重要功能:1、计算、处理数据。2、保存数据
  • fortran 90读取文件的操作有:“顺序读取”和“直接读取”
    • 顺序读取:读写一个文件,只能从开始读取,一步一步向下读取。
    • 直接读取:读写文件时,可以任意跳跃到文件的任何一个位置来读写。
  • 保存文件,有两种格式:“文本文件”和“二进制文件”

文件的操作

open的使用

  • 语法

open(unit=10,file='hello.txt'):unit=10给后面的文件指定一个代码,file='hello.txt'指定所要开启的文件名称

program main
  implicit none
  open(unit=10,file='hello.txt')
  write(10,*) "hello"
  stop
end program main
  • 详论open的参数
open(unit=number,file='···',form='···',status='···',access='···',recl=length,err=label,iostat=iostat,blank='···',position='···',action=action,pad='···',delim='···')
  • number必须是正整数,它可以使用变量或是常量来赋值。number值最好避开1,2,5,6.因为2,6是默认的输出位置,也就是屏幕。1,5则是默认的输入位置,也就是键盘。
  • file='filename':这个字段用来指定所要打开的文件名称,文件名要符合操作系统规定。像是windows下文件名不区分大小写,unix中则会区分大小写,不管使用哪一种操作系统,最好都不要使用中文文件名。
  • form='formatted or unformatted':formatted表示文件使用“文本文件”格式保存,unformatted表示文件使用“二进制文件”格式来保存,默认文本格式
  • status='new' or 'old' or 'replace' or 'scratch' or 'unknown':
    • new表示这个文件原来不存在,是第一次打开
    • old表示这个文件原本就已经存在
    • replace文件若已经存在,会重新创建一次,原本的内容会消失,文件若不存在,会创建新文件
    • scratch表示打开一个暂存盘,这个时候可以不需要指定文件名称,也就是File这一栏可以忽略,因为程序本身会自动取一个文件名,至于文件名是什么也不重要,因为暂存盘会在程序结束后自动删除。
    • unknown由各编译器自定义,通常会同replace的效果,默认设置
  • access='sequential' or 'direct':sequential读写文件的操作会以“顺序”的方法来做读写,也就是“顺序读取文件”,默认设置。direct读写文件的操作可以任意指定位置,这就是“直接读取文件”。
  • recl=length:在顺序读取文件中,recl字段值用来设置一次可以读写多大容量的数据。在打开“直接读取文件”时,recl=length的length值是用来设置文件中每一个模块单元的分区长度。
  • err=label :这个字段用来设置当文件打开发生错误时,程序会跳跃到label所指的行代码处来继续执行程序
  • iostat=var:这个字段会设置一个整数值给后面的整型变量,这是用来说明文件打开的状态,数值有三种情况:1、var>0 读取操作发生错误。2、var=0 读取操作正常。 3、var<0 文件终了
  • blank='null' or 'zero':这用来设置文件输入数字时,当所设置的格式字段中有空格存在时所代表的意义,null,空格代表没有东西,zero,空格部分会自动以0填充
  • Fortran 90添加功能
  • position = 'asis' or 'rewind' or 'append':设置文件打开时的读写位置
    • asis表示文件打开时读写位置,不特别指定,通常就是在文件开头,默认值
    • rewind表示文件打开时的读取位置移到文件的开头
    • append表示文件打开时的读取位置移到文件的结尾
  • action = 'read' or 'write' or 'readwrite':设置所打开文件的读写权限
    • readwrite表示所打开的文件可以用来读取及写入,这是默认值
    • read表示所打开的文件只能用来读取数据
    • write 表示所打开的文件只能用来写入数据
  • pad = 'yes' or 'no'
    • yes在格式化输入时,最前面的不足字段会自动以空格填满,默认值
    • no在格式化输入时,不足的字段不会自动以空格填满
  • delim = 'apostrophe' or 'quote' or 'none'
    • none纯粹输出字符串内容
    • quote输出字符串会在前后加上双引号
    • apostrophe输出字符串会在前后加上单引号

write,read的使用

program main
    implicit none
    character(len=20) :: string
    open(unit=10,file="test.txt")
    write(10,"(A20)") "Good Morning." !写入数据
    rewind(10) !读取位置移到文件的开头
    read(10,"(A20)") string  !读取数据
    write(*,"(A20)") string  !输出数据
    stop
end program main

image-20220323193255671

  • 详论write/read的参数
write/read(unit=number,fmt=format,nml=namelist,rec=record,iostat=stat,err=errlabel,end=endlabel,advance=advance,size=size)
  • unit=number指定write/read所使用的写入/出的位置
  • fml=format指定输入输出格式使用
  • nml=namelist指定读写某个namelist的内容
  • rec=record在直接读取文件中,设置所要读写的文件模块位置
  • iostat=var:这个字段会设置一个整数值给后面的整型变量,这是用来说明文件打开的状态,数值有三种情况:1、var>0 读取操作发生错误。2、var=0 读取操作正常。 3、var<0 文件终了
  • err=errlabel指定在读写过程中发生错误时,会转移到某个行代码来继续执行程序
  • end=endlabel指定在读写到文件末尾时,要转移到某个行代码来继续执行程序
  • Fortran 90添加功能
  • advance = 'yes' or 'no':yes每读写一次会向下移动一行,默认值。no 暂停自动换行的操作。使用这个字段时,一定要设置输出入格式。在屏幕输出时可以使用这个设置来控制write命令是否会自动换行
  • size=countadvance = no时,才可以使用这个字段,它会把这一次输出入的字符数目设置给后面的整型变量

查询文件的状态inquire

  • 编写一个检查某个文件是否存在的程序
program main
    implicit none
    character(len=20) :: filename = "text001.f90" !需要查询的文件名
    logical alive
    inquire(file=filename,exist=alive)
    if(alive) then
        write(*,*) filename,"exist."
    else
        write(*,*) filename,"doesn't exist."
    end if
    stop
end program main
  • 详论inquire的参数
inquire(unit=number,file=filename,iostat=stat,err=label,exist=exist,opened=opened,number=number,named=named,access=access,sequential=sequential,direct=direct,form=form,formatted=formatted,unformatted=unformatted,recl=recl)
  • 基本与open类似
  • exist=sxist检查文件是否存在,会返回一个布尔变量给后面的逻辑变量,返回真值表示文件存在,返回假值表示文件不存在
  • opened=opened检查文件是否已经使用open命令来打开,会返回一个布尔变量给后面的逻辑变量,返回真值表示文件已打开,返回假值表示文件尚未打开
  • number=number由文件名来查询这个文件所给定的代码
  • named=named查询文件是否取了名字,也就是检查文件是否为临时保存盘,返回值为逻辑数
  • sequential=sequential查看文件是否使用顺序格式,会返回一个字符串,字符串可以是yes,no,unknown
  • direct=direct查看文件是否使用直接格式,会返回一个字符串,字符串值可以是yes,no,nuknown
  • form=form查看文件的保存方法,返回一个字符串,formatted,unformatted,undefined未定义
  • formatted=fmt查看文件是否为文本文件,会返回一个字符串
  • unformatted=fmt查看文件是否为二进制文件,会返回一个字符串
  • recl=length返回open文件时recl栏的设置值,nextrec=nr返回下一次文件读写的位置,blank=blank返回值是字符串,用来查看open文件时blank参数所给定的字符串值
  • Fortran 90添加功能
  • position = position返回打开文件时position字段所给定的字符串,字符串值可能是rewind,append,asis,undefined
  • action=action返回打开文件时action字段所赋值的字符串,字符串值可以为read,write,readwrite
  • read=read返回一个字符串,检查文件是否为只读文件
  • write=write返回一个字符串,检查文件是否可以写入
  • readwrite=readwrite返回一个字符串,检查文件是否可以同时读及写
  • delim=delim返回打开文件时,delim字段所设置的字符串,返回值可以是apostrophe,quote,none,undefined
  • pad=pad返回打开文件时pad字段所设置的字符串,返回yesno

其他文件运行命令

  • backspace(unit=number,err=errlabel,iostat=iostat):把文件的读写位置退回一步,其他字段参考上面
  • endfile(unit=number,err=errlabel,iostat=iostat):使用这个命令会把目前文件的读写位置变成文件的结尾,其他字段参考上面
  • rewind(unit=number,err=errlabel,iostat=iostat):把文件的读写位置倒回文件开头
  • close(unit=number,status=string,err=errlabel,iostat=iostat):把文件关闭,不再进行读写操作,status='keep'会在文件关闭后,保留住这个文件,默认状态。status='delete'会在文件关闭后,消除这个文件
  • 编写一个删除文件的程序
program main
    implicit none
    integer,parameter :: fileid=10
    logical alive
    character(len=20) :: filename;
    write(*,*) "filename:"
    read(*,"(A20)") filename
    inquire(file=filename,exist=alive)
    if(alive) then
        open(unit=fileid,file=filename)
        close(fileid,status="delete")
    else
        write(*,*) trim(filename)," doesn't exist."
    end if
    stop
end program main

顺序文件的操作

  • 命令窗口中type可以用来在屏幕上快速地浏览一个文本文件内容,下面实例实现同样的功能
program main
    implicit none
    character(len=79) :: filename
    character(len=79) :: buffer
    integer,parameter :: fileid = 10
    integer :: status = 0
    logical alive
    write(*,*) "filename:"
    read(*,"(A79)") filename
    inquire(file=filename,exist=alive)
    if(alive) then
        open(unit=fileid,file=filename,access = "sequential",status = "old")
        do while(.true.)
            read(unit=fileid,fmt="(A79)",iostat=status) buffer
            if(status/=0) exit !没有数据就跳出循环
            write(*,"(A79)") buffer
        end do
    else
        write(*,*) trim(filename),"doesn't exist."
    end if
    stop
end program main

image-20220324090026579

  • 编写一个写入程序,可以用来记录全班同学的中文、英文及数学的考试成绩
module typedef
    type student
        integer chinese,english,math
    end type
end module

program main
    use typedef
    implicit none
    integer :: students
    type(student),allocatable :: s(:) !声明一个可变大小的一维数组
    character(len=80) :: filename='data.txt'
    integer,parameter :: fileid = 10
    integer :: i
    write(*,*) "班上有多少学生?"
    read(*,*) students
    allocate(s(students),stat=i) !配置内存空间,并判断是否配置成功,i=0成功
    if(i/=0) then
        write(*,*) "Allocate buffer fail."
        stop
    end if
    open(fileid,file=filename)
    do i=1,students
        write(*,"('请输入'I2'号同学的中文、英文及数学成绩')") i
        read(*,*) s(i)%chinese,s(i)%english,s(i)%math
        write(fileid,"('座号:'I2/'中文:'I3'英文:'I3'数学:'I3)") i,s(i) !自定义类型student,包含三个变量,输出可以直接写s(i)
    end do
    close(fileid)
    stop
end program main

image-20220324092958640

  • 读取上一个程序创建的data.txt文件内容

    module typedef
       type student
           integer :: chinese,english,math
       end type
    end module
    
    program main
       use typedef
       implicit none
       type(student) :: s
       character(len=80) :: filename = 'data.txt'
       integer,parameter :: fileid = 10
       logical alive
       integer :: error
       integer :: no
       inquire(file=filename,exist=alive)
       if(.not. alive) then
           write(*,*) trim(filename),"doesn't exist."
           stop
       end if
    
       open(fileid,file=filename)
       do while(.true.)
           read(fileid,"(10XI2,/,10XI3,10XI3,10XI3)",iostat=error) no,s
           if(error/=0) exit
           write(*,"(I2'号 中文:'I3 ' 英文:'I3 ' 数学:' I3 )") no,s
       end do
       close(fileid)
       stop
    end program main

image-20220324095313328

  • 读取一个版面格式比较自由的文件,数据之间可以有任意数量的空行,文件内容如下

  • 数据

姓名  王天才
体重  80.5
身高  195.2
得分  15.8

姓名  李天才
身高  190.3
体重  85.1
得分  10.8




姓名  洪天才
体重  90.8
身高  201.3
得分  19.8

姓名  彭天才
体重  70.2
得分  22.2
身高  185.0

姓名  黄天才
得分  20.1
体重  85.0
身高  190.3
  • 编写一个程序来读入所有选手的数据,并挑选出平均每场可以得20分以上的球员来显示他们的数据
module typedef
    type player
        character(len =80) :: name
        real :: weight,height
        real score
    end type
end module

program main
    use typedef
    implicit none
    character(len=20) :: filename = 'players.txt'
    integer,parameter :: fileid = 20
    logical :: alive
    type(player) :: p
    character(len=10) :: title
    real :: tempnum
    logical,external :: GetNextPlayer !下一个球员数据函数
    integer :: i
    integer :: error

    inquire(file=filename,exist=alive)
    if(.not. alive) then
        write(*,*) trim(filename),"doesn't exist."
        stop
    end if

    open(unit=fileid,file=filename)
    do while(.true.)
        if(GetNextPlayer(fileid,p%name)) then
            do i=1,3
                read(fileid,"(A6,1X,F5.1)",iostat=error) title,tempnum
                if(error/=0) then
                    write(*,*) "文件读取错误!"
                    stop
                end if
                select case(title)
                case('身高')
                    p%height = tempnum
                case('体重')
                    p%weight = tempnum
                case('得分')
                    p%score = tempnum
                case default
                    write(*,*) "出现错误数据!"
                    stop
                end select
            end do
        else
            exit !没有数据,离开循环
        end if
        if(p%score>20.0) then
            write(*,"('姓名:'A10/,'身高:',F5.1,'体重:',F5.1,'得分:',F4.1)") p
        end if
    end do

    stop
end program main

! GetNextPlayer找出下一位球员的数据位置,如果有数据,返回.true.,没有数据,返回.false.
logical function GetNextPlayer(fileid,name)
    implicit none
    integer,intent(in) :: fileid
    character(len=*),intent(out) :: name
    character(len=80) :: title
    integer :: error
    do while(.true.)
        read(fileid,"(A80)",iostat=error) title
        if(error/=0) then
            GetNextPlayer = .false.
            return
        end if
        if(title(1:6) == '姓名') then
            name=title(8:)
            GetNextPlayer = .true.
            return
        end if
    end do
    return
end function GetNextPlayer

image-20220324111136056

  • 编写一个程序,可以同时读写文件,将文件内容前面添加行号后,输出到另一个文件里

    program main
       implicit none
       integer,parameter :: inputfileid = 10,outputfileid = 20
       integer,parameter :: maxbuffer = 200
       character(len=80) :: inputfile,outputfile
       character(len=maxbuffer) buffer
       integer :: count,error
       logical :: alive
    
       write(*,*) "input filename:"
       read(*,"(A80)") inputfile
       write(*,*) "output filename:"
       read(*,"(A80)") outputfile
    
       inquire(file=inputfile,exist=alive)
       if(.not. alive) then
           write(*,*) trim(inputfile),"doesn't exist."
           stop
       end if
    
       open(unit=inputfileid,file=inputfile,status="old")
       open(unit=outputfileid,file=outputfile,status="replace")
       count=1
       do while(.true.)
           read(inputfileid,"(A200)",iostat=error) buffer
           if(error/=0) exit
           write(outputfileid,"(I3,'.',A)") count,trim(buffer)
           count=count+1
       end do
       close(inputfileid)
       close(outputfileid)
       stop
    end program main

image-20220324202919647

直接访问文件的操作

  • 把文件的空间、内容事先分区成好几个同样大小的小模块,这些模块会自动按顺序编号,读写文件时,要先赋值文件读写位置在第几个模块,再来进行读写的工作。

  • 编写一个可以经过棒次来查询选手打击率的程序

image-20220324210336609

 program main
     implicit none
     integer,parameter :: fileid = 10
     character(len=20) :: filename = "list.txt"
     integer :: player,error,length=5
     real :: hit
     logical :: alive

     inquire(file=filename,exist=alive)
     if(.not. alive) then
         write(*,*) trim(filename),"doesn't exist."
         stop
     end if

     open(unit=fileid,file=filename,access="direct",form="formatted",recl=length,status="old")
     do while(.true.)
         write(*,*) "查询第几棒?"
         read(*,*) player
         read(fileid,fmt="(F4.2)",rec=player,iostat=error) hit
         if(error/=0) exit
         write(*,"('打击率:',F4.2)") hit
     end do
     close(fileid)
     stop
 end program main

image-20220324210355925

  • 使用直接访问文件时,要小心使用endfile命令,使用这个命令,会把目前所在的文件位置之后的数据都清除掉,例如使用read/write(···,rec=1)等命令,会把读写位置放在文件头,再使用endfile,这个文件后面的数据都会被清除掉。

二进制文件的操作

  • open命令中的recl字段所设置的整数n会随编译器不同而改变。
  • 输入棒球选手,用二进制文件来记录
program main
    implicit none
    integer,parameter :: fileid=10
    character(len=20) :: filename="list.bin"
    integer :: player,reclnum=4
    real :: hit(9)=(/3.2, 2.8, 3.3, 3.2, 2.9, 2.7, 2.2, 2.3, 1.9/)
    open(unit=fileid,file=filename,form="unformatted",access="direct",recl=reclnum,status="replace")
    do player=1,9
        write(fileid,rec=player) hit(player)
    end do
    close(fileid)
    stop
end program main

image-20220327084954854

  • 查询list.bin中的数据
program main
    implicit none
    integer,parameter :: fileid = 10
    character(len=20) :: filename = "list.bin"
    real :: hit
    integer :: player,error
    logical :: alive
    inquire(file=filename,exist=alive)
    if(.not. alive) then
        write(*,*) trim(filename),"doesn't exist."
        stop
    end if

    open(unit=fileid,file=filename,form="unformatted",access="direct",recl=4,status="old")
    do while(.true.)
        write(*,"('第几棒?')")
        read(*,*) player
        read(fileid,rec=player,iostat=error) hit
        if(error/=0) exit
        write(*,"('打击率:',F5.2)") hit
    end do

    close(fileid)
    stop
end program main

image-20220327085851566

  • 使用二进制数,可以减少内存空间,而且不会造成数据流失,要存放精确及大量的数据时,二进制文件是最优选择

Internal File(内部文件)

  • 内部文件是直接从英文原义译成中文的名词,或许叫它字符串变量文件更为贴切,写入文件的方法是把数据写到一个字符串变量中。
program main
    implicit none
    integer :: a=2,b=3
    character(len=20) :: string
    write(unit=string,fmt="(I2,'+',I2,'=',I2)") a,b,a+b !将数据写入字符串中
    write(*,*) string
    stop
end program main

image-20220327091148254

  • 从字符串中读取数据
program main
    implicit none
    integer :: a
    character(len=20) :: string='123'
    read(string,*) a
    write(*,*) a
    stop
end program main

image-20220327091455088

  • read命令输入时,如果数值型输入为字符串,程序就会死机,为此,可以将输入作为字符串,再判断是否合法,不合法就重新输入
program main
    implicit none
    integer :: i
    integer,external :: GetInteger
    i = GetInteger()
    write(*,*) i
    stop
end program

integer function GetInteger()
    implicit none
    character(len=80) :: string
    logical :: invalid
    integer :: i,code
    invalid=.true.
    do while(invalid)
        write(*,*) "请输入正整数:"
        read(*,"(A80)") string
        invalid=.false.
        do i=1,len_trim(string)
            code = ichar(string(i:i))
            if(code<ichar('0').or.code>ichar('9')) then
                invalid=.true.
                exit
            end if
        end do
    end do

    read(string,*) GetInteger
    return
end function GetInteger

image-20220327092422147

  • 内部文件还可以应用在动态改变输出格式,输出格式可以事先存放在字符串中,程序进行时,动态改变字符串内容就可以改变输出格式。
program main
    implicit none
    integer :: a,b
    character(len=30) :: fmtstring="(I??,'+',I??,'=',I??)"
    integer,external :: GetInteger

    a=GetInteger()
    b=GetInteger()

    !自定义格式输出
    write(fmtstring(3:4),"(I2.2)") int(log10(real(a))+1)
    write(fmtstring(11:12),"(I2.2)") int(log10(real(b))+1)
    write(fmtstring(19:20),"(I2.2)") int(log10(real(a+b))+1)
    write(*,fmtstring) a,b,a+b
    stop
end program main

integer function GetInteger()
    implicit none
    character(len=80) :: string
    logical :: invalid
    integer :: i,code
    invalid=.true.

    do while(invalid)
        write(*,*) "请输入正整数:"
        read(*,"(A80)") string
        invalid=.false.
        do i=1,len_trim(string)
            code = ichar(string(i:i))
            if(code<ichar('0').or.code>ichar('9')) then
                invalid=.true.
                exit
            end if
        end do
    end do

    read(string,*) GetInteger
    return
end function GetInteger

image-20220327095531097

Namelist

  • namelist是很特殊的输入、出方法,把一组相关变量封装在一起,输入、出这一组变量时,只要在write、read中的nml字段赋值要使用哪一个namelist就行了。
  • namelist是声明的一部分,必须编写在程序执行命令前面,且必须要命名
  • namelist /nl_name/ var1,var2,···
program main
    implicit none
    integer :: a=1,b=2,c=3
    namelist /na/ a,b,c !声明
    write(*,nml=na) !调用
    stop
end program main

image-20220327100124204

  • namelist可以用来输入数据,不过通常用来读取文件,不会在键盘输入。先介绍一下键盘输入的实例
program main
    implicit none
    integer :: a,b,c
    namelist /na/ a,b,c
    read(*,nml=na)
    write(*,nml=na)
    stop
end program main

image-20220327100924603

  • 上面可以看出,键盘输入namelist比较麻烦,所以namelist常用在文本文件的输入/输出中,使用read从文件中读取数据时,会自动从目前的位置向下寻找存放namelist的地方

image-20220327115715449

program main
    implicit none
    integer :: a(3)
    logical alive
    character(len=20) :: filename = "namelist1.txt"
    namelist /na/ a

    inquire(file=filename,exist=alive)
    if(.not. alive) then
        write(*,*) trim(filename),"doesn't exist"
    end if
    open(10,file=filename)
    read(10,nml=na)
    write(*,"(3I2)") a
    close(10)

    stop
end program main

文件的应用

  • 编写程序读取成绩单,计算每位同学的总分,及各科的全班平均分数

image-20220327112245312

module typedef
    type student
        integer :: chinese,english,math,natural,social
        integer :: total
    end type
end module typedef

program main
    use typedef
    implicit none
    integer,parameter :: fileid=80
    integer,parameter :: students=5
    character(len=80) :: tempstr
    type(student) :: s(students)
    type(student) :: total
    integer :: i,num,error

    open(fileid,file="grades.txt",status="old",iostat=error)
    if(error/=0) then
        write(*,*) "open grade.txt fail"
        stop
    end if

    read(fileid,"(A80)") tempstr
    total = student(0,0,0,0,0,0)
    do i=1,students,1
        read(fileid,*) num,s(i)%chinese,s(i)%english,s(i)%math,s(i)%natural,s(i)%social
        s(i)%total = s(i)%chinese+s(i)%english+s(i)%math+s(i)%natural+s(i)%social
        total%chinese = total%chinese+s(i)%chinese
        total%english = total%english+s(i)%english
        total%math = total%math+s(i)%math
        total%natural = total%natural+s(i)%natural
        total%social = total%social+s(i)%social
        total%total = total%total+s(i)%total
    end do

    write(*,"(7A7)") "座号","中文","英文","数学","自然","社会","总分"
    do i=1,students
        write(*,"(7I7)") i,s(i)
    end do

    write(*,"(A7,6F7.1)") "平均",real(total%chinese)/real(students),real(total%english)/real(students),&
            real(total%math)/real(students),real(total%natural)/real(students),real(total%social)/real(students),&
            real(total%total)/real(students)
    stop
end program main

image-20220327112306167

  • 编写程序来读取图文件,并且把图片的色彩值反面,并输出
program main
    use typedef
    implicit none
    integer,parameter :: recl_unit = 4
    integer,parameter :: buffer_size = 256*256
    character :: cbuffer(buffer_size)
    integer :: ibuffer(buffer_size)
    integer :: error,i,code,start=1

    open(10,file="photo1.raw",form="unformatted",access="direct",recl=256*256/recl_unit,status="old",iostat=error)
    if(error/=0) then
        write(*,*) "open photo1.raw fail"
        stop
    end if
    !一个像素占一个byte,刚好可以用字符数组把整个文件读入
    read(10,rec=start) cbuffer
    close(10)

    do i=1,buffer_size
        code = ichar(cbuffer(i)) !code返回值会是-128~127,需要转化成0~255之间的数字
        if(code>=0) then
            ibuffer(i)=code
        else
            ibuffer(i)=256+code
        end if
    end do

    !把亮度值反相
    do i=1,buffer_size
        cbuffer(i)=char(255-ibuffer(i))
    end do

    open(10,file="photo1.raw",form="unformatted",access="direct",recl=256*256/recl_unit,status="replace")
    write(10,rec=start) cbuffer
    close(10)

    stop
end program main

习题

  • 1、请改写会显示文本文件内容的实例程序,加入一个功能,在输出行数到一定程序时(占满屏幕的一页),能暂停输出,等用户按enter键后再继续输出。
program main
    implicit none
    character(len=20) :: filename
    character(len=79) :: buffer
    integer,parameter :: fileid=10
    integer :: status = 0,page=0
    logical :: alive
    write(*,*) "please input filename:"
    read(*,"(A20)") filename
    inquire(file=filename,exist=alive)
    if(.not. alive) then
        write(*,*) trim(filename)," doesn't exist"
    else
        open(unit=fileid,file=filename,access="sequential",form="formatted",status="old")
        do while(.true.)
            page=page+1
            read(unit=fileid,fmt="(A79)",iostat=status) buffer
            if(status/=0) then
                exit
            end if
            if(page==6) pause
            write(*,"(A79)") buffer
        end do
    end if
    close(fileid)
    stop
end program main

image-20220327155702519

  • 2、一个经过编码的文本文件,编码方式是把每个字符的ASCII值加上3之后再输出,请编写程序来把这个文件解码。

  • 3、以二进制方式保存了全班20位同学的成绩。包括中文,英文,数学,自然,社会这5个科目的成绩。每个成绩使用长整型方法保存(占4bytes),根据顺序先存放1号同学的5个科目成绩,再存放2号同学的5个科目成绩····,最后存放20号同学的5个科目成绩。请编写程序读出全班同学的成绩,并计算每位同学的总分及全班的各科平均分。

  • 4、一个经过编码的文本文件,编码的方法是根据顺序把每1行第3n-2,2n-1,3n个字符的ASCII值分别加上1,2,3之后再输出,请编写程序把这个文件译码。

    Ex:Hello
    第1个字符 H,编码后char(char("H")+1)="I"
    第2个字符 e,编码后char(char("e")+2)="g"
    ···

    简单地说,编码时会把每个字符按照顺序加上1或2或3之后再输出

  • 5、把文件应用的第一个实例程序加上排名次的功能。

指针

  • 可以用来保存变量,或者是动态使用内存,更进一层则可以应用在特别的“数据结构”上,例如创建“串行结构”、”树状结构“等等

指针基本概念

  • 指针是一种”间接使用数据“的方法。指针变量用来保存一个”内存地址“,当程序要读写指针变量时,实际上会经过两个步骤

    • 1⃣️取出指针中所保存的内存地址
    • 2⃣️到这个内存位置读写数据

    指针变量中所保存的内存地址来源有两种

    • 1⃣️记录其他非指针变量的内存地址
    • 2⃣️程序执行中动态配置一块内存
  • 基本的指针用法,是把指针变量拿来记录另一个目标变量的地址,再经过指针来读写数据

program main
    implicit none
    integer,target :: a=1 !声明一个可以当成目标的变量
    integer,pointer :: p !声明一个可以指向整数的指针
    p=>a !把指针p指向a
    write(*,*) p
    a=2
    write(*,*) p
    p=3
    write(*,*) a
    stop
end program main

image-20220327171533674

  • integer,target :: a=1,加上target后,变量可以把它的内存地址赋值给指针变量
  • 动态配置内存
program main
    implicit none
    integer,pointer :: p
    allocate(p) !配置内存,deallocate(p),释放内存
    p=100
    write(*,*) p
    stop
end program main

使用动态配置内存的方法,对于指针指向问题,需要先释放内存,再重新指向,即先给指针变量一个内存,对于指向哪个变量先不声明

  • associated(pointer,[target]):用来检查指针是否已经设置指向,返回值为布尔变量。如果只放入第一个指针参数,会检查这个指针是否已经赋值好“方向”,如果放入两个参数,则会检查第一个指针变量是否指向第二个变量。

    • 函数null()会返回一个不能使用的内存地址,即可确保associated的正确性。 integer,pointer :: p=>null()

    • 函数nullify(pointer1,[pointer2,···])也可以把指针设置成没有指向的任何内存地址,fortran 90只能使用nullify

      • integer,pointer :: p
        nullify(p)
program main
    implicit none
    integer,pointer :: a=>null()
    integer,target :: b=1,c=2

    write(*,*) associated(a)
    a=>c
    write(*,*) associated(a)
    write(*,*) associated(a,c)
    write(*,*) associated(a,b)
    nullify(a)
    write(*,*) associated(a)

    stop
end program main

image-20220327173853317

  • 指针可以声明成任何数据类型,不管是什么类型,指针变量占用相同的内存空间,只存放内存地址,而非数据。

指针数组

  • 指针也可以声明成数组,声明成数组的指针同样有两种使用方法
    • 1⃣️把指针指到其他数组
    • 2⃣️配置内存空间来使用
program main
    implicit none
    integer,pointer :: a(:)
    integer,target :: b(5)=(/1,2,3,4,5/)
    a=>b !把指针数组a指向数组b
    write(*,*) a
    a=>b(1:3)
    write(*,*) a
    a=>b(1:5:2)
    write(*,*) a
    a=>b(5:1:-1) !指针倒序输出
    write(*,*) a
    stop
end program main

image-20220327174733705

program main
    implicit none
    integer,pointer :: a(:)
    allocated(a(5)) !配置5个整数的内存空间给指针a
    a=(/1,2,3,4,5/)
    write(*,*) a
    deallocate(a) !归还内存空间
    stop
end program main

image-20220327175141168

  • 多维指针数组的使用方法
program main
    implicit none
    integer,pointer :: a(:,:)
    integer,target :: b(3,3,2)
    integer :: i

    forall(i=1:3)
        b(:,i,1) = i
        b(:,i,2) = 2*i
    end forall

    a=>b(:,:,1)
    write(*,"(9I2)") a
    a=>b(1:3:2,1:2,2)
    write(*,"(4I2)") a
    stop
end program main

image-20220327180250287

指针与函数

  • 指针变量一样可以作为参数在函数之间传递,也可以作为函数的返回值。使用时有下面几个策略:
    • 1⃣️要把指针传递给函数时,要声明这个函数的参数使用接口interface
    • 2⃣️指针参数声明时不需要intent这个形容词
    • 3⃣️函数返回值若为指针时,需要定义函数的interface
program main
    implicit none
    integer,target :: a(8)=(/10,15,8,25,9,20,17,19/)
    integer,pointer :: p(:)
    interface
        function getmin(p)
            integer,pointer :: p(:)
            integer,pointer :: getmin
        end function getmin
    end interface
    p=>a(1:8:2)
    write(*,*) getmin(p)

    stop
end program main

function getmin(p)
    implicit none
    integer,pointer :: p(:)
    integer,pointer :: getmin
    integer :: i,s,min

    s=size(p,1) !查询数组的大小
    min=2**30 !先把min设置成很大的值
    do i=1,s
        if(min>p(i)) then
            min=p(i)
            getmin=>p(i)
        end if
    end do
    return
end function

image-20220327182256807

  • 编写interface是一件很麻烦的工作,如果把函数封装在module中,就等于已经编写好使用接口
module func

contains
    function getmin(p)
        implicit none
        integer,pointer :: p(:)
        integer,pointer :: getmin
        integer i,s,min
        s=size(p,1)
        min=2**30
        do i=1,s
            if(min>p(i)) then
                min=p(i)
                getmin=>p(i)
            end if
        end do
        return
    end function getmin
end module

program main
    use func
    implicit none
    integer,target :: a(8)=(/10,15,8,25,9,20,17,19/)
    integer,pointer :: p(:)
    p=>a(1:8:2)
    write(*,*) getmin(p)
    stop
end program main

image-20220327183247205

基本的指针应用

  • 指针指向自定义函数时,指针可以很快地交换数据。
  • 以自定义类型数据来做排序的实例,排序程序常常会需要把两条数据交换,如果不使用指针,交换数据时需要移动很大块的内存空间。
module func
    !自定义类型person
    type person
        character(len=10) :: name
        real :: height,weight
    end type person

    type pperson
        type(person),pointer  :: p
    end type pperson

contains

    subroutine sort(p)
        implicit none
        type(pperson) :: p(:)
        type(pperson) :: temp
        integer :: i,j,s
        s=size(p,1)
        do i=1,s-1
            do j=i+1,s
                if(p(j)%p%height<p(i)%p%height) then
                    temp = p(i)
                    p(i)=p(j)
                    p(j)=temp
                end if
            end do
        end do
        return
    end subroutine sort

end module func

program main
    use func
    implicit none
    type(person),target :: p(5)=(/person("陈同学",180.0,75.0),person("黄同学",170.0,65.0),person("刘同学",175.0,80.0),&
            person("蔡同学",182.0,78.0),person("许同学",178.0,70.0)/)

    type(pperson) :: pt(5)
    integer :: i
    forall(i=1:5)
        pt(i)%p => p(i)
    end forall
    call sort(pt)
    write(*,"(5(A9,F6.1,F5.1/))") (pt(i)%p,i=1,5)
    stop
end program main

image-20220328094734304

  • 指针排序在做数据交换时,只需要交换两条数据的内存地址,不需要去移动这两条数据的内存,当自定义类型中的数据量很大时,执行效率就可以有明显地提升。

指针的高级应用

  • 指针除了可以间接使用变量,以及当成可变大小数组来使用之外,还有一个很重要的用途,它可以创建各种的“串行结构”以及“树状结构”等等
  • 例如,串行结构可以动态规划内存,解决数组过大造成的内存错误,数组过小造成内存浪费

单向串行

  • 数组声明好大小后,发现数组不够用,如果使用deallocate就会丢失数据
  • 在数组中插入一个数据,就需要把数组后面的数据向后移动,数据过多时,极大地影响了执行效率
  • 串行结构可以解决这些问题,来看一个最简单的单向串行结构
module typedef
    implicit none
    type :: datalink
        integer :: i
        type(datalink),pointer :: next
    end type datalink
end module typedef


program main
    use typedef
    implicit none
    type(datalink),target :: node1,node2,node3
    integer :: i
    !单向串行相当于C语言的单链表
    node1%i=1
    node1%next=>node2
    node2%i=2
    node2%next=>node3
    node3%i=3
    nullify(node3%next)

    write(*,*) node1%i
    write(*,*) node1%next%i
    write(*,*) node1%next%next%i
    stop
end program main

image-20220328104328097

  • 使用循环改写程序
module typedef
    implicit none
    type :: datalink
        integer :: i
        type(datalink),pointer :: next
    end type datalink
end module typedef

program main
    use typedef
    implicit none
    type(datalink),target :: node1,node2,node3
    type(datalink),pointer :: p
    integer :: i

    p=>node1
    node1%i=1
    node1%next=>node2
    node2%i=2
    node2%next=>node3
    node3%i=3
    nullify(node3%next)
    do while(.true.)
        write(*,*) p%i
        if(.not. associated(p%next)) exit
        p=>p%next
    end do
    stop
end program main

image-20220328110102740

  • 典型的创建串行方法
module typedef
    implicit none
    type :: datalink
        integer :: i
        type(datalink),pointer :: next
    end type datalink
end module typedef

program main
    use typedef
    implicit none
    type(datalink),pointer :: p,head
    integer :: i,n,err

    !创建头指针
    write(*,*) "Input N:"
    read(*,*) n
    allocate(head)
    head%i=1
    nullify(head%next)

    p=>head
    do i=2,n
        allocate(p%next,stat=err)
        if(err/=0) then
            write(*,*) "Out of memory!"
            stop
        end if
        p=>p%next
        p%i=i
    end do
    nullify(p%next)

    p=>head
    do while(.true.)
        write(*,"(i5)") p%i
        if(.not. associated(p%next)) exit
        p=>p%next
    end do

    stop
end program main

image-20220328111909772

双向串行、环状串行

  • 单向串行只能沿一个方向走,一条接着一条往下读,没有办法往回读,可以回头的串行叫双向串行,首尾相连的叫环状串行。
module typedef
    implicit none
    type :: datalink
        integer :: i
        type(datalink),pointer :: prev !指向上一个数据
        type(datalink),pointer :: next !指向下一个数据
    end type datalink
end module typedef

program main
    use typedef
    implicit none
    type(datalink),target :: node1,node2,node3
    type(datalink),pointer :: p
    integer :: i
    !创建双向串行
    node1=datalink(1,null(),node2)
    node2=datalink(2,node1,node3)
    node3=datalink(3,node2,null())

    write(*,*) "按顺序输出"
    p=>node1
    do while(.true.)
        write(*,*) p%i
        if(.not. associated(p%next)) exit
        p=>p%next
    end do

    write(*,*) "反向输出"
    p=>node3
    do while(.true.)
        write(*,*) p%i
        if(.not. associated(p%prev)) exit
        p=>p%prev
    end do

    stop
end program main

image-20220328151953501

  • 环状串行
module typedef
    implicit none
    type :: datalink
        integer :: i
        type(datalink),pointer :: prev
        type(datalink),pointer :: next
    end type datalink
end module typedef

program main
    use typedef
    implicit none
    type(datalink),target :: node1,node2,node3
    type(datalink),pointer :: p
    integer,parameter :: s=6
    integer :: i

    !创建环状串行
    node1=datalink(1,node3,node2)
    node2=datalink(2,node1,node3)
    node3=datalink(3,node2,node1)

    write(*,*) "从前往后输出"
    p=>node1
    do i=1,s
        write(*,*) p%i
        if(.not. associated(p%next)) exit
        p=>p%next
    end do

    write(*,*) "从后往前输出"
    p=>node3
    do i=1,s
        write(*,*) p%i
        if(.not. associated(p%next)) exit
        p=>p%prev
    end do

    stop
end program main

image-20220328153750283

插入及删除

  • 串行可以快速地插入或删除一条数据,如果是数组就需要移动后面的数据,对于串行,只需要把插入位置的前一个尾指向插入,插入位指向原数据即可

  • 实例,插入数据跟删除数据,独立写成两个函数

module linklist
    implicit none
    type :: datalink
        integer :: i
        type(datalink),pointer :: prev
        type(datalink),pointer :: next
    end type datalink

contains

    subroutine outputlist(list)
        implicit none
        type(datalink),pointer :: list,p
        p=>list
        do while(associated(p))
            write(*,*) p%i
            p=>p%next
        end do
        return
    end subroutine outputlist

    subroutine delitem(item)
        implicit none
        type(datalink),pointer :: item
        type(datalink),pointer :: prev,next

        prev=>item%prev !记录item上一条数据的位置
        next=>item%next !记录item下一条数据的位置
        deallocate(item) !释放item所占内存
        !删除item数据,将上一个数据next指向item的下一个数据prev;将item的下一个数据的prev指向item的prev,即上一个数据的next
        if(associated(prev)) prev%next=>next
        if(associated(next)) next%prev=>prev
        item=>next !item指向next,即删除item
        return
    end subroutine delitem

    subroutine insitem(pos,item,after)
        implicit none
        type(datalink),pointer :: pos,item
        logical :: after
        if(after) then
            !把item插在pos的后面
            item%next=>pos%next
            item%prev=>pos
            if(associated(pos%next)) then  !当pos不是尾时
                pos%next%prev => item
            end if
            pos%next => item
        else
            !item插在pos前面
            item%next=>pos
            item%prev=pos%prev
            if(associated(pos%prev)) then !pos与前一个数据未断开,即不是首数据
                pos%prev%next=>item
            end if
            pos%prev=>item
        end if
        return
    end subroutine insitem

end module linklist

program main
    use linklist
    implicit none
    type(datalink),pointer :: head
    type(datalink),pointer :: item,p
    integer,parameter :: s=5
    integer :: i,n,error

    allocate(head)
    head=datalink(1,null(),null())
    !创建串行
    p=>head
    do i=2,s
        allocate(p%next,stat=error)
        if(error/=0) then
            write(*,*) "out of memory"
            stop
        end if
        p%next = datalink(i,p,null())
        p=>p%next
    end do


    write(*,*) "删除第三个数据"
    call delitem(head%next%next)
    call outputlist(head)

    write(*,*) "插入第三个数据"
    allocate(item)
    item%i=30
    call insitem(head%next,item,.true.)
    call outputlist(head)

    stop
end program main

image-20220328164523256

  • 注意:只有通过allocate函数所配置到的内存才能使用deallocate来释放。

串行的应用

  • 串行可以比较灵活地使用内存,不需要考虑需要记录多少条数据
  • data1.txt和data2.txt是两个班级的段考成绩单,两班的人数不同,请编写一个可以读取成绩的程序,让用户输入文件名来决定要读取哪一个文件,还要提供给用户通过座号来查询成绩的功能。
module linklist
    type student
        integer :: num
        integer :: chinese,english,math,science,social
    end type

    type datalink
        type(student) :: item
        type(datalink),pointer :: next !单向串行
    end type datalink

contains

    function searchList(num,head)
        implicit none
        integer :: num
        type(datalink),pointer :: head,p
        type(datalink),pointer :: searchList

        p=>head
        nullify(searchList)
        do while(associated(p))
            if(p%item%num==num) then
                searchList=>p
                return
            end if
            p=>p%next
        end do
        return
    end function searchList

end module linklist

program main
    use linklist
    implicit none
    character(len=20) :: filename
    character(len=80) :: tempstr
    type(datalink),pointer :: head
    type(datalink),pointer :: p
    integer :: i,error,size

    write(*,*) "please input filename:"
    read(*,*) filename
    open(10,file=filename,status="old",iostat=error)
    if(error/=0) then
        write(*,*) "open file fail!"
        stop
    end if

    allocate(head)
    nullify(head%next)
    p=>head
    size=0
    read(10,"(A80)") tempstr !读入第一行字符串,不需要处理它
    do while(.true.)
        read(10,fmt=*,iostat=error) p%item
        if(error/=0) exit
        size = size+1
        allocate(p%next,stat=error)
        if(error/=0) then
            write(*,*) "out of memory!"
            stop
        end if
        p=>p%next
        nullify(p%next)
    end do

    write(*,"('总共有',I3,'位学生')") size

    do while(.true.)
        write(*,*) "请输入要查询的学生座号:"
        read(*,*) i
        if(i<1.or.i>size) exit
        p=>serchList(i,head)
        if(associated(p)) then
            write(*,"(5(A6,I3))") "中文",p%item%chinese,"英文",p%item%english,"数学",p%item%math,&
            "自然",p%item%science,"社会",p%item%social
        else
            exit !找不到数据,退出循环
        end if
    end do

    write(*,"('座号',I3,'不存在,程序结束。')") i
    stop
end program main
  • 比较好的编写方法,应该是先使用串行来读取文件,读完文件之后,就会知道学生数目,这时候就可以使用另一个可变大小数组来复制串行中的学生成绩,接着再把串行全部删除,查询成绩时直接使用数组来查询就可以了。

习题

  • 1、请问下面的变量,在目前的PC中分别会占用多少内存空间?

integer(kind=4) :: a

real(kind=4) :: b

real(kind=8) :: c

character(len=10) :: str

integer(kind=4),pointer :: pa

real(kind=4),pointer :: pb

real(kind=8),pointer :: pc

character(len=10),pointer :: pstr

type student
integer :: chinese,english,math
end type student
type(student) :: s
type(student),pointer :: ps
  • 2、请说明下面程序段的执行结果。
integer,target :: a=1
integer,target :: b=2
integer,target :: c=3
integer,pointer :: p
p=>a
write(*,*) p
p=>b
write(*,*) p
p=>c
p=5
write(*,*) c

image-20220329103208041

  • 3、想办法把串行的应用的实例改成使用可变大小数组来记录学生数据
  • 4、将单向串行的典型创建串行的实例中加入deallocate函数来释放整个串行的内存

Module及面向对象

结构化与面向对象

  • 结构化与面向对象,是目前设计程序最常使用的两种编写概念。

结构化程序设计概念

  • 结构化程序的特色在于“层次分明”,检查程序代码时,可以把它们分成不同的程序模块
program main
   ······
   ······
   do while(.true.)
      ···
      ···
   end do

   ······
   ······
end program

if(···) then
  ···
  ···
else
  ···
  ···
end if


select case(a)
case(1)
   ····
   ····
case(2)
   ····
   ····
end select

面向对象程序设计概念

  • 新一代的程序语言,除了具有原来的结构化程序设计方法外,还加入了“面向对象”概念。简单地说,面向对象是在做程序代码封装的操作,封装后的程序代码,在使用上会比较安全。
  • 面向对象中很重要的一项工作,就是数据封装。数据经过封装后可以分成两种数据,一种是可以让大家使用的数据,另一种是只能在内部使用的数据,函数亦是如此。
  • 面向对象的两大方向
    • 为了安全,有些数据不能公开,使用封装
    • 经过继承来重复使用程序代码

再论module

module的结构及功能

  • 1、module里面可以声明变量,经常用来声明程序中所需要的常量,或是用来存放全局变量
  • 2、module里面可以定义自定义类型,再经过use命令让程序中的每一个函数都能使用这个类型
  • 3、module里面可以编写函数,通常会把功能相关的函数放在同一个module中,在module外面调用这些函数时,同样要使用use命令
  • 4、module里面的函数,可以直接使用同一个module中所声明的变量,所以module里面的函数,可以经过module里面的变量来互相传递数据

public,private

  • module里面的数据和函数,可以通过publicprivate命令,来区分成分开使用及私下使用。
  • 实例:模拟到银行领钱
module bank
    implicit none
    private money
    public loadMoney,saveMoney,report
    integer :: money=1000000

contains
    subroutine loadMoney(num)  !取钱
        implicit none
        integer :: num
        money=money-num
        return
    end subroutine loadMoney

    subroutine saveMoney(num)  !存钱
        implicit none
        integer :: num
        money=money+num
        return
    end subroutine saveMoney

    subroutine report()  !查询余额
        implicit none
        write(*,"('银行目前库存',I7,'元')") money
        return
    end subroutine report

end module bank

program main
    use bank
    implicit none
    call loadMoney(100)
    call saveMoney(1000)
    call report()
    stop
end program main

image-20220329151942842

  • 程序中money只有在module bank中使用,主函数无法使用,不指明默认public
  • 取款和存款,银行都会记录,实例,比较完整的程序
module bank
    implicit none
    integer :: money =1000000
    integer,parameter :: fileid = 10
    private money,fileid
    private timeLog

contains

    subroutine timeLog()
        implicit none
        integer :: num
        character(len=20) :: date,time
        call date_and_time(date,time)
        write(fileid,"('Data:',A8,'Time:',A2,':',A2,':'A2)") date,time(1:2),time(3:4),time(5:6)
        return
    end subroutine timeLog


    subroutine loadMoney(name,num)
        implicit none
        character(len=*) :: name
        integer :: num
        if(num<=0) then
            write(*,*) "不合理的金额"
            return
        end if
        open(fileid,file="log.txt",position="append")
        if(money>=num) then
            call timeLog()
            write(fileid,"(A10,'领取',I5,'元')") name,num
            money=money-num
        else
            write(fileid,*) "银行目前现金不足"
            write(*,*) "银行目前现金不足"
        end if
        close(fileid)
        return
    end subroutine loadMoney

    subroutine saveMoney(name,num)
        implicit none
        character(len=*) :: name
        integer ::  num
        if(num<=0) then
            write(*,*) "不合理的金额"
            return
        end if
        open(fileid,file="log.txt",position="append")
        call timeLog()
        write(fileid,"(A10,'存入',I5,'元')") name,num
        close(fileid)
        money=money+num
        return
    end subroutine saveMoney
end module bank

program main
    use bank
    implicit none
    call loadMoney("彭先生",100)
    call saveMoney("陈先生",1000)
    stop
end program main

image-20220329154513327

use

  • 编写好module需要use命令才能使用,module中也可以使用另一个module

  • 语法

module A
  implicit none
  integer :: a,b
end module A

module B
  use A
end module B

subroutine sub()
  use B
  implicit none
end subroutine sub
  • 使用module的数量并没有限制,使用多个module会造成变量名称重复的问题,可以在use命令后面,将变量名或函数名临时改名

use A,aa=>va:将module中的变量va临时改名为aa

  • 当重复变量较多时,并且不需要用到module中的所有东西,可以选择module中的一部分东西来使用
module A
implicit none
integer va,vb,vc
end module A

module B
implicit none
integer :: va,vb
end module B

program main
use A,only : vc
use B
implicit none
end program main
  • module A使用module B,可以认为module A继承了module B的数据和函数,不过继承的东西只限制在module B中对外公开的变量及函数,module B私下使用的不会被继承,module A继承module B的原有功能后,可以添加一些函数来扩充功能

    module MA
        implicit none
        real a,b
    contains
        subroutine getx()
            write(*,"('x=',F5.2)") -b/a
            return
        end subroutine getx
    end module MA
    
    module MB
        use MA
        implicit none
        real :: c
    contains
        subroutine getx2()
            real :: a2,d,sqrt_d
            a2=2*a
            d=b*b-4*a*c
            if(d>=0) then
                sqrt_d=sqrt(d)
                write(*,"('x=',F5.2,',',F5.2)") (-b+sqrt_d)/a2,(-b-sqrt_d)/a2
            else
                write(*,*) "无实数解"
            end if
        end subroutine getx2
    end module MB
    
    subroutine sub1()
        use MA
        implicit none
        a=2.0
        b=3.0
        call getx()
        return
    end subroutine sub1
    
    subroutine sub2
        use MB
        implicit none
        a=1.0
        b=4.0
        c=4.0
        call getx2
        return
    end subroutine sub2
    
    program main
        implicit none
        call sub1()
        call sub2()
        stop
    end program main

image-20220329170817585

再论interface

  • interface是用来说明函数的参数及返回值类型,不过当函数封装在module里面时,就不需要再使用interface来做这些说明。

同名函数的重载(overload)

  • overload的意义是:“在程序代码中可以同时拥有多个名称相同,但是参数类型、数目不同的函数,程序会自动根据输入的参数,来决定要调用哪一个函数”
  • module中使用interface,可以用来定义一个虚拟的函数名称
module MA
    implicit none
    interface show !虚拟的函数名称show
        module procedure show_int  !等待选择的函数show_int
        module procedure show_character  !等待选择的函数show_character
    end interface

contains
    subroutine show_int(n)
        implicit none
        integer,intent(in) :: n
        write(*,"('n=',I3)") n
        return
    end subroutine show_int

    subroutine show_character(str)
        implicit none
        character(len=*),intent(in) :: str
        write(*,"('str=',A)") str
        return
    end subroutine show_character
end module MA

program main
    use MA
    implicit none
    call show_int(1) !输入的参数是整数,自动调用show_int
    call show(1)
    call show_character("fortran 95") !输入的参数是字符串,会自动调用show_character
    call show("fortran 95")
    stop
end program main

image-20220329201237335

  • 计算ax+b=0ax^2+bx+c=0这两个式子来看,可以分别显示两个函数来计算,计算ax+b=0的函数需要输入a、b,计算ax^2+bx+c=0的函数需要输入a、b、c的值
module MA
    implicit none
    interface getx
        module procedure getx1
        module procedure getx2
    end interface
contains
    subroutine getx1(a,b)
        real :: a,b
        write(*,"('x=',F5.2)") -b/a
        return
    end subroutine getx1

    subroutine getx2(a,b,c)
        real a,b,c
        real a2,d,sqrt_d
        a2=2*a
        d=b*b-4*a*c
        if(d>=0) then
            sqrt_d=sqrt(d)
            write(*,"('x=',F5.2,',',F5.2)") (-b+sqrt_d)/a2,(-b-sqrt_d)/a2
        else
            write(*,*) "无实数解"
        end if
    end subroutine getx2

end module MA

program main
    use MA
    implicit none
    call getx(1.0,2.0)
    call getx(1.0,3.0,2.0)
    stop
end program main

image-20220329203606724

自定义操作符

  • fortran 基本数值的数据类型,主要有integer、real这两种。使用这两种类型所声明出来的变量,除了可以用来保存数值外,还可以拿来做+,-,*,/数学运算及<,<=,>,>=,==,/=等等的逻辑判断。而使用type所声明的自定义类型,默认时不能拿来做这些运算。不过通过interface的帮忙,可以虚拟出上述的运算符号。
interface operator(+) !在程序代码中,使用a+b时,若a和b的参数符合下面任何函数的两个参数类型,会调用其中一个函数来执行
module procedure add1
module procedure add2
end interface
  • 实例
module MA
    implicit none
    type ta
        integer :: a
    end type ta
    interface operator(+) !这个interface让type(ta)类型变量也能相加
        module procedure add
    end interface
contains
    integer function add(a,b)
        type(ta),intent(in) :: a,b
        add=a%a+b%a
    end function add
end module MA

program main
    use MA
    implicit none
    type(ta) :: a,b
    integer :: c
    a%a=1
    b%a=2
    c=a+b !调用add(a,b)来执行
    write(*,*) c
    stop
end program main

image-20220329215441488

  • 注意:要把运算符号拿来当成虚拟函数名称时,interface后面要先接上operator这个字,再用括号()把运算符号包起来。另外在interface中等待候选的函数,必须明确显示每一个参数属性 intent
  • 实例:黄先生在这个月的5日及20日分别和许律师约谈了1小时45分、2小时18分。请问黄先生这个月花了多少时间和他的律师讨论有关他的遗产分配问题?
module time_util
    implicit none
    type :: time
        integer :: hour,minute
    end type time
    interface operator(+)
        module procedure add
    end interface
contains
    function add(a,b)
        implicit none
        type(time),intent(in) :: a,b
        type(time) :: add
        integer :: minutes,carry
        minutes = a%minute+b%minute
        carry = minutes/60
        add%minute=mod(minutes,60)
        add%hour=a%hour+b%hour+carry
        return
    end function add

    subroutine output(t)
        type(time),intent(in) :: t
        write(*,"(I2,':',I2.2)") t%hour,t%minute
        return
    end subroutine output

end module time_util

program main
    use time_util
    implicit none
    type(time) :: a,b,c
    a=time(1,45)
    b=time(2,18)
    c=a+b
    call output(c)
    stop
end program main

image-20220330093938298

  • 比较完整的实例程序,定义type(time)+real、real+type(time)、小于的判断及两种等号type(time)=real、real=type(time)操作
module time_util
    implicit none
    type :: time
        integer :: hour,minute
    end type time
    interface operator(+)
        module procedure add_time_time
        module procedure add_time_real
        module procedure add_real_time
    end interface

    interface operator(<)
        module procedure time_it_time
    end interface

    interface assignment(=)
        module procedure time_assign_real
        module procedure real_assign_time
    end interface

contains
    function add_time_time(a,b) result(add)
        implicit none
        type(time),intent(in) :: a,b
        type(time) :: add
        integer :: minutes,carry
        minutes=a%minute+b%minute
        carry=minutes/60
        add%minute=mod(minutes,60)
        add%hour=a%hour+b%hour+carry
        return
    end function add_time_time

    function add_time_real(a,b)
        implicit none
        type(time),intent(in) :: a
        real,intent(in) :: b
        type(time) :: add_time_real
        type(time) :: tb
        tb%hour=int(b)
        tb%minute=int((b-tb%hour)*60)
        add_time_real = add_time_time(a,tb)
        return
    end function add_time_real

    function add_real_time(a,b)
        implicit none
        real,intent(in) :: a
        type(time),intent(in) :: b
        type(time) :: add_real_time
        add_real_time = add_time_real(b,a)
        return
    end function add_real_time

    logical function time_it_time(a,b)
        implicit none
        type(time),intent(in) :: a,b
        if(a%hour<b%hour) then
            time_it_time = .true.
            return
        end if
        time_it_time = .false.
        return
    end function time_it_time

    subroutine time_assign_real(a,b)
        implicit none
        type(time),intent(out) :: a
        real,intent(in) :: b
        a%hour=int(b)
        a%minute=int((b-a%hour)*60)
        return
    end subroutine time_assign_real

    subroutine real_assign_time(a,b)
        implicit none
        real,intent(out) :: a
        type(time),intent(in) :: b
        a=b%hour+real(b%minute)/60.0
        return
    end subroutine real_assign_time

    subroutine output(t)
        type(time),intent(in) :: t
        write(*,"(I2,':',I2.2)") t%hour,t%minute
        return
    end subroutine output

end module time_util

program main
    use time_util
    implicit none
    type(time) :: a,b,c
    real :: rt
    a=0.5
    b=0.1+a
    c=a+0.6
    rt=time(1,30)+time(2,30)
    call output(c)
    write(*,*) rt
    write(*,*) a<b
    stop
end program main

image-20220330103327651

实际应用

继承module

  • 编写一个新的module来继承module TextGraphLib,增加把绘图结果输出到文件中的功能
module newGraphLib
    use TextGraphLib
    implicit none
contains
    subroutine outputToFile(filename)
        implicit none
        character(len=*),intent(in) :: filename
        character(len=10) :: fmt="(xxxA)"
        integer :: i
        if(.not. allocated(screen)) return
        open(10,file=filename,status="replace")
        write(fmt(2:4),"(I3.3)") ScreenWidth
        do i=1,ScreenHeight
            write(10,fmt) screen(i,:)
        end do
        close(10)
        return
    end subroutine outputToFile
end module newGraphLib

program main
    use newGraphLib
    implicit none
    call SetScreen(20,20)
    call ClearScreen()
    call DrawCircle(10,10,8)
    call DrawLine(14,6,6,14)
    call outputToFile("text.txt")
    stop
end program main

image-20220330110120167

自定义操作符的应用

  • fortran中没有提供“分数”类型,不过在现实世界中,常常会使用这个类型。分数可以保存更为精确的数值,例如2/3,分数转换成实数后,会变成循环小数,因为有效位的限制,所以为了精度,可以选择分数来存储
module rational_util
    implicit none
    private
    public :: rational,operator(+),operator(-),operator(*),operator(/),assignment(=),output
    type :: rational
        integer :: num !分子
        integer :: denom !分母
    end type rational
    !加法
    interface operator(+)
        module procedure rat_plus_rat
    end interface
    !减法
    interface operator(-)
        module procedure rat_minus_rat
    end interface
    !乘法
    interface operator(*)
        module procedure rat_times_rat
    end interface
    !除法
    interface operator(/)
        module procedure rat_div_rat
    end interface
        !等号
        interface assignment(=)
                module procedure int_eq_rat
                module procedure real_eq_rat
        end interface

                contains
                !整数=分数
                subroutine int_eq_rat(int,rat)
                implicit none
                integer,intent(out) :: int
                type(rational),intent(in) :: rat
        int = rat%num / rat%denom
        return
        end subroutine int_eq_rat

        !浮点数=分数
                subroutine real_eq_rat(float,rat)
                implicit none
                real,intent(out) :: float
                type(rational),intent(in) :: rat
        float = real(rat%num)/real(rat%denom)
                return
                end subroutine real_eq_rat

                !化简分数
                function reduse(a)
                implicit none
                type(rational),intent(in) :: a
        type(rational) :: temp
        integer :: b
        integer :: sign
                type(rational) :: reduse
                if(a%num*a%denom>0) then
                sign=1
                else
                sign=-1
                end if
                temp%num=abs(a%num)
        temp%denom=abs(a%denom)
                b=gcv(temp%num,temp%denom) !找出分子与分母的最大公因子
                reduse%num = temp%num/b*sign
                reduse%denom = temp%denom/b
                return
                end function reduse

                !用辗转相除法找最大公因子
                function gcv(a,b)
                implicit none
                integer,intent(in) :: a,b
        integer :: big,small
                integer :: temp
                integer :: gcv
                big=max(a,b)
                small=min(a,b)
                do while(small>1)
                temp=mod(big,small)
                if(temp==0) exit
                big=small
                small=temp
                        end do
                        gcv=small
                        return
                        end function gcv

                !分数相加
                function rat_plus_rat(rat1,rat2)
                        implicit none
                        type(rational) :: rat_plus_rat
                type(rational),intent(in) :: rat1,rat2
                type(rational) :: act
                act%denom = rat1%denom*rat2%denom
                act%num = rat1%num*rat2%denom+rat2%num*rat1%denom
                rat_plus_rat = reduse(act)
                        return
                        end function rat_plus_rat

                        !分数相减
                        function rat_minus_rat(rat1,rat2)
                        implicit none
                        type(rational) :: rat_minus_rat
                        type(rational),intent(in) :: rat1,rat2
                        type(rational) :: temp
                        temp%denom = rat1%denom*rat2%denom
                temp%num = rat1%num*rat2%denom-rat2%num*rat1%denom
                rat_minus_rat = reduse(temp)
                        return
                        end function rat_minus_rat

                        !分数相乘
                        function rat_times_rat(rat1,rat2)
                        implicit none
                        type(rational) :: rat_times_rat
                        type(rational),intent(in) :: rat1,rat2
                type(rational) :: temp
                temp%denom=rat1%denom*rat2%denom
                temp%num=rat1%num*rat2%num
                rat_times_rat = reduse(temp)
                        return
                        end function rat_times_rat

                        !分数相除
                        function rat_div_rat(rat1,rat2)
                implicit none
                type(rational) :: rat_div_rat
                type(rational),intent(in) :: rat1,rat2
                type(rational) :: temp
                temp%denom=rat1%denom*rat2%num
                temp%num=rat1%num*rat2%denom
                rat_div_rat = reduse(temp)
                        return
                        end function rat_div_rat

                        !输出
                        subroutine output(a)
                        implicit none
                        type(rational),intent(in) :: a
                        if(a%denom/=1) then
                write(*,"(1x,'(',I3,'/',I3,')')") a%num,a%denom
                else
                        write(*,"(1x,I3)") a%num
                        end if
                        return
                        end subroutine output

                        end module rational_util

                        !主程序
                        program main
                        use rational_util
                        implicit none
                        type(rational) :: a,b,c
                        real :: f
        a=rational(1.0,3.0)
                b=rational(2.0,3.0)
                write(*,"(1x,A4)",advance="no") "a="
        call output(a)
                write(*,"(1x,A4)",advance="no") "b="
        call output(b)
                c=a+b
                write(*,"(1x,A4)",advance="no") "a+b="
                call output(c)
                c=a-b
                write(*,"(1x,A4)",advance="no") "a-b="
                call output(c)
        c=a*b
        write(*,"(1x,A4)",advance="no") "a*b="
        call output(c)
                c=a/b
                write(*,"(1x,A4)",advance="no") "a/b="
        call output(c)
                f=c
                write(*,"(f6.2)") f
stop
end program main

image-20220330161116969

习题

  • 1、用interface定义出一个新的虚拟函数area,当调用area只输入一个浮点数时,把它当成是圆的半径值,计算并返回圆面积。当输入两个浮点数时,把它们当成是矩形的两个边长,计算并返回矩形面积。

  • 2、ex1107示例了自定义的时间相加方法,这个程序所能接受的时间单位只有小时及分钟而已。请把它改写成还能接受“秒”这个单位。

  • 3、ex1110中所定义分数类型运算,只定义了加、减、乘、除而已,请试着加入逻辑操作符<,>,<=,>=,==,/=等等。

  • 4、添加一个二维向量类型(x,y),请实现出两个向量的加减法、和实数之间的乘法,以及内积的计算

       ```fortran
       type vector
          real x,y
       end type
       ```

编译器的高级使用

计算机绘图

数值方法

  • 数值方法是fortran语言最主要的应用领域

求解非线性函数

  • 如何求解函数f(x)=0的解,也就是计算函数f(x)的图形和x轴的交点

二分法Bisection

  • 二分法是最简单的解法
    • (一)先猜两个值,使得f(a)*f(b)小于0,也就是f(a)、f(b)必须异号。这样才能保证在a~b之间存在一个c值,使得f(c)=0
    • (二)令c=(a+b)/2,如果f(c)=0就找到了一个解,工作完成
    • (三)f(c)不为0时,如果f(a)、f(c)异号,则以a,c为新的两个猜值来重复步骤二;如果f(b)、 f(c)异号,则以b、c为新的猜值来重复步骤二
module numerical
    implicit none
    real,parameter :: zero=0.00001
contains
    real function bisect(A,B)
        implicit none
        real :: A,B !输入的猜值
        real :: C !用来记录算(A+B)/2
        real :: FA !记录F(A)
        real :: FB  !记录F(B)
        real :: FC  !记录F(C)
        C=(A+B)/2.0
        FC=func(C)
        do while(abs(fc)>zero)
            FA=func(A)
            FB=func(B)
            if(FA*FC<0) then
                B=C
                C=(A+B)/2.0
            else
                A=C
                C=(A+B)/2.0
            end if
            FC=func(c)
        end do
        bisect=C
        return
    end function

    real function func(X)
        implicit none
        real :: X
        func=(X+3)*(X-3)
        return
    end function func
end module numerical

program main
    use numerical
    implicit none
    real :: A,B
    real :: ANS
    do while(.true.)
        write(*,*) "请输入两个猜测值:"
        read(*,*) A,B
        if(func(A)*func(B)<0) exit
        write(*,*) "猜值A,B不正确,请重新输入!"
    end do
    ANS=bisect(A,B)
    write(*,"('x=',F6.3)") ans
    stop
end program main

image-20220401104058309

  • 上面实例只能求解f(x)=(x+3)*(x-3),值为-3和3,但二分法只能求出一个解,不具有求解多个解的方程,且猜值繁琐。
module numerical
    implicit none
    real,parameter :: zero = 0.00001
contains
    real function bisect(A,B,func)
        implicit none
        real :: A,B
        real :: C
        real :: FA,FB,FC
        real,external :: func
        C=(A+B)/2.0
        FC=func(C)
        do while(abs(fc)>zero)
            FA=func(A)
            FB=func(B)
            if(FA*FC<0) then
                B=C
                C=(A+B)/2.0
            else
                A=C
                C=(A+B)/2.0
            end if
            FC=func(C)
        end do
        bisect=C
        return
    end function bisect

    real function f1(X)
        implicit none
        real :: X
        f1=(X+3)*(X-3)
        return
    end function f1

    real function f2(X)
        implicit none
        real :: X
        f2=(X+4)*(X-5)
        return
    end function f2
end module numerical

program main
    use numerical
    implicit none
    real :: A,B
    real :: ans
    do while(.true.)
        write(*,*) "请输入两个猜测值(第一个方程):"
        read(*,*) A,B
        if(f1(A)*f1(B)<0) exit
        write(*,*) "不正确的猜值!"
    end do

    ans=bisect(A,B,f1)
    write(*,"('x=',F6.3)") ans
    do while(.true.)
        write(*,*) "请输入两个猜值(第二个方程):"
        read(*,*) A,B
        if(f2(A)*f2(B)<0) exit
        write(*,*) "不正确的猜值!"
    end do
    ans=bisect(A,B,f2)
    write(*,"('x=',F6.3)") ans
    stop
end program main

image-20220402201000159

割线法Secant

  • 割线法,适合使用图形来解释,它主要是利用线段来逼近结果,过程如下:
    • (一)先选出两个猜测值
    • (二)画一条通过(a,f(a))、(b,f(b))这两点的直线,令这条直线与X轴的交点为c。检查f(c)是否等于0,如果是就找到了一个解。
    • (三)f(c)不为0时,令b、c值为新的两个猜测值a,b,再回到上一个步骤继续。

image-20220403090809996

module numerical
    implicit none
    real,parameter :: zero = 0.00001
contains
    real function secant(a,b,f)
        implicit none
        real :: a,b
        real :: c
        real,external :: f
        real :: fa,fb,fc
        fa=f(a)
        fb=f(b)
        c=a-fa*(b-a)/(fb-fa)
        fc=f(c)
        do while(abs(fc)>zero)
            a=b
            b=c
            fa=f(a)
            fb=f(b)
            c=a-fa*(b-a)/(fb-fa)
            fc=f(c)
        end do
        secant = c
        return
    end function secant

    real function func(x)
        implicit none
        real :: x
        func=sin(x)
        return
    end function func

end module numerical

program main
    use numerical
    implicit none
    real :: a,b
    real :: ans
    write(*,*) "输入两个值:"
    read(*,*) a,b
    ans=secant(a,b,func)
    write(*,"('x=',f8.4)") ans
    stop
end program main

image-20220403090956712

  • 注意:割线法并不一定保证会找到解,也有可能c值会越来越偏离答案,实际使用时,应该随时检查f(c)的值,如果发现f(c)没有向0逼近,表示初始的猜值不能使用,这个时候就没有办法找到正确解。

牛顿法

  • 牛顿法也是利用线段来逼近结果,计算过程如下:
    • (一)先做一个猜值a
    • (二)以f‘(a)为斜率,经过(a,f(a))作一条直线,令这条直线与X轴的交点为b。检查f(b)是否为0,如果是就找到一个解。
    • (三)f(b)不为0时,重新令b为新的猜值a,回到步骤2重复操作

image-20220403091851952

  • 若初始猜值a取的好,f(b)应该会越来越接近0,程序如下
module numerical
implicit none
real,parameter :: zero =0.00001
contains
real function newton(a,f,df)
implicit none
real :: a
real,external :: f  !输入的求值函数
real,external :: df  !f'(x)的函数
real :: b
real :: fb
b=a-f(a)/df(a)
fb=f(b)
do while(abs(f(b))>zero)
a=b
b=a-f(a)/df(a)
fb=f(b)
end do
newton=b
return
end function newton
real function func(x)
implicit none
real :: x
func=sin(x)
return
end function func

real function dfunc(x)
implicit none
real :: x
dfunc=cos(x)
return
end function dfunc
end module numerical

program main
use numerical
implicit none
real :: a
real :: ans
write(*,*) "输入起始初猜值:"
read(*,*) a
ans=newton(a,func,dfunc)
write(*,"('x=',f8.4)") ans
stop
end program main

image-20220403093938079

  • 使用牛顿法跟割线法一样,如果猜值给的不好,会永远无法逼近出结果。实际应用中,应该要检查f(b)是否向0逼近

线性代数

  • 线性代数的数值方法,就是矩阵的应用,二维数组经常被当成矩阵来使用。

矩阵的加、减、乘法

  • 矩阵的加、减只是单纯地把矩阵中相同坐标位置的数字相加、减而已。
real a(m,n),b(m,n),c(m,n)
c=a+b
c=a-b
  • 矩阵的乘法
c=a*b  !这个命令是做c(i,j)=a(i,j)*b(i,j),并不等于矩阵相乘
c=matmul(a,b) !库存函数matmul可以做矩阵相乘

三角矩阵

  • 使用程序来解决这个问题的方法,和用手计算的过程是一样的,同样是把某一行乘上一个系数之后和另外一行相减。做上三角矩阵时,就先把第一行的第一列以下的元素数值都清为0,再把第二行第二列以下的数值都清为0,······如此一直做到第N-1行为止
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,1)
        n=size(matrix,2)
        write(for(2:3),'(I2)') N !完善输出格式
        do i=1,N
            write(*,FMT=for) matrix(i,:)
        end do
        return
    end subroutine output
    !求上三角矩阵的子程序
    subroutine Upper(matrix)
        implicit none
        real :: matrix(:,:)
        integer :: M,N
        integer :: I,J
        real :: E
        M=size(matrix,1)
        N=size(matrix,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 Upper
    !求下三角矩阵的子程序
    subroutine Lower(matrix)
        implicit none
        real :: matrix(:,:)
        integer :: M,N
        real :: I,J,E
        M=size(matrix,1)
        N=size(matrix,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:I)-matrix(I,1:I)*E
            end do
        end do
        return
    end subroutine Lower
end module LinearAlgebra

program main
    use LinearAlgebra
    implicit none
    integer,parameter :: N=3
    real :: A(N,N) = reshape((/1,2,1,3,2,3,2,3,4/),(/N,N/))
    real :: B(N,N)
    write(*,*) "Matrix A:"
    call output(A)
    B=A
    write(*,*) "Upper:"
    call Upper(B)
    call output(B)
    B=A
    write(*,*) "Lower:"
    call Lower(B)
    call output(B)
    stop
end program main

image-20220403193539464

Determinant矩阵的值

  • 求矩阵行列式值的方法,把矩阵转化成上三角或下三角矩阵后,再把对角线上的元素相乘就是其行列式的值。
module LinearAlgebra
    implicit none
contains
    !求矩阵的Determinant值
    real function Determinant(matrix)
        real :: matrix(:,:)
        real,allocatable :: ma(:,:)
        integer :: i,N
        N=size(matrix,1)
        allocate(ma(N,N))
        ma=matrix
        call Upper(ma)
        Determinant=1.0
        do i=1,N
            Determinant = Determinant*ma(i,i)
        end do
        return
    end function Determinant
    subroutine Upper(matrix)
        real :: matrix(:,:)
        integer :: M,N
        integer :: I,J
        real :: E
        M=size(matrix,1)
        N=size(matrix,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 Upper
end module LinearAlgebra

program main
    use LinearAlgebra
    implicit none
    integer,parameter :: N=3
    real :: A(N,N) = reshape((/1,2,1,3,2,3,2,3,4/),(/N,N/))
    write(*,"('det(A)=',f6.2)") Determinant(A)
    stop
end program main

image-20220403205948790

Gauss-Jordan法求联立方程式

image-20220417212809793

module LinearAlgebramodule LinearAlgebra
    implicit none
    contains
    !Guess-jordan法
    subroutine Gauss_Jordan(A,S,ANS)
        implicit none
        real :: A(:,:)
        real :: S(:)
        real :: ANS(:)
        real,allocatable :: B(:,:) !可分配数组B
        integer :: i,N
        N = size(A,1) !行数
        allocate(B(N,N))
        !保存原先矩阵A,S
        B=A
        ANS=S
        call Upper(B,ANS,N) ! 将B化成上三角矩阵
        call Lower(B,ANS,N)  ! 将B化成下三角矩阵
        !求解
        forall(i=1:N)
            ANS(i)=ANS(i)/B(i,i)
        end forall
        return
    end subroutine Gauss_jordan
    !输出等式
    subroutine output(M,S)
        implicit none
        real :: M(:,:),S(:)
        integer :: N,i,j
        N=size(M,1)
        do i=1,N
            write(*,"(1x,f5.2,a1)",advance="NO") M(i,1),'A'
            do j=2,N
                if(M(i,j)<0) then
                    write(*,"('-',f5.2,a1)",advance="NO") -M(i,j),char(64+j)
                else
                    write(*,"('+',f5.2,a1)",advance="NO") M(i,j),char(64+j)
                end if
            end do
            write(*,"('=',f8.4)") S(i)
        end do
        return
    end subroutine output

    !求上三角函数的子程序
    subroutine Upper(M,S,N)
        implicit none
        integer :: N
        real :: M(N,N)
        real :: S(N)
        integer :: i,j
        real :: E
        do i=1,N-1
            do J=i+1,N
                E=M(j,i)/M(i,i)
                M(j,i:N)=M(J,i:N)-M(i,i:N)*E
                S(j)=S(j)-S(i)*E
            end do
        end do
        return
    end subroutine Upper

    !求下三角函数的子程序
    subroutine Lower(M,S,N)
        implicit none
        integer :: N
        real :: M(N,N)
        real :: S(N)
        integer :: i,j
        real :: E
        do i=N,2,-1
            do j=i-1,1,-1
                E=M(i,j)/M(i,i)
                M(j,1:N)=M(j,1:N)-M(i,i:N)*E
                S(j)=S(j)-S(i)*E
            end do
        end do
        return
    end subroutine Lower
end module LinearAlgebra

!0求解联立式
program main
    use LinearAlgebra
    implicit none
    integer,parameter :: N=3 !Size of Matrix
    real :: A(N,N)=reshape((/1,2,3,4,5,6,7,8,8/),(/N,N/))
    real :: S(N) = (/12,15,17/)
    real :: ans(N)
    integer :: i
    write(*,*) "Equation:"
    call output(A,S)
    call Gauss_Jordan(A,S,ANS)
    write(*,*) 'ANS:'
    do i=1,N
        write(*,"(1x,a1,'=',F8.4)") char(64+i),ANS(i)
    end do
    return
end program main

image-20220417224332233

逆矩阵

对角矩阵的运行

积分

梯形法积分

SIMPSON辛普森法积分

插值法与曲线近似

Lagrange Interpolation 多项式插值法

牛顿法Forward Interpolation

最小方差法(Least Square)

曲线近似法(Cubic Spline)

数据结构与算法

排序

冒泡排序法(Bubble Sort)

选择排序法(Selection Sort)

Shell排序法

快速排序法(Quick Sort)

搜索

顺序搜索

二元搜索

散列搜索(Hashing)

堆栈Stack

堆栈的基本范例

堆栈的应用——骑士走棋盘

树状结构

IMSL函数库

线性代数

添加的运算符号

矩阵函数

解线性系统

求解非线性方程

多项式函数

任意函数

求解非线性系统

微积分

积分

多重积分

微分

微分方程

常微分方程(I)

常微分方程(II)

插值与曲线近似

曲线近似

最小方差

附录

fortran库函数

数值运算函数

数学函数

字符函数

数组函数

查询状态函数

二进制运算函数

其他函数

龙格库塔算法

image-20220323092601240

image-20220323101338822

program main
    implicit none
    real :: x,y
    real :: h
    integer :: N,i
    write(*,*) "请输入步长:"
    read(*,*) h
    N = int(3.0*1000)/int(h*1000)
    x=0.0
    y=1.0
    do i=1,N,1
        write(*,*) i
        call rk(x,y,h)
        write(*,*) x,y,"deff=",exp(-0.5*x**2.0)-y
    end do
end program main

function fnf(x,y)
    implicit none
    real :: x,y,fnf
    fnf = -x*y
end function fnf

subroutine rk(x,y,h)
    implicit none
    real :: x,y,h
    real :: k1,k2,k3,k4
    real :: x1,y1
    real,external :: fnf
    x1=x
    y1=y
    k1=h*fnf(x1,y1)
    x1=x+0.5*h
    y1=y+0.5*k1
    k2=h*fnf(x1,y1)
    x1=x+0.5*h
    y1=y+0.5*k2
    k3=h*fnf(x1,y1)
    x1=x+h
    y1=y+k3
    k4=h*fnf(x1,y1)
    y=y+(k1+2.0*k2+2.0*k3+k4)/6.0
    x=x+h
end subroutine rk

image-20220323092957176

image-20220323101048863

program main
    implicit none
    integer(kind=4),parameter :: n=2
    integer(kind=4),parameter :: nm=50,mm=40,m=10
    real(kind=8),parameter :: h=0.001,PI=3.1415926
    real(kind=8) :: x(n)
    integer(kind=4) :: i
    open(20,file="traj.txt")
    x(1)=PI/30.0
    x(2)=2.0
    do i=1,nm*1000
        call rk4(n,h,x)
        if(i>mm*1000) then
            write(20,"(1X,2f10.5)") x
        end if
    end do
    close(20)
end program main

subroutine rk4(n,h,x)
    implicit none
    integer(kind=4) :: n
    real(kind=8) :: h,xx(n),x(n),fx(n)
    real(kind=8) :: kx1(n),kx2(n),kx3(n),kx4(n)
    xx=x
    call fun(xx,fx,n)
    kx1=h*fx
    xx=x+0.5*kx1
    call fun(xx,fx,n)
    kx2=h*fx
    xx=x+0.5*kx2
    call fun(xx,fx,n)
    kx3=h*fx
    xx=x+kx3
    call fun(xx,fx,n)
    kx4=h*fx
    x=x+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
    return
end subroutine rk4

subroutine fun(x,fx,n)
    implicit none
    integer(kind=4) :: n,i,j
    real(kind=8) :: x(n),fx(n)
    real(kind=8) :: omega=1.0d0
    fx(1) = x(2)
    fx(2) = -omega*omega*x(1)
    return
end subroutine fun

image-20220323100934068

image-20220323101417151

program main
    implicit none
    integer(kind=4),parameter :: n=2,nm=20
    real(kind=8),parameter :: h=0.001,PI=3.1415926
    real(kind=8) :: x(n),x0(n)
    integer(kind=4) :: i
    open(20,file="traj1.txt")
    x0(1)=-2*PI-PI/7.0
    do while(x0(1)<2*PI+PI/7.0)
        x0(1)=x0(1)+pi/7.0
        x0(2)=-2.5
        do while(x0(2)<2.5)
            x0(2)=x0(2)+0.5
            x=x0
            do i=1,nm*1000
                call rk4(n,h,x)
                if(mod(i,10)==0) then
                    if(abs(x(1))<2.0*PI) then
                        write(20,"(1x,2f10.5)") x(1),x(2)
                    end if
                end if
            end do
        end do
    end do
    close(20)
end program main

subroutine rk4(n,h,x)
    implicit none
    integer(kind=4) :: n
    real(kind=8) :: h,xx(n),x(n),fx(n)
    real(kind=8) :: kx1(n),kx2(n),kx3(n),kx4(n)
    xx=x
    call fun(xx,fx,n)
    kx1=h*fx
    xx=x+0.5*kx1
    call fun(xx,fx,n)
    kx2=h*fx
    xx=x+0.5*kx2
    call fun(xx,fx,n)
    kx3=h*fx
    xx=x+kx3
    call fun(xx,fx,n)
    kx4=h*fx
    x=x+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
    return
end subroutine rk4

subroutine fun(x,fx,n)
    implicit none
    integer(kind=4) :: n,i,j
    real(kind=8) :: x(n),fx(n)
    real(kind=8) :: omega=1.0d0
    fx(1) = x(2)
    fx(2) = -omega*omega*sin(x(1))
    return
end subroutine fun

image-20220323105007772


文章作者: rep-rebirth
版权声明: 本博客所有文章除特別声明外,均采用 CC BY 4.0 许可协议。转载请注明来源 rep-rebirth !
评论
评论
  目录