8.5. How to setup PIC simulations

This document explains how to setup PIC simulations (statistical weight, PIC parameters, mesh).

8.5.1. How To Choose a Good Statistical Weight for PIC Runs?

One key to creating efficient and accurate particle in cell (PIC) simulations is choosing an appropriate value for the statistical weight of particles per parcel. The user needs to balance desired solids weight fraction, mesh density and solids mass error with computational speed and physics.

The PIC statistical weight in MFIX (IC_PIC_CONST_STATWT or BC_PIC_CONST_STATWT) represents how many particles are represented by each parcel. The parcels themselves must be of discrete number in a simulation, but because statistical weight is a math construct, partial particles may occur mathematically.

To facilitate statistical weight selection, it is desirable to create a spreadsheet that illustrates how its selection can impact simulation conditions, or the following logic may be prescribed.

  1. Collect the following information:

    1. Diameter of an equivalent single spherical particle (m)

    2. Mass density of particles (kg/m 3 )

    3. Max-pack solids fraction of particles for simulation

  2. Choose a mesh for the simulation. Separately, calculate an average cell volume (m 3 ) and a minimum cell volume (m 3 ).

  3. Using max-pack solids fraction (which must be defined for PIC, 1-EP_STAR), calculate an average cellular volume available for solid (m 3) and a minimum cellular volume available for solid (m 3).

  4. Using particle diameter, calculate the number of spherical particles these average cellular volumes available for solid represent. Assuming a uniform mesh, this calculated number, as a rounded-down integer, would represent the maximum statistical weight one can choose in a PIC simulation, although it is highly unlikely one would choose this value. In a variable mesh, the value associated with the minimum cell volume is the maximum statistical weight. (Note: using an average cellular volume for this calculation allows a balance of parcel load in both large and small cells.) Note that a user does not need to know the solids fraction distribution that they expect in a simulation, only the max-pack solids fraction, which is a constant.

  5. Now the balancing act between parcel count, actual solids fraction and mass error begins. First, choose a statistical weight value. An easy way to make a first approximation is to use the maximum statistical weight from (4) and divide by an integer which in your mind represents the number of parcels you want per cell. Statistical weights derived like this will naturally give a solids fraction very close to whatever the original desired value from (1) was, and minimize mass error in a simulation. The statistical weight does not need to be an integer, because partial particles make sense in PIC. While the parcel count in an individual cell will initially begin as an integer, there is no limitation that says an integer number of particles must make up a parcel. Choosing an integer value for statistical weight does make it simple to imagine the parcels as grouped entities of particles, but it may interfere with expectations for representative solids volume fraction per cell. Note, however, that regardless of what the user specifies as a statistical weight, the MFIX program will collect truncated mass values in each cell until enough mass to represent a full parcel accumulates, whereupon a single extra parcel will be placed in the domain.

  6. It is advisable to inspect your fully populated PIC parcel domain to make sure you have not applied a value that does not make visual or physical sense. This is particularly true for very dilute flow and highly packed flow. For dilute flows, choosing a statistical weight that is too large will result in a simulation that cannot possibly capture expected flow dynamics. In a highly packed flow, choosing a statistical weight that is too low may result in non-physical packing, and an overly large local solids fraction. As such, be certain to note the values reported in the MFIX output for desired solids fraction and calculated solids fraction, then assure yourself that the numbers are within some reasonable error band.

8.5.2. Extra info and Examples:

During initialization, the user can choose any region to fill with parcels. If a distinct mass of particles is needed in a region at start-up, the user must pre-conceive an appropriate volume of cells and associated solids fraction to achieve that mass.

For example, if a user wants :

  • 0.1 kg of 0.001 m particles with density 1000 kg/m 3

  • in a volume that is 0.1 x 0.1 x 0.1 m 3 defined by a 10 x 10 x 10 mesh

  • with max-pack of 0.4 (theoretical limit of random packing of spheres is 0.36)

