[pyfr] 45/88: Unify the various shock capturing kernels.

Ghislain Vaillant ghisvail-guest at moszumanska.debian.org
Wed Nov 16 12:05:28 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 25e902df3a9257ce58567efab2d7461b5b47d67b
Author: Freddie Witherden <freddie at witherden.org>
Date:   Wed Jun 1 20:57:16 2016 -0700

    Unify the various shock capturing kernels.
---
 pyfr/backends/mic/generator.py                     |  1 +
 pyfr/backends/openmp/generator.py                  |  1 +
 pyfr/polys.py                                      |  5 ++++
 pyfr/solvers/baseadvecdiff/elements.py             | 34 ++++------------------
 .../solvers/baseadvecdiff/kernels/shocksensor.mako | 22 ++++++++------
 pyfr/solvers/baseadvecdiff/kernels/shockvar.mako   |  9 ------
 pyfr/solvers/baseadvecdiff/system.py               |  4 +--
 7 files changed, 26 insertions(+), 50 deletions(-)

diff --git a/pyfr/backends/mic/generator.py b/pyfr/backends/mic/generator.py
index ade6573..8befc78 100644
--- a/pyfr/backends/mic/generator.py
+++ b/pyfr/backends/mic/generator.py
@@ -17,6 +17,7 @@ class MICKernelGenerator(BaseKernelGenerator):
                           int align = PYFR_ALIGN_BYTES / sizeof(fpdtype_t);
                           int cb, ce;
                           loop_sched_1d(_nx, align, &cb, &ce);
