MPI-AMRVAC 3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
Loading...
Searching...
No Matches
mod_ppm.t
Go to the documentation of this file.
1module mod_ppm
2
3 implicit none
4 private
5
6 public :: ppmlimiter
7 public :: ppmlimitervar
8
9contains
10
11 subroutine ppmlimitervar(ixI^L,ix^L,idims,q,qCT,qLC,qRC)
12
13 ! references:
14 ! Mignone et al 2005, ApJS 160, 199,
15 ! Miller and Colella 2002, JCP 183, 26
16 ! Fryxell et al. 2000 ApJ, 131, 273 (Flash)
17 ! baciotti Phd (http://www.aei.mpg.de/~baiotti/Baiotti_PhD.pdf)
18 ! old version : april 2009
19 ! author: zakaria.meliani@wis.kuleuven.be
20 ! current version : March 2023 by Chun Xia
21
23
24 integer, intent(in) :: ixi^l, ix^l, idims
25 double precision, intent(in) :: q(ixi^s),qct(ixi^s)
26
27 double precision, intent(inout) :: qrc(ixi^s),qlc(ixi^s)
28
29 double precision,dimension(ixI^S) :: dqc,d2qc,ldq
30 double precision,dimension(ixI^S) :: qmin,qmax,tmp
31
32 integer :: lxc^l,lxr^l,lxl^l
33 integer :: ixl^l,ixo^l,ixr^l,ixol^l,ixor^l
34 integer :: hxl^l,hxc^l,hxr^l
35 integer :: kxl^l,kxc^l,kxr^l
36
37 ixomin^d=ixmin^d-kr(idims,^d);ixomax^d=ixmax^d+kr(idims,^d);!ixO[ixMmin-1,ixMmax+1]
38 ixolmin^d=ixomin^d-kr(idims,^d);ixolmax^d=ixomax^d; !ixOL[ixMmin-2,ixMmax+1]
39 ixor^l=ixol^l+kr(idims,^d); !ixOR=[iMmin-1,ixMmax+2]
40 ixl^l=ixo^l-kr(idims,^d); !ixL[ixMmin-2,ixMmax]
41 ixr^l=ixo^l+kr(idims,^d); !ixR=[iMmin,ixMmax+2]
42
43 hxcmin^d=ixomin^d;hxcmax^d=ixmax^d; ! hxC = [ixMmin-1,ixMmax]
44 hxl^l=hxc^l-kr(idims,^d); ! hxL = [ixMmin-2,ixMmax-1]
45 hxr^l=hxc^l+kr(idims,^d); ! hxR = [ixMmin,ixMmax+1]
46
47 kxcmin^d=ixlmin^d-1; kxcmax^d=ixrmax^d; ! kxC=[iMmin-3,ixMmax+2]
48 kxr^l=kxc^l+kr(idims,^d); ! kxR=[iMmin-2,ixMmax+3]
49
50 lxcmin^d=ixomin^d-kr(idims,^d);lxcmax^d=ixomax^d+kr(idims,^d);! lxC=[iMmin-2,ixMmax+2]
51 lxl^l=lxc^l-kr(idims,^d); ! lxL=[iMmin-3,ixMmax+1]
52 lxr^l=lxc^l+kr(idims,^d); ! lxR=[iMmin-1,ixMmax+3]
53
54 dqc(kxc^s)=q(kxr^s)-q(kxc^s)
55 ! Eq. 64, Miller and Colella 2002, JCP 183, 26
56 d2qc(kxc^s)=half*(q(lxr^s)-q(lxl^s))
57 where(dqc(lxc^s)*dqc(lxl^s)>zero)
58 ! Store the sign of dwC in wMin
59 qmin(kxc^s)= sign(one,d2qc(lxc^s))
60 ! Eq. 65, Miller and Colella 2002, JCP 183, 26
61 ldq(lxc^s)= qmin(lxc^s)*min(dabs(d2qc(lxc^s)),&
62 2.0d0*dabs(dqc(lxl^s)),&
63 2.0d0*dabs(dqc(lxc^s)))
64 else where
65 ldq(lxc^s)=zero
66 end where
67
68 ! Eq. 66, Miller and Colella 2002, JCP 183, 26
69 qlc(ixol^s)=qlc(ixol^s)+half*dqc(ixol^s)&
70 +(ldq(ixol^s)-ldq(ixor^s))/6.0d0
71 qrc(ixl^s)=qlc(ixl^s)
72
73 ! make sure that min wCT(i)<wLC(i)<wCT(i+1) same for wRC(i)
74 call extremaq(ixi^l,ixo^l,qct,1,qmax,qmin)
75
76 ! Eq. B8, page 217, Mignone et al 2005, ApJS
77 qrc(ixl^s)=max(qmin(ixo^s),min(qmax(ixo^s),qrc(ixl^s)))
78 qlc(ixo^s)=max(qmin(ixo^s),min(qmax(ixo^s),qlc(ixo^s)))
79
80 ! Eq. B9, page 217, Mignone et al 2005, ApJS
81 where((qrc(ixl^s)-qct(ixo^s))*(qct(ixo^s)-qlc(ixo^s))<=zero)
82 qrc(ixl^s)=qct(ixo^s)
83 qlc(ixo^s)=qct(ixo^s)
84 end where
85
86 qmax(ixo^s)=(qlc(ixo^s)-qrc(ixl^s))&
87 *(qct(ixo^s)-half*(qlc(ixo^s)+qrc(ixl^s)))
88 qmin(ixo^s)=(qlc(ixo^s)-qrc(ixl^s))**2/6.0d0
89
90 tmp(hxl^s)=qrc(hxl^s)
91 ! Eq. B10, page 218, Mignone et al 2005, ApJS
92 where(qmax(hxr^s)>qmin(hxr^s))
93 qrc(hxc^s)= 3.0d0*qct(hxr^s)-2.0d0*qlc(hxr^s)
94 end where
95 ! Eq. B11, page 218, Mignone et al 2005, ApJS
96 where(qmax(hxc^s)<-qmin(hxc^s))
97 qlc(hxc^s)= 3.0d0*qct(hxc^s)-2.0d0*tmp(hxl^s)
98 end where
99
100 end subroutine ppmlimitervar
101
102 subroutine ppmlimiter(ixI^L,ix^L,idims,w,wCT,wLC,wRC)
103
104 ! references:
105 ! Mignone et al 2005, ApJS 160, 199,
106 ! Miller and Colella 2002, JCP 183, 26
107 ! Fryxell et al. 2000 ApJ, 131, 273 (Flash)
108 ! baciotti Phd (http://www.aei.mpg.de/~baiotti/Baiotti_PhD.pdf)
109 ! old version : april 2009
110 ! author: zakaria.meliani@wis.kuleuven.be
111 ! current version : March 2023 by Chun Xia
113
114 integer, intent(in) :: ixi^l, ix^l, idims
115 double precision, intent(in) :: w(ixi^s,1:nw),wct(ixi^s,1:nw)
116
117 double precision, intent(inout) :: wrc(ixi^s,1:nw),wlc(ixi^s,1:nw)
118
119 double precision,dimension(ixI^S,1:nw_recon) :: dwc,d2wc,ldw
120 double precision,dimension(ixI^S,1:nw_recon) :: wmin,wmax,tmp
121 double precision,dimension(ixI^S) :: aa, ab, ac, ki
122
123 double precision, parameter :: betamin=0.75d0, betamax=0.85d0,&
124 zmin=0.25d0, zmax=0.75d0,&
125 eta1=20.0d0,eta2=0.05d0,eps=0.01d0,kappa=0.1d0
126 integer :: lxc^l,lxl^l,lxr^l
127 integer :: ixll^l,ixl^l,ixo^l,ixr^l,ixrr^l,ixol^l,ixor^l
128 integer :: hxl^l,hxc^l,hxr^l
129 integer :: kxll^l,kxl^l,kxc^l,kxr^l,kxrr^l
130 integer :: iw, idimss
131
132 ixomin^d=ixmin^d-kr(idims,^d);ixomax^d=ixmax^d+kr(idims,^d);!ixO[ixMmin-1,ixMmax+1]
133 ixolmin^d=ixomin^d-kr(idims,^d);ixolmax^d=ixomax^d; !ixOL[ixMmin-2,ixMmax+1]
134 ixor^l=ixol^l+kr(idims,^d); !ixOR=[iMmin-1,ixMmax+2]
135 ixl^l=ixo^l-kr(idims,^d); !ixL[ixMmin-2,ixMmax]
136 ixr^l=ixo^l+kr(idims,^d); !ixR=[iMmin,ixMmax+2]
137
138 hxcmin^d=ixomin^d;hxcmax^d=ixmax^d; ! hxC = [ixMmin-1,ixMmax]
139 hxl^l=hxc^l-kr(idims,^d); ! hxL = [ixMmin-2,ixMmax-1]
140 hxr^l=hxc^l+kr(idims,^d); ! hxR = [ixMmin,ixMmax+1]
141
142 kxcmin^d=ixlmin^d-1; kxcmax^d=ixrmax^d; ! kxC=[iMmin-3,ixMmax+2]
143 kxr^l=kxc^l+kr(idims,^d); ! kxR=[iMmin-2,ixMmax+3]
144
145 lxcmin^d=ixomin^d-kr(idims,^d);lxcmax^d=ixomax^d+kr(idims,^d);! lxC=[iMmin-2,ixMmax+2]
146 lxl^l=lxc^l-kr(idims,^d); ! lxL=[iMmin-3,ixMmax+1]
147 lxr^l=lxc^l+kr(idims,^d); ! lxR=[iMmin-1,ixMmax+3]
148
149 dwc(kxc^s,1:nw_recon)=w(kxr^s,1:nw_recon)-w(kxc^s,1:nw_recon)
150 ! Eq. 64, Miller and Colella 2002, JCP 183, 26
151 d2wc(lxc^s,1:nw_recon)=half*(w(lxr^s,1:nw_recon)-w(lxl^s,1:nw_recon))
152 where(dwc(lxc^s,1:nw_recon)*dwc(lxl^s,1:nw_recon)>zero)
153 ! Store the sign of dwC in wMin
154 wmin(lxc^s,1:nw_recon)= sign(one,d2wc(lxc^s,1:nw_recon))
155 ! Eq. 65, Miller and Colella 2002, JCP 183, 26
156 ldw(lxc^s,1:nw_recon)= wmin(lxc^s,1:nw_recon)*min(dabs(d2wc(lxc^s,1:nw_recon)),&
157 2.0d0*dabs(dwc(lxl^s,1:nw_recon)),&
158 2.0d0*dabs(dwc(lxc^s,1:nw_recon)))
159 else where
160 ldw(lxc^s,1:nw_recon)=zero
161 end where
162
163 ! Eq. 66, Miller and Colella 2002, JCP 183, 26
164 wlc(ixol^s,1:nw_recon)=wlc(ixol^s,1:nw_recon)+half*dwc(ixol^s,1:nw_recon)&
165 +(ldw(ixol^s,1:nw_recon)-ldw(ixor^s,1:nw_recon))/6.0d0
166 wrc(ixl^s,1:nw_recon)=wlc(ixl^s,1:nw_recon)
167
168 ! make sure that min w(i)<wLC(i)<w(i+1) same for wRC(i)
169 call extremaw(ixi^l,ixo^l,w,1,wmax,wmin)
170
171 ! Eq. B8, page 217, Mignone et al 2005, ApJS
172 wrc(ixl^s,1:nw_recon)=max(wmin(ixo^s,1:nw_recon)&
173 ,min(wmax(ixo^s,1:nw_recon),wrc(ixl^s,1:nw_recon)))
174 wlc(ixo^s,1:nw_recon)=max(wmin(ixo^s,1:nw_recon)&
175 ,min(wmax(ixo^s,1:nw_recon),wlc(ixo^s,1:nw_recon)))
176
177 ! Eq. B9, page 217, Mignone et al 2005, ApJS
178 where((wrc(ixl^s,1:nw_recon)-wct(ixo^s,1:nw_recon))&
179 *(wct(ixo^s,1:nw_recon)-wlc(ixo^s,1:nw_recon))<=zero)
180 wrc(ixl^s,1:nw_recon)=wct(ixo^s,1:nw_recon)
181 wlc(ixo^s,1:nw_recon)=wct(ixo^s,1:nw_recon)
182 end where
183
184 wmax(ixo^s,1:nw_recon)=(wlc(ixo^s,1:nw_recon)-wrc(ixl^s,1:nw_recon))*&
185 (wct(ixo^s,1:nw_recon)-half*(wlc(ixo^s,1:nw_recon)+wrc(ixl^s,1:nw_recon)))
186 wmin(ixo^s,1:nw_recon)=(wlc(ixo^s,1:nw_recon)-wrc(ixl^s,1:nw_recon))**2/6.0d0
187 tmp(hxl^s,1:nw_recon)=wrc(hxl^s,1:nw_recon)
188 ! Eq. B10, page 218, Mignone et al 2005, ApJS
189 where(wmax(hxr^s,1:nw_recon)>wmin(hxr^s,1:nw_recon))
190 wrc(hxc^s,1:nw_recon)= 3.0d0*wct(hxr^s,1:nw_recon)&
191 -2.0d0*wlc(hxr^s,1:nw_recon)
192 end where
193 ! Eq. B11, page 218, Mignone et al 2005, ApJS
194 where(wmax(hxc^s,1:nw_recon)<-wmin(hxc^s,1:nw_recon))
195 wlc(hxc^s,1:nw_recon)= 3.0d0*wct(hxc^s,1:nw_recon)&
196 -2.0d0*tmp(hxl^s,1:nw_recon)
197 end where
198
199 ! flattening at the contact discontinuities
200 if(flatcd)then
201 ixrr^l=ixr^l+kr(idims,^d); !ixRR=[iMmin+1,ixMmax+3]
202 kxl^l=kxc^l-kr(idims,^d); ! kxL=[iMmin-4,ixMmax+1]
203 call ppm_flatcd(ixi^l,kxc^l,kxl^l,kxr^l,wct,d2wc,aa,ab)
204 if(any(kappa*aa(kxc^s)>=ab(kxc^s)))then
205 do iw=1,nw_recon
206 where(kappa*aa(kxc^s)>=ab(kxc^s).and. dabs(dwc(kxc^s,iw))>smalldouble)
207 wmax(kxc^s,iw) = wct(kxr^s,iw)-2.0d0*wct(kxc^s,iw)+wct(kxl^s,iw)
208 end where
209
210 where(wmax(ixr^s,iw)*wmax(ixl^s,iw)<zero .and.&
211 dabs(wct(ixr^s,iw)-wct(ixl^s,iw))&
212 -eps*min(dabs(wct(ixr^s,iw)),dabs(wct(ixl^s,iw)))>zero &
213 .and. kappa*aa(ixo^s)>=ab(ixo^s)&
214 .and. dabs(dwc(ixo^s,iw))>smalldouble)
215
216 ac(ixo^s)=(wct(ixll^s,iw)-wct(ixrr^s,iw)+4.0d0*dwc(ixo^s,iw))&
217 /(12.0d0*dwc(ixo^s,iw))
218 wmin(ixo^s,iw)=max(zero,min(eta1*(ac(ixo^s)-eta2),one))
219 else where
220 wmin(ixo^s,iw)=zero
221 end where
222
223 where(wmin(hxc^s,iw)>zero)
224 wlc(hxc^s,iw) = wlc(hxc^s,iw)*(one-wmin(hxc^s,iw))&
225 +(wct(hxc^s,iw)+half*ldw(hxc^s,iw))*wmin(hxc^s,iw)
226 end where
227 where(wmin(hxr^s,iw)>zero)
228 wrc(hxc^s,iw) = wrc(hxc^s,iw)*(one-wmin(hxr^s,iw))&
229 +(wct(hxr^s,iw)-half*ldw(hxr^s,iw))*wmin(hxr^s,iw)
230 end where
231 end do
232 end if
233 end if
234
235 ! flattening at the shocks
236 if(flatsh)then
237 ! following MILLER and COLELLA 2002 JCP 183, 26
238 kxcmin^d=ixmin^d-2; kxcmax^d=ixmax^d+2; ! kxC=[ixMmin-2,ixMmax+2]
239 ki=bigdouble
240 do idimss=1,ndim
241 kxl^l=kxc^l-kr(idimss,^d); ! kxL=[ixMmin-3,ixMmax+1]
242 kxr^l=kxc^l+kr(idimss,^d); ! kxR=[ixMmin-1,ixMmax+3]
243 kxll^l=kxl^l-kr(idimss,^d);! kxLL=[ixMmin-4,ixMmax]
244 kxrr^l=kxr^l+kr(idimss,^d);! kxRR=[ixMmin,ixMmax+4]
245
246 call ppm_flatsh(ixi^l,kxc^l,kxll^l,kxl^l,kxr^l,kxrr^l,idimss,wct,aa,ab)
247
248 ! eq. B17, page 218, Mignone et al 2005, ApJS (Xi1 min)
249 ac(kxc^s) = max(zero,min(one,(betamax-aa(kxc^s))/(betamax-betamin)))
250 ! eq. B18, page 218, Mignone et al 2005, ApJS (Xi1)
251 ! recycling aa(ixL^S)
252 where(wct(kxr^s,iw_mom(idimss))<wct(kxl^s,iw_mom(idimss)))
253 aa(kxc^s) = max(ac(kxc^s), min(one,(zmax-ab(kxc^s))/(zmax-zmin)))
254 else where
255 aa(kxc^s) = one
256 end where
257 ixl^l=ixo^l-kr(idimss,^d); ! ixL=[ixMmin-2,ixMmax]
258 ixr^l=ixo^l+kr(idimss,^d); ! ixR=[ixMmin,ixMmax+2]
259 ki(ixo^s)=min(ki(ixo^s),aa(ixl^s),aa(ixo^s),aa(ixr^s))
260 end do
261 ! recycling wMax
262 do iw=1,nw_recon
263 where(dabs(ab(ixo^s)-one)>smalldouble)
264 wmax(ixo^s,iw) = (one-ab(ixo^s))*wct(ixo^s,iw)
265 end where
266
267 where(dabs(ab(hxc^s)-one)>smalldouble)
268 wlc(hxc^s,iw) = ab(hxc^s)*wlc(hxc^s,iw)+wmax(hxc^s,iw)
269 end where
270
271 where(dabs(ab(hxr^s)-one)>smalldouble)
272 wrc(hxc^s,iw) = ab(hxr^s)*wrc(hxc^s,iw)+wmax(hxr^s,iw)
273 end where
274 end do
275 end if
276
277 end subroutine ppmlimiter
278
279 subroutine ppm_flatcd(ixI^L,ixO^L,ixL^L,ixR^L,w,d2w,drho,dp)
281 use mod_physics, only: phys_gamma
282
283 integer, intent(in) :: ixi^l, ixo^l, ixl^l, ixr^l
284 double precision, intent(in) :: w(ixi^s, nw), d2w(ixg^t, 1:nw_recon)
285 double precision, intent(inout) :: drho(ixg^t), dp(ixg^t)
286
287 drho(ixo^s) = phys_gamma*abs(d2w(ixo^s, iw_rho))&
288 /min(w(ixl^s, iw_rho), w(ixr^s, iw_rho))
289 dp(ixo^s) = abs(d2w(ixo^s, iw_e))/min(w(ixl^s, iw_e), w(ixr^s, iw_e))
290
291 end subroutine ppm_flatcd
292
293 subroutine ppm_flatsh(ixI^L,ixO^L,ixLL^L,ixL^L,ixR^L,ixRR^L,idims,w,drho,dp)
295 use mod_physics, only: phys_gamma
296
297 integer, intent(in) :: ixi^l, ixo^l, ixll^l, ixl^l, ixr^l, ixrr^l
298 integer, intent(in) :: idims
299 double precision, intent(in) :: w(ixi^s, nw)
300 double precision, intent(inout) :: drho(ixi^s), dp(ixi^s)
301
302 ! eq. B15, page 218, Mignone and Bodo 2005, ApJS (beta1)
303 where (abs(w(ixrr^s, iw_e)-w(ixll^s, iw_e))>smalldouble)
304 drho(ixo^s) = abs((w(ixr^s, iw_e)-w(ixl^s, iw_e))&
305 /(w(ixrr^s, iw_e)-w(ixll^s, iw_e)))
306 else where
307 drho(ixo^s) = zero
308 end where
309
310 ! eq. 76, page 48, Miller and Colella 2002, JCoPh, adjusted
311 dp(ixo^s) = abs(w(ixr^s, iw_e)-w(ixl^s, iw_e))&
312 /(phys_gamma*(w(ixr^s, iw_e)+w(ixl^s, iw_e)))
313
314 end subroutine ppm_flatsh
315
316 subroutine extremaq(ixI^L,ixO^L,q,nshift,qMax,qMin)
317
319
320 integer,intent(in) :: ixi^l,ixo^l
321 double precision, intent(in) :: q(ixi^s)
322 integer,intent(in) :: nshift
323
324 double precision, intent(out) :: qmax(ixi^s),qmin(ixi^s)
325
326 integer :: ixs^l,ixsr^l,ixsl^l,idims,jdims,kdims,ishift,i,j
327
328 do ishift=1,nshift
329 idims=1
330 ixsr^l=ixo^l+ishift*kr(idims,^d);
331 ixsl^l=ixo^l-ishift*kr(idims,^d);
332 if (ishift==1) then
333 qmax(ixo^s)=max(q(ixo^s),q(ixsr^s),q(ixsl^s))
334 qmin(ixo^s)=min(q(ixo^s),q(ixsr^s),q(ixsl^s))
335 else
336 qmax(ixo^s)=max(qmax(ixo^s),q(ixsr^s),q(ixsl^s))
337 qmin(ixo^s)=min(qmin(ixo^s),q(ixsr^s),q(ixsl^s))
338 end if
339 {^nooned
340 idims=1
341 jdims=idims+1
342 do i=-1,1
343 ixs^l=ixo^l+i*ishift*kr(idims,^d);
344 ixsr^l=ixs^l+ishift*kr(jdims,^d);
345 ixsl^l=ixs^l-ishift*kr(jdims,^d);
346 qmax(ixo^s)=max(qmax(ixo^s),q(ixsr^s),q(ixsl^s))
347 qmin(ixo^s)=min(qmin(ixo^s),q(ixsr^s),q(ixsl^s))
348 end do
349 }
350 {^ifthreed
351 idims=1
352 jdims=idims+1
353 kdims=jdims+1
354 do i=-1,1
355 ixs^l=ixo^l+i*ishift*kr(idims,^d);
356 do j=-1,1
357 ixs^l=ixo^l+j*ishift*kr(jdims,^d);
358 ixsr^l=ixs^l+ishift*kr(kdims,^d);
359 ixsl^l=ixs^l-ishift*kr(kdims,^d);
360 qmax(ixo^s)=max(qmax(ixo^s),q(ixsr^s),q(ixsl^s))
361 qmin(ixo^s)=min(qmin(ixo^s),q(ixsr^s),q(ixsl^s))
362 end do
363 end do
364 }
365 enddo
366
367 end subroutine extremaq
368
369 subroutine extremaw(ixI^L,ixO^L,w,nshift,wMax,wMin)
371
372 integer,intent(in) :: ixi^l,ixo^l
373 double precision, intent(in) :: w(ixi^s,1:nw)
374 integer,intent(in) :: nshift
375
376 double precision, intent(out) :: wmax(ixi^s,1:nw_recon),wmin(ixi^s,1:nw_recon)
377
378 integer :: ixs^l,ixsr^l,ixsl^l,idims,jdims,kdims,ishift,i,j
379
380 do ishift=1,nshift
381 idims=1
382 ixsr^l=ixo^l+ishift*kr(idims,^d);
383 ixsl^l=ixo^l-ishift*kr(idims,^d);
384 if (ishift==1) then
385 wmax(ixo^s,1:nw_recon)= &
386 max(w(ixo^s,1:nw_recon),w(ixsr^s,1:nw_recon),w(ixsl^s,1:nw_recon))
387 wmin(ixo^s,1:nw_recon)= &
388 min(w(ixo^s,1:nw_recon),w(ixsr^s,1:nw_recon),w(ixsl^s,1:nw_recon))
389 else
390 wmax(ixo^s,1:nw_recon)= &
391 max(wmax(ixo^s,1:nw_recon),w(ixsr^s,1:nw_recon),w(ixsl^s,1:nw_recon))
392 wmin(ixo^s,1:nw_recon)= &
393 min(wmin(ixo^s,1:nw_recon),w(ixsr^s,1:nw_recon),w(ixsl^s,1:nw_recon))
394 end if
395 {^nooned
396 idims=1
397 jdims=idims+1
398 do i=-1,1
399 ixs^l=ixo^l+i*ishift*kr(idims,^d);
400 ixsr^l=ixs^l+ishift*kr(jdims,^d);
401 ixsl^l=ixs^l-ishift*kr(jdims,^d);
402 wmax(ixo^s,1:nw_recon)= &
403 max(wmax(ixo^s,1:nw_recon),w(ixsr^s,1:nw_recon),w(ixsl^s,1:nw_recon))
404 wmin(ixo^s,1:nw_recon)= &
405 min(wmin(ixo^s,1:nw_recon),w(ixsr^s,1:nw_recon),w(ixsl^s,1:nw_recon))
406 end do
407 }
408 {^ifthreed
409 idims=1
410 jdims=idims+1
411 kdims=jdims+1
412 do i=-1,1
413 ixs^l=ixo^l+i*ishift*kr(idims,^d);
414 do j=-1,1
415 ixs^l=ixo^l+j*ishift*kr(jdims,^d);
416 ixsr^l=ixs^l+ishift*kr(kdims,^d);
417 ixsl^l=ixs^l-ishift*kr(kdims,^d);
418 wmax(ixo^s,1:nw_recon)= &
419 max(wmax(ixo^s,1:nw_recon),w(ixsr^s,1:nw_recon),w(ixsl^s,1:nw_recon))
420 wmin(ixo^s,1:nw_recon)= &
421 min(wmin(ixo^s,1:nw_recon),w(ixsr^s,1:nw_recon),w(ixsl^s,1:nw_recon))
422 end do
423 end do
424 }
425 enddo
426
427 end subroutine extremaw
428
429end module mod_ppm
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, dimension(3, 3) kr
Kronecker delta tensor.
integer, parameter ndim
Number of spatial dimensions for grid variables.
double precision, dimension(:), allocatable, parameter d
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition mod_physics.t:4
double precision phys_gamma
Definition mod_physics.t:13
subroutine, public ppmlimitervar(ixil, ixl, idims, q, qct, qlc, qrc)
Definition mod_ppm.t:12
subroutine, public ppmlimiter(ixil, ixl, idims, w, wct, wlc, wrc)
Definition mod_ppm.t:103