Skip to content

Commit

Permalink
improve code repeating
Browse files Browse the repository at this point in the history
  • Loading branch information
@Franguuuu
Franguuuu committed Sep 30, 2024
1 parent a80ce96 commit 47e0563
Show file tree
Hide file tree
Showing 5 changed files with 81 additions and 122 deletions.
1 change: 0 additions & 1 deletion src/shared/common/sph_data_containers.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -76,7 +76,6 @@ typedef DataContainerAddressAssemble<SingleVariable> SingleVariables;

/** Generalized mesh data type */
// template <typename DataType>
// using MeshVariable = DiscreteVariable<DataType>;
typedef DataContainerAddressAssemble<MeshVariable> MeshVariableAssemble;

} // namespace SPH
Expand Down
130 changes: 50 additions & 80 deletions src/shared/geometries/level_set.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,101 +11,77 @@ MultilevelLevelSet::MultilevelLevelSet(
Shape &shape, SPHAdaptation &sph_adaptation)
: BaseMeshField("LevelSet_" + shape.getName()), kernel_(*sph_adaptation.getKernel()), shape_(shape), total_levels_(total_levels)
{
mesh_data_set_.push_back(
mesh_data_ptr_vector_keeper_
.template createPtr<MeshWithGridDataPackagesType>(tentative_bounds, reference_data_spacing, 4));

mesh_data_set_[0]->registerMeshVariable<Real>("Levelset");
mesh_data_set_[0]->registerMeshVariable<int>("NearInterfaceID");
mesh_data_set_[0]->registerMeshVariable<Vecd>("LevelsetGradient");
mesh_data_set_[0]->registerMeshVariable<Real>("KernelWeight");
mesh_data_set_[0]->registerMeshVariable<Vecd>("KernelGradient");

Real global_h_ratio = sph_adaptation.ReferenceSpacing() / reference_data_spacing;
global_h_ratio_vec_.push_back(global_h_ratio);
MeshAllDynamics<InitializeDataInACell> initialize_data_in_a_cell{*mesh_data_set_[0], shape_};
FinishDataPackages finish_data_packages{*mesh_data_set_[0], shape_, kernel_, global_h_ratio};
initialize_data_in_a_cell.exec();
finish_data_packages.exec();
probe_signed_distance_set_.push_back(
probe_signed_distance_vector_keeper_
.template createPtr<ProbeSignedDistance>(*mesh_data_set_[0]));
probe_normal_direction_set_.push_back(
probe_normal_direction_vector_keeper_
.template createPtr<ProbeNormalDirection>(*mesh_data_set_[0]));
probe_level_set_gradient_set_.push_back(
probe_level_set_gradient_vector_keeper_
.template createPtr<ProbeLevelSetGradient>(*mesh_data_set_[0]));

for (size_t level = 1; level != total_levels_; ++level)
{
reference_data_spacing *= 0.5;
global_h_ratio *= 2;
initializeLevel(0, reference_data_spacing, global_h_ratio, tentative_bounds);

for (size_t level = 1; level < total_levels_; ++level) {
reference_data_spacing *= 0.5; // Halve the data spacing
global_h_ratio *= 2; // Double the ratio
global_h_ratio_vec_.push_back(global_h_ratio);
/** all mesh levels aligned at the lower bound of tentative_bounds */
mesh_data_set_.push_back(
mesh_data_ptr_vector_keeper_
.template createPtr<MeshWithGridDataPackagesType>(tentative_bounds, reference_data_spacing, 4));
mesh_data_set_[level]->registerMeshVariable<Real>("Levelset");
mesh_data_set_[level]->registerMeshVariable<int>("NearInterfaceID");
mesh_data_set_[level]->registerMeshVariable<Vecd>("LevelsetGradient");
mesh_data_set_[level]->registerMeshVariable<Real>("KernelWeight");
mesh_data_set_[level]->registerMeshVariable<Vecd>("KernelGradient");
MeshAllDynamics<InitializeDataInACellFromCoarse> initialize_data_in_a_cell_from_coarse(*mesh_data_set_[level], *mesh_data_set_[level - 1], shape_);
FinishDataPackages finish_data_packages{*mesh_data_set_[level], shape_, kernel_, global_h_ratio};
initialize_data_in_a_cell_from_coarse.exec();
finish_data_packages.exec();

probe_signed_distance_set_.push_back(
probe_signed_distance_vector_keeper_
.template createPtr<ProbeSignedDistance>(*mesh_data_set_[level]));
probe_normal_direction_set_.push_back(
probe_normal_direction_vector_keeper_
.template createPtr<ProbeNormalDirection>(*mesh_data_set_[level]));
probe_level_set_gradient_set_.push_back(
probe_level_set_gradient_vector_keeper_
.template createPtr<ProbeLevelSetGradient>(*mesh_data_set_[level]));
initializeLevel(level, reference_data_spacing, global_h_ratio, tentative_bounds);
}

clean_interface = new CleanInterface{*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back()};
correct_topology = new CorrectTopology{*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back()};
clean_interface = makeUnique<CleanInterface>(*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back());
correct_topology = makeUnique<CorrectTopology>(*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back());
}
//=================================================================================================//
MultilevelLevelSet::MultilevelLevelSet(
BoundingBox tentative_bounds, MeshWithGridDataPackagesType* coarse_data, Shape &shape, SPHAdaptation &sph_adaptation)
: BaseMeshField("LevelSet_" + shape.getName()), kernel_(*sph_adaptation.getKernel()), shape_(shape)
: BaseMeshField("LevelSet_" + shape.getName()), kernel_(*sph_adaptation.getKernel()), shape_(shape), total_levels_(1)
{
total_levels_ = 1;
Real reference_data_spacing = coarse_data->DataSpacing() * 0.5;
Real global_h_ratio = sph_adaptation.ReferenceSpacing() / reference_data_spacing;
global_h_ratio_vec_.push_back(global_h_ratio);

initializeLevel(0, reference_data_spacing, global_h_ratio, tentative_bounds, coarse_data);

clean_interface = makeUnique<CleanInterface>(*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back());
correct_topology = makeUnique<CorrectTopology>(*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back());
}
//=================================================================================================//
void MultilevelLevelSet::initializeLevel(size_t level, Real reference_data_spacing, Real global_h_ratio, BoundingBox tentative_bounds, MeshWithGridDataPackagesType* coarse_data)
{
mesh_data_set_.push_back(
mesh_data_ptr_vector_keeper_
.template createPtr<MeshWithGridDataPackagesType>(tentative_bounds, reference_data_spacing, 4));
mesh_data_set_[0]->registerMeshVariable<Real>("Levelset");
mesh_data_set_[0]->registerMeshVariable<int>("NearInterfaceID");
mesh_data_set_[0]->registerMeshVariable<Vecd>("LevelsetGradient");
mesh_data_set_[0]->registerMeshVariable<Real>("KernelWeight");
mesh_data_set_[0]->registerMeshVariable<Vecd>("KernelGradient");

MeshAllDynamics<InitializeDataInACellFromCoarse> initialize_data_in_a_cell_from_coarse(*mesh_data_set_[0], *coarse_data, shape_);
FinishDataPackages finish_data_packages{*mesh_data_set_[0], shape_, kernel_, global_h_ratio};
initialize_data_in_a_cell_from_coarse.exec();
RegisterMeshVariable register_mesh_variable;
register_mesh_variable.exec(mesh_data_set_[level]);

if (coarse_data == nullptr) {
MeshAllDynamics<InitializeDataInACell> initialize_data_in_a_cell(*mesh_data_set_[level], shape_);
initialize_data_in_a_cell.exec();
} else {
MeshAllDynamics<InitializeDataInACellFromCoarse> initialize_data_in_a_cell_from_coarse(*mesh_data_set_[level], *coarse_data, shape_);
initialize_data_in_a_cell_from_coarse.exec();
}

FinishDataPackages finish_data_packages(*mesh_data_set_[level], shape_, kernel_, global_h_ratio);
finish_data_packages.exec();

registerProbes(level);
}
//=================================================================================================//
void MultilevelLevelSet::registerProbes(size_t level)
{
probe_signed_distance_set_.push_back(
probe_signed_distance_vector_keeper_
.template createPtr<ProbeSignedDistance>(*mesh_data_set_[0]));
probe_signed_distance_vector_keeper_
.template createPtr<ProbeSignedDistance>(*mesh_data_set_[level]));
probe_normal_direction_set_.push_back(
probe_normal_direction_vector_keeper_
.template createPtr<ProbeNormalDirection>(*mesh_data_set_[0]));
.template createPtr<ProbeNormalDirection>(*mesh_data_set_[level]));
probe_level_set_gradient_set_.push_back(
probe_level_set_gradient_vector_keeper_
.template createPtr<ProbeLevelSetGradient>(*mesh_data_set_[0]));