+                          #pragma omp simd
                           for (int _x = cb; _x < ce; _x++)
                           {{
                               {body}
diff --git a/pyfr/backends/openmp/generator.py b/pyfr/backends/openmp/generator.py
index 4ed7b05..c9677b4 100644
--- a/pyfr/backends/openmp/generator.py
+++ b/pyfr/backends/openmp/generator.py
@@ -17,6 +17,7 @@ class OpenMPKernelGenerator(BaseKernelGenerator):
                           int align = PYFR_ALIGN_BYTES / sizeof(fpdtype_t);
                           int cb, ce;
                           loop_sched_1d(_nx, align, &cb, &ce);
+                          #pragma omp simd
                           for (int _x = cb; _x < ce; _x++)
                           {{
                               {body}
diff --git a/pyfr/polys.py b/pyfr/polys.py
index 2f2b6b0..ce6f830 100644
--- a/pyfr/polys.py
+++ b/pyfr/polys.py
@@ -80,6 +80,11 @@ class BasePolyBasis(object):
     def vdm(self):
         return self.ortho_basis_at(self.pts)
 
+    @lazyprop
+    @chop
+    def invvdm(self):
+        return np.linalg.inv(self.vdm)
+
 
 class LinePolyBasis(BasePolyBasis):
     name = 'line'
diff --git a/pyfr/solvers/baseadvecdiff/elements.py b/pyfr/solvers/baseadvecdiff/elements.py
index b57b7b3..0bd1318 100644
--- a/pyfr/solvers/baseadvecdiff/elements.py
+++ b/pyfr/solvers/baseadvecdiff/elements.py
@@ -1,7 +1,5 @@
 # -*- coding: utf-8 -*-
 
-import numpy as np
-
 from pyfr.backends.base.kernels import ComputeMetaKernel
 from pyfr.solvers.baseadvec import BaseAdvectionElements
 
@@ -80,10 +78,7 @@ class BaseAdvectionDiffusionElements(BaseAdvectionElements):
         if shock_capturing == 'artificial-viscosity':
             tags = {'align'}
 
-            # Register pointwise kernels
-            backend.pointwise.register(
-                'pyfr.solvers.baseadvecdiff.kernels.shockvar'
-            )
+            # Register the kernels
             backend.pointwise.register(
                 'pyfr.solvers.baseadvecdiff.kernels.shocksensor'
             )
@@ -91,41 +86,22 @@ class BaseAdvectionDiffusionElements(BaseAdvectionElements):
             # Obtain the scalar variable to be used for shock sensing
             shockvar = self.convarmap[self.ndims].index(self.shockvar)
 
-            # Common template arguments
+            # Template arguments
             tplargs = dict(
-                nvars=self.nvars, nupts=self.nupts, shockvar=shockvar,
+                nvars=self.nvars, nupts=self.nupts, svar=shockvar,
                 c=self.cfg.items_as('solver-artificial-viscosity', float),
                 order=self.basis.order, ubdegs=self.basis.ubasis.degrees,
+                invvdm=self.basis.ubasis.invvdm.T
             )
 
             # Allocate space for the artificial viscosity
             self.artvisc = backend.matrix((1, self.neles), extent='artvisc',
                                           tags=tags)
 
-            # Scratch space
-            tmp_upts = backend.matrix((2*self.nupts, 1, self.neles),
-                                      aliases=self._vect_upts, tags=tags)
-            svar_upts = tmp_upts.rslice(0, self.nupts)
-            modal_svar_upts = tmp_upts.rslice(self.nupts, 2*self.nupts)
-
-
-            # Extract the scalar variable to be used for shock sensing
-            self.kernels['shockvar'] = lambda: backend.kernel(
-                'shockvar', tplargs=tplargs, dims=[self.nupts, self.neles],
-                u=self.scal_upts_inb, s=svar_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, svar_upts, out=modal_svar_upts
-            )
-
             # Apply the sensor to estimate the required artificial viscosity
             self.kernels['shocksensor'] = lambda: backend.kernel(
                 'shocksensor', tplargs=tplargs, dims=[self.neles],
-                s=modal_svar_upts, artvisc=self.artvisc
+                u=self.scal_upts_inb, artvisc=self.artvisc
             )
         elif shock_capturing == 'none':
             self.artvisc = None
diff --git a/pyfr/solvers/baseadvecdiff/kernels/shocksensor.mako b/pyfr/solvers/baseadvecdiff/kernels/shocksensor.mako
index c677f63..33e853b 100644
--- a/pyfr/solvers/baseadvecdiff/kernels/shocksensor.mako
+++ b/pyfr/solvers/baseadvecdiff/kernels/shocksensor.mako
@@ -2,27 +2,31 @@
 <%inherit file='base'/>
 <%namespace module='pyfr.backends.base.makoutil' name='pyfr'/>
 
+<% se0 = math.log10(c['s0']) %>
+
 <%pyfr:kernel name='shocksensor' ndim='1'
-              s='in fpdtype_t[${str(nupts)}]'
+              u='in fpdtype_t[${str(nupts)}][${str(nvars)}]'
               artvisc='out fpdtype_t'>
     // Smoothness indicator
-    fpdtype_t totEn = 0.0, pnEn = 1e-15;
+    fpdtype_t totEn = 0.0, pnEn = 1e-15, tmp;
 
 % for i, deg in enumerate(ubdegs):
-    totEn += s[${i}]*s[${i}];
+    tmp = ${' + '.join('{jx}*u[{j}][{svar}]'.format(j=j, jx=jx, svar=svar)
+                       for j, jx in enumerate(invvdm[i]) if jx != 0)};
+
+    totEn += tmp*tmp;
 % if deg >= order:
-    pnEn += s[${i}]*s[${i}];
+    pnEn += tmp*tmp;
 % endif
 % endfor
 
-    fpdtype_t se0 = ${math.log10(c['s0'])};
-    fpdtype_t se  = log10(pnEn/totEn);
+    fpdtype_t se  = ${1/math.log(10)}*log(pnEn/totEn);
 
     // Compute cell-wise artificial viscosity
-    fpdtype_t mu = (se < se0 - ${c['kappa']})
+    fpdtype_t mu = (se < ${se0 - c['kappa']})
                  ? 0.0
-                 : ${0.5*c['max-artvisc']}*(1.0 + sin(${0.5*math.pi/c['kappa']}*(se - se0)));
-    mu = (se < se0 + ${c['kappa']}) ? mu : ${c['max-artvisc']};
+                 : ${0.5*c['max-artvisc']}*(1.0 + sin(${0.5*math.pi/c['kappa']}*(se - ${se0})));
+    mu = (se < ${se0 + c['kappa']}) ? mu : ${c['max-artvisc']};
 
     artvisc = mu;
 </%pyfr:kernel>
diff --git a/pyfr/solvers/baseadvecdiff/kernels/shockvar.mako b/pyfr/solvers/baseadvecdiff/kernels/shockvar.mako
deleted file mode 100644
index affb8ed..0000000
--- a/pyfr/solvers/baseadvecdiff/kernels/shockvar.mako
+++ /dev/null
@@ -1,9 +0,0 @@
-# -*- 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 458b5a5..1133615 100644
--- a/pyfr/solvers/baseadvecdiff/system.py
+++ b/pyfr/solvers/baseadvecdiff/system.py
@@ -20,9 +20,7 @@ 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']()
+        if ('eles', 'shocksensor') in kernels:
             q1 << kernels['eles', 'shocksensor']()
             q1 << kernels['mpiint', 'artvisc_fpts_pack']()
         q1 << kernels['eles', 'tgradpcoru_upts']()

-- 
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