Av operator #582

MTCam · 2021-12-16T15:13:07Z

This change set adds artificial viscosity (mostly) as developed by @w-hagen, with the tweaks and improvements we've made along the way including @thomasgibson changes to bring in over-integration capabilities.

Checklist for @lukeolson

inducer

Some initial thoughts below.

inducer · 2022-06-01T21:17:51Z

test/test_av.py

+    Tests that the cells tagging properly tags cells
+    given prescribed solutions.


Expand to describe precisely what is being tested.

Took a stab at this. Want to have a look and sharpen it up @w-hagen ?

inducer · 2022-06-05T22:33:28Z

test/test_av.py

@@ -0,0 +1,265 @@
+"""Test the artificial viscosity functions."""


IMO, this is missing a test for proper handing of boundary conditions, particularly in view of the fact that we had trouble with that.

A test has been added to make sure the AV fluxes returned from the fluid BC constructs match those which we intended for them to be. We collectively share a lack of knowledge about what the BCs should actually be for this operator (i.e., when computing the 2nd part; div(alpha*indicator*grad(Q))). For the first (grad(Q)) part, this operator now mercifully uses grad(Q) from the fluid computation, and is thus using the appropriate fluid BCs directly there - whatever those may be for the given boundary. The added tests check the AV-specific interface fluxes for the 2nd part.

inducer · 2022-06-05T22:33:55Z

test/test_av.py

+    """Test artificial_viscosity.
+
+    Tests the application on a few simple functions
+    to confirm artificial viscosity returns the analytical result.
+    """


Improve the description here.

Took a stab at this in a9743ee.

inducer · 2022-06-05T22:34:29Z

test/test_av.py

+    """Test artificial_viscosity.
+
+    Tests the application on a few simple functions
+    to confirm artificial viscosity returns the analytical result.
+    """


I think this should be replaced by a manufactured-solution test, instead of these two simplistic known-solution tests.

I tried quickly re-use the diffusion operator MMS tests, but turns out not to be as easy of a drop-in as I had imagined.

We'll need some infrastructure backing to make this happen. It looks like the amount of effort to get an MMS infrastructure for AV is about the same as replacing the guts of AV with the diffusion operator. I think we'd prefer using the diffusion operator. Another option would be just beefing up the tests a bit with more analytic functions, perhaps leveraging the symbolic infrastructure.

@lukeolson added a trig test that helps with this a little.

I tried to add a discontinuous coeff test (which I do not expect this operator to pass) - but it turned out to be too much effort for this fluid-specific operator. If this operator was just a (multi)scalar laplacian operator, then it'd be a no-brainer. Alas, now it re-uses the gradient operator and boundary projection routines designed specifically for fluid solutions - and thus requires physically-consistent fluid solutions.

inducer · 2022-06-05T22:35:15Z

test/test_av.py

+    """Test artificial_viscosity.
+
+    Tests the application on a few simple functions
+    to confirm artificial viscosity returns the analytical result.
+    """


This does not test the numerical fluxes, since there are no jumps. As a result, numerical stability of AV is not being tested.

It appears that the diffusion operator has only one test where jumps are tested. We could potentially port that over to AV testing, if that is of interest.

mirgecom/initializers.py

inducer · 2022-06-05T22:44:44Z

mirgecom/artificial_viscosity.py

