46 character(len=*),
intent(in) :: namelim
106 write(*,*)
'Unknown limiter: ', namelim
107 call mpistop(
"No such limiter")
112 integer,
intent(in) :: typelim
114 select case (typelim)
130 subroutine dwlimiter2(dwC,ixI^L,ixC^L,idims,typelim,ldw,rdw,a2max)
134 integer,
intent(in) :: ixI^L, ixC^L, idims
135 double precision,
intent(in) :: dwC(ixi^s)
136 integer,
intent(in) :: typelim
138 double precision,
intent(out),
optional :: ldw(ixi^s)
140 double precision,
intent(out),
optional :: rdw(ixi^s)
141 double precision,
intent(in),
optional :: a2max
143 double precision :: tmp(ixi^s), tmp2(ixi^s)
144 integer :: ixO^L, hxO^L
145 double precision,
parameter :: qsmall=1.
d-12, qsmall2=2.
d-12
146 double precision,
parameter :: eps = sqrt(epsilon(1.0d0))
149 double precision,
parameter :: c_mcbeta=1.4d0
151 double precision,
parameter :: cadalfa=0.5d0, cadbeta=2.0d0, cadgamma=1.6d0
153 double precision :: rdelinv
154 double precision :: ldwA(ixi^s),ldwB(ixi^s),tmpeta(ixi^s)
155 double precision,
parameter :: cadepsilon=1.
d-14, invcadepsilon=1.d14,
cada3_radius=0.1d0
160 ixomin^d=ixcmin^d+
kr(idims,^d); ixomax^d=ixcmax^d;
161 hxo^l=ixo^l-
kr(idims,^d);
173 select case (typelim)
176 tmp(ixo^s)=sign(one,dwc(ixo^s))
177 tmp(ixo^s)=tmp(ixo^s)* &
178 max(zero,min(abs(dwc(ixo^s)),tmp(ixo^s)*dwc(hxo^s)))
179 if (
present(ldw)) ldw(ixo^s) = tmp(ixo^s)
180 if (
present(rdw)) rdw(ixo^s) = tmp(ixo^s)
183 tmp(ixo^s)=sign(one,dwc(ixo^s))
184 tmp(ixo^s)=2*tmp(ixo^s)* &
185 max(zero,min(abs(dwc(ixo^s)),tmp(ixo^s)*dwc(hxo^s),&
186 tmp(ixo^s)*quarter*(dwc(hxo^s)+dwc(ixo^s))))
187 if (
present(ldw)) ldw(ixo^s) = tmp(ixo^s)
188 if (
present(rdw)) rdw(ixo^s) = tmp(ixo^s)
191 tmp(ixo^s)=sign(one,dwc(ixo^s))
192 tmp(ixo^s)=tmp(ixo^s)* &
193 max(zero,min(c_mcbeta*abs(dwc(ixo^s)),c_mcbeta*tmp(ixo^s)*dwc(hxo^s),&
194 tmp(ixo^s)*half*(dwc(hxo^s)+dwc(ixo^s))))
195 if (
present(ldw)) ldw(ixo^s) = tmp(ixo^s)
196 if (
present(rdw)) rdw(ixo^s) = tmp(ixo^s)
199 tmp(ixo^s)=sign(one,dwc(ixo^s))
200 tmp(ixo^s)=tmp(ixo^s)* &
201 max(zero,min(2*abs(dwc(ixo^s)),tmp(ixo^s)*dwc(hxo^s)),&
202 min(abs(dwc(ixo^s)),2*tmp(ixo^s)*dwc(hxo^s)))
203 if (
present(ldw)) ldw(ixo^s) = tmp(ixo^s)
204 if (
present(rdw)) rdw(ixo^s) = tmp(ixo^s)
207 tmp(ixo^s)=2*max(dwc(hxo^s)*dwc(ixo^s),zero) &
208 /(dwc(ixo^s)+dwc(hxo^s)+qsmall)
209 if (
present(ldw)) ldw(ixo^s) = tmp(ixo^s)
210 if (
present(rdw)) rdw(ixo^s) = tmp(ixo^s)
213 tmp(ixo^s)=(dwc(hxo^s)*(dwc(ixo^s)**2+qsmall)&
214 +dwc(ixo^s)*(dwc(hxo^s)**2+qsmall))&
215 /(dwc(ixo^s)**2+dwc(hxo^s)**2+qsmall2)
216 if (
present(ldw)) ldw(ixo^s) = tmp(ixo^s)
217 if (
present(rdw)) rdw(ixo^s) = tmp(ixo^s)
219 tmp(ixo^s)=sign(one,dwc(ixo^s))
220 tmp2(ixo^s)=min(2*abs(dwc(ixo^s)),2*tmp(ixo^s)*dwc(hxo^s))
221 if (
present(ldw))
then 222 ldw(ixo^s)=tmp(ixo^s)* &
223 max(zero,min(tmp2(ixo^s),&
224 (dwc(hxo^s)*tmp(ixo^s)+2*abs(dwc(ixo^s)))*third))
226 if (
present(rdw))
then 227 rdw(ixo^s)=tmp(ixo^s)* &
228 max(zero,min(tmp2(ixo^s),&
229 (2*dwc(hxo^s)*tmp(ixo^s)+abs(dwc(ixo^s)))*third))
234 if (
present(ldw))
then 237 tmp(ixo^s)=dwc(hxo^s)/(dwc(ixo^s) + sign(eps, dwc(ixo^s)))
238 tmp2(ixo^s)=(2+tmp(ixo^s))*third
239 ldw(ixo^s)= max(zero,min(tmp2(ixo^s), &
240 max(-cadalfa*tmp(ixo^s), &
241 min(cadbeta*tmp(ixo^s), tmp2(ixo^s), &
242 cadgamma)))) * dwc(ixo^s)
245 if (
present(rdw))
then 247 tmp(ixo^s)=dwc(ixo^s)/(dwc(hxo^s) + sign(eps, dwc(hxo^s)))
248 tmp2(ixo^s)=(2+tmp(ixo^s))*third
249 rdw(ixo^s)= max(zero,min(tmp2(ixo^s), &
250 max(-cadalfa*tmp(ixo^s), &
251 min(cadbeta*tmp(ixo^s), tmp2(ixo^s), &
252 cadgamma)))) * dwc(hxo^s)
256 tmpeta(ixo^s)=(dwc(ixo^s)**2+dwc(hxo^s)**2)*rdelinv
257 if (
present(ldw))
then 258 tmp(ixo^s)=dwc(hxo^s)/(dwc(ixo^s) + sign(eps, dwc(ixo^s)))
259 ldwa(ixo^s)=(two+tmp(ixo^s))*third
260 where(tmpeta(ixo^s)<=one-cadepsilon)
261 ldw(ixo^s)=ldwa(ixo^s)
262 elsewhere(tmpeta(ixo^s)>=one+cadepsilon)
263 ldwb(ixo^s)= max(zero,min(ldwa(ixo^s), max(-cadalfa*tmp(ixo^s), &
264 min(cadbeta*tmp(ixo^s), ldwa(ixo^s), cadgamma))))
265 ldw(ixo^s)=ldwb(ixo^s)
267 ldwb(ixo^s)= max(zero,min(ldwa(ixo^s), max(-cadalfa*tmp(ixo^s), &
268 min(cadbeta*tmp(ixo^s), ldwa(ixo^s), cadgamma))))
269 tmp2(ixo^s)=(tmpeta(ixo^s)-one)*invcadepsilon
270 ldw(ixo^s)=half*( (one-tmp2(ixo^s))*ldwa(ixo^s) &
271 +(one+tmp2(ixo^s))*ldwb(ixo^s))
273 ldw(ixo^s)=ldw(ixo^s) * dwc(ixo^s)
276 if (
present(rdw))
then 277 tmp(ixo^s)=dwc(ixo^s)/(dwc(hxo^s) + sign(eps, dwc(hxo^s)))
278 ldwa(ixo^s)=(two+tmp(ixo^s))*third
279 where(tmpeta(ixo^s)<=one-cadepsilon)
280 rdw(ixo^s)=ldwa(ixo^s)
281 elsewhere(tmpeta(ixo^s)>=one+cadepsilon)
282 ldwb(ixo^s)= max(zero,min(ldwa(ixo^s), max(-cadalfa*tmp(ixo^s), &
283 min(cadbeta*tmp(ixo^s), ldwa(ixo^s), cadgamma))))
284 rdw(ixo^s)=ldwb(ixo^s)
286 ldwb(ixo^s)= max(zero,min(ldwa(ixo^s), max(-cadalfa*tmp(ixo^s), &
287 min(cadbeta*tmp(ixo^s), ldwa(ixo^s), cadgamma))))
288 tmp2(ixo^s)=(tmpeta(ixo^s)-one)*invcadepsilon
289 rdw(ixo^s)=half*( (one-tmp2(ixo^s))*ldwa(ixo^s) &
290 +(one+tmp2(ixo^s))*ldwb(ixo^s))
292 rdw(ixo^s)=rdw(ixo^s) * dwc(hxo^s)
295 tmpeta(ixo^s)=(sqrt(0.4d0*(dwc(ixo^s)**2+dwc(hxo^s)**2)))&
296 /((a2max+cadepsilon)*
dxlevel(idims)**2)
297 if(
present(ldw))
then 298 tmp(ixo^s)=dwc(hxo^s)/(dwc(ixo^s)+sign(eps,dwc(ixo^s)))
299 ldwa(ixo^s)=(two+tmp(ixo^s))*third
300 where(tmpeta(ixo^s)<=one-cadepsilon)
301 ldw(ixo^s)=ldwa(ixo^s)
302 else where(tmpeta(ixo^s)>=one+cadepsilon)
303 ldwb(ixo^s)=max(zero,min(ldwa(ixo^s),max(-tmp(ixo^s),&
304 min(cadbeta*tmp(ixo^s),ldwa(ixo^s),1.5d0))))
305 ldw(ixo^s)=ldwb(ixo^s)
307 ldwb(ixo^s)=max(zero,min(ldwa(ixo^s),max(-tmp(ixo^s),&
308 min(cadbeta*tmp(ixo^s),ldwa(ixo^s),1.5d0))))
309 tmp2(ixo^s)=(tmpeta(ixo^s)-one)*invcadepsilon
310 ldw(ixo^s)=half*((one-tmp2(ixo^s))*ldwa(ixo^s)&
311 +(one+tmp2(ixo^s))*ldwb(ixo^s))
313 ldw(ixo^s)=ldw(ixo^s)*dwc(ixo^s)
315 if(
present(rdw))
then 316 tmp(ixo^s)=dwc(ixo^s)/(dwc(hxo^s)+sign(eps,dwc(hxo^s)))
317 ldwa(ixo^s)=(two+tmp(ixo^s))*third
318 where(tmpeta(ixo^s)<=one-cadepsilon)
319 rdw(ixo^s)=ldwa(ixo^s)
320 else where(tmpeta(ixo^s)>=one+cadepsilon)
321 ldwb(ixo^s)=max(zero,min(ldwa(ixo^s),max(-tmp(ixo^s),&
322 min(cadbeta*tmp(ixo^s),ldwa(ixo^s),1.5d0))))
323 rdw(ixo^s)=ldwb(ixo^s)
325 ldwb(ixo^s)=max(zero,min(ldwa(ixo^s),max(-tmp(ixo^s),&
326 min(cadbeta*tmp(ixo^s), ldwa(ixo^s), 1.5d0))))
327 tmp2(ixo^s)=(tmpeta(ixo^s)-one)*invcadepsilon
328 rdw(ixo^s)=half*((one-tmp2(ixo^s))*ldwa(ixo^s)&
329 +(one+tmp2(ixo^s))*ldwb(ixo^s))
331 rdw(ixo^s)=rdw(ixo^s)*dwc(hxo^s)
334 call mpistop(
"Error in dwLimiter: unknown limiter")
This module contains definitions of global parameters and variables and some generic functions/subrou...
pure logical function limiter_symmetric(typelim)
integer, parameter limiter_superbee
integer, parameter limiter_cada3
integer, parameter limiter_wenoz5
integer, parameter limiter_wenozp5
integer, parameter limiter_mp5
integer, parameter limiter_albada
integer, parameter limiter_wenoyc3
integer, parameter limiter_mpweno7
Module containing the MP5 (fifth order) flux scheme.
integer, parameter limiter_teno5ad
integer, parameter limiter_woodward
integer, parameter limiter_weno7
integer, parameter limiter_ppm
Module with slope/flux limiters.
integer, parameter limiter_schmid
integer, parameter limiter_wenozp5nm
integer, parameter limiter_weno5nm
integer function limiter_type(namelim)
integer, parameter limiter_koren
double precision cada3_radius
radius of the asymptotic region [0.001, 10], larger means more accurate in smooth region but more ove...
integer, dimension(3, 3) kr
Kronecker delta tensor.
integer, dimension(:), allocatable, parameter d
integer, parameter limiter_weno5cu6
integer, parameter limiter_mcbeta
integer, parameter limiter_venk
subroutine dwlimiter2(dwC, ixIL, ixCL, idims, typelim, ldw, rdw, a2max)
Limit the centered dwC differences within ixC for iw in direction idim. The limiter is chosen accordi...
integer, parameter limiter_exeno7
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
integer, parameter limiter_weno5
double precision, dimension(ndim) dxlevel
double precision schmid_rad
integer, parameter limiter_vanleer
integer, parameter limiter_minmod
integer, parameter limiter_cada
integer, parameter limiter_weno3
integer, parameter limiter_wenoz5nm