Help With Divergence Issue

Hi everyone!

I was hoping to get some help with my model. The run seems to diverge no matter what I change.

I have tried to take out all the complexity including turning off chemical equations, increasing dt_max and decreasing dt_min, turning off heat transfer, and I am still getting the same errors.

I printed out the reaction rate and concentration of my reactant as well. They are reasonably small. I assume divergence happens when something is too large, but I couldn’t find anything too large.

Any help is greatly appreciated, thanks in advance!

I’m attaching all of my files here:

modify_f.py (2.1 KB)
td_jelly_roll_rev.mfx (21.4 KB)
usr_rates.f (2.9 KB)

Hi Sofiya.

I’m not sure what the file modify_f.py is for so I’m ignoring that one.

I’m pretty sure that your reaction rate is several orders of magnitude too high. If I disable reaction 1 (uncheck it), the simulation runs:

Looking at your usr_rates.f I see there is no else statement here:

      IF(X_s(IJK, 1,RM) > c_Limiter) then
	     c_RM  = (RO_s(IJK,1)*X_s(IJK,1,RM)/MW_s(1,RM))
      ENDIF

so c_RM could have an undetermined value. Adding ELSE c_RM = 0 did not fix the problem.

I changed

! k for reaction rate
      DOUBLE PRECISION, parameter :: k=1d-16

and was able to get the simulation to run:

I set up a monitor for max. reaction rate across the entire domain, that’s the plot that’s being shown here. If I increase k I get a divergence.

– Charles

Thank you Charles! I’m confused about why the else statement did not work. I thought that c_RM was linear related to (X_s(IJK, 1,RM), so in my mind as long as (X_s(IJK, 1,RM) is not too small, c_RM will be kept in check. Did I misunderstand something?

So the if statement gives a lower bound of c_RM. We also checked the upper bound. We did print out reaction rate, which is linearly proportional to c_RM. The number was reasonable so we thought c_RM was fine. So I’m very confused about what actually is blowing up. Do you have any suggestions on what to check?

Also we ran the simulation without chemical reaction and it was fine, so the problem lies within the reaction. We tried k=1d-16 and it was running without diverging. But we can’t just turn the k values to 1d-16 because it’s definitely not physical. There must be something else that we need to change in order to get the code running without a small k. Any input is appreciated!

Thanks,
Jack

I only took a quick look, but I can see a few issues:

  1. You need to save data in vtk files so you can visualize the result. This will generally help troubleshoot. For example, verify that the case without reaction provides expected hydrodynamics.
  2. You set a solids pressure in the Initial Condition regions to 101325Pa. The solids pressure is not the same as the gas phase pressure. I don’t think this is the cause of non-convergence because you have a fixed bed, but better to leave it undefined.
  3. In the usr_rate.f file, you do a test on `X_s(IJK,1,RM). This is not sufficient for TFM. Each cell will hold a value for RM specie fraction even if there is no solids in that cell. It is better to test first if there are solids in the cell, say ep_s(IJK,1)>some small number. Similarly, it is better to use ROP_S(IJK,1) i.e. the bulk density, rather that the solids density RO_S(IJK,1).

Now the convergence issue seems to come form the solids fraction equation. I was able to get it to run by checking “Disable close pack” in the Solids>Material pane, Advance section (at the bottom). Should be fine since you have a fixed bed that is already packed as per your packed be void fraction of 0.7.

1 Like

Thanks a lot Jeff! I turned off the close packing, and it’s working now.

1 Like