+    pass
+
+
+def av_laplacian_operator(discr, boundaries, fluid_state, alpha, gas_model=None,


Let me ask the elephant-in-the-room question here: Given that we have variable-coefficient diffusion, why do we want this? As @majosm can possibly attest, getting diffusion with elementwise-discontinuous coefficients right is not entirely trivial. And the already-implemented diffusion operator literally implements the same operator. It doesn't help that this operator is not well-tested.

I'd be OK with this becoming a deprecated wrapper of the diffusion operator.

inducer · 2022-06-05T22:49:27Z

mirgecom/artificial_viscosity.py

+    return -div_operator(discr, dd_vol, dd_faces, interp_to_vol_quad(r), r_bnd)
+
+
+def smoothness_indicator(discr, u, kappa=1.0, s0=-6.0):


I think this can be built simply out of the modal discretization facility in meshmode. There's definitely no need for this to make use of a loopy kernel. @kaushikcfd had also worked on a more lazy-compatible version of this. I am surprised to see that that is not part of this.

I believe the version you are thinking of is added by @kaushikcfd in #666. It requires the fusion array context and does not work without it afaik.

The tagging bit won't work, but the einsum bit should be fine. I don't have a preference for which PR the einsum is introduced in.

What part of the tagging is the issue? I'm looking over the code in #666, and it mostly seems like it ought to go even with arraycontext/meshmode main?

Tag types like DiscretizationDOFAxisTag haven't made it to meshmode main.

Do we have a path forward for this item?

fwiw, here's the error i get when trying to run with @kaushikcfd indicator:

indicator = smoothness_indicator(discr, fluid_state.mass_density, File "/Users/mtcampbe/CEESD/devel/tau/main/emirge/mirgecom/mirgecom/artificial_viscosity.py", line 318, in smoothness_indicator from meshmode.transform_metadata import DiscretizationDOFAxisTag ImportError: cannot import name 'DiscretizationDOFAxisTag' from 'meshmode.transform_metadata' (/Users/mtcampbe/CEESD/devel/tau/

The actx.tag_axis part would be illegal, but without that the call to actx.einsum should be legal on meshmode main. I'm ok with the rewrite to einsum either here or in #666.

I just added @kaushikcfd version of this in (732f005). The tests passed locally with main (:tada:).

inducer · 2022-06-05T22:50:53Z

It might make sense to squash this so that the Github UI does not bog down.

lukeolson · 2022-06-06T19:12:00Z

examples/doublemach-mpi.py

+            health_error = True
+            logger.info(f"{rank=}: NANs/Infs in pressure data.")
+
+        from mirgecom.simutil import allsync


I believe these should be global_reduce now

lukeolson · 2022-06-06T19:18:32Z

examples/doublemach-mpi.py

+        default_simplex_group_factory, QuadratureSimplexGroupFactory
+
+    order = 3
+    discr = EagerDGDiscretization(


These should be updated to the new collection. But this is done in #360

Let's make a dedicated issue about updating mirgecom to use discretization collections. This cannot be done example-by-example. It requires a large collection of supporting updates throughout mirgecom.

Working to centralize and then replace the discretization collection in a dedicated effort was precisely the reason we added this code:
https://github.com/illinois-ceesd/mirgecom/blob/main/mirgecom/discretization.py

7b9ee93 uses the wrapper, and gets rid of the deprecation warning.

lukeolson · 2022-06-06T19:24:01Z

examples/doublemach-mpi.py

+from mirgecom.logging_quantities import (
+    initialize_logmgr,
+    logmgr_add_many_discretization_quantities,
+    logmgr_add_device_name,


use logmgr_cl_device_name

lukeolson · 2022-06-06T19:27:06Z

examples/doublemach-mpi.py

+    current_state = make_fluid_state(cv=current_cv, gas_model=gas_model)
+
+    visualizer = make_visualizer(discr,
+                                 discr.order if discr.dim == 2 else discr.order)


.order is deprecated (it appears)

Also, what is the point of this line.

I think that was used for a bit to plot at a higher discretization order but looks like I changed that and never got rid of the if statement. That part can be removed.

Oh sure - but this particular line says something like "if sthing then discr.order else discr.order"

Co-authored-by: Andreas Klöckner <[email protected]>

…ard-dcoll

lukeolson · 2022-06-10T22:36:42Z

mirgecom/artificial_viscosity.py

+        0.5 * (1.0 + actx.np.sin(np.pi * (indicator - s0) / (2.0 * kappa))),
+        0.0 * indicator,
+    )
+    indicator = actx.np.where(yesnou, 1.0 + 0.0 * indicator, sin_indicator)


np.where() works with a scalar field... this would avoid an extra dof array by just using 1.0

I'm not sure if that's extended in the array context though worth checking

I do not know if it has changed (i doubt it), but I definitely needed each "rhs" args to be arrays of equal size when using this function in the past.

heh, whadda ya know.

mirgecom/test/test_various.py

Line 34 in 532549a

def test_actx_np_where(actx_factory):

Went ahead and fixed this in (732f005)

MTCam added Production Feeder Feeds the production branch enhancement New feature or request labels Dec 16, 2021

MTCam mentioned this pull request Feb 18, 2022

Production #551

Draft

inducer mentioned this pull request Apr 15, 2022

Lazy compilation fails for isolator driver #640

Open

MTCam requested review from lukeolson and w-hagen May 27, 2022 15:55

MTCam marked this pull request as ready for review May 27, 2022 15:55

MTCam changed the base branch from ns-update to main June 1, 2022 21:20

inducer reviewed Jun 5, 2022

View reviewed changes

MTCam mentioned this pull request Jun 6, 2022

Refactor AV to use the diffusion operator #684

Open

Add @w-hagen artificial viscosity operator.

b9bb3f4

MTCam force-pushed the av-operator branch from 71f43be to b9bb3f4 Compare June 6, 2022 14:25

lukeolson reviewed Jun 6, 2022

View reviewed changes

MTCam removed the request for review from w-hagen June 6, 2022 20:09

MTCam and others added 9 commits June 6, 2022 16:48

eliminate some deprecated constructs

7b9ee93

Sharpen the test comments and simplify just a bit.

a9743ee

Replace 'an exact' --> 'a'

2911563

Co-authored-by: Andreas Klöckner <[email protected]>

Move toward EagerDiscretization --> Discretization Collection

e741a24

Fix op --> oper instances.

ce8e9e0

Fix up weak_{grad, div} --> weak_local_{grad, div}

658183f

Merge branch 'main' into move-toward-dcoll

2f39e05

Use proper dofdescs

884ebcd

Merge remote-tracking branch 'origin/move-toward-dcoll' into move-tow…

fa1ad16

…ard-dcoll

lukeolson reviewed Jun 10, 2022

View reviewed changes

MTCam added 2 commits June 10, 2022 19:32

Use base domain, not quad domain

c7a5677

Merge branch 'main' into av-operator

cb753a3

MTCam and others added 13 commits June 12, 2022 09:24

Excise unneeded array value arg to actx.np.where

732f005

Customize production env.

2f9c0d3

Notation fix in div/grad

3721471

Sharpen and clean up some superfluous patterns.

1bf3fe3

add sin test for convergence

82f617c

Add trig test from @lukeo slightly twisted.

68b70d2

Sharpen docs, note limitations and AV replacement.

cb23f1f

debug latex

bfa1143

Add boundary handling test

c101096

Move toward dcoll

7047c72

remove superfluous 1.0

5f6e1e3

Undo production customization

c156836

Merge upstream

3fa86c4

lukeolson approved these changes Jun 15, 2022

View reviewed changes

Merge branch 'main' into av-operator

2d953c1

MTCam enabled auto-merge (squash) June 15, 2022 16:48

MTCam merged commit 117714c into main Jun 15, 2022

MTCam deleted the av-operator branch June 15, 2022 17:19

MTCam mentioned this pull request Jun 25, 2022

Artificial viscosity #213

Closed

6 tasks

		Tests that the cells tagging properly tags cells
		given prescribed solutions.

		@@ -0,0 +1,265 @@
		"""Test the artificial viscosity functions."""

		pass


		def av_laplacian_operator(discr, boundaries, fluid_state, alpha, gas_model=None,

		return -div_operator(discr, dd_vol, dd_faces, interp_to_vol_quad(r), r_bnd)


		def smoothness_indicator(discr, u, kappa=1.0, s0=-6.0):

Av operator #582

Av operator #582

Conversation

MTCam commented Dec 16, 2021 • edited by lukeolson Loading

inducer left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

MTCam Jun 6, 2022 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

inducer commented Jun 5, 2022

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

w-hagen Jun 6, 2022 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

lukeolson Jun 10, 2022 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

MTCam commented Dec 16, 2021 •

edited by lukeolson

Loading

MTCam Jun 6, 2022 •

edited

Loading

w-hagen Jun 6, 2022 •

edited

Loading

lukeolson Jun 10, 2022 •

edited

Loading