1. Cremora Pot Operation:
A fine powder, such as Coffeemate, is layered over the propellant at about a 1:40 ratio of lifting charge to creamer. Upon ignition, the fine powder is lifted up and ignited into a fireball. Classic dust explosion; particles have very large surface area to mass ratio. This allows them to catch on fire with much less energy than would a bulk material.
2. Smoothed-particle Hydrodynamics (SPH) Overview:
Gridless lagrangian method for fluid flows. Discretises and performs calculations on individual particles rather than volume elements. I chose this method because it seemed to naturally fit the philosophy of a dust explosion, and a meshfree method seemed like a more natural way to program.
Particle elements represented as integrals, delta function approximated by a Kernel Approximation, where h is the smoothing lenght.
Now you can plug this in to physical relationships, such as Navier-Stokes.
3. Math and Code Representation:
void compute_density(flow* s, constants* args) {
int i = 0;
int j = 0;
int n = s->n;
float* __restrict__ rho = s->rho;
const float* __restrict__ x = s->x;
float h = args->h;
float h3 = h*h*h;
float h12 = ( h3*h3 )*( h3*h3 );
float C = 4 * s->mass / M_PI / h12;
memset(rho, 0, n*sizeof(float));
for (i = 0; i < n; ++i) {
rho[i] += 4 * s->mass / M_PI / h3;
for (j = i+1; j < n; ++j) {
float dx = x[2*i+0]-x[2*j+0];
float dy = x[2*i+1]-x[2*j+1];
float dz = x[2*i+2]-x[2*j+2];
float r3 = dx*dx + dy*dy + dz*dz;
float sub = h3-r3;
if (sub > 0) {
float rho_ij = C*sub*sub*sub;
rho[i] += rho_ij;
rho[j] += rho_ij;
}
}
}
}
Acceleration:
void compute_accel(flow* currFlow, constants* args) {
int i;
int j;
const float h = args->h;
const float rho0 = args->rho0;
const float k = args->k;
const float mu = args->mu;
const float g = args->g;
const float mass = currFlow->mass;
const float h3 = h*h*h;
const float* __restrict__ rho = currFlow->rho;
const float* __restrict__ x = currFlow->x;
const float* __restrict__ v = currFlow->v;
float* __restrict__ a = currFlow->a;
int n = currFlow->n;
compute_density(currFlow, args);
for (i = 0; i < n; ++i) {
a[2*i+0] = 0;
a[2*i+1] = 0;
a[2*i+2] = -g; //gravity
}
float C0 = mass / M_PI / ( (h3)*(h3) );
float Cp = 15*k;
float Cv = -40*mu;
for (i = 0; i < n; ++i) {
const float rhoi = rho[i];
for (j = i+1; j < n; ++j) {
float dx = x[2*i+0]-x[2*j+0];
float dy = x[2*i+1]-x[2*j+1];
float dz = x[2*i+2]-x[2*j+2];
float r3 = dx*dx + dy*dy + dz*dz;
if (r3 < h3) {
const float rhoj = rho[j];
float q = cbrt(r3)/h;
float u = 1-q;
float w0 = C0 * u/rhoi/rhoj;
float wp = w0 * Cp * (rhoi+rhoj-2*rho0) * u/q;
float wv = w0 * Cv;
float dvx = v[2*i+0]-v[2*j+0];
float dvy = v[2*i+1]-v[2*j+1];
float dvz = v[2*i+2]-v[2*j+2];
a[2*i+0] += (wp*dx + wv*dvx);
a[2*i+1] += (wp*dy + wv*dvy);
a[2*i+2] += (wp*dz + wv*dvz);
a[2*j+0] -= (wp*dx + wv*dvx);
a[2*j+1] -= (wp*dy + wv*dvy);
a[2*i+2] -= (wp*dz + wv*dvz);
}
}
}
}
Other cool functions:
First Code Version: Defining a Particle System
Written in C, with GLUT. In this version pressing space creates a simple explosion of particles from the center. I had trouble getting the math to connect to the particle system.
Second Code Revision: Implementing SPH
In this version an SPH-governed group of particles interacts with each other in 3d.
4. Conclusions and Future Work:
I learned all about SPH! It's really intuitive, moreso that grid solutions I think. Also, I learned that C is really terrible for math when you have numbers that might blow up, or become imaginary. Also, I learned to start final projects earlier.
In the future, I'm going to fix the math problems with my code, and seperate some threads so I can run with more particles on my GPU.