Extract public version of SweepLB

1 files changed, 89 insertions, 0 deletions

diff --git a/src/LLBM/collide.h b/src/LLBM/collide.h
new file mode 100644
index 0000000..c08aade
--- /dev/null
+++ b/src/LLBM/collide.h
@@ -0,0 +1,89 @@
+#pragma once
+
+#include "concepts.h"
+
+struct BgkCollideO {
+
+template <concepts::Arithmetic V, concepts::Arithmetic T>
+static void apply(V f_curr[19], V f_next[19], T tau) {
+	V rho;
+	V u[3];
+	V m0 = f_curr[10] + f_curr[13] + f_curr[17];
+	V m1 = f_curr[14] + f_curr[5] + f_curr[6];
+	V m2 = f_curr[1] + f_curr[2] + f_curr[4];
+	V m3 = m0 + m1 + m2 + f_curr[0] + f_curr[11] + f_curr[12] + f_curr[15] + f_curr[16] + f_curr[18] + f_curr[3] + f_curr[7] + f_curr[8] + f_curr[9];
+	V m4 = -f_curr[11] + f_curr[7];
+	V m5 = -f_curr[15] + f_curr[3];
+	V m6 = T{1}/m3;
+	V m7 = f_curr[0] - f_curr[18];
+	rho = m3;
+	u[0] = m6*(m0 + m4 + m5 - f_curr[1] - f_curr[5] - f_curr[8]);
+	u[1] = m6*(m1 + m4 + m7 - f_curr[12] - f_curr[13] - f_curr[4]);
+	u[2] = m6*(m2 + m5 + m7 - f_curr[14] - f_curr[16] - f_curr[17]);
+	V x0 = T{1}/tau;
+	V x1 = T{0.0138888888888889}*x0;
+	V x2 = u[1] + u[2];
+	V x3 = T{9.00000000000000}*(x2*x2);
+	V x4 = T{6.00000000000000}*u[2];
+	V x5 = u[1]*u[1];
+	V x6 = T{3.00000000000000}*x5;
+	V x7 = -x6;
+	V x8 = u[0]*u[0];
+	V x9 = T{3.00000000000000}*x8;
+	V x10 = T{2.00000000000000} - x9;
+	V x11 = x10 + x7;
+	V x12 = x11 + x4;
+	V x13 = u[2]*u[2];
+	V x14 = T{3.00000000000000}*x13;
+	V x15 = -x14;
+	V x16 = T{6.00000000000000}*u[1];
+	V x17 = x15 + x16;
+	V x18 = T{6.00000000000000}*u[0];
+	V x19 = -u[0];
+	V x20 = x19 + u[2];
+	V x21 = x14 + x6 + T{-2.00000000000000};
+	V x22 = x21 + x9;
+	V x23 = x22 - x4;
+	V x24 = T{0.0277777777777778}*x0;
+	V x25 = T{6.00000000000000}*x13;
+	V x26 = u[0] + u[2];
+	V x27 = T{9.00000000000000}*(x26*x26);
+	V x28 = x15 + x18;
+	V x29 = -u[1];
+	V x30 = x29 + u[2];
+	V x31 = x19 + u[1];
+	V x32 = -x16;
+	V x33 = x18 + x22;
+	V x34 = T{6.00000000000000}*x5;
+	V x35 = u[0] + u[1];
+	V x36 = T{9.00000000000000}*(x35*x35);
+	V x37 = T{6.00000000000000}*x8;
+	V x38 = x9 + T{-2.00000000000000};
+	V x39 = x29 + u[0];
+	V x40 = -x18 + x22;
+	V x41 = -u[2];
+	V x42 = x41 + u[1];
+	V x43 = x22 + x4;
+	V x44 = x41 + u[0];
+	f_next[0] = x1*(rho*(x12 + x17 + x3) - T{72.0000000000000}*f_curr[0]) + f_curr[0];
+	f_next[1] = -x1*(rho*(x18 + x23 - T{9.00000000000000}*x20*x20) + T{72.0000000000000}*f_curr[1]) + f_curr[1];
+	f_next[2] = x24*(rho*(x12 + x25) - T{36.0000000000000}*f_curr[2]) + f_curr[2];
+	f_next[3] = x1*(rho*(x12 + x27 + x28) - T{72.0000000000000}*f_curr[3]) + f_curr[3];
+	f_next[4] = -x1*(rho*(x16 + x23 - T{9.00000000000000}*x30*x30) + T{72.0000000000000}*f_curr[4]) + f_curr[4];
+	f_next[5] = -x1*(rho*(x32 + x33 - T{9.00000000000000}*x31*x31) + T{72.0000000000000}*f_curr[5]) + f_curr[5];
+	f_next[6] = x24*(rho*(x10 + x17 + x34) - T{36.0000000000000}*f_curr[6]) + f_curr[6];
+	f_next[7] = x1*(rho*(x11 + x17 + x18 + x36) - T{72.0000000000000}*f_curr[7]) + f_curr[7];
+	f_next[8] = -x24*(rho*(x18 + x21 - x37) + T{36.0000000000000}*f_curr[8]) + f_curr[8];
+	f_next[9] = -T{0.166666666666667}*x0*(rho*x22 + T{6.00000000000000}*f_curr[9]) + f_curr[9];
+	f_next[10] = x24*(rho*(x28 + x37 + x7 + T{2.00000000000000}) - T{36.0000000000000}*f_curr[10]) + f_curr[10];
+	f_next[11] = -x1*(rho*(x16 + x33 - x36) + T{72.0000000000000}*f_curr[11]) + f_curr[11];
+	f_next[12] = -x24*(rho*(x14 + x16 - x34 + x38) + T{36.0000000000000}*f_curr[12]) + f_curr[12];
+	f_next[13] = -x1*(rho*(x16 + x40 - T{9.00000000000000}*x39*x39) + T{72.0000000000000}*f_curr[13]) + f_curr[13];
+	f_next[14] = -x1*(rho*(x32 + x43 - T{9.00000000000000}*x42*x42) + T{72.0000000000000}*f_curr[14]) + f_curr[14];
+	f_next[15] = -x1*(rho*(-x27 + x33 + x4) + T{72.0000000000000}*f_curr[15]) + f_curr[15];
+	f_next[16] = -x24*(rho*(-x25 + x38 + x4 + x6) + T{36.0000000000000}*f_curr[16]) + f_curr[16];
+	f_next[17] = -x1*(rho*(x4 + x40 - T{9.00000000000000}*x44*x44) + T{72.0000000000000}*f_curr[17]) + f_curr[17];
+	f_next[18] = -x1*(rho*(x16 - x3 + x43) + T{72.0000000000000}*f_curr[18]) + f_curr[18];
+}
+
+};