If
写法
if (条件) then
满足条件时执行的部分
else
不满足条件时执行的部分
end if
不满足条件时不执行,则
if (条件) then
满足条件时执行的部分
end if
再精简
if (条件) 满足条件时执行的部分
判断条件增加
if (条件1) then
执行部分1
else if (条件2) then
执行部分2
else if (条件3) then
执行部分3
else
以上条件皆不满足时执行的部分
end if
条件式中不等式的不等号写法
写法 | 意味 |
---|---|
< | 小于 |
<= | 小于等于 |
== | 等于 |
>= | 大于等于 |
> | 大于 |
/= | 不等于 |
用例
素数筛(埃拉托斯特尼筛法)
找出[2, n]之间的素数
第一步:写n个0
第二步:将2的倍数(除2自身)全部标记为1
第三步:将3的倍数(除3自身)全部标记为1
......
第n步:将n的倍数(除n自身)全部标记为1
剩下还是0的即为素数
program abc
implicit none
integer :: i, j, k = 1
integer, parameter :: n = 2000
integer, dimension(n) :: a, b
a = 0
do i = 2, n
if (a(i) == 0) then
b(k) = i
k = k + 1
do j = i * i, n, i
a(j) = 1
end do
end if
end do
print '(1x, 15i6)', b(1: k - 1)
print *,'press any key to continue'
read *,
stop
end
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47
53 59 61 67 71 73 79 83 89 97 101 103 107 109 113
127 131 137 139 149 151 157 163 167 173 179 181 191 193 197
199 211 223 227 229 233 239 241 251 257 263 269 271 277 281
283 293 307 311 313 317 331 337 347 349 353 359 367 373 379
383 389 397 401 409 419 421 431 433 439 443 449 457 461 463
467 479 487 491 499 503 509 521 523 541 547 557 563 569 571
577 587 593 599 601 607 613 617 619 631 641 643 647 653 659
661 673 677 683 691 701 709 719 727 733 739 743 751 757 761
769 773 787 797 809 811 821 823 827 829 839 853 857 859 863
877 881 883 887 907 911 919 929 937 941 947 953 967 971 977
983 991 997 1009 1013 1019 1021 1031 1033 1039 1049 1051 1061 1063 1069
1087 1091 1093 1097 1103 1109 1117 1123 1129 1151 1153 1163 1171 1181 1187
1193 1201 1213 1217 1223 1229 1231 1237 1249 1259 1277 1279 1283 1289 1291
1297 1301 1303 1307 1319 1321 1327 1361 1367 1373 1381 1399 1409 1423 1427
1429 1433 1439 1447 1451 1453 1459 1471 1481 1483 1487 1489 1493 1499 1511
1523 1531 1543 1549 1553 1559 1567 1571 1579 1583 1597 1601 1607 1609 1613
1619 1621 1627 1637 1657 1663 1667 1669 1693 1697 1699 1709 1721 1723 1733
1741 1747 1753 1759 1777 1783 1787 1789 1801 1811 1823 1831 1847 1861 1867
1871 1873 1877 1879 1889 1901 1907 1913 1931 1933 1949 1951 1973 1979 1987
1993 1997 1999
press any key to continue
利用迭代法求数组的最大 / 小值
先将数组中的第一个数同时指定为最大值和最小值
对第2个数,如果其大 / 小于前次指定的最大 / 小值,则将其指定为新最大 / 小值
......
对第n个数,如果其大 / 小于前次指定的最大 / 小值,则将其指定为新最大 / 小值
比较完成,输出最新一轮比较出的最大 / 小值
program abc
implicit none
integer :: i, n
integer, parameter :: nn = 2000
real, dimension(nn) :: a
real :: maax, miin
print *, 'n = '
read *, n
do i = 1, n
print '(1x, 2ha(, i1, 2h)= )', i
read *, a(i)
end do
maax = a(1)
miin = a(1)
do i = 2, n
if(a(i) > maax) maax = a(i)
if(a(i) < miin) miin = a(i)
end do
print *, 'maximum', maax
print *, 'minimum', miin
print *,'press any key to continue'
read *,
stop
end
n =
6
a(1)=
54515.231
a(2)=
45345.11
a(3)=
5143515012
a(4)=
1231.011
a(5)=
1513.01515
a(6)=
1631194.001144
maximum 5.14351514E+09
minimum 1231.01099
press any key to continue
交换法(冒泡排序)
对第1到第n - 1项
将每一项与后项比较,如果该项较后项大,则将两项交换位置
此时最后一项应是最大的,故对第1到第n - 2 项
将每一项与后项比较,如果该项较后项大,则将两项交换位置
......
对第1项
与后项比较,大则换位
排序完成
program abc
implicit none
integer :: i, n, ter
integer, parameter :: nn = 2000
real, dimension(nn) :: a
real :: bigger
print *, 'n = '
read *, n
print *, 'a = '
read *, a(1: n)
do ter = n - 1, 1, -1
do i = 1, ter
if (a(i) > a(i + 1)) then
bigger= a(i)
a(i) = a(i + 1)
a(i + 1) = bigger
end if
end do
end do
print *, a(1: n)
print *,'press any key to continue'
read *,
stop
end
n =
6
a =
5115
6161
52152
121
3513
1315
121.000000 1315.00000 3513.00000 5115.00000 6161.00000 52152.0000
press any key to continue
JUMP
Exit
出现于循环节(do)之中,用来退出循环
do 循环条件
if (循环执行到某次总会满足的一个条件) then
循环执行的部分
exit
end if
end do
多重循环的场合
命名1: do
命名2: do
if (...) exit 命名1
if (...) exit 命名2
end do
end do
不加指定退出哪一个循环节则只退出最内侧的那个
Go to
跳到指定序号的那一行
if (...) go to 10
...
10 print *, '10'
Cycle
出现于循环节(do)之中,用来初始化循环
do ...
if (...) cycle
循环执行的部分
end do
Do while
写法
满足条件时一直重复
do while (条件)
循环执行的部分
end do
用例
二分查找算法(整数)
确定搜索的上下限
求上下限的平均数并取整
如果要查找的数正好是平均数,则过程结束
如果要查找的数在平均数与上 / 下限之间,则将平均数指定为新的最下 / 上限
......
由于上限 - 下限 = 2时再操作一次平均数必定是最终结果
故循环条件是上限 > 下限 + 1
例如求分段函数的取值
$$\begin{array}\ y=\begin{cases}300\quad x\in (0,2]\\\ 550\quad\ x\in \ (2,3]\\\ 880\quad\ x\in \left(3,5\right]\\\ 1100\quad\ x\in \left(5,7\right]\\\ 1600+500\times \left(\frac{x-7.000001}{5}\right)\quad\ x\in \left(7,\infty \right)\ \end{cases}\end{array}\begin{array}\ x\in N*\end{array}$$
program abc
implicit none
integer :: d, u, m, t
integer, parameter :: n = 4
integer, parameter :: q = 3
integer, dimension(n) :: x, y
integer, dimension(q) :: a, b
x = (/2, 3, 5, 7/)
y = (/300, 550, 880, 1100/)
do t = 1, 3
print '(1x, a2, i2, a2)', 'x(', t, ')='
read *, a(t)
if(a(t) <= x(1)) then
b(t) = y(1)
else if(a(t) > x(n)) then
b(t) = 1600 + 500 * int((a(t) - 7.000001) / 5)
else
d = 1
u = n
do while(d + 1 < u)
m = (u + d) / 2
if (a(t) == x(m)) then
u = m
exit
else if(a(t) < x(m)) then
u = m
else
d = m
end if
end do
b(t) = y(u)
end if
end do
do t = 1, 3
print '(1x, a2, i2, a2, i5)', 'y(', t, ')=', b(t)
end do
print *,'press any key to continue'
read *,
stop
end
x( 1)=
1
x( 2)=
6
x( 3)=
15
y( 1)= 300
y( 2)= 1100
y( 3)= 2100
press any key to continue
当然,将输入的值与有界部分每段的上界比较,然后找出该值所在的区间的方法也行得通
program abc
implicit none
integer :: i, t
integer, parameter :: n = 4
integer, parameter :: q = 3
integer, dimension(n) :: x, y
integer, dimension(q) :: a, b
x = (/2, 3, 5, 7/)
y = (/300, 550, 880, 1100/)
do t = 1, 3
print '(1x, a2, i2, a2)', 'x(', t, ')='
read *, a(t)
if(a(t) > x(n)) then
b(t) = 1600 + 500 * int((a(t) - 7.000001) / 5)
else
do i = 1, 4
if(a(t) <= x(i)) then
b(t) = y(i)
exit
end if
end do
end if
end do
do t = 1, 3
print '(1x, a2, i2, a2, i5)', 'y(', t, ')=', b(t)
end do
print *,'press any key to continue'
read *,
stop
end
x( 1)=
1
x( 2)=
6
x( 3)=
13
y( 1)= 300
y( 2)= 1100
y( 3)= 2100
press any key to continue
但是设想极端情况,比如要查找的数在有界部分的很后端,那么用逐个对比的方法就会比不断二分来得慢了
副程序
程序较长时,分主次更合适。
副程序一般用来摞subroutine(子程序)和function(函数)。
副程序是独立于主程序的,所以咱一般把它写在
end program
的后一行。
也是因此,需要在主程序的宣言中宣言副程序要用到的各参数,以便调动副程序中的参数到主程序中。
如果是想要在副程序中再创建副程序,步骤与上述相同,也就是写在
end subroutine
或者
end function
的后行。
你会注意到这个时候同一个宣言要在主副里头都要写一道,很麻烦。如果想要省略一道,可以用
contains
把external(外部的)subroutine或function变成internal(内部的)subroutine或function。
具体参照以下内容。
写法
Subroutine / function ( external )
在主程序的宣言后追加
interface
subroutine 子程序名(用到的数)
implicit none
型, intent (in) :: 输入用数扔这
型, intent (out) :: 输出用数扔这
型, intent (inout) :: 输入输出两用数扔这
end subroutine 子程序名
end interface
interface
function 函数名(用到的数)
implicit none
型, intent (in) :: 输入用数扔这
型, intent (out) :: 输出用数扔这
型, intent (inout) :: 输入输出两用数扔这
end function 函数名
end interface
在主程序结束以后写
subroutine 子程序名(用到的数)
implicit none
型, intent (in) :: 列出输入用数
型, intent (out) :: 列出输出用数
型, intent (inout) :: 列出输入输出两用数
对subroutine中要用到的其他进行型宣言
处理
结果代入输出用数或输入输出两用数
return
end subroutine 子程序名
function 函数名(用到的数)
implicit none
型, intent (in) :: 列出输入用数
型, intent (out) :: 列出输出用数
型, intent (inout) :: 列出输入输出两用数
对function中要用到的其他进行型宣言
处理
函数名 = 处理结果
return
end function 函数名
Subroutine / function ( internal )
主程序stop以后用写包含的若干副程序,然后end掉主程序
program 主程序名
implicit none
型宣言
处理
stop
contains !包含
subroutine 子程序名(用到的数)
implicit none
型, intent (in) :: 列出输入用数
型, intent (out) :: 列出输出用数
型, intent (inout) :: 列出输入输出两用数
对subroutine中要用到的其他进行型宣言
处理
结果代入输出用数或输入输出两用数
return
end subroutine 子程序名
function 函数名(用到的数)
implicit none
型, intent (in) :: 列出输入用数
型, intent (out) :: 列出输出用数
型, intent (inout) :: 列出输入输出两用数
对function中要用到的其他进行型宣言
处理
函数名 = 处理结果
return
end function 函数名
end program
用例
分数的加法
呼出子程序时,用以下表达
call 子程序名(用到的数)
如计算
$$\frac{a}{b} + \frac{c}{d} = \frac{e}{f}$$
program abc
implicit none
integer :: a, b, c, d, e, f
interface
subroutine fs(a, b, c, d, e, f) !外部子程序与主程序相互独立,可以重复使用主程序中已经使用到的字母
implicit none
integer, intent(in) :: a, b, c, d
integer, intent(out) :: e, f
end subroutine
end interface
print *, 'a, b, c, d = '
read *, a, b, c, d
call fs(a, b, c, d, e, f) !呼出子程序
print *, 'e = ', e
print *, 'f = ', f
print *,'press any key to continue'
read *,
stop
end
!以下是通分运算再约分的子程序
subroutine fs(a, b, c, d, e, f)
implicit none
integer, intent(in) :: a, b, c, d
integer, intent(out) :: e, f
integer :: g
interface
subroutine gcd(a, b, c)
implicit none
integer, intent(in) :: a, b
integer, intent(out) :: c
end subroutine
end interface
e = a * d + b * c
f = b * d
call gcd(e, f, g) !子程序中呼出子程序
e = e / g
f = f / g
return
end subroutine
!以下是辗转相除法求最大公约数的子程序
subroutine gcd(a, b, c)
implicit none
integer, intent(in) :: a, b
integer, intent(out) :: c
integer :: d, e, f
d = a
e = b
do
f = mod(d, e)
if(f == 0) then
c = e
exit
else
d = e
e = f
end if
end do
return
end subroutine
a, b, c, d =
64641,65465,622,8485
e = 117839623
f = 111094105
press any key to continue
向量的夹角
根据函数的定义,函数名在此是唯一的输出值,故在宣言时只用intent(in)
,函数名不用宣言为intent(out)
另外函数的呼出也相对简单,形如笔算的f(x)型就可以
如计算3维向量a, b的夹角
program abc
implicit none
real, dimension(3) :: a, b
real :: s, t
interface
function dp(a, b, n)
implicit none
real :: dp
real, dimension(*), intent(in) :: a, b
integer, intent(in) :: n
end function
end interface
print *, 'Enter a 3 dimensional vector a'
read *, a
print *, 'Enter a 3 dimensional vector b'
read *, b
s = dp(a, b, 3) / sqrt(dp(a, a, 3) * dp(b, b, 3)) !直接把定义好的函数代入运算
t = acos(s) * 180 / 3.1415926535
print *, 'The angle between a and b is', t, 'degree(s).'
print *,'press any key to continue'
read *,
stop
end
!以下时n维向量的内积计算函数
function dp(a, b, n)
implicit none
real :: dp
real, dimension(*), intent(in) :: a, b
integer, intent(in) :: n
integer :: i
real :: s
s = 0.0
do i = 1, n
s = a(i) * b(i) + s
end do
dp = s
return
end function
Enter a 3 dimensional vector a
44,55,66
Enter a 3 dimensional vector b
11,77,99
The angle between a and b is 22.1832886 degree(s).
press any key to continue
二项式系数
二分法求根
确定初始的x1,x2,使f(x1)<0,f(x2)>0
求x1,x2的平均数x
若f(x)<0,则将x作为新的x1,否则将其作为新的x2重复操作
......
达到指定回数时,停止操作
如求
$$\cos\ x-x=0$$
的根
program abc
implicit none
real :: x, y, x1, x2
integer :: r
x1 = 0.0
x2 = 1.0
do r = 1, 10000
x = (x1 + x2) / 2.0
y = f(x)
if(y < 0) then
x1 = x
else
x2 = x
end if
end do
print *, 'x = ', x
print *,'press any key to continue'
read *,
stop
!主程序中包含以下函数
contains
function f(x)
implicit none
real :: f
real, intent(in) :: x
f = x - cos(x)
return
end function
end
x = 0.739085078
press any key to continue
牛顿法求根
先看一个简单的例子
$$x^{2}-99=0$$
这就相当于是求面积是99的正方形的边长
于是可以有以下思路
首先找到一组长宽,使其乘积等于给定的正方形面积
接下来把长宽的平均数指定为新长,面积除以新长指定为新宽
......
直到长宽之差几乎为0,结束过程
program abc
implicit none
real :: l, w, x2
real, parameter :: p = 0.0000001
print *, 'square = '
read *, x2
l = 1.0
w = x2
do while(abs((l - w) / l) > p)
l = (l + w) * 0.5
w = x2 / l
end do
print *, 'root = ', l
print *,'press any key to continue'
read *,
stop
end
square =
99
root = 9.94987488
press any key to continue
如果要求任意函数与x轴的交点,如
$$\cos\ x-x=0$$
则使用牛顿法
通过不断作切线,直到f(x)几乎为0
$$x_{n+1} = x_{n} - \frac{f(x_{n})}{f^{\prime}(x_{n})}$$
逼近思想与与上述长方化正方是一致的
数值求导可以利用导数的定义求,为避免麻烦此处直接使用笔算的表达式
program abc
implicit none
real :: x, y, dy, p = 0.000001
integer :: r
x = 1.0
do r = 1, 100
y = f(x)
dy = df(x)
if(abs(dy) < p) then
print *, 'Stationary Point'
stop
end if
x = x - y / dy
if(abs(y) < p) exit !y值足够小
if(abs(y / dy) < p) exit !x(n+1)和x(n)的距离足够小
end do
print *, 'x = ', x
print *,'press any key to continue'
read *,
stop
contains
function f(x)
implicit none
real :: f
real, intent(in) :: x
f = cos(x) - x
return
end function
function df(x)
implicit none
real :: df
real, intent(in) :: x
df = -sin(x) - 1
return
end function
end
x = 0.739085138
press any key to continue
辛普森积分法
$$\int_{a}^{b} f(x) d x \approx \frac{h}{3}\left[f\left(x_{0}\right)+2 \sum_{j=1}^{n / 2-1} f\left(x_{2 j}\right)+4 \sum_{j=1}^{n / 2} f\left(x_{2 j-1}\right)+f\left(x_{n}\right)\right]$$
其中
$$h=(b-a) / n \\x_{i}=a+i h(i=0, \ldots, n)$$
program abc
implicit none
integer :: j, n
real :: a, b, h, s1, s2, s
interface
function f(x)
implicit none
real :: f
real, intent(in) :: x
end function
end interface
print *, 'a, b = '
read *, a, b
print *, 'n = '
read *, n
h = (b - a) / real(n)
s1 = 0.0
do j = 2, n - 2, 2
s1 = s1 + f(a + j * h)
end do
s2 = 0.0
do j = 1, n - 1, 2
s2 = s2 + f(a + j * h)
end do
s = (f(a) + 2.0 * s1 + 4.0 * s2 + f(b)) * h / 3.0
print *, 's = ', s
print *,'press any key to continue'
read *,
stop
end
function f(x)
implicit none
real :: f
real, intent(in) :: x
f = sin(x)**cos(x)
return
end function
a, b =
0,1
n =
1000000
s = 0.501324713
press any key to continue