XML guide: definition of attributes of rigid objects

Would you please consider the following questions about the attributes of floating objects given in the XML guide (version April 2019, section 2.6)?

These attributes can be computed by GenCase or be custom-defined in the code. In order to specify the values correctly, I would need to know more about their definitions.

The attribute center is defined as the 'gravity center of the rigid object'.
  • Q1: Is this the 'centre of mass'? I presume so, this is just about matching terminologies.
  • Q2: Which reference system should the coordinates be given in? The global one (the domain, as from inside definition) or the local one (the object, as from inside a setmkbound for example)?
  • Q3: Reformulation of Q2. Say that I draw an object after rotating the global reference system (rotate x=...), and do not revert to the original system (no matrixreset) before moving on to the floatings specifications. Will the coordinates of "center" live in the original or in the modified reference system?
  • Q4: Linked to Q1. Does the code takes into account eccentricity when the centre of mass is not the centroid of the figure?
The attribute inertia is defined as 'momentum of inertia of the rigid object'. There is a slip of the pen since it should be 'moment of inertia'.
  • Q5: With respect to which axes are these moments of inertia?
  • Q6: Is the code informed about the principal moments of inertia of objects with a predefined shape such as drawbox, drawsphere?
  • Q7: Similar to Q3. If the object I drew is no longer aligned with the global reference system, which axes are considered for the moments of inertia?
Although I anticipate that the Project Chrono implementation clarifies these question, coding the old way may still be attractive for simple enough cases. So, in your answers, references to what the basic code does do and to what it does not do are equally valued.

Thanks in advance for dealing with this.