clean_interface = new CleanInterface{*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back()};
correct_topology = new CorrectTopology{*mesh_data_set_.back(), kernel_, global_h_ratio_vec_.back()};
.template createPtr<ProbeLevelSetGradient>(*mesh_data_set_[level]));
probe_kernel_integral_set_.push_back(
probe_kernel_integral_vector_keeper_
.template createPtr<ProbeKernelIntegral>(*mesh_data_set_[level]));
probe_kernel_gradient_integral_set_.push_back(
probe_kernel_gradient_integral_vector_keeper_
.template createPtr<ProbeKernelGradientIntegral>(*mesh_data_set_[level]));
}
//=================================================================================================//
size_t MultilevelLevelSet::getCoarseLevel(Real h_ratio)
Expand Down Expand Up @@ -158,33 +134,27 @@ size_t MultilevelLevelSet::getProbeLevel(const Vecd &position)
Real MultilevelLevelSet::probeKernelIntegral(const Vecd &position, Real h_ratio)
{
if(mesh_data_set_.size() == 1){
ProbeKernelIntegral refine_probe{*mesh_data_set_[0]};
return refine_probe.update(position);
return probe_kernel_integral_set_[0]->update(position);
}
size_t coarse_level = getCoarseLevel(h_ratio);
Real alpha = (global_h_ratio_vec_[coarse_level + 1] - h_ratio) /
(global_h_ratio_vec_[coarse_level + 1] - global_h_ratio_vec_[coarse_level]);
ProbeKernelIntegral coarse_probe{*mesh_data_set_[coarse_level]};
ProbeKernelIntegral fine_probe{*mesh_data_set_[coarse_level + 1]};
Real coarse_level_value = coarse_probe.update(position);
Real fine_level_value = fine_probe.update(position);
Real coarse_level_value = probe_kernel_integral_set_[coarse_level]->update(position);
Real fine_level_value = probe_kernel_integral_set_[coarse_level + 1]->update(position);

return alpha * coarse_level_value + (1.0 - alpha) * fine_level_value;
}
//=================================================================================================//
Vecd MultilevelLevelSet::probeKernelGradientIntegral(const Vecd &position, Real h_ratio)
{
if(mesh_data_set_.size() == 1){
ProbeKernelGradientIntegral refine_probe{*mesh_data_set_[0]};
return refine_probe.update(position);
return probe_kernel_gradient_integral_set_[0]->update(position);
}
size_t coarse_level = getCoarseLevel(h_ratio);
Real alpha = (global_h_ratio_vec_[coarse_level + 1] - h_ratio) /
(global_h_ratio_vec_[coarse_level + 1] - global_h_ratio_vec_[coarse_level]);
ProbeKernelGradientIntegral coarse_probe{*mesh_data_set_[coarse_level]};
ProbeKernelGradientIntegral fine_probe{*mesh_data_set_[coarse_level + 1]};
Vecd coarse_level_value = coarse_probe.update(position);
Vecd fine_level_value = fine_probe.update(position);
Vecd coarse_level_value = probe_kernel_gradient_integral_set_[coarse_level]->update(position);
Vecd fine_level_value = probe_kernel_gradient_integral_set_[coarse_level + 1]->update(position);

return alpha * coarse_level_value + (1.0 - alpha) * fine_level_value;
}
Expand Down
18 changes: 11 additions & 7 deletions src/shared/geometries/level_set.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -63,17 +63,17 @@ class MultilevelLevelSet : public BaseMeshField

