Skip to content

Commit

Permalink
restore examples to version found at production
Browse files Browse the repository at this point in the history
  • Loading branch information
a-alveyblanc committed Feb 24, 2025
1 parent ae29c3e commit ced5d6a
Show file tree
Hide file tree
Showing 3 changed files with 395 additions and 128 deletions.
5 changes: 2 additions & 3 deletions examples/euler/acoustic_pulse.py
Original file line number Diff line number Diff line change
Expand Up @@ -146,7 +146,7 @@ def run_acoustic_pulse(actx,

from meshmode.mesh.generation import generate_regular_rect_mesh

dim = 3
dim = 2
box_ll = -0.5
box_ur = 0.5
group_cls = TensorProductElementGroup if tpe else None
Expand Down Expand Up @@ -257,8 +257,7 @@ def main(ctx_factory, order=3, final_time=1, resolution=16,
queue = cl.CommandQueue(cl_ctx)

if lazy:
from grudge.array_context import FusionContractorArrayContext
actx = FusionContractorArrayContext(
actx = PytatoPyOpenCLArrayContext(
queue,
allocator=cl_tools.MemoryPool(cl_tools.ImmediateAllocator(queue)),
)
Expand Down
100 changes: 32 additions & 68 deletions examples/euler/vortex.py
Original file line number Diff line number Diff line change
Expand Up @@ -32,16 +32,11 @@
from grudge.models.euler import (
vortex_initial_condition,
EulerOperator,
# EntropyStableEulerOperator
EntropyStableEulerOperator
)
import grudge.op as op
from grudge.array_context import (
FusionContractorArrayContext,
PyOpenCLArrayContext,
PytatoPyOpenCLArrayContext
)
from grudge.array_context import PyOpenCLArrayContext, PytatoPyOpenCLArrayContext
from grudge.shortcuts import rk4_step
from meshmode.mesh import SimplexElementGroup


logger = logging.getLogger(__name__)
Expand All @@ -51,8 +46,7 @@ def run_vortex(actx, order=3, resolution=8, final_time=5,
esdg=False,
overintegration=False,
flux_type="central",
visualize=False,
group_cls=SimplexElementGroup):
visualize=False):

logger.info(
"""
Expand All @@ -63,11 +57,10 @@ def run_vortex(actx, order=3, resolution=8, final_time=5,
entropy stable: %s\n
overintegration: %s\n
flux_type: %s\n
visualize: %s\n
element type: %s
visualize: %s
""",
order, final_time, resolution, esdg,
overintegration, flux_type, visualize, group_cls
overintegration, flux_type, visualize
)

# eos-related parameters
Expand All @@ -81,20 +74,19 @@ def run_vortex(actx, order=3, resolution=8, final_time=5,
a=(0, -5),
b=(20, 5),
nelements_per_axis=(2*resolution, resolution),
periodic=(True, True),
group_cls=group_cls)
periodic=(True, True))

from meshmode.discretization.poly_element import (
InterpolatoryEdgeClusteredGroupFactory,
QuadratureGroupFactory
QuadratureSimplexGroupFactory,
default_simplex_group_factory,
)

from grudge.discretization import make_discretization_collection
from grudge.dof_desc import DISCR_TAG_BASE, DISCR_TAG_QUAD

if esdg:
case = "esdg-vortex"
operator_cls = EulerOperator
operator_cls = EntropyStableEulerOperator
else:
case = "vortex"
operator_cls = EulerOperator
Expand All @@ -110,8 +102,9 @@ def run_vortex(actx, order=3, resolution=8, final_time=5,
dcoll = make_discretization_collection(
actx, mesh,
discr_tag_to_group_factory={
DISCR_TAG_BASE: InterpolatoryEdgeClusteredGroupFactory(order),
DISCR_TAG_QUAD: QuadratureGroupFactory(2*order)
DISCR_TAG_BASE: default_simplex_group_factory(
base_dim=mesh.dim, order=order),
DISCR_TAG_QUAD: QuadratureSimplexGroupFactory(2*order)
}
)

Expand All @@ -129,6 +122,8 @@ def run_vortex(actx, order=3, resolution=8, final_time=5,
def rhs(t, q):
return euler_operator.operator(actx, t, q)

compiled_rhs = actx.compile(rhs)

fields = vortex_initial_condition(actx.thaw(dcoll.nodes()))

from grudge.dt_utils import h_min_from_volume
Expand All @@ -149,42 +144,24 @@ def rhs(t, q):

step = 0
t = 0.0

compiled_rhs = actx.compile(rhs)

import time
start = time.time()
fields = rk4_step(fields, t, dt, compiled_rhs)
elapsed = time.time()

norm_q = actx.to_numpy(op.norm(dcoll, fields, 2))
logger.info("[%04d] t = %.5f |q| = %.5e, walltime (s) = %.5f",
step, t, norm_q, abs(start - elapsed))

t += dt
step += 1

while t < final_time:
if visualize:
vis.write_vtk_file(
f"{exp_name}-{step:04d}.vtu",
[
("rho", fields.mass),
("energy", fields.energy),
("momentum", fields.momentum)
]
)
assert norm_q < 200
if step % 10 == 0:
norm_q = actx.to_numpy(op.norm(dcoll, fields, 2))
logger.info("[%04d] t = %.5f |q| = %.5e", step, t, norm_q)

if visualize:
vis.write_vtk_file(
f"{exp_name}-{step:04d}.vtu",
[
("rho", fields.mass),
("energy", fields.energy),
("momentum", fields.momentum)
]
)
assert norm_q < 200

fields = actx.thaw(actx.freeze(fields))
start = time.time()
fields = rk4_step(fields, t, dt, compiled_rhs)
elapsed = time.time()

norm_q = actx.to_numpy(op.norm(dcoll, fields, 2))
logger.info("[%04d] t = %.5f |q| = %.5e, walltime (s) = %.5f",
step, t, norm_q, abs(start - elapsed))

t += dt
step += 1

Expand All @@ -196,23 +173,14 @@ def main(ctx_factory, order=3, final_time=5, resolution=8,
overintegration=False,
lf_stabilization=False,
visualize=False,
lazy=False,
use_tpe=False):
lazy=False):
cl_ctx = ctx_factory()
queue = cl.CommandQueue(cl_ctx)

if use_tpe:
from meshmode.mesh import TensorProductElementGroup
group_cls = TensorProductElementGroup
else:
from meshmode.mesh import SimplexElementGroup
group_cls = SimplexElementGroup

if lazy:
actx = FusionContractorArrayContext(
actx = PytatoPyOpenCLArrayContext(
queue,
allocator=cl_tools.MemoryPool(cl_tools.ImmediateAllocator(queue)),
use_tp_transforms=use_tpe
)
else:
actx = PyOpenCLArrayContext(
Expand Down Expand Up @@ -240,8 +208,7 @@ def main(ctx_factory, order=3, final_time=5, resolution=8,
esdg=esdg,
overintegration=overintegration,
flux_type=flux_type,
visualize=visualize,
group_cls=group_cls
visualize=visualize
)


Expand All @@ -262,8 +229,6 @@ def main(ctx_factory, order=3, final_time=5, resolution=8,
help="write out vtk output")
parser.add_argument("--lazy", action="store_true",
help="switch to a lazy computation mode")
parser.add_argument("--tpe", action="store_true",
help="use a tensor-product discretization")
args = parser.parse_args()

logging.basicConfig(level=logging.INFO)
Expand All @@ -275,5 +240,4 @@ def main(ctx_factory, order=3, final_time=5, resolution=8,
overintegration=args.oi,
lf_stabilization=args.lf,
visualize=args.visualize,
lazy=args.lazy,
use_tpe=args.tpe)
lazy=args.lazy)
Loading

0 comments on commit ced5d6a

Please sign in to comment.