From 6da0d42e698312968d5a90d1ab5faac2377f35bf Mon Sep 17 00:00:00 2001
From: Behzad Safaei <iwia103h@a0522.nhr.fau.de>
Date: Sun, 4 May 2025 01:04:57 +0200
Subject: [PATCH] Add option for scaling the ballast app
---
examples/modular/ballast/ballast.cpp | 79 ++++++++++++++++++++--------
examples/modular/ballast/ballast.py | 3 +-
runtime/utility/read_from_file.cpp | 14 +++--
runtime/utility/read_from_file.hpp | 2 +-
4 files changed, 70 insertions(+), 28 deletions(-)
diff --git a/examples/modular/ballast/ballast.cpp b/examples/modular/ballast/ballast.cpp
index 20e55af..fbbf79c 100644
--- a/examples/modular/ballast/ballast.cpp
+++ b/examples/modular/ballast/ballast.cpp
@@ -33,17 +33,26 @@ void set_feature_properties(std::shared_ptr<PairsAccessor> &ac){
ac->syncTypeDampingTan();
}
+void print_usage(){
+ std::cerr << "Required args: num_beds_x num_beds_y num_timesteps" << std::endl;
+ exit(-1);
+}
+
int main(int argc, char **argv) {
+ if (argc<4) print_usage();
+
auto pairs_sim = std::make_shared<PairsSimulation>();
pairs_sim->initialize();
auto ac = std::make_shared<PairsAccessor>(pairs_sim.get());
auto pairs_runtime = pairs_sim->getPairsRuntime();
// ===================================================================================================
- // Domain
+ // Domain
// ===================================================================================================
- double domain_size[3] = {1.5, 1.0, 0.8};
- pairs_runtime->initDomain(&argc, &argv, 0, 0, 0, domain_size[0], domain_size[1], domain_size[2], false);
+ int nbeds[2] = {std::stoi(argv[1]), std::stoi(argv[2])};
+ double bed_size[3] = {1.5, 1.0, 0.65};
+ double domain_size[3] = {bed_size[0]*nbeds[0], bed_size[1]*nbeds[1], bed_size[2]};
+ pairs_runtime->initDomain(&argc, &argv, 0, 0, 0, domain_size[0], domain_size[1], domain_size[2], false);
// ===================================================================================================
// Spheres
@@ -61,7 +70,14 @@ int main(int argc, char **argv) {
// Settled bed with box in place
// ------------------------------------------------------------------
- pairs::read_spheres(pairs_runtime, "../examples/modular/ballast/spheres_bed_with_box.txt");
+ for(int i=0; i<nbeds[0]; ++i){
+ for(int j=0; j<nbeds[1]; ++j){
+ pairs::read_spheres(pairs_runtime,
+ "../examples/modular/ballast/spheres_bed_with_box.txt",
+ {i*bed_size[0], j*bed_size[1], 0.0});
+
+ }
+ }
// ===================================================================================================
// Box (sleeper)
@@ -70,10 +86,21 @@ int main(int argc, char **argv) {
double box_density = 8000;
double box_size[3] = {0.4, 0.3, 0.15};
double box_pos[3] = {0.779191, 0.497293, 0.484325};
- auto box_uid = pairs::create_box(pairs_runtime, box_pos[0], box_pos[1], box_pos[2],
- 0, 0, 0,
- box_size[0], box_size[1], box_size[2],
- box_density, 1, pairs::flags::GLOBAL);
+ int num_boxes = nbeds[0] * nbeds[1];
+ pairs::id_t box_uid[num_boxes];
+ for(int i=0; i<nbeds[0]; ++i){
+ for(int j=0; j<nbeds[1]; ++j){
+ box_uid[i*nbeds[1] + j] = pairs::create_box(pairs_runtime,
+ box_pos[0] + i*bed_size[0],
+ box_pos[1] + j*bed_size[1],
+ box_pos[2],
+ 0, 0, 0,
+ box_size[0], box_size[1], box_size[2],
+ box_density, 1, pairs::flags::GLOBAL);
+ }
+ }
+
+ pairs::Vector3<double> external_force(100, 0, 0);
// ===================================================================================================
// Halfspace (bottom plate)
@@ -92,7 +119,7 @@ int main(int argc, char **argv) {
pairs_sim->update_mass_and_inertia();
pairs_sim->updateDomain();
- int num_timesteps = 500000;
+ int num_timesteps = std::stoi(argv[3]);
int vtk_freq = 5000;
double dt = 1e-5;
set_feature_properties(ac);
@@ -127,14 +154,15 @@ int main(int argc, char **argv) {
// Apply external force to Box
// ------------------------------------------------------------------
- auto box_idx = ac->uidToIdxLocal(box_uid);
- ac->syncForce(PairsAccessor::Host);
- ac->syncTorque(PairsAccessor::Host);
- pairs::Vector3<double> external_force(500, 0, 0);
- ac->setForce(box_idx, ac->getForce(box_idx) + external_force);
- ac->setTorque(box_idx, {0.0, 0.0, 0.0}); // Reset torque to prevent rotation
- ac->syncForce(PairsAccessor::Host);
- ac->syncTorque(PairsAccessor::Host);
+ for(int i=0; i<num_boxes; ++i){
+ auto box_idx = ac->uidToIdxLocal(box_uid[i]);
+ ac->syncForce(PairsAccessor::Host);
+ ac->syncTorque(PairsAccessor::Host);
+ ac->setForce(box_idx, ac->getForce(box_idx) + external_force);
+ ac->setTorque(box_idx, {0.0, 0.0, 0.0}); // Reset torque to prevent rotation
+ ac->syncForce(PairsAccessor::Host);
+ ac->syncTorque(PairsAccessor::Host);
+ }
// Update positions and reneighbor
// ------------------------------------------------------------------
@@ -147,9 +175,15 @@ int main(int argc, char **argv) {
pairs::vtk_write_data(pairs_runtime, "output/local_spheres", 0, pairs_sim->nlocal(), t);
if(rank==0){
+ auto runtime = std::chrono::duration<double>(std::chrono::high_resolution_clock::now() - start).count();
+ std::cout << "Timestep: " << t << " - Time = " << dt*t << " - Runtime = " << runtime << std::endl;
ac->syncPosition(PairsAccessor::Host);
- std::cout << "Timestep: " << t << " - Time = " << dt*t << " - Box position = " << ac->getPosition(box_idx) << std::endl;
+ for(int i=0; i<num_boxes; ++i){
+ auto box_idx = ac->uidToIdxLocal(box_uid[i]);
+ std::cout << "\t Box (" << i << ") position = " << ac->getPosition(box_idx) << std::endl;
+ }
pairs::vtk_with_rotation(pairs_runtime, pairs::Shapes::Box, "output/local_boxes", 0, pairs_sim->nlocal(), t);
+
}
}
}
@@ -172,12 +206,15 @@ int main(int argc, char **argv) {
std::cout << "PUPS: " << pups << std::endl;
ac->syncUid(PairsAccessor::Host);
- auto plate_idx = ac->uidToIdxLocal(plate_uid);
- auto box_idx = ac->uidToIdxLocal(box_uid);
+ for(int i=0; i<num_boxes; ++i){
+ auto box_idx = ac->uidToIdxLocal(box_uid[i]);
+ std::cout << "box (" << i << ") final position = " << ac->getPosition(box_idx) << std::endl;
+ }
- std::cout << "box final position = " << ac->getPosition(box_idx) << std::endl;
+ auto plate_idx = ac->uidToIdxLocal(plate_uid);
std::cout << "plate final position = " << ac->getPosition(plate_idx) << std::endl;
}
+ pairs::log_timers(pairs_runtime);
ac->end();
pairs_sim->end();
diff --git a/examples/modular/ballast/ballast.py b/examples/modular/ballast/ballast.py
index 2950b80..0fb1571 100644
--- a/examples/modular/ballast/ballast.py
+++ b/examples/modular/ballast/ballast.py
@@ -99,7 +99,8 @@ psim.add_feature_property('type', 'damping_tan', pairs.real())
psim.add_feature_property('type', 'friction', pairs.real())
# psim.set_domain_partitioner(pairs.block_forest())
-psim.set_domain_partitioner(pairs.regular_domain_partitioner())
+# psim.set_domain_partitioner(pairs.regular_domain_partitioner())
+psim.set_domain_partitioner(pairs.regular_domain_partitioner_xy())
psim.pbc([True, True, False])
psim.build_cell_lists()
diff --git a/runtime/utility/read_from_file.cpp b/runtime/utility/read_from_file.cpp
index c19695f..7e325a8 100644
--- a/runtime/utility/read_from_file.cpp
+++ b/runtime/utility/read_from_file.cpp
@@ -139,7 +139,7 @@ void write_spheres(PairsRuntime *pr, const char *filename){
MPI_Barrier(MPI_COMM_WORLD);
}
-void read_spheres(PairsRuntime *pr, const char *filename){
+void read_spheres(PairsRuntime *pr, const char *filename, std::array<double, 3> offset){
int n = pr->getTrackedVariableAsInteger("nlocal");
auto shapes = pr->getAsIntegerProperty(pr->getPropertyByName("shape"));
auto flags = pr->getAsIntegerProperty(pr->getPropertyByName("flags"));
@@ -168,16 +168,20 @@ void read_spheres(PairsRuntime *pr, const char *filename){
double mass;
in_stream >> shape >> flag >> type >> x >> y >> z >> radius >> mass;
+ double shifted_posx = x + offset[0];
+ double shifted_posy = y + offset[1];
+ double shifted_posz = z + offset[2];
+
if(shape == Shapes::Sphere){
- bool within_subdom = pr->getDomainPartitioner()->isWithinSubdomain(x, y, z);
+ bool within_subdom = pr->getDomainPartitioner()->isWithinSubdomain(shifted_posx, shifted_posy, shifted_posz);
if(within_subdom || (flag & (flags::INFINITE | flags::GLOBAL))){
shapes(n) = shape;
flags(n) = flag;
types(n) = type;
uids(n) = (flag & (flags::INFINITE |flags::GLOBAL)) ? UniqueID::createGlobal(pr) : UniqueID::create(pr);
- positions(n,0) = x;
- positions(n,1) = y;
- positions(n,2) = z;
+ positions(n,0) = shifted_posx;
+ positions(n,1) = shifted_posy;
+ positions(n,2) = shifted_posz;
radii(n) = radius;
masses(n) = mass;
++n;
diff --git a/runtime/utility/read_from_file.hpp b/runtime/utility/read_from_file.hpp
index f1cb633..ba7d681 100644
--- a/runtime/utility/read_from_file.hpp
+++ b/runtime/utility/read_from_file.hpp
@@ -14,7 +14,7 @@ namespace pairs {
void write_boxes(PairsRuntime *pr, const char *filename);
void read_boxes(PairsRuntime *pr, const char *filename);
void write_spheres(PairsRuntime *pr, const char *filename);
-void read_spheres(PairsRuntime *pr, const char *filename);
+void read_spheres(PairsRuntime *pr, const char *filename, std::array<double, 3> offset = {0.0, 0.0, 0.0});
void read_grid_data(PairsRuntime *ps, const char *filename, real_t *grid_buffer);
--
GitLab