[pyfr] 02/88: Refactor integrators to support dual time stepping.
Ghislain Vaillant
ghisvail-guest at moszumanska.debian.org
Wed Nov 16 12:05:23 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 3170f1744ef90bcf441a1b625570e06654720a9c
Author: Freddie Witherden <freddie at witherden.org>
Date: Tue Feb 9 14:17:49 2016 +0000
Refactor integrators to support dual time stepping.
---
pyfr/integrators/__init__.py | 37 ++++++++----
pyfr/integrators/base.py | 93 +++++++++++++++++++++++++------
pyfr/integrators/dual/__init__.py | 5 ++
pyfr/integrators/dual/base.py | 12 ++++
pyfr/integrators/dual/controllers.py | 7 +++
pyfr/integrators/dual/pseudosteppers.py | 7 +++
pyfr/integrators/dual/steppers.py | 7 +++
pyfr/integrators/std/__init__.py | 4 ++
pyfr/integrators/std/base.py | 24 ++++++++
pyfr/integrators/{ => std}/controllers.py | 53 ++----------------
pyfr/integrators/{ => std}/steppers.py | 31 ++++-------
pyfr/solvers/euler/elements.py | 2 +
12 files changed, 186 insertions(+), 96 deletions(-)
diff --git a/pyfr/integrators/__init__.py b/pyfr/integrators/__init__.py
index db145f0..689ae2c 100644
--- a/pyfr/integrators/__init__.py
+++ b/pyfr/integrators/__init__.py
@@ -2,25 +2,42 @@
import re
-from pyfr.integrators.controllers import BaseController
-from pyfr.integrators.steppers import BaseStepper
+from pyfr.integrators.dual import (BaseDualController, BaseDualPseudoStepper,
+ BaseDualStepper)
+from pyfr.integrators.std import BaseStdController, BaseStdStepper
from pyfr.util import subclass_where
def get_integrator(backend, systemcls, rallocs, mesh, initsoln, cfg):
- # Look-up the controller, stepper and writer names
- c = cfg.get('solver-time-integrator', 'controller')
- s = cfg.get('solver-time-integrator', 'scheme')
+ form = cfg.get('solver-time-integrator', 'formulation', 'std')
- controller = subclass_where(BaseController, controller_name=c)
- stepper = subclass_where(BaseStepper, stepper_name=s)
+ if form == 'std':
+ cn = cfg.get('solver-time-integrator', 'controller')
+ sn = cfg.get('solver-time-integrator', 'scheme')
+
+ cc = subclass_where(BaseStdController, controller_name=cn)
+ sc = subclass_where(BaseStdStepper, stepper_name=sn)
+
+ bases = [(cn, cc), (sn, sc)]
+ elif form == 'dual':
+ cn = cfg.get('solver-time-integrator', 'controller')
+ pn = cfg.get('solver-time-integrator', 'pseudo-scheme')
+ sn = cfg.get('solver-time-integrator', 'scheme')
+
+ cc = subclass_where(BaseDualController, controller_name=cn)
+ pc = subclass_where(BaseDualPseudoStepper, pseudo_stepper_name=pn)
+ sn = subclass_where(BaseDualStepper, stepper_name=sn)
+
+ bases = [(cn, cc), (pn, pc), (sn, sc)]
+ else:
+ raise ValueError('Invalid integrator formulation')
# Determine the integrator name
- name = '_'.join([c, s])
+ name = '_'.join([form] + list(bn for bn, bc in bases) + ['integrator'])
name = re.sub('(?:^|_|-)([a-z])', lambda m: m.group(1).upper(), name)
- # Composite the classes together to form a new type
- integrator = type(name, (controller, stepper), dict(name=name))
+ # Composite the base classes together to form a new type
+ integrator = type(name, tuple(bc for bn, bc in bases), dict(name=name))
# Construct and return an instance of this new integrator class
return integrator(backend, systemcls, rallocs, mesh, initsoln, cfg)
diff --git a/pyfr/integrators/base.py b/pyfr/integrators/base.py
index 4074bfc..5bc3ccc 100644
--- a/pyfr/integrators/base.py
+++ b/pyfr/integrators/base.py
@@ -2,11 +2,13 @@
from abc import ABCMeta, abstractmethod, abstractproperty
from collections import deque
+import re
import numpy as np
from pyfr.inifile import Inifile
from pyfr.mpiutil import get_comm_rank_root, get_mpi
+from pyfr.plugins import get_plugin
from pyfr.util import proxylist
@@ -17,9 +19,12 @@ class BaseIntegrator(object, metaclass=ABCMeta):
self.cfg = cfg
self.isrestart = initsoln is not None
- # Sanity checks
- if self._controller_needs_errest and not self._stepper_has_errest:
- raise TypeError('Incompatible stepper/controller combination')
+ # Ensure the system is compatible with our formulation
+ if self.formulation not in systemcls.elementscls.formulations:
+ raise RuntimeError(
+ 'System {0} does not support time stepping formulation {1}'
+ .format(systemcls.name, self.formulation)
+ )
# Start time
self.tstart = cfg.getfloat('solver-time-integrator', 'tstart', 0.0)
@@ -32,36 +37,90 @@ class BaseIntegrator(object, metaclass=ABCMeta):
else:
self.tcurr = self.tstart
+ # List of target times to advance to
self.tlist = deque([self.tend])
- # Determine the amount of temp storage required by thus method
+ # Accepted and rejected step counters
+ self.nacptsteps = 0
+ self.nrjctsteps = 0
+ self.nacptchain = 0
+
+ # Current and minimum time steps
+ self._dt = self.cfg.getfloat('solver-time-integrator', 'dt')
+ self.dtmin = 1.0e-12
+
+ # Determine the amount of temp storage required by this method
nreg = self._stepper_nregs
# Construct the relevant mesh partition
self.system = systemcls(backend, rallocs, mesh, initsoln, nreg, cfg)
+ # Storage register banks
+ self._regs, self._regidx = self._get_reg_banks(nreg)
+
# Extract the UUID of the mesh (to be saved with solutions)
self.mesh_uuid = mesh['mesh_uuid']
# Get a queue for subclasses to use
self._queue = backend.queue()
- # Get the number of degrees of freedom in this partition
- ndofs = sum(self.system.ele_ndofs)
+ # Global degree of freedom count
+ self._gnofs = self._get_gndofs()
+
+ # Bank index of solution
+ self._idxcurr = 0
+
+ # Solution cache
+ self._curr_soln = None
+ # Event handlers for advance_to
+ self.completed_step_handlers = proxylist(self._get_plugins())
+
+ # Delete the memory-intensive elements map from the system
+ del self.system.ele_map
+
+ def _get_reg_banks(self, nreg):
+ regs, regidx = [], list(range(nreg))
+
+ # Create a proxylist of matrix-banks for each storage register
+ for i in regidx:
+ regs.append(
+ proxylist([self.backend.matrix_bank(em, i)
+ for em in self.system.ele_banks])
+ )
+
+ return regs, regidx
+
+ def _get_gndofs(self):
comm, rank, root = get_comm_rank_root()
+ # Get the number of degrees of freedom in this partition
+ ndofs = sum(self.system.ele_ndofs)
+
# Sum to get the global number over all partitions
- self._gndofs = comm.allreduce(ndofs, op=get_mpi('sum'))
+ return comm.allreduce(ndofs, op=get_mpi('sum'))
+
+ def _get_plugins(self):
+ plugins = []
+
+ for s in self.cfg.sections():
+ m = re.match('soln-plugin-(.+?)(?:-(.+))?$', s)
+ if m:
+ cfgsect, name, suffix = m.group(0), m.group(1), m.group(2)
+
+ # Instantiate
+ plugins.append(get_plugin(name, self, cfgsect, suffix))
- def _kernel(self, name, nargs):
+ return plugins
+
+ def _get_kernels(self, name, nargs, **kwargs):
# Transpose from [nregs][neletypes] to [neletypes][nregs]
transregs = zip(*self._regs)
# Generate an kernel for each element type
kerns = proxylist([])
for tr in transregs:
- kerns.append(self.backend.kernel(name, *tr[:nargs]))
+ kerns.append(self.backend.kernel(name, *tr[:nargs], **kwargs))
return kerns
@@ -84,6 +143,14 @@ class BaseIntegrator(object, metaclass=ABCMeta):
tlist.extend(ta)
tlist.extend(tb)
+ @property
+ def soln(self):
+ # If we do not have the solution cached then fetch it
+ if not self._curr_soln:
+ self._curr_soln = self.system.ele_scal_upts(self._idxcurr)
+
+ return self._curr_soln
+
@abstractmethod
def step(self, t, dt):
pass
@@ -93,14 +160,6 @@ class BaseIntegrator(object, metaclass=ABCMeta):
pass
@abstractproperty
- def _controller_needs_errest(self):
- pass
-
- @abstractproperty
- def _stepper_has_errest(self):
- pass
-
- @abstractproperty
def _stepper_nfevals(self):
pass
diff --git a/pyfr/integrators/dual/__init__.py b/pyfr/integrators/dual/__init__.py
new file mode 100644
index 0000000..c6de43c
--- /dev/null
+++ b/pyfr/integrators/dual/__init__.py
@@ -0,0 +1,5 @@
+# -*- coding: utf-8 -*-
+
+from pyfr.integrators.dual.controllers import BaseDualController
+from pyfr.integrators.dual.pseudosteppers import BaseDualPseudoStepper
+from pyfr.integrators.dual.steppers import BaseDualStepper
diff --git a/pyfr/integrators/dual/base.py b/pyfr/integrators/dual/base.py
new file mode 100644
index 0000000..6c1036f
--- /dev/null
+++ b/pyfr/integrators/dual/base.py
@@ -0,0 +1,12 @@
+# -*- coding: utf-8 -*-
+
+from pyfr.integrators.base import BaseIntegrator
+
+
+class BaseDualIntegrator(BaseIntegrator):
+ formulation = 'dual'
+
+ def __init__(self, *args, **kwargs):
+ super().__init__(*args, **kwargs)
+
+ pass
diff --git a/pyfr/integrators/dual/controllers.py b/pyfr/integrators/dual/controllers.py
new file mode 100644
index 0000000..19bd585
--- /dev/null
+++ b/pyfr/integrators/dual/controllers.py
@@ -0,0 +1,7 @@
+# -*- coding: utf-8 -*-
+
+from pyfr.integrators.dual.base import BaseDualIntegrator
+
+
+class BaseDualController(BaseDualIntegrator):
+ pass
diff --git a/pyfr/integrators/dual/pseudosteppers.py b/pyfr/integrators/dual/pseudosteppers.py
new file mode 100644
index 0000000..92667fc
--- /dev/null
+++ b/pyfr/integrators/dual/pseudosteppers.py
@@ -0,0 +1,7 @@
+# -*- coding: utf-8 -*-
+
+from pyfr.integrators.dual.base import BaseDualIntegrator
+
+
+class BaseDualPseudoStepper(BaseDualIntegrator):
+ pass
diff --git a/pyfr/integrators/dual/steppers.py b/pyfr/integrators/dual/steppers.py
new file mode 100644
index 0000000..03a0761
--- /dev/null
+++ b/pyfr/integrators/dual/steppers.py
@@ -0,0 +1,7 @@
+# -*- coding: utf-8 -*-
+
+from pyfr.integrators.dual.base import BaseDualIntegrator
+
+
+class BaseDualStepper(BaseDualIntegrator):
+ pass
diff --git a/pyfr/integrators/std/__init__.py b/pyfr/integrators/std/__init__.py
new file mode 100644
index 0000000..07f064f
--- /dev/null
+++ b/pyfr/integrators/std/__init__.py
@@ -0,0 +1,4 @@
+# -*- coding: utf-8 -*-
+
+from pyfr.integrators.std.controllers import BaseStdController
+from pyfr.integrators.std.steppers import BaseStdStepper
diff --git a/pyfr/integrators/std/base.py b/pyfr/integrators/std/base.py
new file mode 100644
index 0000000..c799b7b
--- /dev/null
+++ b/pyfr/integrators/std/base.py
@@ -0,0 +1,24 @@
+# -*- coding: utf-8 -*-
+
+from abc import abstractproperty
+
+from pyfr.integrators.base import BaseIntegrator
+
+
+class BaseStdIntegrator(BaseIntegrator):
+ formulation = 'std'
+
+ def __init__(self, *args, **kwargs):
+ super().__init__(*args, **kwargs)
+
+ # Sanity checks
+ if self._controller_needs_errest and not self._stepper_has_errest:
+ raise TypeError('Incompatible stepper/controller combination')
+
+ @abstractproperty
+ def _controller_needs_errest(self):
+ pass
+
+ @abstractproperty
+ def _stepper_has_errest(self):
+ pass
diff --git a/pyfr/integrators/controllers.py b/pyfr/integrators/std/controllers.py
similarity index 74%
rename from pyfr/integrators/controllers.py
rename to pyfr/integrators/std/controllers.py
index 286e3f1..cefc7da 100644
--- a/pyfr/integrators/controllers.py
+++ b/pyfr/integrators/std/controllers.py
@@ -1,57 +1,22 @@
# -*- coding: utf-8 -*-
import math
-import re
-from pyfr.integrators.base import BaseIntegrator
+from pyfr.integrators.std.base import BaseStdIntegrator
from pyfr.mpiutil import get_comm_rank_root, get_mpi
-from pyfr.plugins import get_plugin
from pyfr.util import memoize, proxylist
-class BaseController(BaseIntegrator):
+class BaseStdController(BaseStdIntegrator):
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
- # Current and minimum time steps
- self._dt = self.cfg.getfloat('solver-time-integrator', 'dt')
- self.dtmin = 1.0e-12
-
# Solution filtering frequency
self._fnsteps = self.cfg.getint('soln-filter', 'nsteps', '0')
- # Bank index of solution
- self._idxcurr = 0
-
- # Solution cache
- self._curr_soln = None
-
- # Accepted and rejected step counters
- self.nacptsteps = 0
- self.nrjctsteps = 0
- self.nacptchain = 0
-
# Stats on the most recent step
self.stepinfo = []
- # Event handlers for advance_to
- self.completed_step_handlers = proxylist([])
-
- # Load any plugins specified in the config file
- for s in self.cfg.sections():
- m = re.match('soln-plugin-(.+?)(?:-(.+))?$', s)
- if m:
- cfgsect, name, suffix = m.group(0), m.group(1), m.group(2)
-
- # Instantiate
- plugin = get_plugin(name, self, cfgsect, suffix)
-
- # Register as an event handler
- self.completed_step_handlers.append(plugin)
-
- # Delete the memory-intensive elements map from the system
- del self.system.ele_map
-
def _accept_step(self, dt, idxcurr, err=None):
self.tcurr += dt
self.nacptsteps += 1
@@ -87,16 +52,8 @@ class BaseController(BaseIntegrator):
def nsteps(self):
return self.nacptsteps + self.nrjctsteps
- @property
- def soln(self):
- # If we do not have the solution cached then fetch it
- if not self._curr_soln:
- self._curr_soln = self.system.ele_scal_upts(self._idxcurr)
-
- return self._curr_soln
-
-class NoneController(BaseController):
+class StdNoneController(BaseStdController):
controller_name = 'none'
@property
@@ -118,7 +75,7 @@ class NoneController(BaseController):
self._accept_step(dt, idxcurr)
-class PIController(BaseController):
+class StdPIController(BaseStdController):
controller_name = 'pi'
def __init__(self, *args, **kwargs):
@@ -148,7 +105,7 @@ class PIController(BaseController):
@memoize
def _get_errest_kerns(self):
- return self._kernel('errest', nargs=3)
+ return self._get_kernels('errest', nargs=3)
def _errest(self, x, y, z):
comm, rank, root = get_comm_rank_root()
diff --git a/pyfr/integrators/steppers.py b/pyfr/integrators/std/steppers.py
similarity index 90%
rename from pyfr/integrators/steppers.py
rename to pyfr/integrators/std/steppers.py
index 5a7e86e..c28c6f2 100644
--- a/pyfr/integrators/steppers.py
+++ b/pyfr/integrators/std/steppers.py
@@ -1,24 +1,13 @@
# -*- coding: utf-8 -*-
-from pyfr.integrators.base import BaseIntegrator
-from pyfr.util import memoize, proxylist
+from pyfr.integrators.std.base import BaseStdIntegrator
+from pyfr.util import memoize
-class BaseStepper(BaseIntegrator):
+class BaseStdStepper(BaseStdIntegrator):
def __init__(self, *args, **kwargs):
super().__init__(*args, **kwargs)
- backend = self.backend
- elemats = self.system.ele_banks
-
- # Create a proxylist of matrix-banks for each storage register
- self._regs = regs = []
- self._regidx = regidx = []
- for i in range(self._stepper_nregs):
- b = proxylist([backend.matrix_bank(em, i) for em in elemats])
- regs.append(b)
- regidx.append(i)
-
# Add kernel cache
self._axnpby_kerns = {}
@@ -30,7 +19,7 @@ class BaseStepper(BaseIntegrator):
@memoize
def _get_axnpby_kerns(self, n):
- return self._kernel('axnpby', nargs=n)
+ return self._get_kernels('axnpby', nargs=n)
def _add(self, *args):
# Get a suitable set of axnpby kernels
@@ -43,7 +32,7 @@ class BaseStepper(BaseIntegrator):
self._queue % axnpby(*args[::2])
-class EulerStepper(BaseStepper):
+class StdEulerStepper(BaseStdStepper):
stepper_name = 'euler'
@property
@@ -72,7 +61,7 @@ class EulerStepper(BaseStepper):
return ut
-class TVDRK3Stepper(BaseStepper):
+class StdTVDRK3Stepper(BaseStdStepper):
stepper_name = 'tvd-rk3'
@property
@@ -117,7 +106,7 @@ class TVDRK3Stepper(BaseStepper):
return r1
-class RK4Stepper(BaseStepper):
+class StdRK4StdStepper(BaseStdStepper):
stepper_name = 'rk4'
@property
@@ -180,7 +169,7 @@ class RK4Stepper(BaseStepper):
return r1
-class RKVdH2RStepper(BaseStepper):
+class StdRKVdH2RStepper(BaseStdStepper):
# Coefficients
a = []
b = []
@@ -244,7 +233,7 @@ class RKVdH2RStepper(BaseStepper):
return (r2, rold, rerr) if errest else r2
-class RK34Stepper(RKVdH2RStepper):
+class StdRK34Stepper(StdRKVdH2RStepper):
stepper_name = 'rk34'
a = [
@@ -272,7 +261,7 @@ class RK34Stepper(RKVdH2RStepper):
return 3
-class RK45Stepper(RKVdH2RStepper):
+class StdRK45Stepper(StdRKVdH2RStepper):
stepper_name = 'rk45'
a = [
diff --git a/pyfr/solvers/euler/elements.py b/pyfr/solvers/euler/elements.py
index 4c89e16..fd94b1b 100644
--- a/pyfr/solvers/euler/elements.py
+++ b/pyfr/solvers/euler/elements.py
@@ -4,6 +4,8 @@ from pyfr.solvers.baseadvec import BaseAdvectionElements
class BaseFluidElements(object):
+ formulations = ['std', 'dual']
+
privarmap = {2: ['rho', 'u', 'v', 'p'],
3: ['rho', 'u', 'v', 'w', 'p']}
--
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