Error in running SQP simulation

Hai Developer,
I have been trying to simulate the flow over pine needles bed (of dimension 0.5 m x 0.08 m x 0.022 m) using SQP model without considering any chemical reactions. The cylindrical particle of diameter 0.879 mm and height 5 mm is used. The minimum cell size is chosen such that it is more than 1.63 times the particle diameter (bounding diameter of 10 mm). The X_min and X_max are given as pressure outlet boundary conditions and all other faces of the enclosure box are given as wall. I could not run the simulation successfully. I am facing the following two issues while running the simulation.

(1) Pine needle bed is placed on top of the asbestos wall (of dimensions 0.68 m x 0.08 m x 0.011 m). I created the asbestos wall (impermeable) as six planes and include them using the internal surface option. When I start the simulation with this, I am getting error like Dt >Dt_min.

(2) When I removed the asbestos wall and considered pine needles bed is placed on the bottom wall of the domain, then the simulation starts to run but suddenly stops and throwing the floating point exception- erroneous arithmetic operation error. The error report is also attached below along with the mfix input file for your kind reference.

Error: Unexpected solver message:
Program received signal SIGFPE: Floating-point exception - erroneous arithmetic operation.
Backtrace for this error:
#0 sq_rotation_mod_MOD_qincrotate2
** at des/sq_rotation.f:624**
#1 cfnewvalues_mod_MOD_superdem_cfnewvalues
** at des/cfnewvalues.f:388**
#2 des_time_march_MOD_des_time_step
** at des/des_time_march.f:219**
#3 run_dem
** at mfix.f:214**
#4 run_mfix
** at mfix.f:149**
#5 main_MOD_run_mfix0
** at main.f:85**
Floating point exception (core dumped)

I request you to kindly help me in solving these issues and please explain me why these errors are coming and where I am going wrong.

Thanking you,

Regards,
Shanmugasundaram D
Pine_bed.mfx (17.4 KB)

A couple of things you could try:

  1. Turn off the floating point exception traps (Numerics>Advanced pane)
  2. Use a rounder corner (set the shape exponent n=4)

Hai Jeff Dietiker,
Thanks for your response. I have tried both the things you suggested. Now, it is not showing error immediately and it is running for few time steps. After 0.025 s, the console printed " overlap between two superquadric particles is too large!". Also, the velocity of the particles are nonphysical (like it varies between 16 to 60 m/s). Moreover, particles are flying and scattered and finally simulation stopped running at 0.05 s.
Why is it showing as “overlap between two superquadric particles is too large!” ? Where am I going wrong?.
Can I atleast arrest the rotational component of particle velocity alone and making other two components to be present?
My Mfix file is also attached below for your kind reference.
Pine_bed.mfx (20.1 KB)

Thanks,

A few additional comments on this case:

  1. Please upgrade to the latest MFiX, version 24.1 A few SQP-related bugs have been fixed.

  2. Your case is using the SYAM_OBRIEN drag model which is not compatible with SQP. In MFiX 24.1 this leads to a warning:

Drag type SYAM_OBRIEN not compatible with MFiX-SQP solver, setting to SQP_DIFELICE_GANSER
  1. Your particle seeding is not succeeding. This also creates an error popup in MFiX 24.1, but is easy to miss in earlier versions. Currently, achieving a high volume fraction with such elongated particles is not feasable, because the seeding algorithm uses the bounding sphere, rather than the actual particle shape, to determine possible overlaps. You may have to either use an external script to generate the particle input file, or else seed the particles over the entire region and let them settle down into the bed.

  2. Regarding the numerical overflows - we know that the SQP contact detection algorithm tends to produce very large values, which is why we added the option to turn off overflow checking. But the reported overflow was in a rotation algorithm - sq_rotation_mod_MOD_qincrotate2 - which I’ve never seen before. I’m not sure that turning off overflow checking is advisable in this case - I’d like to know how we’re getting invalid operation out of a rotation. However, I was not able to reproduce the reported bug with MFiX 24.1 (it may have been fixed already)

  3. How long did you have to run the simulation before the error happened? It’s helpful if you submit a full bug report (Main menu → Submit bug report) which includes, as well as your input files, all the job logs. This is very useful to us in debugging.

I suggest the following course of action:

  • Upgrade to MFiX 24.1
  • Re-enable overflow checking
  • Fix particle seeding issue
  • If model still fails, submit full bug report.

Thanks!

– Charles

Hai Charles,
Thanks for your detailed response. I followed your suggestions and now the particle seeding is successful (after reducing the particle length). The same message of “overlapping between two SQP particles are too large” has been shown in the console. This time, I am not getting the same error related to rotation algorithm but rather some other error I am getting with MFIX 24.1. The bug report of the same is attached below for your kind reference.

  1. This particle seeding error makes me to change either the initial bed volume fraction or the particle size, which I don’t want to do. It’s affecting the initial condition of the simulation.

  2. In the previous mail, you mentioned to “seed the particle in the entire domain allow to settle down in the bed”. I don’t understand this statement. If I seed the particle in the entire domain, how can I get my bed initial conditions correctly (for eg. volume fraction of particles in the bed though I maintain the mass in the domain constant). I request you to elaborate more about this procedure.

  3. In my pine bed simulation, I need to have correct initial bed volume fraction. If I write a particle input file, can I avoid this overlapping error and have the particle size and initial volume fraction according to my need? If so, could you please suggest me a suitable starting point (example files or documents related to it) for beginning to write the particle input file?
    Kindly suggest a suitable way to seed the particle with the required volume fraction and particle size.

