diff --git a/lbmpy/creationfunctions.py b/lbmpy/creationfunctions.py
index 936ec0dd4ddd88c55655385b041552683a0c24e8..3e7474dc3732d1ab2341d6a19a6108bc2b6406e4 100644
--- a/lbmpy/creationfunctions.py
+++ b/lbmpy/creationfunctions.py
@@ -75,9 +75,8 @@ LES methods:
Fluctuating LB:
-- ``fluctuating=(variance1, variance2, )``: enables fluctuating lattice Boltzmann by randomizing collision process.
- Pass sequence of variances for each moment, or `True` to use symbolic variances
-
+- ``fluctuating``: enables fluctuating lattice Boltzmann by randomizing collision process.
+ Pass dictionary with parameters to ``lbmpy.fluctuatinglb.add_fluctuations_to_collision_rule``
Optimization Parameters
@@ -175,11 +174,12 @@ import sympy as sp
import lbmpy.forcemodels as forcemodels
from lbmpy.fieldaccess import (
- AAEvenTimeStepAccessor, AAOddTimeStepAccessor, CollideOnlyInplaceAccessor, EsoTwistEvenTimeStepAccessor,
- EsoTwistOddTimeStepAccessor, PdfFieldAccessor, PeriodicTwoFieldsAccessor, StreamPullTwoFieldsAccessor,
- StreamPushTwoFieldsAccessor)
-from lbmpy.fluctuatinglb import fluctuation_correction, fluctuating_variance_equations
-from lbmpy.methods import create_mrt3, create_mrt_orthogonal, create_mrt_raw, create_srt, create_trt, create_trt_kbc
+ AAEvenTimeStepAccessor, AAOddTimeStepAccessor, CollideOnlyInplaceAccessor,
+ EsoTwistEvenTimeStepAccessor, EsoTwistOddTimeStepAccessor, PdfFieldAccessor,
+ PeriodicTwoFieldsAccessor, StreamPullTwoFieldsAccessor, StreamPushTwoFieldsAccessor)
+from lbmpy.fluctuatinglb import add_fluctuations_to_collision_rule, method_with_rescaled_equilibrium_values
+from lbmpy.methods import (
+ create_mrt3, create_mrt_orthogonal, create_mrt_raw, create_srt, create_trt, create_trt_kbc)
from lbmpy.methods.creationfunctions import create_generic_mrt
from lbmpy.methods.cumulantbased import CumulantBasedLbMethod
from lbmpy.methods.entropic import add_entropy_condition, add_iterative_entropy_condition
@@ -193,8 +193,6 @@ from pystencils import Assignment, AssignmentCollection, create_kernel
from pystencils.cache import disk_cache_no_fallback
from pystencils.data_types import collate_types
from pystencils.field import Field, get_layout_of_array
-from pystencils.rng import random_symbol
-from pystencils.simp.assignment_collection import SymbolGen
from pystencils.stencil import have_same_entries
@@ -309,6 +307,9 @@ def create_lb_collision_rule(lb_method=None, optimization={}, **kwargs):
if rho_in is not None and isinstance(rho_in, Field):
rho_in = rho_in.center
+ if params['fluctuating']:
+ lb_method = method_with_rescaled_equilibrium_values(lb_method)
+
if u_in is not None:
density_rhs = sum(lb_method.pre_collision_pdf_symbols) if rho_in is None else rho_in
eqs = [Assignment(cqc.zeroth_order_moment_symbol, density_rhs)]
@@ -323,25 +324,7 @@ def create_lb_collision_rule(lb_method=None, optimization={}, **kwargs):
collision_rule = simplification(collision_rule)
if params['fluctuating']:
- variances = params["fluctuating"]
- if params["fluctuating"] is True:
- variances = [v for v, _ in zip(iter(SymbolGen("variance")), lb_method.moments)]
-
- correction = fluctuation_correction(lb_method, random_symbol(collision_rule.subexpressions, dim=lb_method.dim),
- variances)
-
- for i, corr in enumerate(correction):
- collision_rule.main_assignments[i] = Assignment(collision_rule.main_assignments[i].lhs,
- collision_rule.main_assignments[i].rhs + corr)
-
- if params["temperature"] is not None:
- temperature = sp.Symbol("temperature") if params["temperature"] is True else params["temperature"]
- variance_equations = fluctuating_variance_equations(method=lb_method,
- temperature=temperature,
- c_s_sq=params["c_s_sq"],
- variances=variances)
- collision_rule.subexpressions += variance_equations
- collision_rule.topological_sort(sort_subexpressions=True, sort_main_assignments=False)
+ add_fluctuations_to_collision_rule(collision_rule, **params['fluctuating'])
if params['entropic']:
if params['smagorinsky']:
diff --git a/lbmpy/fieldaccess.py b/lbmpy/fieldaccess.py
index a82384968ec3a7802a7db0f607598211c8b736e6..571db5845ac328c6bdb10c185875baecab7f436f 100644
--- a/lbmpy/fieldaccess.py
+++ b/lbmpy/fieldaccess.py
@@ -9,6 +9,7 @@ from pystencils.stencil import inverse_direction
__all__ = ['PdfFieldAccessor', 'CollideOnlyInplaceAccessor', 'StreamPullTwoFieldsAccessor',
'AAEvenTimeStepAccessor', 'AAOddTimeStepAccessor',
+ 'PeriodicTwoFieldsAccessor', 'StreamPushTwoFieldsAccessor',
'EsoTwistEvenTimeStepAccessor', 'EsoTwistOddTimeStepAccessor',
'visualize_pdf_field_accessor', 'visualize_field_mapping']
diff --git a/lbmpy/fluctuatinglb.py b/lbmpy/fluctuatinglb.py
index 076e063948c446657c59c7fe2d19005aecb3581b..caa7009983984bcce70e5f67671e99a31060ca99 100644
--- a/lbmpy/fluctuatinglb.py
+++ b/lbmpy/fluctuatinglb.py
@@ -2,41 +2,45 @@
to generate a fluctuating LBM the equilibrium moment values have to be scaled and an additive (random)
correction term is added to the collision rule
-
-
-Usage example:
-
- >>> from lbmpy.session import *
- >>> from pystencils.rng import random_symbol
- >>> method = create_lb_method(stencil='D2Q9', method='MRT')
- >>> rescaled_method = method_with_rescaled_equilibrium_values(method)
- >>> cr = create_lb_collision_rule(lb_method=rescaled_method)
- >>> correction = fluctuation_correction(rescaled_method,
- ... rng_generator=random_symbol(cr.subexpressions, dim=method.dim))
- >>> for i, corr in enumerate(correction):
- ... cr.main_assignments[i] = ps.Assignment(cr.main_assignments[i].lhs,
- ... cr.main_assignments[i].rhs + corr)
- >>>
-
"""
+import numpy as np
import sympy as sp
from lbmpy.moments import MOMENT_SYMBOLS
+from pystencils import Assignment, TypedSymbol
+from pystencils.rng import PhiloxFourFloats, random_symbol
from pystencils.simp.assignment_collection import SymbolGen
-from pystencils import Assignment
-def fluctuating_variance_equations(method, temperature=sp.Symbol("fluctuating_temperature"),
- c_s_sq=sp.Symbol("c_s") ** 2,
- variances=SymbolGen("variance")):
+def add_fluctuations_to_collision_rule(collision_rule, temperature=None, variances=(),
+ block_offsets=(0, 0, 0), seed=TypedSymbol("seed", np.uint32),
+ rng_node=PhiloxFourFloats, c_s_sq=sp.Rational(1, 3)):
+ """"""
+ if not (temperature and not variances) or (temperature and variances):
+ raise ValueError("Fluctuating LBM: Pass either 'temperature' or 'variances'.")
+
+ method = collision_rule.method
+ if not variances:
+ variances = fluctuating_variance_from_temperature(method, temperature, c_s_sq)
+
+ rng_symbol_gen = random_symbol(collision_rule.subexpressions, seed,
+ rng_node=rng_node, dim=method.dim, offsets=block_offsets)
+ correction = fluctuation_correction(method, rng_symbol_gen, variances)
+
+ for i, corr in enumerate(correction):
+ collision_rule.main_assignments[i] = Assignment(collision_rule.main_assignments[i].lhs,
+ collision_rule.main_assignments[i].rhs + corr)
+
+
+def fluctuating_variance_from_temperature(method, temperature, c_s_sq=sp.Symbol("c_s") ** 2):
"""Produces variance equations according to (3.54) in Schiller08"""
normalization_factors = abs(method.moment_matrix) * sp.Matrix(method.weights)
- relaxation_rates_sqr = [r ** 2 for r in method.relaxation_rates]
density = method.zeroth_order_equilibrium_moment_symbol
+ if method.conserved_quantity_computation.zero_centered_pdfs:
+ density += 1
mu = temperature * density / c_s_sq
-
- return [Assignment(v, sp.sqrt(mu * norm * (1 - rr)))
- for v, norm, rr in zip(iter(variances), normalization_factors, relaxation_rates_sqr)]
+ return [sp.sqrt(mu * norm * (1 - (1 - rr) ** 2))
+ for norm, rr in zip(normalization_factors, method.relaxation_rates)]
def method_with_rescaled_equilibrium_values(base_method):
diff --git a/lbmpy/methods/conservedquantitycomputation.py b/lbmpy/methods/conservedquantitycomputation.py
index 180a2f6b170204322ec2984887d0485e99db9a97..d34939b7e89477ebe2b2f4f531978a9bc2acb184 100644
--- a/lbmpy/methods/conservedquantitycomputation.py
+++ b/lbmpy/methods/conservedquantitycomputation.py
@@ -243,10 +243,11 @@ def get_equations_for_zeroth_and_first_order_moment(stencil, symbolic_pdfs,
\rho = \sum_{d \in S} f_d
u_j = \sum_{d \in S} f_d u_jd
- :param stencil: called :math:`S` above
- :param symbolic_pdfs: called :math:`f` above
- :param symbolic_zeroth_moment: called :math:`\rho` above
- :param symbolic_first_moments: called :math:`u` above
+ Args:
+ stencil: called :math:`S` above
+ symbolic_pdfs: called :math:`f` above
+ symbolic_zeroth_moment: called :math:`\rho` above
+ symbolic_first_moments: called :math:`u` above
"""
def filter_out_plus_terms(expr):
result = 0
diff --git a/lbmpy_tests/test_fluctuation.py b/lbmpy_tests/test_fluctuation.py
index 88b01cad2b3b6b29237a39ab46748490d38b37b3..be886c280189da7fb460891501cf88322eed20c6 100644
--- a/lbmpy_tests/test_fluctuation.py
+++ b/lbmpy_tests/test_fluctuation.py
@@ -1,9 +1,13 @@
+import numpy as np
+
from lbmpy.scenarios import create_channel
+
def test_fluctuating_generation_pipeline():
ch = create_channel((40, 10), method='mrt3', relaxation_rates=[1.5, 1, 1], force=1e-5,
- fluctuating=True,
- temperature=1e-9,
- kernel_params={'time_step': 1})
+ fluctuating={'temperature': 1e-9},
+ kernel_params={'time_step': 1, 'seed': 312},
+ optimization={'cse_global': True})
ch.run(10)
+ assert np.max(ch.velocity[:, :]) < 0.1