Hi Ruochun:
The script is attached below.
Thank you,
Weigang
On Thursday, June 20, 2024 at 7:27:35 PM UTC+8 Ruochun Zhang wrote:
> I assume you are doing it via a custom force model. In that case, this
> problem is very similar to the exact demo we provided (which runs). You
> have to tell us what you changed in order for us to give suggestions.
>
> Thank you,
> Ruochun
>
> On Thursday, June 20, 2024 at 6:55:10 PM UTC+8 [email protected] wrote:
>
>> Hi Ruochun,
>>
>> I am currently modelling the impact of a rock block against a plate.
>> However, as shown below, the spheres penetrate through the plate without
>> any contact with the plate. Could you give some suggestions?
>> [image: 微信图片_20240620185350.jpg]
>>
>> Thank you,
>> Weigang
>>
>> On Thursday, June 20, 2024 at 5:02:38 PM UTC+8 Ruochun Zhang wrote:
>>
>>> Hi Weigang,
>>>
>>> Yes. And that's exactly what we did in the fracture demo to begin with,
>>> is it not?
>>>
>>> Thank you,
>>> Ruochun
>>>
>>> On Thursday, June 20, 2024 at 3:07:07 PM UTC+8 [email protected]
>>> wrote:
>>>
>>>> Hi Ruochun,
>>>>
>>>> Do you mean that we can add a if statements in the
>>>> ForceModelWithFracture.cu to differentiate types of contacts?
>>>>
>>>> Thank you,
>>>> Weigang
>>>>
>>>>
>>>>
>>>> On Thursday, June 13, 2024 at 8:30:45 PM UTC+8 Ruochun Zhang wrote:
>>>>
>>>>> Currently, DEM-Engine custom model functionality permits only one
>>>>> script, meaning that in principal, you are writing a unified force model
>>>>> that goes for sphere-sphere, sphere-boundary, *and *sphere-mesh. Just
>>>>> think of the sphere-compressor contact as a special case, which is the
>>>>> contact between a normal sphere and an infinitely-large sphere.
>>>>>
>>>>> About enforcing different models for different objects, you can
>>>>> associate different objects with respective family numbers, then use if
>>>>> statements to differentiate types of contacts. I talked about this in one
>>>>> of the earlier replies in this thread, please go back and have a look.
>>>>>
>>>>> Thank you,
>>>>> Ruochun
>>>>>
>>>>>
>>>>> On Thursday, June 13, 2024 at 1:31:37 PM UTC+8 [email protected]
>>>>> wrote:
>>>>>
>>>>>> Hi Ruochun,
>>>>>>
>>>>>> Thank you for your reply. In fact, I am puzzled that there is not any
>>>>>> sentence in the Fracture-Box demo to declare the contact model between
>>>>>> spheres and the compressing plate. So, how to declare the contact model,
>>>>>> especially for the case in which there has been another contact model,
>>>>>> such
>>>>>> as the ForceModelWithFracture.cu?
>>>>>>
>>>>>> Thank you,
>>>>>> Weigang
>>>>>>
>>>>>> On Thursday, June 13, 2024 at 12:03:52 AM UTC+8 Ruochun Zhang wrote:
>>>>>>
>>>>>>> Hi Weigang,
>>>>>>>
>>>>>>> Note that in general, if you use a custom force model, then that
>>>>>>> model alone is used, there is no fallback so I wouldn't call that
>>>>>>> anything
>>>>>>> has "defaulted" to Hertz–Mindlin. Specific to the Fracture-Box demo,
>>>>>>> it's
>>>>>>> just that the force model (ForceModelWithFracture.cu) is written in a
>>>>>>> way
>>>>>>> that when the bond between two particles breaks, the "else" branch of
>>>>>>> that
>>>>>>> big if statement is executed, and something similar to Hertz–Mindlin is
>>>>>>> enforced.
>>>>>>>
>>>>>>> I said it's similar since it's not a pure spring–damper model now,
>>>>>>> it's a spring–damper model with the resting location (zero-force
>>>>>>> location)
>>>>>>> being the penetration depth the moment when the bond broke. This
>>>>>>> implementation avoids a potential problem with a pure spring–damper
>>>>>>> model:
>>>>>>> Without the adjusted resting location, a newly broken-away pair of
>>>>>>> particles would see the bond suddenly gone but there were still
>>>>>>> left-over
>>>>>>> penetration, therefore immediately pushing each other away with high
>>>>>>> velocity, resulting in non-physical results. Shout-out to Bona for
>>>>>>> implementing this.
>>>>>>>
>>>>>>> Thank you!
>>>>>>> Ruochun
>>>>>>>
>>>>>>>
>>>>>>> On Wednesday, June 12, 2024 at 8:33:55 PM UTC+8 [email protected]
>>>>>>> wrote:
>>>>>>>
>>>>>>>> Hi Ruochun,
>>>>>>>>
>>>>>>>> Thank you for your great work. The demo now is runs well. I have a
>>>>>>>> question about the Fracture-Box demo. After the bond between two
>>>>>>>> particles
>>>>>>>> breaks, is it right that the contact model between them change
>>>>>>>> defaultly
>>>>>>>> from ForceModelWithFracture to Hertz–Mindlin model?
>>>>>>>>
>>>>>>>> Thank you,
>>>>>>>> Weigang
>>>>>>>>
>>>>>>>> On Friday, May 31, 2024 at 8:05:24 PM UTC+8 Ruochun Zhang wrote:
>>>>>>>>
>>>>>>>>> Hi Weigang,
>>>>>>>>>
>>>>>>>>> The problem appears to be fixed now. Please pull the newest main
>>>>>>>>> branch and try the fracture demo again. I suggest you do a clean
>>>>>>>>> rebuild
>>>>>>>>> from an empty directory, because only a force model file is changed
>>>>>>>>> and if
>>>>>>>>> you just "make", cmake probably thinks nothing needs to be done. Or
>>>>>>>>> you can
>>>>>>>>> simply manually copy the modified cu file over to the appropriate
>>>>>>>>> location
>>>>>>>>> in the build folder.
>>>>>>>>>
>>>>>>>>> The problem was that we left an uninitialized variable in the
>>>>>>>>> force model file. On data center cards, the compiler zeros it out
>>>>>>>>> automatically, but on consumer cards it does not. We therefore missed
>>>>>>>>> it
>>>>>>>>> initially. Should be all good now.
>>>>>>>>>
>>>>>>>>> Thank you,
>>>>>>>>> Ruochun
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> On Tue, May 28, 2024, 11:00 PM Ruochun Zhang <[email protected]>
>>>>>>>>> wrote:
>>>>>>>>>
>>>>>>>>>> Hi Weigang,
>>>>>>>>>>
>>>>>>>>>> It's not like recommending a card; rather, if you happen to have
>>>>>>>>>> a data center card available (A100, A5000 etc.), you can circumvent
>>>>>>>>>> this
>>>>>>>>>> problem for now by using them. Otherwise, let's hope we find a
>>>>>>>>>> solution
>>>>>>>>>> very soon, or at least find a workaround. We'll keep you posted.
>>>>>>>>>>
>>>>>>>>>> Thank you,
>>>>>>>>>> Ruochun
>>>>>>>>>>
>>>>>>>>>> On Monday, May 27, 2024 at 7:27:10 PM UTC+8 [email protected]
>>>>>>>>>> wrote:
>>>>>>>>>>
>>>>>>>>>>> Hi Ruochun,
>>>>>>>>>>>
>>>>>>>>>>> My card is indeed a gaming card NVIDIA GeForce RTX 2070 Super.
>>>>>>>>>>> Could you recommend a card?
>>>>>>>>>>>
>>>>>>>>>>> Thanks,
>>>>>>>>>>> Weigang
>>>>>>>>>>>
>>>>>>>>>>> On Saturday, May 25, 2024 at 10:44:17 PM UTC+8 Ruochun Zhang
>>>>>>>>>>> wrote:
>>>>>>>>>>>
>>>>>>>>>>>> Hi Weigang,
>>>>>>>>>>>>
>>>>>>>>>>>> I would like to answer Question 2 first. You can use variables
>>>>>>>>>>>> *AOwnerFamily
>>>>>>>>>>>> *and *BOwnerFamily *to refer to clump (or owner) A's and clump
>>>>>>>>>>>> B's family numbers directly in the force model. So you can write
>>>>>>>>>>>> *if
>>>>>>>>>>>> *statements that treat each case accordingly. This is
>>>>>>>>>>>> mentioned in the DEME paper: You can check out Table 2.
>>>>>>>>>>>>
>>>>>>>>>>>> For Question 1, the situation is more complicated. We are
>>>>>>>>>>>> actually able to reproduce the problem you encountered (simulation
>>>>>>>>>>>> freezing), but on consumer/gaming cards only. It runs perfectly on
>>>>>>>>>>>> data
>>>>>>>>>>>> center cards and that's what we used for generating this
>>>>>>>>>>>> experiment,
>>>>>>>>>>>> therefore it did not appear a problem for us. It's likely due to
>>>>>>>>>>>> different
>>>>>>>>>>>> GPU architectures running the same code differently, we'll try to
>>>>>>>>>>>> find a
>>>>>>>>>>>> solution and let you know in the following days, so stay tuned. By
>>>>>>>>>>>> the way,
>>>>>>>>>>>> you are using a gaming card to run the simulations, correct?
>>>>>>>>>>>>
>>>>>>>>>>>> Thank you,
>>>>>>>>>>>> Ruochun
>>>>>>>>>>>>
>>>>>>>>>>>> On Friday, May 24, 2024 at 8:47:00 AM UTC+8 [email protected]
>>>>>>>>>>>> wrote:
>>>>>>>>>>>>
>>>>>>>>>>>>> Thank you for your reply! I have two questions about the demo
>>>>>>>>>>>>> "DEMdemo-Fracture-Box".
>>>>>>>>>>>>>
>>>>>>>>>>>>> (1) After i have run this demo for one day, there is still
>>>>>>>>>>>>> nothing outputed. Is it fact that the DEME will run a long time
>>>>>>>>>>>>> once a
>>>>>>>>>>>>> fracture model is used?
>>>>>>>>>>>>>
>>>>>>>>>>>>> (2) How to tell DEME the contact force models between
>>>>>>>>>>>>> different familys? I did not see any snippet which undertakes
>>>>>>>>>>>>> this work,
>>>>>>>>>>>>> although the demo has at least two kinds of contact force model.
>>>>>>>>>>>>> The
>>>>>>>>>>>>> question i really want to ask is that, if we have more than three
>>>>>>>>>>>>> familys
>>>>>>>>>>>>> in our model and there are more than two familys consist of
>>>>>>>>>>>>> spherical
>>>>>>>>>>>>> particles, how to tell DEME the contact force model between
>>>>>>>>>>>>> different
>>>>>>>>>>>>> familys, and tell DEME the contact force model between the
>>>>>>>>>>>>> particles of a
>>>>>>>>>>>>> family in the same time.
>>>>>>>>>>>>>
>>>>>>>>>>>>> Thank you,
>>>>>>>>>>>>> Weigang
>>>>>>>>>>>>>
>>>>>>>>>>>>> On Thursday, May 23, 2024 at 10:10:14 PM UTC+8 Ruochun Zhang
>>>>>>>>>>>>> wrote:
>>>>>>>>>>>>>
>>>>>>>>>>>>>> No. That's because the *sphere components* in a clump are
>>>>>>>>>>>>>> just shape representations, and they have no physics, only
>>>>>>>>>>>>>> material
>>>>>>>>>>>>>> properties. If they did have interaction with each other then
>>>>>>>>>>>>>> every clump
>>>>>>>>>>>>>> would've immediately exploded upon simulation starting.
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Thank you,
>>>>>>>>>>>>>> Ruochun
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> On Thursday, May 23, 2024 at 12:16:32 PM UTC+8
>>>>>>>>>>>>>> [email protected] wrote:
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Hi Ruochun,
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Thank you for your reply! I have a question about the clump.
>>>>>>>>>>>>>>> Do the particles in a clump interact with each other via any
>>>>>>>>>>>>>>> force model?
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Thank you,
>>>>>>>>>>>>>>> Weigang
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> On Monday, May 20, 2024 at 10:26:03 PM UTC+8 Ruochun Zhang
>>>>>>>>>>>>>>> wrote:
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Hi Weigang,
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Yes. Again, DEME only cares about the position and radius
>>>>>>>>>>>>>>>> of the initial spheres, so you can supply it however you want.
>>>>>>>>>>>>>>>> If you have
>>>>>>>>>>>>>>>> a file that records this information using your own format,
>>>>>>>>>>>>>>>> then the
>>>>>>>>>>>>>>>> easiest thing is to have your own C++ or Python script that
>>>>>>>>>>>>>>>> reads it into
>>>>>>>>>>>>>>>> arrays and then feeds them to the solver. If your file
>>>>>>>>>>>>>>>> shares the format with DEME's standard clump output file, or
>>>>>>>>>>>>>>>> your file is
>>>>>>>>>>>>>>>> simply the output of *WriteClumpFile*, then you can
>>>>>>>>>>>>>>>> conveniently load clump positions and orientations by
>>>>>>>>>>>>>>>> *ReadClumpXyzFromCsv
>>>>>>>>>>>>>>>> *and *ReadClumpQuatFromCsv *(examples are in GRCPrep_Part1
>>>>>>>>>>>>>>>> and 2). Note that currently, sphere radius cannot be read from
>>>>>>>>>>>>>>>> standard
>>>>>>>>>>>>>>>> clump output files since it's part of the clump template
>>>>>>>>>>>>>>>> information thus
>>>>>>>>>>>>>>>> not in the clump output.
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Thank you,
>>>>>>>>>>>>>>>> Ruochun
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> On Monday, May 20, 2024 at 12:56:56 PM UTC+8
>>>>>>>>>>>>>>>> [email protected] wrote:
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>> Hi Everyone,
>>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>> I would like to generate a block which consists of a
>>>>>>>>>>>>>>>>> assembly spheres. Is it possible to generate the block by
>>>>>>>>>>>>>>>>> other code, and
>>>>>>>>>>>>>>>>> then load it into DEM-Engine via a file containg the
>>>>>>>>>>>>>>>>> informations (such as
>>>>>>>>>>>>>>>>> position and radius) of the spheres?
>>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>>> Thank you,
>>>>>>>>>>>>>>>>> Weigang
>>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> --
>>>>>>>>>> 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 on the web visit
>>>>>>>>>> https://groups.google.com/d/msgid/projectchrono/1910af21-d1e3-43dd-a779-0c229e6150a6n%40googlegroups.com
>>>>>>>>>>
>>>>>>>>>> <https://groups.google.com/d/msgid/projectchrono/1910af21-d1e3-43dd-a779-0c229e6150a6n%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 on the web visit
https://groups.google.com/d/msgid/projectchrono/f0d53f93-3d0d-439e-b1c8-4c02030144a5n%40googlegroups.com.
// Copyright (c) 2021, SBEL GPU Development Team
// Copyright (c) 2021, University of Wisconsin - Madison
//
// SPDX-License-Identifier: BSD-3-Clause
// =============================================================================
// Fracture
// =============================================================================
#include <DEM/API.h>
#include <DEM/HostSideHelpers.hpp>
#include <DEM/utils/Samplers.hpp>
#include <chrono>
#include <cmath>
#include <cstdio>
#include <filesystem>
#include <iostream>
#include <fstream>
using namespace deme;
const double math_PI = 3.1415927;
int main() {
DEMSolver DEMSim;
DEMSim.SetVerbosity(INFO);
DEMSim.SetOutputFormat(OUTPUT_FORMAT::CSV);
DEMSim.SetOutputContent(OUTPUT_CONTENT::VEL);
DEMSim.SetMeshOutputFormat(MESH_FORMAT::VTK);
DEMSim.SetContactOutputContent(OWNER | FORCE | CNT_WILDCARD);
// This demo could lead to large numbers of per-sphere contacts, so to be safe...
DEMSim.SetErrorOutAvgContacts(500);
// E, nu, CoR, mu, Crr...
auto mat_type_container =
DEMSim.LoadMaterial({{"E", 1e8}, {"nu", 0.3}, {"CoR", 0.7}, {"mu", 0.80}, {"Crr", 0.10}});
auto mat_type_particle = DEMSim.LoadMaterial({{"E", 1e8}, {"nu", 0.20}, {"CoR", 0.5}, {"mu", 0.5}, {"Crr", 0.05}});
// If you don't have this line, then values will take average between 2 materials, when they are in contact
DEMSim.SetMaterialPropertyPair("CoR", mat_type_container, mat_type_particle, 0.2);
DEMSim.SetMaterialPropertyPair("Crr", mat_type_container, mat_type_particle, 0.5);
DEMSim.SetMaterialPropertyPair("mu", mat_type_container, mat_type_particle, 0.5);
// We can specify the force model using a file.
auto my_force_model = DEMSim.ReadContactForceModel("ForceModelWithFractureModel.cu");
// Those following lines are needed. We must let the solver know that those var names are history variable etc.
my_force_model->SetMustHaveMatProp({"E", "nu", "CoR", "mu", "Crr"});
my_force_model->SetMustPairwiseMatProp({"CoR", "mu", "Crr"});
// Pay attention to the extra per-contact wildcard `unbroken' here.
my_force_model->SetPerContactWildcards(
{"delta_time", "delta_tan_x", "delta_tan_y", "delta_tan_z", "unbroken", "initialLength"});
float world_size = 1.5;
float container_diameter = 0.05;
float terrain_density = 2.80e3;
float sphere_rad = 0.0020;
float step_size = 2e-7;
float fact_radius = 0.9;
DEMSim.InstructBoxDomainDimension(world_size, world_size, world_size);
// No need to add simulation `world' boundaries, b/c we'll add a cylinderical container manually
DEMSim.InstructBoxDomainBoundingBC("all", mat_type_container);
// DEMSim.SetInitBinSize(sphere_rad * 5);
// Now add a cylinderical boundary along with a bottom plane
double bottom = -0;
double top = 0.10;
// auto walls = DEMSim.AddExternalObject();
// std::string shake_pattern_xz = " 0.00 * sin( 500 * 2 * deme::PI * t)";
// walls->AddPlane(make_float3(0, 0, bottom), make_float3(0, 0, 1), mat_type_container);
// walls->AddPlane(make_float3(0, 0, world_size / 2. - world_size / 20.), make_float3(0, 0, -1),
// mat_type_container); walls->SetFamily(2);
auto walls = DEMSim.AddWavefrontMeshObject("../data/mesh/funnel_left.obj", mat_type_container);
float3 move = make_float3(0.05, 0.00, 0 - 2 * sphere_rad); // z
float4 rot = make_float4(0.7071, 0, 0, 0.7071);
walls->Scale(make_float3(0.8, 0.1, 2.0));
walls->Move(move, rot);
walls->SetFamily(2);
DEMSim.SetFamilyFixed(2);
auto cylinder = DEMSim.AddExternalObject();
cylinder->AddCylinder(make_float3(0), make_float3(0, 0, 1), 1.8 * container_diameter / 2., mat_type_container, 0);
cylinder->SetFamily(10);
DEMSim.SetFamilyFixed(10);
// DEMSim.SetFamilyPrescribedLinVel(10, shake_pattern_xz, shake_pattern_xz, shake_pattern_xz);
//auto plate = DEMSim.AddWavefrontMeshObject("../data/mesh/funnel_left.obj", mat_type_container);
//move = make_float3(0.05, 0.00, top + 2.0 * sphere_rad); // z
//rot = make_float4(0.7071, 0, 0, 0.7071);
//plate->Scale(make_float3(0.8, 0.1, 2.0));
//plate->Move(move, rot);
//plate->SetFamily(20);
//DEMSim.SetFamilyFixed(20);
// Track it
auto plate_tracker = DEMSim.Track(walls);
//DEMSim.SetFamilyPrescribedLinVel(22, "0", "0", to_string_with_precision(-0.0050));
// auto total_mass_finder = DEMSim.CreateInspector("clump_mass");
// auto max_v_finder = DEMSim.CreateInspector("clump_max_absv");
// Define the terrain particle templates
// Calculate its mass and MOI
float sphere_vol = 4. / 3. * math_PI * sphere_rad * sphere_rad * sphere_rad;
float mass = terrain_density * sphere_vol;
// Then load it to system
std::shared_ptr<DEMClumpTemplate> my_template = DEMSim.LoadSphereType(mass, sphere_rad, mat_type_particle);
// Sampler to sample
HCPSampler sampler(2.0 * sphere_rad);
float fill_height = top;
float3 fill_center = make_float3(0, 0, bottom + fill_height / 2);
const float fill_radius = container_diameter / 2. - sphere_rad * 0.;
float3 fill = make_float3(fill_radius, fill_radius, bottom + fill_height / 2);
auto input_xyz = sampler.SampleBox(fill_center, fill);
// auto input_xyz = sampler(fill_center, fill_radius, fill_height / 2 - sphere_rad * 0.);
auto particles = DEMSim.AddClumps(my_template, input_xyz);
particles->SetFamily(1);
std::cout << "Total num of particles: " << particles->GetNumClumps() << std::endl;
DEMSim.SetFamilyPrescribedLinVel(1, "0", "0", to_string_with_precision(-1.0));
std::filesystem::path out_dir = std::filesystem::current_path();
std::string nameOutFolder = "R" + std::to_string(sphere_rad) + "_Int" + std::to_string(fact_radius) + "";
out_dir += "/BlockImpact_" + nameOutFolder;
remove_all(out_dir);
create_directory(out_dir);
// Some inspectors
auto max_z_finder = DEMSim.CreateInspector("clump_max_z");
auto min_z_finder = DEMSim.CreateInspector("clump_min_z");
DEMSim.SetFamilyExtraMargin(1, fact_radius * sphere_rad);
DEMSim.SetInitTimeStep(step_size);
DEMSim.SetGravitationalAcceleration(make_float3(0, 0.00, 1 * -9.81));
DEMSim.Initialize();
// DEMSim.DisableContactBetweenFamilies(20, 1);
std::cout << "Initial number of contacts: " << DEMSim.GetNumContacts() << std::endl;
float sim_end = 5;
unsigned int fps = 1000;
unsigned int datafps = 1000;
unsigned int modfpsGeo = datafps / fps;
float frame_time = 1.0 / datafps;
std::cout << "Output at " << fps << " FPS" << std::endl;
unsigned int out_steps = (unsigned int)(1.0 / (datafps * step_size));
unsigned int frame_count = 0;
unsigned int step_count = 0;
bool status_1 = true;
bool status_2 = true;
// DEMSim.DisableContactBetweenFamilies(10, 1);
double L0;
double stress;
std::string nameOutFile = "data_R" + std::to_string(sphere_rad) + "_Int" + std::to_string(fact_radius) + ".csv";
std::ofstream csvFile(nameOutFile);
DEMSim.SetFamilyContactWildcardValueAll(1, "initialLength", 0.0);
// DEMSim.SetFamilyContactWildcardValueAll(1, "damage", 0.0);
DEMSim.SetFamilyContactWildcardValueAll(1, "unbroken", 0.0);
// Simulation loop
for (float t = 0; t < sim_end; t += frame_time) {
// DEMSim.ShowThreadCollaborationStats();
std::cout << "Initial number of contacts: " << DEMSim.GetNumContacts() << std::endl;
if (frame_count % 1 == 0) {
float3 forces = plate_tracker->ContactAcc();
double L = max_z_finder->GetValue() - min_z_finder->GetValue() + 2 * sphere_rad;
std::cout << "Time: " << t << std::endl;
std::cout << "Force [N]: " << stress << std::endl;
if (frame_count % modfpsGeo == 0) {
char filename[200];
char meshname[200];
char cnt_filename[200];
std::cout << "Outputting frame: " << frame_count / modfpsGeo << std::endl;
sprintf(filename, "%s/DEMdemo_output_%04d.csv", out_dir.c_str(), frame_count / modfpsGeo);
sprintf(meshname, "%s/DEMdemo_mesh_%04d.vtk", out_dir.c_str(), frame_count / modfpsGeo);
sprintf(cnt_filename, "%s/DEMdemo_contact_%04d.csv", out_dir.c_str(), frame_count / modfpsGeo);
DEMSim.WriteSphereFile(std::string(filename));
DEMSim.WriteMeshFile(std::string(meshname));
DEMSim.WriteContactFile(std::string(cnt_filename));
}
}
DEMSim.DoDynamics(frame_time);
frame_count++;
}
csvFile.close();
DEMSim.ShowTimingStats();
std::cout << "Fracture demo exiting..." << std::endl;
return 0;
}