void writeMeshFieldToPlt(std::ofstream &output_file) override
{
for (size_t l = 0; l != total_levels_; ++l)
{
WriteMeshFieldToPlt write_mesh_field_to_plt(*mesh_data_set_[l]);
write_mesh_field_to_plt.update(output_file);
}
for(size_t l = 0; l != total_levels_; ++l)
WriteMeshFieldToPlt(*mesh_data_set_[l]).update(output_file);
}

protected:
inline size_t getProbeLevel(const Vecd &position);
inline size_t getCoarseLevel(Real h_ratio);

void initializeLevel(size_t level, Real reference_data_spacing, Real global_h_ratio, BoundingBox tentative_bounds, MeshWithGridDataPackagesType* coarse_data = nullptr);
void registerProbes(size_t level);

Kernel &kernel_;
Shape &shape_; /**< the geometry is described by the level set. */
size_t total_levels_; /**< level 0 is the coarsest */
Expand All @@ -82,13 +82,17 @@ class MultilevelLevelSet : public BaseMeshField
StdVec<ProbeSignedDistance *> probe_signed_distance_set_;
StdVec<ProbeNormalDirection *> probe_normal_direction_set_;
StdVec<ProbeLevelSetGradient *> probe_level_set_gradient_set_;
StdVec<ProbeKernelIntegral *> probe_kernel_integral_set_;
StdVec<ProbeKernelGradientIntegral *> probe_kernel_gradient_integral_set_;
UniquePtrsKeeper<MeshWithGridDataPackagesType> mesh_data_ptr_vector_keeper_;
UniquePtrsKeeper<ProbeSignedDistance> probe_signed_distance_vector_keeper_;
UniquePtrsKeeper<ProbeNormalDirection> probe_normal_direction_vector_keeper_;
UniquePtrsKeeper<ProbeLevelSetGradient> probe_level_set_gradient_vector_keeper_;
UniquePtrsKeeper<ProbeKernelIntegral> probe_kernel_integral_vector_keeper_;
UniquePtrsKeeper<ProbeKernelGradientIntegral> probe_kernel_gradient_integral_vector_keeper_;

