LBPM/cpu/D3Q7.cpp
2023-10-23 04:18:20 -04:00

227 lines
9.2 KiB
C++

/*
Copyright 2013--2018 James E. McClure, Virginia Polytechnic & State University
Copyright Equnior ASA
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
// CPU Functions for D3Q7 Lattice Boltzmann Methods
extern "C" void ScaLBL_Scalar_Pack(int *list, int count, double *sendbuf,
double *Data, int N) {
//....................................................................................
// Pack distribution q into the send buffer for the listed lattice sites
// dist may be even or odd distributions stored by stream layout
//....................................................................................
int idx, n;
for (idx = 0; idx < count; idx++) {
n = list[idx];
sendbuf[idx] = Data[n];
}
}
extern "C" void ScaLBL_Scalar_Unpack(int *list, int count, double *recvbuf,
double *Data, int N) {
//....................................................................................
// Pack distribution q into the send buffer for the listed lattice sites
// dist may be even or odd distributions stored by stream layout
//....................................................................................
int idx, n;
for (idx = 0; idx < count; idx++) {
n = list[idx];
Data[n] = recvbuf[idx];
}
}
extern "C" void ScaLBL_D3Q7_Unpack(int q, int *list, int start, int count,
double *recvbuf, double *dist, int N) {
//....................................................................................
// Unack distribution from the recv buffer
// Distribution q matche Cqx, Cqy, Cqz
// swap rule means that the distributions in recvbuf are OPPOSITE of q
// dist may be even or odd distributions stored by stream layout
//....................................................................................
int n, idx;
for (idx = 0; idx < count; idx++) {
// Get the value from the list -- note that n is the index is from the send (non-local) process
n = list[idx];
// unpack the distribution to the proper location
if (!(n < 0))
dist[q * N + n] = recvbuf[start + idx];
//dist[q*N+n] = recvbuf[start+idx];
}
}
extern "C" void ScaLBL_PackDenD3Q7(int *list, int count, double *sendbuf,
int number, double *Data, int N) {
//....................................................................................
// Pack distribution into the send buffer for the listed lattice sites
//....................................................................................
int idx, n, component;
for (idx = 0; idx < count; idx++) {
for (component = 0; component < number; component++) {
n = list[idx];
sendbuf[idx * number + component] = Data[number * n + component];
Data[number * n + component] =
0.0; // Set the data value to zero once it's in the buffer!
}
}
}
extern "C" void ScaLBL_UnpackDenD3Q7(int *list, int count, double *recvbuf,
int number, double *Data, int N) {
//....................................................................................
// Unack distribution from the recv buffer
// Sum to the existing density value
//....................................................................................
int idx, n, component;
for (idx = 0; idx < count; idx++) {
for (component = 0; component < number; component++) {
n = list[idx];
Data[number * n + component] += recvbuf[idx * number + component];
}
}
}
extern "C" void ScaLBL_D3Q7_Reflection_BC_z(int *list, double *dist, int count,
int Np) {
int n;
for (int idx = 0; idx < count; idx++) {
n = list[idx];
double f5 = 0.222222222222222222222222 - dist[6 * Np + n];
dist[6 * Np + n] = f5;
}
}
extern "C" void ScaLBL_D3Q7_Reflection_BC_Z(int *list, double *dist, int count,
int Np) {
int n;
for (int idx = 0; idx < count; idx++) {
n = list[idx];
double f6 = 0.222222222222222222222222 - dist[5 * Np + n];
dist[5 * Np + n] = f6;
}
}
extern "C" void ScaLBL_D3Q7_Init(char *ID, double *f_even, double *f_odd,
double *Den, int Nx, int Ny, int Nz) {
int n, N;
N = Nx * Ny * Nz;
double value;
for (n = 0; n < N; n++) {
if (ID[n] > 0) {
value = Den[n];
f_even[n] = 0.3333333333333333 * value;
f_odd[n] = 0.1111111111111111 * value; //double(100*n)+1.f;
f_even[N + n] = 0.1111111111111111 * value; //double(100*n)+2.f;
f_odd[N + n] = 0.1111111111111111 * value; //double(100*n)+3.f;
f_even[2 * N + n] = 0.1111111111111111 * value; //double(100*n)+4.f;
f_odd[2 * N + n] = 0.1111111111111111 * value; //double(100*n)+5.f;
f_even[3 * N + n] = 0.1111111111111111 * value; //double(100*n)+6.f;
} else {
for (int q = 0; q < 3; q++) {
f_even[q * N + n] = -1.0;
f_odd[q * N + n] = -1.0;
}
f_even[3 * N + n] = -1.0;
}
}
}
//*************************************************************************
extern "C" void ScaLBL_D3Q7_Swap(char *ID, double *disteven, double *distodd,
int Nx, int Ny, int Nz) {
int i, j, k, n, nn, N;
// distributions
double f1, f2, f3, f4, f5, f6;
N = Nx * Ny * Nz;
for (n = 0; n < N; n++) {
//.......Back out the 3-D indices for node n..............
k = n / (Nx * Ny);
j = (n - Nx * Ny * k) / Nx;
i = n - Nx * Ny * k - Nz * j;
if (ID[n] > 0) {
//........................................................................
// Retrieve even distributions from the local node (swap convention)
// f0 = disteven[n]; // Does not particupate in streaming
f1 = distodd[n];
f3 = distodd[N + n];
f5 = distodd[2 * N + n];
//........................................................................
//........................................................................
// Retrieve odd distributions from neighboring nodes (swap convention)
//........................................................................
nn = n + 1; // neighbor index (pull convention)
if (!(i + 1 < Nx))
nn -= Nx; // periodic BC along the x-boundary
//if (i+1<Nx){
f2 = disteven[N + nn]; // pull neighbor for distribution 2
if (!(f2 < 0.0)) {
distodd[n] = f2;
disteven[N + nn] = f1;
}
//}
//........................................................................
nn = n + Nx; // neighbor index (pull convention)
if (!(j + 1 < Ny))
nn -= Nx * Ny; // Perioidic BC along the y-boundary
//if (j+1<Ny){
f4 = disteven[2 * N + nn]; // pull neighbor for distribution 4
if (!(f4 < 0.0)) {
distodd[N + n] = f4;
disteven[2 * N + nn] = f3;
// }
}
//........................................................................
nn = n + Nx * Ny; // neighbor index (pull convention)
if (!(k + 1 < Nz))
nn -= Nx * Ny * Nz; // Perioidic BC along the z-boundary
//if (k+1<Nz){
f6 = disteven[3 * N + nn]; // pull neighbor for distribution 6
if (!(f6 < 0.0)) {
distodd[2 * N + n] = f6;
disteven[3 * N + nn] = f5;
// }
}
}
}
}
//*************************************************************************
extern "C" void ScaLBL_D3Q7_Density(char *ID, double *disteven, double *distodd,
double *Den, int Nx, int Ny, int Nz) {
char id;
int n;
double f0, f1, f2, f3, f4, f5, f6;
int N = Nx * Ny * Nz;
for (n = 0; n < N; n++) {
id = ID[n];
if (id > 0) {
// Read the distributions
f0 = disteven[n];
f2 = disteven[N + n];
f4 = disteven[2 * N + n];
f6 = disteven[3 * N + n];
f1 = distodd[n];
f3 = distodd[N + n];
f5 = distodd[2 * N + n];
// Compute the density
Den[n] = f0 + f1 + f2 + f3 + f4 + f5 + f6;
}
}
}