MFC
MFC copied to clipboard
All physical quantities become "infinity".
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.

These are my input file. input.zip
What I am doing wrong? Could you advice, please?
Interesting, thank you for pointing this out. We will look shortly! @anandrdbz @mrodrig6 would you mind having a look as well?
@boombiubla did you consider a smaller time step?
@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
@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 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?
@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 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.🤔

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.
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?
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.
Would you mind sending me some input files for learning? Maybe files about JC Meng's works are better.
There are example input files available in the repository. I do not have Meng's input files on hand.
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
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.
Hmmmm, I tried once again like below. But it seems not successful. Is there any wrong in my input file? norm.zip
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.
OK, thanks for your attention and suggestions!
@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.
@boombiubla we may have found the issue with your configuration. Will update soon.
Thank you!
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)
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?
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.
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.
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
Closing this issue as the fixed