From 5174572fb869927b79243a69ece32fa884a699a0 Mon Sep 17 00:00:00 2001 From: austin Date: Sat, 30 Nov 2024 14:32:35 -0500 Subject: [PATCH] Migrate C code --- .gitignore | 5 +- .gitmodules | 3 + .vscode/extensions.json | 3 +- README.md | 4 + src/fluid/includes/chunk.h | 19 + src/fluid/includes/chunkmask.h | 48 + src/fluid/includes/electrosphere_FluidSim.h | 31 + src/fluid/includes/mainFunctions.h | 101 ++ src/fluid/includes/simulation.h | 8 + src/fluid/includes/utilities.h | 12 + src/fluid/lib/stb | 1 + src/fluid/src/chunkmask.c | 61 ++ src/fluid/src/densitystep.c | 230 +++++ src/fluid/src/fluidsim.c | 542 ++++++++++ src/fluid/src/javainterface.c | 170 ++++ src/fluid/src/velocitystep.c | 936 ++++++++++++++++++ .../fluid/manager/ServerFluidChunk.java | 214 +++- .../simulator/FluidAcceleratedSimulator.java | 48 + 18 files changed, 2405 insertions(+), 31 deletions(-) create mode 100644 .gitmodules create mode 100644 src/fluid/includes/chunk.h create mode 100644 src/fluid/includes/chunkmask.h create mode 100644 src/fluid/includes/electrosphere_FluidSim.h create mode 100644 src/fluid/includes/mainFunctions.h create mode 100644 src/fluid/includes/simulation.h create mode 100644 src/fluid/includes/utilities.h create mode 160000 src/fluid/lib/stb create mode 100644 src/fluid/src/chunkmask.c create mode 100644 src/fluid/src/densitystep.c create mode 100644 src/fluid/src/fluidsim.c create mode 100644 src/fluid/src/javainterface.c create mode 100644 src/fluid/src/velocitystep.c create mode 100644 src/main/java/electrosphere/server/fluid/simulator/FluidAcceleratedSimulator.java diff --git a/.gitignore b/.gitignore index e3db7342..04525c1c 100644 --- a/.gitignore +++ b/.gitignore @@ -68,4 +68,7 @@ #Memory debug files /heap.dump /heap.dump.hwcache -/tmp \ No newline at end of file +/tmp + +#native libraries +/shared-folder diff --git a/.gitmodules b/.gitmodules new file mode 100644 index 00000000..9ae1a725 --- /dev/null +++ b/.gitmodules @@ -0,0 +1,3 @@ +[submodule "src/fluid/lib/stb"] + path = src/fluid/lib/stb + url = https://github.com/nothings/stb.git diff --git a/.vscode/extensions.json b/.vscode/extensions.json index 81a24224..9471da5e 100644 --- a/.vscode/extensions.json +++ b/.vscode/extensions.json @@ -4,6 +4,7 @@ "vscjava.vscode-java-debug", "slevesque.shader", "vscjava.vscode-java-test", - "dtoplak.vscode-glsllint" + "dtoplak.vscode-glsllint", + "ms-vscode.cpptools-extension-pack" ] } \ No newline at end of file diff --git a/README.md b/README.md index 9c1d0382..6e7547bc 100644 --- a/README.md +++ b/README.md @@ -7,6 +7,10 @@ A multiplayer-focused game engine ## Building + +### Cloning +When cloning the repo, make sure to grab all submodules with `git clone --recurse-submodules git@git.austinwhoover.com:studiorailgun/Renderer.git` + ### Windows 1. Install - [gitbash](https://git-scm.com/downloads) diff --git a/src/fluid/includes/chunk.h b/src/fluid/includes/chunk.h new file mode 100644 index 00000000..586f37ee --- /dev/null +++ b/src/fluid/includes/chunk.h @@ -0,0 +1,19 @@ +#ifndef CHUNK_H +#define CHUNK_H + +/** + * A chunk +*/ +typedef struct { + float * d[27]; + float * d0[27]; + float * u[27]; + float * v[27]; + float * w[27]; + float * u0[27]; + float * v0[27]; + float * w0[27]; + int chunkMask; +} Chunk; + +#endif \ No newline at end of file diff --git a/src/fluid/includes/chunkmask.h b/src/fluid/includes/chunkmask.h new file mode 100644 index 00000000..a0003dc2 --- /dev/null +++ b/src/fluid/includes/chunkmask.h @@ -0,0 +1,48 @@ +#include + +#ifndef CHUNKMASK_H +#define CHUNKMASK_H + +#define CENTER_LOC 13 + +#define CHUNK_222 1 +#define CHUNK_122 2 +#define CHUNK_022 4 +#define CHUNK_212 8 +#define CHUNK_112 16 +#define CHUNK_012 32 +#define CHUNK_202 64 +#define CHUNK_102 128 +#define CHUNK_002 256 + +#define CHUNK_221 512 +#define CHUNK_121 1024 +#define CHUNK_021 2048 +#define CHUNK_211 4096 +#define CHUNK_111 8192 +#define CHUNK_011 16384 +#define CHUNK_201 32768 +#define CHUNK_101 65536 +#define CHUNK_001 131072 + +#define CHUNK_220 262144 +#define CHUNK_120 524288 +#define CHUNK_020 1048576 +#define CHUNK_210 2097152 +#define CHUNK_110 4194304 +#define CHUNK_010 8388608 +#define CHUNK_200 16777216 +#define CHUNK_100 33554432 +#define CHUNK_000 67108864 + +extern const uint32_t CHUNK_INDEX_ARR[]; + + +//control offsetting the advect sampler location if a valid neighbor chunk is hit +extern const char CHUNK_NORMALIZE_U[]; + +extern const char CHUNK_NORMALIZE_V[]; + +extern const char CHUNK_NORMALIZE_W[]; + +#endif \ No newline at end of file diff --git a/src/fluid/includes/electrosphere_FluidSim.h b/src/fluid/includes/electrosphere_FluidSim.h new file mode 100644 index 00000000..84ade9d9 --- /dev/null +++ b/src/fluid/includes/electrosphere_FluidSim.h @@ -0,0 +1,31 @@ +/* DO NOT EDIT THIS FILE - it is machine generated */ +#include +/* Header for class electrosphere_FluidSim */ + +#ifndef _Included_electrosphere_FluidSim +#define _Included_electrosphere_FluidSim +#ifdef __cplusplus +extern "C" { +#endif +#undef electrosphere_FluidSim_DIM +#define electrosphere_FluidSim_DIM 18L +#undef electrosphere_FluidSim_DIFFUSION_CONSTANT +#define electrosphere_FluidSim_DIFFUSION_CONSTANT 1.0E-5f +#undef electrosphere_FluidSim_VISCOSITY_CONSTANT +#define electrosphere_FluidSim_VISCOSITY_CONSTANT 1.0E-5f +#undef electrosphere_FluidSim_LINEARSOLVERTIMES +#define electrosphere_FluidSim_LINEARSOLVERTIMES 20L +#undef electrosphere_FluidSim_GRAVITY +#define electrosphere_FluidSim_GRAVITY -1000.0f +/* + * Class: electrosphere_FluidSim + * Method: simulate + * Signature: (Ljava/util/List;F)V + */ +JNIEXPORT void JNICALL Java_electrosphere_FluidSim_simulate + (JNIEnv *, jclass, jobject, jfloat); + +#ifdef __cplusplus +} +#endif +#endif diff --git a/src/fluid/includes/mainFunctions.h b/src/fluid/includes/mainFunctions.h new file mode 100644 index 00000000..0d6040b3 --- /dev/null +++ b/src/fluid/includes/mainFunctions.h @@ -0,0 +1,101 @@ +#ifndef MAINFUNC +#define MAINFUNC + + +/* + * Class: electrosphere_FluidSim + * Method: addSourceToVectors + * Signature: (II[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;FFF)V + */ +static inline void addSourceToVectors + (int, int, float **, float **, float **, float **, float **, float **, float, float, float); + +/* + * Class: electrosphere_FluidSim + * Method: solveVectorDiffuse + * Signature: (II[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;FFF)V + */ +static inline void solveVectorDiffuse + (int, int, float **, float **, float **, float **, float **, float **, float, float, float); + +/* + * Class: electrosphere_FluidSim + * Method: setupProjection + * Signature: (II[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;FFF)V + */ +static inline void setupProjection + ( + int N, + int chunk_mask, + float ** ur, + float ** vr, + float ** wr, + float ** pr, + float ** divr, + float DIFFUSION_CONST, + float VISCOSITY_CONST, + float dt); + +/* + * Class: electrosphere_FluidSim + * Method: solveProjection + * Signature: (II[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;FFF)V + */ +static inline void solveProjection + (int, int, float **, float **, float **, float **, float **, float **, float, float, float); + +/* + * Class: electrosphere_FluidSim + * Method: finalizeProjection + * Signature: (II[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;FFF)V + */ +static inline void finalizeProjection + (int, int, float **, float **, float **, float **, float **, float **, float, float, float); + +/* + * Class: electrosphere_FluidSim + * Method: advectVectors + * Signature: (II[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;FFF)V + */ +static inline void advectVectors + (int, int, float **, float **, float **, float **, float **, float **, float, float, float); + +/* + * Class: electrosphere_FluidSim + * Method: addDensity + * Signature: (II[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;F)V + */ +static inline void addDensity + (int, int, float **, float **, float); + +/* + * Class: electrosphere_FluidSim + * Method: solveDiffuseDensity + * Signature: (II[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;FFF)V + */ +static inline void solveDiffuseDensity + (int, int, float **, float **, float **, float **, float **, float, float, float); + +/* + * Class: electrosphere_FluidSim + * Method: advectDensity + * Signature: (II[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;[Ljava/nio/ByteBuffer;FFF)V + */ +static inline void advectDensity(uint32_t chunk_mask, int N, float ** d, float ** d0, float ** ur, float ** vr, float ** wr, float dt); + +static inline void setBoundsToNeighborsRaw + ( + int N, + int chunk_mask, + int vector_dir, + float ** neighborArray); + +static inline void copyNeighborsRaw + ( + int N, + int chunk_mask, + int cx, + int vector_dir, + float ** neighborArray); + +#endif \ No newline at end of file diff --git a/src/fluid/includes/simulation.h b/src/fluid/includes/simulation.h new file mode 100644 index 00000000..4f86ca57 --- /dev/null +++ b/src/fluid/includes/simulation.h @@ -0,0 +1,8 @@ +#ifndef SIMULATION_H +#define SIMULATION_H + +#include "./chunk.h" + +void simulate(int numChunks, Chunk ** passedInChunks, float timestep); + +#endif \ No newline at end of file diff --git a/src/fluid/includes/utilities.h b/src/fluid/includes/utilities.h new file mode 100644 index 00000000..ab128a5d --- /dev/null +++ b/src/fluid/includes/utilities.h @@ -0,0 +1,12 @@ +#include + +#ifndef UTILITIES_H +#define UTILITIES_H + +#define SWAP(x0,x) {float *tmp=x0;x0=x;x=tmp;} +#define IX(i,j,k) ((i)+(N)*(j)+(N*N)*(k)) +#define CK(m,n,o) ((m)+(n)*(3)+(o)*(3)*(3)) +#define GET_ARR_RAW(src,i) src[i] +#define ARR_EXISTS(chunk_mask,m,n,o) (chunk_mask & CHUNK_INDEX_ARR[CK(m,n,o)]) > 0 + +#endif \ No newline at end of file diff --git a/src/fluid/lib/stb b/src/fluid/lib/stb new file mode 160000 index 00000000..5c205738 --- /dev/null +++ b/src/fluid/lib/stb @@ -0,0 +1 @@ +Subproject commit 5c205738c191bcb0abc65c4febfa9bd25ff35234 diff --git a/src/fluid/src/chunkmask.c b/src/fluid/src/chunkmask.c new file mode 100644 index 00000000..0b2d0775 --- /dev/null +++ b/src/fluid/src/chunkmask.c @@ -0,0 +1,61 @@ +#include +#include "../includes/utilities.h" +#include "../includes/chunkmask.h" + +const uint32_t CHUNK_INDEX_ARR[] = { + CHUNK_000, CHUNK_100, CHUNK_200, + CHUNK_010, CHUNK_110, CHUNK_210, + CHUNK_020, CHUNK_120, CHUNK_220, + + CHUNK_001, CHUNK_101, CHUNK_201, + CHUNK_011, CHUNK_111, CHUNK_211, + CHUNK_021, CHUNK_121, CHUNK_221, + + CHUNK_002, CHUNK_102, CHUNK_202, + CHUNK_012, CHUNK_112, CHUNK_212, + CHUNK_022, CHUNK_122, CHUNK_222, +}; + + +//control offsetting the advect sampler location if a valid neighbor chunk is hit +const char CHUNK_NORMALIZE_U[] = { + 1, 0, -1, + 1, 0, -1, + 1, 0, -1, + + 1, 0, -1, + 1, 0, -1, + 1, 0, -1, + + 1, 0, -1, + 1, 0, -1, + 1, 0, -1, +}; + +const char CHUNK_NORMALIZE_V[] = { + 1, 1, 1, + 0, 0, 0, + -1, -1, -1, + + 1, 1, 1, + 0, 0, 0, + -1, -1, -1, + + 1, 1, 1, + 0, 0, 0, + -1, -1, -1, +}; + +const char CHUNK_NORMALIZE_W[] = { + 1, 1, 1, + 1, 1, 1, + 1, 1, 1, + + 0, 0, 0, + 0, 0, 0, + 0, 0, 0, + + -1, -1, -1, + -1, -1, -1, + -1, -1, -1, +}; \ No newline at end of file diff --git a/src/fluid/src/densitystep.c b/src/fluid/src/densitystep.c new file mode 100644 index 00000000..0cd1266e --- /dev/null +++ b/src/fluid/src/densitystep.c @@ -0,0 +1,230 @@ +#include +#include +#include + +#include "../includes/utilities.h" +#include "../includes/chunkmask.h" + + +/** + * Adds density to the density array +*/ +static inline void addDensity( + int N, + int chunk_mask, + float ** d, + float ** d0, + float dt +){ + int i; + int size=N*N*N; + float * x = GET_ARR_RAW(d,CENTER_LOC); + float * s = GET_ARR_RAW(d0,CENTER_LOC); + for(i=0; iN-1){ + for(i=i-8; i < N-1; i++){ + x[IX(i,j,k)] = (x0[IX(i,j,k)] + a*(x[IX(i-1,j,k)]+x[IX(i+1,j,k)]+x[IX(i,j-1,k)]+x[IX(i,j+1,k)]+x[IX(i,j,k-1)]+x[IX(i,j,k+1)]))/c; + } + } + } + } +} + +/** + * Advects the density based on the vectors +*/ +static inline void advectDensity(uint32_t chunk_mask, int N, float ** d, float ** d0, float ** ur, float ** vr, float ** wr, float dt){ + int i, j, k, i0, j0, k0, i1, j1, k1; + int m,n,o; + float x, y, z, s0, t0, s1, t1, u1, u0, dtx,dty,dtz; + + dtx=dty=dtz=dt*N; + + float * center_d = GET_ARR_RAW(d,CENTER_LOC); + float * center_d0 = GET_ARR_RAW(d0,CENTER_LOC); + + float * u = GET_ARR_RAW(ur,CENTER_LOC); + float * v = GET_ARR_RAW(vr,CENTER_LOC); + float * w = GET_ARR_RAW(wr,CENTER_LOC); + + for(k=1; k= N-1){ m += 1; } + if(y < 1){ n -= 1; } + if(y >= N-1){ n += 1; } + if(z < 1){ o -= 1; } + if(z >= N-1){ o += 1; } + + //If the out of bounds coordinate is in bounds for a neighbor chunk, use that chunk as source instead + // if(CK(m,n,o) != CENTER_LOC){ + // printf("Looking in border chunk\n"); + // } + // if(x > 16){ + // printf("%f %d %d %d\n",m,n,o); + // } + // if(CK(m,n,o) != CENTER_LOC && ARR_EXISTS(chunk_mask,m,n,o)){ + // // printf("Hit other chunk\n"); + // d0 = GET_ARR(env,jrd0,CK(m,n,o)); + // x = x + CHUNK_NORMALIZE_U[CK(m,n,o)] * (N-1); + // y = y + CHUNK_NORMALIZE_V[CK(m,n,o)] * (N-1); + // z = z + CHUNK_NORMALIZE_W[CK(m,n,o)] * (N-1); + // } + + if(x < 0.001f){ + //cases to consider: + //m = 0, x = -10 + //m = 2, x = 0.01 + x=0.001f; + i0=(int)0; + i1=1; + s0 = 0.999f; + s1 = 0.001f; + } else if(x >= N - 1){ + //cases to consider: + //m = 0, x = 17.01 + //m = 2, x = 20 + x = N-1; + i0=(int)N-2; + i1=N-1; + s0 = 0.001f; + s1 = 0.999f; + } else { + i0=(int)x; + i1=i0+1; + s1 = x-i0; + s0 = 1-s1; + } + + //clamp location within chunk + // if (x<0.5f) x=0.5f; + // if (x>N+0.5f) x=N+0.5f; + if (y<0.5f) y=0.5f; + if (y>N+0.5f) y=N+0.5f; + if (z<0.5f) z=0.5f; + if (z>N+0.5f) z=N+0.5f; + + //get actual indices + // i0=(int)x; + // i1=i0+1; + j0=(int)y; + j1=j0+1; + k0=(int)z; + k1=k0+1; + + //calculate percentage of each index + // s1 = x-i0; + // s0 = 1-s1; + t1 = y-j0; + t0 = 1-t1; + u1 = z-k0; + u0 = 1-u1; + + if(i0 >= N){ + i0 = N - 1; + } + // if(i0 < 0){ + // i0 = 0; + // } + if(j0 >= N){ + j0 = N - 1; + } + // if(j0 < 0){ + // j0 = 0; + // } + if(k0 >= N){ + k0 = N - 1; + } + // if(k0 < 0){ + // k0 = 0; + // } + if(i1 >= N){ + i1 = N - 1; + } + // if(i1 < 0){ + // i1 = 0; + // } + if(j1 >= N){ + j1 = N - 1; + } + // if(j1 < 0){ + // j1 = 0; + // } + if(k1 >= N){ + k1 = N - 1; + } + // if(k1 < 0){ + // k1 = 0; + // } + center_d[IX(i,j,k)] = + s0*( + t0*u0*center_d0[IX(i0,j0,k0)]+ + t1*u0*center_d0[IX(i0,j1,k0)]+ + t0*u1*center_d0[IX(i0,j0,k1)]+ + t1*u1*center_d0[IX(i0,j1,k1)] + )+ + s1*( + t0*u0*center_d0[IX(i1,j0,k0)]+ + t1*u0*center_d0[IX(i1,j1,k0)]+ + t0*u1*center_d0[IX(i1,j0,k1)]+ + t1*u1*center_d0[IX(i1,j1,k1)] + ); + } + } + } +} \ No newline at end of file diff --git a/src/fluid/src/fluidsim.c b/src/fluid/src/fluidsim.c new file mode 100644 index 00000000..06e21628 --- /dev/null +++ b/src/fluid/src/fluidsim.c @@ -0,0 +1,542 @@ +#include +#include "../includes/utilities.h" +#include "../includes/chunkmask.h" +#include "../includes/electrosphere_FluidSim.h" +#include "../includes/mainFunctions.h" +#include "../includes/chunk.h" +#include "../includes/simulation.h" + +#include "./chunkmask.c" +#include "./velocitystep.c" +#include "./densitystep.c" + +#ifndef SAVE_STEPS +#define SAVE_STEPS 0 +#endif + + +#define DIM 18 +#define LINEARSOLVERTIMES 20 +#define REALLY_SMALL_VALUE 0.00001 + +#define DIFFUSION_CONSTANT 0.00001 +#define VISCOSITY_CONSTANT 0.00001 + +char fileNameBuff[50]; + +//all chunks +Chunk ** chunks = NULL; + +//jni help: +//https://stackoverflow.com/questions/39823375/clarification-about-getfieldid + +static inline void saveStep(float * values, const char * name); + +void simulate( + int numChunks, + Chunk ** passedInChunks, + jfloat timestep +){ + chunks = passedInChunks; + + // printf("%p\n",chunks[0].d); + saveStep(chunks[0]->u[CENTER_LOC], "./chunks/beginU"); + saveStep(chunks[0]->v[CENTER_LOC], "./chunks/beginV"); + saveStep(chunks[0]->w[CENTER_LOC], "./chunks/beginW"); + saveStep(chunks[0]->u0[CENTER_LOC], "./chunks/beginU0"); + saveStep(chunks[0]->v0[CENTER_LOC], "./chunks/beginV0"); + saveStep(chunks[0]->w0[CENTER_LOC], "./chunks/beginW0"); + + //solve chunk mask + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + addSourceToVectors( + DIM, + currentChunk->chunkMask, + currentChunk->u, + currentChunk->v, + currentChunk->w, + currentChunk->u0, + currentChunk->v0, + currentChunk->w0, + DIFFUSION_CONSTANT, + VISCOSITY_CONSTANT, + timestep + ); + saveStep(currentChunk->u[CENTER_LOC], "./chunks/addSrcU"); + saveStep(currentChunk->v[CENTER_LOC], "./chunks/addSrcV"); + saveStep(currentChunk->w[CENTER_LOC], "./chunks/addSrcW"); + saveStep(currentChunk->u0[CENTER_LOC], "./chunks/addSrcU0"); + saveStep(currentChunk->v0[CENTER_LOC], "./chunks/addSrcV0"); + saveStep(currentChunk->w0[CENTER_LOC], "./chunks/addSrcW0"); + } + //swap all vector fields + { + //swap vector fields + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + + float * tmpArr; + for(int j = 0; j < 27; j++){ + tmpArr = currentChunk->u[j]; + currentChunk->u[j] = currentChunk->u0[j]; + currentChunk->u0[j] = tmpArr; + } + for(int j = 0; j < 27; j++){ + tmpArr = currentChunk->v[j]; + currentChunk->v[j] = currentChunk->v0[j]; + currentChunk->v0[j] = tmpArr; + } + for(int j = 0; j < 27; j++){ + tmpArr = currentChunk->w[j]; + currentChunk->w[j] = currentChunk->w0[j]; + currentChunk->w0[j] = tmpArr; + } + } + //copy neighbors + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w0); + } + } + saveStep(chunks[0]->u[CENTER_LOC], "./chunks/swapU"); + saveStep(chunks[0]->v[CENTER_LOC], "./chunks/swapV"); + saveStep(chunks[0]->w[CENTER_LOC], "./chunks/swapW"); + saveStep(chunks[0]->u0[CENTER_LOC], "./chunks/swapU0"); + saveStep(chunks[0]->v0[CENTER_LOC], "./chunks/swapV0"); + saveStep(chunks[0]->w0[CENTER_LOC], "./chunks/swapW0"); + // printf("after swap vecs u\n"); + // printLayer(chunks[0]->u[CENTER_LOC],targetLayer); + // printf("after swap vecs u0\n"); + // printLayer(chunks[0]->u0[CENTER_LOC],targetLayer); + //solve vector diffusion + { + for(int l = 0; l < LINEARSOLVERTIMES; l++){ + //solve vector diffusion + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + solveVectorDiffuse(DIM,currentChunk->chunkMask,currentChunk->u,currentChunk->v,currentChunk->w,currentChunk->u0,currentChunk->v0,currentChunk->w0,DIFFUSION_CONSTANT,VISCOSITY_CONSTANT,timestep); + } + if(SAVE_STEPS){ + sprintf(fileNameBuff, "./chunks/diffuseUStep%dx", l); + saveStep(chunks[0]->u[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseUStep%dx0", l); + saveStep(chunks[0]->u0[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseVStep%dx", l); + saveStep(chunks[0]->v[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseVStep%dx0", l); + saveStep(chunks[0]->v0[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseWStep%dx", l); + saveStep(chunks[0]->w[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseWStep%dx0", l); + saveStep(chunks[0]->w0[CENTER_LOC], fileNameBuff); + } + //update array for vectors + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,3,currentChunk->w); + // setBoundsToNeighborsRaw(DIM,chunkMask,1,currentChunk->u0); + // setBoundsToNeighborsRaw(DIM,chunkMask,2,currentChunk->v0); + // setBoundsToNeighborsRaw(DIM,chunkMask,3,currentChunk->w0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->u); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->v); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->w); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->w0); + } + if(SAVE_STEPS){ + sprintf(fileNameBuff, "./chunks/diffuseUStep%dxBnd", l); + saveStep(chunks[0]->u[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseUStep%dx0Bnd", l); + saveStep(chunks[0]->u0[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseVStep%dxBnd", l); + saveStep(chunks[0]->v[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseVStep%dx0Bnd", l); + saveStep(chunks[0]->v0[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseWStep%dxBnd", l); + saveStep(chunks[0]->w[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/diffuseWStep%dx0Bnd", l); + saveStep(chunks[0]->w0[CENTER_LOC], fileNameBuff); + } + } + } + saveStep(chunks[0]->u[CENTER_LOC], "./chunks/diffuseU"); + saveStep(chunks[0]->v[CENTER_LOC], "./chunks/diffuseV"); + saveStep(chunks[0]->w[CENTER_LOC], "./chunks/diffuseW"); + saveStep(chunks[0]->u0[CENTER_LOC], "./chunks/diffuseU0"); + saveStep(chunks[0]->v0[CENTER_LOC], "./chunks/diffuseV0"); + saveStep(chunks[0]->w0[CENTER_LOC], "./chunks/diffuseW0"); + //solve projection + { + //update array for vectors + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + // setBoundsToNeighborsRaw(DIM,chunkMask,1,currentChunk->u); + // setBoundsToNeighborsRaw(DIM,chunkMask,2,currentChunk->v); + // setBoundsToNeighborsRaw(DIM,chunkMask,3,currentChunk->w); + // setBoundsToNeighborsRaw(DIM,chunkMask,1,currentChunk->u0); + // setBoundsToNeighborsRaw(DIM,chunkMask,2,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v0); + } + //setup projection + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setupProjection(DIM,currentChunk->chunkMask,currentChunk->u,currentChunk->v,currentChunk->w,currentChunk->u0,currentChunk->v0,DIFFUSION_CONSTANT,VISCOSITY_CONSTANT,timestep); + } + + saveStep(chunks[0]->v0[CENTER_LOC], "./chunks/setupProj1Div"); + saveStep(chunks[0]->u0[CENTER_LOC], "./chunks/setupProj1P"); + + //update array for vectors + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,0,currentChunk->u0); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,0,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->v0); + } + + saveStep(chunks[0]->v0[CENTER_LOC], "./chunks/setupProj1DivBnd"); + saveStep(chunks[0]->u0[CENTER_LOC], "./chunks/setupProj1PBnd"); + + //samples u0, v0 + //sets u0 + //these should have just been mirrored in the above + // + //Perform main projection solver + for(int l = 0; l < LINEARSOLVERTIMES; l++){ + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + solveProjection(DIM,currentChunk->chunkMask,currentChunk->u,currentChunk->v,currentChunk->w,currentChunk->u0,currentChunk->v0,currentChunk->w0,DIFFUSION_CONSTANT,VISCOSITY_CONSTANT,timestep); + } + if(SAVE_STEPS){ + sprintf(fileNameBuff, "./chunks/proj1Step%dx", l); + saveStep(chunks[0]->u0[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/proj1Step%dx0", l); + saveStep(chunks[0]->v0[CENTER_LOC], fileNameBuff); + } + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,0,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->u0); + } + if(SAVE_STEPS){ + sprintf(fileNameBuff, "./chunks/proj1Step%dxBnd", l); + saveStep(chunks[0]->u0[CENTER_LOC], fileNameBuff); + sprintf(fileNameBuff, "./chunks/proj1Step%dx0Bnd", l); + saveStep(chunks[0]->v0[CENTER_LOC], fileNameBuff); + } + } + //samples u,v,w,u0 + //sets u,v,w + //Finalize projection + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + finalizeProjection(DIM,currentChunk->chunkMask,currentChunk->u,currentChunk->v,currentChunk->w,currentChunk->u0,currentChunk->v0,currentChunk->w0,DIFFUSION_CONSTANT,VISCOSITY_CONSTANT,timestep); + } + + saveStep(chunks[0]->u[CENTER_LOC], "./chunks/finalizeProj1U"); + saveStep(chunks[0]->v[CENTER_LOC], "./chunks/finalizeProj1V"); + saveStep(chunks[0]->w[CENTER_LOC], "./chunks/finalizeProj1W"); + // exit(0); + //set boundaries a final time for u,v,w + //... + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,3,currentChunk->w); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u0); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v0); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,3,currentChunk->w0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w0); + } + } + saveStep(chunks[0]->u[CENTER_LOC], "./chunks/projU"); + saveStep(chunks[0]->v[CENTER_LOC], "./chunks/projV"); + saveStep(chunks[0]->w[CENTER_LOC], "./chunks/projW"); + saveStep(chunks[0]->u0[CENTER_LOC], "./chunks/projU0"); + saveStep(chunks[0]->v0[CENTER_LOC], "./chunks/projV0"); + saveStep(chunks[0]->w0[CENTER_LOC], "./chunks/projW0"); + // exit(0); + //swap all vector fields + { + //swap vector fields + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + + float * tmpArr; + for(int j = 0; j < 27; j++){ + tmpArr = currentChunk->u[j]; + currentChunk->u[j] = currentChunk->u0[j]; + currentChunk->u0[j] = tmpArr; + } + for(int j = 0; j < 27; j++){ + tmpArr = currentChunk->v[j]; + currentChunk->v[j] = currentChunk->v0[j]; + currentChunk->v0[j] = tmpArr; + } + for(int j = 0; j < 27; j++){ + tmpArr = currentChunk->w[j]; + currentChunk->w[j] = currentChunk->w0[j]; + currentChunk->w0[j] = tmpArr; + } + } + //copy neighbors + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w0); + } + } + saveStep(chunks[0]->u[CENTER_LOC], "./chunks/swap2U"); + saveStep(chunks[0]->v[CENTER_LOC], "./chunks/swap2V"); + saveStep(chunks[0]->w[CENTER_LOC], "./chunks/swap2W"); + saveStep(chunks[0]->u0[CENTER_LOC], "./chunks/swap2U0"); + saveStep(chunks[0]->v0[CENTER_LOC], "./chunks/swap2V0"); + saveStep(chunks[0]->w0[CENTER_LOC], "./chunks/swap2W0"); + //advect vectors across boundaries + { + //update border arrs + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,3,currentChunk->w); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u0); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v0); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,3,currentChunk->w0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w0); + } + //advect + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + advectVectors(DIM,currentChunk->chunkMask,currentChunk->u,currentChunk->v,currentChunk->w,currentChunk->u0,currentChunk->v0,currentChunk->w0,DIFFUSION_CONSTANT,VISCOSITY_CONSTANT,timestep); + } + //update neighbor arr + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,3,currentChunk->w); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w); + } + } + saveStep(chunks[0]->u[CENTER_LOC], "./chunks/advectU"); + saveStep(chunks[0]->v[CENTER_LOC], "./chunks/advectV"); + saveStep(chunks[0]->w[CENTER_LOC], "./chunks/advectW"); + saveStep(chunks[0]->u0[CENTER_LOC], "./chunks/advectU0"); + saveStep(chunks[0]->v0[CENTER_LOC], "./chunks/advectV0"); + saveStep(chunks[0]->w0[CENTER_LOC], "./chunks/advectW0"); + //solve projection + { + //update array for vectors + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,3,currentChunk->w); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u0); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v0); + } + //setup projection + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setupProjection(DIM,currentChunk->chunkMask,currentChunk->u,currentChunk->v,currentChunk->w,currentChunk->u0,currentChunk->v0,DIFFUSION_CONSTANT,VISCOSITY_CONSTANT,timestep); + } + //update array for vectors + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u0); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v0); + } + //samples u0, v0 + //sets u0 + //these should have just been mirrored in the above + // + //Perform main projection solver + for(int l = 0; l < LINEARSOLVERTIMES; l++){ + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + solveProjection(DIM,currentChunk->chunkMask,currentChunk->u,currentChunk->v,currentChunk->w,currentChunk->u0,currentChunk->v0,currentChunk->w0,DIFFUSION_CONSTANT,VISCOSITY_CONSTANT,timestep); + } + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,0,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->u0); + } + } + //samples u,v,w,u0 + //sets u,v,w + //Finalize projection + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + finalizeProjection(DIM,currentChunk->chunkMask,currentChunk->u,currentChunk->v,currentChunk->w,currentChunk->u0,currentChunk->v0,currentChunk->w0,DIFFUSION_CONSTANT,VISCOSITY_CONSTANT,timestep); + } + //set boundaries a final time for u,v,w + //... + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,3,currentChunk->w); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,1,currentChunk->u0); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,2,currentChunk->v0); + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,3,currentChunk->w0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,1,currentChunk->u0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,2,currentChunk->v0); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,3,currentChunk->w0); + } + } + + + + + ///------------------------------------------------------------------------------------------------------------------------------------------------------------------------ + ///------------------------------------------------------------------------------------------------------------------------------------------------------------------------ + ///------------------------------------------------------------------------------------------------------------------------------------------------------------------------ + ///------------------------------------------------------------------------------------------------------------------------------------------------------------------------ + ///------------------------------------------------------------------------------------------------------------------------------------------------------------------------ + + + + + //add density + { + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + addDensity(DIM,currentChunk->chunkMask,currentChunk->d,currentChunk->d0,timestep); + } + } + //swap all density arrays + { + //swap vector fields + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + + float * tmpArr; + for(int j = 0; j < 27; j++){ + tmpArr = currentChunk->d[j]; + currentChunk->d[j] = currentChunk->d0[j]; + currentChunk->d0[j] = tmpArr; + } + } + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->d); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->d0); + } + } + //diffuse density + { + for(int l = 0; l < LINEARSOLVERTIMES; l++){ + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + solveDiffuseDensity(DIM,currentChunk->chunkMask,currentChunk->d,currentChunk->d0,currentChunk->u,currentChunk->v,currentChunk->w,DIFFUSION_CONSTANT,VISCOSITY_CONSTANT,timestep); + } + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,0,currentChunk->d); + } + } + } + //swap all density arrays + { + //swap vector fields + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + float * tmpArr; + for(int j = 0; j < 27; j++){ + tmpArr = currentChunk->d[j]; + currentChunk->d[j] = currentChunk->d0[j]; + currentChunk->d0[j] = tmpArr; + } + } + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->d); + copyNeighborsRaw(DIM,currentChunk->chunkMask,0,0,currentChunk->d0); + } + } + //advect density + { + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + advectDensity(currentChunk->chunkMask,DIM,currentChunk->d,currentChunk->d0,currentChunk->u,currentChunk->v,currentChunk->w,timestep); + } + } + //mirror densities + { + for(int i = 0; i < numChunks; i++){ + Chunk * currentChunk = chunks[i]; + setBoundsToNeighborsRaw(DIM,currentChunk->chunkMask,0,currentChunk->d); + } + } +} + + +static inline void saveStep(float * values, const char * name){ + if(SAVE_STEPS){ + FILE *fp; + int N = DIM; + + // ... fill the array somehow ... + + fp = fopen(name, "w"); + // check for error here + + for(int x = 0; x < DIM; x++){ + for(int y = 0; y < DIM; y++){ + for(int z = 0; z < DIM; z++){ + float val = values[IX(x,y,z)]; + if(val < REALLY_SMALL_VALUE && val > -REALLY_SMALL_VALUE){ + val = 0; + } + fprintf(fp, "%f\t", val); + } + fprintf(fp, "\n"); + } + fprintf(fp, "\n"); + } + + fclose(fp); + } +} \ No newline at end of file diff --git a/src/fluid/src/javainterface.c b/src/fluid/src/javainterface.c new file mode 100644 index 00000000..ef4e5581 --- /dev/null +++ b/src/fluid/src/javainterface.c @@ -0,0 +1,170 @@ +#include + +//library includes +//include stb ds +#define STB_DS_IMPLEMENTATION +#include "../lib/stb/stb_ds.h" + +//local includes +#include "../includes/chunk.h" +#include "../includes/chunkmask.h" +#include "../includes/utilities.h" +#include "../includes/simulation.h" + + +//defines + + +//function defines +#define getChunk(i) (*env)->CallObjectMethod(env,chunkList,jListGet,i) +#define getBuffArr(buffId) (*env)->GetObjectField(env,chunkJRaw,buffId) +#define setBuffArr(buffId,value) (*env)->SetObjectField(env,chunkJRaw,buffId,value) +#define GET_ARR(env,src,i) (*env)->GetDirectBufferAddress(env,(*env)->GetObjectArrayElement(env,src,i)) + + + +//declarations +void readInChunks(JNIEnv * env, jobject chunkList); +int calculateChunkMask(JNIEnv * env, jobjectArray jrx); +uint32_t matrix_transform(JNIEnv * env, jobjectArray jrx); + + + +//the list of chunks +Chunk ** javaChunkView = NULL; +//the number of chunks +int numChunks = 0; + + +JNIEXPORT void JNICALL Java_electrosphere_FluidSim_simulate( + JNIEnv * env, + jclass fluidSimClass, + jobject chunkList, + jfloat dt + ){ + readInChunks(env,chunkList); + simulate(numChunks,javaChunkView,dt); +} + + +/** + * Reads chunks into the dynamic array +*/ +void readInChunks(JNIEnv * env, jobject chunkList){ + jclass listClass = (*env)->FindClass(env,"java/util/List"); + jclass fluidSimClass = (*env)->FindClass(env,"electrosphere/FluidSim"); + //JNIEnv *env, jclass clazz, const char *name, const char *sig + jmethodID jListSize = (*env)->GetMethodID(env, listClass, "size", "()I"); + jmethodID jListGet = (*env)->GetMethodID(env, listClass, "get", "(I)Ljava/lang/Object;"); + jmethodID jListAdd = (*env)->GetMethodID(env, listClass, "add", "(Ljava/lang/Object;)Z"); + + //ByteBuffer[] + jfieldID dJId = (*env)->GetFieldID(env,fluidSimClass,"density","[Ljava/nio/ByteBuffer;"); + jfieldID d0JId = (*env)->GetFieldID(env,fluidSimClass,"densityAddition","[Ljava/nio/ByteBuffer;"); + jfieldID uJId = (*env)->GetFieldID(env,fluidSimClass,"uVector","[Ljava/nio/ByteBuffer;"); + jfieldID vJId = (*env)->GetFieldID(env,fluidSimClass,"vVector","[Ljava/nio/ByteBuffer;"); + jfieldID wJId = (*env)->GetFieldID(env,fluidSimClass,"wVector","[Ljava/nio/ByteBuffer;"); + jfieldID u0JId = (*env)->GetFieldID(env,fluidSimClass,"uAdditionVector","[Ljava/nio/ByteBuffer;"); + jfieldID v0JId = (*env)->GetFieldID(env,fluidSimClass,"vAdditionVector","[Ljava/nio/ByteBuffer;"); + jfieldID w0JId = (*env)->GetFieldID(env,fluidSimClass,"wAdditionVector","[Ljava/nio/ByteBuffer;"); + jfieldID chunkmaskJId = (*env)->GetFieldID(env,fluidSimClass,"chunkMask","I"); + + //the number of chunks + numChunks = (*env)->CallIntMethod(env,chunkList,jListSize); + + //current chunk (this) + jobject chunkJRaw; + //current chunk fields + jobjectArray jd; + jobjectArray jd0; + jobjectArray u; + jobjectArray v; + jobjectArray w; + jobjectArray u0; + jobjectArray v0; + jobjectArray w0; + int chunkMask; + + //solve chunk mask + for(int i = 0; i < numChunks; i++){ + chunkJRaw = getChunk(i); + chunkMask = calculateChunkMask(env,getBuffArr(dJId)); + (*env)->SetIntField(env,chunkJRaw,chunkmaskJId,chunkMask); + + Chunk * newChunk; + if(i >= stbds_arrlen(javaChunkView)){ + // printf("allocate chunk %d\n",i); + // fflush(stdout); + newChunk = (Chunk *)malloc(sizeof(Chunk)); + // printf("new chunk %p\n",newChunk); + // fflush(stdout); + stbds_arrput(javaChunkView,newChunk); + // printf("new chunk %p\n",chunks[i]); + // fflush(stdout); + } else { + newChunk = javaChunkView[i]; + // printf("get chunk %d: %p\n",i,newChunk); + // fflush(stdout); + } + jd = (*env)->GetObjectField(env,chunkJRaw,dJId); + jd0 = (*env)->GetObjectField(env,chunkJRaw,d0JId); + u = (*env)->GetObjectField(env,chunkJRaw,uJId); + v = (*env)->GetObjectField(env,chunkJRaw,vJId); + w = (*env)->GetObjectField(env,chunkJRaw,wJId); + u0 = (*env)->GetObjectField(env,chunkJRaw,u0JId); + v0 = (*env)->GetObjectField(env,chunkJRaw,v0JId); + w0 = (*env)->GetObjectField(env,chunkJRaw,w0JId); + newChunk->chunkMask = chunkMask; + for(int j = 0; j < 27; j++){ + if((chunkMask & CHUNK_INDEX_ARR[j]) > 0){ + newChunk->d[j] = GET_ARR(env,jd,j); + newChunk->d0[j] = GET_ARR(env,jd0,j); + newChunk->u[j] = GET_ARR(env,u,j); + newChunk->v[j] = GET_ARR(env,v,j); + newChunk->w[j] = GET_ARR(env,w,j); + newChunk->u0[j] = GET_ARR(env,u0,j); + newChunk->v0[j] = GET_ARR(env,v0,j); + newChunk->w0[j] = GET_ARR(env,w0,j); + } + } + } +} + + +/** + * Calculates the bitmask for available chunks for the provided chunk's neighbor array +*/ +int calculateChunkMask(JNIEnv * env, jobjectArray jrx){ + return matrix_transform(env,jrx); +} + +/** + * Calculates a mask that represents all nearby chunks that are actually accessible and exist +*/ +uint32_t matrix_transform(JNIEnv * env, jobjectArray jrx){ + + //The returned value, an availability mask that contains the availability of each neighbor chunk + uint32_t rVal = 0; + + //Add to maks for initial chunks + for(int i = 0; i < CENTER_LOC; i++){ + if((*env)->GetObjectArrayElement(env,jrx,i)!=NULL){ + rVal = rVal + 1; + } + rVal = rVal << 1; + } + //add 1 for center chunk because we already have that + rVal = rVal + 1; + rVal = rVal << 1; + //continue on for remaining chunks + for(int i = CENTER_LOC+1; i < 27; i++){ + if((*env)->GetObjectArrayElement(env,jrx,i)!=NULL){ + rVal = rVal + 1; + } + if(i < 26){ + rVal = rVal << 1; + } + } + + return rVal; +} \ No newline at end of file diff --git a/src/fluid/src/velocitystep.c b/src/fluid/src/velocitystep.c new file mode 100644 index 00000000..66665347 --- /dev/null +++ b/src/fluid/src/velocitystep.c @@ -0,0 +1,936 @@ +#include +#include +#include + +#include "../includes/utilities.h" +#include "../includes/chunkmask.h" + + +#define BOUND_NO_DIR 0 +#define BOUND_DIR_U 1 +#define BOUND_DIR_V 2 +#define BOUND_DIR_W 3 + +#define SET_BOUND_IGNORE 0 +#define SET_BOUND_USE_NEIGHBOR 1 + +#define LINEARSOLVERTIMES 20 + +static inline void add_source(int N, float * x, float * s, float dt); +static inline void advect(uint32_t chunk_mask, int N, int b, float ** jrd, float ** jrd0, float * u, float * v, float * w, float dt); + + +/* + * Adds force to all vectors + */ +static inline void addSourceToVectors + ( + int N, + int chunk_mask, + float ** jru, + float ** jrv, + float ** jrw, + float ** jru0, + float ** jrv0, + float ** jrw0, + float DIFFUSION_CONST, + float VISCOSITY_CONST, + float dt){ + add_source(N,GET_ARR_RAW(jru,CENTER_LOC),GET_ARR_RAW(jru0,CENTER_LOC),dt); + add_source(N,GET_ARR_RAW(jrv,CENTER_LOC),GET_ARR_RAW(jrv0,CENTER_LOC),dt); + add_source(N,GET_ARR_RAW(jrw,CENTER_LOC),GET_ARR_RAW(jrw0,CENTER_LOC),dt); +} + +/** + * Adds from a source array to a destination array +*/ +static inline void add_source(int N, float * x, float * s, float dt){ + int i; + int size=N*N*N; + for(i=0; iN-1){ + for(i=i-8; i < N-1; i++){ + u[IX(i,j,k)] = (u0[IX(i,j,k)] + a*(u[IX(i-1,j,k)]+u[IX(i+1,j,k)]+u[IX(i,j-1,k)]+u[IX(i,j+1,k)]+u[IX(i,j,k-1)]+u[IX(i,j,k+1)]))/c; + } + } + } + } + + //transform v direction + for(k=1; kN-1){ + for(i=i-8; i < N-1; i++){ + v[IX(i,j,k)] = (v0[IX(i,j,k)] + a*(v[IX(i-1,j,k)]+v[IX(i+1,j,k)]+v[IX(i,j-1,k)]+v[IX(i,j+1,k)]+v[IX(i,j,k-1)]+v[IX(i,j,k+1)]))/c; + } + } + } + } + + //transform w direction + for(k=1; kN-1){ + for(i=i-8; i < N-1; i++){ + w[IX(i,j,k)] = (w0[IX(i,j,k)] + a*(w[IX(i-1,j,k)]+w[IX(i+1,j,k)]+w[IX(i,j-1,k)]+w[IX(i,j+1,k)]+w[IX(i,j,k-1)]+w[IX(i,j,k+1)]))/c; + } + } + } + } +} + +/* + * Sets up a projection system of equations + */ +static inline void setupProjection( + int N, + int chunk_mask, + float ** ur, + float ** vr, + float ** wr, + float ** pr, + float ** divr, + float DIFFUSION_CONST, + float VISCOSITY_CONST, + float dt +){ + int i, j, k; + + __m256 nVector = _mm256_set1_ps(N); + __m256 constScalar = _mm256_set1_ps(-1.0/3.0); + __m256 zeroVec = _mm256_set1_ps(0); + __m256 vector, vector2, vector3; + + float * u = GET_ARR_RAW(ur,CENTER_LOC); + float * v = GET_ARR_RAW(vr,CENTER_LOC); + float * w = GET_ARR_RAW(wr,CENTER_LOC); + + float * p = GET_ARR_RAW(pr,CENTER_LOC); + float * div = GET_ARR_RAW(divr,CENTER_LOC); + + float scalar = 1.0/3.0; + float h = 1.0/N; + + for(k=1; kN-1){ + for(i=i-8; i < N-1; i++){ + p[IX(i,j,k)] = (div[IX(i,j,k)] + a*(p[IX(i-1,j,k)]+p[IX(i+1,j,k)]+p[IX(i,j-1,k)]+p[IX(i,j+1,k)]+p[IX(i,j,k-1)]+p[IX(i,j,k+1)]))/c; + } + } + } + } +} + +/* + * Finalizes a projection (subtract curl, set bounds, etc) + */ +static inline void finalizeProjection( + int N, + int chunk_mask, + float ** jru, + float ** jrv, + float ** jrw, + float ** jru0, + float ** jrv0, + float ** jrw0, + float DIFFUSION_CONST, + float VISCOSITY_CONST, + float dt +){ + int i, j, k; + __m256 constScalar = _mm256_set1_ps(0.5f*N); + __m256 vector, vector2, vector3; + + float * u = GET_ARR_RAW(jru,CENTER_LOC); + float * v = GET_ARR_RAW(jrv,CENTER_LOC); + float * w = GET_ARR_RAW(jrw,CENTER_LOC); + + float * p = GET_ARR_RAW(jru0,CENTER_LOC); + float * div = GET_ARR_RAW(jrv0,CENTER_LOC); + + for ( k=1 ; k= N){ m -= 1; } + if(y < 0){ n += 1; } + else if(y >= N){ n -= 1; } + if(z < 0){ o += 1; } + else if(z >= N){ o -= 1; } + + //If the out of bounds coordinate is in bounds for a neighbor chunk, use that chunk as source instead + if(CK(m,n,o) != CENTER_LOC && ARR_EXISTS(chunk_mask,m,n,o)){ + + // if(i == 1 && j == 1 && k == 1){ + // printf("\narr indices: %d %d %d\n\n",m,n,o); + // } + + //cases: + //if x = 17.01, m = 2 + // 17 in current array is 1 in neighbor + // 18 in current array is 2 in neighbor + // 19 in current array is 3 in neighbor + //want to sample neighbor array at 1 & 2 + //x becomes 1.01, sampling new array (keep in mind that 0 in the new array should contain the current array values) + //modification: subtract 16 + + //cases: + //if x = 16.99, m = 2 + // 16 in current array is 0 in neighbor + // 17 in current array is 1 in neighbor + // 18 in current array is 2 in neighbor + // 19 in current array is 3 in neighbor + //want to sample current array still + //x becomes 1.01, sampling new array (keep in mind that 0 in the new array should contain the current array values) + //modification: no modification + + //if x = 0.01, m = 0 + // 0 in current array is 16 in neighbor + //-1 in current array is 15 in neighbor + //-2 in current array is 14 in neighbor + //want to sample current array still + //x becomes 15.01, sampling new array (keep in mind that 17 in the new array should contain the current array) + //modification: no modification + + //if x = -0.01, m = 0 + // 0 in current array is 16 in neighbor + //-1 in current array is 15 in neighbor + //-2 in current array is 14 in neighbor + //want to sample -1 & 0, so i0 becomes 15 + //x becomes 15.99, sampling new array (keep in mind that 17 in the new array should contain the current array) + //modification: add 16 + + //if x = -2, m = 0 + // 0 in current array is 16 in neighbor + //-1 in current array is 15 in neighbor + //-2 in current array is 14 in neighbor + //x becomes 14, sampling new array (keep in mind that 17 in the new array should contain the current array) + //modification: add 16 + + + // printf("Hit other chunk\n"); + d0 = GET_ARR_RAW(jrd0,CK(m,n,o)); + x = x + CHUNK_NORMALIZE_U[CK(m,n,o)] * (N-2); + // printf("%d => %f\n",m,x); + y = y + CHUNK_NORMALIZE_V[CK(m,n,o)] * (N-2); + z = z + CHUNK_NORMALIZE_W[CK(m,n,o)] * (N-2); + } + + //clamp location within chunk + //get indices, and calculate percentage to pull from each index + if(x < 0.001f){ + //cases to consider: + //m = 0, x = -10 + //m = 2, x = 0.01 + x=0.001f; + i0=(int)0; + i1=1; + s0 = 0.999f; + s1 = 0.001f; + } else if(x > N - 1){ + //cases to consider: + //m = 0, x = 17.01 + //m = 2, x = 20 + x = N-1; + i0=(int)N-2; + i1=N-1; + s0 = 0.001f; + s1 = 0.999f; + } else { + i0=(int)x; + i1=i0+1; + s1 = x-i0; + s0 = 1-s1; + } + + if(y < 0.001f){ + //cases to consider: + //m = 0, x = -10 + //m = 2, x = 0.01 + y=0.001f; + j0=(int)0; + j1=1; + t0 = 0.999f; + t1 = 0.001f; + } else if(y > N - 1){ + //cases to consider: + //m = 0, x = 17.01 + //m = 2, x = 20 + y = N-1; + j0=(int)N-2; + j1=N-1; + t0 = 0.001f; + t1 = 0.999f; + } else { + j0=(int)y; + j1=j0+1; + t1 = y-j0; + t0 = 1-t1; + } + + + if(z < 0.001f){ + //cases to consider: + //m = 0, x = -10 + //m = 2, x = 0.01 + z=0.001f; + k0=(int)0; + k1=1; + u0 = 0.999f; + u1 = 0.001f; + } else if(z > N - 1){ + //cases to consider: + //m = 0, x = 17.01 + //m = 2, x = 20 + z = N-1; + k0=(int)N-2; + k1=N-1; + u0 = 0.001f; + u1 = 0.999f; + } else { + k0=(int)z; + k1=k0+1; + u1 = z-k0; + u0 = 1-u1; + } + + // if (x<0.001f) x=0.001f; + // if (x>N+0.5f) x=N+0.5f; + // if (y<0.001f) y=0.001f; + // if (y>N+0.5f) y=N+0.5f; + // if (z<0.001f) z=0.001f; + // if (z>N+0.5f) z=N+0.5f; + + //get actual indices + // i0=(int)x; + // i1=i0+1; + // j0=(int)y; + // j1=j0+1; + // k0=(int)z; + // k1=k0+1; + + //calculate percentage of each index + // s1 = x-i0; + // s0 = 1-s1; + // t1 = y-j0; + // t0 = 1-t1; + // u1 = z-k0; + // u0 = 1-u1; + + // if(i0 >= N){ + // i0 = N - 1; + // } + // if(i0 < 0){ + // i0 = 0; + // } + if(j0 >= N){ + j0 = N - 1; + } + // if(j0 < 0){ + // j0 = 0; + // } + if(k0 >= N){ + k0 = N - 1; + } + // if(k0 < 0){ + // k0 = 0; + // } + // if(i1 >= N){ + // i1 = N - 1; + // } + // if(i1 < 0){ + // i1 = 0; + // } + if(j1 >= N){ + j1 = N - 1; + } + // if(j1 < 0){ + // j1 = 0; + // } + if(k1 >= N){ + k1 = N - 1; + } + // if(k1 < 0){ + // k1 = 0; + // } + d[IX(i,j,k)] = + s0*( + t0*u0*d0[IX(i0,j0,k0)]+ + t1*u0*d0[IX(i0,j1,k0)]+ + t0*u1*d0[IX(i0,j0,k1)]+ + t1*u1*d0[IX(i0,j1,k1)] + )+ + s1*( + t0*u0*d0[IX(i1,j0,k0)]+ + t1*u0*d0[IX(i1,j1,k0)]+ + t0*u1*d0[IX(i1,j0,k1)]+ + t1*u1*d0[IX(i1,j1,k1)] + ); + } + } + } +} + +/** + * Sets the bounds of this cube to those of its neighbor +*/ +static inline void setBoundsToNeighborsRaw( + int N, + int chunk_mask, + int vector_dir, + float ** neighborArray +){ + int DIM = N; + float * target = GET_ARR_RAW(neighborArray,CENTER_LOC); + float * source; + //set the faces bounds + for(int x=1; x < DIM-1; x++){ + for(int y = 1; y < DIM-1; y++){ + //((x)+(DIM)*(y) + (DIM)*(DIM)*(z)) + target[IX(0,x,y)] = vector_dir==BOUND_DIR_U ? -target[IX(1,x,y)] : target[IX(1,x,y)]; + target[IX(DIM-1,x,y)] = vector_dir==BOUND_DIR_U ? -target[IX(DIM-2,x,y)] : target[IX(DIM-2,x,y)]; + target[IX(x,0,y)] = vector_dir==BOUND_DIR_V ? -target[IX(x,1,y)] : target[IX(x,1,y)]; + target[IX(x,DIM-1,y)] = vector_dir==BOUND_DIR_V ? -target[IX(x,DIM-2,y)] : target[IX(x,DIM-2,y)]; + target[IX(x,y,0)] = vector_dir==BOUND_DIR_W ? -target[IX(x,y,1)] : target[IX(x,y,1)]; + target[IX(x,y,DIM-1)] = vector_dir==BOUND_DIR_W ? -target[IX(x,y,DIM-2)] : target[IX(x,y,DIM-2)]; + } + } + //sets the edges of the chunk + for(int x = 1; x < DIM-1; x++){ + target[IX(x,0,0)] = (float)(0.5f * (target[IX(x,1,0)] + target[IX(x,0,1)])); + target[IX(x,DIM-1,0)] = (float)(0.5f * (target[IX(x,DIM-2,0)] + target[IX(x,DIM-1,1)])); + target[IX(x,0,DIM-1)] = (float)(0.5f * (target[IX(x,1,DIM-1)] + target[IX(x,0,DIM-2)])); + target[IX(x,DIM-1,DIM-1)] = (float)(0.5f * (target[IX(x,DIM-2,DIM-1)] + target[IX(x,DIM-1,DIM-2)])); + + target[IX(0,x,0)] = (float)(0.5f * (target[IX(1,x,0)] + target[IX(0,x,1)])); + target[IX(DIM-1,x,0)] = (float)(0.5f * (target[IX(DIM-2,x,0)] + target[IX(DIM-1,x,1)])); + target[IX(0,x,DIM-1)] = (float)(0.5f * (target[IX(1,x,DIM-1)] + target[IX(0,x,DIM-2)])); + target[IX(DIM-1,x,DIM-1)] = (float)(0.5f * (target[IX(DIM-2,x,DIM-1)] + target[IX(DIM-1,x,DIM-2)])); + + + target[IX(0,0,x)] = (float)(0.5f * (target[IX(1,0,x)] + target[IX(0,1,x)])); + target[IX(DIM-1,0,x)] = (float)(0.5f * (target[IX(DIM-2,0,x)] + target[IX(DIM-1,1,x)])); + target[IX(0,DIM-1,x)] = (float)(0.5f * (target[IX(1,DIM-1,x)] + target[IX(0,DIM-2,x)])); + target[IX(DIM-1,DIM-1,x)] = (float)(0.5f * (target[IX(DIM-2,DIM-1,x)] + target[IX(DIM-1,DIM-2,x)])); + + } + //sets the corners of the chunk + target[IX(0,0,0)] = (float)((target[IX(1,0,0)]+target[IX(0,1,0)]+target[IX(0,0,1)])/3.0); + target[IX(DIM-1,0,0)] = (float)((target[IX(DIM-2,0,0)]+target[IX(DIM-1,1,0)]+target[IX(DIM-1,0,1)])/3.0); + target[IX(0,DIM-1,0)] = (float)((target[IX(1,DIM-1,0)]+target[IX(0,DIM-2,0)]+target[IX(0,DIM-1,1)])/3.0); + target[IX(0,0,DIM-1)] = (float)((target[IX(0,0,DIM-2)]+target[IX(1,0,DIM-1)]+target[IX(0,1,DIM-1)])/3.0); + target[IX(DIM-1,DIM-1,0)] = (float)((target[IX(DIM-2,DIM-1,0)]+target[IX(DIM-1,DIM-2,0)]+target[IX(DIM-1,DIM-1,1)])/3.0); + target[IX(0,DIM-1,DIM-1)] = (float)((target[IX(1,DIM-1,DIM-1)]+target[IX(0,DIM-2,DIM-1)]+target[IX(0,DIM-1,DIM-2)])/3.0); + target[IX(DIM-1,0,DIM-1)] = (float)((target[IX(DIM-1,0,DIM-2)]+target[IX(DIM-2,0,DIM-1)]+target[IX(DIM-1,1,DIM-1)])/3.0); + target[IX(DIM-1,DIM-1,DIM-1)] = (float)((target[IX(DIM-1,DIM-1,DIM-2)]+target[IX(DIM-1,DIM-2,DIM-1)]+target[IX(DIM-1,DIM-1,DIM-2)])/3.0); +} + + +/** + * This exclusively copies neighbors to make sure zeroing out stuff doesn't break sim +*/ +static inline void copyNeighborsRaw( + int N, + int chunk_mask, + int cx, + int vector_dir, + float ** neighborArray +){ + int DIM = N; + float * target = GET_ARR_RAW(neighborArray,CENTER_LOC); + float * source; + + + // + // + // PLANES + // + // + // __m512 transferVector;// = _mm512_set1_ps(0.5*N); + + //__m256 vector = _mm256_loadu_ps(&p[IX(i-1,j,k)]); + //vector = _mm256_add_ps(vector,_mm256_loadu_ps(&p[IX(i+1,j,k)])); + //vector = _mm256_add_ps(vector,_mm256_loadu_ps(&p[IX(i,j-1,k)])); + //_mm256_storeu_ps(&p[IX(i,j,k)],vector); + //__m256 + //_mm256_loadu_ps + //_mm256_storeu_ps + if(ARR_EXISTS(chunk_mask,0,1,1)){ + source = GET_ARR_RAW(neighborArray,CK(0,1,1)); + for(int x=1; x < DIM-1; x++){ + // transferVector = _mm512_loadu_ps(&source[IX(DIM-2,x,1)]); + // _mm512_storeu_ps(&target[IX(0,x,1)],_mm512_loadu_ps(&source[IX(DIM-2,x,1)])); + for(int y = 1; y < DIM-1; y++){ + target[IX(0,x,y)] = source[IX(DIM-2,x,y)]; + } + } + } + + if(ARR_EXISTS(chunk_mask,2,1,1)){ + source = GET_ARR_RAW(neighborArray,CK(2,1,1)); + for(int x=1; x < DIM-1; x++){ + // _mm512_storeu_ps(&target[IX(DIM-1,x,1)],_mm512_loadu_ps(&source[IX(1,x,1)])); + for(int y = 1; y < DIM-1; y++){ + target[IX(DIM-1,x,y)] = source[IX(1,x,y)]; + } + } + } + + if(ARR_EXISTS(chunk_mask,1,0,1)){ + source = GET_ARR_RAW(neighborArray,CK(1,0,1)); + for(int x=1; x < DIM-1; x++){ + for(int y = 1; y < DIM-1; y++){ + target[IX(x,0,y)] = source[IX(x,DIM-2,y)]; + } + } + } + + if(ARR_EXISTS(chunk_mask,1,2,1)){ + source = GET_ARR_RAW(neighborArray,CK(1,2,1)); + for(int x=1; x < DIM-1; x++){ + for(int y = 1; y < DIM-1; y++){ + target[IX(x,DIM-1,y)] = source[IX(x,1,y)]; + } + } + } + + if(ARR_EXISTS(chunk_mask,1,1,0)){ + source = GET_ARR_RAW(neighborArray,CK(1,1,0)); + for(int x=1; x < DIM-1; x++){ + for(int y = 1; y < DIM-1; y++){ + target[IX(x,y,0)] = source[IX(x,y,DIM-2)]; + } + } + } + + if(ARR_EXISTS(chunk_mask,1,1,2)){ + source = GET_ARR_RAW(neighborArray,CK(1,1,2)); + for(int x=1; x < DIM-1; x++){ + for(int y = 1; y < DIM-1; y++){ + target[IX(x,y,DIM-1)] = source[IX(x,y,1)]; + } + } + } + + + // + // + // EDGES + // + // + if(ARR_EXISTS(chunk_mask,0,0,1)){ + source = GET_ARR_RAW(neighborArray,CK(0,0,1)); + for(int x=1; x < DIM-1; x++){ + target[IX(0,0,x)] = source[IX(DIM-2,DIM-2,x)]; + } + } + + if(ARR_EXISTS(chunk_mask,2,0,1)){ + source = GET_ARR_RAW(neighborArray,CK(2,0,1)); + for(int x=1; x < DIM-1; x++){ + target[IX(DIM-1,0,x)] = source[IX(1,DIM-2,x)]; + } + } + + if(ARR_EXISTS(chunk_mask,0,2,1)){ + source = GET_ARR_RAW(neighborArray,CK(0,2,1)); + for(int x=1; x < DIM-1; x++){ + target[IX(0,DIM-1,x)] = source[IX(DIM-2,1,x)]; + } + } + + if(ARR_EXISTS(chunk_mask,2,2,1)){ + source = GET_ARR_RAW(neighborArray,CK(2,2,1)); + for(int x=1; x < DIM-1; x++){ + target[IX(DIM-1,DIM-1,x)] = source[IX(1,1,x)]; + } + } + + // + // + + if(ARR_EXISTS(chunk_mask,0,1,0)){ + source = GET_ARR_RAW(neighborArray,CK(0,1,0)); + for(int x=1; x < DIM-1; x++){ + target[IX(0,x,0)] = source[IX(DIM-2,x,DIM-2)]; + } + } + + if(ARR_EXISTS(chunk_mask,2,1,0)){ + source = GET_ARR_RAW(neighborArray,CK(2,1,0)); + for(int x=1; x < DIM-1; x++){ + target[IX(DIM-1,x,0)] = source[IX(1,x,DIM-2)]; + } + } + + if(ARR_EXISTS(chunk_mask,0,1,2)){ + source = GET_ARR_RAW(neighborArray,CK(0,1,2)); + for(int x=1; x < DIM-1; x++){ + target[IX(0,x,DIM-1)] = source[IX(DIM-2,x,1)]; + } + } + + if(ARR_EXISTS(chunk_mask,2,1,2)){ + source = GET_ARR_RAW(neighborArray,CK(2,1,2)); + for(int x=1; x < DIM-1; x++){ + target[IX(DIM-1,x,DIM-1)] = source[IX(1,x,1)]; + } + } + + // + // + + if(ARR_EXISTS(chunk_mask,1,0,0)){ + source = GET_ARR_RAW(neighborArray,CK(1,0,0)); + for(int x=1; x < DIM-1; x++){ + target[IX(x,0,0)] = source[IX(x,DIM-2,DIM-2)]; + } + } + + if(ARR_EXISTS(chunk_mask,1,2,0)){ + source = GET_ARR_RAW(neighborArray,CK(1,2,0)); + for(int x=1; x < DIM-1; x++){ + target[IX(x,DIM-1,0)] = source[IX(x,1,DIM-2)]; + } + } + + if(ARR_EXISTS(chunk_mask,1,0,2)){ + source = GET_ARR_RAW(neighborArray,CK(1,0,2)); + for(int x=1; x < DIM-1; x++){ + target[IX(x,0,DIM-1)] = source[IX(x,DIM-2,1)]; + } + } + + if(ARR_EXISTS(chunk_mask,1,2,2)){ + source = GET_ARR_RAW(neighborArray,CK(1,2,2)); + for(int x=1; x < DIM-1; x++){ + target[IX(x,DIM-1,DIM-1)] = source[IX(x,1,1)]; + } + } + + + // + // + // CORNERS + // + // + + if(ARR_EXISTS(chunk_mask,0,0,0)){ + source = GET_ARR_RAW(neighborArray,CK(0,0,0)); + target[IX(0,0,0)] = source[IX(DIM-2,DIM-2,DIM-2)]; + } + + if(ARR_EXISTS(chunk_mask,2,0,0)){ + source = GET_ARR_RAW(neighborArray,CK(2,0,0)); + target[IX(DIM-1,0,0)] = source[IX(1,DIM-2,DIM-2)]; + } + + if(ARR_EXISTS(chunk_mask,0,2,0)){ + source = GET_ARR_RAW(neighborArray,CK(0,2,0)); + target[IX(0,DIM-1,0)] = source[IX(DIM-2,1,DIM-2)]; + } + + if(ARR_EXISTS(chunk_mask,2,2,0)){ + source = GET_ARR_RAW(neighborArray,CK(2,2,0)); + target[IX(DIM-1,DIM-1,0)] = source[IX(1,1,DIM-2)]; + } + + // + // + + if(ARR_EXISTS(chunk_mask,0,0,2)){ + source = GET_ARR_RAW(neighborArray,CK(0,0,2)); + target[IX(0,0,DIM-1)] = source[IX(DIM-2,DIM-2,1)]; + } + + if(ARR_EXISTS(chunk_mask,2,0,2)){ + source = GET_ARR_RAW(neighborArray,CK(2,0,2)); + target[IX(DIM-1,0,DIM-1)] = source[IX(1,DIM-2,1)]; + } + + if(ARR_EXISTS(chunk_mask,0,2,2)){ + source = GET_ARR_RAW(neighborArray,CK(0,2,2)); + target[IX(0,DIM-1,DIM-1)] = source[IX(DIM-2,1,1)]; + } + + if(ARR_EXISTS(chunk_mask,2,2,2)){ + source = GET_ARR_RAW(neighborArray,CK(2,2,2)); + target[IX(DIM-1,DIM-1,DIM-1)] = source[IX(1,1,1)]; + } + + + +} diff --git a/src/main/java/electrosphere/server/fluid/manager/ServerFluidChunk.java b/src/main/java/electrosphere/server/fluid/manager/ServerFluidChunk.java index 847ffa2b..c7f9bc45 100644 --- a/src/main/java/electrosphere/server/fluid/manager/ServerFluidChunk.java +++ b/src/main/java/electrosphere/server/fluid/manager/ServerFluidChunk.java @@ -1,6 +1,8 @@ package electrosphere.server.fluid.manager; +import java.nio.ByteBuffer; +import java.nio.ByteOrder; import java.nio.FloatBuffer; import org.joml.Vector3f; @@ -16,30 +18,79 @@ import electrosphere.server.terrain.manager.ServerTerrainChunk; public class ServerFluidChunk { - int worldX, worldY, worldZ; + /** + * Number of adjacent arrays + */ + static final int ARRAY_CT = 27; + + /** + * Index of the center buffer + */ + static final int CENTER_BUFF = 13; + + + /** + * The world x coordinate of this chunk + */ + int worldX; + + /** + * The world y coordinate of this chunk + */ + int worldY; + + /** + * The world z coordinate of this chunk + */ + int worldZ; + + + /** + * The float view of the center weight buffer + */ FloatBuffer weights; + + /** + * The float view of the center velocity x buffer + */ FloatBuffer velocityX; + + /** + * The float view of the center velocity y buffer + */ FloatBuffer velocityY; + + /** + * The float view of the center velocity z buffer + */ FloatBuffer velocityZ; - public ServerFluidChunk( - int worldX, - int worldY, - int worldZ, - FloatBuffer weights, - FloatBuffer velocityX, - FloatBuffer velocityY, - FloatBuffer velocityZ - ) { - this.worldX = worldX; - this.worldY = worldY; - this.worldZ = worldZ; - this.weights = weights; - this.velocityX = velocityX; - this.velocityY = velocityY; - this.velocityZ = velocityZ; - } + /** + * The array of all adjacent weight buffers for the fluid sim + */ + ByteBuffer[] bWeights = new ByteBuffer[ARRAY_CT]; + /** + * The array of all adjacent velocity x buffers for the fluid sim + */ + ByteBuffer[] bVelocityX = new ByteBuffer[ARRAY_CT]; + + /** + * The array of all adjacent velocity y buffers for the fluid sim + */ + ByteBuffer[] bVelocityY = new ByteBuffer[ARRAY_CT]; + + /** + * The array of all adjacent velocity z buffers for the fluid sim + */ + ByteBuffer[] bVelocityZ = new ByteBuffer[ARRAY_CT]; + + /** + * Constructor + * @param worldX The world x coordinate + * @param worldY The world y coordinate + * @param worldZ The world z coordinate + */ public ServerFluidChunk( int worldX, int worldY, @@ -48,20 +99,46 @@ public class ServerFluidChunk { this.worldX = worldX; this.worldY = worldY; this.worldZ = worldZ; - this.weights = BufferUtils.createFloatBuffer(ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION); - this.velocityX = BufferUtils.createFloatBuffer(ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION); - this.velocityY = BufferUtils.createFloatBuffer(ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION); - this.velocityZ = BufferUtils.createFloatBuffer(ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION); + + //allocate + this.bWeights[CENTER_BUFF] = BufferUtils.createByteBuffer(ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * 4); + this.bVelocityX[CENTER_BUFF] = BufferUtils.createByteBuffer(ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * 4); + this.bVelocityY[CENTER_BUFF] = BufferUtils.createByteBuffer(ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * 4); + this.bVelocityZ[CENTER_BUFF] = BufferUtils.createByteBuffer(ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * ServerTerrainChunk.CHUNK_DIMENSION * 4); + + //order + this.bWeights[CENTER_BUFF].order(ByteOrder.LITTLE_ENDIAN); + this.bVelocityX[CENTER_BUFF].order(ByteOrder.LITTLE_ENDIAN); + this.bVelocityY[CENTER_BUFF].order(ByteOrder.LITTLE_ENDIAN); + this.bVelocityZ[CENTER_BUFF].order(ByteOrder.LITTLE_ENDIAN); + + //get float view + this.weights = this.bWeights[CENTER_BUFF].asFloatBuffer(); + this.velocityX = this.bVelocityX[CENTER_BUFF].asFloatBuffer(); + this.velocityY = this.bVelocityY[CENTER_BUFF].asFloatBuffer(); + this.velocityZ = this.bVelocityZ[CENTER_BUFF].asFloatBuffer(); } + /** + * Gets the world x coordinate + * @return The world x coordinate + */ public int getWorldX() { return worldX; } + /** + * Gets the world y coordinate + * @return The world y coordinate + */ public int getWorldY() { return worldY; } + /** + * Gets the world z coordinate + * @return The world z coordinate + */ public int getWorldZ() { return worldZ; } @@ -74,6 +151,10 @@ public class ServerFluidChunk { return new Vector3i(worldX,worldY,worldZ); } + /** + * Gets the weights buffer + * @return The weight buffer + */ public FloatBuffer getWeights() { return weights; } @@ -109,63 +190,129 @@ public class ServerFluidChunk { weights.put(this.IX(x,y,z),weight); } - //get velocity x + /** + * Gets the velocity x buffer + * @return The velocity x buffer + */ public FloatBuffer getVelocityX() { return velocityX; } + /** + * Gets the x velocity at the point + * @param x The x coordinate + * @param y The y coordinate + * @param z The z coordinate + * @return The x velocity at the point + */ public float getVelocityX(int x, int y, int z){ return velocityX.get(this.IX(x, y, z)); } - //set velocity x + /** + * Sets the velocity x buffer + * @param velocityX The velocity x buffer + */ public void setVelocityX(FloatBuffer velocityX) { this.velocityX = velocityX; } + /** + * Sets the x velocity at the point + * @param x The x coordinate + * @param y The y coordinate + * @param z The z coordinate + * @param velocity The x velocity + */ public void setVelocityX(int x, int y, int z, float velocity){ this.velocityX.put(this.IX(x,y,z),velocity); } - //get velocity y + /** + * Gets the velocity y buffer + * @return The velocity y buffer + */ public FloatBuffer getVelocityY() { return velocityY; } + /** + * Gets the y velocity at the point + * @param x The x coordinate + * @param y The y coordinate + * @param z The z coordinate + * @return The y velocity at the point + */ public float getVelocityY(int x, int y, int z){ return velocityY.get(this.IX(x, y, z)); } - //set velocity y + /** + * Sets the velocity y buffer + * @param velocityY The velocity y buffer + */ public void setVelocityY(FloatBuffer velocityY) { this.velocityY = velocityY; } + /** + * Sets the y velocity at the point + * @param x The x coordinate + * @param y The y coordinate + * @param z The z coordinate + * @param velocity The y velocity + */ public void setVelocityY(int x, int y, int z, float velocity){ this.velocityY.put(this.IX(x,y,z),velocity); } - //get velocity z + /** + * Gets the velocity z buffer + * @return The velocity z buffer + */ public FloatBuffer getVelocityZ() { return velocityZ; } + /** + * Gets the z velocity at the point + * @param x The x coordinate + * @param y The y coordinate + * @param z The z coordinate + * @return The z velocity at the point + */ public float getVelocityZ(int x, int y, int z){ return velocityZ.get(this.IX(x, y, z)); } - //set velocity z + /** + * Sets the velocity z buffer + * @param velocityZ The velocity z buffer + */ public void setVelocityZ(FloatBuffer velocityZ) { this.velocityZ = velocityZ; } + /** + * Sets the z velocity at the point + * @param x The x coordinate + * @param y The y coordinate + * @param z The z coordinate + * @param velocity The z velocity + */ public void setVelocityZ(int x, int y, int z, float velocity){ this.velocityZ.put(this.IX(x,y,z),velocity); } - //get a velocity at a given x, y and z as a Vector3f + /** + * Gets the velocity at a given point as a vector3f + * @param x The x coordinate + * @param y The y coordinate + * @param z The z coordinate + * @return The velocity + */ public Vector3f getVelocity(int x, int y, int z){ int index = this.IX(x,y,z); return new Vector3f( @@ -176,6 +323,15 @@ public class ServerFluidChunk { } //set a velocity at a given x, y, and z given three ints + /** + * Sets the full velocity at a given point + * @param x The x coordinate + * @param y The y coordinate + * @param z The z coordinate + * @param velX The x component of the velocity + * @param velY The y component of the velocity + * @param velZ The z component of the velocity + */ public void setVelocity(int x, int y, int z, float velX, float velY, float velZ){ int index = this.IX(x,y,z); velocityX.put(index, velX); diff --git a/src/main/java/electrosphere/server/fluid/simulator/FluidAcceleratedSimulator.java b/src/main/java/electrosphere/server/fluid/simulator/FluidAcceleratedSimulator.java new file mode 100644 index 00000000..fdb9bf33 --- /dev/null +++ b/src/main/java/electrosphere/server/fluid/simulator/FluidAcceleratedSimulator.java @@ -0,0 +1,48 @@ +package electrosphere.server.fluid.simulator; + +import java.io.File; +import java.util.List; + +import electrosphere.server.fluid.manager.ServerFluidChunk; +import electrosphere.server.terrain.manager.ServerTerrainChunk; + +/** + * A c-accelerated fluid simulator + */ +public class FluidAcceleratedSimulator implements ServerFluidSimulator { + + /** + * Load fluid sim library + */ + static { + String osName = System.getProperty("os.name").toLowerCase(); + if(osName.contains("win")){ + System.load(new File("./shared-folder/libfluidsim.dll").toPath().toAbsolutePath().toString()); + } else { + System.load(new File("./shared-folder/libfluidsim.so").toPath().toAbsolutePath().toString()); + } + } + + @Override + public boolean simulate(ServerFluidChunk fluidChunk, ServerTerrainChunk terrainChunk, int worldX, int worldY, int worldZ){ + // TODO Auto-generated method stub + throw new UnsupportedOperationException("Unimplemented method 'simulate'"); + } + + /** + * Main simulation function + * @param chunks The list of chunks to simulate with + * @param timestep The timestep to simulate + */ + private static void simulateWrapper(List chunks, float timestep){ + FluidAcceleratedSimulator.simulate(chunks,timestep); + } + + /** + * Main native simulation function + * @param chunks The list of chunks to simulate with + * @param timestep The timestep to simulate + */ + private static native void simulate(List chunks, float timestep); + +}