PySPH Attempting to implement the Mixture Particle model - AlexanderPuckhaber/FluidSimulationThesis GitHub Wiki
The Mixture Particle Model
The volume fraction technique outlined in [insert source] uses a Mixture Particle, a SPH particle that includes multiple phases/materials within itself using the Volume Fraction method.
It has all the standard SPH properties such as position, velocity, acceleration, mass, density, volume, etc. However, the different phases/materials can have different properties. For example, PySPH keeps track of spring force data for solids in a dozen variables. The Mixture Particle must contain all data for the phases/materials, data on the volume fractions of the phases/materials, as well as mixing velocities of the phases/materials for modeling diffusion. Finally, it must also condense the diverse properties into a total volume, mass, density, etc. for the ordinary SPH interactions.
Using the MixtureParticle class to make custom Mixture Particles
The MixtureParticle class is designed to hold all the data for the template of a Mixture Particle. It is constructed with base properties, and different phases/materials are added to it.
It generates a dict of properties that looks like this example [water and oil]:
{'x':0, 'y':0, 'z': 0, 'h':5, other base values, 'water__frac':0.8, 'water_m': 0.1, 'oil__frac':0.2, 'oil_m':0.09}
These properties and default values can be used in the create_particles method to make a bunch of Mixture Particles.
Equations of State for Mixture Particles
[put more here, link stuff]
The Volume Fraction method requires that all particles maintain near-constant volume. Depending on the properties of the component phases/materials, properties such as overall mass, density, etc. can change; however for the Volume Fractions to have meaning, the Mixture Particle must have a constant volume [find source].
I wanted to simulate a mutiphase interaction between fluids of different densities, so I divided the fluid in hydrostatic_tank.py into two halves, and made the top half more dense than the bottom one. I multiplied the mass of the particles in the bottom half by 0.5, assuming that it would change the density. This is the result:

It appears successful - the high-density "water" begins to switch places with the low-density "oil" to become more stable. However, it seems a bit bouncy, so I wondered if the method could scale up. Here is an animation of oil and water where the oil particles are 0.05 times the mass of the water particles.

The simulation is visibly unstable. It was clear that this approach could not scale and was not accurate because it changed the Volume and Pressure of the particles, and left the density unchanged. I needed to find a way to change the density of the particles. Unfortunately, directly changing the density only worked for the first frame of the simulation: afterwards, the density of the fluid became uniform and equivalent to rho0, or the water density. When oil was set to 1.0rho0 and the density of water was increased to 1.1rho0, the water particles jumped up and out of the tank. Clearly, there was some Equation that was changing other variables to satisfy a constant rho0 value for other particles.
Unfortunately, PySPH doesn't seem to have good support for weakly compressible liquids with different densities (It probably has good support for gases because it has programs to handle shockwaves and such, but I haven't investigated). The rayleigh_taylor.py program simulates the interaction between fluids with different densities, where the less dense fluid is below the fluid with greater density. The potential energy leads to a Rayleigh-Taylor instability. The example program creates a process_term that differentiates the target densities (rho0) and pressures (p0) of the liquids and injects it into the equations.
Requiring the user to create a process_term to handle different density phases/materials is OK to do occasionally, but it is a required feature of all Multiphase (Mixture Particle) simulations. The author of PySPH has noted "this is an ugly hack" about injecting parameters for different rho0 values and suggests that in the future, the equations would use the unique rho0 value for each particle. So I looked at all the Equations in (hydrostatic_tank.py)[] and modified an equation involving rho0, TaitEOSHGCorrection, and created TaitEOSHGCorrectionVariableRho to handle unique rho0 values. But it didn't work. So I printed out every Equation used in rayleigh_taylor.py and found that an Equation called StateEquation was automatically generated in the Scheme TVFScheme.
So I modified the StateEquation, making it use the individual rho0 of each particle for computing pressure and renamed it StateEquationVariableRho0. I put it in a modification of TVFScheme and tested the new method.