[pyfr] 13/88: Rework artificial viscosity.
Ghislain Vaillant
ghisvail-guest at moszumanska.debian.org
Wed Nov 16 12:05:25 UTC 2016
This is an automated email from the git hooks/post-receive script.
ghisvail-guest pushed a commit to branch master
in repository pyfr.
commit 0138e8d5bd07a28f3bfb43d2938cbd4a25c730f1
Author: Freddie Witherden <freddie at witherden.org>
Date: Sun Apr 17 17:15:01 2016 -0700
Rework artificial viscosity.
Support is now present at the BaseAdvectionDiffusion level as
opposed to just Navier-Stokes. Further, sensing in the case of
Navier-Stokes is now based on the density field as opposed to the
entropy. Additionally, artificial viscosity is treated
separately from the physical viscosity (if any). For Navier-Stokes
this means that the density field now has viscosity applied to it.
---
pyfr/solvers/base/elements.py | 8 --
pyfr/solvers/baseadvecdiff/elements.py | 92 ++++++++++++++++-
pyfr/solvers/baseadvecdiff/inters.py | 44 ++++++++
pyfr/solvers/baseadvecdiff/kernels/artvisc.mako | 10 ++
.../kernels/shocksensor.mako} | 3 +-
pyfr/solvers/baseadvecdiff/kernels/shockvar.mako | 9 ++
pyfr/solvers/baseadvecdiff/system.py | 13 ++-
pyfr/solvers/navstokes/elements.py | 87 ++--------------
pyfr/solvers/navstokes/inters.py | 115 ++++++---------------
pyfr/solvers/navstokes/kernels/bcs/ghost.mako | 4 +-
pyfr/solvers/navstokes/kernels/entropy.mako | 15 ---
pyfr/solvers/navstokes/kernels/flux.mako | 12 +--
pyfr/solvers/navstokes/kernels/intcflux.mako | 7 +-
pyfr/solvers/navstokes/kernels/mpicflux.mako | 7 +-
pyfr/solvers/navstokes/kernels/tflux.mako | 5 +-
pyfr/solvers/navstokes/system.py | 64 ------------
16 files changed, 225 insertions(+), 270 deletions(-)
diff --git a/pyfr/solvers/base/elements.py b/pyfr/solvers/base/elements.py
index bc06cde..27cf3c8 100644
--- a/pyfr/solvers/base/elements.py
+++ b/pyfr/solvers/base/elements.py
@@ -360,14 +360,6 @@ class BaseElements(object, metaclass=ABCMeta):
return (self._vect_fpts.mid,)*nfp, rcmap, rcstri
- def get_avis_fpts_for_inter(self, eidx, fidx):
- nfp = self.nfacefpts[fidx]
-
- rcmap = [(fpidx, eidx) for fpidx in self._srtd_face_fpts[fidx][eidx]]
- cstri = ((self._avis_fpts.leadsubdim,),)*nfp
-
- return (self._avis_fpts.mid,)*nfp, rcmap, cstri
-
def get_ploc_for_inter(self, eidx, fidx):
fpts_idx = self._srtd_face_fpts[fidx][eidx]
return self.plocfpts[fpts_idx,eidx]
diff --git a/pyfr/solvers/baseadvecdiff/elements.py b/pyfr/solvers/baseadvecdiff/elements.py
index afb3da6..437a37b 100644
--- a/pyfr/solvers/baseadvecdiff/elements.py
+++ b/pyfr/solvers/baseadvecdiff/elements.py
@@ -1,7 +1,9 @@
# -*- coding: utf-8 -*-
-from pyfr.solvers.baseadvec import BaseAdvectionElements
+import numpy as np
+
from pyfr.backends.base.kernels import ComputeMetaKernel
+from pyfr.solvers.baseadvec import BaseAdvectionElements
class BaseAdvectionDiffusionElements(BaseAdvectionElements):
@@ -22,6 +24,79 @@ class BaseAdvectionDiffusionElements(BaseAdvectionElements):
return bufs
+ def _set_backend_art_visc(self, backend):
+ nvars, neles = self.nvars, self.neles
+ nupts, nfpts = self.nupts, self.nfpts
+ tags = {'align'}
+
+ # Register pointwise kernels
+ backend.pointwise.register(
+ 'pyfr.solvers.baseadvecdiff.kernels.shockvar'
+ )
+ backend.pointwise.register(
+ 'pyfr.solvers.baseadvecdiff.kernels.shocksensor'
+ )
+
+ # Obtain the scalar variable to be used for shock sensing
+ shockvar = self.convarmap[self.ndims].index(self.shockvar)
+
+ # Common template arguments
+ tplargs = dict(
+ nvars=nvars, nupts=nupts, nfpts=nfpts,
+ c=self.cfg.items_as('solver-artificial-viscosity', float),
+ order=self.basis.order, ubdegs=self.basis.ubasis.degrees,
+ shockvar=shockvar
+ )
+
+ # Allocate required scratch space for artificial viscosity
+ if 'flux' in self.antialias:
+ self.avis_qpts = backend.matrix((self.nqpts, 1, neles),
+ extent='avis_qpts', tags=tags)
+ self.avis_upts = backend.matrix((nupts, 1, neles),
+ aliases=self.avis_qpts, tags=tags)
+ else:
+ self.avis_upts = backend.matrix((nupts, 1, neles),
+ extent='avis_upts', tags=tags)
+
+ if nfpts >= nupts:
+ self._avis_fpts = backend.matrix((nfpts, 1, neles),
+ extent='avis_fpts', tags=tags)
+ avis_upts_temp = backend.matrix(
+ (nupts, 1, neles), aliases=self._avis_fpts, tags=tags
+ )
+ else:
+ avis_upts_temp = backend.matrix((nupts, 1, neles),
+ extent='avis_fpts', tags=tags)
+ self._avis_fpts = backend.matrix(
+ (nfpts, 1, neles), aliases=avis_upts_temp, tags=tags
+ )
+
+ # Extract the scalar variable to be used for shock sensing
+ self.kernels['shockvar'] = lambda: backend.kernel(
+ 'shockvar', tplargs=tplargs, dims=[nupts, neles],
+ u=self.scal_upts_inb, s=self.avis_upts
+ )
+
+ # Obtain the modal coefficients
+ rcpvdm = np.linalg.inv(self.basis.ubasis.vdm.T)
+ rcpvdm = backend.const_matrix(rcpvdm, tags={'align'})
+ self.kernels['shockvar_modal'] = lambda: backend.kernel(
+ 'mul', rcpvdm, self.avis_upts, out=avis_upts_temp
+ )
+
+ if 'flux' in self.antialias:
+ ame_e = self.avis_qpts
+ tplargs['nrow_amu'] = self.nqpts
+ else:
+ ame_e = self.avis_upts
+ tplargs['nrow_amu'] = nupts
+
+ # Apply the sensor to estimate the required artificial viscosity
+ self.kernels['shocksensor'] = lambda: backend.kernel(
+ 'shocksensor', tplargs=tplargs, dims=[neles],
+ s=avis_upts_temp, amu_e=ame_e, amu_f=self._avis_fpts
+ )
+
def set_backend(self, backend, nscal_upts):
super().set_backend(backend, nscal_upts)
@@ -72,3 +147,18 @@ class BaseAdvectionDiffusionElements(BaseAdvectionElements):
return ComputeMetaKernel(muls)
self.kernels['gradcoru_qpts'] = gradcoru_qpts
+
+ # Shock capturing
+ shock_capturing = self.cfg.get('solver', 'shock-capturing', 'none')
+ if shock_capturing == 'artificial-viscosity':
+ self._set_backend_art_visc(backend)
+ elif shock_capturing == 'none':
+ self.avis_qpts = self.avis_upts = None
+ else:
+ raise ValueError('Invalid shock capturing scheme')
+
+ def get_avis_fpts_for_inter(self, eidx, fidx):
+ nfp = self.nfacefpts[fidx]
+
+ rcmap = [(fpidx, eidx) for fpidx in self._srtd_face_fpts[fidx][eidx]]
+ return (self._avis_fpts.mid,)*nfp, rcmap
\ No newline at end of file
diff --git a/pyfr/solvers/baseadvecdiff/inters.py b/pyfr/solvers/baseadvecdiff/inters.py
index 873fa64..f84ac2f 100644
--- a/pyfr/solvers/baseadvecdiff/inters.py
+++ b/pyfr/solvers/baseadvecdiff/inters.py
@@ -15,6 +15,13 @@ class BaseAdvectionDiffusionIntInters(BaseAdvectionIntInters):
self._vect0_lhs = self._vect_view(lhs, 'get_vect_fpts_for_inter')
self._vect0_rhs = self._vect_view(rhs, 'get_vect_fpts_for_inter')
+ # Generate the additional view matrices for artificial viscosity
+ if cfg.get('solver', 'shock-capturing') == 'artificial-viscosity':
+ self._avis0_lhs = self._view(lhs, 'get_avis_fpts_for_inter')
+ self._avis0_rhs = self._view(rhs, 'get_avis_fpts_for_inter')
+ else:
+ self._avis0_lhs = self._avis0_rhs = None
+
# Additional kernel constants
self._tpl_c.update(cfg.items_as('solver-interfaces', float))
@@ -76,6 +83,37 @@ class BaseAdvectionDiffusionMPIInters(BaseAdvectionMPIInters):
self.kernels['vect_fpts_recv'] = lambda: NullMPIKernel()
self.kernels['vect_fpts_unpack'] = lambda: NullComputeKernel()
+ # Generate the additional kernels/views for artificial viscosity
+ if cfg.get('solver', 'shock-capturing') == 'artificial-viscosity':
+ self._avis0_lhs = self._xchg_view(lhs, 'get_avis_fpts_for_inter')
+ self._avis0_rhs = be.xchg_matrix_for_view(self._avis0_lhs)
+
+ # If we need to send our artificial viscosity to the RHS
+ if self._tpl_c['ldg-beta'] != -0.5:
+ self.kernels['avis_fpts_pack'] = lambda: be.kernel(
+ 'pack', self._avis0_lhs
+ )
+ self.kernels['avis_fpts_send'] = lambda: be.kernel(
+ 'send_pack', self._avis0_lhs, self._rhsrank, self.MPI_TAG
+ )
+ else:
+ self.kernels['avis_fpts_pack'] = lambda: NullComputeKernel()
+ self.kernels['avis_fpts_send'] = lambda: NullMPIKernel()
+
+ # If we need to recv artificial viscosity from the RHS
+ if self._tpl_c['ldg-beta'] != 0.5:
+ self.kernels['avis_fpts_recv'] = lambda: be.kernel(
+ 'recv_pack', self._avis0_rhs, self._rhsrank, self.MPI_TAG
+ )
+ self.kernels['avis_fpts_unpack'] = lambda: be.kernel(
+ 'unpack', self._avis0_rhs
+ )
+ else:
+ self.kernels['avis_fpts_recv'] = lambda: NullMPIKernel()
+ self.kernels['avis_fpts_unpack'] = lambda: NullComputeKernel()
+ else:
+ self._avis0_lhs = self._avis0_rhs = None
+
class BaseAdvectionDiffusionBCInters(BaseAdvectionBCInters):
def __init__(self, be, lhs, elemap, cfgsect, cfg):
@@ -86,3 +124,9 @@ class BaseAdvectionDiffusionBCInters(BaseAdvectionBCInters):
# Additional kernel constants
self._tpl_c.update(cfg.items_as('solver-interfaces', float))
+
+ # Generate the additional view matrices for artificial viscosity
+ if cfg.get('solver', 'shock-capturing') == 'artificial-viscosity':
+ self._avis0_lhs = self._view(lhs, 'get_avis_fpts_for_inter')
+ else:
+ self._avis0_lhs = None
diff --git a/pyfr/solvers/baseadvecdiff/kernels/artvisc.mako b/pyfr/solvers/baseadvecdiff/kernels/artvisc.mako
new file mode 100644
index 0000000..1aa8e4a
--- /dev/null
+++ b/pyfr/solvers/baseadvecdiff/kernels/artvisc.mako
@@ -0,0 +1,10 @@
+# -*- coding: utf-8 -*-
+<%namespace module='pyfr.backends.base.makoutil' name='pyfr'/>
+
+<%pyfr:macro name='artificial_viscosity_add' params='grad_uin, fout, artvisc'>
+% if shock_capturing == 'artificial-viscosity':
+% for i, j in pyfr.ndrange(ndims, nvars):
+ fout[${i}][${j}] -= artvisc*grad_uin[${i}][${j}];
+% endfor
+% endif
+</%pyfr:macro>
diff --git a/pyfr/solvers/navstokes/kernels/avis.mako b/pyfr/solvers/baseadvecdiff/kernels/shocksensor.mako
similarity index 96%
rename from pyfr/solvers/navstokes/kernels/avis.mako
rename to pyfr/solvers/baseadvecdiff/kernels/shocksensor.mako
index b9bc178..676c162 100644
--- a/pyfr/solvers/navstokes/kernels/avis.mako
+++ b/pyfr/solvers/baseadvecdiff/kernels/shocksensor.mako
@@ -7,8 +7,7 @@
pi = math.pi
%>
-
-<%pyfr:kernel name='avis' ndim='1'
+<%pyfr:kernel name='shocksensor' ndim='1'
s='in fpdtype_t[${str(nupts)}]'
amu_e='out fpdtype_t[${str(nrow_amu)}]'
amu_f='out fpdtype_t[${str(nfpts)}]'>
diff --git a/pyfr/solvers/baseadvecdiff/kernels/shockvar.mako b/pyfr/solvers/baseadvecdiff/kernels/shockvar.mako
new file mode 100644
index 0000000..affb8ed
--- /dev/null
+++ b/pyfr/solvers/baseadvecdiff/kernels/shockvar.mako
@@ -0,0 +1,9 @@
+# -*- coding: utf-8 -*-
+<%inherit file='base'/>
+<%namespace module='pyfr.backends.base.makoutil' name='pyfr'/>
+
+<%pyfr:kernel name='shockvar' ndim='2'
+ u='in fpdtype_t[${str(nvars)}]'
+ s='out fpdtype_t'>
+ s = u[${shockvar}];
+</%pyfr:kernel>
diff --git a/pyfr/solvers/baseadvecdiff/system.py b/pyfr/solvers/baseadvecdiff/system.py
index ca952de..dee45d3 100644
--- a/pyfr/solvers/baseadvecdiff/system.py
+++ b/pyfr/solvers/baseadvecdiff/system.py
@@ -20,8 +20,12 @@ class BaseAdvectionDiffusionSystem(BaseAdvectionSystem):
q1 << kernels['eles', 'copy_soln']()
q1 << kernels['iint', 'con_u']()
q1 << kernels['bcint', 'con_u'](t=t)
+ if ('eles', 'shockvar') in kernels:
+ q1 << kernels['eles', 'shockvar']()
+ q1 << kernels['eles', 'shockvar_modal']()
+ q1 << kernels['eles', 'shocksensor']()
+ q1 << kernels['mpiint', 'avis_fpts_pack']()
q1 << kernels['eles', 'tgradpcoru_upts']()
-
q2 << kernels['mpiint', 'scal_fpts_send']()
q2 << kernels['mpiint', 'scal_fpts_recv']()
q2 << kernels['mpiint', 'scal_fpts_unpack']()
@@ -33,7 +37,12 @@ class BaseAdvectionDiffusionSystem(BaseAdvectionSystem):
q1 << kernels['eles', 'gradcoru_upts']()
q1 << kernels['eles', 'gradcoru_fpts']()
q1 << kernels['mpiint', 'vect_fpts_pack']()
- runall([q1])
+ if ('eles', 'shockvar') in kernels:
+ q2 << kernels['mpiint', 'avis_fpts_send']()
+ q2 << kernels['mpiint', 'avis_fpts_recv']()
+ q2 << kernels['mpiint', 'avis_fpts_unpack']()
+
+ runall([q1, q2])
if ('eles', 'gradcoru_qpts') in kernels:
q1 << kernels['eles', 'gradcoru_qpts']()
diff --git a/pyfr/solvers/navstokes/elements.py b/pyfr/solvers/navstokes/elements.py
index 12a6135..b5994cc 100644
--- a/pyfr/solvers/navstokes/elements.py
+++ b/pyfr/solvers/navstokes/elements.py
@@ -1,108 +1,35 @@
# -*- coding: utf-8 -*-
-import numpy as np
-
from pyfr.solvers.baseadvecdiff import BaseAdvectionDiffusionElements
from pyfr.solvers.euler.elements import BaseFluidElements
class NavierStokesElements(BaseFluidElements, BaseAdvectionDiffusionElements):
+ # Use the density field for shock sensing
+ shockvar = 'rho'
+
def set_backend(self, backend, nscalupts):
super().set_backend(backend, nscalupts)
backend.pointwise.register('pyfr.solvers.navstokes.kernels.tflux')
+ shock_capturing = self.cfg.get('solver', 'shock-capturing')
visc_corr = self.cfg.get('solver', 'viscosity-correction', 'none')
if visc_corr not in {'sutherland', 'none'}:
raise ValueError('Invalid viscosity-correction option')
tplargs = dict(ndims=self.ndims, nvars=self.nvars,
- visc_corr=visc_corr,
+ shock_capturing=shock_capturing, visc_corr=visc_corr,
c=self.cfg.items_as('constants', float))
- shock_capturing = self.cfg.get('solver', 'shock-capturing', 'none')
- if shock_capturing == 'artificial-viscosity':
- # Allocate required scratch space for artificial viscosity
- nupts, nfpts, neles = self.nupts, self.nfpts, self.neles
- tags = {'align'}
-
- if 'flux' in self.antialias:
- avis_qpts = backend.matrix((self.nqpts, 1, neles),
- extent='avis_qpts', tags=tags)
- avis_upts = backend.matrix((nupts, 1, neles),
- aliases=avis_qpts, tags=tags)
- else:
- avis_upts = backend.matrix((nupts, 1, neles),
- extent='avis_upts', tags=tags)
-
- if nfpts >= nupts:
- self._avis_fpts = backend.matrix((nfpts, 1, neles),
- extent='avis_fpts', tags=tags)
- avis_upts_temp = backend.matrix(
- (nupts, 1, neles), aliases=self._avis_fpts, tags=tags
- )
- else:
- avis_upts_temp = backend.matrix((nupts, 1, neles),
- extent='avis_fpts', tags=tags)
- self._avis_fpts = backend.matrix(
- (nfpts, 1, neles), aliases=avis_upts_temp, tags=tags
- )
-
- # Entropy kernel
- backend.pointwise.register(
- 'pyfr.solvers.navstokes.kernels.entropy'
- )
-
- self.kernels['entropy'] = lambda: backend.kernel(
- 'entropy', tplargs=tplargs, dims=[nupts, neles],
- u=self.scal_upts_inb, s=avis_upts
- )
-
- # Modal entropy coefficient kernel
- rcpvdm = np.linalg.inv(self.basis.ubasis.vdm.T)
- rcpvdm = self._be.const_matrix(rcpvdm, tags={'align'})
- self.kernels['modal_entropy'] = lambda: backend.kernel(
- 'mul', rcpvdm, avis_upts, out=avis_upts_temp
- )
-
- # Artificial viscosity kernel
- backend.pointwise.register('pyfr.solvers.navstokes.kernels.avis')
- art_visc_tplargs = dict(
- tplargs, nupts=nupts, nfpts=nfpts, order=self.basis.order,
- ubdegs=self.basis.ubasis.degrees
- )
- art_visc_tplargs['c'].update(
- self.cfg.items_as('solver-artificial-viscosity', float)
- )
-
- if 'flux' in self.antialias:
- ame_e = avis_qpts
- art_visc_tplargs['nrow_amu'] = self.nqpts
- else:
- ame_e = avis_upts
- art_visc_tplargs['nrow_amu'] = nupts
-
- # Element-wise artificial viscosity kernel
- self.kernels['art_visc'] = lambda: backend.kernel(
- 'avis', tplargs=art_visc_tplargs, dims=[neles],
- s=avis_upts_temp, amu_e=ame_e, amu_f=self._avis_fpts
- )
-
- tplargs['art_vis'] = 'mu'
- elif shock_capturing == 'none':
- avis_upts = avis_qpts = None
- tplargs['art_vis'] = 'none'
- else:
- raise ValueError('Invalid shock-capturing option')
-
if 'flux' in self.antialias:
self.kernels['tdisf'] = lambda: backend.kernel(
'tflux', tplargs=tplargs, dims=[self.nqpts, self.neles],
u=self._scal_qpts, smats=self.smat_at('qpts'),
- f=self._vect_qpts, amu=avis_qpts
+ f=self._vect_qpts, amu=self.avis_qpts
)
else:
self.kernels['tdisf'] = lambda: backend.kernel(
'tflux', tplargs=tplargs, dims=[self.nupts, self.neles],
u=self.scal_upts_inb, smats=self.smat_at('upts'),
- f=self._vect_upts, amu=avis_upts
+ f=self._vect_upts, amu=self.avis_upts
)
diff --git a/pyfr/solvers/navstokes/inters.py b/pyfr/solvers/navstokes/inters.py
index 7cad8ee..d86a438 100644
--- a/pyfr/solvers/navstokes/inters.py
+++ b/pyfr/solvers/navstokes/inters.py
@@ -2,7 +2,6 @@
import numpy as np
-from pyfr.backends.base import NullComputeKernel, NullMPIKernel
from pyfr.solvers.baseadvecdiff import (BaseAdvectionDiffusionBCInters,
BaseAdvectionDiffusionIntInters,
BaseAdvectionDiffusionMPIInters)
@@ -14,33 +13,25 @@ class NavierStokesIntInters(BaseAdvectionDiffusionIntInters):
# Pointwise template arguments
rsolver = self.cfg.get('solver-interfaces', 'riemann-solver')
- visc_corr = self.cfg.get('solver', 'viscosity-correction', 'none')
+ visc_corr = self.cfg.get('solver', 'viscosity-correction')
+ shock_capturing = self.cfg.get('solver', 'shock-capturing')
tplargs = dict(ndims=self.ndims, nvars=self.nvars, rsolver=rsolver,
- visc_corr=visc_corr, c=self._tpl_c)
-
- # Generate the additional view matrices for artificial viscosity
- shock_capturing = self.cfg.get('solver', 'shock-capturing', 'none')
- if shock_capturing == 'artificial-viscosity':
- avis0_lhs = self._view(lhs, 'get_avis_fpts_for_inter')
- avis0_rhs = self._view(rhs, 'get_avis_fpts_for_inter')
- tplargs['art_vis'] = 'mu'
- else:
- avis0_lhs = avis0_rhs = None
- tplargs['art_vis'] = 'none'
-
- self._be.pointwise.register('pyfr.solvers.navstokes.kernels.intconu')
- self._be.pointwise.register('pyfr.solvers.navstokes.kernels.intcflux')
-
- self.kernels['con_u'] = lambda: self._be.kernel(
+ visc_corr=visc_corr, shock_capturing=shock_capturing,
+ c=self._tpl_c)
+
+ be.pointwise.register('pyfr.solvers.navstokes.kernels.intconu')
+ be.pointwise.register('pyfr.solvers.navstokes.kernels.intcflux')
+
+ self.kernels['con_u'] = lambda: be.kernel(
'intconu', tplargs=tplargs, dims=[self.ninterfpts],
ulin=self._scal0_lhs, urin=self._scal0_rhs,
ulout=self._vect0_lhs, urout=self._vect0_rhs
)
- self.kernels['comm_flux'] = lambda: self._be.kernel(
+ self.kernels['comm_flux'] = lambda: be.kernel(
'intcflux', tplargs=tplargs, dims=[self.ninterfpts],
ul=self._scal0_lhs, ur=self._scal0_rhs,
gradul=self._vect0_lhs, gradur=self._vect0_rhs,
- amul=avis0_lhs, amur=avis0_rhs,
+ amul=self._avis0_lhs, amur=self._avis0_rhs,
magnl=self._mag_pnorm_lhs, nl=self._norm_pnorm_lhs
)
@@ -50,58 +41,25 @@ class NavierStokesMPIInters(BaseAdvectionDiffusionMPIInters):
super().__init__(be, lhs, rhsrank, rallocs, elemap, cfg)
# Pointwise template arguments
- rsolver = self.cfg.get('solver-interfaces', 'riemann-solver')
- visc_corr = self.cfg.get('solver', 'viscosity-correction', 'none')
+ rsolver = cfg.get('solver-interfaces', 'riemann-solver')
+ visc_corr = cfg.get('solver', 'viscosity-correction')
+ shock_capturing = cfg.get('solver', 'shock-capturing')
tplargs = dict(ndims=self.ndims, nvars=self.nvars, rsolver=rsolver,
- visc_corr=visc_corr, c=self._tpl_c)
-
- # Generate the additional kernels/views for artificial viscosity
- shock_capturing = self.cfg.get('solver', 'shock-capturing', 'none')
- if shock_capturing == 'artificial-viscosity':
- avis0_lhs = self._xchg_view(lhs, 'get_avis_fpts_for_inter')
- avis0_rhs = be.xchg_matrix_for_view(avis0_lhs)
-
- # If we need to send our artificial viscosity to the RHS
- if self._tpl_c['ldg-beta'] != -0.5:
- self.kernels['avis_fpts_pack'] = lambda: be.kernel(
- 'pack', avis0_lhs
- )
- self.kernels['avis_fpts_send'] = lambda: be.kernel(
- 'send_pack', avis0_lhs, self._rhsrank, self.MPI_TAG
- )
- else:
- self.kernels['avis_fpts_pack'] = lambda: NullComputeKernel()
- self.kernels['avis_fpts_send'] = lambda: NullMPIKernel()
-
- # If we need to recv artificial viscosity from the RHS
- if self._tpl_c['ldg-beta'] != 0.5:
- self.kernels['avis_fpts_recv'] = lambda: be.kernel(
- 'recv_pack', avis0_rhs, self._rhsrank, self.MPI_TAG
- )
- self.kernels['avis_fpts_unpack'] = lambda: be.kernel(
- 'unpack', avis0_rhs
- )
- else:
- self.kernels['avis_fpts_recv'] = lambda: NullMPIKernel()
- self.kernels['avis_fpts_unpack'] = lambda: NullComputeKernel()
-
- tplargs['art_vis'] = 'mu'
- else:
- avis0_lhs = avis0_rhs = None
- tplargs['art_vis'] = 'none'
-
- self._be.pointwise.register('pyfr.solvers.navstokes.kernels.mpiconu')
- self._be.pointwise.register('pyfr.solvers.navstokes.kernels.mpicflux')
-
- self.kernels['con_u'] = lambda: self._be.kernel(
+ visc_corr=visc_corr, shock_capturing=shock_capturing,
+ c=self._tpl_c)
+
+ be.pointwise.register('pyfr.solvers.navstokes.kernels.mpiconu')
+ be.pointwise.register('pyfr.solvers.navstokes.kernels.mpicflux')
+
+ self.kernels['con_u'] = lambda: be.kernel(
'mpiconu', tplargs=tplargs, dims=[self.ninterfpts],
ulin=self._scal0_lhs, urin=self._scal0_rhs, ulout=self._vect0_lhs
)
- self.kernels['comm_flux'] = lambda: self._be.kernel(
+ self.kernels['comm_flux'] = lambda: be.kernel(
'mpicflux', tplargs=tplargs, dims=[self.ninterfpts],
ul=self._scal0_lhs, ur=self._scal0_rhs,
gradul=self._vect0_lhs, gradur=self._vect0_rhs,
- amul=avis0_lhs, amur=avis0_rhs,
+ amul=self._avis0_lhs, amur=self._avis0_rhs,
magnl=self._mag_pnorm_lhs, nl=self._norm_pnorm_lhs
)
@@ -113,32 +71,25 @@ class NavierStokesBaseBCInters(BaseAdvectionDiffusionBCInters):
super().__init__(be, lhs, elemap, cfgsect, cfg)
# Pointwise template arguments
- rsolver = self.cfg.get('solver-interfaces', 'riemann-solver')
- visc_corr = self.cfg.get('solver', 'viscosity-correction', 'none')
+ rsolver = cfg.get('solver-interfaces', 'riemann-solver')
+ visc_corr = cfg.get('solver', 'viscosity-correction', 'none')
+ shock_capturing = cfg.get('solver', 'shock-capturing', 'none')
tplargs = dict(ndims=self.ndims, nvars=self.nvars, rsolver=rsolver,
- visc_corr=visc_corr, c=self._tpl_c, bctype=self.type,
+ visc_corr=visc_corr, shock_capturing=shock_capturing,
+ c=self._tpl_c, bctype=self.type,
bccfluxstate=self.cflux_state)
- # Generate the additional view matrices for artificial viscosity
- shock_capturing = self.cfg.get('solver', 'shock-capturing', 'none')
- if shock_capturing == 'artificial-viscosity':
- avis0_lhs = self._view(lhs, 'get_avis_fpts_for_inter')
- tplargs['art_vis'] = 'mu'
- else:
- avis0_lhs = None
- tplargs['art_vis'] = 'none'
-
- self._be.pointwise.register('pyfr.solvers.navstokes.kernels.bcconu')
- self._be.pointwise.register('pyfr.solvers.navstokes.kernels.bccflux')
+ be.pointwise.register('pyfr.solvers.navstokes.kernels.bcconu')
+ be.pointwise.register('pyfr.solvers.navstokes.kernels.bccflux')
- self.kernels['con_u'] = lambda: self._be.kernel(
+ self.kernels['con_u'] = lambda: be.kernel(
'bcconu', tplargs=tplargs, dims=[self.ninterfpts],
ulin=self._scal0_lhs, ulout=self._vect0_lhs,
nlin=self._norm_pnorm_lhs, ploc=self._ploc
)
- self.kernels['comm_flux'] = lambda: self._be.kernel(
+ self.kernels['comm_flux'] = lambda: be.kernel(
'bccflux', tplargs=tplargs, dims=[self.ninterfpts],
- ul=self._scal0_lhs, gradul=self._vect0_lhs, amul=avis0_lhs,
+ ul=self._scal0_lhs, gradul=self._vect0_lhs, amul=self._avis0_lhs,
magnl=self._mag_pnorm_lhs, nl=self._norm_pnorm_lhs,
ploc=self._ploc
)
diff --git a/pyfr/solvers/navstokes/kernels/bcs/ghost.mako b/pyfr/solvers/navstokes/kernels/bcs/ghost.mako
index 24c69d2..38ccd64 100644
--- a/pyfr/solvers/navstokes/kernels/bcs/ghost.mako
+++ b/pyfr/solvers/navstokes/kernels/bcs/ghost.mako
@@ -1,5 +1,6 @@
<%namespace module='pyfr.backends.base.makoutil' name='pyfr'/>
+<%include file='pyfr.solvers.baseadvecdiff.kernels.artvisc'/>
<%include file='pyfr.solvers.euler.kernels.rsolvers.${rsolver}'/>
<%include file='pyfr.solvers.navstokes.kernels.flux'/>
@@ -12,7 +13,8 @@
${pyfr.expand('bc_ldg_grad_state', 'ul', 'nl', 'gradul', 'gradur')};
fpdtype_t fvr[${ndims}][${nvars}] = {{0}};
- ${pyfr.expand('viscous_flux_add', 'ur', 'gradur', 'amul', 'fvr')};
+ ${pyfr.expand('viscous_flux_add', 'ur', 'gradur', 'fvr')};
+ ${pyfr.expand('artificial_viscosity_add', 'gradur', 'fvr', 'amul')};
// Inviscid (Riemann solve) state
${pyfr.expand('bc_rsolve_state', 'ul', 'nl', 'ur', 'ploc', 't')};
diff --git a/pyfr/solvers/navstokes/kernels/entropy.mako b/pyfr/solvers/navstokes/kernels/entropy.mako
deleted file mode 100644
index 7c876f9..0000000
--- a/pyfr/solvers/navstokes/kernels/entropy.mako
+++ /dev/null
@@ -1,15 +0,0 @@
-# -*- coding: utf-8 -*-
-<%inherit file='base'/>
-<%namespace module='pyfr.backends.base.makoutil' name='pyfr'/>
-
-<%pyfr:kernel name='entropy' ndim='2'
- u='in fpdtype_t[${str(nvars)}]'
- s='out fpdtype_t'>
- fpdtype_t invrho = 1.0/u[0], E = u[${nvars - 1}];
-
- // Compute the pressure
- fpdtype_t p = ${c['gamma'] - 1}*(E - 0.5*invrho*${pyfr.dot('u[{i}]', i=(1, ndims + 1))});
-
- // Compute Entropy
- s = p*pow(invrho, ${c['gamma']});
-</%pyfr:kernel>
diff --git a/pyfr/solvers/navstokes/kernels/flux.mako b/pyfr/solvers/navstokes/kernels/flux.mako
index f4a0792..fd7122f 100644
--- a/pyfr/solvers/navstokes/kernels/flux.mako
+++ b/pyfr/solvers/navstokes/kernels/flux.mako
@@ -2,7 +2,7 @@
<%namespace module='pyfr.backends.base.makoutil' name='pyfr'/>
% if ndims == 2:
-<%pyfr:macro name='viscous_flux_add' params='uin, grad_uin, amuin, fout'>
+<%pyfr:macro name='viscous_flux_add' params='uin, grad_uin, fout'>
fpdtype_t rho = uin[0], rhou = uin[1], rhov = uin[2], E = uin[3];
fpdtype_t rcprho = 1.0/rho;
@@ -30,10 +30,6 @@
fpdtype_t mu_c = ${c['mu']};
% endif
-% if art_vis == 'mu':
- mu_c += amuin;
-% endif
-
// Compute temperature derivatives (c_v*dT/d[x,y])
fpdtype_t T_x = rcprho*(E_x - (rcprho*rho_x*E + u*u_x + v*v_x));
fpdtype_t T_y = rcprho*(E_y - (rcprho*rho_y*E + u*u_y + v*v_y));
@@ -50,7 +46,7 @@
fout[1][3] += u*t_xy + v*t_yy + -mu_c*${c['gamma']/c['Pr']}*T_y;
</%pyfr:macro>
% elif ndims == 3:
-<%pyfr:macro name='viscous_flux_add' params='uin, grad_uin, amuin, fout'>
+<%pyfr:macro name='viscous_flux_add' params='uin, grad_uin, fout'>
fpdtype_t rho = uin[0];
fpdtype_t rhou = uin[1], rhov = uin[2], rhow = uin[3];
fpdtype_t E = uin[4];
@@ -87,10 +83,6 @@
fpdtype_t mu_c = ${c['mu']};
% endif
-% if art_vis == 'mu':
- mu_c += amuin;
-% endif
-
// Compute temperature derivatives (c_v*dT/d[x,y,z])
fpdtype_t T_x = rcprho*(E_x - (rcprho*rho_x*E + u*u_x + v*v_x + w*w_x));
fpdtype_t T_y = rcprho*(E_y - (rcprho*rho_y*E + u*u_y + v*v_y + w*w_y));
diff --git a/pyfr/solvers/navstokes/kernels/intcflux.mako b/pyfr/solvers/navstokes/kernels/intcflux.mako
index 1c07f20..cc377a6 100644
--- a/pyfr/solvers/navstokes/kernels/intcflux.mako
+++ b/pyfr/solvers/navstokes/kernels/intcflux.mako
@@ -2,6 +2,7 @@
<%inherit file='base'/>
<%namespace module='pyfr.backends.base.makoutil' name='pyfr'/>
+<%include file='pyfr.solvers.baseadvecdiff.kernels.artvisc'/>
<%include file='pyfr.solvers.euler.kernels.rsolvers.${rsolver}'/>
<%include file='pyfr.solvers.navstokes.kernels.flux'/>
@@ -22,12 +23,14 @@
% if beta != -0.5:
fpdtype_t fvl[${ndims}][${nvars}] = {{0}};
- ${pyfr.expand('viscous_flux_add', 'ul', 'gradul', 'amul', 'fvl')};
+ ${pyfr.expand('viscous_flux_add', 'ul', 'gradul', 'fvl')};
+ ${pyfr.expand('artificial_viscosity_add', 'gradul', 'fvl', 'amul')};
% endif
% if beta != 0.5:
fpdtype_t fvr[${ndims}][${nvars}] = {{0}};
- ${pyfr.expand('viscous_flux_add', 'ur', 'gradur', 'amur', 'fvr')};
+ ${pyfr.expand('viscous_flux_add', 'ur', 'gradur', 'fvr')};
+ ${pyfr.expand('artificial_viscosity_add', 'gradur', 'fvr', 'amur')};
% endif
% for i in range(nvars):
diff --git a/pyfr/solvers/navstokes/kernels/mpicflux.mako b/pyfr/solvers/navstokes/kernels/mpicflux.mako
index 34ddb12..10a6be0 100644
--- a/pyfr/solvers/navstokes/kernels/mpicflux.mako
+++ b/pyfr/solvers/navstokes/kernels/mpicflux.mako
@@ -2,6 +2,7 @@
<%inherit file='base'/>
<%namespace module='pyfr.backends.base.makoutil' name='pyfr'/>
+<%include file='pyfr.solvers.baseadvecdiff.kernels.artvisc'/>
<%include file='pyfr.solvers.euler.kernels.rsolvers.${rsolver}'/>
<%include file='pyfr.solvers.navstokes.kernels.flux'/>
@@ -22,12 +23,14 @@
% if beta != -0.5:
fpdtype_t fvl[${ndims}][${nvars}] = {{0}};
- ${pyfr.expand('viscous_flux_add', 'ul', 'gradul', 'amul', 'fvl')};
+ ${pyfr.expand('viscous_flux_add', 'ul', 'gradul', 'fvl')};
+ ${pyfr.expand('artificial_viscosity_add', 'gradul', 'fvl', 'amul')};
% endif
% if beta != 0.5:
fpdtype_t fvr[${ndims}][${nvars}] = {{0}};
- ${pyfr.expand('viscous_flux_add', 'ur', 'gradur', 'amur', 'fvr')};
+ ${pyfr.expand('viscous_flux_add', 'ur', 'gradur', 'fvr')};
+ ${pyfr.expand('artificial_viscosity_add', 'gradur', 'fvr', 'amur')};
% endif
% for i in range(nvars):
diff --git a/pyfr/solvers/navstokes/kernels/tflux.mako b/pyfr/solvers/navstokes/kernels/tflux.mako
index 1989b85..813e63d 100644
--- a/pyfr/solvers/navstokes/kernels/tflux.mako
+++ b/pyfr/solvers/navstokes/kernels/tflux.mako
@@ -1,6 +1,8 @@
# -*- coding: utf-8 -*-
<%inherit file='base'/>
<%namespace module='pyfr.backends.base.makoutil' name='pyfr'/>
+
+<%include file='pyfr.solvers.baseadvecdiff.kernels.artvisc'/>
<%include file='pyfr.solvers.euler.kernels.flux'/>
<%include file='pyfr.solvers.navstokes.kernels.flux'/>
@@ -13,7 +15,8 @@
fpdtype_t ftemp[${ndims}][${nvars}];
fpdtype_t p, v[${ndims}];
${pyfr.expand('inviscid_flux', 'u', 'ftemp', 'p', 'v')};
- ${pyfr.expand('viscous_flux_add', 'u', 'f', 'amu', 'ftemp')};
+ ${pyfr.expand('viscous_flux_add', 'u', 'f', 'ftemp')};
+ ${pyfr.expand('artificial_viscosity_add', 'f', 'ftemp', 'amu')};
// Transform the fluxes
% for i, j in pyfr.ndrange(ndims, nvars):
diff --git a/pyfr/solvers/navstokes/system.py b/pyfr/solvers/navstokes/system.py
index 7d82b89..2428c6e 100644
--- a/pyfr/solvers/navstokes/system.py
+++ b/pyfr/solvers/navstokes/system.py
@@ -14,67 +14,3 @@ class NavierStokesSystem(BaseAdvectionDiffusionSystem):
intinterscls = NavierStokesIntInters
mpiinterscls = NavierStokesMPIInters
bbcinterscls = NavierStokesBaseBCInters
-
- def rhs(self, t, uinbank, foutbank):
- runall = self.backend.runall
- q1, q2 = self._queues
- kernels = self._kernels
-
- self.eles_scal_upts_inb.active = uinbank
- self.eles_scal_upts_outb.active = foutbank
-
- q1 << kernels['eles', 'disu']()
- q1 << kernels['mpiint', 'scal_fpts_pack']()
- runall([q1])
-
- if ('eles', 'copy_soln') in kernels:
- q1 << kernels['eles', 'copy_soln']()
- q1 << kernels['iint', 'con_u']()
- q1 << kernels['bcint', 'con_u'](t=t)
- q1 << kernels['eles', 'tgradpcoru_upts']()
- if ('eles', 'art_visc') in kernels:
- q1 << kernels['eles', 'entropy']()
- q1 << kernels['eles', 'modal_entropy']()
- q1 << kernels['eles', 'art_visc']()
- q1 << kernels['mpiint', 'avis_fpts_pack']()
-
- q2 << kernels['mpiint', 'scal_fpts_send']()
- q2 << kernels['mpiint', 'scal_fpts_recv']()
- q2 << kernels['mpiint', 'scal_fpts_unpack']()
-
- runall([q1, q2])
-
- q1 << kernels['mpiint', 'con_u']()
- q1 << kernels['eles', 'tgradcoru_upts']()
- q1 << kernels['eles', 'gradcoru_upts']()
- q1 << kernels['eles', 'gradcoru_fpts']()
- q1 << kernels['mpiint', 'vect_fpts_pack']()
- if ('eles', 'avis') in kernels:
- q2 << kernels['mpiint', 'avis_fpts_send']()
- q2 << kernels['mpiint', 'avis_fpts_recv']()
- q2 << kernels['mpiint', 'avis_fpts_unpack']()
-
- runall([q1, q2])
-
- if ('eles', 'gradcoru_qpts') in kernels:
- q1 << kernels['eles', 'gradcoru_qpts']()
- q1 << kernels['eles', 'tdisf']()
- q1 << kernels['eles', 'tdivtpcorf']()
- q1 << kernels['iint', 'comm_flux']()
- q1 << kernels['bcint', 'comm_flux'](t=t)
-
- q2 << kernels['mpiint', 'vect_fpts_send']()
- q2 << kernels['mpiint', 'vect_fpts_recv']()
- q2 << kernels['mpiint', 'vect_fpts_unpack']()
-
- runall([q1, q2])
-
- q1 << kernels['mpiint', 'comm_flux']()
- q1 << kernels['eles', 'tdivtconf']()
- if ('eles', 'tdivf_qpts') in kernels:
- q1 << kernels['eles', 'tdivf_qpts']()
- q1 << kernels['eles', 'negdivconf'](t=t)
- q1 << kernels['eles', 'divf_upts']()
- else:
- q1 << kernels['eles', 'negdivconf'](t=t)
- runall([q1])
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/pyfr.git
More information about the debian-science-commits
mailing list