CleanInterface *clean_interface;
CorrectTopology *correct_topology;
UniquePtr<CleanInterface> clean_interface;
UniquePtr<CorrectTopology> correct_topology;
};
} // namespace SPH
#endif // LEVEL_SET_H
14 changes: 14 additions & 0 deletions src/shared/meshes/all_mesh_dynamics.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,20 @@

namespace SPH
{
class RegisterMeshVariable
{
public:
explicit RegisterMeshVariable(){};
~RegisterMeshVariable(){};

void exec(MeshWithGridDataPackagesType *mesh_data){
mesh_data->registerMeshVariable<Real>("Levelset");
mesh_data->registerMeshVariable<int>("NearInterfaceID");
mesh_data->registerMeshVariable<Vecd>("LevelsetGradient");
mesh_data->registerMeshVariable<Real>("KernelWeight");
mesh_data->registerMeshVariable<Vecd>("KernelGradient");
}
};
class FinishDataPackages : public BaseMeshDynamics
{
public:
Expand Down
40 changes: 6 additions & 34 deletions src/shared/meshes/mesh_with_data_packages.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -43,32 +43,6 @@ using namespace std::placeholders;

namespace SPH
{
/**
* @class BaseDataPackage
* @brief Abstract base class for a data package,
* by which the data in a derived class can be on- or off-grid.
* The data package can be defined in a cell of a background mesh so the pkg_index is
* the cell location on the mesh.
*/
class BaseDataPackage
{
public:
BaseDataPackage() : cell_index_on_mesh_(Arrayi::Zero()), state_indicator_(0){};
virtual ~BaseDataPackage(){};
void setInnerPackage() { state_indicator_ = 1; };
bool isInnerPackage() { return state_indicator_ != 0; };
void setCorePackage() { state_indicator_ = 2; };
bool isCorePackage() { return state_indicator_ == 2; };
void setCellIndexOnMesh(const Arrayi &cell_index) { cell_index_on_mesh_ = cell_index; }
Arrayi CellIndexOnMesh() const { return cell_index_on_mesh_; }

protected:
Arrayi cell_index_on_mesh_; /**< index of this data package on the background mesh, zero if it is not on the mesh. */
/** reserved value: 0 not occupying background mesh, 1 occupying.
* guide to use: larger for high priority of the data package. */
int state_indicator_;
};

/**
* @class MeshWithGridDataPackages
* @brief Abstract class for mesh with grid-based data packages.
Expand All @@ -85,7 +59,6 @@ class BaseDataPackage
* All these data packages are indexed by a concurrent vector inner_data_pkgs_.
* Note that a data package should be not near the mesh bound, otherwise one will encounter the error "out of range".
*/
// template <class GridDataPackageType>
template <int PKG_SIZE>
class MeshWithGridDataPackages : public Mesh
{
Expand Down Expand Up @@ -136,21 +109,20 @@ class MeshWithGridDataPackages : public Mesh
void deleteIndexDataMatrix(); /**< delete memories for metadata of data packages. */

/** resize all mesh variable data field with `num_grid_pkgs_` size(initially only singular data) */
template <typename DataType>
struct ResizeMeshVariableData
{
void operator()(MeshVariableAssemble &all_mesh_variables_,
template <typename DataType>
void operator()(DataContainerAddressKeeper<MeshVariable<DataType>> &all_mesh_variables_,
const size_t num_grid_pkgs_)
{
constexpr int type_index = DataTypeIndex<DataType>::value;
for (size_t l = 0; l != std::get<type_index>(all_mesh_variables_).size(); ++l)
for (size_t l = 0; l != all_mesh_variables_.size(); ++l)
{
MeshVariable<DataType> *variable = std::get<type_index>(all_mesh_variables_)[l];
MeshVariable<DataType> *variable = all_mesh_variables_[l];
variable->allocateAllMeshVariableData(num_grid_pkgs_);
}
}
};
DataAssembleOperation<ResizeMeshVariableData> resize_mesh_variable_data_;
OperationOnDataAssemble<MeshVariableAssemble, ResizeMeshVariableData> resize_mesh_variable_data_{all_mesh_variables_};

/** probe by applying bi and tri-linear interpolation within the package. */
template <class DataType>
Expand Down Expand Up @@ -204,7 +176,7 @@ class MeshWithGridDataPackages : public Mesh
void for_each_cell_data(const FunctionOnData &function);
void resizeMeshVariableData()
{
resize_mesh_variable_data_(all_mesh_variables_, num_grid_pkgs_);
resize_mesh_variable_data_(num_grid_pkgs_);
}

template <typename DataType>
Expand Down

0 comments on commit 47e0563

Please sign in to comment.