Thanking you,

Regards,
Shanmugasundaram D
pine_bed_2024-04-02T230711.536862.zip (1.9 MB)

A few comments:

  1. The seeding algorithm in MFiX is lattice-based and won’t do well with elongated particles, since the seeding lattice uses the diameter of the bounding sphere … but your pine needles of course can get much closer than that.

  2. By “seed the entire domain” what I meant was to use random seeding in the area above the desired bed location, and to let the particles settle due to gravity. Since you can control the number of particles, you should be able to achieve the desired density in the bed, once it settles.

  3. You should be able to write a Python script to emit particle locations and orientations. I am also working on some code for seeding SQP particles, which should be available in a few days or so.

  4. I noticed that you have changed the shape of the pine needles and they are much less elongated now… and it looks like you also fixed the problem where one of your inlets was set up as an outlet :slight_smile: I hope we can get this running with the original particle shape!

– Charles

Hai Charles,
Thanks for your valuable comments. I followed your suggestion and attempted few things accordingly.
(1) I reduced the particle size to 2 mm and carried out the simulation and very quickly (less than 0.001 s) simulation ended with floating point exception and erroneous arithmetic operation. The detailed bug report of the same is attached below for your kind reference.
pine_bed_2024-04-05T105150.188720.zip (1.7 MB)

(2) As per your another suggestion, I increased my bed height to 150 mm (from 20 mm) and seeded the particles in that entire region by specifying the number of pine particles. In this case, the simulation starts to run nicely until 0.016 s. But, again it stopped with the errroneous arithmetic operation. The detailed bug report for this case is also attached below for your kind reference.
pine_bed_2024-04-05T112212.996761.zip (2.1 MB)

(3) In extended bed height case, When I change the DES neighbor grid search option to grid based (from N-square), the simulation ran for few time steps and stopped quickly by throwing the below shown error.


pine_bed_2024-04-05T120222.259976.zip (1.6 MB)

I hope that I have set all the conditions physically correct and followed your suggestions. But could you please let me know why these errors are coming.

I request you to kindly let me know where I am going wrong and suggest me some ways to avoid these errors.

Thanking you,

Regards,
Shanmugasundaram D

@Shanmugasundaram
I want to thank you for these detailed and complete bug reports. I have been running your cases and some SQP cases of my own and I can reproduce these failures pretty reliably. As you know, there is an option to disable checks for floating-point exceptions, you can try enabling this as a workaround, but I’m not sure if it will affect the validity of your results. It would really be best to get to the bottom of what’s causing the exceptions and fix the code, and I am determined to to so! - however the code is fairly complex and the original author is no longer at NETL, so this will take a bit of time. I will post more here soon, I just wanted to let you know that we are looking into it and hope to get this resolved before long.

Thanks for your patience,

– Charles

“Just for fun”, here’s a successful SQP run which shows settling of a randomly-seeded bed - I was experimenting with different approaches to generating the particle configuration, this particular approach only was able to reach a density of about 0.15 - however you can see that after it settles the height is approximately 1/4 of the original, so the density is closer to 0.60. (These are superquadric particles with parameters 0.01, 0.01, 0.05, 2, 8)

1 Like

Hai Charles,
It is really nice to see the settling of SQP particles in the above video. I am wondering to know about the approach you used to run this case successfully. Is it the same geometry what I have used for my simulation? If so, please let me know where I am going wrong?
I am really thankful to you and NETL community for putting sincere attempts in resolving the issue.
Though my objective is to simulate the pine needles movement during its combustion, due to the problem I am facing currently , I am trying to simulate the case with static particles as the initial step. By doing so, I have the following doubts,

(1) At present, In order to achieve ignition, I specified some portion of pine bed with a higher temperature of 1700 K. But, due to air flow from the mass inlet boundary reduces the temperature and ignition is not achieved properly. In avoiding this situation, how can we specify the igniter (constant heat flux for specified time) to initially ignite the particles? The point source option does not have option to specify heat flux. Please let me know the way to achieve it.

(2) In order to specify the high temperature to particular portion of the bed, I splitted the bed into two regions (ignition bed and pine bed). The ignition bed dimension extends (along x-axis) from -0.25 to -0.2 m and the pine bed starts immediately from same -0.2 m and runs upto +0.25 m. In doing so, I see (in VTU) these two regions get splitted with some gap in between them (the gap is roughly about one cell size, I guess). Why is this gap coming? Due to this, heat transfer is not happening properly and thus I could not achieve the ignition. In order to avoid this gap, do I need to activate any advanced option keywords?? The image of the same is attached below for your kind notice.

(2) Due to some heat transfer, the temperature of some portion of the cold pine needles gets increased to around 400 K. I see the pyrolysis started to happen from the top layer of the pine needles. But, the pyrolysis rate is showing zero for the particles (bottom layer) which are next to the wall though their temperatures are around 400 K. Am I specifying any wrong input and please let me know how to correct it? The complete case files and the corresponding screenshots of the temperature and pyrolysis rate of the pine bed are attached below for your kind reference.


usr_rates_des.f (8.4 KB)
usr_rates.f (4.6 KB)
usr1_des.f (1.5 KB)
Pine_bed_old.mfx (35.7 KB)

Thanking you,

Regards,
Shanmugasundaram D