First, calculate the volume of 1 particle: V_particle= \(\pi\)/6 (0.001)3 = 5.236E-10 m 3

Calculate the mass of 1 particle: m_particle=(1000 kg/m 3 )(5.236E-10 m 3) = 5.236E-7 kg

Calculate single cellular volume: V_cell=(0.1m) 3/((10×10×10))=1.0E-6 m 3

Calculate available volume for parcels: (1-0.4)( 1.0E-6 m 3) = 6.0E-7 m 3

For simplicity, we assume that the particles are evenly distributed among the 1000 cells. That means, each cell needs 0.100 kg/1000cells = 0.0001 kg of solid. In turn, that means each cell needs 0.0001 kg / 5.236E-7 kg = 190.985 particles. From this particle loading, you can choose a statistical weight based on the number of particles. For example, you could choose 190.985 as the statistical weight and get 1 big parcel (unwise), or you could choose, say, a statistical weight of 10 and get 19 parcels with a small local mass error. The small errors are accumulated, and when large enough, an extra parcel will be placed locally in the mesh to account for the local loss. Note that if the user does not want to fill an initial area at maximum packing, simply choose a different value for solids fraction and recalculate. Having a preconception of how “loose” a particle bed should be is very helpful for reducing the transient in a simulation start-up.

Likewise, if a mass inflow boundary is prescribed in the PIC model, a similar logic must be followed. The key is to make sure that the statistical weight used for the boundary flow condition is reasonable (as in, the number of particles per parcel will fit into the cell(s) that abut the boundary condition specification). MFIX will control the mass itself by accounting individual particle mass and randomly assigning parcels through the BC window. Note that smaller statistical weights produce less local mass error, but the user must choose their own balance between statistical weight and computational speed.

8.5.3. Meshing Basics for PIC Runs

There are three common mistakes a user may make when setting up a PIC run. The first represents an entire misunderstanding of MFIX-PIC which is necessarily a 3-D model. That means, you cannot specify a 2-D simulation with MFIX-PIC. Second, the interpolation algorithm which undergirds the MFIX-PIC model requires a minimum of 3-cells in each direction through a geometry to give a robust solution. It is important to check small geometric regions of the mesh that might pass through an orifice or channel to assure that there are at least 3 cells available for calculations to work in a hydrodynamically stable way. Smaller cell distributions will not crash the model, but you may get unexpected flow behavior. Finally, choosing statistical weights that in comparison to mesh give less than one parcel per mesh cell will give errant solutions. (See section on choosing good statistical weights.) Note that one of the benefits of PIC is that it allows for a somewhat coarse mesh with a large particle count, but a user needs to balance mesh and parcel count in a reasonable way.

8.5.4. Choosing Primary PIC Parameters

PIC parcels do not behave like DEM particles. There are no Newtonian physical laws that govern collisions. Instead, PIC utilizes a solids stress ( \(\tau_p\)) field that influences solids motion through a coupling with the momentum equation. Specifically, the solids stress is calculated as:

\(\tau_p=\frac{P_p \epsilon_p^\gamma}{\max[\epsilon_{cp} - \epsilon_p , \delta(1-\epsilon_p)]}\)

