518 lines
28 KiB
Fortran
518 lines
28 KiB
Fortran
! ==============================================================================
|
||
! 通过 gfortran ./test.f90 -o ./run && ./run 运行
|
||
! 程序名：
|
||
! 习题
|
||
! 目的:
|
||
! 练习用
|
||
! 修订记录:
|
||
! 日期 编程者 改动描述
|
||
! =================== ============= =====================================
|
||
! 2021-04-26 16:53:03 Sola 测验8-1 7 测试数组调用
|
||
! 2021-04-30 18:03:39 Sola 习题8-7 输出数组中正负数和0的个数
|
||
! 2021-04-30 18:49:30 Sola 习题8-8 可分配数组输入数据定义并计算行列和
|
||
! 2021-04-30 19:26:59 Sola 习题8-10 没看懂啥意思, 跳过
|
||
! 2021-04-30 19:27:16 Sola 习题8-11 在赋值中修改数组范围，太长了，跳过
|
||
! 2021-04-30 19:30:55 Sola 习题8-12 DO循环和WHERE结构比较
|
||
! 2021-04-30 20:10:23 Sola 习题8-13 计算年平均温度
|
||
! 2021-04-30 20:46:31 Sola 习题8-14 矩阵乘法 8-15 8-16
|
||
! 2021-05-01 00:32:18 Sola 习题8-17 相对极大(鞍点)
|
||
! 2021-05-01 02:01:58 Sola 习题8-18 金属盘温度分布
|
||
! 程序结构:
|
||
!
|
||
! ==============================================================================
|
||
! 模块：
|
||
module Module_8
|
||
implicit none
|
||
! 数据字典
|
||
! 声明常数
|
||
REAL, PARAMETER :: PI=3.14159265 ! PI值
|
||
REAL, PARAMETER :: e=2.718281828459 ! 自然对数
|
||
INTEGER, PARAMETER :: arrayLength=20 ! 数组基准长度
|
||
! 声明变量
|
||
! 创建显式接口
|
||
contains
|
||
! subroutine SubName(varName1,varName2)
|
||
! implicit none
|
||
! ! 数据字典
|
||
! end subroutine SubName
|
||
end module Module_8
|
||
! ==============================================================================
|
||
! 主程序：
|
||
program ProName
|
||
implicit none
|
||
! call Exercises8_1_7
|
||
! call Exercises8_7
|
||
! call Exercises8_8
|
||
! call Exercises8_9
|
||
! call Exercises8_10
|
||
! call Exercises8_11
|
||
! call Exercises8_12
|
||
! call Exercises8_13
|
||
! call Exercises8_14
|
||
! call Exercises8_17
|
||
call Exercises8_18
|
||
! 数据字典
|
||
! 声明常量
|
||
! 声明变量
|
||
! 变量初始化
|
||
! 数据输入
|
||
! 运算过程
|
||
! 结果输出
|
||
end program ProName
|
||
! ==============================================================================
|
||
! 子程序
|
||
! subroutine SubName(varName1,varName2)
|
||
! use MouName
|
||
! implicit none
|
||
! Type, intent(inout) :: varName
|
||
! end subroutine SubName
|
||
! ==============================================================================
|
||
! 函数
|
||
! function FunName(varName1,varName2)
|
||
! use MouName
|
||
! implicit none
|
||
! end function FunName
|
||
! ==============================================================================
|
||
! 测试8-1.7 使用数组从数组中调取元素
|
||
subroutine Exercises8_1_7
|
||
implicit none
|
||
integer, dimension(4) :: list1=(/1,4,2,2/) ! 定义数组1
|
||
integer, dimension(3) :: list2=(/1,2,3/) ! 定义数组2
|
||
integer, dimension(5, 5) :: array ! 定义数组3
|
||
integer :: i, j ! 循环参数
|
||
do i = 1, 5
|
||
do j = 1, 5
|
||
array(i, j) = i + 10*j ! 数组3赋值
|
||
end do
|
||
end do
|
||
write (*,'(1X, 4I4)') array(list1, list2) ! 测试通过数组调用数组
|
||
end subroutine Exercises8_1_7
|
||
! 习题8-7 统计数组中正负数和0的个数
|
||
subroutine Exercises8_7
|
||
implicit none
|
||
real, dimension(-50:50, 0:100) :: values ! 定义随机数组
|
||
integer :: i, j ! 循环参数
|
||
integer :: numPositive, numNegative, numZero ! 统计参数
|
||
! 初始化变量
|
||
numPositive = 0
|
||
numNegative = 0
|
||
numZero = 0
|
||
! FORALL方式随机赋值(好像不能这么用, 所以注释)
|
||
! forall (i=-50:50, j=0:100)
|
||
! call random_seed()
|
||
! call random_number(values(i, j))
|
||
! end forall
|
||
! 嵌套循环赋值
|
||
do i = -50, 50
|
||
do j = 0, 100
|
||
call random_seed() ! 根据日期, 时间获取随机数种子
|
||
call random_number(values(i, j)) ! 调用随机数子程序给数组赋值
|
||
end do
|
||
end do
|
||
values = values * 200. - 100. ! 调整数组范围
|
||
! 循环计算正负数和0的个数
|
||
do i = -50, 50
|
||
do j = 0, 100
|
||
if (values(i, j) > 0) then
|
||
numPositive = numPositive + 1
|
||
elseif (values(i, j) < 0) then
|
||
numNegative = numNegative + 1
|
||
else
|
||
numZero = numZero + 1
|
||
end if
|
||
end do
|
||
end do
|
||
! 输出结果
|
||
write(*,1) ' 测试组', numPositive, numNegative, numZero
|
||
write(*,1) ' 对照组', count(values > 0.), count(values < 0.), count(values == 0.)
|
||
1 format(1X, A/&
|
||
&1X, '数组中有正数 ', I5, '个'/&
|
||
&1X, '数组中有负数 ', I5, '个'/&
|
||
&1X, '数组中有0 ', I5, '个')
|
||
end subroutine Exercises8_7
|
||
! 习题8-8 可分配数组输入数据定义并计算行列和, 注意IO操作的返回值参数和数组分配不同
|
||
subroutine Exercises8_8
|
||
implicit none
|
||
real, dimension(:, :), allocatable :: array ! 定义可分配数组
|
||
integer :: errorLevel ! 错误码
|
||
integer :: a, b ! 数组范围
|
||
integer :: i ! 循环参数
|
||
character(len=20) :: inforError ! 错误信息
|
||
! 打开一个临时文件
|
||
open(unit=1, status='scratch', iostat=errorLevel, iomsg=inforError)
|
||
if (errorLevel /= 0) stop inforError ! 检测打开文件是否成功
|
||
! 输入参数 8-8
|
||
! write(1, *) 2, 4
|
||
! write(1, *) -24., -1121., 812.1, 11.1
|
||
! write(1, *) 35.6, 8.1E3, 135.23, -17.3
|
||
! 输入参数 8-9
|
||
write(1, *) 4, 5
|
||
write(1, *) 33., -12., 16., 0.5, -1.9
|
||
write(1, *) -6., -14., 3.5, 11., 2.1
|
||
write(1, *) 4.4, 1.1, -7.1, 9.3, -16.1
|
||
write(1, *) 0.3, 6.2, -9.9, -12., 6.8
|
||
! write(1, *) , , , ,
|
||
rewind(unit=1) ! 回到文件开头
|
||
read(1,*) a, b ! 读取数组范围
|
||
! 检查数组范围
|
||
if (a > 100) a = 100
|
||
if (b > 100) b = 100
|
||
allocate(array(a, b), stat=errorLevel, errmsg=inforError) ! 分配数组
|
||
if (errorLevel /= 0) stop inforError ! 检测数组分配是否成功
|
||
do i = 1, a
|
||
read(1, *) array(i, :) ! 数组赋值
|
||
end do
|
||
! 输出统计结果
|
||
1 format(1X, 'Sum of ', A, ' ', I3, ' = ', F10.3)
|
||
do i = 1, a
|
||
write(*, 1) 'row', i, sum(array(i, :))
|
||
end do
|
||
do i = 1, b
|
||
write(*, 1) 'col', i, sum(array(:, i))
|
||
end do
|
||
if (allocated(array)) deallocate(array, stat=errorLevel)
|
||
if (errorLevel /= 0) print *, "array: Deallocation request denied"
|
||
close(unit=1)
|
||
end subroutine Exercises8_8
|
||
! 习题8-10 没看懂啥意思, 跳过
|
||
! 习题8-11 在赋值中修改数组范围，太长了，跳过, 写了个简化版的，利用F03的特性实现
|
||
! 不过有一点需要注意，好像不可以用allocate语句多次定义，有点理解无缝是什么意思了
|
||
subroutine Exercises8_11
|
||
implicit none
|
||
real, dimension(:), allocatable :: array ! 定义可分配数组
|
||
integer :: i ! 循环变量
|
||
integer :: errorLevel ! 错误码
|
||
character(len=20) :: inforError ! 错误信息
|
||
real :: temp ! 临时变量
|
||
allocate(array(1), stat=errorLevel, errmsg=inforError) ! 分配数组
|
||
if (errorLevel /= 0) stop inforError ! 检测数组分配状态
|
||
call random_seed() ! 根据日期时间产生随机数种子
|
||
call random_number(array(1)) ! 用随机数给数组赋值
|
||
do i = 2, 50 ! 循环产生随机数, 并给数组赋值
|
||
call random_seed()
|
||
call random_number(temp)
|
||
array = (/array, temp/) ! 利用F03特性扩展数组
|
||
end do
|
||
write(*,'(1X, 10F6.3, " ")') array ! 10个一行输出数组内容
|
||
if (allocated(array)) deallocate(array, stat=errorLevel)
|
||
if (errorLevel /= 0) print *, "array: Deallocation request denied"
|
||
end subroutine Exercises8_11
|
||
! 习题8-12 DO循环和WHERE结构比较
|
||
subroutine Exercises8_12
|
||
implicit none
|
||
integer, parameter :: arrayLength=10 ! 设置数组基本长度
|
||
real, dimension(100*arrayLength, arrayLength, 3*arrayLength) :: arr, arrayTemp ! 定义数组
|
||
integer :: i, j, k ! 循环参数
|
||
! Where结构中也不可调用子程序
|
||
! real :: temp
|
||
! arr = 0.
|
||
! where(arr == 0.)
|
||
! call random_seed()
|
||
! call random_number(temp)
|
||
! arr = temp
|
||
! end where
|
||
do i = 1, 100*arrayLength
|
||
do j = 1, arrayLength
|
||
do k = 1, 3*arrayLength
|
||
call random_seed()
|
||
call random_number(arr(i, j, k)) ! 对数组进行循环赋值
|
||
end do
|
||
end do
|
||
end do
|
||
arr = arr*1300 ! 调整数组范围
|
||
arrayTemp = arr ! 拷贝数组
|
||
write(*, *) 'Do循环判断开始, 数组总和为', sum(arr) ! 参照
|
||
do i = 1, 100*arrayLength
|
||
do j = 1, arrayLength
|
||
do k = 1, 3*arrayLength
|
||
if (arr(i, j, k) > 1000.) arr(i, j, k) = 1000. ! DO判断
|
||
end do
|
||
end do
|
||
end do
|
||
write(*, *) 'Do循环判断结束, 数组总和为', sum(arr) ! 输出结果
|
||
write(*, *) 'Where结构开始' ! 参照
|
||
where(arrayTemp > 1000.)
|
||
arrayTemp = 1000. ! WHERE判断, Where结构真的方便...
|
||
end where
|
||
write(*, *) 'Where结构结束, 数组总和为', sum(arrayTemp) ! 输出结果
|
||
end subroutine Exercises8_12
|
||
! 习题8-13 计算年平均温度
|
||
subroutine Exercises8_13
|
||
implicit none
|
||
real, dimension(6, 6) :: array ! 定义数组
|
||
integer :: i ! 循环参数
|
||
integer :: errorLevel ! 错误码
|
||
character(len=20) :: inforError ! 错误信息
|
||
open(unit=1, status='scratch', iostat=errorLevel, iomsg=inforError) ! 打开临时文件
|
||
if (errorLevel /= 0) stop inforError ! 判断打开成功与否
|
||
write(1, 1) 68.2, 72.1, 72.5, 74.1, 74.4, 74.2, & ! 输入统计数据
|
||
&69.4, 71.1, 71.1, 71.9, 73.1, 73.6, &
|
||
&68.9, 70.5, 70.9, 71.5, 72.8, 73.0, &
|
||
&68.6, 69.9, 70.4, 70.8, 71.5, 72.2, &
|
||
&68.1, 69.3, 69.8, 70.2, 70.9, 71.2, &
|
||
&68.3, 68.8, 69.6, 70.0, 70.5, 70.9
|
||
1 format(6F5.1) ! 设定统计数据输入格式
|
||
rewind(unit=1) ! 回到文件开头
|
||
read(1,*) (array(i, :), i = 1, 6) ! 给数组赋值
|
||
2 format(1X, F4.1, A8, ' 的平均气温为 ', F4.1, ' ℃') ! 设定结果输出格式
|
||
do i = 1, 6
|
||
write(*, 2) 29.5 + real(i)*0.5, '°N lat ', sum(array(i, :))/6. ! 输出纬度平均气温
|
||
end do
|
||
do i = 1, 6
|
||
write(*, 2) 89.5 + real(i)*0.5, '°W long', sum(array(:, i))/6. ! 输出经度平均气温
|
||
end do
|
||
write(*, 3) '所有统计地点年平均温度为 ', sum(array)/36., ' ℃' ! 输出所有地区平均气温
|
||
3 format(1X, A, F4.1, A) ! 设定结果输出格式
|
||
close(unit=1) ! 关闭文件, 此句可忽略
|
||
end subroutine Exercises8_13
|
||
! 习题8-14 矩阵乘法, 注意, F03虽然可以自动分配大小, 但是还是得先给他分配一个大小
|
||
subroutine Exercises8_14
|
||
implicit none
|
||
! real, dimension(:,:), allocatable :: MatrixMultplication
|
||
real, dimension(:, :), allocatable :: array1, array2, array3
|
||
integer :: i, j, k
|
||
integer :: a, b
|
||
integer :: errorLevel ! 错误码
|
||
character(len=20) :: inforError ! 错误信息
|
||
open(unit=1, status='scratch', iostat=errorLevel, iomsg=inforError) ! 打开临时文件
|
||
if (errorLevel /= 0) stop inforError ! 判断打开成功与否
|
||
open(unit=2, status='scratch', iostat=errorLevel, iomsg=inforError) ! 打开临时文件
|
||
if (errorLevel /= 0) stop inforError ! 判断打开成功与否
|
||
! 习题8-14使用
|
||
! 匹配的, 输入矩阵1
|
||
! write (1, '(A)') " 2, 2",&
|
||
! &" 3.0, -1.0",&
|
||
! &" 1.0, 2.0"
|
||
! ! 不匹配的, 输入矩阵1
|
||
! ! write (1, '(A)') " 2, 2, 3",&
|
||
! ! &" 3.0, -1.0, 1.0",&
|
||
! ! &" 1.0, 2.0, 1.0"
|
||
! ! 输入矩阵2
|
||
! write (2, '(A)') " 2, 2",&
|
||
! &" 1.0, 4.0",&
|
||
! &" 2.0, -3.0"
|
||
! 习题8-15使用
|
||
write (1, '(A)') " 2, 4",&
|
||
&" 1.0, -5.0, 4.0, 2.0",&
|
||
&" -6.0, -4.0, 2.0, 2.0"
|
||
write (2, '(A)') " 4, 3",&
|
||
&" 1.0, -2.0, -1.0",&
|
||
&" 2.0, 3.0, 4.0",&
|
||
&" 0.0, -1.0, 2.0",&
|
||
&" 0.0, -3.0, 1.0"
|
||
rewind(unit=1) ! 回到文件1第一行
|
||
rewind(unit=2) ! 回到文件2第一行
|
||
read(1, *) a, b ! 读取矩阵1大小
|
||
allocate(array1(a, b), stat=errorLevel) ! 分配矩阵1
|
||
if (errorLevel /= 0) print *, "array1: Allocation request denied"
|
||
read(1, *) (array1(i, :), i = 1, a) ! 给矩阵1赋值
|
||
read(2, *) a, b ! 读取矩阵2大小
|
||
allocate(array2(a, b), stat=errorLevel) ! 分配矩阵2
|
||
if (errorLevel /= 0) print *, "array2: Allocation request denied"
|
||
read(2, *) (array2(i, :), i = 1, a) ! 给矩阵2赋值
|
||
! 关闭打开的临时文件(实际上一般可以不管)
|
||
close(unit=1)
|
||
close(unit=2)
|
||
if (size(array1, 2) /= size(array2, 1)) stop "警告! 数组维度不匹配!"! 检测是否可进行矩阵乘法
|
||
allocate(array3(size(array1, 1), size(array2, 2)), stat=errorLevel) ! 分配矩阵3
|
||
if (errorLevel /= 0) print *, "array3: Allocation request denied"
|
||
! array3 = MatrixMultplication(array1, array2)
|
||
array3 = 0. ! 初始化矩阵3
|
||
! 习题8-14, 8-15使用
|
||
! FORALL结构给矩阵3赋值, 注意这里会警告重复计算, 别管他, 设计就是这样
|
||
! forall(i = 1:size(array1, 1), j = 1:size(array1, 2), k = 1:size(array2, 2))
|
||
! array3(i, k) = array1(i, j)*array2(j, k) + array3(i, k)
|
||
! end forall
|
||
! 习题8-16 使用(改写后的)
|
||
array3 = matmul(array1, array2)
|
||
! 关闭数组, 释放资源
|
||
if (allocated(array2)) deallocate(array2, stat=errorLevel)
|
||
if (errorLevel /= 0) print *, "array2: Deallocation request denied"
|
||
if (allocated(array1)) deallocate(array1, stat=errorLevel)
|
||
if (errorLevel /= 0) print *, "array1: Deallocation request denied"
|
||
write(*, 2) (array3(i, :), i = 1, size(array3, 1)) ! 输出矩阵乘法计算结果
|
||
2 format(2(1X, F7.2, 1X))
|
||
if (allocated(array3)) deallocate(array3, stat=errorLevel)
|
||
if (errorLevel /= 0) print *, "array3: Deallocation request denied"
|
||
! contains
|
||
! function MatrixMultplication(arrayInput1, arrayInput2)
|
||
! implicit none
|
||
! real, dimension(:,:), allocatable, intent(in) :: arrayInput1, arrayInput2
|
||
! real, dimension(:,:), allocatable :: MatrixMultplication
|
||
! integer :: i, j, k
|
||
! if (size(arrayInput1, 2) /= size(arrayInput2, 1)) stop "警告! 数组维度不匹配!"
|
||
! MatrixMultplication = 0.
|
||
! forall(i = 1:size(arrayInput1, 1), j = 1:size(arrayInput1, 2), k = 1:size(arrayInput2, 2))
|
||
! MatrixMultplication(i, k) = arrayInput1(i, j)*arrayInput2(j, k) + MatrixMultplication(i, k)
|
||
! end forall
|
||
! end function MatrixMultplication
|
||
end subroutine Exercises8_14
|
||
! 习题8-17 相对极大(鞍点)
|
||
! 根据鞍点的定义, 可以直接遍历程序中每个点,,, 不过这样好像效率不太高的样子
|
||
! 1. 鞍点比上一行的值都大, 比下一行的值都大, 比左一列的值都大, 比右一列的值都大
|
||
! 2. 检测筛选之后的值
|
||
! 所以, 需要有一个输入矩阵, 有一个标定矩阵, 确定其位置的, i, j, 以及坐标数组(可扩展)
|
||
subroutine Exercises8_17
|
||
implicit none
|
||
real, dimension(:, :), allocatable :: arrayInput ! 输入矩阵
|
||
logical, dimension(:, :), allocatable :: arrayCheck ! 检测矩阵
|
||
integer :: i, j ! 循环参数
|
||
integer, dimension(:, :), allocatable :: arrayPoint ! 鞍点坐标
|
||
integer :: sizeX, sizeY ! 矩阵大小
|
||
integer :: errorLevel ! 错误码
|
||
integer :: num ! 鞍点个数
|
||
character(len=20) :: inforError ! 错误信息
|
||
! 变量初始化
|
||
num = 0
|
||
open(unit=1, status='scratch', iostat=errorLevel, iomsg=inforError) ! 打开临时文件
|
||
if (errorLevel /= 0) stop inforError ! 确认临时文件成功打开
|
||
write (1, '(A)') " 8, 8",& ! 向临时文件写入矩阵
|
||
&" 2.0, -1.0, -2.0, 1.0, 3.0, -5.0, 2.0, 1.0",&
|
||
&" -2.0, 0.0, -2.5, 5.0, -2.0, 2.0, 1.0, 0.0",&
|
||
&" -3.0, -3.0, -3.0, 3.0, 0.0, 0.0, -1.0, -2.0",&
|
||
&" -4.5, -4.0, -7.0, 6.0, 1.0, -3.0, 0.0, 5.0",&
|
||
&" -3.5, -3.0, -5.0, 0.0, 4.0, 17.0, 11.0, 5.0",&
|
||
&" -9.0, -6.0, -5.0, -3.0, 1.0, 2.0, 0.0, 0.5",&
|
||
&" -7.0, -4.0, -5.0, -3.0, 2.0, 4.0, 3.0, -1.0",&
|
||
&" -6.0, -5.0, -5.0, -2.0, 0.0, 1.0, 2.0, 5.0"
|
||
rewind(unit=1) ! 返回初始行
|
||
read(1, *) sizeY, sizeX ! 获得矩阵大小
|
||
allocate(arrayInput(sizeX, sizeY), stat=errorLevel) ! 定义输入矩阵
|
||
if (errorLevel /= 0) print *, "arrayInput: Allocation request denied"
|
||
allocate(arrayCheck(sizeX, sizeY), stat=errorLevel) ! 定义检测矩阵
|
||
if (errorLevel /= 0) print *, "arrayCheck: Allocation request denied"
|
||
allocate(arrayPoint(2, 1), stat=errorLevel) ! 定义鞍点坐标数组
|
||
if (errorLevel /= 0) print *, "arrayPoint: Allocation request denied"
|
||
read(1, *) (arrayInput(i, :), i = 1, sizeY) ! 给输入矩阵赋值
|
||
close(unit=1, iostat=errorLevel, iomsg=inforError) ! 关闭临时文件
|
||
if (errorLevel /= 0) stop inforError ! 检测是否关闭成功
|
||
arrayCheck = .FALSE. ! 初始化检测矩阵
|
||
do i = 2, sizeX - 1 ! 标记行间最大值
|
||
arrayCheck(i, 2:sizeY-1) = (arrayInput(i, 2:sizeY-1) > arrayInput(i-1, 2:sizeY-1)) .and. &
|
||
&(arrayInput(i, 2:sizeY-1) > arrayInput(i+1, 2:sizeY-1))
|
||
end do
|
||
do j = 2, sizeY - 1 ! 标记行列间最大值
|
||
arrayCheck(2:sizeX-1, j) = (arrayInput(2:sizeX-1, j) > arrayInput(2:sizeX-1, j-1)) .and. &
|
||
&(arrayInput(2:sizeX-1, j) > arrayInput(2:sizeX-1, j+1)) .and. &
|
||
& arrayCheck(2:sizeX-1, j)
|
||
end do
|
||
forall (i=2:sizeX-1, j=2:sizeY-1, arrayCheck(i, j)) ! 标记点与左上比较并迭代
|
||
arrayCheck(i, j) = (arrayInput(i, j) > arrayInput(i-1, j-1))
|
||
end forall
|
||
forall (i=2:sizeX-1, j=2:sizeY-1, arrayCheck(i, j)) ! 标记点与左下比较并迭代
|
||
arrayCheck(i, j) = (arrayInput(i, j) > arrayInput(i+1, j-1))
|
||
end forall
|
||
forall (i=2:sizeX-1, j=2:sizeY-1, arrayCheck(i, j)) ! 标记点与右下比较并迭代
|
||
arrayCheck(i, j) = (arrayInput(i, j) > arrayInput(i+1, j+1))
|
||
end forall
|
||
forall (i=2:sizeX-1, j=2:sizeY-1, arrayCheck(i, j)) ! 标记点与右上得到鞍点
|
||
arrayCheck(i, j) = (arrayInput(i, j) > arrayInput(i-1, j+1))
|
||
end forall
|
||
if (allocated(arrayInput)) deallocate(arrayInput, stat=errorLevel) ! 关闭输入矩阵
|
||
if (errorLevel /= 0) print *, "arrayInput: Deallocation request denied"
|
||
do i = 2, sizeX - 1
|
||
do j = 2, sizeY - 1
|
||
if (arrayCheck(j, i) .eqv. .TRUE.) then
|
||
num = num + 1
|
||
arrayPoint(:, num) = (/i, j/) ! 遍历检测矩阵, 得到鞍点坐标
|
||
end if
|
||
end do
|
||
end do
|
||
if (allocated(arrayCheck)) deallocate(arrayCheck, stat=errorLevel) ! 关闭检测矩阵
|
||
if (errorLevel /= 0) print *, "arrayCheck: Deallocation request denied"
|
||
write(*, 1) '矩阵中的鞍点坐标如下:', (arrayPoint(:, i), i = 1, num) ! 输出结果
|
||
1 format(1X, A/&
|
||
&6('(', I2, ',', I2, ') ')) ! 定义输出格式
|
||
if (allocated(arrayPoint)) deallocate(arrayPoint, stat=errorLevel) ! 关闭坐标数组
|
||
if (errorLevel /= 0) print *, "arrayPoint: Deallocation request denied"
|
||
end subroutine Exercises8_17
|
||
! 习题8-18 金属盘温度分布, 总感觉这题好像有点难度(如果要优化算法的话, 直接穷举确实容易倒是, 就是浪费计算资源, 模仿自然温度扩散写个循环吧)
|
||
! 基本路线原则如下:
|
||
! 1. 计算前进方向一格周围的格子情况, 以确定下一次前进方向, 并累计前进次数
|
||
! 1.1 当前进次数为0时, 朝方向1前进1格, 前进次数+1,
|
||
! 1.2 当前进次数大于1时:
|
||
! 1.2.1 检测前进方向格子状态
|
||
! 1.2.1.1 如果为0(未计算), 前进一格, 并检测 优先侧方向 格子状态
|
||
! 1.2.1.1.1 如果为0(未计算), 则下次前进方向为 优先侧方向
|
||
! 1.2.1.1.2 如果为1(已计算)或-1(固定值), 检测 前进方向 格子状态
|
||
! 1.2.1.1.2.1 如果为0(未计算), 则下次前进方向为 当前方向
|
||
! 1.2.1.1.2.2 如果为1(已计算)或-1(固定值). 检测 非优先侧方向
|
||
! 1.2.1.1.2.2.1 如果为0(未计算), 则下次前进方向为 非优先侧方向
|
||
! 1.2.1.1.2.2.2 如果为1(已计算)或-1(固定值), 这时候无路可走, 检测所有格子状态
|
||
! 1.2.1.1.2.2.2.1 所有格子为1, 检测此次迭代计算结果
|
||
! 1.2.1.1.2.2.2.1.1 符合条件, 退出
|
||
! 1.2.1.1.2.2.2.1.1 不符合条件, 初始化检测状态, 进行下一次循环
|
||
! 1.2.1.1.2.2.2.2 存在格子不为1, 选定可选位置中
|
||
! 1.2.1.2 如果为1(已计算或固定值)
|
||
! 突然发现越来越复杂了, 舍弃这种思路
|
||
! 按照距离优先计算, 优先计算离热源近的方块, 但问题是这样只能适应单热源的部分, 无法计算多热源(说得好像之前那种就能计算多热源了)
|
||
! 为了处理多热源的问题, 引入时序步进的想法, 假定热量传递速度相同
|
||
! 这样就是多热源画同样的同心圆, 在此范围内未被计算过的方块会被标记为2
|
||
! 对标记为2的方块进行统计, 周围被标记为0的格子数为2的, 被标记为3
|
||
! 优先计算标记为2的方块, 再计算标记为3的方块
|
||
! 步进计算, 直到无方块被标记
|
||
! 从距离的角度, 选取欧式距离似乎最为合适, 选取曼哈顿距离则计算最为简单, 最终还是选取了欧式距离
|
||
! 有时候感觉自己就是个ZZ, 直接算出来每个点的距离, 然后排序计算不就好了, 干嘛费那么大劲还去考虑步进每次包括的点位和点位计算顺序
|
||
subroutine Exercises8_18
|
||
implicit none
|
||
real, dimension(10, 10) :: metalPlate, metalPlateOld ! 金属盘表面温度, 及对照组
|
||
integer, dimension(63, 2) :: pointList ! 优化后的计算顺序
|
||
integer :: i, j, k ! 循环参数(什么都干的那种)
|
||
real :: temp, distense ! 用于判断是否符合要求
|
||
! 初始化计算顺序用数组
|
||
metalPlate = 0.
|
||
metalPlate(2:9, 2:9) = 1.
|
||
metalPlate(3, 8) = 0.
|
||
! 计算需要计算温度处距离热源的距离
|
||
forall (i = 2: 9, j = 2: 9, metalPlate(i, j) == 1.)
|
||
metalPlate(i, j) = (i-3)**2 + (j-8)**2
|
||
end forall
|
||
i = 0
|
||
! 输出计算顺序
|
||
do
|
||
i = i + 1
|
||
pointList(i, :) = minloc(metalPlate, mask=metalPlate>0.)
|
||
metalPlate(pointList(i, 1), pointList(i, 2)) = 0.
|
||
if (all(metalPlate == 0.)) exit
|
||
end do
|
||
! 初始化金属盘表面温度
|
||
metalPlate = 20.
|
||
metalPlate(2:9, 2:9) = 50.
|
||
metalPlate(3, 8) = 100.
|
||
j = 0
|
||
! 开始迭代计算金属盘表面温度
|
||
do
|
||
distense = 0.
|
||
do i = 1, size(pointList, 1)
|
||
temp = metalPlate(pointList(i, 1), pointList(i, 2))
|
||
metalPlate(pointList(i, 1), pointList(i, 2)) = (metalPlate(pointList(i, 1)+1, pointList(i, 2)) + &
|
||
&metalPlate(pointList(i, 1), pointList(i, 2)+1) + metalPlate(pointList(i, 1)-1, pointList(i, 2)) &
|
||
&+ metalPlate(pointList(i, 1), pointList(i, 2)-1))/4.
|
||
if (abs(temp-metalPlate(pointList(i, 1), pointList(i, 2))) > distense)&
|
||
& distense = abs(temp-metalPlate(pointList(i, 1), pointList(i, 2)))
|
||
end do
|
||
if (distense < 0.01) exit
|
||
j = j + 1
|
||
end do
|
||
write(*, '(" 共经历", I3, "次迭代计算")') j
|
||
write(*, *) '节点(5, 5)处稳定状态下温度为', metalPlate(5, 5)
|
||
write(*, '(" ", 10F7.2)') (metalPlate(i, :), i = 1, 10)
|
||
! 初始化金属盘表明温度
|
||
metalPlate = 20.
|
||
metalPlate(2:9, 2:9) = 50.
|
||
metalPlate(3, 8) = 100.
|
||
k = 0
|
||
! 开始迭代计算金属盘表面温度
|
||
do
|
||
metalPlateOld = metalPlate
|
||
forall(i = 1:10, j = 1:10, metalPlate(i, j) /= 20. .and. metalPlate(i, j) /= 100.)
|
||
metalPlate(i, j) = (metalPlate(i+1, j) + metalPlate(i, j+1) + metalPlate(i-1, j) + metalPlate(i, j-1))/4.
|
||
end forall
|
||
k = k + 1
|
||
if (all(abs(metalPlate-metalPlateOld) < 0.01)) exit
|
||
end do
|
||
write(*, '(" 共经历", I3, "次迭代计算")') k
|
||
write(*, *) '节点(5, 5)处稳定状态下温度为', metalPlate(5, 5)
|
||
write(*, '(" ", 10F7.2)') (metalPlate(i, :), i = 1, 10)
|
||
end subroutine Exercises8_18
|