DISPATCH
solver_mod.f90
1 !===============================================================================
2 !> This module contains all experiment specific information necessary to solve
3 !> the heat diffusion problem in DISPATCH
4 !===============================================================================
5 MODULE solver_mod
6  USE io_mod
7  USE trace_mod
8  USE extras_mod
9  implicit none
10  private
11  type, public, extends(extras_t):: solver_t
12  real, pointer:: radii(:,:,:)
13  real, pointer:: temperature(:,:,:)
14  real:: initial_temperature=1.0
15  real:: outside_temperature=2.0
16  real:: radius=0.45
17  contains
18  procedure:: init
19  procedure:: update
20  procedure:: boundary_condition
21  end type
22 CONTAINS
23 
24 !===============================================================================
25 !> Setup an initial state. The box size is 1, with the origin at 0,0,0, so
26 !> the max radius is 0.5
27 !===============================================================================
28 SUBROUTINE init (self)
29  class(solver_t):: self
30  integer:: m(6), i1, i2, i3, it
31  !-----------------------------------------------------------------------------
32  call trace%begin('solver_t%init')
33  !-----------------------------------------------------------------------------
34  ! This initializers a default do-nothing solver, and allocates the
35  ! self%mem() 6-dimenstional array (nx,ny,nz,nv,nt,nw), in which the 4th
36  ! index (= 1) here is the variable index, and the 5th index is the time slot.
37  !-----------------------------------------------------------------------------
38  self%nv = 1
39  self%nt = 5
40  call self%patch_t%init
41  self%temperature => self%mem(:,:,:,1,1,1)
42  self%temperature = self%initial_temperature
43  !-----------------------------------------------------------------------------
44  ! For convenience, we also allocate the radius (distance from origin) for
45  ! each cell. The array should have the same size as the temperature:
46  !-----------------------------------------------------------------------------
47  m = shape(self%mem)
48  allocate (self%radii(m(1),m(2),m(3)))
49  do i3=1,m(3)
50  do i2=1,m(2)
51  do i1=1,m(1)
52  !---------------------------------------------------------------------------
53  ! Each mesh (x,y,z) has properties stored in a data type (object), where
54  ! %p is the position (of the center), and %r are the relative cell positions
55  !---------------------------------------------------------------------------
56  self%radii(i1,i2,i3) = sqrt((self%mesh(1)%p + self%mesh(1)%r(i1))**2 + &
57  (self%mesh(2)%p + self%mesh(2)%r(i2))**2 + &
58  (self%mesh(3)%p + self%mesh(3)%r(i3))**2)
59  end do
60  end do
61  end do
62  !-----------------------------------------------------------------------------
63  ! Set temperature inside and then apply boundary condition
64  !-----------------------------------------------------------------------------
65  call self%boundary_condition (1)
66  call trace%end
67 END SUBROUTINE init
68 
69 !===============================================================================
70 !> Apply the boundary condition T=temperature(1) outside the given radius
71 !===============================================================================
72 SUBROUTINE boundary_condition (self, it)
73  class(solver_t):: self
74  integer:: it
75  real, pointer, contiguous:: f(:,:,:)
76  !-----------------------------------------------------------------------------
77  call trace%begin('solver_t%boundary_condition')
78  !-----------------------------------------------------------------------------
79  f => self%mem(:,:,:,1,it,1)
80  where (self%radii > self%radius)
81  f = self%outside_temperature
82  end where
83  call trace%end
84 END SUBROUTINE boundary_condition
85 
86 !===============================================================================
87 !> Update the solution, by adding a fraction times the diffusion operator
88 !===============================================================================
89 SUBROUTINE update (self)
90  class(solver_t):: self
91  real, pointer, contiguous:: v(:,:,:)
92  real, allocatable:: d2f(:,:,:)
93  integer:: m(3), i1, i2, i3
94  integer, save:: itimer=0
95  !----------------------------------------------------------------------------
96  call trace%begin('solver_t%update', itimer=itimer)
97  call self%output
98  !----------------------------------------------------------------------------
99  ! Copy the variable (temperature) from the current slot (%it) to the %new slot
100  !----------------------------------------------------------------------------
101  v => self%mem(:,:,:,1,self%new,1)
102  v = self%mem(:,:,:,1,self%it ,1)
103  call self%boundary_condition (self%new)
104  !----------------------------------------------------------------------------
105  ! Allocate a scratch variable for the diffusion (Laplace) operator
106  !----------------------------------------------------------------------------
107  m = shape(v)
108  allocate (d2f(m(1),m(2),m(3)))
109  !----------------------------------------------------------------------------
110  ! Evaluate a 6-point diffusion operator over the internal patch region
111  ! (no need to do it in the guard zones)
112  !----------------------------------------------------------------------------
113  do i3=self%mesh(3)%li,self%mesh(3)%ui
114  do i2=self%mesh(2)%li,self%mesh(2)%ui
115  do i1=self%mesh(1)%li,self%mesh(1)%ui
116  d2f(i1,i2,i3) = v(i1+1,i2 ,i3 ) &
117  + v(i1-1,i2 ,i3 ) &
118  + v(i1 ,i2+1,i3 ) &
119  + v(i1 ,i2-1,i3 ) &
120  + v(i1 ,i2 ,i3+1) &
121  + v(i1 ,i2 ,i3-1) &
122  - v(i1 ,i3 ,i3 )*6.0
123  end do
124  end do
125  end do
126  !----------------------------------------------------------------------------
127  ! Set the time step, update the temperature, and apply boundary condition
128  !----------------------------------------------------------------------------
129  self%dtime = self%courant
130  v = v + self%dtime*d2f
131  call self%boundary_condition (self%new)
132  !----------------------------------------------------------------------------
133  ! Clean up and count operations
134  !----------------------------------------------------------------------------
135  deallocate (d2f)
136  call self%counter_update
137  call trace%end (itimer)
138 END SUBROUTINE update
139 
140 END MODULE solver_mod
Template version of a module intended for adding extra features, at a level between the basic task an...
Definition: extras_mod.f90:54
Definition: io_mod.f90:4
This module contains all experiment specific information necessary to solve the heat diffusion proble...
Definition: solver_mod.f90:5