In MFIX, the parameters that the user can control to define this solids stress are an empirical pressure constant (PSFAC_FRIC_PIC ( \(P_p\) ), a volume fraction exponential scale factor (FRIC_EXP_PIC ( \(\gamma\)), void fraction at maximal close packing(EP_STAR ( \(\epsilon_{cp}\)) and a non-singularity constant (FRIC_NON_SING_FAC ( \(\delta\)). The only other variable at play within the solids stress model is solids void fraction ( \(\epsilon_p\) ) which is calculated dynamically as the simulation progresses.

The parameters used in PIC are very robust. The default values of \(P_p\) =100, \(\gamma\) =3, \(\delta\) =1.0E-7 should work well for most systems. There should be no reason to change \(\delta\). The value for \(\epsilon_{cp}\) will vary depending on your application, although it is worth noting that there is a theoretical limit of ~0.36 for the random packing of equal spheres. It is important to note that \(\epsilon_{cp}\) represents the void or gas fraction per cell, and not the solids fraction. The effects of changing the values of \(P_p\) and \(\gamma\) are very problem specific. For example, sedimentation cases (on which the solids model was originally built) will model best with \(P_p\) =1, \(\gamma\) =2 so that gravity has the most effect on the system through the momentum equation. But, very dynamic cases (U_mf>10) where we expect the solids to more greatly affect the hydrodynamics of the system may respond better with selections like \(P_p\) =1000, \(\gamma\) =3. Essentially, like all CFD models, it may be necessary to tune the PIC model to your specific case. As a rule of thumb, larger \(P_p\) at constant \(\gamma\) will result in particles that move more dynamically (as the selection will linearly increase the solids stress). But, larger \(\gamma\) at constant \(P_p\) will reduce the solids stress exponentially and dampen overall solids motion. The largest values of solid stress will occur when \(P_p\) is large and \(\gamma\) is small.

8.5.5. Choosing Secondary PIC parameters

Like DEM, the user can control restitution factors that affect parcel motion. MPPIC_COEFF_EN_WALL applies a scale factor to the normal component of velocity after a parcel interacts with a wall. Likewise, MPPIC_COEFF_ET_WALL applies a scale factor to the tangential components of velocity after a parcel interacts with a wall. There are two other empirical factors available in MFIX-PIC. The first is MPPIC_COEFF_EN1 which is an overall scale factor that can dampen the effect of the entire solids stress model, and MPPIC_VELFAC_COEFF which is a scale factor that effects the value of bulk solids velocity through an evaluation of gas-solids slip velocity. These last two factors should remain at default values unless you are working as an advanced user and understand the open source PIC code implicitly.

8.5.6. Collision Damping

An empirical collision damping feature has been added to the MFIX-PIC model. The feature includes a calculation for collision frequency based on statistical likelihood. From this frequency, a small change in velocity is calculated and applied to parcels locally. The effect is an overall damping of solids velocity, particularly in areas where dense particle flow is apparent. To invoke the feature, the user must select “Collision Damping” in the GUI (Solids>PIC pane) or set the keyword is PIC_COLLISION_DAMPING=.TRUE. in the .mfx file.

8.5.7. Frequently Asked Questions

FAQ: Why can’t I define parcels per cell instead of particles per parcel?

When creating the MFIX-PIC model, programmers chose to define particles per parcel for two reasons.

  1. The user has full control over the statistical weights of the PIC parcels and is not reliant on the computer code to calculate those values without knowledge of the simulation being run.

  2. In most CFD simulations, mesh can be highly variable. When a user defines parcels/cell, the smallest cell will limit the number of parcels that can be placed in any cell at initialization. By defining particles/parcel, the check we suggest for examining the smallest cell only limits the overall size of the parcel, and a more even solids void fraction can be achieved at start-up. An argument can be made that multiple initialization regions can be created to smooth solids void fraction using the parcel/cell method, but this only complicates set-up.

FAQ: What is a good number (or range) for the number of parcels per cells

It depends on the kind of simulation you are trying to run, how many processors you have available to you, and what kind of fidelity you expect. Avoid going below 3 parcels/cell in a serial run where there is very little dynamic (like settling).

Typically, around 10 parcels/cell works smoothly. If the simulation looks chunky, reduce the number of particles/parcel, thus upping the parcels/cell count.

FAQ: What solids model parameter values work best for the PIC model when running 2 densities that are roughly 10x different in magnitude?

A density difference like this (say 300 kg/m 3 and 3,000 kg/m 3) will cause natural physical separation in most systems. To capture the best physics in a PIC model where density drives strong gravitational effects, choose \(P_p\) =1 and \(\gamma\) =5 to begin. You may still need to adjust these parameters but these choices will give you a good place to start.