So I'd like to note first that the methods in *HostSideHelpers.hpp* are a collection of various utilities that DEME needs. It's fine that you use them, but they are not necessarily a complete MBD toolbox.
To answer your question, *RotateQuat *is used to rotate a quaternion. As a quaternion is the representation of an orientation and has no spatial coordinate information, when you rotate it the "rotation point" is irrelevant. You are probably thinking about rotating a point and more specifically, the origin of the local coordinate system you assign to the object in question. If that's what you are talking about, then to rotate it, translate everything so the rotation point goes through the origin, then use Rodrigues' formula to rotate the point (or vector), then translate everything back (undo the first step). This is again, MBD stuff and not specific to DEME. Thank you, Ruochun On Mon, Jul 8, 2024 at 5:08 PM simboys <[email protected]> wrote: > Hi, ruochun > > The RotateQuat(const float4& quat, const float3& axis, const > float& theta) just need to set the rotating axis, so does the rotating > point is (0, 0, 0) default? If I want to modify the rotating point, how to > implement it ? > > Best regards, > Wenxuan > > 在2024年7月1日星期一 UTC+8 16:18:50<Ruochun Zhang> 写道： > >> Hi Wenxuan, >> >> Regarding debugging on Linux, I recommend seeking professional >> assistance. However, I find that VSCode's GUI for stepping through code is >> sufficient for my needs. >> >> For enforcing motion on an entity step-by-step, use >> *DoDynamics(step_size)* or *DoStepDynamics()* instead of >> *DoDynamics(frame_duration)*. This way, you can count the number of >> steps and use if statements to control outputs. The RotatingDrum demo >> provides a good example of this. >> >> The performance impact should be minimal, as DEME is designed for >> fine-grain operations like these to enable more flexible co-simulation. The >> two GPU streams synchronize only with *DoDynamicsThenSync()*, not >> *DoDynamics()*. Also, transferring the explicit position and velocity of >> one or two entities to GPUs at each step introduces only a small and >> manageable overhead, too. >> >> Regarding your script, it doesn't fully implement my specified formula. I >> didn't use *fabs*, and the linear velocity should have a scaling factor >> of 5, which seems to be missing. Even with my specified formula, I can't >> guarantee its correctness (without knowing details like what the coordinate >> system is, and what the "expected results" is). But more realistically, >> it's that I cannot do your research for you. You see, this is becoming more >> of a consulting task than troubleshooting. >> >> I suggest starting simple, using only one of the implemented motions, >> verifying it works as expected, and then gradually adding components of the >> second one to approach your desired outcome. >> >> Thank you, >> >> Ruochun >> >> On Saturday, June 29, 2024 at 2:16:35 AM UTC+8 [email protected] wrote: >> >>> Hi, Ruochun >>> >>> I have tried it, but the geometric motion doesn't seem to achieve >>> the expected result. Can you help me troubleshoot the error? Ths mesh obj >>> file are attached. >>> >>> #include <core/ApiVersion.h> >>> #include <core/utils/ThreadManager.h> >>> #include <DEM/API.h> >>> #include <DEM/HostSideHelpers.hpp> >>> #include <DEM/utils/Samplers.hpp> >>> >>> #include <cstdio> >>> #include <chrono> >>> #include <filesystem> >>> >>> using namespace deme; >>> using namespace std::filesystem; >>> >>> int main() >>> { >>> DEMSolver DEMSim; >>> DEMSim.SetVerbosity(STEP_METRIC); >>> // For general use cases, you want to set the verbosity to INFO: It's >>> also a bit faster than STEP_METRIC. >>> // DEMSim.SetVerbosity(INFO); >>> DEMSim.SetOutputFormat(OUTPUT_FORMAT::CSV); >>> DEMSim.SetOutputContent(OUTPUT_CONTENT::ABSV); >>> DEMSim.SetMeshOutputFormat(MESH_FORMAT::VTK); >>> >>> // If you don't need individual force information, then this option >>> makes the solver run a bit faster. >>> DEMSim.SetNoForceRecord(); >>> >>> // E, nu, CoR, mu, Crr... >>> auto mat_type_mixer = DEMSim.LoadMaterial({{"E", 1e8}, {"nu", 0.3}, { >>> "CoR", 0.6}, {"mu", 0.5}, {"Crr", 0.0}}); >>> auto mat_type_granular = DEMSim.LoadMaterial({{"E", 1e8}, {"nu", 0.3}, { >>> "CoR", 0.6}, {"mu", 0.2}, {"Crr", 0.0}}); >>> // If you don't have this line, then mu between mixer material and >>> granular material will be 0.35 (average of the >>> // two). >>> DEMSim.SetMaterialPropertyPair("mu", mat_type_mixer, mat_type_granular, >>> 0.5); >>> >>> float step_size = 5e-5; >>> const double world_size = 10; >>> const float chamber_height = world_size / 3.; >>> const float fill_height = chamber_height; >>> const float chamber_bottom = -world_size / 2.; >>> const float fill_bottom = chamber_bottom + chamber_height; >>> >>> DEMSim.InstructBoxDomainDimension(world_size, world_size, world_size); >>> DEMSim.InstructBoxDomainBoundingBC("all", mat_type_granular); >>> >>> >>> auto mixer = DEMSim.AddWavefrontMeshObject((GET_DATA_PATH() / >>> "mesh/Chute_Rotating.obj").string(), mat_type_mixer); >>> std::cout << "Total num of triangles: " << mixer->GetNumTriangles() << >>> std::endl; >>> mixer->SetFamily(10); >>> mixer->InformCentroidPrincipal(make_float3(0, 0, 6.793), make_float4(0, >>> 0, 0, 1)); >>> >>> >>> double w1 = 3.14159; >>> double w2 = 1.0; >>> std::string angVel_x = "fabs(1.0 * cos(deme::PI * t))"; >>> std::string angVel_y = "fabs(1.0 * sin(deme::PI * t))"; >>> std::string angVel_z = "deme::PI"; >>> >>> >>> std::string vel_x = "(1.0 * cos(1.0 * t) * sin(deme::PI * t) + sin(1.0 >>> * t) * deme::PI * cos(deme::PI * t))"; >>> std::string vel_y = "(-1.0 * cos(1.0 * t) * cos(deme::PI * t) + sin(1.0 >>> * t) * deme::PI * sin(deme::PI * t))"; >>> std::string vel_z = "0"; >>> >>> DEMSim.SetFamilyPrescribedAngVel(10, angVel_x, angVel_y, angVel_z); >>> DEMSim.SetFamilyPrescribedLinVel(10, vel_x, vel_y, vel_z); >>> >>> >>> >>> >>> float granular_rad = 0.005; >>> // Calculate its mass and MOI >>> float mass = 2.6e3 * 5.5886717; // in kg or g >>> float3 MOI = make_float3(2.928, 2.6029, 3.9908) * 2.6e3; >>> std::shared_ptr<DEMClumpTemplate> template_granular = DEMSim. >>> LoadClumpType(mass, MOI, GetDEMEDataFile("clumps/3_clump.csv"), >>> mat_type_granular); >>> template_granular->Scale(granular_rad); >>> >>> // Track the mixer >>> auto mixer_tracker = DEMSim.Track(mixer); >>> >>> >>> >>> DEMSim.SetInitTimeStep(step_size); >>> DEMSim.SetGravitationalAcceleration(make_float3(0, 0, -9.81)); >>> >>> // Initialize the simulation system >>> DEMSim.Initialize(); >>> >>> path out_dir = current_path(); >>> out_dir += "/DEMdemo_RotatingChute"; >>> create_directory(out_dir); >>> >>> float sim_end = 2.0; >>> unsigned int fps = 20; >>> float frame_time = 1.0 / fps; >>> >>> std::cout << "Output at " << fps << " FPS" << std::endl; >>> unsigned int currframe = 0; >>> >>> >>> // mixer_tracker->SetOriQ(make_float4(0, 0, 0, 1)); >>> for (float t = 0; t < sim_end; t += frame_time) // 0.05迭代一次 >>> { >>> >>> std::cout << "Position: " << mixer_tracker->GetPos()[0] << " " << >>> mixer_tracker->GetPos()[1] << " " << mixer_tracker->GetPos()[2] << std:: >>> endl; >>> std::vector<float> Qua = mixer_tracker->GetOriQ(); // x, y, z, w >>> // std::cout << Qua[0] << " " << Qua[1] << " " << Qua[2] << " " << >>> Qua[3] << std::endl; >>> >>> >>> std::cout << "Frame: " << currframe << std::endl; >>> char filename[200], meshfilename[200], cnt_filename[200]; >>> sprintf(filename, "%s/DEMdemo_output_%04d.csv", out_dir.c_str(), >>> currframe); >>> sprintf(meshfilename, "%s/DEMdemo_mesh_%04d.vtk", out_dir.c_str(), >>> currframe); >>> sprintf(cnt_filename, "%s/Contact_pairs_%04d.csv", out_dir.c_str(), >>> currframe++); >>> >>> DEMSim.WriteMeshFile(std::string(meshfilename)); >>> DEMSim.DoDynamics(frame_time); >>> >>> } >>> >>> std::cout << "DEMdemo_RotatingChute simulating complete!" << std::endl; >>> return 0; >>> } >>> >>> 在2024年6月29日星期六 UTC+8 00:33:30<simboys> 写道： >>> >>>> Hi, Ruochun >>>> >>>> Thank you for your reply, I wiil try it. By the way, can you >>>> share me the method how to debug the code step-by-step in linux? In >>>> addition, if I use the DEMTracker to set the location and rotation of >>>> geometry at each time step, will this affect computational efficiency? And >>>> how to set the location of geometry at each time step? From the >>>> demo, DoDynamics(frameTime) is many times longer than timestep, and if we >>>> set the frameTime by timestep, it will result in a huge data file. >>>> >>>> Best regards, >>>> Wenxuan >>>> >>>> 在2024年6月27日星期四 UTC+8 23:06:58<Ruochun Zhang> 写道： >>>> >>>>> Hi Wenxuan, >>>>> >>>>> 1. When you load a mesh or an analytical object, it's by default in >>>>> family 255, a reserved and automatically fixed family. You probably >>>>> changed >>>>> your mesh into a family that does not have any motion prescription, then >>>>> in >>>>> that case it follows the "normal physics" like all the particles you load >>>>> in, including falling under gravity. If you wish to fix that family, just >>>>> use *SetFamilyFixed*. >>>>> >>>>> 2. This is more complicated. First, we have to understand that in >>>>> DEME, strictly speaking, objects work with their principal frames of >>>>> inertia (centered at CoM, having a diagonal inertia matrix). This is a >>>>> current limitation due to how the inertia matrix-defining APIs are >>>>> designed: you can only input the principal moments of inertia. But this >>>>> also means that setting the frame correctly is not necessary if you >>>>> completely prescribe the motion of that object (a.k.a. this object does >>>>> not >>>>> respond to physics), and you don't care about getting a torque reading. >>>>> However, in your case, I suggest you do so, since you'll prescribe a >>>>> complex motion, and standardizing the frame to use at least the CoM frame >>>>> can save you some head-scratching. I suggest that you make it such that in >>>>> the mesh file you load in, the (0,0,0) point is right at the CoM of this >>>>> mesh; but if not, you can adjust the CoM and principal axes with >>>>> *InformCentroidPrincipal* after loading in the mesh. >>>>> >>>>> After that is taken care of, you need to combine the two rotations >>>>> into a rotation about this object's own CoM and a linear translation of >>>>> this object's CoM (because family motion prescription allows for setting >>>>> these two motions only). Use your MBD knowledge to do that. I don't >>>>> guarantee the correctness of this, but some back-of-the-envelop >>>>> calculation >>>>> seems to put the angular velocity to be (w2×cos(w1×t), w2×sin(w1×t), w1), >>>>> and the linear velocity to be >>>>> 5×(w2×cos(w2×t)×sin(w1×t)+sin(w2×t)×w1×cos(w1×t), >>>>> -w2×cos(w2×t)×cos(w1×t)+sin(w2×t)×w1×sin(w1×t), 0). Again, please verify >>>>> yourself. After that, prescribe the object's family's motion using it. >>>>> >>>>> But if this is too complicated and you can write your code to >>>>> calculate the location and quaternion of this object at any given *t >>>>> *(which >>>>> for some people is easier than mathematically deriving the analytical >>>>> expression for the linear and rotational motions), then you can set this >>>>> object fixed, then use the tracker methods *SetPos *and *SetOriQ *to >>>>> enforce the location and rotation of this object *at each time step*. >>>>> It should achieve the same assuming this object is not moving at a very >>>>> high velocity. >>>>> >>>>> Thank you, >>>>> Ruochun >>>>> >>>>> On Thursday, June 27, 2024 at 2:49:43 PM UTC+8 [email protected] >>>>> wrote: >>>>> >>>>>> Hi, >>>>>> I want to achieve a composite rotational motion of a geometric >>>>>> body, as follows >>>>>> Motion 1: Rotate at a speed of w1 rad/s, with the rotation axis being >>>>>> the z-axis >>>>>> Motion 2: Rotate at a speed of w2 rad/s, with a rotation vector of [1 >>>>>> * cos (w1 * time), 1 * sin (w1 * time), 0], and an action point of [0, 0, >>>>>> 5], as shown in attached figure. >>>>>> [image: rotating_Chute.png] >>>>>> >>>>>> >>>>>> However, the following issues were encountered: >>>>>> 1. After using DEMTracker to track geometry, even if I don't add any >>>>>> motion, the geometry will continue to drop down. How can I fix the >>>>>> geometry? >>>>>> 2. How to achieve the composite rotational motion mentioned above? >>>>>> >>>>>> Best regards, >>>>>> Wenxuan Xu >>>>>> >>>>> -- > 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/358e5ffd-7b54-49a1-88f4-8b230e898af6n%40googlegroups.com > <https://groups.google.com/d/msgid/projectchrono/358e5ffd-7b54-49a1-88f4-8b230e898af6n%40googlegroups.com?utm_medium=email&utm_source=footer> > . > -- Ruochun Zhang Email: [email protected] Email: [email protected] <[email protected]> Tel: 832-353-5111 -- 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/CAHvQpOucTEwH%2BACU53uM3A4vdr8jrOBURzOV4_Eh-vayjXmwYg%40mail.gmail.com.
