24#ifndef HYPERPLANE_LATTICE_3D_HH
25#define HYPERPLANE_LATTICE_3D_HH
33int HyperplaneLattice3D<T>::computeMaxLatticeDistance(Cuboid3D<T>&& cuboid)
const
35 const Vector<T,3> origin = cuboid.getOrigin();
36 const Vector<int,3> extend = cuboid.getExtent();
37 const T deltaR = cuboid.getDeltaR();
39 T maxPhysDistance = T();
44 for (
int iDim=0; iDim<3; ++iDim) {
45 tmpO[iDim] = origin[iDim] - _origin[iDim];
46 tmpE[iDim] = origin[iDim] + extend[iDim]*deltaR - _origin[iDim];
48 tmp =
util::sqrt(tmpO[0]*tmpO[0] + tmpO[1]*tmpO[1] + tmpO[2]*tmpO[2]);
49 if (maxPhysDistance < tmp) {
50 maxPhysDistance = tmp;
52 tmp =
util::sqrt((tmpE[0]*tmpE[0] + tmpO[1]*tmpO[1] + tmpO[2]*tmpO[2]));
53 if (maxPhysDistance < tmp) {
54 maxPhysDistance = tmp;
56 tmp =
util::sqrt(tmpO[0]*tmpO[0] + tmpE[1]*tmpE[1] + tmpO[2]*tmpO[2]);
57 if (maxPhysDistance < tmp) {
58 maxPhysDistance = tmp;
60 tmp =
util::sqrt(tmpO[0]*tmpO[0] + tmpO[1]*tmpO[1] + tmpE[2]*tmpE[2]);
61 if (maxPhysDistance < tmp) {
62 maxPhysDistance = tmp;
64 tmp =
util::sqrt(tmpO[0]*tmpO[0] + tmpE[1]*tmpE[1] + tmpE[2]*tmpE[2]);
65 if (maxPhysDistance < tmp) {
66 maxPhysDistance = tmp;
68 tmp =
util::sqrt(tmpE[0]*tmpE[0] + tmpO[1]*tmpO[1] + tmpE[2]*tmpE[2]);
69 if (maxPhysDistance < tmp) {
70 maxPhysDistance = tmp;
72 tmp =
util::sqrt(tmpE[0]*tmpE[0] + tmpE[1]*tmpE[1] + tmpO[2]*tmpO[2]);
73 if (maxPhysDistance < tmp) {
74 maxPhysDistance = tmp;
76 tmp =
util::sqrt(tmpE[0]*tmpE[0] + tmpE[1]*tmpE[1] + tmpE[2]*tmpE[2]);
77 if (maxPhysDistance < tmp) {
78 maxPhysDistance = tmp;
81 return int(maxPhysDistance/_h) + 1;
85void HyperplaneLattice3D<T>::constructCuboid(
86 CuboidGeometry3D<T>& geometry,
int maxLatticeDistance)
89 int minX = -maxLatticeDistance;
90 int maxX = maxLatticeDistance;
91 int minY = -maxLatticeDistance;
92 int maxY = maxLatticeDistance;
95 for (
int iX = -maxLatticeDistance; iX < maxLatticeDistance; ++iX ) {
96 for (
int iY = -maxLatticeDistance; iY < maxLatticeDistance; ++iY ) {
97 if ( geometry.getC(getPhysR(iX, iY), iC) ) {
108 for (
int iX = maxLatticeDistance; iX > -maxLatticeDistance; --iX ) {
109 for (
int iY = -maxLatticeDistance; iY < maxLatticeDistance; ++iY ) {
110 if ( geometry.getC(getPhysR(iX, iY), iC) ) {
121 for (
int iY = -maxLatticeDistance; iY < maxLatticeDistance; ++iY ) {
122 for (
int iX = -maxLatticeDistance; iX < maxLatticeDistance; ++iX ) {
123 if ( geometry.getC(getPhysR(iX, iY), iC) ) {
134 for (
int iY = maxLatticeDistance; iY > -maxLatticeDistance; --iY ) {
135 for (
int iX = -maxLatticeDistance; iX < maxLatticeDistance; ++iX ) {
136 if ( geometry.getC(getPhysR(iX, iY), iC) ) {
147 _nx = maxX - minX + 1;
148 _ny = maxY - minY + 1;
150 _origin = _origin + T(minX)*_u + T(minY)*_v;
154void HyperplaneLattice3D<T>::setToResolution(
int resolution)
157 T newH = _nx*_h/(T) resolution;
159 _ny = (int) (_ny*_h/newH) + 1;
163 T newH = _ny*_h/(T) resolution;
165 _nx = (int) (_nx*_h/newH) + 1;
175 : _hyperplane(hyperplane),
176 _origin(hyperplane.origin),
179 _h(geometry.getMinDeltaR())
184 const int maxLatticeDistance = computeMaxLatticeDistance(geometry.
getMotherCuboid());
186 constructCuboid(geometry, maxLatticeDistance);
192 : _hyperplane(hyperplane),
193 _origin(hyperplane.origin),
196 _h(geometry.getMinDeltaR())
201 const int maxLatticeDistance = computeMaxLatticeDistance(geometry.
getMotherCuboid());
203 constructCuboid(geometry, maxLatticeDistance);
205 if ( resolution > 0 ) {
206 setToResolution(resolution);
213 : _hyperplane(hyperplane),
214 _origin(hyperplane.origin),
226 const int maxLatticeDistance = computeMaxLatticeDistance(geometry.
getMotherCuboid());
228 constructCuboid(geometry, maxLatticeDistance);
235 : _hyperplane(hyperplane),
236 _origin(hyperplane.origin),
257 _origin[0] + T(planeX)*_u[0] + T(planeY)*_v[0],
258 _origin[1] + T(planeX)*_u[1] + T(planeY)*_v[1],
259 _origin[2] + T(planeX)*_u[2] + T(planeY)*_v[2]
A cuboid geometry represents a voxel mesh.
Cuboid3D< T > getMotherCuboid()
Returns the smallest cuboid that includes all cuboids of the structure.
T getMinDeltaR() const
Returns the minimum delta in the structure.
Vector< T, 3 > getVectorV() const
Vector< T, 3 > getPhysR(const int &planeX, const int &planeY) const
Transform 2d lattice coordinates to their physical 3d location.
Vector< T, 3 > _v
Span vector of the lattice, normalized to grid width _h.
Vector< T, 3 > getPhysOrigin() const
Vector< T, 3 > getVectorU() const
Vector< T, 3 > _u
Span vector of the lattice, normalized to grid width _h.
const Hyperplane3D< T > & getHyperplane() const
HyperplaneLattice3D(CuboidGeometry3D< T > &geometry, Hyperplane3D< T > hyperplane)
Constructor for automatic discretization.
T _h
Distance between discrete lattice points.
cpu::simd::Pack< T > sqrt(cpu::simd::Pack< T > value)
bool nearZero(const ADf< T, DIM > &a)
Top level namespace for all of OpenLB.
constexpr Vector< T, D > normalize(const ScalarVector< T, D, IMPL > &a, T scale=T{1})
Definition of a analytical 2D plane embedded in 3D space.