Sure. I modified the SingleWheelTest demo to run the wheel under different slip ratios. I set a for loop in main function. But only the first slip ratio ran.
在2023年4月3日星期一 UTC+8 21:08:18<Dan Negrut> 写道： > Perhaps we can help if you drop here the file that runs the loop to carry > out the analysis. We need more info to try to help you. > > > > Dan > > --------------------------------------------- > > Bernard A. and Frances M. Weideman Professor > > NVIDIA CUDA Fellow > > Department of Mechanical Engineering > > Department of Computer Science > > University of Wisconsin - Madison > > 4150ME, 1513 University Avenue > > Madison, WI 53706-1572 > > 608 772 0914 <(608)%20772-0914> > > http://sbel.wisc.edu/ > > http://projectchrono.org/ > > --------------------------------------------- > > > > *From:* [email protected] <[email protected]> *On > Behalf Of *Chiny Lan > *Sent:* Monday, April 3, 2023 4:45 AM > *To:* ProjectChrono <[email protected]> > *Subject:* [chrono] Simulation with a set of parameters in one program > > > > Hi, > > > > I want to simulate sand with different parameters using FSI. I used a for > loop in main function to implement the series of simulation. But when the > first simulation done, the program also closed. > > > > How can I do seris simulation with a set of parameters in one program? > > > > Thanks and regards. > > -- > 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/46b75e69-cef3-4c4f-9b00-26603e5bfff7n%40googlegroups.com > > <https://groups.google.com/d/msgid/projectchrono/46b75e69-cef3-4c4f-9b00-26603e5bfff7n%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/0503b868-9cbd-48f5-895e-a8da1ab431e4n%40googlegroups.com.
// ============================================================================= // PROJECT CHRONO - http://projectchrono.org // // Copyright (c) 2014 projectchrono.org // All rights reserved. // // Use of this source code is governed by a BSD-style license that can be found // in the LICENSE file at the top level of the distribution and at // http://projectchrono.org/license-chrono.txt. // // ============================================================================= // Author: Wei Hu, Radu Serban // ============================================================================= // 2023.04.01: Qingning Lan // Modified the code to run single wheel test with for loop in different slip // Slip ratio = -1.0, -0.6, -0.3, 0.0, 0.3, 0.6, 0.8 // ============================================================================= #include <cassert> #include <cstdlib> #include <iostream> #include <vector> #include <string> #include <sstream> #include <cmath> #include "chrono/ChConfig.h" #include "chrono/physics/ChSystemSMC.h" #include "chrono/physics/ChBody.h" #include "chrono/physics/ChInertiaUtils.h" #include "chrono/physics/ChLinkMotorRotationAngle.h" #include "chrono/utils/ChUtilsCreators.h" #include "chrono/utils/ChUtilsGenerators.h" #include "chrono/utils/ChUtilsGeometry.h" #include "chrono/assets/ChBoxShape.h" #include "chrono/utils/ChUtilsInputOutput.h" #include "chrono/assets/ChTriangleMeshShape.h" #include "chrono/core/ChTimer.h" #include "chrono_fsi/ChSystemFsi.h" #include "chrono_fsi/ChVisualizationFsi.h" #include "chrono_thirdparty/filesystem/path.h" // Chrono namespaces using namespace chrono; using namespace chrono::fsi; using namespace chrono::geometry; using namespace std; // Set Slip ratio vector std::vector<double> slip_ratios = {-1.0, -0.6, -0.3, 0.0, 0.3, 0.6, 0.8}; std::vector<string> slip_ratios_str = {"-1.0", "-0.6", "-0.3", "0.0", "0.3", "0.6", "0.8"}; // Physical properties of terrain particles double iniSpacing = 0.004; double kernelLength = 0.004; double density = 1700.0; // Dimension of the terrain container double smalldis = 1.0e-9; double bxDim = 1.0 + smalldis; double byDim = 0.4 + smalldis; double bzDim = 0.18 + smalldis; // Size of the wheel double wheel_radius = 0.135; double wheel_width = 0.165; double velocity = 0.04; // Gravity is 80N // double total_mass = 8.1549439348; double total_mass = 13.25; // Simulation time and stepsize // Simulate 20sec for each slip ratio, displacement = 0.6m double total_time = 16.0; double dT = 1.0e-4; // linear actuator and angular actuator auto actuator = chrono_types::make_shared<ChLinkLinActuator>(); auto motor = chrono_types::make_shared<ChLinkMotorRotationAngle>(); // Save data as csv files to see the results off-line using Paraview bool output = true; int out_fps = 10; // Enable/disable run-time visualization (if Chrono::OpenGL is available) bool render = true; float render_fps = 100; // Verbose terminal output bool verbose_fsi = true; bool verbose = true; //------------------------------------------------------------------ // Function to save wheel to Paraview VTK files //------------------------------------------------------------------ void WriteCylinderVTK(const std::string& filename, double radius, double length, const ChFrame<>& frame, unsigned int res) { std::ofstream outf; outf.open(filename, std::ios::app); outf << "# vtk DataFile Version 2.0\nUnstructured Grid Example\nASCII\n\n" << std::endl; outf << "DATASET UNSTRUCTURED_GRID\nPOINTS " << 2 * res << " float\n"; for (int i = 0; i < res; i++) { auto w = frame.TransformPointLocalToParent( ChVector<>(radius * cos(2 * i * 3.1415 / res), -1 * length / 2, radius * sin(2 * i * 3.1415 / res))); outf << w.x() << " " << w.y() << " " << w.z() << "\n"; } for (int i = 0; i < res; i++) { auto w = frame.TransformPointLocalToParent( ChVector<>(radius * cos(2 * i * 3.1415 / res), +1 * length / 2, radius * sin(2 * i * 3.1415 / res))); outf << w.x() << " " << w.y() << " " << w.z() << "\n"; } outf << "CELLS " << res + res << "\t" << 5 * (res + res) << "\n"; for (int i = 0; i < res - 1; i++) { outf << "4 " << i << " " << i + 1 << " " << i + res + 1 << " " << i + res << "\n"; } outf << "4 " << res - 1 << " " << 0 << " " << res << " " << 2 * res - 1 << "\n"; for (int i = 0; i < res / 4; i++) { outf << "4 " << i << " " << i + 1 << " " << +res / 2 - i - 1 << " " << +res / 2 - i << "\n"; } for (int i = 0; i < res / 4; i++) { outf << "4 " << i + res << " " << i + 1 + res << " " << +res / 2 - i - 1 + res << " " << +res / 2 - i + res << "\n"; } outf << "4 " << +res / 2 << " " << 1 + res / 2 << " " << +res - 1 << " " << 0 << "\n"; for (int i = 1; i < res / 4; i++) { outf << "4 " << i + res / 2 << " " << i + 1 + res / 2 << " " << +res / 2 - i - 1 + res / 2 << " " << +res / 2 - i + res / 2 << "\n"; } outf << "4 " << 3 * res / 2 << " " << 1 + 3 * res / 2 << " " << +2 * res - 1 << " " << +res << "\n"; for (int i = 1; i < res / 4; i++) { outf << "4 " << i + 3 * res / 2 << " " << i + 1 + 3 * res / 2 << " " << +2 * res - i - 1 << " " << +2 * res - i << "\n"; } outf << "\nCELL_TYPES " << res + res << "\n"; for (int iele = 0; iele < (res + res); iele++) { outf << "9\n"; } } //------------------------------------------------------------------ // Create the objects of the MBD system. Rigid bodies, and if FSI, // their BCE representation are created and added to the systems //------------------------------------------------------------------ void CreateSolidPhase(ChSystemSMC& sysMBS, ChSystemFsi& sysFSI, double wheel_slip) { // Wheel angular velocity double wheel_AngVel; if (wheel_slip >= 0) { wheel_AngVel = velocity / wheel_radius / (1 - wheel_slip); } else { wheel_AngVel = velocity / wheel_radius * (1 + wheel_slip); } // Common contact material auto cmaterial = chrono_types::make_shared<ChMaterialSurfaceSMC>(); cmaterial->SetYoungModulus(1e8); cmaterial->SetFriction(0.9f); cmaterial->SetRestitution(0.4f); cmaterial->SetAdhesion(0); // Create a container -- always FIRST body in the system auto ground = chrono_types::make_shared<ChBody>(); ground->SetIdentifier(-1); ground->SetBodyFixed(true); sysMBS.AddBody(ground); ground->GetCollisionModel()->ClearModel(); chrono::utils::AddBoxContainer(ground, cmaterial, ChFrame<>(), ChVector<>(bxDim, byDim, bzDim), 0.1, ChVector<int>(0, 0, -1), false); ground->GetCollisionModel()->BuildModel(); ground->SetCollide(true); // Add BCE particles attached on the walls into FSI system sysFSI.AddContainerBCE(ground, ChFrame<>(), ChVector<>(bxDim, byDim, 2 * bzDim), ChVector<int>(2, 0, -1)); // Create the wheel -- always SECOND body in the system auto wheel = chrono_types::make_shared<ChBody>(); // Set the general wheel properties ChVector<> gyration = chrono::utils::CalcCylinderGyration(wheel_radius, wheel_width / 2).diagonal(); ChVector<> wheel_IniPos(-bxDim / 2 + 1.2 * wheel_radius, 0.0, wheel_radius + bzDim + 1.0 * iniSpacing); ChVector<> wheel_IniVel(0.0, 0.0, 0.0); wheel->SetPos(wheel_IniPos); wheel->SetPos_dt(wheel_IniVel); wheel->SetWvel_loc(ChVector<>(0.0, 0.0, 0.0)); wheel->SetMass(total_mass * 1.0 / 2.0); wheel->SetInertiaXX(total_mass * gyration); // Set the collision wheel->SetCollide(false); wheel->SetBodyFixed(false); wheel->GetCollisionModel()->ClearModel(); // wheel->GetCollisionModel()->SetSafeMargin(iniSpacing); // chrono::utils::AddCylinderGeometry(wheel.get(), cmaterial, wheel_radius, wheel_width / 2, VNULL, QUNIT); wheel->GetCollisionModel()->AddCylinder(cmaterial, wheel_radius, wheel_radius, wheel_width / 2, VNULL, QUNIT); // wheel->GetCollisionModel()->SetEnvelope(0.001); wheel->GetCollisionModel()->BuildModel(); // Add the wheel to the MBD system sysMBS.AddBody(wheel); // Add the wheel to the FSI system sysFSI.AddFsiBody(wheel); // Add BCE particles attached on the wheel into FSI system sysFSI.AddCylinderBCE(wheel, ChFrame<>(VNULL, Q_from_AngX(CH_C_PI_2)), wheel_radius, wheel_width, true); // Create the chassis -- always THIRD body in the system auto chassis = chrono_types::make_shared<ChBody>(); chassis->SetMass(total_mass * 1.0 / 2.0); chassis->SetPos(wheel->GetPos()); chassis->SetCollide(false); chassis->SetBodyFixed(false); // Add geometry of the chassis. chassis->GetCollisionModel()->ClearModel(); chrono::utils::AddBoxGeometry(chassis.get(), cmaterial, ChVector<>(0.1, 0.01, 0.01), ChVector<>(0, 0, 0)); chassis->GetCollisionModel()->BuildModel(); sysMBS.AddBody(chassis); // Create the axle -- always FOURTH body in the system auto axle = chrono_types::make_shared<ChBody>(); axle->SetMass(total_mass * 1.0 / 2.0); axle->SetPos(wheel->GetPos()); axle->SetCollide(false); axle->SetBodyFixed(false); // Add geometry of the axle. axle->GetCollisionModel()->ClearModel(); chrono::utils::AddSphereGeometry(axle.get(), cmaterial, 0.01, ChVector<>(0, 0, 0)); axle->GetCollisionModel()->BuildModel(); sysMBS.AddBody(axle); // Connect the chassis to the containing bin (ground) through a translational joint and create a linear actuator. auto prismatic1 = chrono_types::make_shared<ChLinkLockPrismatic>(); prismatic1->Initialize(ground, chassis, ChCoordsys<>(chassis->GetPos(), Q_from_AngY(CH_C_PI_2))); prismatic1->SetName("prismatic_chassis_ground"); sysMBS.AddLink(prismatic1); auto actuator_fun = chrono_types::make_shared<ChFunction_Ramp>(0.0, velocity); actuator->Initialize(ground, chassis, false, ChCoordsys<>(chassis->GetPos(), QUNIT), ChCoordsys<>(chassis->GetPos() + ChVector<>(1, 0, 0), QUNIT)); actuator->SetName("actuator"); actuator->SetDistanceOffset(1); actuator->SetActuatorFunction(actuator_fun); sysMBS.AddLink(actuator); // Connect the axle to the chassis through a vertical translational joint. auto prismatic2 = chrono_types::make_shared<ChLinkLockPrismatic>(); prismatic2->Initialize(chassis, axle, ChCoordsys<>(chassis->GetPos(), QUNIT)); prismatic2->SetName("prismatic_axle_chassis"); sysMBS.AddLink(prismatic2); // Connect the wheel to the axle through a engine joint. motor->SetName("engine_wheel_axle"); motor->Initialize(wheel, axle, ChFrame<>(wheel->GetPos(), chrono::Q_from_AngAxis(-CH_C_PI / 2.0, ChVector<>(1, 0, 0)))); motor->SetAngleFunction(chrono_types::make_shared<ChFunction_Ramp>(0, wheel_AngVel)); sysMBS.AddLink(motor); } // ============================================================================= int main(int argc, char* argv[]) { // Output directories and settings, set the Release directory as the working directory const std::string out_dir_base = "output/"; if (!filesystem::create_directory(filesystem::path(out_dir_base))) std::cout << "Error creating directory " << out_dir_base << std::endl; // A loop to run the simulation for different slip ratios for (int i = 0; i < slip_ratios.size(); i++) { // Create the output directory std::string out_dir = out_dir_base + "S" + slip_ratios_str[i]; // Create oputput directories if (!filesystem::create_directory(filesystem::path(out_dir))) { std::cerr << "Error creating directory " << out_dir << std::endl; return 1; } if (!filesystem::create_directory(filesystem::path(out_dir + "/particles"))) { std::cerr << "Error creating directory " << out_dir + "/particles" << std::endl; return 1; } if (!filesystem::create_directory(filesystem::path(out_dir + "/fsi"))) { std::cerr << "Error creating directory " << out_dir + "/fsi" << std::endl; return 1; } if (!filesystem::create_directory(filesystem::path(out_dir + "/vtk"))) { std::cerr << "Error creating directory " << out_dir + "/vtk" << std::endl; return 1; } // Create the MBS and FSI systems ChSystemSMC sysMBS; ChSystemFsi sysFSI(&sysMBS); ChVector<> gravity = ChVector<>(0, 0, -9.81); sysMBS.Set_G_acc(gravity); sysFSI.Set_G_acc(gravity); sysFSI.SetVerbose(verbose_fsi); // Use the default input file or you may enter your input parameters as a command line argument std::string inputJson = "demo_FSI_SingleWheelTest.json"; sysFSI.ReadParametersFromFile(inputJson); // Set the terrain container size sysFSI.SetContainerDim(ChVector<>(bxDim, byDim, bzDim)); // Set SPH discretization type, consistent or inconsistent sysFSI.SetDiscreType(false, false); // Set wall boundary condition sysFSI.SetWallBC(BceVersion::ADAMI); // Set rigid body boundary condition sysFSI.SetRigidBodyBC(BceVersion::ADAMI); // Set cohsion of the granular material sysFSI.SetCohesionForce(0); // Setup the SPH method sysFSI.SetSPHMethod(FluidDynamics::WCSPH); // Set up the periodic boundary condition (if not, set relative larger values) ChVector<> cMin(-bxDim / 2 * 10, -byDim / 2 - 0.5 * iniSpacing, -bzDim * 3); ChVector<> cMax(bxDim / 2 * 10, byDim / 2 + 0.5 * iniSpacing, bzDim * 10); sysFSI.SetBoundaries(cMin, cMax); // Initialize the SPH particles chrono::utils::GridSampler<> sampler(iniSpacing); ChVector<> boxCenter(0.0, 0.0, bzDim / 2); ChVector<> boxHalfDim(bxDim / 2, byDim / 2, bzDim / 2); sysFSI.AddBoxSPH(boxCenter, boxHalfDim); // Create Solid region and attach BCE SPH particles CreateSolidPhase(sysMBS, sysFSI, slip_ratios[i]); // Set simulation data output length sysFSI.SetOutputLength(0); // Construction of the FSI system must be finalized before running sysFSI.Initialize(); auto wheel = sysMBS.Get_bodylist()[1]; ChVector<> force = actuator->Get_react_force(); ChVector<> torque = motor->Get_react_torque(); ChVector<> w_pos = wheel->GetPos(); ChVector<> w_vel = wheel->GetPos_dt(); ChVector<> angvel = wheel->GetWvel_loc(); // Write the information into a txt file std::ofstream myFile; std::ofstream myDBP_Torque; if (output) { myFile.open(out_dir + "/results.txt", std::ios::trunc); myDBP_Torque.open(out_dir + "/DBP_Torque.txt", std::ios::trunc); } // Create a run-tme visualizer ChVisualizationFsi fsi_vis(&sysFSI); if (render) { fsi_vis.SetTitle("Chrono::FSI single wheel demo"); fsi_vis.SetCameraPosition(ChVector<>(0.8 * bxDim, -3 * byDim, 3 * bzDim), ChVector<>(0, 1, -0.5)); fsi_vis.SetCameraMoveScale(0.05f); fsi_vis.EnableBoundaryMarkers(true); fsi_vis.Initialize(); } // Start the simulation unsigned int output_steps = (unsigned int)round(1 / (out_fps * dT)); unsigned int render_steps = (unsigned int)round(1 / (render_fps * dT)); double time = 0.0; int current_step = 0; ChTimer<> timer; timer.start(); double sinkage; while (time < total_time) { // Get the infomation of the wheel force = actuator->Get_react_force(); torque = motor->Get_react_torque(); w_pos = wheel->GetPos(); w_vel = wheel->GetPos_dt(); angvel = wheel->GetWvel_loc(); sinkage = bzDim + wheel_radius - w_pos.z() + 0.5 * iniSpacing; if (verbose) { std::cout << "slip ratio: " << slip_ratios_str[i] << std::endl; std::cout << "time: " << time << std::endl; std::cout << " wheel position: " << w_pos << std::endl; std::cout << " wheel linear velocity: " << w_vel << std::endl; std::cout << " wheel angular velocity: " << angvel << std::endl; std::cout << " drawbar pull: " << force << std::endl; std::cout << " wheel torque: " << torque << std::endl; std::cout << " sinkage: " << sinkage << std::endl; } if (output) { myFile << time << "\t" << w_pos.x() << "\t" << w_pos.y() << "\t" << w_pos.z() << "\t" << w_vel.x() << "\t" << w_vel.y() << "\t" << w_vel.z() << "\t" << angvel.x() << "\t" << angvel.y() << "\t" << angvel.z() << "\t" << force.x() << "\t" << force.y() << "\t" << force.z() << "\t" << torque.x() << "\t" << torque.y() << "\t" << torque.z() << "\n"; myDBP_Torque << time << "\t" << force.x() << "\t" << torque.z() << "\n"; } if (output && current_step % output_steps == 0) { std::cout << "-------- Output" << std::endl; sysFSI.PrintParticleToFile(out_dir + "/particles"); sysFSI.PrintFsiInfoToFile(out_dir + "/fsi", time); static int counter = 0; std::string filename = out_dir + "/vtk/wheel." + std::to_string(counter++) + ".vtk"; WriteCylinderVTK(filename, wheel_radius, wheel_width, sysFSI.GetFsiBodies()[0]->GetFrame_REF_to_abs(), 100); } // Render SPH particles if (render && current_step % render_steps == 0) { if (!fsi_vis.Render()) break; } // Call the FSI solver sysFSI.DoStepDynamics_FSI(); time += dT; current_step++; } timer.stop(); std::cout << "\nSimulation time: " << timer() << " seconds\n" << std::endl; if (output) { myFile.close(); myDBP_Torque.close(); } return 0; } }
