1 // =============================================================================
2 // PROJECT CHRONO - http://projectchrono.org
3 //
4 // Copyright (c) 2019 projectchrono.org
5 // All rights reserved.
6 //
7 // Use of this source code is governed by a BSD-style license that can be found
8 // in the LICENSE file at the top level of the distribution and at
9 // http://projectchrono.org/license-chrono.txt.
10 //
11 // =============================================================================
12 // Authors: Conlain Kelly
13 // =============================================================================
14 // Chrono::Gpu simulation of a rectangular bed of granular material which is
15 // first let to settle and then compressed by advancing one of the box walls
16 // into the material.
17 // =============================================================================
18 
19 #include <iostream>
20 #include <string>
21 #include <cmath>
22 
23 #include "chrono/core/ChGlobal.h"
24 #include "chrono/utils/ChUtilsSamplers.h"
25 
26 #include "chrono_gpu/ChGpuData.h"
27 #include "chrono_gpu/physics/ChSystemGpu.h"
28 #include "chrono_gpu/utils/ChGpuJsonParser.h"
29 
30 #include "chrono_thirdparty/filesystem/path.h"
31 
32 using namespace chrono;
33 using namespace chrono::gpu;
34 
main(int argc,char * argv[])35 int main(int argc, char* argv[]) {
36     ChGpuSimulationParameters params;
37 
38     // Some of the default values might be overwritten by user via command line
39     if (argc != 2 || ParseJSON(gpu::GetDataFile(argv[1]), params) == false) {
40         std::cout << "Usage:\n./demo_GPU_movingBoundary <json_file>" << std::endl;
41         return 1;
42     }
43 
44     // Setup simulation. A convenient thing we can do is to move the Big Box Domain by (X/2, Y/2, Z/2) which is done
45     // with the fourth constructor param, so the coordinate range we are now working with is (0,0,0) to (X,Y,Z), instead
46     // of (-X/2,-Y/2,-Z/2) to (X/2, Y/2, Z/2).
47     ChSystemGpu gpu_sys(params.sphere_radius, params.sphere_density,
48                         ChVector<float>(params.box_X, params.box_Y, params.box_Z),
49                         ChVector<float>(params.box_X / 2, params.box_Y / 2, params.box_Z / 2));
50 
51     gpu_sys.SetPsiFactors(params.psi_T, params.psi_L);
52 
53     gpu_sys.SetKn_SPH2SPH(params.normalStiffS2S);
54     gpu_sys.SetKn_SPH2WALL(params.normalStiffS2W);
55     gpu_sys.SetGn_SPH2SPH(params.normalDampS2S);
56     gpu_sys.SetGn_SPH2WALL(params.normalDampS2W);
57 
58     gpu_sys.SetKt_SPH2SPH(params.tangentStiffS2S);
59     gpu_sys.SetKt_SPH2WALL(params.tangentStiffS2W);
60     gpu_sys.SetGt_SPH2SPH(params.tangentDampS2S);
61     gpu_sys.SetGt_SPH2WALL(params.tangentDampS2W);
62     gpu_sys.SetStaticFrictionCoeff_SPH2SPH(params.static_friction_coeffS2S);
63     gpu_sys.SetStaticFrictionCoeff_SPH2WALL(params.static_friction_coeffS2W);
64 
65     gpu_sys.SetCohesionRatio(params.cohesion_ratio);
66     gpu_sys.SetAdhesionRatio_SPH2WALL(params.adhesion_ratio_s2w);
67 
68     gpu_sys.SetGravitationalAcceleration(ChVector<float>(params.grav_X, params.grav_Y, params.grav_Z));
69     gpu_sys.SetParticleOutputMode(params.write_mode);
70 
71     std::string out_dir = GetChronoOutputPath() + "GPU/";
72     filesystem::create_directory(filesystem::path(out_dir));
73     out_dir = out_dir + params.output_dir;
74     filesystem::create_directory(filesystem::path(out_dir));
75 
76     // The box to fill with particles
77     ChVector<float> hdims((float)(params.box_X / 2.0 - 1.2), (float)(params.box_Y / 2.0 - 1.2),
78                           (float)(params.box_Z / 10.0 - 1.2));
79     ChVector<float> center((float)(params.box_X / 2), (float)(params.box_Y / 2), (float)(params.box_Z / 10.0));
80 
81     // Fill box with bodies
82     std::vector<ChVector<float>> body_points =
83         utils::PDLayerSampler_BOX<float>(center, hdims, 2.f * params.sphere_radius, 1.05f);
84 
85     gpu_sys.SetParticles(body_points);
86 
87     // Set the position of the BD
88     gpu_sys.SetBDFixed(true);
89 
90     gpu_sys.SetTimeIntegrator(CHGPU_TIME_INTEGRATOR::FORWARD_EULER);
91     gpu_sys.SetFrictionMode(CHGPU_FRICTION_MODE::MULTI_STEP);
92     gpu_sys.SetFixedStepSize(params.step_size);
93 
94     gpu_sys.SetVerbosity(params.verbose);
95 
96     // start outside BD by 10 cm
97     ChVector<float> plane_pos(-10, 0, 0);
98     ChVector<float> plane_normal(1, 0, 0);
99 
100     size_t plane_bc_id = gpu_sys.CreateBCPlane(plane_pos, plane_normal, false);
101 
102     // Function prescibing the motion of the advancing plane.
103     // Begins outside of the domain.
104     std::function<double3(float)> plane_pos_func = [&params](float t) {
105         double3 pos = {0, 0, 0};
106 
107         // move at 10 cm/s
108         constexpr float vel = 10;
109 
110         // after 1 second the plane will be at the edge of the BD, and will continue in thereafter
111         pos.x = vel * t;
112 
113         return pos;
114     };
115 
116     gpu_sys.Initialize();
117 
118     gpu_sys.SetBCOffsetFunction(plane_bc_id, plane_pos_func);
119 
120     int fps = 50;
121     // assume we run for at least one frame
122     float frame_step = 1.0f / fps;
123     float curr_time = 0;
124     int currframe = 0;
125     unsigned int total_frames = (unsigned int)((float)params.time_end * fps);
126 
127     char filename[100];
128     sprintf(filename, "%s/step%06d.csv", out_dir.c_str(), currframe++);
129     gpu_sys.WriteParticleFile(std::string(filename));
130 
131     std::cout << "frame step is " << frame_step << std::endl;
132 
133     // Run settling experiments
134     while (curr_time < params.time_end) {
135         gpu_sys.AdvanceSimulation(frame_step);
136         curr_time += frame_step;
137         printf("rendering frame %u of %u\n", currframe, total_frames);
138         sprintf(filename, "%s/step%06d.csv", out_dir.c_str(), currframe++);
139         gpu_sys.WriteParticleFile(std::string(filename));
140     }
141 
142     return 0;
143 }
144