diff --git a/boundaries/boundaryconditions.py b/boundaries/boundaryconditions.py
index f551f3aeede52ee6929d7af118902eb3055fe30e..16f2dd8519956255386c42d63aa996cc9330657b 100644
--- a/boundaries/boundaryconditions.py
+++ b/boundaries/boundaryconditions.py
@@ -1,10 +1,10 @@
import sympy as sp
-from lbmpy.simplificationfactory import createSimplificationStrategy
+from pystencils import Field, Assignment
from pystencils.astnodes import SympyAssignment
from pystencils.sympyextensions import getSymmetricPart
-from pystencils import Field
-from lbmpy.boundaries.boundaryhandling import BoundaryOffsetInfo, LbmWeightInfo
from pystencils.data_types import createType
+from lbmpy.simplificationfactory import createSimplificationStrategy
+from lbmpy.boundaries.boundaryhandling import BoundaryOffsetInfo, LbmWeightInfo
class Boundary(object):
@@ -64,7 +64,7 @@ class NoSlip(Boundary):
def __call__(self, pdfField, directionSymbol, lbMethod, **kwargs):
neighbor = BoundaryOffsetInfo.offsetFromDir(directionSymbol, lbMethod.dim)
inverseDir = BoundaryOffsetInfo.invDir(directionSymbol)
- return [sp.Eq(pdfField[neighbor](inverseDir), pdfField(directionSymbol))]
+ return [Assignment(pdfField[neighbor](inverseDir), pdfField(directionSymbol))]
def __hash__(self):
# All boundaries of these class behave equal -> should also be equal (as long as name is equal)
@@ -126,7 +126,7 @@ class UBB(Boundary):
if self._adaptVelocityToForce:
cqc = lbMethod.conservedQuantityComputation
shiftedVelEqs = cqc.equilibriumInputEquationsFromInitValues(velocity=velocity)
- velocity = [eq.rhs for eq in shiftedVelEqs.extract(cqc.firstOrderMomentSymbols).mainEquations]
+ velocity = [eq.rhs for eq in shiftedVelEqs.extract(cqc.firstOrderMomentSymbols).mainAssignments]
c_s_sq = sp.Rational(1, 3)
weightOfDirection = LbmWeightInfo.weightOfDirection
@@ -142,12 +142,12 @@ class UBB(Boundary):
densityEquations = cqc.outputEquationsFromPdfs(pdfFieldAccesses, {'density': densitySymbol})
densitySymbol = lbMethod.conservedQuantityComputation.definedSymbols()['density']
result = densityEquations.allEquations
- result += [sp.Eq(pdfField[neighbor](inverseDir),
- pdfField(direction) - velTerm * densitySymbol)]
+ result += [Assignment(pdfField[neighbor](inverseDir),
+ pdfField(direction) - velTerm * densitySymbol)]
return result
else:
- return [sp.Eq(pdfField[neighbor](inverseDir),
- pdfField(direction) - velTerm)]
+ return [Assignment(pdfField[neighbor](inverseDir),
+ pdfField(direction) - velTerm)]
class FixedDensity(Boundary):
@@ -159,16 +159,16 @@ class FixedDensity(Boundary):
def __call__(self, pdfField, directionSymbol, lbMethod, **kwargs):
"""Boundary condition that fixes the density/pressure at the obstacle"""
- def removeAsymmetricPartOfMainEquations(eqColl, dofs):
- newMainEquations = [sp.Eq(e.lhs, getSymmetricPart(e.rhs, dofs)) for e in eqColl.mainEquations]
- return eqColl.copy(newMainEquations)
+ def removeAsymmetricPartOfmainAssignments(eqColl, dofs):
+ newmainAssignments = [Assignment(e.lhs, getSymmetricPart(e.rhs, dofs)) for e in eqColl.mainAssignments]
+ return eqColl.copy(newmainAssignments)
neighbor = BoundaryOffsetInfo.offsetFromDir(directionSymbol, lbMethod.dim)
inverseDir = BoundaryOffsetInfo.invDir(directionSymbol)
cqc = lbMethod.conservedQuantityComputation
velocity = cqc.definedSymbols()['velocity']
- symmetricEq = removeAsymmetricPartOfMainEquations(lbMethod.getEquilibrium(), dofs=velocity)
+ symmetricEq = removeAsymmetricPartOfmainAssignments(lbMethod.getEquilibrium(), dofs=velocity)
substitutions = {sym: pdfField(i) for i, sym in enumerate(lbMethod.preCollisionPdfSymbols)}
symmetricEq = symmetricEq.copyWithSubstitutionsApplied(substitutions)
@@ -178,15 +178,15 @@ class FixedDensity(Boundary):
densitySymbol = cqc.definedSymbols()['density']
density = self._density
- densityEq = cqc.equilibriumInputEquationsFromInitValues(density=density).insertSubexpressions().mainEquations[0]
+ densityEq = cqc.equilibriumInputEquationsFromInitValues(density=density).insertSubexpressions().mainAssignments[0]
assert densityEq.lhs == densitySymbol
transformedDensity = densityEq.rhs
conditions = [(eq_i.rhs, sp.Equality(directionSymbol, i))
- for i, eq_i in enumerate(symmetricEq.mainEquations)] + [(0, True)]
+ for i, eq_i in enumerate(symmetricEq.mainAssignments)] + [(0, True)]
eq_component = sp.Piecewise(*conditions)
- subExprs = [sp.Eq(eq.lhs, transformedDensity if eq.lhs == densitySymbol else eq.rhs)
+ subExprs = [Assignment(eq.lhs, transformedDensity if eq.lhs == densitySymbol else eq.rhs)
for eq in symmetricEq.subexpressions]
return subExprs + [SympyAssignment(pdfField[neighbor](inverseDir),
@@ -197,8 +197,8 @@ class NeumannByCopy(Boundary):
def __call__(self, pdfField, directionSymbol, lbMethod, **kwargs):
neighbor = BoundaryOffsetInfo.offsetFromDir(directionSymbol, lbMethod.dim)
inverseDir = BoundaryOffsetInfo.invDir(directionSymbol)
- return [sp.Eq(pdfField[neighbor](inverseDir), pdfField(inverseDir)),
- sp.Eq(pdfField[neighbor](directionSymbol), pdfField(directionSymbol))]
+ return [Assignment(pdfField[neighbor](inverseDir), pdfField(inverseDir)),
+ Assignment(pdfField[neighbor](directionSymbol), pdfField(directionSymbol))]
def __hash__(self):
# All boundaries of these class behave equal -> should also be equal
@@ -212,8 +212,8 @@ class StreamInZero(Boundary):
def __call__(self, pdfField, directionSymbol, lbMethod, **kwargs):
neighbor = BoundaryOffsetInfo.offsetFromDir(directionSymbol, lbMethod.dim)
inverseDir = BoundaryOffsetInfo.invDir(directionSymbol)
- return [sp.Eq(pdfField[neighbor](inverseDir), 0),
- sp.Eq(pdfField[neighbor](directionSymbol), 0)]
+ return [Assignment(pdfField[neighbor](inverseDir), 0),
+ Assignment(pdfField[neighbor](directionSymbol), 0)]
def __hash__(self):
# All boundaries of these class behave equal -> should also be equal
diff --git a/boundaries/boundaryhandling.py b/boundaries/boundaryhandling.py
index 990b17e2c1bd30d09e7129f2ad2504b6b4821664..2e79359201b7db6da745dcadbb1ccd62ae9986b7 100644
--- a/boundaries/boundaryhandling.py
+++ b/boundaries/boundaryhandling.py
@@ -1,10 +1,10 @@
import numpy as np
import sympy as sp
-from lbmpy.stencils import inverseDirection
-from pystencils import TypedSymbol, createIndexedKernel
+from pystencils import TypedSymbol, createIndexedKernel, Assignment
from pystencils.backends.cbackend import CustomCppCode
from pystencils.boundaries import BoundaryHandling
from pystencils.boundaries.boundaryhandling import BoundaryOffsetInfo
+from lbmpy.stencils import inverseDirection
class LatticeBoltzmannBoundaryHandling(BoundaryHandling):
@@ -99,6 +99,6 @@ def createLatticeBoltzmannBoundaryKernel(pdfField, indexField, lbMethod, boundar
elements = [BoundaryOffsetInfo(lbMethod.stencil), LbmWeightInfo(lbMethod)]
indexArrDtype = indexField.dtype.numpyDtype
dirSymbol = TypedSymbol("dir", indexArrDtype.fields['dir'][0])
- elements += [sp.Eq(dirSymbol, indexField[0]('dir'))]
+ elements += [Assignment(dirSymbol, indexField[0]('dir'))]
elements += boundaryFunctor(pdfField=pdfField, directionSymbol=dirSymbol, lbMethod=lbMethod, indexField=indexField)
return createIndexedKernel(elements, [indexField], target=target, cpuOpenMP=openMP)
diff --git a/chapman_enskog/chapman_enskog.py b/chapman_enskog/chapman_enskog.py
index 286cbfff166fc675c943c341c2727cc4a18b64e1..9ab2ed46539861d5593aa8f8a6ac4f927ee6d081 100644
--- a/chapman_enskog/chapman_enskog.py
+++ b/chapman_enskog/chapman_enskog.py
@@ -77,7 +77,7 @@ class CeMoment(sp.Symbol):
class LbMethodEqMoments:
def __init__(self, lbMethod):
- self._eq = tuple(e.rhs for e in lbMethod.getEquilibrium().mainEquations)
+ self._eq = tuple(e.rhs for e in lbMethod.getEquilibrium().mainAssignments)
self._momentCache = dict()
self._postCollisionMomentCache = dict()
self._stencil = lbMethod.stencil
diff --git a/creationfunctions.py b/creationfunctions.py
index b18cb2c8c2a6b37e8d253b927a1070bcf84f6e4c..5772ce3067176a64a2868431956ca77bec3e6496 100644
--- a/creationfunctions.py
+++ b/creationfunctions.py
@@ -152,8 +152,12 @@ For example, to modify the AST one can run::
import sympy as sp
from copy import copy
-from lbmpy.turbulence_models import addSmagorinskyModel
from pystencils.cache import diskcacheNoFallback
+from pystencils.data_types import collateTypes
+from pystencils.assignment_collection.assignment_collection import AssignmentCollection
+from pystencils.field import getLayoutOfArray, Field
+from pystencils import createKernel, Assignment
+from lbmpy.turbulence_models import addSmagorinskyModel
from lbmpy.methods import createSRT, createTRT, createOrthogonalMRT, createKBCTypeTRT, \
createRawMRT, createThreeRelaxationRateMRT
from lbmpy.methods.entropic import addIterativeEntropyCondition, addEntropyCondition
@@ -164,10 +168,6 @@ import lbmpy.forcemodels as forcemodels
from lbmpy.simplificationfactory import createSimplificationStrategy
from lbmpy.updatekernels import StreamPullTwoFieldsAccessor, PeriodicTwoFieldsAccessor, CollideOnlyInplaceAccessor, \
createLBMKernel, createStreamPullWithOutputKernel
-from pystencils.data_types import collateTypes
-from pystencils.equationcollection.equationcollection import EquationCollection
-from pystencils.field import getLayoutOfArray, Field
-from pystencils import createKernel
def createLatticeBoltzmannFunction(ast=None, optimizationParams={}, **kwargs):
@@ -275,10 +275,10 @@ def createLatticeBoltzmannCollisionRule(lbMethod=None, optimizationParams={}, **
cqc = lbMethod.conservedQuantityComputation
if params['velocityInput'] is not None:
- eqs = [sp.Eq(cqc.zerothOrderMomentSymbol, sum(lbMethod.preCollisionPdfSymbols))]
+ eqs = [Assignment(cqc.zerothOrderMomentSymbol, sum(lbMethod.preCollisionPdfSymbols))]
velocityField = params['velocityInput']
- eqs += [sp.Eq(uSym, velocityField(i)) for i, uSym in enumerate(cqc.firstOrderMomentSymbols)]
- eqs = EquationCollection(eqs, [])
+ eqs += [Assignment(uSym, velocityField(i)) for i, uSym in enumerate(cqc.firstOrderMomentSymbols)]
+ eqs = AssignmentCollection(eqs, [])
collisionRule = lbMethod.getCollisionRule(conservedQuantityEquations=eqs)
else:
collisionRule = lbMethod.getCollisionRule()
diff --git a/innerloopsplit.py b/innerloopsplit.py
index f3b6fa8370381e490ab5164ad2aa1cc848d05717..10417c607140d6ea4c7cdf3f8504c0bd3c22805c 100644
--- a/innerloopsplit.py
+++ b/innerloopsplit.py
@@ -29,7 +29,7 @@ def createLbmSplitGroups(lbmCollisionEqs):
if preCollisionSymbols.intersection(lbmCollisionEqs.getDependentSymbols([e.lhs]))}
otherWrittenFields = []
- for eq in lbmCollisionEqs.mainEquations:
+ for eq in lbmCollisionEqs.mainAssignments:
if eq.lhs not in postCollisionSymbols and isinstance(eq.lhs, Field.Access):
otherWrittenFields.append(eq.lhs)
if eq.lhs not in nonCenterPostCollisionSymbols:
@@ -44,7 +44,7 @@ def createLbmSplitGroups(lbmCollisionEqs):
directionGroups = defaultdict(list)
dim = len(stencil[0])
- for direction, eq in zip(stencil, lbmCollisionEqs.mainEquations):
+ for direction, eq in zip(stencil, lbmCollisionEqs.mainAssignments):
if direction == tuple([0]*dim):
splitGroups[0].append(eq.lhs)
continue
diff --git a/methods/abstractlbmethod.py b/methods/abstractlbmethod.py
index 00cbdc5ab539d9932b2b1db3fe6a63daf0bd21b6..cf9345c26fa29a82a07067eee77dcd6d5516859e 100644
--- a/methods/abstractlbmethod.py
+++ b/methods/abstractlbmethod.py
@@ -1,13 +1,13 @@
import abc
import sympy as sp
from collections import namedtuple
-from pystencils.equationcollection import EquationCollection
+from pystencils.assignment_collection import AssignmentCollection
RelaxationInfo = namedtuple('RelaxationInfo', ['equilibriumValue', 'relaxationRate'])
-class LbmCollisionRule(EquationCollection):
+class LbmCollisionRule(AssignmentCollection):
def __init__(self, lbmMethod, *args, **kwargs):
super(LbmCollisionRule, self).__init__(*args, **kwargs)
self.method = lbmMethod
diff --git a/methods/conservedquantitycomputation.py b/methods/conservedquantitycomputation.py
index 10aa1c49cab1974ac12c873fcaae3d510ba47e8a..37fddc87af0fa608ceb55a3774c9e6d29fca94f7 100644
--- a/methods/conservedquantitycomputation.py
+++ b/methods/conservedquantitycomputation.py
@@ -1,9 +1,8 @@
import abc
-from collections import OrderedDict
-
import sympy as sp
-from pystencils.equationcollection import EquationCollection
-from pystencils.field import Field
+from collections import OrderedDict
+from pystencils.assignment_collection import AssignmentCollection
+from pystencils.field import Field, Assignment
class AbstractConservedQuantityComputation(abc.ABCMeta('ABC', (object,), {})):
@@ -23,7 +22,8 @@ class AbstractConservedQuantityComputation(abc.ABCMeta('ABC', (object,), {})):
"""
- @abc.abstractproperty
+ @property
+ @abc.abstractmethod
def conservedQuantities(self):
"""
Dict, mapping names (symbol) to dimensionality (int)
@@ -37,7 +37,8 @@ class AbstractConservedQuantityComputation(abc.ABCMeta('ABC', (object,), {})):
Returns a dict, mapping names of conserved quantities to their symbols
"""
- @abc.abstractproperty
+ @property
+ @abc.abstractmethod
def defaultValues(self):
"""
Returns a dict of symbol to default value, where "default" means that
@@ -129,84 +130,84 @@ class DensityVelocityComputation(AbstractConservedQuantityComputation):
def equilibriumInputEquationsFromPdfs(self, pdfs):
dim = len(self._stencil[0])
- eqColl = getEquationsForZerothAndFirstOrderMoment(self._stencil, pdfs, self._symbolOrder0,
+ eq_coll = getEquationsForZerothAndFirstOrderMoment(self._stencil, pdfs, self._symbolOrder0,
self._symbolsOrder1[:dim])
if self._compressible:
- eqColl = divideFirstOrderMomentsByRho(eqColl, dim)
+ eq_coll = divideFirstOrderMomentsByRho(eq_coll, dim)
- eqColl = applyForceModelShift('equilibriumVelocityShift', dim, eqColl, self._forceModel, self._compressible)
- return eqColl
+ eq_coll = applyForceModelShift('equilibriumVelocityShift', dim, eq_coll, self._forceModel, self._compressible)
+ return eq_coll
def equilibriumInputEquationsFromInitValues(self, density=1, velocity=(0, 0, 0)):
dim = len(self._stencil[0])
zerothOrderMoment = density
- firstOrderMoments = velocity[:dim]
- velOffset = [0] * dim
+ first_order_moments = velocity[:dim]
+ vel_offset = [0] * dim
if self._compressible:
if self._forceModel and hasattr(self._forceModel, 'macroscopicVelocityShift'):
- velOffset = self._forceModel.macroscopicVelocityShift(zerothOrderMoment)
+ vel_offset = self._forceModel.macroscopicVelocityShift(zerothOrderMoment)
else:
if self._forceModel and hasattr(self._forceModel, 'macroscopicVelocityShift'):
- velOffset = self._forceModel.macroscopicVelocityShift(sp.Rational(1, 1))
+ vel_offset = self._forceModel.macroscopicVelocityShift(sp.Rational(1, 1))
zerothOrderMoment -= sp.Rational(1, 1)
- eqs = [sp.Eq(self._symbolOrder0, zerothOrderMoment)]
+ eqs = [Assignment(self._symbolOrder0, zerothOrderMoment)]
- firstOrderMoments = [a - b for a, b in zip(firstOrderMoments, velOffset)]
- eqs += [sp.Eq(l, r) for l, r in zip(self._symbolsOrder1, firstOrderMoments)]
+ first_order_moments = [a - b for a, b in zip(first_order_moments, vel_offset)]
+ eqs += [Assignment(l, r) for l, r in zip(self._symbolsOrder1, first_order_moments)]
- return EquationCollection(eqs, [])
+ return AssignmentCollection(eqs, [])
def outputEquationsFromPdfs(self, pdfs, outputQuantityNamesToSymbols):
dim = len(self._stencil[0])
- eqColl = getEquationsForZerothAndFirstOrderMoment(self._stencil, pdfs, self._symbolOrder0, self._symbolsOrder1)
+ ac = getEquationsForZerothAndFirstOrderMoment(self._stencil, pdfs, self._symbolOrder0, self._symbolsOrder1)
if self._compressible:
- eqColl = divideFirstOrderMomentsByRho(eqColl, dim)
+ ac = divideFirstOrderMomentsByRho(ac, dim)
else:
- eqColl = addDensityOffset(eqColl)
+ ac = addDensityOffset(ac)
- eqColl = applyForceModelShift('macroscopicVelocityShift', dim, eqColl, self._forceModel, self._compressible)
+ ac = applyForceModelShift('macroscopicVelocityShift', dim, ac, self._forceModel, self._compressible)
- mainEquations = []
- eqs = OrderedDict([(eq.lhs, eq.rhs) for eq in eqColl.allEquations])
+ main_assignments = []
+ eqs = OrderedDict([(eq.lhs, eq.rhs) for eq in ac.allEquations])
if 'density' in outputQuantityNamesToSymbols:
- densityOutputSymbol = outputQuantityNamesToSymbols['density']
- if isinstance(densityOutputSymbol, Field):
- densityOutputSymbol = densityOutputSymbol()
- if densityOutputSymbol != self._symbolOrder0:
- mainEquations.append(sp.Eq(densityOutputSymbol, self._symbolOrder0))
+ density_output_symbol = outputQuantityNamesToSymbols['density']
+ if isinstance(density_output_symbol, Field):
+ density_output_symbol = density_output_symbol()
+ if density_output_symbol != self._symbolOrder0:
+ main_assignments.append(Assignment(density_output_symbol, self._symbolOrder0))
else:
- mainEquations.append(sp.Eq(self._symbolOrder0, eqs[self._symbolOrder0]))
+ main_assignments.append(Assignment(self._symbolOrder0, eqs[self._symbolOrder0]))
del eqs[self._symbolOrder0]
if 'velocity' in outputQuantityNamesToSymbols:
- velOutputSymbols = outputQuantityNamesToSymbols['velocity']
- if isinstance(velOutputSymbols, Field):
- field = velOutputSymbols
- velOutputSymbols = [field(i) for i in range(len(self._symbolsOrder1))]
- if tuple(velOutputSymbols) != tuple(self._symbolsOrder1):
- mainEquations += [sp.Eq(a, b) for a, b in zip(velOutputSymbols, self._symbolsOrder1)]
+ vel_output_symbols = outputQuantityNamesToSymbols['velocity']
+ if isinstance(vel_output_symbols, Field):
+ field = vel_output_symbols
+ vel_output_symbols = [field(i) for i in range(len(self._symbolsOrder1))]
+ if tuple(vel_output_symbols) != tuple(self._symbolsOrder1):
+ main_assignments += [Assignment(a, b) for a, b in zip(vel_output_symbols, self._symbolsOrder1)]
else:
for u_i in self._symbolsOrder1:
- mainEquations.append(sp.Eq(u_i, eqs[u_i]))
+ main_assignments.append(Assignment(u_i, eqs[u_i]))
del eqs[u_i]
if 'momentumDensity' in outputQuantityNamesToSymbols:
# get zeroth and first moments again - force-shift them if necessary
# and add their values directly to the main equations assuming that subexpressions are already in
# main equation collection
# Is not optimal when velocity and momentumDensity are calculated together, but this is usually not the case
- momentumDensityOutputSymbols = outputQuantityNamesToSymbols['momentumDensity']
- momDensityEqColl = getEquationsForZerothAndFirstOrderMoment(self._stencil, pdfs,
+ momentum_density_output_symbols = outputQuantityNamesToSymbols['momentumDensity']
+ mom_density_eq_coll = getEquationsForZerothAndFirstOrderMoment(self._stencil, pdfs,
self._symbolOrder0, self._symbolsOrder1)
- momDensityEqColl = applyForceModelShift('macroscopicVelocityShift', dim, momDensityEqColl,
+ mom_density_eq_coll = applyForceModelShift('macroscopicVelocityShift', dim, mom_density_eq_coll,
self._forceModel, self._compressible)
- for sym, val in zip(momentumDensityOutputSymbols, momDensityEqColl.mainEquations[1:]):
- mainEquations.append(sp.Eq(sym, val.rhs))
+ for sym, val in zip(momentum_density_output_symbols, mom_density_eq_coll.mainAssignments[1:]):
+ main_assignments.append(Assignment(sym, val.rhs))
- eqColl = eqColl.copy(mainEquations, [sp.Eq(a, b) for a, b in eqs.items()])
- return eqColl.newWithoutUnusedSubexpressions()
+ ac = ac.copy(main_assignments, [Assignment(a, b) for a, b in eqs.items()])
+ return ac.newWithoutUnusedSubexpressions()
def __repr__(self):
return "ConservedValueComputation for %s" % (", " .join(self.conservedQuantities.keys()),)
@@ -231,7 +232,7 @@ def getEquationsForZerothAndFirstOrderMoment(stencil, symbolicPdfs, symbolicZero
:param symbolicZerothMoment: called :math:`\rho` above
:param symbolicFirstMoments: called :math:`u` above
"""
- def filterOutPlusTerms(expr):
+ def filter_out_plus_terms(expr):
result = 0
for term in expr.args:
if not type(term) is sp.Mul:
@@ -241,34 +242,34 @@ def getEquationsForZerothAndFirstOrderMoment(stencil, symbolicPdfs, symbolicZero
dim = len(stencil[0])
subexpressions = []
- pdfSum = sum(symbolicPdfs)
+ pdf_sum = sum(symbolicPdfs)
u = [0] * dim
for f, offset in zip(symbolicPdfs, stencil):
for i in range(dim):
u[i] += f * int(offset[i])
- plusTerms = [set(filterOutPlusTerms(u_i).args) for u_i in u]
+ plus_terms = [set(filter_out_plus_terms(u_i).args) for u_i in u]
for i in range(dim):
- rhs = plusTerms[i]
+ rhs = plus_terms[i]
for j in range(i):
- rhs -= plusTerms[j]
- eq = sp.Eq(sp.Symbol("vel%dTerm" % (i,)), sum(rhs))
+ rhs -= plus_terms[j]
+ eq = Assignment(sp.Symbol("vel%dTerm" % (i,)), sum(rhs))
subexpressions.append(eq)
for subexpression in subexpressions:
- pdfSum = pdfSum.subs(subexpression.rhs, subexpression.lhs)
+ pdf_sum = pdf_sum.subs(subexpression.rhs, subexpression.lhs)
for i in range(dim):
u[i] = u[i].subs(subexpressions[i].rhs, subexpressions[i].lhs)
equations = []
- equations += [sp.Eq(symbolicZerothMoment, pdfSum)]
- equations += [sp.Eq(u_i_sym, u_i) for u_i_sym, u_i in zip(symbolicFirstMoments, u)]
+ equations += [Assignment(symbolicZerothMoment, pdf_sum)]
+ equations += [Assignment(u_i_sym, u_i) for u_i_sym, u_i in zip(symbolicFirstMoments, u)]
- return EquationCollection(equations, subexpressions)
+ return AssignmentCollection(equations, subexpressions)
-def divideFirstOrderMomentsByRho(equationCollection, dim):
+def divideFirstOrderMomentsByRho(assignment_collection, dim):
"""
Assumes that the equations of the passed equation collection are the following
- rho = f_0 + f_1 + ...
@@ -277,38 +278,39 @@ def divideFirstOrderMomentsByRho(equationCollection, dim):
Returns a new equation collection where the u terms (first order moments) are divided by rho.
The dim parameter specifies the number of first order moments. All subsequent equations are just copied over.
"""
- oldEqs = equationCollection.mainEquations
+ oldEqs = assignment_collection.mainAssignments
rho = oldEqs[0].lhs
- newFirstOrderMomentEq = [sp.Eq(eq.lhs, eq.rhs / rho) for eq in oldEqs[1:dim+1]]
- newEqs = [oldEqs[0]] + newFirstOrderMomentEq + oldEqs[dim+1:]
- return equationCollection.copy(newEqs)
+ new_first_order_moment_eq = [Assignment(eq.lhs, eq.rhs / rho) for eq in oldEqs[1:dim+1]]
+ new_eqs = [oldEqs[0]] + new_first_order_moment_eq + oldEqs[dim+1:]
+ return assignment_collection.copy(new_eqs)
-def addDensityOffset(equationCollection, offset=sp.Rational(1, 1)):
+def addDensityOffset(assignment_collection, offset=sp.Rational(1, 1)):
"""
Assumes that first equation is the density (zeroth moment). Changes the density equations by adding offset to it.
"""
- oldEqs = equationCollection.mainEquations
- newDensity = sp.Eq(oldEqs[0].lhs, oldEqs[0].rhs + offset)
- return equationCollection.copy([newDensity] + oldEqs[1:])
+ oldEqs = assignment_collection.mainAssignments
+ newDensity = Assignment(oldEqs[0].lhs, oldEqs[0].rhs + offset)
+ return assignment_collection.copy([newDensity] + oldEqs[1:])
-def applyForceModelShift(shiftMemberName, dim, equationCollection, forceModel, compressible, reverse=False):
+def applyForceModelShift(shiftMemberName, dim, assignment_collection, forceModel, compressible, reverse=False):
"""
- Modifies the first order moment equations in equationCollection according to the force model shift.
+ Modifies the first order moment equations in assignment collection according to the force model shift.
It is applied if force model has a method named shiftMemberName. The equations 1: dim+1 of the passed
equation collection are assumed to be the velocity equations.
"""
if forceModel is not None and hasattr(forceModel, shiftMemberName):
- oldEqs = equationCollection.mainEquations
- density = oldEqs[0].lhs if compressible else sp.Rational(1, 1)
- oldVelEqs = oldEqs[1:dim + 1]
- shiftFunc = getattr(forceModel, shiftMemberName)
- velOffsets = shiftFunc(density)
+ old_eqs = assignment_collection.mainAssignments
+ density = old_eqs[0].lhs if compressible else sp.Rational(1, 1)
+ old_vel_eqs = old_eqs[1:dim + 1]
+ shift_func = getattr(forceModel, shiftMemberName)
+ vel_offsets = shift_func(density)
if reverse:
- velOffsets = [-v for v in velOffsets]
- shiftedVelocityEqs = [sp.Eq(oldEq.lhs, oldEq.rhs + offset) for oldEq, offset in zip(oldVelEqs, velOffsets)]
- newEqs = [oldEqs[0]] + shiftedVelocityEqs + oldEqs[dim + 1:]
- return equationCollection.copy(newEqs)
+ vel_offsets = [-v for v in vel_offsets]
+ shifted_velocity_eqs = [Assignment(oldEq.lhs, oldEq.rhs + offset)
+ for oldEq, offset in zip(old_vel_eqs, vel_offsets)]
+ new_eqs = [old_eqs[0]] + shifted_velocity_eqs + old_eqs[dim + 1:]
+ return assignment_collection.copy(new_eqs)
else:
- return equationCollection
+ return assignment_collection
diff --git a/methods/cumulantbased.py b/methods/cumulantbased.py
index 4ed1d9acca259f5e8ceddc1108dd5a694a160797..0a991677bcbdcf6c0a6e4b8634345ce91c2dfaa2 100644
--- a/methods/cumulantbased.py
+++ b/methods/cumulantbased.py
@@ -1,10 +1,11 @@
import sympy as sp
from collections import OrderedDict
+from pystencils import Assignment
+from pystencils.sympyextensions import fastSubs, replaceAdditive
from lbmpy.methods.abstractlbmethod import AbstractLbMethod, LbmCollisionRule, RelaxationInfo
from lbmpy.methods.conservedquantitycomputation import AbstractConservedQuantityComputation
from lbmpy.moments import MOMENT_SYMBOLS, extractMonomials, momentMatrix, monomialToPolynomialTransformationMatrix
from lbmpy.cumulants import cumulantAsFunctionOfRawMoments, rawMomentAsFunctionOfCumulants
-from pystencils.sympyextensions import fastSubs, replaceAdditive
class CumulantBasedLbMethod(AbstractLbMethod):
@@ -100,7 +101,7 @@ class CumulantBasedLbMethod(AbstractLbMethod):
def getEquilibriumTerms(self):
equilibrium = self.getEquilibrium()
- return sp.Matrix([eq.rhs for eq in equilibrium.mainEquations])
+ return sp.Matrix([eq.rhs for eq in equilibrium.mainAssignments])
def getCollisionRule(self, conservedQuantityEquations=None, momentSubexpressions=False,
preCollisionSubexpressions=True, postCollisionSubexpressions=False):
@@ -110,14 +111,16 @@ class CumulantBasedLbMethod(AbstractLbMethod):
def _computeWeights(self):
replacements = self._conservedQuantityComputation.defaultValues
- eqColl = self.getEquilibrium().copyWithSubstitutionsApplied(replacements).insertSubexpressions()
- newMainEqs = [sp.Eq(e.lhs,
- replaceAdditive(e.rhs, 1, sum(self.preCollisionPdfSymbols), requiredMatchReplacement=1.0))
- for e in eqColl.mainEquations]
+ eq = self.getEquilibrium()
+ eqColl = eq.copyWithSubstitutionsApplied(replacements, substituteOnLhs=False).insertSubexpressions()
+ newMainEqs = [Assignment(e.lhs,
+ replaceAdditive(e.rhs, 1, sum(self.preCollisionPdfSymbols),
+ requiredMatchReplacement=1.0))
+ for e in eqColl.mainAssignments]
eqColl = eqColl.copy(newMainEqs)
weights = []
- for eq in eqColl.mainEquations:
+ for eq in eqColl.mainAssignments:
value = eq.rhs.expand()
assert len(value.atoms(sp.Symbol)) == 0, "Failed to compute weights"
weights.append(value)
@@ -133,9 +136,9 @@ class CumulantBasedLbMethod(AbstractLbMethod):
def substituteConservedQuantities(expressions, cqe):
cqe = cqe.insertSubexpressions()
- substitutionDict = {eq.rhs: eq.lhs for eq in cqe.mainEquations}
- density = cqe.mainEquations[0].lhs
- substitutionDict.update({density * eq.rhs: density * eq.lhs for eq in cqe.mainEquations[1:]})
+ substitutionDict = {eq.rhs: eq.lhs for eq in cqe.mainAssignments}
+ density = cqe.mainAssignments[0].lhs
+ substitutionDict.update({density * eq.rhs: density * eq.lhs for eq in cqe.mainAssignments[1:]})
return [fastSubs(e, substitutionDict) for e in expressions]
f = self.preCollisionPdfSymbols
@@ -161,7 +164,7 @@ class CumulantBasedLbMethod(AbstractLbMethod):
moments = substituteConservedQuantities(moments, conservedQuantityEquations)
if momentSubexpressions:
symbols = [tupleToSymbol(t, "m") for t in higherOrderIndices]
- subexpressions += [sp.Eq(sym, moment) for sym, moment in zip(symbols, moments[numLowerOrderIndices:])]
+ subexpressions += [Assignment(sym, moment) for sym, moment in zip(symbols, moments[numLowerOrderIndices:])]
moments = moments[:numLowerOrderIndices] + symbols
# 3) Transform moments to monomial cumulants
@@ -170,7 +173,7 @@ class CumulantBasedLbMethod(AbstractLbMethod):
if preCollisionSubexpressions:
symbols = [tupleToSymbol(t, "preC") for t in higherOrderIndices]
- subexpressions += [sp.Eq(sym, c)
+ subexpressions += [Assignment(sym, c)
for sym, c in zip(symbols, monomialCumulants[numLowerOrderIndices:])]
monomialCumulants = monomialCumulants[:numLowerOrderIndices] + symbols
@@ -183,7 +186,7 @@ class CumulantBasedLbMethod(AbstractLbMethod):
if postCollisionSubexpressions:
symbols = [tupleToSymbol(t, "postC") for t in higherOrderIndices]
- subexpressions += [sp.Eq(sym, c)
+ subexpressions += [Assignment(sym, c)
for sym, c in zip(symbols, relaxedMonomialCumulants[numLowerOrderIndices:])]
relaxedMonomialCumulants = relaxedMonomialCumulants[:numLowerOrderIndices] + symbols
@@ -191,20 +194,20 @@ class CumulantBasedLbMethod(AbstractLbMethod):
cumulantDict = {idx: value for idx, value in zip(indices, relaxedMonomialCumulants)}
collidedMoments = [rawMomentAsFunctionOfCumulants(idx, cumulantDict) for idx in indices]
result = momentTransformationMatrix.inv() * sp.Matrix(collidedMoments)
- mainEquations = [sp.Eq(sym, val) for sym, val in zip(self.postCollisionPdfSymbols, result)]
+ mainAssignments = [Assignment(sym, val) for sym, val in zip(self.postCollisionPdfSymbols, result)]
# 6) Add forcing terms
if self._forceModel is not None and includeForceTerms:
forceModelTerms = self._forceModel(self)
forceTermSymbols = sp.symbols("forceTerm_:%d" % (len(forceModelTerms,)))
- forceSubexpressions = [sp.Eq(sym, forceModelTerm)
+ forceSubexpressions = [Assignment(sym, forceModelTerm)
for sym, forceModelTerm in zip(forceTermSymbols, forceModelTerms)]
subexpressions += forceSubexpressions
- mainEquations = [sp.Eq(eq.lhs, eq.rhs + forceTermSymbol)
- for eq, forceTermSymbol in zip(mainEquations, forceTermSymbols)]
+ mainAssignments = [Assignment(eq.lhs, eq.rhs + forceTermSymbol)
+ for eq, forceTermSymbol in zip(mainAssignments, forceTermSymbols)]
sh = {'relaxationRates': list(self.relaxationRates)}
- return LbmCollisionRule(self, mainEquations, subexpressions, simplificationHints=sh)
+ return LbmCollisionRule(self, mainAssignments, subexpressions, simplificationHints=sh)
diff --git a/methods/entropic.py b/methods/entropic.py
index daa98c59ec9f38db0c2da7d8dedf18fd5e936111..4703b0c64d4f5d67db0db89b9963452881242dc5 100644
--- a/methods/entropic.py
+++ b/methods/entropic.py
@@ -1,4 +1,5 @@
import sympy as sp
+from pystencils import Assignment
from pystencils.transformations import fastSubs
from lbmpy.relaxationrates import getShearRelaxationRate
@@ -45,20 +46,21 @@ def addEntropyCondition(collisionRule, omegaOutputField=None):
dhSymbols = [sp.Symbol("entropicDh_%d" % (i,)) for i in range(Q)]
feqSymbols = [sp.Symbol("entropicFeq_%d" % (i,)) for i in range(Q)]
- subexprs = [sp.Eq(a, b) for a, b in zip(dsSymbols, ds)] + \
- [sp.Eq(a, b) for a, b in zip(dhSymbols, dh)] + \
- [sp.Eq(a, f_i + ds_i + dh_i) for a, f_i, ds_i, dh_i in zip(feqSymbols, fSymbols, dsSymbols, dhSymbols)]
+ subexprs = [Assignment(a, b) for a, b in zip(dsSymbols, ds)] + \
+ [Assignment(a, b) for a, b in zip(dhSymbols, dh)] + \
+ [Assignment(a, f_i + ds_i + dh_i) for a, f_i, ds_i, dh_i in
+ zip(feqSymbols, fSymbols, dsSymbols, dhSymbols)]
optimalOmegaH = _getEntropyMaximizingOmega(omega_s, feqSymbols, dsSymbols, dhSymbols)
- subexprs += [sp.Eq(omega_h, optimalOmegaH)]
+ subexprs += [Assignment(omega_h, optimalOmegaH)]
newUpdateEquations = []
constPart = decomp.constantExprs()
- for updateEq in collisionRule.mainEquations:
+ for updateEq in collisionRule.mainAssignments:
index = collisionRule.method.postCollisionPdfSymbols.index(updateEq.lhs)
- newEq = sp.Eq(updateEq.lhs, constPart[index] + omega_s * dsSymbols[index] + omega_h * dhSymbols[index])
+ newEq = Assignment(updateEq.lhs, constPart[index] + omega_s * dsSymbols[index] + omega_h * dhSymbols[index])
newUpdateEquations.append(newEq)
newCollisionRule = collisionRule.copy(newUpdateEquations, collisionRule.subexpressions + subexprs)
newCollisionRule.simplificationHints['entropic'] = True
@@ -101,7 +103,7 @@ def addIterativeEntropyCondition(collisionRule, freeOmega=None, newtonIterations
polynomialSubexpressions = []
rrPolynomials = []
for i, constantExpr in enumerate(decomp.constantExprs()):
- constantExprEq = sp.Eq(decomp.symbolicConstantExpr(i), constantExpr)
+ constantExprEq = Assignment(decomp.symbolicConstantExpr(i), constantExpr)
polynomialSubexpressions.append(constantExprEq)
rrPolynomial = constantExprEq.lhs
@@ -109,21 +111,21 @@ def addIterativeEntropyCondition(collisionRule, freeOmega=None, newtonIterations
for idx, f in enumerate(factors[i]):
power = idx + 1
symbolicFactor = decomp.symbolicRelaxationRateFactors(freeOmega, power)[i]
- polynomialSubexpressions.append(sp.Eq(symbolicFactor, f))
+ polynomialSubexpressions.append(Assignment(symbolicFactor, f))
rrPolynomial += freeOmega ** power * symbolicFactor
rrPolynomials.append(rrPolynomial)
- newUpdateEquations.append(sp.Eq(collisionRule.method.postCollisionPdfSymbols[i], rrPolynomial))
+ newUpdateEquations.append(Assignment(collisionRule.method.postCollisionPdfSymbols[i], rrPolynomial))
# 2) get equilibrium from method and define subexpressions for it
- eqTerms = [eq.rhs for eq in collisionRule.method.getEquilibrium().mainEquations]
+ eqTerms = [eq.rhs for eq in collisionRule.method.getEquilibrium().mainAssignments]
eqSymbols = sp.symbols("entropicFeq_:%d" % (len(eqTerms,)))
- eqSubexpressions = [sp.Eq(a, b) for a, b in zip(eqSymbols, eqTerms)]
+ eqSubexpressions = [Assignment(a, b) for a, b in zip(eqSymbols, eqTerms)]
# 3) find coefficients of entropy derivatives
entropyDiff = sp.diff(discreteApproxEntropy(rrPolynomials, eqSymbols), freeOmega)
coefficientsFirstDiff = [c.expand() for c in reversed(sp.poly(entropyDiff, freeOmega).all_coeffs())]
symCoeffDiff1 = sp.symbols("entropicDiffCoeff_:%d" % (len(coefficientsFirstDiff,)))
- coefficientEqs = [sp.Eq(a, b) for a, b in zip(symCoeffDiff1, coefficientsFirstDiff)]
+ coefficientEqs = [Assignment(a, b) for a, b in zip(symCoeffDiff1, coefficientsFirstDiff)]
symCoeffDiff2 = [(i+1) * coeff for i, coeff in enumerate(symCoeffDiff1[1:])]
# 4) define Newtons method update iterations
@@ -136,7 +138,7 @@ def addIterativeEntropyCondition(collisionRule, freeOmega=None, newtonIterations
lhsOmega = intermediateOmegas[omega_idx+1]
diff1Poly = sum([coeff * rhsOmega**i for i, coeff in enumerate(symCoeffDiff1)])
diff2Poly = sum([coeff * rhsOmega**i for i, coeff in enumerate(symCoeffDiff2)])
- newtonEq = sp.Eq(lhsOmega, rhsOmega - diff1Poly / diff2Poly)
+ newtonEq = Assignment(lhsOmega, rhsOmega - diff1Poly / diff2Poly)
newtonIterationEquations.append(newtonEq)
# 5) final update equations
@@ -204,7 +206,7 @@ class RelaxationRatePolynomialDecomposition(object):
return [sp.Symbol("entFacOmega_%d_%d_%d" % (i, omegaIdx, power)) for i in range(Q)]
def relaxationRateFactors(self, relaxationRate):
- updateEquations = self._collisionRule.mainEquations
+ updateEquations = self._collisionRule.mainAssignments
result = []
for updateEquation in updateEquations:
@@ -231,13 +233,13 @@ class RelaxationRatePolynomialDecomposition(object):
def constantExprs(self):
subsDict = {rr: 0 for rr in self._freeRelaxationRates}
subsDict.update({rr: 0 for rr in self._fixedRelaxationRates})
- updateEquations = self._collisionRule.mainEquations
+ updateEquations = self._collisionRule.mainAssignments
return [fastSubs(eq.rhs, subsDict) for eq in updateEquations]
def equilibriumExprs(self):
subsDict = {rr: 1 for rr in self._freeRelaxationRates}
subsDict.update({rr: 1 for rr in self._fixedRelaxationRates})
- updateEquations = self._collisionRule.mainEquations
+ updateEquations = self._collisionRule.mainAssignments
return [fastSubs(eq.rhs, subsDict) for eq in updateEquations]
def symbolicEquilibrium(self):
diff --git a/methods/entropic_eq_srt.py b/methods/entropic_eq_srt.py
index e3f3dbef26801e68ff2f0dd848142ec4a4105fc0..1b90716873e656863a8b5daa9ee6982a33a7ba8e 100644
--- a/methods/entropic_eq_srt.py
+++ b/methods/entropic_eq_srt.py
@@ -1,4 +1,5 @@
import sympy as sp
+from pystencils import Assignment
from lbmpy.maxwellian_equilibrium import getWeights
from lbmpy.methods.abstractlbmethod import AbstractLbMethod, LbmCollisionRule
from lbmpy.methods.conservedquantitycomputation import DensityVelocityComputation
@@ -52,16 +53,16 @@ class EntropicEquilibriumSRT(AbstractLbMethod):
f_i *= (2 - B) * ((2 * u_a + B) / (1 - u_a)) ** e_ia
eq.append(f_i)
- collisionEqs = [sp.Eq(lhs, (1 - relaxationRate) * f_i + relaxationRate * eq_i)
+ collisionEqs = [Assignment(lhs, (1 - relaxationRate) * f_i + relaxationRate * eq_i)
for lhs, f_i, eq_i in zip(self.postCollisionPdfSymbols, self.preCollisionPdfSymbols, eq)]
if (self._forceModel is not None) and includeForceTerms:
forceModelTerms = self._forceModel(self)
forceTermSymbols = sp.symbols("forceTerm_:%d" % (len(forceModelTerms, )))
- forceSubexpressions = [sp.Eq(sym, forceModelTerm)
+ forceSubexpressions = [Assignment(sym, forceModelTerm)
for sym, forceModelTerm in zip(forceTermSymbols, forceModelTerms)]
allSubexpressions += forceSubexpressions
- collisionEqs = [sp.Eq(eq.lhs, eq.rhs + forceTermSymbol)
+ collisionEqs = [Assignment(eq.lhs, eq.rhs + forceTermSymbol)
for eq, forceTermSymbol in zip(collisionEqs, forceTermSymbols)]
return LbmCollisionRule(self, collisionEqs, allSubexpressions)
diff --git a/methods/momentbased.py b/methods/momentbased.py
index 407403dca58f4d991f711ebdfdf22c2b8f11be7e..f277400a29e2fe95d714d99bf1b41ae15b1ec9ec 100644
--- a/methods/momentbased.py
+++ b/methods/momentbased.py
@@ -1,10 +1,10 @@
import sympy as sp
from collections import OrderedDict
-
+from pystencils import Assignment
+from pystencils.sympyextensions import replaceAdditive
from lbmpy.methods.abstractlbmethod import AbstractLbMethod, LbmCollisionRule, RelaxationInfo
from lbmpy.methods.conservedquantitycomputation import AbstractConservedQuantityComputation
from lbmpy.moments import MOMENT_SYMBOLS, momentMatrix
-from pystencils.sympyextensions import replaceAdditive
class MomentBasedLbMethod(AbstractLbMethod):
@@ -81,7 +81,7 @@ class MomentBasedLbMethod(AbstractLbMethod):
def getEquilibriumTerms(self):
equilibrium = self.getEquilibrium()
- return sp.Matrix([eq.rhs for eq in equilibrium.mainEquations])
+ return sp.Matrix([eq.rhs for eq in equilibrium.mainAssignments])
def getCollisionRule(self, conservedQuantityEquations=None):
D = sp.diag(*self.relaxationRates)
@@ -156,13 +156,13 @@ class MomentBasedLbMethod(AbstractLbMethod):
eqColl = self.getEquilibrium(includeForceTerms=False)
eqColl = eqColl.copyWithSubstitutionsApplied(replacements, substituteOnLhs=False).insertSubexpressions()
- newMainEqs = [sp.Eq(e.lhs,
+ newMainEqs = [Assignment(e.lhs,
replaceAdditive(e.rhs, 1, sum(self.preCollisionPdfSymbols), requiredMatchReplacement=1.0))
- for e in eqColl.mainEquations]
+ for e in eqColl.mainAssignments]
eqColl = eqColl.copy(newMainEqs)
weights = []
- for eq in eqColl.mainEquations:
+ for eq in eqColl.mainAssignments:
value = eq.rhs.expand()
assert len(value.atoms(sp.Symbol)) == 0, "Failed to compute weights " + str(value)
weights.append(value)
@@ -175,7 +175,7 @@ class MomentBasedLbMethod(AbstractLbMethod):
m_eq = sp.Matrix(self.momentEquilibriumValues)
collisionRule = f + M.inv() * D * (m_eq - M * f)
- collisionEqs = [sp.Eq(lhs, rhs) for lhs, rhs in zip(self.postCollisionPdfSymbols, collisionRule)]
+ collisionEqs = [Assignment(lhs, rhs) for lhs, rhs in zip(self.postCollisionPdfSymbols, collisionRule)]
if conservedQuantityEquations is None:
conservedQuantityEquations = self._conservedQuantityComputation.equilibriumInputEquationsFromPdfs(f)
@@ -189,10 +189,10 @@ class MomentBasedLbMethod(AbstractLbMethod):
if self._forceModel is not None and includeForceTerms:
forceModelTerms = self._forceModel(self)
forceTermSymbols = sp.symbols("forceTerm_:%d" % (len(forceModelTerms,)))
- forceSubexpressions = [sp.Eq(sym, forceModelTerm)
+ forceSubexpressions = [Assignment(sym, forceModelTerm)
for sym, forceModelTerm in zip(forceTermSymbols, forceModelTerms)]
allSubexpressions += forceSubexpressions
- collisionEqs = [sp.Eq(eq.lhs, eq.rhs + forceTermSymbol)
+ collisionEqs = [Assignment(eq.lhs, eq.rhs + forceTermSymbol)
for eq, forceTermSymbol in zip(collisionEqs, forceTermSymbols)]
simplificationHints['forceTerms'] = forceTermSymbols
@@ -219,7 +219,7 @@ class MomentBasedLbMethod(AbstractLbMethod):
newRR = [subexpressions[sp.sympify(e)] if sp.sympify(e) in subexpressions else e
for e in rr]
- substitutions = [sp.Eq(e[1], e[0]) for e in subexpressions.items()]
+ substitutions = [Assignment(e[1], e[0]) for e in subexpressions.items()]
return substitutions, sp.diag(*newRR)
else:
return [], relaxationMatrix
diff --git a/methods/momentbasedsimplifications.py b/methods/momentbasedsimplifications.py
index fcf3ee14f957719f7aa42488eacd90232a96c1cf..9d727d7a3473dab72f73f8bd0921eec8369ad130 100644
--- a/methods/momentbasedsimplifications.py
+++ b/methods/momentbasedsimplifications.py
@@ -1,9 +1,10 @@
"""
This module holds special transformations for simplifying the collision equations of moment-based methods.
-All of these transformations operate on :class:`pystencils.equationcollection.EquationCollection` and need special
+All of these transformations operate on :class:`pystencils.AssignmentCollection` and need special
simplification hints, which are set by the MomentBasedLbMethod.
"""
import sympy as sp
+from pystencils import Assignment
from pystencils.sympyextensions import replaceAdditive, replaceSecondOrderProducts, extractMostCommonFactor
@@ -19,9 +20,9 @@ def replaceSecondOrderVelocityProducts(lbmCollisionEqs):
result = []
substitutions = []
u = sh['velocity']
- for i, s in enumerate(lbmCollisionEqs.mainEquations):
+ for i, s in enumerate(lbmCollisionEqs.mainAssignments):
newRhs = replaceSecondOrderProducts(s.rhs, u, positive=None, replaceMixed=substitutions)
- result.append(sp.Eq(s.lhs, newRhs))
+ result.append(Assignment(s.lhs, newRhs))
res = lbmCollisionEqs.copy(result)
res.subexpressions += substitutions
return res
@@ -41,11 +42,11 @@ def factorRelaxationRates(lbmCollisionEqs):
relaxationRates = sp.Matrix(sh['relaxationRates']).atoms(sp.Symbol)
result = []
- for s in lbmCollisionEqs.mainEquations:
+ for s in lbmCollisionEqs.mainAssignments:
newRhs = s.rhs
for rp in relaxationRates:
newRhs = newRhs.collect(rp)
- result.append(sp.Eq(s.lhs, newRhs))
+ result.append(Assignment(s.lhs, newRhs))
return lbmCollisionEqs.copy(result)
@@ -67,12 +68,12 @@ def factorDensityAfterFactoringRelaxationTimes(lbmCollisionEqs):
relaxationRates = sp.Matrix(sh['relaxationRates']).atoms(sp.Symbol)
result = []
rho = sh['density']
- for s in lbmCollisionEqs.mainEquations:
+ for s in lbmCollisionEqs.mainAssignments:
newRhs = s.rhs
for rp in relaxationRates:
coeff = newRhs.coeff(rp)
newRhs = newRhs.subs(coeff, coeff.collect(rho))
- result.append(sp.Eq(s.lhs, newRhs))
+ result.append(Assignment(s.lhs, newRhs))
return lbmCollisionEqs.copy(result)
@@ -89,13 +90,13 @@ def replaceDensityAndVelocity(lbmCollisionEqs):
rho = sh['density']
u = sh['velocity']
- substitutions = lbmCollisionEqs.extract([rho] + list(u)).insertSubexpressions().mainEquations
+ substitutions = lbmCollisionEqs.extract([rho] + list(u)).insertSubexpressions().mainAssignments
result = []
- for s in lbmCollisionEqs.mainEquations:
+ for s in lbmCollisionEqs.mainAssignments:
newRhs = s.rhs
for replacement in substitutions:
newRhs = replaceAdditive(newRhs, replacement.lhs, replacement.rhs, requiredMatchReplacement=0.5)
- result.append(sp.Eq(s.lhs, newRhs))
+ result.append(Assignment(s.lhs, newRhs))
return lbmCollisionEqs.copy(result)
@@ -122,11 +123,11 @@ def replaceCommonQuadraticAndConstantTerm(lbmCollisionEqs):
return lbmCollisionEqs
if len(f_eq_common.args) > 1:
- f_eq_common = sp.Eq(sp.Symbol('f_eq_common'), f_eq_common)
+ f_eq_common = Assignment(sp.Symbol('f_eq_common'), f_eq_common)
result = []
- for s in lbmCollisionEqs.mainEquations:
+ for s in lbmCollisionEqs.mainAssignments:
newRhs = replaceAdditive(s.rhs, f_eq_common.lhs, f_eq_common.rhs, requiredMatchReplacement=0.5)
- result.append(sp.Eq(s.lhs, newRhs))
+ result.append(Assignment(s.lhs, newRhs))
res = lbmCollisionEqs.copy(result)
res.subexpressions.append(f_eq_common)
return res
@@ -145,7 +146,7 @@ def cseInOpposingDirections(lbmCollisionEqs):
sh = lbmCollisionEqs.simplificationHints
assert 'relaxationRates' in sh, "Needs simplification hint 'relaxationRates': Sequence of relaxation rates"
- updateRules = lbmCollisionEqs.mainEquations
+ updateRules = lbmCollisionEqs.mainAssignments
stencil = lbmCollisionEqs.method.stencil
relaxationRates = sp.Matrix(sh['relaxationRates']).atoms(sp.Symbol)
@@ -171,7 +172,7 @@ def cseInOpposingDirections(lbmCollisionEqs):
if len(relaxationRates) == 0:
foundSubexpressions, newTerms = sp.cse(updateRules, symbols=replacementSymbolGenerator,
order='None', optimizations=[])
- substitutions += [sp.Eq(f[0], f[1]) for f in foundSubexpressions]
+ substitutions += [Assignment(f[0], f[1]) for f in foundSubexpressions]
updateRules = newTerms
else:
@@ -193,15 +194,15 @@ def cseInOpposingDirections(lbmCollisionEqs):
foundSubexpressions, newTerms = sp.cse(handledTerms, symbols=replacementSymbolGenerator,
order='None', optimizations=[])
- substitutions += [sp.Eq(f[0], f[1]) for f in foundSubexpressions]
+ substitutions += [Assignment(f[0], f[1]) for f in foundSubexpressions]
- updateRules = [sp.Eq(ur.lhs, ur.rhs.subs(relaxationRate * oldTerm, newCoefficient * newTerm))
+ updateRules = [Assignment(ur.lhs, ur.rhs.subs(relaxationRate * oldTerm, newCoefficient * newTerm))
for ur, newTerm, oldTerm in zip(updateRules, newTerms, terms)]
result += updateRules
for term, substitutedVar in newCoefficientSubstitutions.items():
- substitutions.append(sp.Eq(substitutedVar, term))
+ substitutions.append(Assignment(substitutedVar, term))
result.sort(key=lambda e: lbmCollisionEqs.method.postCollisionPdfSymbols.index(e.lhs))
res = lbmCollisionEqs.copy(result)
@@ -223,7 +224,7 @@ def __getCommonQuadraticAndConstantTerms(lbmCollisionEqs):
pdfSymbols = lbmCollisionEqs.freeSymbols - relaxationRates
center = tuple([0] * dim)
- t = lbmCollisionEqs.mainEquations[stencil.index(center)].rhs
+ t = lbmCollisionEqs.mainAssignments[stencil.index(center)].rhs
for rp in relaxationRates:
t = t.subs(rp, 1)
diff --git a/phasefield/kerneleqs.py b/phasefield/kerneleqs.py
index 008dda4d3b677c91fd653c802f956f2936503be1..7b7184d30a5e22f1fb8866a36a915e5bb9bf97fc 100644
--- a/phasefield/kerneleqs.py
+++ b/phasefield/kerneleqs.py
@@ -1,4 +1,5 @@
import sympy as sp
+from pystencils import Assignment
from pystencils.finitedifferences import Discretization2ndOrder
from lbmpy.phasefield.analytical import chemicalPotentialsFromFreeEnergy, substituteLaplacianBySum, \
forceFromPhiAndMu, symmetricTensorLinearization, pressureTensorFromFreeEnergy, forceFromPressureTensor
@@ -15,7 +16,7 @@ def muKernel(freeEnergy, orderParameters, phiField, muField, dx=1):
chemicalPotential = substituteLaplacianBySum(chemicalPotential, dim)
chemicalPotential = chemicalPotential.subs({op: phiField(i) for i, op in enumerate(orderParameters)})
discretize = Discretization2ndOrder(dx=dx)
- return [sp.Eq(muField(i), discretize(mu_i)) for i, mu_i in enumerate(chemicalPotential)]
+ return [Assignment(muField(i), discretize(mu_i)) for i, mu_i in enumerate(chemicalPotential)]
def forceKernelUsingMu(forceField, phiField, muField, dx=1):
@@ -23,7 +24,7 @@ def forceKernelUsingMu(forceField, phiField, muField, dx=1):
assert muField.indexDimensions == 1
force = forceFromPhiAndMu(phiField.vecCenter, mu=muField.vecCenter, dim=muField.spatialDimensions)
discretize = Discretization2ndOrder(dx=dx)
- return [sp.Eq(forceField(i),
+ return [Assignment(forceField(i),
discretize(f_i)).expand() for i, f_i in enumerate(force)]
@@ -35,8 +36,7 @@ def pressureTensorKernel(freeEnergy, orderParameters, phiField, pressureTensorFi
discretize = Discretization2ndOrder(dx=dx)
eqs = []
for index, linIndex in indexMap.items():
- eq = sp.Eq(pressureTensorField(linIndex),
- discretize(p[index]).expand())
+ eq = Assignment(pressureTensorField(linIndex), discretize(p[index]).expand())
eqs.append(eq)
return eqs
@@ -50,7 +50,7 @@ def forceKernelUsingPressureTensor(forceField, pressureTensorField, extraForce=N
if extraForce:
f += extraForce
discretize = Discretization2ndOrder(dx=dx)
- return [sp.Eq(forceField(i), discretize(f_i).expand())
+ return [Assignment(forceField(i), discretize(f_i).expand())
for i, f_i in enumerate(f)]
@@ -81,7 +81,7 @@ class CahnHilliardFDStep:
updateEqs = []
for i in range(numPhases):
rhs = cahnHilliardFdEq(i, self.phiField, muField, velField, mobilities[i], dx, dt)
- updateEqs.append(sp.Eq(self.tmpField(i), rhs))
+ updateEqs.append(Assignment(self.tmpField(i), rhs))
self.updateEqs = updateEqs
self.updateEqs = equationModifier(updateEqs)
self.kernel = createKernel(self.updateEqs, target=target).compile()
diff --git a/phasefield/phasefieldstep.py b/phasefield/phasefieldstep.py
index dc52d57cce75b9d6c94ad96df8a55c4ebbc3dc46..167f5320e8b522629b30ad0d513bd585a4b82647 100644
--- a/phasefield/phasefieldstep.py
+++ b/phasefield/phasefieldstep.py
@@ -10,7 +10,7 @@ from lbmpy.phasefield.analytical import chemicalPotentialsFromFreeEnergy, symmet
from pystencils.boundaries.boundaryhandling import FlagInterface
from pystencils.boundaries.inkernel import addNeumannBoundary
from pystencils.datahandling import SerialDataHandling
-from pystencils.equationcollection.simplifications import sympyCseOnEquationList
+from pystencils.assignment_collection.simplifications import sympyCseOnEquationList
from pystencils.slicing import makeSlice, SlicedGetter
diff --git a/simplificationfactory.py b/simplificationfactory.py
index 64eb7e9d770fd51be1ca3ec889bd70efee78d814..d93f6f5f99f85406539b40f4d064c75fd7da4f74 100644
--- a/simplificationfactory.py
+++ b/simplificationfactory.py
@@ -3,12 +3,12 @@ import sympy as sp
from lbmpy.innerloopsplit import createLbmSplitGroups
from lbmpy.methods.cumulantbased import CumulantBasedLbMethod
-from pystencils.equationcollection.simplifications import applyOnAllEquations, \
- subexpressionSubstitutionInMainEquations, sympyCSE, addSubexpressionsForDivisions
+from pystencils.assignment_collection.simplifications import applyOnAllEquations, \
+ subexpressionSubstitutionInmainAssignments, sympyCSE, addSubexpressionsForDivisions
def createSimplificationStrategy(lbmMethod, doCseInOpposingDirections=False, doOverallCse=False, splitInnerLoop=False):
- from pystencils.equationcollection import SimplificationStrategy
+ from pystencils.assignment_collection import SimplificationStrategy
from lbmpy.methods import MomentBasedLbMethod
from lbmpy.methods.momentbasedsimplifications import replaceSecondOrderVelocityProducts, \
factorDensityAfterFactoringRelaxationTimes, factorRelaxationRates, cseInOpposingDirections, \
@@ -27,7 +27,7 @@ def createSimplificationStrategy(lbmMethod, doCseInOpposingDirections=False, doO
s.add(replaceDensityAndVelocity)
s.add(replaceCommonQuadraticAndConstantTerm)
s.add(factorDensityAfterFactoringRelaxationTimes)
- s.add(subexpressionSubstitutionInMainEquations)
+ s.add(subexpressionSubstitutionInmainAssignments)
if splitInnerLoop:
s.add(createLbmSplitGroups)
elif isinstance(lbmMethod, CumulantBasedLbMethod):
diff --git a/turbulence_models.py b/turbulence_models.py
index 22038f2c779bd4b7c226307d2ce13179d8b9b15a..379fc5f1f50a1c80c79dd25905c4ccfa01171292 100644
--- a/turbulence_models.py
+++ b/turbulence_models.py
@@ -1,5 +1,6 @@
import sympy as sp
+from pystencils import Assignment
from lbmpy.relaxationrates import getShearRelaxationRate
@@ -27,13 +28,13 @@ def addSmagorinskyModel(collisionRule, smagorinskyConstant, omegaOutputField=Non
adaptedOmega = sp.Symbol("smagorinskyOmega")
# for derivation see notebook demo_custom_LES_model.pynb
- eqs = [sp.Eq(tau_0, 1 / omega_s),
- sp.Eq(secondOrderNeqMoments, frobeniusNorm(secondOrderMomentTensor(fNeq, method.stencil), factor=2) / rho),
- sp.Eq(adaptedOmega, 2 / (tau_0 + sp.sqrt(18 * smagorinskyConstant**2 * secondOrderNeqMoments + tau_0**2)))]
+ eqs = [Assignment(tau_0, 1 / omega_s),
+ Assignment(secondOrderNeqMoments, frobeniusNorm(secondOrderMomentTensor(fNeq, method.stencil), factor=2) / rho),
+ Assignment(adaptedOmega, 2 / (tau_0 + sp.sqrt(18 * smagorinskyConstant**2 * secondOrderNeqMoments + tau_0**2)))]
collisionRule.subexpressions += eqs
- collisionRule.topologicalSort(subexpressions=True, mainEquations=False)
+ collisionRule.topologicalSort(subexpressions=True, mainAssignments=False)
if omegaOutputField:
- collisionRule.mainEquations.append(sp.Eq(omegaOutputField.center, adaptedOmega))
+ collisionRule.mainAssignments.append(Assignment(omegaOutputField.center, adaptedOmega))
return collisionRule
diff --git a/updatekernels.py b/updatekernels.py
index 84ca738f64b775a31c65bf72bca428f0d069194c..dba48166891dbbe33f6a3e4089785ae8d6faeb0e 100644
--- a/updatekernels.py
+++ b/updatekernels.py
@@ -1,13 +1,12 @@
import numpy as np
import sympy as sp
-from lbmpy.methods.abstractlbmethod import LbmCollisionRule
-from pystencils import Field
-from pystencils.equationcollection.equationcollection import EquationCollection
+from pystencils import Field, Assignment
+from pystencils.assignment_collection.assignment_collection import AssignmentCollection
from pystencils.field import createNumpyArrayWithLayout, layoutStringToTuple
from pystencils.sympyextensions import fastSubs
-from lbmpy.fieldaccess import StreamPullTwoFieldsAccessor, Pseudo2DTwoFieldsAccessor, PeriodicTwoFieldsAccessor, \
- CollideOnlyInplaceAccessor
+from lbmpy.methods.abstractlbmethod import LbmCollisionRule
+from lbmpy.fieldaccess import StreamPullTwoFieldsAccessor, PeriodicTwoFieldsAccessor, CollideOnlyInplaceAccessor
# -------------------------------------------- LBM Kernel Creation -----------------------------------------------------
@@ -117,9 +116,9 @@ def createStreamPullOnlyKernel(stencil, numpyField=None, srcFieldName="src", dst
dst = Field.createFromNumpyArray(dstFieldName, numpyField, indexDimensions=1)
accessor = StreamPullTwoFieldsAccessor()
- eqs = [sp.Eq(a, b) for a, b in zip(accessor.write(dst, stencil), accessor.read(src, stencil))
- if sp.Eq(a, b) != True]
- return EquationCollection(eqs, [])
+ eqs = [Assignment(a, b) for a, b in zip(accessor.write(dst, stencil), accessor.read(src, stencil))
+ if Assignment(a, b) != True]
+ return AssignmentCollection(eqs, [])
def createStreamPullWithOutputKernel(lbMethod, srcField, dstField, output):
@@ -127,13 +126,13 @@ def createStreamPullWithOutputKernel(lbMethod, srcField, dstField, output):
cqc = lbMethod.conservedQuantityComputation
streamed = sp.symbols("streamed_:%d" % (len(stencil),))
accessor = StreamPullTwoFieldsAccessor()
- streamEqs = [sp.Eq(a, b) for a, b in zip(streamed, accessor.read(srcField, stencil))
- if sp.Eq(a, b) != True]
+ streamEqs = [Assignment(a, b) for a, b in zip(streamed, accessor.read(srcField, stencil))
+ if Assignment(a, b) != True]
outputEqCollection = cqc.outputEquationsFromPdfs(streamed, output)
- writeEqs = [sp.Eq(a, b) for a, b in zip(accessor.write(dstField, stencil), streamed)]
+ writeEqs = [Assignment(a, b) for a, b in zip(accessor.write(dstField, stencil), streamed)]
subExprs = streamEqs + outputEqCollection.subexpressions
- mainEqs = outputEqCollection.mainEquations + writeEqs
+ mainEqs = outputEqCollection.mainAssignments + writeEqs
return LbmCollisionRule(lbMethod, mainEqs, subExprs, simplificationHints=outputEqCollection.simplificationHints)
# ---------------------------------- Pdf array creation for various layouts --------------------------------------------
@@ -173,5 +172,5 @@ def addOutputFieldForConservedQuantities(collisionRule, conservedQuantitiesToOut
def writeQuantitiesToField(collisionRule, symbols, outputField):
if not hasattr(symbols, "__len__"):
symbols = [symbols]
- eqs = [sp.Eq(outputField(i), s) for i, s in enumerate(symbols)]
- return collisionRule.copy(collisionRule.mainEquations + eqs)
+ eqs = [Assignment(outputField(i), s) for i, s in enumerate(symbols)]
+ return collisionRule.copy(collisionRule.mainAssignments + eqs)