Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Multiphase stabilization and mass-momentum consistency updates #1265

Open
wants to merge 7 commits into
base: master
Choose a base branch
from

Conversation

wjhorne
Copy link
Contributor

@wjhorne wjhorne commented Jun 10, 2024

This PR contains all of the multiphase work that was used to successfully run hybrid cases of waves passing over solid bodies. The main changes include

  • mass-momentum consistency via modifying the mass flux when vof is used. This was found to be a substantial improvement for moving interfaces over long runs.
  • upwinding stabilization with switching in the momentum equation at interfaces. This was a critical addition to ensure stability for high CFLs and harsh flow events that are encountered when waves pass over bodies.

@wjhorne
Copy link
Contributor Author

wjhorne commented Jun 20, 2024

@psakievich Do you know anyone that could take a look at this to start getting the ball rolling on getting it in?

@psakievich
Copy link
Contributor

I will try to take a look but I'm not sure how timely I can be. @BumseokLee would you also be willing to take a look? @mbkuhn your input would also be helpful. Also @rcknaus if you have some spare time your input is always appreciated. We have a P/T you can charge for your time if needed.

@psakievich psakievich self-requested a review June 25, 2024 20:30
@BumseokLee
Copy link

Hi Phil,
I am on a business trip until July 6. Also, I don’t have experiences with the multiphase stabilization or mass-momentum consistency in Nalu-Wind. I think it would be better to assign someone else to the work.
Bumseok

@psakievich
Copy link
Contributor

@wjhorne will you add some regression tests that utilize the user functions you've added?

Copy link
Contributor

@marchdf marchdf left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Great work! A bunch of minor comments.

include/edge_kernels/VOFAdvectionEdgeAlg.h Show resolved Hide resolved
virtual ~DropletVelocityAuxFunction() {}

using AuxFunction::do_evaluate;
virtual void do_evaluate(
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Kind of unfortunate we call this with 2 verbs ;) Probably not relevant to this PR.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It is strange for sure. There are a few places where patterns like this show up. There is also stuff like EquationSystem.C vs EquationSystems.C which can be confusing when working in the code. Likely worth a pass through the code to fix them down the road.

@@ -165,6 +165,9 @@
#include <user_functions/KovasznayVelocityAuxFunction.h>
#include <user_functions/KovasznayPressureAuxFunction.h>

#include <user_functions/DropletVelocityAuxFunction.h>
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Any reason for the ordering here? Want to group the multiphase aux functions together?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This got re-arranged during a conflict. I'll clean it up.

@@ -270,12 +281,12 @@ VolumeOfFluidEquationSystem::register_inflow_bc(

stk::mesh::MetaData& meta_data = realm_.meta_data();

// register boundary data; gamma_bc
// register boundary data; vof_bc
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Should we just remove these comments (and the following)? The code is pretty explicit I think.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I purged a bunch of these comments. We should probably do a pass through the code base to get more of them since this is quite common throughout the code.

ScalarFieldType* density_ = realm_.meta_data().get_field<double>(
stk::topology::NODE_RANK, "density");
std::vector<double> userSpec(1);
std::vector<double> userSpec(3);
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

we could initialize this with an initializer list instead of like this.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I agree and changed this over

// widths that are ~2 cells thick
DblType density_upwinding_factor = 1.0;
DblType alphaUpw_w_vof = alphaUpw;
DblType om_alphaUpw_w_vof = 1.0 - alphaUpw_w_vof;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We get so many questions about upwinding in nalu-wind and all the settings. It feels like we should document these extra modifications?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Along those lines. I would be in favor of adding a warning or throwing an error if the upwinding is put in a state that you have discovered is unstable. I'm assuming pure central will not work with VOF? Or is you modification below sufficient to ensure that the local upwinding is sufficient? Not knowing is also an acceptable answer, but I would like a clear guide line for the best practice, preferably in the code that is checked at runtime.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The documentation should definitely reflect the changes here. Our general finding is that pure central at the multiphase interface does not work well for VOF when you have high density ratios and strong flow dynamics.

The modifications here respect the users settings except explicitly at the multiphase interface when VOF is turned on. I think, at best, an appropriate warning would be to inform the user that by enabling VOF you have implicitly turned on upwinding at the interface.

? (alphaUpw * uIpL[i] + om_alphaUpw * uiIp)
: (alphaUpw * uIpR[i] + om_alphaUpw * uiIp);
const DblType uiUpw =
(mdot > 0.0) ? (alphaUpw_w_vof * uIpL[i] + om_alphaUpw_w_vof * uiIp)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

the suffix _vof is confusing to me. Can we keep it as it is and just modify the value if vof is on? My thinking is that if I read this code I am thinking that it is vof specific.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yeah I agree. It seems like the vof changes are to modifying the general upwinding? If there is not a second set of behavior that requires these variables for vof then I think it makes sense to just remove these variables and modify the original as needed when VOF is active.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It ends up being identical to alphaUpw and om_alphaUpw if a non-VOF case is run, or you are not at the interface with VOF. I'll change over the naming.

const double z = coords[2];
const double interface_thickness = 0.015;

fieldPtr[0] = 0.0;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

probably don't need this one, Just make the line after than a plain =?

fieldPtr[0] += -0.5 * (std::erf(y / interface_thickness) + 1.0) + 1.0;

// air-water
fieldPtr[0] = 1000.0 * fieldPtr[0] + (1.0 - fieldPtr[0]);
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

this is confusing? Maybe make the first assignement a variable and then use it here?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I will change this over. There was a lot of fast iteration on this function since it was used for a lot of our method development resulting in the cruft you see. I'll scan through some of the other user functions to make sure they don't have similar issues.

const double z0 = coords[2];
const double z1 = domain_height_;

// These need to come from elsewhere
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

yeah it should. Before we get into a mess like amr-wind and density.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Constants we should stash during construction. There is more freedom on the API there and we can glob them into a type specific struct.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This user function was added during some of the testing of the pressure behavior of some fully AMR-wind/Nalu-Wind coupled cases. Specifying the pressure is not necessary for this case. I will remove it completely

@wjhorne
Copy link
Contributor Author

wjhorne commented Jul 18, 2024

@wjhorne will you add some regression tests that utilize the user functions you've added?

We have fully-coupled tests that utilize the user functions. The only exception is the droplet case, which is also covered by VOFDroplet in our regression tests.

@wjhorne
Copy link
Contributor Author

wjhorne commented Jul 18, 2024

I pushed in changes addressing all of the review here. Let me know how it looks and I'll keep iterating

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

4 participants