search

AMReX-Combustion / PeleLMeX

An adaptive mesh hydrodynamics simulation code for low Mach number reacting flows without level sub-cycling.
https://amrex-combustion.github.io/PeleLMeX/
BSD 3-Clause "New" or "Revised" License
33 stars 41 forks source link

Soret isothermal wall fix #391

Closed ThomasHowarth closed 5 days ago

ThomasHowarth commented 5 months ago

This PR addresses #389 . Happy to discuss any of the implementation, and let me know if you can see spots for improvement/clarity. Currently, I've not implemented the necessary changes required for EB.

drummerdoc commented 5 months ago

Is it possible to check in the test case?

ThomasHowarth commented 3 months ago

Sorry for the delay, busy month or so...

I've added in the changes so it is computed when computing divu, I've added in the same boundary input to the diffuse_scalar call when running with EB. Code also aborts on Isothermal+Soret+Efield. I am also still seeing some species flux in the species balance files, on the order of 3e-7 or so - any thoughts on this?

baperry2 commented 2 months ago

Sorry @ThomasHowarth, I know you've got a few big PRs in the queue waiting for review. I'm going to try to get to them this week.

I've dug into this one a bit and think I figured out the where the remaining conservation error is coming from. With the lagged correction approach, I was only thinking of the implicitly solved $\tilde{D}^{n+1,k+1}$ term in the update equation below, but, as you've now implemented the corrections must also be applied when computing $D^n$ and ${D}^{n+1,k}$, but those are explicitly computed and can't be lagged. Basically, the Soret and Wbar fluxes should be computed first in computeDifferentialDiffusionFluxes and those values should be used to set the inhomogeneous neumann conditions for the diffusive fluxes. For the Soret fluxes that should be easy with some code rearrangement. For the Wbar fluxes it's a little tricky because the Wbar fluxes at the isothermal boundary depend on the species gradients at the boundary, but also affect what the boundary condition needs to be for the species gradients, leading to an implicit definition of the boundary condition. I think it should be okay to compute the Soret terms first, use the boundary species gradients from the Soret terms in computing Wbar fluxes, then use the combined boundary gradients from both Soret and Wbar for the diffusive fluxes. Iterating around the Wbar flux computation would be the fully correct thing though.

Screenshot 2024-09-04 at 10 14 21 AM

Let me know if this makes sense. This is definitely a tricky problem and I'm not sure if I fully understand it. But I did hack a test for the Soret part of the fix (while turning Wbar fluxes off) and that did lead to correct mass conservation.

BTW, I also submitted #410 to make using temporals to track conservation a bit more user friendly.

ThomasHowarth commented 2 months ago

Without the Wbar correction, I have also managed to eliminate the conservation error (Up to 10^-23). However, iterating around the wbar fluxes seems to blow up the gradients, and I'm not entirely sure why. I'll keep working on it, but I've also pushed the latest version so you can take a look if you like.

baperry2 commented 2 months ago

Nice progress a @ThomasHowarth. This has turned out to be much more involved than expected!

If a solution doesn't jump right out at you for the exploding wbar iterations, do you think we can merge this without the wbar iterations (or maybe hardcoding to 1 "iteration") and then have a separate PR for the final correction? There's a lot here that I don't want to risk falling out of date and even without the last bit it's more correct than what's in there now.

ThomasHowarth commented 2 months ago

I think perhaps the best course is to just do the Soret correction for the conservation, disable the use_wbar for consistency (in these cases) and leave a warning/todo note. If this seems reasonable I can tidy things up and leave it like this for now? As you say, at least it will work, if not perfect yet.

ThomasHowarth commented 1 month ago

So, I've done some more thinking about why things blow up. Since we have to do some portion of the solve explicitly rather than lagged, I think we revert back to the point I made in the discussion, and the Y coefficient of the Wbar term becomes a problem. As I say, I can either disable the Wbar corrections in this case and everything works, otherwise this is going to become even more complicated...

ThomasHowarth commented 1 week ago

@baperry2 I have resorted to disabling wbar for this particular case - the boundary condition becomes very complicated otherwise. Some of the machinery is still left (e.g. filling boundaries for wbar and computing gradients of wbar using gradients of Y), but it should be bypassed. I'm leaving it there in case someone wants to try and figure this out down the line - I can also remove it. This way the solver is both conservative and consistent, even if the diffusion model isn't perfect.

Given the number of commits, and switch arounds, it is probably best to squash this merge.

baperry2 commented 5 days ago

@jrood-nrel any idea why the Docs workflow is not running? It seems to be not running since you updated the workflow ubuntu version.

This doesn't touch Docs, so i'm going to merge, but we should figure that out since it appears to be affecting all PRs

jrood-nrel commented 5 days ago

It requires certain files to have changed. I assume that's it. https://github.com/AMReX-Combustion/PeleLMeX/blob/b3d27d4d40943772bdd47a9be8335d91ff9f73a5/.github/workflows/docs.yml#L7-L9 . I would also be fine having it always run.

baperry2 commented 5 days ago

I think this is fine, I'll just remove it from the list of "required" checks then, otherwise it will show up as "expected" and prevent auto-merging for cases where it does not run.