Added periodic boundconds to utils

This commit is contained in:
jpekkila
2019-10-15 17:37:57 +03:00
parent ff1ad37047
commit fd9dc7ca98
2 changed files with 67 additions and 0 deletions

View File

@@ -74,6 +74,71 @@ acMeshRandomize(AcMesh* mesh)
return AC_SUCCESS;
}
AcResult
acMeshApplyPeriodicBounds(AcMesh* mesh)
{
const AcMeshInfo info = mesh->info;
for (int w = 0; w < NUM_VTXBUF_HANDLES; ++w) {
const int3 start = (int3){0, 0, 0};
const int3 end = (int3){info.int_params[AC_mx], info.int_params[AC_my],
info.int_params[AC_mz]};
const int nx = info.int_params[AC_nx];
const int ny = info.int_params[AC_ny];
const int nz = info.int_params[AC_nz];
const int nx_min = info.int_params[AC_nx_min];
const int ny_min = info.int_params[AC_ny_min];
const int nz_min = info.int_params[AC_nz_min];
// The old kxt was inclusive, but our mx_max is exclusive
const int nx_max = info.int_params[AC_nx_max];
const int ny_max = info.int_params[AC_ny_max];
const int nz_max = info.int_params[AC_nz_max];
for (int k_dst = start.z; k_dst < end.z; ++k_dst) {
for (int j_dst = start.y; j_dst < end.y; ++j_dst) {
for (int i_dst = start.x; i_dst < end.x; ++i_dst) {
// If destination index is inside the computational domain, return since
// the boundary conditions are only applied to the ghost zones
if (i_dst >= nx_min && i_dst < nx_max && j_dst >= ny_min && j_dst < ny_max &&
k_dst >= nz_min && k_dst < nz_max)
continue;
// Find the source index
// Map to nx, ny, nz coordinates
int i_src = i_dst - nx_min;
int j_src = j_dst - ny_min;
int k_src = k_dst - nz_min;
// Translate (s.t. the index is always positive)
i_src += nx;
j_src += ny;
k_src += nz;
// Wrap
i_src %= nx;
j_src %= ny;
k_src %= nz;
// Map to mx, my, mz coordinates
i_src += nx_min;
j_src += ny_min;
k_src += nz_min;
const size_t src_idx = acVertexBufferIdx(i_src, j_src, k_src, info);
const size_t dst_idx = acVertexBufferIdx(i_dst, j_dst, k_dst, info);
ERRCHK(src_idx < acVertexBufferSize(info));
ERRCHK(dst_idx < acVertexBufferSize(info));
mesh->vertex_buffer[w][dst_idx] = mesh->vertex_buffer[w][src_idx];
}
}
}
}
return AC_SUCCESS;
}
AcResult
acMeshClear(AcMesh* mesh)
{

View File

@@ -39,6 +39,8 @@ AcResult acMeshSet(const AcReal value, AcMesh* mesh);
AcResult acMeshRandomize(AcMesh* mesh);
AcResult acMeshApplyPeriodicBounds(AcMesh* mesh);
AcResult acMeshClear(AcMesh* mesh);
#ifdef __cplusplus