Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
! Heat equation solver in 2D.
program heat_solve
use heat
use core
use io
use setup
use utilities
use mpi
implicit none
real(dp), parameter :: a = 0.5 ! Diffusion constant
type(field) :: current, previous ! Current and previus temperature fields
real(dp) :: dt ! Time step
integer :: nsteps ! Number of time steps
integer, parameter :: image_interval = 500 ! Image output interval
integer, parameter :: checkpoint_interval = 200 ! restart interval
type(parallel_data) :: parallelization
integer :: ierr
integer :: iter, iter0
real(kind=dp) :: start, stop ! Timers
call mpi_init(ierr)
call initialize(current, previous, nsteps, parallelization, iter0)
! Draw the picture of the initial state
call write_field(current, iter0, parallelization)
iter0 = iter0 + 1
! Largest stable time step
dt = current%dx**2 * current%dy**2 / &
& (2.0 * a * (current%dx**2 + current%dy**2))
! Main iteration loop, save a picture every
! image_interval steps
start = mpi_wtime()
do iter = iter0, iter0 + nsteps
call exchange_init(previous, parallelization)
call evolve_interior(current, previous, a, dt)
call exchange_finalize(parallelization)
call evolve_edges(current,previous, a, dt)
if (mod(iter, image_interval) == 0) then
call write_field(current, iter, parallelization)
if (mod(iter, checkpoint_interval) == 0) then
call write_restart(current, parallelization, iter)
end if
end if
call swap_fields(current, previous)
end do
stop = mpi_wtime()
if (parallelization % rank == 0) then
write(*,'(A,F7.3,A)') 'Iteration took ', stop - start, ' seconds.'
write(*,'(A,G0)') 'Reference value at 5,5: ', previous % data(5,5)
end if
call finalize(current, previous)
call mpi_finalize(ierr)
end program heat_solve