Bug in simulations involving superquadratic particles with zero gravity

MFiX Bug report
#1

Dear MFIX developers,

I recently encountered a bug in the MFIX version 23.1 when I ran a DEM simulation involving superquadratic particles and a gravity vector set to zero in all directions.
The bug occurred in the subroutines Sq_angle_with_gravity in the file des/sq_properties.f.
Using zero gravity led to a division by 0.

I solved the problem by modifying the line:
“gravity=gravity/gravity_norm”
to:
“IF(gravity_norm == 0.0) THEN
gravity=(/0,1,0/)
ELSE
gravity=gravity/gravity_norm
ENDIF”

I don’t know if it’s the most appropriate solution for everyone, but I wanted to inform you of this little bug.

Thanks for your great work.

Best,

Alex

#2

Thanks for the suggestion. Your change should be safe - when the magnitute of the g-vector is zero, the angle is undefined, so it shouldn’t matter what value the function Sq_angle_with_gravity returns in that case.

I would probably write:

IF(gravity_norm == 0.0) THEN
    Sq_angle_with_gravity = 0
ELSE

i.e. skip the rest of the calculation and just return 0 if |g| is 0.

But, I took a closer look at the code to see where Sq_angle_with_gravity is used, and it seems that this value is not used at all. The call sites look like this:

! Orentation angle
     call Sq_angle_with_gravity(Qc, axi,m,n, cos_angle)
!     des_usr_var(10,NP)=cos_angle

There is a commented-out line storing cos_angle into des_usr_var but other than this, the variable cos_angle is not used. So I think the right fix is to remove these unneeded calls completely. I will do this for the next MFiX release.

Thanks again for your report.