Hi Liam, You know your needs the best, how you do it depends on you.
But again, to get a densely packed configuration, you'll usually likely need to run the simulation and let it settle. Looks like from the picture you still have initial mesh-particle penetration, so running it probably won't go well. For me, it won't be difficult to carefully select a couple of regions in the space so after spawning the particles, there won't be initial penetration; or, I feel the easiest way is still sampling a big box region without the mesh, then add the mesh into the simulation scene, even if it involves simulating sticking the mesh into the pile of particles, which is not too difficult. Thank you, Ruochun On Tuesday, November 5, 2024 at 1:51:21 PM UTC+8 [email protected] wrote: > Thanks for the response! > > I updated the code to this as reccomended and ran a lot of different > scenarios and began to think that maybe I don't need the basin and could > have a particle domain that just surrounds the shell as the image attached > shows. I've been really trying to show some decent progress since we've > last talked. Now I just need to densely pack the particles and decrease the > size eventually to start conveying some stuff! Then maybe trying to write > the code for preventing it spawning in mesh. Let me know what you think > > #include <DEM/API.h> > > #include <DEM/HostSideHelpers.hpp> > > #include <DEM/utils/Samplers.hpp> > > #include <chrono> > > #include <iostream> > > #include <filesystem> > > #include <cmath> // For sqrt function > > > using namespace deme; > > using namespace std::chrono; > > namespace fs = std::filesystem; > > > int main() { > > // Initialize the DEM solver > > DEMSolver DEMSim; > > DEMSim.SetVerbosity(INFO); > > DEMSim.SetOutputFormat(OUTPUT_FORMAT::CSV); > > DEMSim.SetOutputContent(OUTPUT_CONTENT::ABSV); > > DEMSim.SetMeshOutputFormat(MESH_FORMAT::VTK); > > > // Material properties with reduced stiffness and increased damping > > auto mat_type_mixer = DEMSim.LoadMaterial({ > > {"E", 1e5f}, // Reduced stiffness to prevent high contact forces > > {"nu", 0.3f}, > > {"CoR", 0.2f}, > > {"mu", 0.3f}, > > {"Crr", 0.2f} // Increased rolling resistance > > }); > > > // Adjusted simulation domain dimensions > > float domain_x = 3.0f; // Increased domain to accommodate spawning area > > float domain_y = 3.0f; > > float domain_z = 3.0f; > > DEMSim.InstructBoxDomainDimension(domain_x, domain_y, domain_z); > > DEMSim.InstructBoxDomainBoundingBC("none", mat_type_mixer); // We'll > add custom walls instead > > > // Load mesh objects (Auger and Drum) without scaling > > // Commented out the basin as per user request > > /* > > auto basin = DEMSim.AddWavefrontMeshObject( > > > "/home/butterco/DEM-Engine/data/mesh/augertest_shapes/body_1_1_collision.obj", > > mat_type_mixer > > ); > > basin->SetFamily(0); > > DEMSim.SetFamilyFixed(0); // Basin is fixed > > */ > > > auto auger = DEMSim.AddWavefrontMeshObject( > > > "/home/butterco/DEM-Engine/data/mesh/augertest_shapes/body_2_1_collision.obj", > > mat_type_mixer > > ); > > auto drum = DEMSim.AddWavefrontMeshObject( > > > "/home/butterco/DEM-Engine/data/mesh/augertest_shapes/body_3_1_collision.obj", > > mat_type_mixer > > ); > > > > > drum->Move(make_float3(0, 0, -0.6f), make_float4(0, 0.7071f, 0.7071f, > 0)); > > auger->Move(make_float3(0, 0, 0.85f), make_float4(1, 0, 0, 0)); > > > // Set family numbers for simulation properties > > // basin->SetFamily(0); // Basin is commented out > > auger->SetFamily(1); > > drum->SetFamily(2); > > > // Fix the auger > > DEMSim.SetFamilyFixed(1); // Auger is completely fixed > > > > > // DEMSim.SetFamilyFixed(2); // Removed to allow drum rotation > > > // Assign prescribed angular velocity to the drum > > DEMSim.SetFamilyPrescribedAngVel(2, "0", "0", "0.5"); // Drum rotates > around Z-axis at 0.5 rad/s > > > // Define particle template with decreased radius > > float particle_radius = 0.01f; // Particle radius decreased to 10 mm > > float mass = 2.6e3f * (4.0f / 3.0f) * M_PI * pow(particle_radius, 3); > // Mass in kg (~0.1088 kg) > > float3 MOI = (2.0f / 5.0f) * mass * particle_radius * particle_radius > * make_float3(1.0f, 1.0f, 1.0f); // MOI for a sphere (~4.352e-5 kg·m²) > > auto particle_template = DEMSim.LoadSphereType(particle_radius, mass, > mat_type_mixer); > > > // Sampler setup to spawn particles around the drum and auger > > float drum_radius = 0.02f; // Drum radius in meters > > > // Reverted inner and outer radii > > float inner_radius = 0.036f; // Original inner radius (3.6 cm) > > float outer_radius = 0.15f; // Reduced outer radius to 0.15 meters to > match particle domain > > > // Adjusted fill height and fill bottom based on spawning volume > calculation > > float fill_height = 0.2f; // Fill height in meters (20 cm) > > float fill_bottom = -0.7f; // Bottom of the fill region in meters > > float3 fill_center = make_float3(0.0f, 0.0f, fill_bottom + fill_height > / 2.0f); // (0, 0, -0.6) > > > // Debugging: Print inner and outer radii and fill height > > std::cout << "Inner Radius: " << inner_radius << " meters" << > std::endl; > > std::cout << "Outer Radius: " << outer_radius << " meters" << > std::endl; > > std::cout << "Fill Height: " << fill_height << " meters" << std::endl; > > > // Adjust sampler spacing to prevent overlaps and slightly increase > particle count > > HCPSampler sampler(2.8f * particle_radius); // Sampler spacing set to > 0.028m (2.8x particle_radius) > > > // Sample particles in the cylinder with outer radius > > auto input_xyz_all = sampler.SampleCylinderZ(fill_center, > outer_radius, fill_height); > > > // Filter out particles within the inner radius to create a donut shape > > std::vector<float3> input_xyz; > > for (const auto& pos : input_xyz_all) { > > float dx = pos.x - fill_center.x; > > float dy = pos.y - fill_center.y; > > float dist = sqrt(dx * dx + dy * dy); > > if (dist >= inner_radius && dist <= outer_radius) { > > input_xyz.push_back(pos); > > } > > } > > > std::cout << "Total number of particles generated: " << > input_xyz.size() << std::endl; > > > if (!input_xyz.empty()) { > > // Add clumps and get the clump batch > > auto clump_batch = DEMSim.AddClumps(particle_template, input_xyz); > > clump_batch->SetFamily(3); // Assign particles to family 3, which > is not fixed > > } else { > > std::cerr << "Error: No particles were generated. Adjust sampling > parameters." << std::endl; > > return -1; > > } > > > // Adjust the initial time step size > > DEMSim.SetInitTimeStep(0.5e-6f); // Decrease time step for higher > stability > > > // Set gravity and other simulation parameters > > DEMSim.SetGravitationalAcceleration(make_float3(0.0f, 0.0f, -9.81f)); > > DEMSim.SetCDUpdateFreq(40); > > DEMSim.SetExpandSafetyAdder(2.0f); > > > // Increase the maximum allowed average contacts per sphere > > DEMSim.SetErrorOutAvgContacts(300); // Increased from 250 to 300 > > > // Increase the allowed system max velocity > > DEMSim.SetErrorOutVelocity(3000.0f); // Increased from 2000 to 3000 > > > DEMSim.SetAdaptiveBinSizeUpperProactivity(0.5f); > > DEMSim.SetAdaptiveBinSizeLowerProactivity(0.15f); > > > // Set initial bin size to prevent bins from becoming overpopulated > > DEMSim.SetInitBinSize(25.0f * particle_radius); // 25 * 0.01 = 0.25m > > > // Increase maximum spheres per bin if necessary > > DEMSim.SetMaxSphereInBin(50000); > > > // Add walls to prevent particles from falling over the edge > > auto walls = DEMSim.AddExternalObject(); > > walls->AddPlane(make_float3(0, 0, -0.7f), make_float3(0, 0, 1), > mat_type_mixer); // Bottom plane at Z=-0.7m > > walls->AddCylinder( > > make_float3(0, 0, 0), > > make_float3(0, 0, 1), > > outer_radius, // 0.15m > > mat_type_mixer > > ); // Outer cylinder > > walls->SetFamily(4); // Assign walls to a new family > > DEMSim.SetFamilyFixed(4); // Fix the walls > > > // Initialize the simulation > > DEMSim.Initialize(); > > > // Create output directory > > std::string output_dir = "/home/butterco/DEM-Engine/data/output/"; > > fs::create_directories(output_dir); > > > auto max_velocity_finder = DEMSim.CreateInspector("clump_max_absv"); > > > // Extend the simulation duration > > float sim_end = 20.0f; // Simulation end time in seconds > > unsigned int fps = 20; // Frames per second > > float frame_time = 1.0f / fps; > > unsigned int curr_frame = 0; > > > // Simulation loop > > for (float t = 0.0f; t < sim_end; t += frame_time) { > > std::cout << "Simulating frame " << curr_frame << std::endl; > > > char filename[200], mesh_filename[200]; > > sprintf(filename, "%sDEM_output_%04d.csv", output_dir.c_str(), > curr_frame); > > sprintf(mesh_filename, "%sDEM_mesh_output_%04d.vtk", > output_dir.c_str(), curr_frame++); > > DEMSim.WriteSphereFile(std::string(filename)); > > DEMSim.WriteMeshFile(std::string(mesh_filename)); > > > DEMSim.DoDynamics(frame_time); > > > float max_velocity = max_velocity_finder->GetValue(); > > std::cout << "Max velocity of any particle: " << max_velocity << > std::endl; > > } > > > DEMSim.ShowTimingStats(); > > std::cout << "Simulation complete." << std::endl; > > > return 0; > > } > best, > > Liam > > On Mon, Nov 4, 2024, 9:11 PM Ruochun Zhang <[email protected]> wrote: > >> Hi Liam, >> >> We see your simulation crashed at the first step, which shows problems >> with the initial setup. >> >> In your script, you didn't scale the clump info in the 3_clump.csv file, >> and left a note saying there's no need. Did you modify the file or is it >> still the same file bundled with DEM-Engine? If it's the same, then the >> clump has a size ≈ 2, as large as your simulation domain. This can explain >> the error you had, since if 70000 such particles overlap together, no >> matter how you divide the domain, a bin will intersect with most of them. >> >> In the picture you did not render the particles, but why? You need to >> visualize them to understand the problem. >> >> About the utility that prevents particles from spawning in the mesh, I >> think it's a nice-to-have and won't be too difficult to implement, but it >> is currently not our priority to make it a part of the package. You could >> easily write it yourself using C++ or Python though. DEM-Engine has support >> for smartly disabling these particles before the simulation starts, but >> it's an advanced usage. For now I suggest sampling the particles at safe >> locations and then letting them settle. >> >> Thank you, >> Ruochun >> >> On Tuesday, November 5, 2024 at 4:15:08 AM UTC+8 [email protected] >> wrote: >> >>> the thing im having trouble seeing possible is the mixer demo spawns the >>> particles above the mesh when the demo initializes and im not sure compared >>> to programs like StarCCM+ or ansys Rocky where there are options for >>> assisting in it not spawning within the mesh. would it be something of the >>> sort of making a solver within the code given the collision file to have >>> the particles not contact the mesh, but im not sure how one would code that >>> exactly, im just putting my thoughts out there. or is it better to spawn a >>> lot of particles spawn in above the mesh and wait for all of them to >>> settle >>> >>> On Monday, November 4, 2024 at 1:04:40 PM UTC-6 Liam Murray wrote: >>> >>>> I have tried to make progress on this simulation from the last time i >>>> posted on here. I got the shell to spin successfully and the meshes lined >>>> up thanks to feedback. What im trying to do is fill up the basin and shell >>>> meshes sort of like the mixer demo so that it is filled all of the way >>>> when >>>> the simulation initializes to start getting data and not having to wait >>>> for >>>> the particles to convey to the top. Ive tried messing with a lot of these >>>> parameters, but im just continuing to get errors like " outputs this >>>> error: >>>> I know setting max sphere in bin isnt the greatest practice or max >>>> contacts, but i keep running into this again and again and i wish i could >>>> spawn it in cleanly into the basin. the image attached is what the setup >>>> looks like without particles currently. >>>> >>>> >>>> butterco@NASALaptop:~/DEM-Engine/build/bin$ ./augertest Total number >>>> of particles generated: 71996 WARNING! One GPU device is detected. On >>>> consumer cards, DEME's performance edge is limited with only one GPU. Try >>>> allocating 2 GPU devices if possible. WARNING! Some clumps do not have >>>> their family numbers specified, so defaulted to 0 Number of total active >>>> devices: 1 User-specified X-dimension range: [-0.75, 0.75] User-specified >>>> Y-dimension range: [-0.75, 0.75] User-specified Z-dimension range: [-1, 1] >>>> User-specified dimensions should NOT be larger than the following >>>> simulation world. The dimension of the simulation world: >>>> 1.7999999523162842, 1.7999999523162842, 3.5999999046325684 Simulation >>>> world >>>> X range: [-0.9, 0.9] Simulation world Y range: [-0.9, 0.9] Simulation >>>> world >>>> Z range: [-1.2, 2.4] The length unit in this simulation is: >>>> 1.3096723358585471e-11 The edge length of a voxel: 8.5830686202825746e-07 >>>> The initial time step size: 5e-06 The initial edge length of a bin: >>>> 0.021110625075200021 The initial number of bins: 1264716 The total number >>>> of clumps: 71996 The combined number of component spheres: 215988 The >>>> total >>>> number of analytical objects: 1 The total number of meshes: 3 Grand total >>>> number of owners: 72000 The number of material types: 1 History-based >>>> Hertzian contact model is in use. The solver to set to adaptively change >>>> the contact margin size. Bin 884 contains 35373 sphere components, >>>> exceeding maximum allowance (32768). If you want the solver to run despite >>>> this, set allowance higher via SetMaxSphereInBin before simulation starts. >>>> -------- Simulation crashed "potentially" due to too many geometries in a >>>> bin -------- The dT reported max velocity is 0 >>>> ------------------------------------ If the velocity is huge, then the >>>> simulation probably diverged due to encountering large particle >>>> velocities. >>>> Decreasing the step size could help, and remember to check if your >>>> simulation objects are initially within the domain you specified. >>>> ------------------------------------ If the velocity is fair, and you are >>>> using a custom force model, one thing to do is to >>>> SetForceCalcThreadsPerBlock to a small number like 128 (see README.md >>>> troubleshooting for details). If you are going to discuss this on forum >>>> https://groups.google.com/g/projectchrono, please include a visual >>>> rendering of the simulation before crash. terminate called after throwing >>>> an instance of 'std::runtime_error' what(): GPU Assertion: unspecified >>>> launch failure. This happened in >>>> /home/butterco/DEM-Engine/src/algorithms/DEMCubContactDetection.cu:384 >>>> Aborted (core dumped) butterco@NASALaptop:~/DEM-Engine/build/bin$" >>>> >>>> " >>>> #include <DEM/API.h> >>>> #include <DEM/HostSideHelpers.hpp> >>>> #include <DEM/utils/Samplers.hpp> >>>> #include <chrono> >>>> #include <iostream> >>>> #include <filesystem> >>>> >>>> using namespace deme; >>>> using namespace std::chrono; >>>> namespace fs = std::filesystem; >>>> >>>> int main() { >>>> DEMSolver DEMSim; >>>> DEMSim.SetVerbosity(INFO); >>>> DEMSim.SetOutputFormat(OUTPUT_FORMAT::CSV); >>>> DEMSim.SetOutputContent(OUTPUT_CONTENT::ABSV); >>>> DEMSim.SetMeshOutputFormat(MESH_FORMAT::VTK); >>>> >>>> // Material properties >>>> auto mat_type_mixer = DEMSim.LoadMaterial({{"E", 5e5}, {"nu", 0.3}, >>>> {"CoR", 0.2}, {"mu", 0.3}, {"Crr", 0.0}}); >>>> >>>> // Set simulation domain dimensions >>>> DEMSim.InstructBoxDomainDimension(1.5, 1.5, 2.0); >>>> DEMSim.InstructBoxDomainBoundingBC("all", mat_type_mixer); >>>> >>>> // Load mesh objects (Basin, Auger, Drum) >>>> auto basin = >>>> DEMSim.AddWavefrontMeshObject("/home/butterco/DEM-Engine/data/mesh/augertest_shapes/body_1_1_collision.obj", >>>> >>>> mat_type_mixer); >>>> auto auger = >>>> DEMSim.AddWavefrontMeshObject("/home/butterco/DEM-Engine/data/mesh/augertest_shapes/body_2_1_collision.obj", >>>> >>>> mat_type_mixer); >>>> auto drum = >>>> DEMSim.AddWavefrontMeshObject("/home/butterco/DEM-Engine/data/mesh/augertest_shapes/body_3_1_collision.obj", >>>> >>>> mat_type_mixer); >>>> >>>> // Set initial positions and orientations >>>> basin->Move(make_float3(-0.1, -0.1, -0.70), make_float4(0, 0, 0, >>>> 1)); >>>> drum->Move(make_float3(0, 0, -0.6), make_float4(0, 0.7071, 0.7071, >>>> 0)); >>>> auger->Move(make_float3(0, 0, 0.85), make_float4(1, 0, 0, 0)); >>>> >>>> // Set family numbers for simulation properties >>>> basin->SetFamily(0); >>>> auger->SetFamily(1); >>>> drum->SetFamily(2); >>>> DEMSim.SetFamilyFixed(0); >>>> DEMSim.SetFamilyFixed(1); >>>> DEMSim.SetFamilyPrescribedAngVel(2, "0", "0", "0.1"); >>>> >>>> // Define particle template with increased radius >>>> float particle_radius = 0.005f; // Increased radius to reduce >>>> particle count >>>> float particle_volume = (4.0f / 3.0f) * 3.14159f * >>>> std::pow(particle_radius, 3); >>>> float particle_mass = 2.6e3f * particle_volume; // Density * >>>> volume (kg) >>>> float3 particle_moi = make_float3(1.8327927f, 2.1580013f, >>>> 0.77010059f) * (std::pow(particle_radius, 2) / std::pow(0.005f, 2)); >>>> >>>> auto particle_template = DEMSim.LoadClumpType(particle_mass, >>>> particle_moi, "/home/butterco/DEM-Engine/data/clumps/3_clump.csv", >>>> mat_type_mixer); >>>> // No scaling needed as the particle radius matches the clump file >>>> >>>> // Sampler setup with increased spacing >>>> HCPSampler basin_sampler(3.0f * particle_radius); // Increased >>>> spacing >>>> float3 basin_center = make_float3(0.0f, 0.0f, -0.65f); >>>> float basin_radius = 0.07f; >>>> float basin_height = 0.15f; >>>> auto basin_positions = basin_sampler.SampleCylinderZ(basin_center, >>>> basin_radius, basin_height); >>>> >>>> std::cout << "Total number of particles generated: " << >>>> basin_positions.size() << std::endl; >>>> >>>> if (!basin_positions.empty()) { >>>> DEMSim.AddClumps(particle_template, basin_positions); >>>> } else { >>>> std::cerr << "Error: No particles were generated. Adjust >>>> sampling parameters." << std::endl; >>>> return -1; >>>> } >>>> >>>> DEMSim.SetInitTimeStep(5e-6f); >>>> DEMSim.SetGravitationalAcceleration(make_float3(0.0f, 0.0f, >>>> -9.81f)); >>>> DEMSim.SetCDUpdateFreq(40); >>>> DEMSim.SetExpandSafetyAdder(2.0f); >>>> DEMSim.SetErrorOutAvgContacts(50); >>>> DEMSim.SetAdaptiveBinSizeUpperProactivity(0.5f); >>>> DEMSim.SetAdaptiveBinSizeLowerProactivity(0.15f); >>>> >>>> // Set initial bin size to prevent bins from becoming overpopulated >>>> DEMSim.SetInitBinSize(25.0f * particle_radius); >>>> >>>> // Increase maximum spheres per bin if necessary >>>> DEMSim.SetMaxSphereInBin(50000); >>>> >>>> // Initialize the simulation >>>> DEMSim.Initialize(); >>>> >>>> // Create output directory >>>> std::string output_dir = "/home/butterco/DEM-Engine/data/output/"; >>>> fs::create_directories(output_dir); >>>> >>>> auto max_velocity_finder = DEMSim.CreateInspector("clump_max_absv"); >>>> >>>> // Set simulation time and frame rate >>>> float sim_end = 1.0f; >>>> unsigned int fps = 10; >>>> float frame_time = 1.0f / fps; >>>> unsigned int curr_frame = 0; >>>> >>>> // Simulation loop >>>> for (float t = 0.0f; t < sim_end; t += frame_time) { >>>> std::cout << "Simulating frame " << curr_frame << std::endl; >>>> >>>> char filename[200], mesh_filename[200]; >>>> sprintf(filename, "%sDEM_output_%04d.csv", output_dir.c_str(), >>>> curr_frame); >>>> sprintf(mesh_filename, "%sDEM_mesh_output_%04d.vtk", >>>> output_dir.c_str(), curr_frame++); >>>> DEMSim.WriteSphereFile(std::string(filename)); >>>> DEMSim.WriteMeshFile(std::string(mesh_filename)); >>>> >>>> DEMSim.DoDynamics(frame_time); >>>> >>>> float max_velocity = max_velocity_finder->GetValue(); >>>> std::cout << "Max velocity of any particle: " << max_velocity >>>> << std::endl; >>>> } >>>> >>>> DEMSim.ShowTimingStats(); >>>> std::cout << "Simulation complete." << std::endl; >>>> >>>> return 0; >>>> } >>>> " >>>> >>>> >>>> >>>> -- >> You received this message because you are subscribed to a topic in the >> Google Groups "ProjectChrono" group. >> To unsubscribe from this topic, visit >> https://groups.google.com/d/topic/projectchrono/ewhC4AYLF-Y/unsubscribe. >> To unsubscribe from this group and all its topics, send an email to >> [email protected]. >> To view this discussion visit >> https://groups.google.com/d/msgid/projectchrono/c7552937-946a-4488-a32c-27cc429f24c9n%40googlegroups.com >> >> <https://groups.google.com/d/msgid/projectchrono/c7552937-946a-4488-a32c-27cc429f24c9n%40googlegroups.com?utm_medium=email&utm_source=footer> >> . >> > -- You received this message because you are subscribed to the Google Groups "ProjectChrono" group. To unsubscribe from this group and stop receiving emails from it, send an email to [email protected]. To view this discussion visit https://groups.google.com/d/msgid/projectchrono/33ae610f-5c39-4263-8f59-24b464cae99en%40googlegroups.com.
