test_entropic.py 1.71 KB
Newer Older
1
2
import numpy as np
import sympy as sp
Martin Bauer's avatar
Martin Bauer committed
3

4
from lbmpy.forcemodels import Guo
Frederik Hennig's avatar
Frederik Hennig committed
5
from lbmpy.methods.momentbased.entropic_eq_srt import create_srt_entropic
6
7
8
9
10
11
12
from lbmpy.scenarios import create_lid_driven_cavity
from lbmpy.stencils import get_stencil


def test_entropic_methods():
    sc_kbc = create_lid_driven_cavity((20, 20), method='trt-kbc-n4', relaxation_rate=1.9999,
                                      entropic_newton_iterations=3, entropic=True, compressible=True,
13
                                      force=(-1e-10, 0), force_model="luo")
14
15

    sc_srt = create_lid_driven_cavity((40, 40), relaxation_rate=1.9999, lid_velocity=0.05, compressible=True,
16
                                      force=(-1e-10, 0), force_model="luo")
17
18

    sc_entropic = create_lid_driven_cavity((40, 40), method='entropic-srt', relaxation_rate=1.9999,
19
                                           lid_velocity=0.05, compressible=True, force=(-1e-10, 0), force_model="luo")
20
21
22
23
24
25
26
27
28
29
30

    sc_srt.run(1000)
    sc_kbc.run(1000)
    sc_entropic.run(1000)
    assert np.isnan(np.max(sc_srt.velocity[:, :]))
    assert np.isfinite(np.max(sc_kbc.velocity[:, :]))
    assert np.isfinite(np.max(sc_entropic.velocity[:, :]))


def test_entropic_srt():
    stencil = get_stencil("D2Q9")
Markus Holzer's avatar
Markus Holzer committed
31
32
    relaxation_rate = 1.8
    method = create_srt_entropic(stencil, relaxation_rate, Guo((0, 1e-6)), True)
33
34
    assert method.zeroth_order_equilibrium_moment_symbol == sp.symbols("rho")
    assert method.first_order_equilibrium_moment_symbols == sp.symbols("u_:2")
Markus Holzer's avatar
Markus Holzer committed
35
36
37
38
39
40
41
42
43

    eq = method.get_equilibrium()
    terms = method.get_equilibrium_terms()
    rel = method.relaxation_rates

    for i in range(len(terms)):
        assert sp.simplify(eq.main_assignments[i].rhs - terms[i]) == 0
        assert rel[i] == 1.8