MFC icon indicating copy to clipboard operation
MFC copied to clipboard

All physical quantities become "infinity".

Open HIT-CFD-group opened this issue 3 years ago • 22 comments

Hi,

I tried to reproduce the case of S.Sembian(2016), but icfl condition in the run_time.inf became infinity . And all the physical quantities including alpha, pressure, velocity also became infinity.

The infinity firstly comes in several grids, and then grows into the whole grids. image image image image image

These are my input file. input.zip

What I am doing wrong? Could you advice, please?

HIT-CFD-group avatar Apr 22 '22 13:04 HIT-CFD-group

Interesting, thank you for pointing this out. We will look shortly! @anandrdbz @mrodrig6 would you mind having a look as well?

sbryngelson avatar Apr 22 '22 14:04 sbryngelson

@boombiubla did you consider a smaller time step?

sbryngelson avatar Apr 22 '22 14:04 sbryngelson

@sbryngelson Yes, and this is the run_time.inf of the case of time step equal to 1.5e-8. The time step of the file i sent before is 2.3e-8. run_time_dt=1.55e-8.txt

HIT-CFD-group avatar Apr 23 '22 06:04 HIT-CFD-group

@boombiubla I was playing around with your test case here, using the usual numerical tricks available in MFC. I didn't have much success. I think part of the problem could be the lack of normalization. For example, your pressures are very large compared to other quantities in the code (like volume fractions). This isn't always a problem, but could be causing one here. My suggestion is to non-dimensionalize.

sbryngelson avatar Apr 23 '22 13:04 sbryngelson

@sbryngelson Thanks for your suggestion. Does it mean that I should change the form of equation in the MFC code? And where is the specific location of corresponding code?

HIT-CFD-group avatar Apr 24 '22 08:04 HIT-CFD-group

@boombiubla Oh no, I just meant adjusting your input variables (initial conditions) inside of the input script. So you could normalize everything such that the max. initial pressure is 1, for example. The code is agnostic to units.

sbryngelson avatar Apr 24 '22 18:04 sbryngelson

@sbryngelson Yes, I agree with you that I should take a non-dimentional ratio to prevent the large value, but I am still confused about that you told me not to change the equation form.

When the variables like velocity and pressure are normalized, I think the equations should take the same form of normalization like below pictures about the change from form one to form two.🤔

image image

HIT-CFD-group avatar Apr 26 '22 05:04 HIT-CFD-group

I see what you mean. Yes, for the non-dimensional coefficients that show up, like Re, we do use those. In your case you have no body forces or viscosity (at least for now!), so you don't have to input any such coefficients. Instead, just normalize your state variables to the appropriate renormalization group.

sbryngelson avatar Apr 26 '22 12:04 sbryngelson

Thanks! I sent the assignment to do the calculation just now. And maybe I will study the case with body forces and viscosity soon. What should I do at that time?

HIT-CFD-group avatar Apr 27 '22 11:04 HIT-CFD-group

The viscosity is determined via the input Reynolds number. This should be in the user documentation. There is no support for body forces right now, though they're trivial to add. They were also added by another user and could be merged in if needed.

sbryngelson avatar Apr 27 '22 12:04 sbryngelson

Would you mind sending me some input files for learning? Maybe files about JC Meng's works are better.

HIT-CFD-group avatar May 01 '22 05:05 HIT-CFD-group

There are example input files available in the repository. I do not have Meng's input files on hand.

sbryngelson avatar May 01 '22 19:05 sbryngelson

I tried to non-dimensionalize the pressure like 2.py, but it was still not successful. At the beginning, I forgot to normalize the pi_inf, and it became infinity in step 1828. Then I fixed it, and it became infinity in step 1827. norm_pre.zip

HIT-CFD-group avatar May 02 '22 06:05 HIT-CFD-group

It looks like you didn't nondimensionalize velocity and density. You cannot nondimensionalize one variable by a reference quantity (p0 in your case) and then leave the rest in usual units.

