Working fluid sim
This commit is contained in:
parent
ec055379d6
commit
7ab4168a5c
2247
hs_err_pid559270.log
2247
hs_err_pid559270.log
File diff suppressed because it is too large
Load Diff
@ -8,6 +8,7 @@ import java.util.LinkedList;
|
|||||||
import java.util.List;
|
import java.util.List;
|
||||||
import java.util.Random;
|
import java.util.Random;
|
||||||
import java.util.Set;
|
import java.util.Set;
|
||||||
|
// import jdk.incubator.vector.FloatVector;
|
||||||
|
|
||||||
import org.joml.Vector2i;
|
import org.joml.Vector2i;
|
||||||
import org.lwjgl.BufferUtils;
|
import org.lwjgl.BufferUtils;
|
||||||
@ -35,31 +36,39 @@ public class FluidSim {
|
|||||||
|
|
||||||
public static final int DIM = 18;
|
public static final int DIM = 18;
|
||||||
|
|
||||||
float[][][][] data = new float[4][DIM][DIM][DIM];
|
// float[][][][] data = new float[4][DIM][DIM][DIM];
|
||||||
float[][][] x = new float[DIM][DIM][DIM];
|
float[][][] x = new float[DIM][DIM][DIM];
|
||||||
float[][][] x0 = new float[DIM][DIM][DIM];
|
float[][][] x0 = new float[DIM][DIM][DIM];
|
||||||
float[][][] u = new float[DIM][DIM][DIM];
|
float[][][] u = new float[DIM][DIM][DIM];
|
||||||
float[][][] v = new float[DIM][DIM][DIM];
|
float[][][] v = new float[DIM][DIM][DIM];
|
||||||
float[][][] w = new float[DIM][DIM][DIM];
|
float[][][] w = new float[DIM][DIM][DIM];
|
||||||
float[][][][] oldData = new float[4][DIM][DIM][DIM];
|
float[][][] u0 = new float[DIM][DIM][DIM];
|
||||||
|
float[][][] v0 = new float[DIM][DIM][DIM];
|
||||||
|
float[][][] w0 = new float[DIM][DIM][DIM];
|
||||||
|
// float[][][][] oldData = new float[4][DIM][DIM][DIM];
|
||||||
|
|
||||||
static final float DIFFUSION_CONSTANT = 0.01f;
|
static final float DIFFUSION_CONSTANT = 0.0001f;
|
||||||
|
static final float VISCOSITY_CONSTANT = 0.0001f;
|
||||||
|
|
||||||
|
static final int LINEARSOLVERTIMES = 20;
|
||||||
|
|
||||||
|
static final float GRAVITY = -10000f;
|
||||||
|
|
||||||
public void setup(){
|
public void setup(){
|
||||||
Random rand = new Random(0);
|
Random rand = new Random(1);
|
||||||
//make a cube of water in the center
|
//make a cube of water in the center
|
||||||
for(int x = 0; x < DIM; x++){
|
for(int i = 0; i < DIM; i++){
|
||||||
for(int y = 0; y < DIM; y++){
|
for(int j = 0; j < DIM; j++){
|
||||||
for(int z = 0; z < DIM; z++){
|
for(int k = 0; k < DIM; k++){
|
||||||
if(
|
if(
|
||||||
Math.abs(8 - x) < 4 &&
|
Math.abs(16 - i) < 4 &&
|
||||||
Math.abs(8 - y) < 4 &&
|
Math.abs(8 - j) < 4 &&
|
||||||
Math.abs(8 - z) < 4
|
Math.abs(10 - k) < 4
|
||||||
){
|
){
|
||||||
data[0][x][y][z] = 1;
|
x[i][j][k] = 1;
|
||||||
data[1][x][y][z] = rand.nextFloat() * 0.01f;
|
u[i][j][k] = rand.nextFloat() * 0.1f;
|
||||||
data[2][x][y][z] = -1f;
|
v[i][j][k] = -1f;
|
||||||
data[3][x][y][z] = rand.nextFloat() * 0.01f;
|
w[i][j][k] = rand.nextFloat() * 0.1f;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
@ -69,63 +78,72 @@ public class FluidSim {
|
|||||||
/**
|
/**
|
||||||
* Runs a frame of the fluid simulation
|
* Runs a frame of the fluid simulation
|
||||||
*/
|
*/
|
||||||
public void simulate(float timestep){
|
public void simulate(int step, float timestep){
|
||||||
//diffuse vectors
|
//diffuse vectors
|
||||||
moveData(false, 1);
|
// moveData(oldData[1], data[1]);
|
||||||
diffuse(timestep, 1, 1);
|
// diffuse(timestep, data[1], oldData[1], 1);
|
||||||
moveData(false, 2);
|
// moveData(oldData[2], data[2]);
|
||||||
diffuse(timestep, 2, 2);
|
// diffuse(timestep, data[2], oldData[2], 2);
|
||||||
moveData(false, 2);
|
// moveData(oldData[3], data[3]);
|
||||||
diffuse(timestep, 3, 3);
|
// diffuse(timestep, data[3], oldData[3], 3);
|
||||||
//project
|
// //project
|
||||||
project();
|
// project(data[1],data[2],data[3],oldData[1],oldData[2]);
|
||||||
moveData(false, 1);
|
// moveData(oldData[1], data[1]);
|
||||||
moveData(false, 2);
|
// moveData(oldData[2], data[2]);
|
||||||
moveData(false, 3);
|
// moveData(oldData[3], data[3]);
|
||||||
//advect vectors
|
// //advect vectors
|
||||||
advect(timestep, 1, 1);
|
// advect(timestep, 1, 1);
|
||||||
advect(timestep, 2, 2);
|
// advect(timestep, 2, 2);
|
||||||
advect(timestep, 3, 3);
|
// advect(timestep, 3, 3);
|
||||||
project();
|
// project(data[1],data[2],data[3],oldData[1],oldData[2]);
|
||||||
//diffuse density
|
|
||||||
moveData(false, 0);
|
|
||||||
diffuse(timestep, 0, 0);
|
// //diffuse density
|
||||||
//advect density
|
// moveData(oldData[0], data[0]);
|
||||||
moveData(false, 0);
|
// diffuse(timestep, data[0], oldData[0], 0);
|
||||||
advect(timestep, 0, 0);
|
// //advect density
|
||||||
|
// moveData(oldData[0], data[0]);
|
||||||
|
// advect(timestep, 0, 0);
|
||||||
|
|
||||||
|
if(step % 1500 == 0){
|
||||||
|
for(int i = 0; i < 4; i++){
|
||||||
|
for(int j = 0; j < 4; j++){
|
||||||
|
x[15 - i][15][15 - j] = 1;
|
||||||
|
v[15 - i][15][15 - j] = -1;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
// x[16][16][16] = 1;
|
||||||
|
|
||||||
|
vel_step(DIM, DIM, DIM, u, v, w, u0, v0, w0, VISCOSITY_CONSTANT, timestep);
|
||||||
|
|
||||||
|
dens_step(DIM, DIM, DIM, x, x0, u, v, w, DIFFUSION_CONSTANT, timestep);
|
||||||
|
|
||||||
sum();
|
sum();
|
||||||
}
|
}
|
||||||
|
|
||||||
/**
|
|
||||||
* Diffuses liquid particles
|
void lin_solve ( int M, int N, int O, int b, float[][][] x, float[][][] x0, float a, float c ){
|
||||||
* @param timestep
|
int i, j, k, l;
|
||||||
*/
|
|
||||||
void diffuse(float timestep, int map, int b){
|
// iterate the solver
|
||||||
float diff = DIFFUSION_CONSTANT;
|
for ( l=0 ; l<LINEARSOLVERTIMES ; l++ ) {
|
||||||
float a = timestep * diff * DIM * DIM * DIM;
|
// update for each cell
|
||||||
int approximationIterations = 20;
|
for ( i=1 ; i<M-1 ; i++ ) { for ( j=1 ; j<N-1 ; j++ ) { for ( k=1 ; k<O-1 ; k++ ) {
|
||||||
for(int i = 0; i < approximationIterations; i++){
|
x[i][j][k] = (x0[i][j][k] + a*(x[i-1][j][k]+x[i+1][j][k]+x[i][j-1][k]+x[i][j+1][k]+x[i][j][k-1]+x[i][j][k+1]))/c;
|
||||||
//TODO: parallelize w/ vector api
|
}}}
|
||||||
for(int x = 1; x < DIM-1; x++){
|
setBounds(x, b);
|
||||||
for(int y = 1; y < DIM-1; y++){
|
|
||||||
for(int z = 1; z < DIM-1; z++){
|
|
||||||
data[map][x][y][z] = (oldData[map][x][y][z] +
|
|
||||||
a * (
|
|
||||||
data[map][x-1][y+0][z+0] +
|
|
||||||
data[map][x+1][y+0][z+0] +
|
|
||||||
data[map][x+0][y-1][z+0] +
|
|
||||||
data[map][x+0][y+1][z+0] +
|
|
||||||
data[map][x+0][y+0][z-1] +
|
|
||||||
data[map][x+0][y+0][z+1]
|
|
||||||
))/ (1.0f + 6.0f * a);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
setBounds(data[map],b);
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
|
void diffuse(int M, int N, int O, int b, float[][][] x, float[][][] x0, float diff, float dt){
|
||||||
|
int max = Math.max(Math.max(M, N), Math.max(N, O));
|
||||||
|
float a=dt*diff*max*max*max;
|
||||||
|
lin_solve ( M, N, O, b, x, x0, a, 1+6*a );
|
||||||
|
}
|
||||||
|
|
||||||
void setBounds(float[][][] target, int b){
|
void setBounds(float[][][] target, int b){
|
||||||
for(int x=1; x < DIM-1; x++){
|
for(int x=1; x < DIM-1; x++){
|
||||||
for(int y = 1; y < DIM-1; y++){
|
for(int y = 1; y < DIM-1; y++){
|
||||||
@ -137,161 +155,155 @@ public class FluidSim {
|
|||||||
target[x][y][DIM-1] = b==3 ? -target[x][y][DIM-2] : target[x][y][DIM-2];
|
target[x][y][DIM-1] = b==3 ? -target[x][y][DIM-2] : target[x][y][DIM-2];
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
target[0][0][0] = (float)(0.333*(target[1][0][0]+target[0][1][0]+target[0][0][1]));
|
|
||||||
target[DIM-1][0][0] = (float)(0.333*(target[DIM-2][0][0]+target[DIM-1][1][0]+target[DIM-1][0][1]));
|
|
||||||
target[0][DIM-1][0] = (float)(0.333*(target[1][DIM-1][0]+target[0][DIM-2][0]+target[0][DIM-1][1]));
|
|
||||||
target[0][0][DIM-1] = (float)(0.333*(target[0][0][DIM-2]+target[1][0][DIM-1]+target[0][1][DIM-1]));
|
|
||||||
target[DIM-1][DIM-1][0] = (float)(0.333*(target[DIM-2][DIM-1][0]+target[DIM-1][DIM-2][0]+target[DIM-1][DIM-1][1]));
|
|
||||||
target[0][DIM-1][DIM-1] = (float)(0.333*(target[1][DIM-1][DIM-1]+target[0][DIM-2][DIM-1]+target[0][DIM-1][DIM-2]));
|
|
||||||
target[DIM-1][0][DIM-1] = (float)(0.333*(target[DIM-1][0][DIM-2]+target[DIM-2][0][DIM-1]+target[DIM-1][1][DIM-1]));
|
|
||||||
target[DIM-1][DIM-1][DIM-1] = (float)(0.333*(target[DIM-1][DIM-1][DIM-2]+target[DIM-1][DIM-2][DIM-1]+target[DIM-1][DIM-1][DIM-2]));
|
|
||||||
}
|
|
||||||
|
|
||||||
void advect(float timestep, int map, int b){
|
|
||||||
float dt0 = timestep * DIM;
|
|
||||||
for(int x = 1; x < DIM-1; x++){
|
for(int x = 1; x < DIM-1; x++){
|
||||||
for(int y = 1; y < DIM-1; y++){
|
target[x][0][0] = (float)(0.5f * (target[x][1][0] + target[x][0][1]));
|
||||||
for(int z = 1; z < DIM-1; z++){
|
target[x][DIM-1][0] = (float)(0.5f * (target[x][DIM-2][0] + target[x][DIM-1][1]));
|
||||||
float vx = x * dt0 * oldData[1][x][y][z];
|
target[x][0][DIM-1] = (float)(0.5f * (target[x][1][DIM-1] + target[x][0][DIM-2]));
|
||||||
float vy = y * dt0 * oldData[2][x][y][z];
|
target[x][DIM-1][DIM-1] = (float)(0.5f * (target[x][DIM-2][DIM-1] + target[x][DIM-1][DIM-2]));
|
||||||
float vz = z * dt0 * oldData[3][x][y][z];
|
|
||||||
if(vx < 0.5f){
|
target[0][x][0] = (float)(0.5f * (target[1][x][0] + target[0][x][1]));
|
||||||
vx = 0.5f;
|
target[DIM-1][x][0] = (float)(0.5f * (target[DIM-2][x][0] + target[DIM-1][x][1]));
|
||||||
}
|
target[0][x][DIM-1] = (float)(0.5f * (target[1][x][DIM-1] + target[0][x][DIM-2]));
|
||||||
if(vx > DIM + 0.5f){
|
target[DIM-1][x][DIM-1] = (float)(0.5f * (target[DIM-2][x][DIM-1] + target[DIM-1][x][DIM-2]));
|
||||||
vx = DIM + 0.5f;
|
|
||||||
}
|
|
||||||
if(vy < 0.5f){
|
target[0][0][x] = (float)(0.5f * (target[1][0][x] + target[0][1][x]));
|
||||||
vy = 0.5f;
|
target[DIM-1][0][x] = (float)(0.5f * (target[DIM-2][0][x] + target[DIM-1][1][x]));
|
||||||
}
|
target[0][DIM-1][x] = (float)(0.5f * (target[1][DIM-1][x] + target[0][DIM-2][x]));
|
||||||
if(vy > DIM + 0.5f){
|
target[DIM-1][DIM-1][x] = (float)(0.5f * (target[DIM-2][DIM-1][x] + target[DIM-1][DIM-2][x]));
|
||||||
vy = DIM + 0.5f;
|
|
||||||
}
|
|
||||||
if(vz < 0.5f){
|
|
||||||
vz = 0.5f;
|
|
||||||
}
|
|
||||||
if(vz > DIM + 0.5f){
|
|
||||||
vz = DIM + 0.5f;
|
|
||||||
}
|
|
||||||
int i0 = (int)x;
|
|
||||||
int i1 = (int)x+1;
|
|
||||||
int j0 = (int)y;
|
|
||||||
int j1 = (int)y+1;
|
|
||||||
int k0 = (int)x;
|
|
||||||
int k1 = (int)x+1;
|
|
||||||
float r1 = vx - i0;
|
|
||||||
float r0 = 1 - r1;
|
|
||||||
float s1 = vy - j0;
|
|
||||||
float s0 = 1 - s1;
|
|
||||||
float t1 = vz - k0;
|
|
||||||
float t0 = 1 - t1;
|
|
||||||
|
|
||||||
data[map][x][y][z] =
|
|
||||||
t0 * (
|
|
||||||
r0*(s0*oldData[map][i0][j0][k0] + s1*oldData[map][i0][j1][k0])+
|
|
||||||
r1*(s0*oldData[map][i1][j0][k0] + s1*oldData[map][i1][j1][k0])
|
|
||||||
) +
|
|
||||||
t1 * (
|
|
||||||
r0*(s0*oldData[map][i0][j0][k1] + s1*oldData[map][i0][j1][k1])+
|
|
||||||
r1*(s0*oldData[map][i1][j0][k1] + s1*oldData[map][i1][j1][k1])
|
|
||||||
)
|
|
||||||
;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
}
|
||||||
setBounds(data[map], b);
|
target[0][0][0] = (float)((target[1][0][0]+target[0][1][0]+target[0][0][1])/3.0);
|
||||||
|
target[DIM-1][0][0] = (float)((target[DIM-2][0][0]+target[DIM-1][1][0]+target[DIM-1][0][1])/3.0);
|
||||||
|
target[0][DIM-1][0] = (float)((target[1][DIM-1][0]+target[0][DIM-2][0]+target[0][DIM-1][1])/3.0);
|
||||||
|
target[0][0][DIM-1] = (float)((target[0][0][DIM-2]+target[1][0][DIM-1]+target[0][1][DIM-1])/3.0);
|
||||||
|
target[DIM-1][DIM-1][0] = (float)((target[DIM-2][DIM-1][0]+target[DIM-1][DIM-2][0]+target[DIM-1][DIM-1][1])/3.0);
|
||||||
|
target[0][DIM-1][DIM-1] = (float)((target[1][DIM-1][DIM-1]+target[0][DIM-2][DIM-1]+target[0][DIM-1][DIM-2])/3.0);
|
||||||
|
target[DIM-1][0][DIM-1] = (float)((target[DIM-1][0][DIM-2]+target[DIM-2][0][DIM-1]+target[DIM-1][1][DIM-1])/3.0);
|
||||||
|
target[DIM-1][DIM-1][DIM-1] = (float)((target[DIM-1][DIM-1][DIM-2]+target[DIM-1][DIM-2][DIM-1]+target[DIM-1][DIM-1][DIM-2])/3.0);
|
||||||
}
|
}
|
||||||
|
|
||||||
void moveData(boolean newDataIsSrc, int dim){
|
void advect(int M, int N, int O, int b, float[][][] d, float[][][] d0, float[][][] u, float[][][] v, float[][][] w, float dt){
|
||||||
|
int i, j, k, i0, j0, k0, i1, j1, k1;
|
||||||
|
float x, y, z, s0, t0, s1, t1, u1, u0, dtx,dty,dtz;
|
||||||
|
|
||||||
|
dtx=dty=dtz=dt*Math.max(Math.max(M, N), Math.max(N, O));
|
||||||
|
|
||||||
|
for ( i=1 ; i<M-1 ; i++ ) { for ( j=1 ; j<N-1 ; j++ ) { for ( k=1 ; k<O-1 ; k++ ) {
|
||||||
|
x = i-dtx*u[i][j][k]; y = j-dty*v[i][j][k]; z = k-dtz*w[i][j][k];
|
||||||
|
if (x<0.5f) x=0.5f; if (x>M+0.5f) x=M+0.5f; i0=(int)x; i1=i0+1;
|
||||||
|
if (y<0.5f) y=0.5f; if (y>N+0.5f) y=N+0.5f; j0=(int)y; j1=j0+1;
|
||||||
|
if (z<0.5f) z=0.5f; if (z>O+0.5f) z=O+0.5f; k0=(int)z; k1=k0+1;
|
||||||
|
|
||||||
|
s1 = x-i0; s0 = 1-s1; t1 = y-j0; t0 = 1-t1; u1 = z-k0; u0 = 1-u1;
|
||||||
|
if(i0 >= DIM){
|
||||||
|
i0 = DIM - 1;
|
||||||
|
}
|
||||||
|
if(j0 >= DIM){
|
||||||
|
j0 = DIM - 1;
|
||||||
|
}
|
||||||
|
if(k0 >= DIM){
|
||||||
|
k0 = DIM - 1;
|
||||||
|
}
|
||||||
|
if(i1 >= DIM){
|
||||||
|
i1 = DIM - 1;
|
||||||
|
}
|
||||||
|
if(j1 >= DIM){
|
||||||
|
j1 = DIM - 1;
|
||||||
|
}
|
||||||
|
if(k1 >= DIM){
|
||||||
|
k1 = DIM - 1;
|
||||||
|
}
|
||||||
|
d[i][j][k] = s0*(t0*u0*d0[i0][j0][k0]+t1*u0*d0[i0][j1][k0]+t0*u1*d0[i0][j0][k1]+t1*u1*d0[i0][j1][k1])+
|
||||||
|
s1*(t0*u0*d0[i1][j0][k0]+t1*u0*d0[i1][j1][k0]+t0*u1*d0[i1][j0][k1]+t1*u1*d0[i1][j1][k1]);
|
||||||
|
}}}
|
||||||
|
|
||||||
|
setBounds(d, b);
|
||||||
|
}
|
||||||
|
|
||||||
|
void moveData(float[][][] set1, float[][][] set2){
|
||||||
float tmp = 0;
|
float tmp = 0;
|
||||||
for(int x = 0; x < DIM; x++){
|
for(int x = 0; x < DIM; x++){
|
||||||
for(int y = 0; y < DIM; y++){
|
for(int y = 0; y < DIM; y++){
|
||||||
for(int z = 0; z < DIM; z++){
|
for(int z = 0; z < DIM; z++){
|
||||||
if(newDataIsSrc){
|
tmp = set1[x][y][z];
|
||||||
tmp = data[dim][x][y][z];
|
set1[x][y][z] = set2[x][y][z];
|
||||||
data[dim][x][y][z] = oldData[dim][x][y][z];
|
set2[x][y][z] = tmp;
|
||||||
oldData[dim][x][y][z] = tmp;
|
|
||||||
} else {
|
|
||||||
tmp = oldData[dim][x][y][z];
|
|
||||||
oldData[dim][x][y][z] = data[dim][x][y][z];
|
|
||||||
data[dim][x][y][z] = tmp;
|
|
||||||
}
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
void project(){
|
void project(int M, int N, int O, float[][][] u, float[][][] v, float[][][] w, float[][][] p, float[][][] div){
|
||||||
float h = 1.0f/DIM;
|
int i, j, k;
|
||||||
float[][][] p = oldData[1];
|
|
||||||
float[][][] div = oldData[2];
|
for ( i=1 ; i<M-1 ; i++ ) { for ( j=1 ; j<N-1 ; j++ ) { for ( k=1 ; k<O-1 ; k++ ) {
|
||||||
for(int x=1; x < DIM - 1; x++ ){
|
div[i][j][k] = (float)(-1.0/3.0*((u[i+1][j][k]-u[i-1][j][k])/M+(v[i][j+1][k]-v[i][j-1][k])/M+(w[i][j][k+1]-w[i][j][k-1])/M));
|
||||||
for(int y = 1 ; y < DIM - 1; y++ ){
|
p[i][j][k] = 0;
|
||||||
for(int z = 1; z < DIM - 1; z++){
|
}}}
|
||||||
div[x][y][z] = (float)(-0.5*h*(
|
|
||||||
data[1][x+1][y][z]-
|
|
||||||
data[1][x-1][y][z]+
|
|
||||||
data[2][x][y+1][z]-
|
|
||||||
data[2][x][y-1][z]+
|
|
||||||
data[3][x][y][z+1]-
|
|
||||||
data[3][x][y][z-1]
|
|
||||||
));
|
|
||||||
p[x][y][z] = 0;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
setBounds(div, 0);
|
setBounds(div, 0);
|
||||||
setBounds(p, 0);
|
setBounds(p, 0);
|
||||||
//TODO: parallelize with vectorized jacobi
|
|
||||||
for(int k=0; k<20; k++){
|
lin_solve ( M, N, O, 0, p, div, 1, 6 );
|
||||||
for(int x=1; x < DIM - 1; x++ ){
|
|
||||||
for(int y = 1 ; y < DIM - 1; y++ ){
|
for ( i=1 ; i<M-1 ; i++ ) { for ( j=1 ; j<N-1 ; j++ ) { for ( k=1 ; k<O-1 ; k++ ) {
|
||||||
for(int z = 1; z < DIM - 1; z++){
|
u[i][j][k] -= 0.5f*M*(p[i+1][j][k]-p[i-1][j][k]);
|
||||||
p[x][y][z] = (
|
v[i][j][k] -= 0.5f*M*(p[i][j+1][k]-p[i][j-1][k]);
|
||||||
div[x][y][z]+
|
w[i][j][k] -= 0.5f*M*(p[i][j][k+1]-p[i][j][k-1]);
|
||||||
p[x-1][y][z]+
|
}}}
|
||||||
p[x+1][y][z]+
|
|
||||||
p[x][y-1][z]+
|
setBounds(u, 1);
|
||||||
p[x][y+1][z]+
|
setBounds(v, 2);
|
||||||
p[x][y][z-1]+
|
setBounds(w, 3);
|
||||||
p[x][y][z+1]
|
}
|
||||||
)/6;
|
|
||||||
}
|
void dens_step(int M, int N, int O, float[][][] x, float[][][] x0, float[][][] u, float[][][] v, float[][][] w, float diff, float dt){
|
||||||
}
|
// add_source ( M, N, O, x, x0, dt );
|
||||||
|
moveData ( x0, x ); diffuse ( M, N, O, 0, x, x0, diff, dt );
|
||||||
|
moveData ( x0, x ); advect ( M, N, O, 0, x, x0, u, v, w, dt );
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
|
void add_source(int M, int N, int O, float[][][] x, float[][][] density, float dt){
|
||||||
|
int i, j, k;
|
||||||
|
for ( i=1 ; i<M-1 ; i++ ) { for ( j=1 ; j<N-1 ; j++ ) { for ( k=1 ; k<O-1 ; k++ ) {
|
||||||
|
if(density[i][j][k] > 0.0f){
|
||||||
|
x[i][j][k] += dt * GRAVITY * density[i][j][k];
|
||||||
}
|
}
|
||||||
setBounds(p, 0);
|
}}}
|
||||||
}
|
}
|
||||||
for(int x=1; x < DIM - 1; x++ ){
|
|
||||||
for(int y = 1 ; y < DIM - 1; y++ ){
|
void vel_step(int M, int N, int O, float[][][] u, float[][][] v, float[][][] w, float[][][] u0, float[][][] v0, float[][][] w0, float visc, float dt ){
|
||||||
for(int z = 1; z < DIM - 1; z++){
|
// add_source ( M, N, O, u, u0, dt ); add_source ( M, N, O, v, v0, dt );add_source ( M, N, O, w, w0, dt );
|
||||||
data[1][x][y][z] -= 0.5*(p[x+1][y][z]-p[x-1][y][z])/h;
|
add_source ( M, N, O, v, x, dt );
|
||||||
data[2][x][y][z] -= 0.5*(p[x][y+1][z]-p[x][y-1][z])/h;
|
moveData ( u0, u ); diffuse ( M, N, O, 1, u, u0, visc, dt );
|
||||||
data[3][x][y][z] -= 0.5*(p[x][y][z+1]-p[x][y][z-1])/h;
|
moveData ( v0, v ); diffuse ( M, N, O, 2, v, v0, visc, dt );
|
||||||
}
|
moveData ( w0, w ); diffuse ( M, N, O, 3, w, w0, visc, dt );
|
||||||
}
|
project ( M, N, O, u, v, w, u0, v0 );
|
||||||
}
|
moveData ( u0, u ); moveData ( v0, v );moveData ( w0, w );
|
||||||
setBounds(data[1], 1);
|
advect ( M, N, O, 1, u, u0, u0, v0, w0, dt ); advect ( M, N, O, 2, v, v0, u0, v0, w0, dt );advect ( M, N, O, 3, w, w0, u0, v0, w0, dt );
|
||||||
setBounds(data[2], 2);
|
project ( M, N, O, u, v, w, u0, v0 );
|
||||||
setBounds(data[3], 3);
|
|
||||||
}
|
}
|
||||||
|
|
||||||
void sum(){
|
void sum(){
|
||||||
double[] sum = new double[8];
|
double[] sum = new double[8];
|
||||||
for(int x = 0; x < DIM; x++){
|
for(int i = 0; i < DIM; i++){
|
||||||
for(int y = 0; y < DIM; y++){
|
for(int j = 0; j < DIM; j++){
|
||||||
for(int z = 0; z < DIM; z++){
|
for(int k = 0; k < DIM; k++){
|
||||||
sum[0] = sum[0] + data[0][x][y][z];
|
sum[0] = sum[0] + x[i][j][k];
|
||||||
sum[1] = sum[1] + data[1][x][y][z];
|
sum[1] = sum[1] + u[i][j][k];
|
||||||
sum[2] = sum[2] + data[2][x][y][z];
|
sum[2] = sum[2] + v[i][j][k];
|
||||||
sum[3] = sum[3] + data[3][x][y][z];
|
sum[3] = sum[3] + w[i][j][k];
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
System.out.println(sum[0] + " " + sum[1] + " " + sum[2] + " " + sum[3] );
|
// System.out.println(sum[0] + " " + sum[1] + " " + sum[2] + " " + sum[3]);
|
||||||
}
|
}
|
||||||
|
|
||||||
public float[][][][] getData(){
|
public float[][][] getData(){
|
||||||
return data;
|
return x;
|
||||||
}
|
}
|
||||||
|
|
||||||
}
|
}
|
||||||
|
|||||||
@ -29,15 +29,17 @@ public class Main {
|
|||||||
FluidSim sim = new FluidSim();
|
FluidSim sim = new FluidSim();
|
||||||
sim.setup();
|
sim.setup();
|
||||||
Mesh.meshInitially(sim);
|
Mesh.meshInitially(sim);
|
||||||
|
int i = 0;
|
||||||
while(true){
|
while(true){
|
||||||
try {
|
try {
|
||||||
TimeUnit.MILLISECONDS.sleep(5);
|
TimeUnit.MILLISECONDS.sleep(1);
|
||||||
} catch (InterruptedException e) {
|
} catch (InterruptedException e) {
|
||||||
e.printStackTrace();
|
e.printStackTrace();
|
||||||
}
|
}
|
||||||
sim.simulate(0.0001f);
|
sim.simulate(i,0.001f);
|
||||||
Mesh.remesh(sim);
|
Mesh.remesh(sim);
|
||||||
GLFWContext.redraw();
|
GLFWContext.redraw();
|
||||||
|
i++;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|||||||
@ -70,7 +70,7 @@ public class MyTest {
|
|||||||
int iterations = 0;
|
int iterations = 0;
|
||||||
while(true){
|
while(true){
|
||||||
|
|
||||||
sim.simulate(0.08f);
|
sim.simulate(1, 0.2f);
|
||||||
|
|
||||||
//framerate calculations
|
//framerate calculations
|
||||||
double currentTime = glfwGetTime();
|
double currentTime = glfwGetTime();
|
||||||
@ -79,7 +79,7 @@ public class MyTest {
|
|||||||
lastFrame = currentTime;
|
lastFrame = currentTime;
|
||||||
|
|
||||||
iterations++;
|
iterations++;
|
||||||
System.out.print("\r" + currentTime + " " + frameRate + " _ " + iterations + " _ " + Runtime.getRuntime().totalMemory() / 1024 / 1024);
|
// System.out.print("\r" + currentTime + " " + frameRate + " _ " + iterations + " _ " + Runtime.getRuntime().totalMemory() / 1024 / 1024);
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
|
|||||||
@ -188,7 +188,7 @@ public class Mesh {
|
|||||||
GL45.glUseProgram(shaderProgram);
|
GL45.glUseProgram(shaderProgram);
|
||||||
|
|
||||||
viewMatrix = new Matrix4f();
|
viewMatrix = new Matrix4f();
|
||||||
viewMatrix.setLookAt(new Vector3f(-10,10,0), new Vector3f(8,8,8), new Vector3f(0,1,0));
|
viewMatrix.setLookAt(new Vector3f(-8,20,0), new Vector3f(8,8,8), new Vector3f(0,1,0));
|
||||||
projectionMatrix = new Matrix4f();
|
projectionMatrix = new Matrix4f();
|
||||||
projectionMatrix.setPerspective(90,1920.0f/1080.0f,0.0001f,1000f);
|
projectionMatrix.setPerspective(90,1920.0f/1080.0f,0.0001f,1000f);
|
||||||
model = new Matrix4f().identity();
|
model = new Matrix4f().identity();
|
||||||
@ -820,7 +820,7 @@ public class Mesh {
|
|||||||
return new Vector3f(x,y,z);
|
return new Vector3f(x,y,z);
|
||||||
}
|
}
|
||||||
|
|
||||||
public static TerrainChunkData generateTerrainChunkData(float[][][][] data){
|
public static TerrainChunkData generateTerrainChunkData(float[][][] data){
|
||||||
|
|
||||||
// 5 6
|
// 5 6
|
||||||
// +-------------+ +-----5-------+ ^ Y
|
// +-------------+ +-----5-------+ ^ Y
|
||||||
@ -849,15 +849,15 @@ public class Mesh {
|
|||||||
|
|
||||||
|
|
||||||
|
|
||||||
for(int x = 0; x < data[0].length - 1; x++){
|
for(int x = 0; x < data.length - 1; x++){
|
||||||
for(int y = 0; y < data[0][0].length - 1; y++){
|
for(int y = 0; y < data[0].length - 1; y++){
|
||||||
for(int z = 0; z < data[0][0][0].length - 1; z++){
|
for(int z = 0; z < data[0][0].length - 1; z++){
|
||||||
//push the current cell's values into the gridcell
|
//push the current cell's values into the gridcell
|
||||||
currentCell.setValues(
|
currentCell.setValues(
|
||||||
new Vector3f(x+0,y+0,z+0), new Vector3f(x+0,y+0,z+1), new Vector3f(x+1,y+0,z+1), new Vector3f(x+1,y+0,z+0),
|
new Vector3f(x+0,y+0,z+0), new Vector3f(x+0,y+0,z+1), new Vector3f(x+1,y+0,z+1), new Vector3f(x+1,y+0,z+0),
|
||||||
new Vector3f(x+0,y+1,z+0), new Vector3f(x+0,y+1,z+1), new Vector3f(x+1,y+1,z+1), new Vector3f(x+1,y+1,z+0),
|
new Vector3f(x+0,y+1,z+0), new Vector3f(x+0,y+1,z+1), new Vector3f(x+1,y+1,z+1), new Vector3f(x+1,y+1,z+0),
|
||||||
data[0][x+0][y+0][z+0], data[0][x+0][y+0][z+1], data[0][x+1][y+0][z+1], data[0][x+1][y+0][z+0],
|
data[x+0][y+0][z+0], data[x+0][y+0][z+1], data[x+1][y+0][z+1], data[x+1][y+0][z+0],
|
||||||
data[0][x+0][y+1][z+0], data[0][x+0][y+1][z+1], data[0][x+1][y+1][z+1], data[0][x+1][y+1][z+0]
|
data[x+0][y+1][z+0], data[x+0][y+1][z+1], data[x+1][y+1][z+1], data[x+1][y+1][z+0]
|
||||||
);
|
);
|
||||||
//polygonize the current gridcell
|
//polygonize the current gridcell
|
||||||
polygonize(currentCell, 0.01, triangles, vertMap, verts, normals, trianglesSharingVert);
|
polygonize(currentCell, 0.01, triangles, vertMap, verts, normals, trianglesSharingVert);
|
||||||
|
|||||||
@ -5,38 +5,10 @@
|
|||||||
out vec4 FragColor;
|
out vec4 FragColor;
|
||||||
|
|
||||||
|
|
||||||
layout (std140) uniform Lights {
|
vec3 color = vec3(0.3,0.7,0.9);
|
||||||
// this is how many because we have to align
|
|
||||||
// bytes it SHOULD in multiples of 16, this
|
|
||||||
// take it where it ACTUALLY is
|
|
||||||
//
|
|
||||||
//refer: https://learnopengl.com/Advanced-OpenGL/Advanced-GLSL
|
|
||||||
//
|
|
||||||
// base alignment aligned offset
|
|
||||||
//direct light
|
|
||||||
vec3 dLDirection; // 16 0
|
|
||||||
vec3 dLAmbient; // 16 16
|
|
||||||
vec3 dLDiffuse; // 16 32
|
|
||||||
vec3 dLSpecular; // 16 48
|
|
||||||
|
|
||||||
//point light
|
vec3 dLDirection = vec3(0.1,-0.9,0.1);
|
||||||
vec3 pLposition[NR_POINT_LIGHTS]; // 16*10 64
|
vec3 dLDiffuse = vec3(0.5,0.5,0.5);
|
||||||
float pLconstant[NR_POINT_LIGHTS]; // 16*10 224
|
|
||||||
float pLlinear[NR_POINT_LIGHTS]; // 16*10 384
|
|
||||||
float pLquadratic[NR_POINT_LIGHTS]; // 16*10 544
|
|
||||||
vec3 pLambient[NR_POINT_LIGHTS]; // 16*10 704
|
|
||||||
vec3 pLdiffuse[NR_POINT_LIGHTS]; // 16*10 864
|
|
||||||
vec3 pLspecular[NR_POINT_LIGHTS]; // 16*10 1024
|
|
||||||
|
|
||||||
//for a total size of 1184
|
|
||||||
|
|
||||||
};
|
|
||||||
|
|
||||||
struct Material {
|
|
||||||
sampler2D diffuse;
|
|
||||||
sampler2D specular;
|
|
||||||
float shininess;
|
|
||||||
};
|
|
||||||
|
|
||||||
in vec3 FragPos;
|
in vec3 FragPos;
|
||||||
in vec3 Normal;
|
in vec3 Normal;
|
||||||
@ -46,23 +18,9 @@ uniform vec3 viewPos;
|
|||||||
// uniform DirLight dirLight;
|
// uniform DirLight dirLight;
|
||||||
// uniform PointLight pointLights[NR_POINT_LIGHTS];
|
// uniform PointLight pointLights[NR_POINT_LIGHTS];
|
||||||
// uniform SpotLight spotLight;
|
// uniform SpotLight spotLight;
|
||||||
uniform Material material;
|
|
||||||
|
|
||||||
//texture stuff
|
|
||||||
// uniform sampler2D ourTexture;
|
|
||||||
uniform int hasTransparency;
|
|
||||||
// uniform sampler2D specularTexture;
|
|
||||||
|
|
||||||
//light depth map
|
|
||||||
uniform sampler2D shadowMap;
|
|
||||||
|
|
||||||
|
|
||||||
// function prototypes
|
// function prototypes
|
||||||
// vec3 CalcDirLight(vec3 normal, vec3 viewDir);
|
vec3 CalcDirLight(vec3 normal, vec3 viewDir);
|
||||||
// vec3 CalcPointLight(int i, vec3 normal, vec3 fragPos, vec3 viewDir);
|
|
||||||
// vec3 CalcSpotLight(vec3 normal, vec3 fragPos, vec3 viewDir);
|
|
||||||
// float calcLightIntensityTotal(vec3 normal);
|
|
||||||
// float ShadowCalculation(vec4 fragPosLightSpace, vec3 lightDir, vec3 normal);
|
|
||||||
|
|
||||||
void main(){
|
void main(){
|
||||||
vec3 norm = normalize(Normal);
|
vec3 norm = normalize(Normal);
|
||||||
@ -77,125 +35,29 @@ void main(){
|
|||||||
|
|
||||||
// //calculate final color
|
// //calculate final color
|
||||||
// vec3 finalColor = textureColor * lightIntensity;
|
// vec3 finalColor = textureColor * lightIntensity;
|
||||||
// vec3 lightAmount = CalcDirLight(norm, viewDir);
|
vec3 lightAmount = CalcDirLight(norm, viewDir);
|
||||||
// for(int i = 0; i < NR_POINT_LIGHTS; i++){
|
// for(int i = 0; i < NR_POINT_LIGHTS; i++){
|
||||||
// lightAmount += CalcPointLight(i, norm, FragPos, viewDir);
|
// lightAmount += CalcPointLight(i, norm, FragPos, viewDir);
|
||||||
// }
|
// }
|
||||||
|
|
||||||
|
vec3 color = vec3(0.3,0.7,0.9) * lightAmount;
|
||||||
|
|
||||||
//this final calculation is for transparency
|
//this final calculation is for transparency
|
||||||
FragColor = vec4(0.3,0.7,0.9,1.0);
|
FragColor = vec4(color,1.0);
|
||||||
}
|
}
|
||||||
|
|
||||||
|
vec3 CalcDirLight(vec3 normal, vec3 viewDir){
|
||||||
|
vec3 lightDir = normalize(-dLDirection);
|
||||||
|
// diffuse shading
|
||||||
|
float diff = max(dot(normal, lightDir), 0.0);
|
||||||
|
// specular shading
|
||||||
|
vec3 reflectDir = reflect(-lightDir, normal);
|
||||||
|
float spec = pow(max(dot(viewDir, reflectDir), 0.0), 0.6);
|
||||||
|
// combine results
|
||||||
|
// vec3 texColor = texture(material.diffuse, TexCoord).rgb;
|
||||||
|
vec3 diffuse = dLDiffuse * diff;
|
||||||
|
|
||||||
// vec3 getColor(vec2 texPlane1, vec2 texPlane2, vec2 texPlane3, vec3 normal, Material material){
|
vec3 specular = spec * color;
|
||||||
|
|
||||||
// vec3 weights = abs(normal);
|
|
||||||
|
|
||||||
// vec3 albedoX = texture(material.diffuse, texPlane1).rgb;
|
|
||||||
// vec3 albedoY = texture(material.diffuse, texPlane2).rgb;
|
|
||||||
// vec3 albedoZ = texture(material.diffuse, texPlane3).rgb;
|
|
||||||
|
|
||||||
|
|
||||||
// return (albedoX * weights.x + albedoY * weights.y + albedoZ * weights.z);
|
|
||||||
// }
|
|
||||||
|
|
||||||
// //
|
|
||||||
// float calcLightIntensityAmbient(){
|
|
||||||
// //calculate average of ambient light
|
|
||||||
// float avg = (dLAmbient.x + dLAmbient.y + dLAmbient.z)/3.0;
|
|
||||||
// return avg;
|
|
||||||
// }
|
|
||||||
|
|
||||||
// //
|
|
||||||
// float calcLightIntensityDir(vec3 normal){
|
|
||||||
// vec3 lightDir = normalize(-dLDirection);
|
|
||||||
// // diffuse shading
|
|
||||||
// float diff = max(dot(normal, lightDir), 0.0);
|
|
||||||
|
|
||||||
// return diff;
|
return diffuse + specular;
|
||||||
// }
|
}
|
||||||
|
|
||||||
// //
|
|
||||||
// float calcLightIntensityTotal(vec3 normal){
|
|
||||||
// //ambient intensity
|
|
||||||
// float ambientLightIntensity = calcLightIntensityAmbient();
|
|
||||||
|
|
||||||
// //get direct intensity
|
|
||||||
// float directLightIntensity = calcLightIntensityDir(normal);
|
|
||||||
|
|
||||||
// //sum
|
|
||||||
// float total = ambientLightIntensity + directLightIntensity;
|
|
||||||
// return total;
|
|
||||||
// }
|
|
||||||
|
|
||||||
// //
|
|
||||||
// vec3 getTotalLightColor(vec3 normal){
|
|
||||||
// //get the direct light color adjusted for intensity
|
|
||||||
// vec3 diffuseLightColor = dLDiffuse * calcLightIntensityDir(normal);
|
|
||||||
|
|
||||||
// //sum light colors
|
|
||||||
// vec3 totalLightColor = diffuseLightColor;
|
|
||||||
// return totalLightColor;
|
|
||||||
// }
|
|
||||||
|
|
||||||
// vec3 CalcPointLight(int i, vec3 normal, vec3 fragPos, vec3 viewDir){
|
|
||||||
// vec3 lightDir = normalize(pLposition[i] - fragPos);
|
|
||||||
// // diffuse shading
|
|
||||||
// float diff = max(dot(normal, lightDir), 0.0);
|
|
||||||
// // specular shading
|
|
||||||
// // vec3 reflectDir = reflect(-lightDir, normal);
|
|
||||||
// // float spec = pow(max(dot(viewDir, reflectDir), 0.0), material.shininess);
|
|
||||||
// // attenuation
|
|
||||||
// float distance = length(pLposition[i] - fragPos);
|
|
||||||
// float attenuation = 1.0 / (pLconstant[i] + pLlinear[i] * distance + pLquadratic[i] * (distance * distance));
|
|
||||||
// // combine results
|
|
||||||
// vec3 ambient = pLambient[i];
|
|
||||||
// vec3 diffuse = pLdiffuse[i] * diff;
|
|
||||||
// ambient *= attenuation;
|
|
||||||
// diffuse *= attenuation;
|
|
||||||
// // specular *= attenuation;
|
|
||||||
// vec3 specular = vec3(0,0,0);
|
|
||||||
|
|
||||||
// vec3 finalValue = (ambient + diffuse + specular);
|
|
||||||
// finalValue = vec3(max(finalValue.x,0),max(finalValue.y,0),max(finalValue.z,0));
|
|
||||||
|
|
||||||
// return finalValue;
|
|
||||||
// }
|
|
||||||
|
|
||||||
|
|
||||||
// float ShadowCalculation(vec4 fragPosLightSpace, vec3 lightDir, vec3 normal){
|
|
||||||
|
|
||||||
// // perform perspective divide
|
|
||||||
// vec3 projCoords = fragPosLightSpace.xyz / fragPosLightSpace.w;
|
|
||||||
|
|
||||||
// //transform to NDC
|
|
||||||
// projCoords = projCoords * 0.5 + 0.5;
|
|
||||||
|
|
||||||
// //get closest depth from light's POV
|
|
||||||
// float closestDepth = texture(shadowMap, projCoords.xy).r;
|
|
||||||
|
|
||||||
// //get depth of current fragment
|
|
||||||
// float currentDepth = projCoords.z;
|
|
||||||
|
|
||||||
// //calculate bias
|
|
||||||
// float bias = max(0.05 * (1.0 - dot(normal, lightDir)), 0.005);
|
|
||||||
|
|
||||||
// //calculate shadow value
|
|
||||||
// float shadow = currentDepth - bias > closestDepth ? 1.0 : 0.0;
|
|
||||||
|
|
||||||
// if(projCoords.z > 1.0){
|
|
||||||
// shadow = 0.0;
|
|
||||||
// }
|
|
||||||
|
|
||||||
// //calculate dot product, if it is >0 we know they're parallel-ish therefore should disregard the shadow mapping
|
|
||||||
// //ie the fragment is already facing away from the light source
|
|
||||||
// float dotprod = dot(normalize(lightDir),normalize(normal));
|
|
||||||
|
|
||||||
// if(dotprod > 0.0){
|
|
||||||
// shadow = 0.0;
|
|
||||||
// }
|
|
||||||
|
|
||||||
// // shadow = currentDepth;
|
|
||||||
|
|
||||||
// return shadow;
|
|
||||||
// }
|
|
||||||
|
|||||||
Loading…
Reference in New Issue
Block a user