22 integer,
parameter:: nslice = 4
25 integer:: time_order=-1
30 procedure:: update_single
38 SUBROUTINE init (self)
40 integer,
save:: verbose=0, time_order=-1
41 namelist /timestep_params/ verbose, time_order
43 call trace_begin (
'timestep_t%init', 1)
45 if (time_order==-1)
then 48 read (io%input,timestep_params)
49 if (io%master)
write (*,timestep_params)
52 self%verbose = verbose
53 self%time_order = time_order
55 select case (time_order)
58 self%wt = [ 1.0, 0.0, 0.0, 0.0]
60 self%wt = [ 1.5, -.5, 0.0, 0.0]
62 self%wt = [23./12., -4./3., 5./12., 0.0]
64 self%wt = [55./24., -59./24., 37/24., -3./8.]
66 if (io%master) print*,
'ERROR: illegal timestep time_order =', time_order
78 SUBROUTINE update (self, id, iit, dt, mem, mesh, lock)
83 integer,
dimension(:):: iit
84 real(8),
dimension(:):: dt
85 real(kind=KindScalarVar),
dimension(:,:,:,:,:,:):: mem
87 integer,
save:: itimer=0
89 call trace_begin (
'timestep_t%update', itimer=itimer)
90 if (self%time_order==-1)
call io%abort(
'timestep_mod not initialized!')
94 where (mesh%lower_boundary)
97 where (mesh%upper_boundary)
100 if (
size(mem,6)==1.and.self%time_order>0)
then 101 if (io%master) print *,id,
'WARNING: setting time_order=0 in timestep%update' 104 call lock%set (
'timestep')
105 select case (self%time_order)
108 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt ),1) = &
109 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-1),1) + &
110 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-1),2)*(self%wt(1)*dt(iit(nt-1)))
112 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt ),1) = &
113 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-1),1) + &
114 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-1),2)*(self%wt(1)*dt(iit(nt-1))) + &
115 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-2),2)*(self%wt(2)*dt(iit(nt-2)))
117 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt ),1) = &
118 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-1),1) + &
119 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-1),2)*(self%wt(1)*dt(iit(nt-1))) + &
120 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-2),2)*(self%wt(2)*dt(iit(nt-2))) + &
121 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-3),2)*(self%wt(3)*dt(iit(nt-3)))
123 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt ),1) = &
124 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-1),1) + &
125 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-1),2)*(self%wt(1)*dt(iit(nt-1))) + &
126 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-2),2)*(self%wt(2)*dt(iit(nt-2))) + &
127 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-3),2)*(self%wt(3)*dt(iit(nt-3))) + &
128 mem(l(1):u(1),l(2):u(2),l(3):u(3),:,iit(nt-4),2)*(self%wt(4)*dt(iit(nt-4)))
130 if (self%verbose>3)
then 131 print
'(a,i4,2x,a,1p,6e10.2,2x,6i5)',
'timestep, id:', &
133 minval(mem(l(1):u(1),l(2):u(2),l(3):u(3),1,iit(nt-1),2)), &
134 maxval(mem(l(1):u(1),l(2):u(2),l(3):u(3),1,iit(nt-1),2)), &
135 minval(mem(l(1):u(1),l(2):u(2),l(3):u(3),1,iit(nt-1),1)), &
136 maxval(mem(l(1):u(1),l(2):u(2),l(3):u(3),1,iit(nt-1),1)), &
137 minval(mem(l(1):u(1),l(2):u(2),l(3):u(3),1,iit(nt ),1)), &
138 maxval(mem(l(1):u(1),l(2):u(2),l(3):u(3),1,iit(nt ),1)), shape(mem)
140 call lock%unset (
'timestep')
141 call trace_end (itimer)
142 END SUBROUTINE update
145 SUBROUTINE update_single (self, id, iit, dt, mem, dmem, mesh, lock, time_order)
150 integer,
dimension(:):: iit
151 real(8),
dimension(:):: dt
152 real(kind=KindScalarVar),
dimension(:,:,:,:):: mem, dmem
153 integer,
intent(in):: time_order
156 integer,
save:: itimer=0
158 call trace_begin (
'timestep_t%update_single', itimer=itimer)
159 if (time_order <= 0)
return 164 where (mesh%lower_boundary)
167 where (mesh%upper_boundary)
171 select case (time_order)
174 wt = [ 1.0, 0.0, 0.0, 0.0]
175 mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt )) = &
176 mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1)) + &
177 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1))*(wt(1)*dt(iit(nt-1)))
179 wt = [ 1.5, -.5, 0.0, 0.0]
180 mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt )) = &
181 mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1)) + &
182 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1))*(wt(1)*dt(iit(nt-1))) + &
183 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-2))*(wt(2)*dt(iit(nt-2)))
185 wt = [23./12., -4./3., 5./12., 0.0]
186 mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt )) = &
187 mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1)) + &
188 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1))*(wt(1)*dt(iit(nt-1))) + &
189 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-2))*(wt(2)*dt(iit(nt-2))) + &
190 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-3))*(wt(3)*dt(iit(nt-3)))
192 wt = [55./24., -59./24., 37/24., -3./8.]
193 mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt )) = &
194 mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1)) + &
195 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1))*(wt(1)*dt(iit(nt-1))) + &
196 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-2))*(wt(2)*dt(iit(nt-2))) + &
197 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-3))*(wt(3)*dt(iit(nt-3))) + &
198 dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-4))*(wt(4)*dt(iit(nt-4)))
200 if (self%verbose>3)
then 201 print
'(a,i4,2x,a,1p,6e10.2,2x,6i5)',
'timestep, id:', &
203 minval(dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1))), &
204 maxval(dmem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1))), &
205 minval(mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1))), &
206 maxval(mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt-1))), &
207 minval(mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt ))), &
208 maxval(mem(l(1):u(1),l(2):u(2),l(3):u(3),iit(nt ))), shape(mem)
210 call trace_end (itimer)
211 END SUBROUTINE update_single
Template module for mesh.
The lock module uses nested locks, to allow versatile use of locks, where a procedure may want to mak...