sbryngelson avatar May 02 '22 12:05 sbryngelson

Hmmmm, I tried once again like below. But it seems not successful. Is there any wrong in my input file? norm.zip

HIT-CFD-group avatar May 05 '22 02:05 HIT-CFD-group

Hi @boombiubla -- We are still looking at your case. Your normalizations look fine. I suspect the problem is with the problem setup itself. Perhaps it's simply too strong a shock wave for the numerical method we use. Right now it looks like we are computing negative sound speeds, which are associated with the pressure becoming more negative than the liquid stiffness (pi_inf). Also, notice that the paper you reference is doing a 1D simulation, which should be much easier from the numerical point of view.

In either case, I hope to have a better answer for you soon.

sbryngelson avatar May 05 '22 17:05 sbryngelson

OK, thanks for your attention and suggestions!

HIT-CFD-group avatar May 06 '22 01:05 HIT-CFD-group

@boombiubla I thought I would circle back to this. I have not had a chance to do too much more testing on this problem. My suspicion is that the problem is either too challenging for the numerical method or must be recast, for example by a different non-dimensionalization, to be stable. I am leaving the issue open because I want to understand the details of where things are going awry, but cannot promise a solution.

sbryngelson avatar May 27 '22 14:05 sbryngelson

@boombiubla we may have found the issue with your configuration. Will update soon.

sbryngelson avatar Sep 27 '22 19:09 sbryngelson

Thank you!

HIT-CFD-group avatar Oct 02 '22 08:10 HIT-CFD-group

Hi @boombiubla, I was able to run your test case for 2500 time steps with dt = 2.3e-7 (I used 10 times coarser grid so that it runs fast) which is equivalent to 25,000 time steps with your original dt. NaNs do not show up !

The issue in this case seems to be the location of the boundary condition in the transverse direction, which is too close to the water bubble. This causes the negative pressures inside the liquid to escape into air, which in turn causes a NaN as pi_inf for air is 0

The fix is simply to move the location of bc_y%end to be far away from the bubble. For this purpose I set bc_y%end to be 0.1 instead of 0.037 in your case. However this will increase the cost of the simulation, but MFC also provides a fix for that by stretching the grid so that it is finer near the bubble and coarse away from it. I shall send you a modified input file with sufficient grid stretching later tonight to demonstrate this (this should also be available in the documentation)

anandrdbz avatar Oct 04 '22 22:10 anandrdbz

Thanks for your answer@anandrdbz. I tried another case, and found the pres_inf in the region of negative pressure is big enough to avoid pres+pres_inf less than 0 in the questioned timestep. And I tried to make the boundary far away from the reaction region, it helps but not enough. Decreasing the weno order helps but the low accuracy is not accepted for my need. 图片 The region with the infinity ccfl is around the shock wave instead of the air-water interface last time. 图片

By the way, is MFC-4.0.0 the updated version you said@sbryngelson? Is there any user-guide for it?

HIT-CFD-group avatar Oct 17 '22 02:10 HIT-CFD-group

We are looking into this @boombiubla -- can you comment on what else you tried? We believe the method is "doing its job", so for this problem it may be the the case that a smaller time step or larger domain is required. Also, consider using mp_weno and mapped_weno. The updated version of MFC is this repository. The user guide has a link on the README.

sbryngelson avatar Nov 09 '22 13:11 sbryngelson

The problem has been identified and has to do with a bug in the boundary conditions. @anandrdbz has fixed this and will push it soon.

sbryngelson avatar Dec 30 '22 04:12 sbryngelson

Hi @boombiubla, I have managed to make a 3D shock droplet case that works for this case. You can adjust the Mach number by changing the value of Min in case.py (currently set to 2.4 in the case file). It employs a stretched grid to move the boundaries far enough away from the droplet to prevent the NaNs

case.py.zip

anandrdbz avatar Feb 10 '23 17:02 anandrdbz

Closing this issue as the fixed

sbryngelson avatar Feb 14 '23 14:02 sbryngelson