Comments

  • Some responses to your questions here:

    Q1: Is this the 'centre of mass'? I presume so, this is just about matching terminologies.
    YES

    Q2: Which reference system should the coordinates be given in? The global one (the domain, as from inside definition) or the local one (the object, as from inside a setmkbound for example)?
    Global ones of course, however this is easy to test to get a faster response ;)

    Q3: Reformulation of Q2. Say that I draw an object after rotating the global reference system (rotate x=...), and do not revert to the original system (no matrixreset) before moving on to the floatings specifications. Will the coordinates of "center" live in the original or in the modified reference system?
    Response in https://dual.sphysics.org/vanilla/discussion/1631/centre-of-mass-of-eccentric-objects-rotated-with-rotateaxis#latest

    Q4: Linked to Q1. Does the code takes into account eccentricity when the centre of mass is not the centroid of the figure?
    Not sure, let us check...

    The attribute inertia is defined as 'momentum of inertia of the rigid object'. There is a slip of the pen since it should be 'moment of inertia'.
    Yes, we have to correct this. Thanks

    Q5: With respect to which axes are these moments of inertia?
    X-axis, Y-axis and Z-axis

    Q6: Is the code informed about the principal moments of inertia of objects with a predefined shape such as drawbox, drawsphere?
    No need, since we use the general formulation to compute any figure, so that you can check that results for a box, for example, will converge to theoretical value with accuracy.

    Q7: Similar to Q3. If the object I drew is no longer aligned with the global reference system, which axes are considered for the moments of inertia?
    X-axis, Y-axis and Z-axis
  • edited June 13
    Thanks for following this up, first of all.

    Re Q2, you have been fast enough to make it worth waiting :-)

    Re Q5 and Q7, the givens for one object I have are the moments of inertia with respect to the principal axes of inertia, which always have the centre of mass as origin.

    Do I then understand it well that the correct procedure to input this object is
    1. to convert principal moments of inertia into those with respect to the lattice/global reference system;
    2. to give these latter moments as an attributes within ' floating'
    3. to let DualSPHysics calculate the moments of inertia wrt the global axes at each time step, in the line of the answer to Q6.
    This link to moments of inertia on wikipedia is for the sake of a common terminology. Corrections and integrations always welcome.
  • In order to compute the initial moment of inertia:
    a) you can allow GenCase to compute the value according to the moment of inertia of a set of points with mass
    b) you can define the values you know in the XML

    however, after this time=0, the code will compute the rotations using that moment of inertia, not sure what you mean with "DualSPHysics calculate the moments of inertia wrt the global axes at each time step"?
  • edited June 14
    Thanks Alex. I try to explain.

    From Q5/Q7 I understood that we impose the conservation of angular momentum with respect to the lattice/global axes at each time step. Logical.

    If the moments of inertia are referred to the principal axes of inertia, these are properties/attributes of the object and need not to be recomputed at each time step. However, as time goes by, the solver would need to convert these object attributes into the moments with respect of the global axes, based on the motion and rotation of the object, again to enforce the conservation of angular momentum around those global axes.

    If the moments of inertia are given relative to the global reference axes, and the object roams and rotates, the moments of inertia change at each time step. This because the distribution of masses with respect to those global axes has changed. And the code should take care of this. Now I understand that is the way DSPH works.

    This was my (perhaps confused/confusing) thinking behind "DualSPHysics calculate the moments of inertia wrt the global axes at each time step". But this topic comes closer to the other companion post https://dual.sphysics.org/vanilla/discussion/1631/centre-of-mass-of-eccentric-objects-rotated-with-rotateaxis , so I will continue it there.
  • edited June 21

    After some testing around, I felt I have to resume with this thread. Please refer to the cubes in the figure.

    At the design stage, they have the same mass of 1kg, the edge length is 1m and let's suppose they have uniform density. Out of these properties, the principal axes of inertia are aligned with the edges and pass for the centre of mass, which is also the centre of the figure in this case of uniform density. The principal moments of inertia are also equal because of the object's symmetry. From the books I got that the value of the principal moments of inertia is (mL^2)/6=0.1667 kgm^2.

    In the input stage, I am not setting these moments of inertia myself, but leave it to DualSPHysics to compute their values. I can set the properties of these objects as floatings in the input xml file (Case_Def.xml) and read the matrix of moments of inertia in the output xml file (Case.xml).

    The white cube at the origin is the test case. I expect the inertia matrix to be [I 0 0; 0 I 0; 0 0 I] where I is a particle-based approximation of the moment of the inertia, a single value in this case. Thanks for the answer to Q6 for clarifying the approximate nature of the calculation. With the particle resolution I chose, I=0.2 kg/m^2 instead of 0.1667. I can indeed recognize that the computed value tends to the expected one with the reduced particle distance. So this 0.2 is OK. The computed inertia matrix is [0.2 0 0; 0 0.2 0; 0 0 0.2]. Take it as ground truth.

    Look now at the blue cubes. They have been translated away from the origin without rotation. The inertia matrix from the output xml is again [0.2 0 0; 0 0.2 0; 0 0 0.2] for all cubes. Hence I conclude that the moments of inertia reported by DualSPHysics are those with respect to the principal axes of inertia, thus the local ones, in contrast with answer to Q5/Q6. If they were the moments of inertia with respect to the lattice reference system, some of all of the moments should have increased because of the greater distance between the object and the axes. Also we should have some cross-moments off the diagonal, for some axes are no longer axes of symmetry for the cubes.

    Let's look now at the green cubes. They have been both translated and rotated by different angles, namely 30, 45, 75 around the vertical axis through the centre of mass. If the conclusion above still held --- that the moments of inertia are with respect to the principal axes of inertia ---, I should get the same inertia matrix [0.2 0 0; 0 0.2 0; 0 0 0.2]. This is not the case though. For example, the value for the light green cube is [0.28 0 0 ; 0 .28 0 ; 0 0.28 0].

    This then lead me to think that the moments of inertia are perhaps computed with respect to a system of axes parallel to the lattice axis and passing through the centre of mass/figure. Then I would have expected either the same or different inertia matrices depending on different angles. However, they are the same for 30 and 45 degrees [0.28 0 0 ; 0 .28 0 ; 0 0.28 0], while for 75 degrees the matrix is [0.2 0 0; 0 0.2 0; 0 0 0.2], the same as the non-rotated cubes.

    Also, this situation holds even when I halve the interparticle distance. It might well be a presentational glitch. At this very moment I am not quite sure what I am looking at, and risks of great confusion follow. Beside possible flaws in my thinking, there is an evident lack of consistency in the results. I am available to send the input and output xml file if this is useful.

    On the positive note, the calculation of the centres of mass is fine across these experiments.

    Thanks for any advice and clarification.

  • Good point: "with respect to a system of axes parallel to the lattice axis and passing through the centre of mass/figure"

    Then I was wrong in the responses I wrote in the other thread, sorry.

    Perhaps @jmdalonso and @RCanelas can help here


    Regards

    Alex

  • edited June 23

    Thanks Alex for the clarification.

    Some additional information: I used interparticle distances of 0.1 and 0.05, and even with dx=0.01 I obtain results of the same sort.

    I also thought that, if the moments of inertia are computed for a local system of axes parallel to the global ones and passing through the centre of mass, the inertia tensors of all cubes in the figure should be the same.

    If I am not wrong, this should be so because, given a direction r passing from the centre of mass, the moment of inertia with respect to that direction Ir is given by Ir = Ia cos² alpha + Ib cos² beta + Ic cos² gamma. There, a,b,c are the principal axes of inertia; and alpha, beta, gamma are the direction cosines of r with respect to the axes a,b,c.

    In this case of cubes rotated around the z axis, thus on the xy plane, it holds cos beta = sin alpha, and cos gamma=0 for any angle of rotation around z, say alpha. Then, since Ia=Ib=Ic=I because of a cube's symmetry, it follows Ir=I too. If that' s correct, this may help check the goodness of the results too.

  • Hi

    I have tested the GenCase for calculation of canter of mass and the inertia matrix.

    I have defined a cube with the "massbody" of 60.

    <drawbox>

                <boxfill>solid</boxfill>

                <point x="1.0" y="2.0" z="1.5" />

                <size x="5" y="1.5" z="1.8" />

              </drawbox>


    The GenCase accurately computes the center of Mass, that is :

    <center x="3.5" y="2.75" z="2.4" units_comment="metres (m)" />


    By hand calculation, the inertia matrix is:

    v11=Ixx=826.8 , v12=-Ixy=-577.5 , v13=-Ixz=-8.4

    v21=Iyx=-577.5 , v2=Iyy=1221.8 , v23=-Iyz=-396

    v31=Izx=-8.4 , v32=-Izy=-396 , v33=Izz=1325


    But, the GenCase computes as follow:

    <inertia units_comment="kg*m^2">

              <values v11="29.1" v12="-2.42655e-16" v13="4.18149e-16" />

              <values v21="-2.42655e-16" v22="144.6" v23="-9.34394e-17" />

              <values v31="4.18149e-16" v32="-9.34394e-17" v33="139.5" />

            </inertia>



    What is wrong?


    (The xml file is attached.)

  • I tried several times with a Box ... and results with GenCase are good

    Try please with a rectangular box and compare with theoretical data

  • In the previous comment, a rectangular box is used.

  • You have not computed by hand in the proper way.....

    Note that theoretically Izz=M(a^2+b^2)/12.... so that 136.25 kg/m^2

    and Izz(GenCase)=139.5 not bad if using individual particles

  • I computed the inertia with respect to the xyz axis, and not to the center of the mass of the object.

    So, we can conclude that the GenCase computes the inertia matrix with respect to the center of the mass of each object.

    Since for solving the rotation motion of an object, how does the code deal with the inertia values with respect to the xyz axis?

  • The code computes motions and rotations respect to the center of the mass of each object.

    You can see the source code in function RunFloating in JSphCpuSingle.cpp

  • edited July 7

    To summarise the progress of this long thread:

    • in DualSPHysics without Project Chrono the moments of inertia of a floating object are referred to a reference frame: 1. with origin in the centre of mass of the object; 2. with axes parallel to the global reference system, that of the lattice in fact;
    • These moments of inertia can be computed by GenCase within an approximation using the particle distribution that approximates the floating; alternatively, they can be customized with the attribute fullinertia of floating in the xml file.
    • During the simulations DualSPHysics uses the customized values, if given; else the computed ones.

    Please @Alex check whether we are on the same page and correct, for benefit of future readers.

    The open question of this thread is the following. The principal matrix of inertia of a cube is an isotropic tensor. The inertia matrix of a cube takes a single value with respect to any reference system passing from the centre of mass, regardless of how you rotate the cube. Computationally, though, for some rotations the moments of inertia of the same cube turn out to take different values, at the same particle resolution. I am aware that this issue is to be triaged in due course, so let's wait for confirmation/confutation.

  • Thanks for the clarification!

    What is the difference when we want to implement "Chrono"?

    Regards

  • edited July 8

    @amir It is a very valid question indeed. I would possibly open another thread, for future readers to find relevant answers to your question easily. Thanks for sharing your interest in this topic.

Sign In or Register to comment.