CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM > OpenFOAM Pre-Processing

Inertia in rigidBodyMotion

Register Blogs Community New Posts Updated Threads Search

Like Tree14Likes

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   October 13, 2019, 20:39
Default Inertia in rigidBodyMotion
  #1
Member
 
Morteza
Join Date: Jan 2018
Posts: 30
Rep Power: 8
mortezahdr is on a distinguished road
Hello Foamers,

I have a question about the inertia entry in dynamicMeshDict for a rigidBodyMotion solver:

type rigidBody;
parent root;

// Cuboid mass -- very simple setup. For more complicated shapes
// see utility surfaceInertia
mass 0.0005;
rho 0.1;
inertia (0.0002 0 0 0.01 0 0.01);
centreOfMass (-0.11 -0.1 0.008);

How can I calculate inertia for a complex geometry? The surfaceInertial tool does not give such an information.

Thanks.
mortezahdr is offline   Reply With Quote

Old   October 14, 2019, 16:30
Default
  #2
Member
 
Morteza
Join Date: Jan 2018
Posts: 30
Rep Power: 8
mortezahdr is on a distinguished road
The answer to my question:
The surfaceInertia tool gives the following matrix as "Inertia tensor around centre of mass":
(Ixx Ixy Ixz Iyx Iyy Iyz Izx Izy Izz)
The entry "inertia" in rigidBodyMotion dynamicMeshDict requires six elements of the aforementioned matrix:
(Ixx Ixy Ixz Iyy Iyz Izz)

The surfaceInertia tool works like:
surfaceInertia surface.stl

In order to extract a STL file out of a patch existing in your grid you may use:
surfaceMeshTriangulate -patches '(patchName)' output.stl
mortezahdr is offline   Reply With Quote

Old   March 4, 2021, 18:02
Unhappy
  #3
Member
 
Deutschland
Join Date: Sep 2020
Posts: 69
Rep Power: 5
vava10 is on a distinguished road
Hey,

I am simulating a ship motion using dynamicMEsh and I tried surfaceInertia and following result

mass 0.208208201605;
centreOfMass (2.63037231059 -1.17234355678e-05 0.0244737162072);
momentOfInertia (0.00199990142845 0.225304731178 0.2257076848);

But I am given a mass of 77 kg
I another thread posted in cfd forum it was said the moment of inertia is mass moment of inertia. Id that so?

If yes, is it possible to find it from the momentOfInertia (0.00199990142845 0.225304731178 0.2257076848); or using surfaceInertia

thanks in advance
vava10 is offline   Reply With Quote

Old   March 4, 2021, 18:11
Default
  #4
Member
 
Morteza
Join Date: Jan 2018
Posts: 30
Rep Power: 8
mortezahdr is on a distinguished road
Hello,

You should give the density value when you run surface inertia. For example:
surfaceInertia -density 5 file.stl
This way, you will get the correct mass, if your density*volume (of the STL file) results in 77.
vava10 likes this.
mortezahdr is offline   Reply With Quote

Old   March 13, 2021, 12:28
Default Mass property
  #5
Member
 
Deutschland
Join Date: Sep 2020
Posts: 69
Rep Power: 5
vava10 is on a distinguished road
Hey,


THANK YOU FOR THE REPLY

I am trying to simulate a ship flow. I have few questions regarding dynamic mesh. Hope someone can help.

1. in DTCHull case entire ship is considered in the snappyHexDict even though one half is simulated. Does this not affect the parameters in the entire case setup? So is it fine if I use my entire geometry rather than cutting the geometry into half and then using it? If I use the entire geometry wall all think should I pay attention to?

2. why is the rhoInf in the dynamic mesh 1. should it not be the density of the ship (mass/volume). I have the weight for my geometry (77 kg) and I got the volume using the Solidwork. and used to find density for my case.

4. The moment of inerta. there are 3 moment of inertia
a. Principal axes of inertia and principal moments of inertia: ( kilograms *
square meters ) Taken at the center of mass.
b.Moments of inertia: ( kilograms * square meters )Taken at the center of
mass and aligned with the output coordinate system
c. Moments of inertia: ( kilograms * square meters )
Taken at the output coordinate system.


My guess for the moment of inertia was about the center of mass since that is what asked in the dynamic mesh. Bot I try to clarify it using Solidwork, I am not getting the same moment of inertia nor mass using the rhoInf as 1
so I tried to find the density using the mass in the dynamicMeah and the volume I got from the software. But I still don't have answers


Please please help. I am lost and out of ideas

Kind regards
vava10
vava10 is offline   Reply With Quote

Old   March 16, 2021, 04:30
Default
  #6
Senior Member
 
Join Date: Sep 2013
Posts: 353
Rep Power: 20
Bloerb will become famous soon enough
1) snappyHexMesh uses the stl file of the ship and cuts this out of the mesh block. Since the mesh block covers only half the rest isn't used. You could use half, but it is easier to use the full ship, since this makes it easier to detect edges in the meshing run. E.g if your half ship stl file has one point at x=1e-15 and your mesh goes to x=0 you might get artifacts during meshing.


2) rhoInf is there because you are using basically two solvers. One for the rigid body and one for the flow simulation. Coupled together on a patch. Incompressible openfoam solvers do not use a density. They instead solve for p/rho. Hence you need to add rhoInf so that the rigid body solver can calculate the fluid forces correctly. rhoInf is hence the density of the fluid. I haven't looked in the code, but this is probably not read for solvers which use a density field. Like interFoam etc.


3) no 3 question asked


4) There is one moment of inertia. You can however transform it into a different coordinate system. Which basically means rotating your part until x y z align with principle axis. You obviously need to calculate it here using the coordinate system you are calculating in. And around the centre of mass.



surfaceInertia -density 1 file.stl calculates the moment of inertia of the stl file around it's centre of mass with density 1.

surfaceInertia -density 1 -referencePoint '(0 0 1)' file.stl calculates the moment of inertia of the stl file around this point
tfuwa, marxioxyz, zyfsoton and 1 others like this.
Bloerb is offline   Reply With Quote

Old   March 16, 2021, 04:52
Default
  #7
Senior Member
 
Join Date: Sep 2013
Posts: 353
Rep Power: 20
Bloerb will become famous soon enough
I should also add, that the inertia of the stl file is not necessarily the one used in the calculation. This would mean that the ship is of uniform density. But in reality the engine/cargo etc are heavier. Hence you might not get the exact value just from the stl file since this is only used for meshing.
vava10 likes this.
Bloerb is offline   Reply With Quote

Old   March 16, 2021, 05:32
Default inertia
  #8
Member
 
Deutschland
Join Date: Sep 2020
Posts: 69
Rep Power: 5
vava10 is on a distinguished road
Thank you so much for the reply I have been stuck soooooo hard

I was finding the moment of inertia of the ship from the stl (MeshLab, Solidworks and surfaceInertia) and using it for the calculation. My kayak was sinking

what I have is a kayak and it does not have any additions cargo or engines.
and I have 77 kg as the weight of the kayak. I obtained the volume of the mesh using MeshLAb and used this volume and given mass to find the density.

This density was used for the calculation of the moment of inertia

If the inertia from the STL file is not used for the calculation, that probably can be the reason for sinking?

How exactly should I proceed for finding the moment of inertia for the calculation purpose?

Should I just use the moment of inertia for the bounding box? if so I need mass for that. should I use 77kg for the mass??



Again thank you
vava10 is offline   Reply With Quote

Old   March 16, 2021, 08:56
Default
  #9
Senior Member
 
Join Date: Sep 2013
Posts: 353
Rep Power: 20
Bloerb will become famous soon enough
Just use the inertia of your cad program or the openfoam surfaceInertia as demonstrated above.
And no the inertia is just a measure for how difficult it is to change a systems velocity. It has however no influence whatsover on a boat sinking. (well it makes it easier for the boat to roll over and fill up)


It will sink if the integral of the pressure around the boat (it's lift force) is lower than it's weight force (mass*gravity). The mass of your rigidBody is specified in the dynamicmeshdict and the lift force can be calculated using the forces function object. Or you can simply calculate how much water it displaces. Since this is what it boils down to and should equal.
vava10 likes this.
Bloerb is offline   Reply With Quote

Old   March 16, 2021, 09:51
Post Mass property
  #10
Member
 
Deutschland
Join Date: Sep 2020
Posts: 69
Rep Power: 5
vava10 is on a distinguished road
Thank you for the help

1.So Since the mesh block covers only half the rest isn't used, does that mean the mass properties of the complete stl obtained using surfaceInertia -density <scalar> should be divided by 2 to use it in the dynamicMeshDict?

P:S I am obtaining density by dividing given mass/volume obtained using MeshLab

2.Incompressible openfoam solvers do not use a density. They instead solve for p/rho. Hence you need to add rhoInf so that the rigid body solver can calculate the fluid forces correctly. rhoInf is hence the density of the fluid.. In the controlDict for finding the forces i am using rhoInf of 1000.

But the rhoInf in dynamicMEsh seems different because I checked the DTCHull case setup provided by openFoam. In dynamicMeshDict rhoInf is 1


Code:
/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  v2012                                 |
|   \\  /    A nd           | Website:  www.openfoam.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      dynamicMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dynamicFvMesh   dynamicMotionSolverFvMesh;

motionSolverLibs (sixDoFRigidBodyMotion);

motionSolver    sixDoFRigidBodyMotion;

patches         (hull);
innerDistance   0.3;
outerDistance   1;

centreOfMass    (2.929541 0 0.2);
mass            412.73;
momentOfInertia (40 921 921);
rhoInf          1;
report          on;

value           uniform (0 0 0);

accelerationRelaxation 0.4;

solver
{
    type Newmark;
}

constraints
{
    zAxis
    {
        sixDoFRigidBodyMotionConstraint line;
        direction     (0 0 1);
    }
    yPlane
    {
        sixDoFRigidBodyMotionConstraint axis;
        axis          (0 1 0);
    }
}

restraints
{
    translationDamper
    {
        sixDoFRigidBodyMotionRestraint linearDamper;
        coeff         8596;
    }
    rotationDamper
    {
        sixDoFRigidBodyMotionRestraint sphericalAngularDamper;
        coeff         11586;
    }
}


// ************************************************************************* //
vava10 is offline   Reply With Quote

Old   March 16, 2021, 10:07
Default
  #11
Senior Member
 
Join Date: Sep 2013
Posts: 353
Rep Power: 20
Bloerb will become famous soon enough
The rigid motion solver basically solves F=m*a. With m the mass inside your dynamicMeshDict. a the acceleration you want to solve for. The same applies for it's velocity v and the distance s. Forces can be either from constraints/restraints in your dynamicMeshDict or fluid forces from the flow calculation. In other words it does nothing else than a forces function object. Or even more precise it even uses the forces function object library internally to calculate the force. rhoInf is hence the same one used in the actual function object. And this is only read if you are using an incompressible solver without a density field (e.g simpleFoam pisoFoam pimpleFoam icoFoam ). If you are running interFoam this entry is not used but rather the solvers density field. It calculates the resulting force vector and point of attack (the moment) and does so on the specified patch.

This is the implementation
Code:
        dictionary forcesDict;

        forcesDict.add("type", functionObjects::forces::typeName);
        forcesDict.add("patches", patches_);
        forcesDict.add("rhoInf", rhoInf_);
        forcesDict.add("rho", rhoName_);
        forcesDict.add("CofR", motion_.centreOfRotation());

        functionObjects::forces f("forces", db(), forcesDict);
        f.calcForcesMoment();

        motion_.update
        (
            firstIter,
            ramp*(f.forceEff() + motion_.mass()*g.value()),
            ramp
           *(
               f.momentEff()
             + motion_.mass()*(motion_.momentArm() ^ g.value())
            ),
            t.deltaTValue(),
            t.deltaT0Value()
        );
It constructs a forces function object uses the function objects calcForcesMoment routine and supplies this to the motion solver for the rigid body update.


Hence symmetry is not baked into the rigid motion solvers. You need to half the weight.
vava10 likes this.
Bloerb is offline   Reply With Quote

Old   March 16, 2021, 10:23
Default
  #12
Member
 
Deutschland
Join Date: Sep 2020
Posts: 69
Rep Power: 5
vava10 is on a distinguished road
I am using interFoam solver with sixDoFRigidBodyMotion so rhoInf 1000?

(I am sorry if the questions are stupid but I am completely new to CFD. I appreciate your help GREATLY)



And I also have a point mass at the centre of gravity representing the passenger weight. How will this include into dynamic mesh?


77 is the combined weight of boat and passenger. I was using this value for the calculation and this could be the reason why it was sinking
vava10 is offline   Reply With Quote

Old   March 16, 2021, 10:36
Default
  #13
Senior Member
 
Join Date: Sep 2013
Posts: 353
Rep Power: 20
Bloerb will become famous soon enough
Start reading my posts above. InterFOAM does not care about rhoInf. InterFoam uses an internal density field. Hence pressure is in Pa not divided by density. You can enter any value you want. 1, 1000 100000. It does not matter.
And your second question isn't even worth answering since you clearly do not even try. You are using a single rigid body. It hence has ONE mass. Hence mass in your dynamicMeshDict is 77. But you need to divide it by half if you are using a symmetryPlane.

And if you have a passenger you clearly need to add that to the moment of inertia.
vava10 likes this.
Bloerb is offline   Reply With Quote

Old   March 16, 2021, 10:57
Default
  #14
Member
 
Deutschland
Join Date: Sep 2020
Posts: 69
Rep Power: 5
vava10 is on a distinguished road
Quote:
Originally Posted by Bloerb View Post
Start reading my posts above. InterFOAM does not care about rhoInf. InterFoam uses an internal density field. Hence pressure is in Pa not divided by density. You can enter any value you want. 1, 1000 100000. It does not matter.
And your second question isn't even worth answering since you clearly do not even try. You are using a single rigid body. It hence has ONE mass. Hence mass in your dynamicMeshDict is 77. But you need to divide it by half if you are using a symmetryPlane.

And if you have a passenger you clearly need to add that to the moment of inertia.
Thank you for your help. I was trying to find the reason for my ship was sinking. I did put the moment of Inertia considering for the half geometry and weight 77/2. But it was sinking. I thought the error could be in the mass properties or rhoInf. I just wanted to clarify.
vava10 is offline   Reply With Quote

Old   March 12, 2022, 15:47
Default
  #15
Member
 
Tony Zhang
Join Date: Nov 2019
Location: soton
Posts: 45
Rep Power: 6
zyfsoton is on a distinguished road
Quote:
Originally Posted by Bloerb View Post
Start reading my posts above. InterFOAM does not care about rhoInf. InterFoam uses an internal density field. Hence pressure is in Pa not divided by density. You can enter any value you want. 1, 1000 100000. It does not matter.
And your second question isn't even worth answering since you clearly do not even try. You are using a single rigid body. It hence has ONE mass. Hence mass in your dynamicMeshDict is 77. But you need to divide it by half if you are using a symmetryPlane.

And if you have a passenger you clearly need to add that to the moment of inertia.
Dear All,

I just came across the same problem. I am simulating a ship in calm water condition using OpenFOAMv7. I have successfully simulated the hull in fixed condition and now want to make hull free to heave and pitch.

I have changed the dynamicMeshDict and also the centre of rotation located in controlDict file. When I enter interFoam, it seems ok from the beginning but then the time step drops to a very small value, and linear & angular velocity became very big, then the error is reported.

Code:
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  7
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
Build  : 7
Exec   : interFoam
Date   : Mar 12 2022
Time   : 20:22:11
Host   : "Tony"
PID    : 29051
I/O    : uncollated
Case   : /home/OpenFOAM/tony-7/run/soton-snappDone
nProcs : 1
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster (fileModificationSkew 10)
allowSystemOperations : Allowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 0

Selecting dynamicFvMesh dynamicMotionSolverFvMesh
Selecting motion solver: rigidBodyMotion
Selecting rigidBodySolver Newmark

PIMPLE: No convergence criteria found


PIMPLE: No corrector convergence criteria found
        Calculations will do 3 corrections


PIMPLE: Operating solver in transient mode with 3 outer correctors


Reading field p_rgh

Reading field U

Reading/calculating face flux field phi

Reading transportProperties

Selecting incompressible transport model Newtonian
Selecting incompressible transport model Newtonian
Selecting turbulence model type RAS
Selecting RAS turbulence model kOmegaSST
Selecting patchDistMethod meshWave
RAS
{
    RASModel        kOmegaSST;
    turbulence      on;
    printCoeffs     on;
    alphaK1         0.85;
    alphaK2         1;
    alphaOmega1     0.5;
    alphaOmega2     0.856;
    gamma1          0.555556;
    gamma2          0.44;
    beta1           0.075;
    beta2           0.0828;
    betaStar        0.09;
    a1              0.31;
    b1              1;
    c1              10;
    F3              false;
}


Reading g

Reading hRef
Calculating field g.h

No MRF models present

No finite volume options present
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000693716, No Iterations 14
time step continuity errors : sum local = 6.13717e-11, global = 1.60616e-11, cumulative = 1.60616e-11
Constructing face velocity Uf

Courant Number mean: 0.000801996 max: 0.151765
forces totForce:
    Not including porosity effects

Starting time loop

Courant Number mean: 0.000801996 max: 0.151765
Interface Courant Number mean: 0 max: 0
deltaT = 0.000119999
Time = 0.000119999

PIMPLE: Iteration 1
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.84309e-07 0.170527)
    Orientation: (1 0 -1.50333e-06 0 1 0 1.50333e-06 0 1)
    Linear velocity: (0 0 0.0468714)
    Angular velocity: (0 0.0250556 0)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000752343, No Iterations 12
time step continuity errors : sum local = 8.96779e-11, global = 2.75192e-11, cumulative = 4.35808e-11
smoothSolver:  Solving for alpha.water, Initial residual = 9.59931e-08, Final residual = 7.56999e-12, No Iterations 1
Phase-1 volume fraction = 0.536621  Min(alpha.water) = 0  Max(alpha.water) = 1
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.536621  Min(alpha.water) = -2.60186e-10  Max(alpha.water) = 1
GAMG:  Solving for p_rgh, Initial residual = 1, Final residual = 4.30904e-08, No Iterations 16
time step continuity errors : sum local = 9.00489e-12, global = -5.89621e-12, cumulative = 3.76846e-11
PIMPLE: Iteration 2
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.84309e-07 0.170258)
    Orientation: (1 0 0.000142199 0 1 0 -0.000142199 0 1)
    Linear velocity: (0 0 -4.43959)
    Angular velocity: (0 -2.37 0)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000771057, No Iterations 12
time step continuity errors : sum local = 1.17222e-09, global = 8.5936e-10, cumulative = 8.97045e-10
smoothSolver:  Solving for alpha.water, Initial residual = 3.61728e-06, Final residual = 1.1771e-11, No Iterations 4
Phase-1 volume fraction = 0.536621  Min(alpha.water) = -1.78302e-34  Max(alpha.water) = 1.00001
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.536621  Min(alpha.water) = -6.09734e-08  Max(alpha.water) = 1.00001
GAMG:  Solving for p_rgh, Initial residual = 0.00425161, Final residual = 4.07793e-08, No Iterations 18
time step continuity errors : sum local = 2.38151e-11, global = 1.71943e-11, cumulative = 9.14239e-10
PIMPLE: Iteration 3
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.84309e-07 0.177328)
    Orientation: (0.999993 0 -0.0036235 0 1 0 0.0036235 0 0.999993)
    Linear velocity: (-5.55112e-17 0 113.406)
    Angular velocity: (0 60.3923 0)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000625751, No Iterations 12
time step continuity errors : sum local = 2.49812e-08, global = -1.83191e-08, cumulative = -1.74049e-08
smoothSolver:  Solving for alpha.water, Initial residual = 9.46965e-05, Final residual = 9.63334e-11, No Iterations 8
Phase-1 volume fraction = 0.536621  Min(alpha.water) = -3.11261e-43  Max(alpha.water) = 1.00002
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.536621  Min(alpha.water) = -0.000972877  Max(alpha.water) = 1.00002
GAMG:  Solving for p_rgh, Initial residual = 0.00824394, Final residual = 4.5727e-08, No Iterations 19
time step continuity errors : sum local = 3.4876e-10, global = 2.68154e-10, cumulative = -1.71367e-08
smoothSolver:  Solving for omega, Initial residual = 0.000320324, Final residual = 4.87026e-08, No Iterations 7
smoothSolver:  Solving for k, Initial residual = 1, Final residual = 4.00958e-08, No Iterations 12
ExecutionTime = 97.03 s  ClockTime = 97 s

forces totForce write:
    sum of forces:
        pressure : (-2.24726e+07 -64526.1 -2.77865e+08)
        viscous  : (-4.05731 -0.459928 -83.9502)
        porous   : (0 0 0)
    sum of moments:
        pressure : (27959.1 -1.82698e+08 184772)
        viscous  : (-0.210328 -23.3017 0.395045)
        porous   : (0 0 0)

Courant Number mean: 0.00589213 max: 163.91
Interface Courant Number mean: 0.000496804 max: 96.5416
deltaT = 1.83026e-05
Time = 0.000138302

PIMPLE: Iteration 1
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.84309e-07 0.171931)
    Orientation: (1 0 -0.000746945 0 1 0 0.000746945 0 1)
    Linear velocity: (0 0 -703.147)
    Angular velocity: (0 -374.727 0)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000949863, No Iterations 12
time step continuity errors : sum local = 3.45495e-08, global = 2.62545e-08, cumulative = 9.11781e-09
smoothSolver:  Solving for alpha.water, Initial residual = 7.18828e-05, Final residual = 2.83216e-11, No Iterations 13
Phase-1 volume fraction = 0.536621  Min(alpha.water) = -1.27194e-18  Max(alpha.water) = 1.00002
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.536621  Min(alpha.water) = -7.54874e-14  Max(alpha.water) = 1.00002
GAMG:  Solving for p_rgh, Initial residual = 0.00601612, Final residual = 2.78072e-08, No Iterations 20
time step continuity errors : sum local = 2.05845e-10, global = 1.51615e-10, cumulative = 9.26942e-09
PIMPLE: Iteration 2
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.84309e-07 0.472629)
    Orientation: (0.98707 0 -0.16029 0 1 0 0.16029 0 0.98707)
    Linear velocity: (0 0 32155.3)
    Angular velocity: (0 17135.1 0)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000999593, No Iterations 11
time step continuity errors : sum local = 1.67722e-06, global = 1.34222e-06, cumulative = 1.35149e-06
smoothSolver:  Solving for alpha.water, Initial residual = 0.00402176, Final residual = 6.77955e-11, No Iterations 35
Phase-1 volume fraction = 0.536621  Min(alpha.water) = -3.62896e-06  Max(alpha.water) = 1.00009
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.536621  Min(alpha.water) = -3.79914  Max(alpha.water) = 1.08871
#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::sigFpe::sigHandler(int) at ??:?
#2  ? in "/lib/x86_64-linux-gnu/libc.so.6"
#3  Foam::GAMGSolver::scale(Foam::Field<double>&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, Foam::Field<double> const&, unsigned char) const at ??:?
#4  Foam::GAMGSolver::Vcycle(Foam::PtrList<Foam::lduMatrix::smoother> const&, Foam::Field<double>&, Foam::Field<double> const&, Foam::Field<double>&, Foam::Field<double>&, Foam::Field<double>&, Foam::Field<double>&, Foam::Field<double>&, Foam::PtrList<Foam::Field<double> >&, Foam::PtrList<Foam::Field<double> >&, unsigned char) const at ??:?
#5  Foam::GAMGSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const at ??:?
#6  Foam::fvMatrix<double>::solveSegregated(Foam::dictionary const&) at ??:?
#7  Foam::fvMatrix<double>::solve(Foam::dictionary const&) in "/home/OpenFOAM/OpenFOAM-7/platforms/linux64GccDPInt32Opt/bin/interFoam"
#8  Foam::fvMatrix<double>::solve() in "/home/OpenFOAM/OpenFOAM-7/platforms/linux64GccDPInt32Opt/bin/interFoam"
#9  ? in "/home/OpenFOAM/OpenFOAM-7/platforms/linux64GccDPInt32Opt/bin/interFoam"
#10  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
#11  ? in "/home/OpenFOAM/OpenFOAM-7/platforms/linux64GccDPInt32Opt/bin/interFoam"
Floating point exception (core dumped)
The dynamicMeshDict is shown:


Code:
/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  7
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      dynamicMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dynamicFvMesh   dynamicMotionSolverFvMesh;

motionSolverLibs   ("librigidBodyMeshMotion.so");

motionSolver    rigidBodyMotion;

report          on;

solver
{
    type Newmark;
}

accelerationRelaxation 0.4;

bodies
{
    hull
    {
        type            rigidBody;
        parent          root;

        centreOfMass    (1.86617058305 -4.84309040657e-07 0.170524016331);
        mass            229.5;
        inertia         (-0.0122617483413 -3.05290279692e-07 0.00566631344423 -0.430050980166 -3.19195024184e-08 -0.435116910148);
        transform       (1 0 0 0 1 0 0 0 1) (1.86617058305 -4.84309040657e-07 0.170524016331);

        joint
        {
            type    composite;
            joints
            (
                {
                    type Pz;
                }
                {
                    type Ry;
                }
            );
        }

        patches         (hull);
        innerDistance   0.3;
        outerDistance   1;
    }
}

restraints
{
    translationDamper
    {
        type linearDamper;
        body hull;
        coeff 8596;
    }

    rotationDamper
    {
        type sphericalAngularDamper;
        body hull;
        coeff 11586;
    }
}


// ************************************************************************* //
Also, when using surfaceInertia utility, the mass is negative and also reported: Non-unique eigenvectors.

Code:
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  7
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
Build  : 7
Exec   : surfaceInertia /home/OpenFOAM/Tony-7/run/soton-snappDone/constant/triSurface/KCS_SOTON_ModelScale.stl
Date   : Mar 12 2022
Time   : 20:03:08
Host   : "Tony"
PID    : 27323
I/O    : uncollated
Case   : /home/
nProcs : 1
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster (fileModificationSkew 10)
allowSystemOperations : Allowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
--> FOAM Warning : 
    From function int main(int, char**)
    in file surfaceInertia.C at line 107
    Negative mass detected, the surface may be inside-out.
--> FOAM Warning : 
    From function int main(int, char**)
    in file surfaceInertia.C at line 324
    Non-unique eigenvectors, cannot compute transformation from Cartesian axes

Density: 1
Mass: -0.4534574004
Centre of mass: (1.86617058305 -4.84309040657e-07 0.170524016331)
Surface area: 5.76446631858
Inertia tensor around centre of mass: 
(-0.0122617483413 -3.05290279692e-07 0.00566631344423 -3.05290279692e-07 -0.430050980166 -3.19195024184e-08 0.00566631344423 -3.19195024184e-08 -0.435116910148)
eigenValues (principal moments): (-0.435192825846 -0.430050980166 -0.0121858326428)
eigenVectors (principal axes): 
(-0.0133965207656 5.41183306631e-06 0.999910262575)
(8.03987066635e-07 0.999999999985 -5.40154715773e-06)
(0.999910262589 -7.31552980241e-07 0.0133965207697)

Writing scaled principal axes at centre of mass of "/home/OpenFOAM/Tony-7/run/soton-snappDone/constant/triSurface/KCS_SOTON_ModelScale.stl" to "axes.obj"

End
The controlDict file is as follows:
Code:
/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.0                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application     interFoam;

startFrom       startTime;

startTime       0;

stopAt          endTime;

endTime         50;

deltaT          0.0001;

writeControl    adjustableRunTime;

writeInterval   5;

purgeWrite      0;

writeFormat     binary;

writePrecision  6;

writeCompression off;

timeFormat      general;

timePrecision   6;

runTimeModifiable yes;

adjustTimeStep  yes;

maxCo           25;
maxAlphaCo      15;
maxDeltaT       0.01;


functions
{
    ///ALLFORCE_SCOMM    
    totForce
    {
        type forces;
        functionObjectLibs ( "libforces.so" );
        patches (hull);
        rhoInf  999.10260000000000000000; 
        rho rho; 
        UName   U;
        log     on;
        writeControl   timeStep;
        writeInterval  1;
        CofR (1.86617058305 -4.84309040657e-07 0.170524016331); 
    }

   
}
Please let me know how this error can be resolved. Any hint will be really appreciated!

Many thanks,
Tony
zyfsoton is offline   Reply With Quote

Old   March 14, 2022, 11:17
Default
  #16
Senior Member
 
Join Date: Sep 2013
Posts: 353
Rep Power: 20
Bloerb will become famous soon enough
Your surfaceInertia results in a negative value, because your stl is inside out. I.e you need to flip the normals of the stl file. For example:

Code:
surfaceOrient input.stl '(10 10 10)' output.stl
With a point outside the geometry or with the -inside option inside as well. But that is more of a cosmetic thing....the inertia might be wrong though.
From the limited info I have another real problem is this
Code:
 maxCo           25;
maxAlphaCo      15;
Are you really trying to solve a coupled problem with these settings? Reduce those below 1. Atleast until it is solving. And from your log file....
Code:
.......
Courant Number mean: 0.000801996 max: 0.151765
Courant Number mean: 0.000801996 max: 0.151765
Interface Courant Number mean: 0 max: 0
deltaT = 0.000119999
Time = 0.000119999

PIMPLE: Iteration 1
    Centre of rotation: (1.86617 -4.84309e-07 0.170527)
    Orientation: (1 0 -1.50333e-06 0 1 0 1.50333e-06 0 1)
    Linear velocity: (0 0 0.0468714)
    Angular velocity: (0 0.0250556 0)
    
    Centre of rotation: (1.86617 -4.84309e-07 0.170258)
    Orientation: (1 0 0.000142199 0 1 0 -0.000142199 0 1)
    Linear velocity: (0 0 -4.43959)
    Angular velocity: (0 -2.37 0)
    
    Centre of rotation: (1.86617 -4.84309e-07 0.177328)
    Orientation: (0.999993 0 -0.0036235 0 1 0 0.0036235 0 0.999993)
    Linear velocity: (-5.55112e-17 0 113.406)
    Angular velocity: (0 60.3923 0)

forces totForce write:
    sum of forces:
        pressure : (-2.24726e+07 -64526.1 -2.77865e+08)
        viscous  : (-4.05731 -0.459928 -83.9502)
        porous   : (0 0 0)
    sum of moments:
        pressure : (27959.1 -1.82698e+08 184772)
        viscous  : (-0.210328 -23.3017 0.395045)
        porous   : (0 0 0)

Courant Number mean: 0.00589213 max: 163.91
Interface Courant Number mean: 0.000496804 max: 96.5416
deltaT = 1.83026e-05
Time = 0.000138302
Your ship is moving at 113 m/s in z direction with an immense angular velocity in the first time step. This is due to on of three things or a mixture of all 3.
1 more outer iterations needed: 3 outer loops are not even close to enough to converge the solution...as you can see your velocity is drastically increasing from outer loop to outer loop...hence increase those drastically and add a convergence criterion
2 Your initial rigid body definition is wrong: check the definition of your dampers, center of rotation and gravity/ the mass etc. Keep in mind the translation you are adding to those with the transform. Your ships weight is really 200 kg?

3 your flow solution is wrong: it could also be that your initial force on the hull is extremely inaccurate because your problem is ill initialized. I.e your first pressure and flow velocity solution is immensely inaccurate which results in a force on the boat which is wrong. Write out the first time step and look at the pressure and drag forces on the hull




Also use the rigid body state function object instead of forces
Code:
 {                                    

 type           rigidBodyState;                                    

 libs           ("librigidBodyState.so");                                    

 angleFormat    degrees;                                

 }
You might also want to add gravity (g) to your rigid body calculation. I.e add a g line to your dict.
Bloerb is offline   Reply With Quote

Old   March 15, 2022, 09:57
Default
  #17
Member
 
Tony Zhang
Join Date: Nov 2019
Location: soton
Posts: 45
Rep Power: 6
zyfsoton is on a distinguished road
[QUOTE=Bloerb;824081]


Dear Bloerb,

Many thanks for your detailed reply!

1) I used the meshLab to invert my STL file and now surfaceInertia gives me a positive value. Since this is a model-scale ship hull (Lpp=3.773m), the mass is 229.5kg.

2) In terms of the controlDict file, I use the tutorial of DTCHullMoving (https://develop.openfoam.com/Develop.../DTCHullMoving) as a base and then change the corresponding parameters. In this tutorial, the maxCo is 25 and maxAlphaCo is 15. But I have changed to 1 and 1 respectively. But I am a little bit confused about this tutorial, it has restraints coefficients in the dynamicDict file, is the DTCHull free to heave and pitch?

3)I have modified the dynamicMeshDict file and controlDict file according to your suggestions and also, add one more out iterations (from 3 to 4). The details of dynamicMeshDict file and controlDict are attached:

Code:
/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  7
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      dynamicMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dynamicFvMesh   dynamicMotionSolverFvMesh;

motionSolverLibs   ("librigidBodyMeshMotion.so");

motionSolver    rigidBodyMotion;

report          on;

solver
{
    type Newmark;
}

accelerationRelaxation 0.4;

bodies
{
    hull
    {
        type            rigidBody;
        parent          root;

        centreOfMass    (1.86617058305 -4.84309040657e-07 0.170524016331);
        mass            229.5;
        inertia         (6.20581170707 0.00015451091794 -2.8677863327 217.653741789 1.61548268875e-05 220.217667153);
        transform       (0.999910262589 -7.31552980241e-07 0.0133965207697 8.03987066635e-07 0.999999999985 -5.40154715773e-06 -0.0133965207656 5.41183306631e-06 0.999910262575) (1.86617058305 -4.84309040657e-07 0.170524016331);

        joint
        {
            type    composite;
            joints
            (
                {
                    type Pz;
                }
                {
                    type Ry;
                }
            );
        }

        patches         (hull);
        innerDistance   0.3;
        outerDistance   1;
    }
}

restraints
{
    translationDamper
    {
        type linearDamper;
        body hull;
        coeff 7217;
    }

    rotationDamper
    {
        type sphericalAngularDamper;
        body hull;
        coeff 9727;
    }
}


// ************************************************************************* //

Code:
/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.0                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application     interFoam; 

startTime 	0; 

startFrom	startTime;

stopAt          endTime;

endTime         30; 

deltaT          0.01;    

writeControl    adjustableRunTime; 

writeInterval   1; 

purgeWrite      0;

writeFormat	binary;

writePrecision  6;

writeCompression uncompressed;

timeFormat      general;

timePrecision   6;

runTimeModifiable yes;

adjustTimeStep  yes;

maxCo           1;
maxAlphaCo      1;
maxDeltaT       0.1;


functions
{
        
    totForce
    {
        type forces;
        functionObjectLibs ( "libforces.so" );
        patches (hull);
        rhoInf  999.1; 
        rho rho; 
        UName   U;
        log     on;
        writeControl   timeStep;
        writeInterval  1;
        CofR (1.86617058305 -4.84309040657e-07 0.170524016331); 
    }

}


// ************************************************************************* //
This time, the simulation ran about 0.85s and then crashed:

Code:
PIMPLE: Iteration 1
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.86048e-07 0.170203)
    Orientation: (0.976646 -1.94577e-06 -0.214856 8.03987e-07 1 -5.40155e-06 0.214856 5.10266e-06 0.976646)
    Linear velocity: (-3.03881e+07 12276 2.26815e+09)
    Angular velocity: (799.724 9.94697e+08 -5372.9)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.00099521, No Iterations 11
time step continuity errors : sum local = 4.88642e-15, global = 8.62449e-16, cumulative = 1.00854e-08
smoothSolver:  Solving for alpha.water, Initial residual = 1.97572e-12, Final residual = 3.06837e-17, No Iterations 1
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
GAMG:  Solving for p_rgh, Initial residual = 0.984401, Final residual = 4.9944e-08, No Iterations 184
time step continuity errors : sum local = 4.04897e-16, global = 4.39492e-18, cumulative = 1.00854e-08
PIMPLE: Iteration 2
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.86048e-07 0.170203)
    Orientation: (0.976646 -1.94577e-06 -0.214856 8.03987e-07 1 -5.40155e-06 0.214856 5.10266e-06 0.976646)
    Linear velocity: (-9.72907e+07 39302.8 7.26173e+09)
    Angular velocity: (2791.19 3.47168e+09 -18752.5)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000900381, No Iterations 11
time step continuity errors : sum local = 2.35324e-15, global = 4.49763e-16, cumulative = 1.00854e-08
smoothSolver:  Solving for alpha.water, Initial residual = 1.42279e-12, Final residual = 3.00484e-17, No Iterations 1
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
GAMG:  Solving for p_rgh, Initial residual = 0.604041, Final residual = 4.93915e-08, No Iterations 108
time step continuity errors : sum local = 2.10258e-15, global = 2.55039e-17, cumulative = 1.00854e-08
PIMPLE: Iteration 3
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.86048e-07 0.170203)
    Orientation: (0.976646 -1.94577e-06 -0.214856 8.03987e-07 1 -5.40155e-06 0.214856 5.10266e-06 0.976646)
    Linear velocity: (-4.93233e+08 199253 3.68147e+10)
    Angular velocity: (13820.6 1.71901e+10 -92853)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000842764, No Iterations 13
time step continuity errors : sum local = 2.97691e-15, global = 2.53722e-16, cumulative = 1.00854e-08
smoothSolver:  Solving for alpha.water, Initial residual = 4.0255e-12, Final residual = 2.98851e-17, No Iterations 1
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
GAMG:  Solving for p_rgh, Initial residual = 0.815495, Final residual = 4.87706e-08, No Iterations 115
time step continuity errors : sum local = 3.47524e-15, global = -4.23773e-17, cumulative = 1.00854e-08
PIMPLE: Iteration 4
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.86048e-07 0.170203)
    Orientation: (0.976646 -1.94577e-06 -0.214856 8.03987e-07 1 -5.40155e-06 0.214856 5.10266e-06 0.976646)
    Linear velocity: (4.16159e+08 -168117 -3.10619e+10)
    Angular velocity: (-10656.3 -1.32544e+10 71594.1)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000775915, No Iterations 11
time step continuity errors : sum local = 4.08494e-15, global = 1.85092e-15, cumulative = 1.00854e-08
smoothSolver:  Solving for alpha.water, Initial residual = 5.98422e-12, Final residual = 4.25305e-17, No Iterations 1
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
GAMG:  Solving for p_rgh, Initial residual = 0.854974, Final residual = 4.8514e-08, No Iterations 107
time step continuity errors : sum local = 6.66147e-15, global = 8.05259e-17, cumulative = 1.00854e-08
smoothSolver:  Solving for omega, Initial residual = 1.57746e-08, Final residual = 1.60001e-11, No Iterations 1
smoothSolver:  Solving for k, Initial residual = 4.55715e-09, Final residual = 7.22209e-13, No Iterations 1
ExecutionTime = 5819.41 s  ClockTime = 5838 s

forces totForce write:
    sum of forces:
        pressure : (2.77316e+33 -1.35636e+33 1.28499e+34)
        viscous  : (6.50696e+11 -1.01531e+11 1.85448e+11)
        porous   : (0 0 0)
    sum of moments:
        pressure : (1.94071e+32 1.41346e+34 1.34365e+33)
        viscous  : (2.87414e+10 3.03959e+11 1.41163e+11)
        porous   : (0 0 0)

Courant Number mean: 7.31397e-08 max: 2.0195
Interface Courant Number mean: 8.31808e-11 max: 2.34747e-05
deltaT = 7.22564e-22
Time = 0.861924

PIMPLE: Iteration 1
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.86048e-07 0.170203)
    Orientation: (0.976646 -1.94577e-06 -0.214856 8.03987e-07 1 -5.40155e-06 0.214856 5.10266e-06 0.976646)
    Linear velocity: (7.71943e+07 -31184.4 -5.76175e+09)
    Angular velocity: (-827.559 -1.02932e+09 5559.92)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000902617, No Iterations 12
time step continuity errors : sum local = 2.92597e-15, global = -8.46416e-17, cumulative = 1.00854e-08
smoothSolver:  Solving for alpha.water, Initial residual = 4.48082e-12, Final residual = 3.09593e-17, No Iterations 1
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50709e-13  Max(alpha.water) = 1
GAMG:  Solving for p_rgh, Initial residual = 0.816907, Final residual = 4.85215e-08, No Iterations 115
time step continuity errors : sum local = 4.28464e-15, global = 5.24605e-17, cumulative = 1.00854e-08
PIMPLE: Iteration 2
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.86048e-07 0.170203)
    Orientation: (0.976646 -1.94577e-06 -0.214856 8.03987e-07 1 -5.40155e-06 0.214856 5.10266e-06 0.976646)
    Linear velocity: (-1.69724e+08 68563.8 1.26681e+10)
    Angular velocity: (5025.58 6.25083e+09 -33764.1)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000996818, No Iterations 7
time step continuity errors : sum local = 7.95507e-15, global = 2.225e-15, cumulative = 1.00854e-08
smoothSolver:  Solving for alpha.water, Initial residual = 1.49025e-11, Final residual = 5.13629e-17, No Iterations 1
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50708e-13  Max(alpha.water) = 1
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.50708e-13  Max(alpha.water) = 1
GAMG:  Solving for p_rgh, Initial residual = 0.811872, Final residual = 4.90291e-08, No Iterations 119
time step continuity errors : sum local = 1.41397e-14, global = -1.73317e-16, cumulative = 1.00854e-08
PIMPLE: Iteration 3
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.86048e-07 0.170203)
    Orientation: (0.976646 -1.94577e-06 -0.214856 8.03987e-07 1 -5.40155e-06 0.214856 5.10266e-06 0.976646)
    Linear velocity: (8.29019e+09 -3.34901e+06 -6.18776e+11)
    Angular velocity: (-227874 -2.8343e+11 1.53096e+06)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.00099994, No Iterations 69
time step continuity errors : sum local = 5.62137e-07, global = -5.51393e-09, cumulative = 4.57149e-09
smoothSolver:  Solving for alpha.water, Initial residual = 5.64671e-07, Final residual = 2.0355e-11, No Iterations 1
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.29505e-13  Max(alpha.water) = 1
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -6.84271e-11  Max(alpha.water) = 1
GAMG:  Solving for p_rgh, Initial residual = 0.999722, Final residual = 4.9867e-08, No Iterations 339
time step continuity errors : sum local = 2.80599e-11, global = -4.54197e-15, cumulative = 4.57148e-09
PIMPLE: Iteration 4
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.86051e-07 0.170202)
    Orientation: (0.976646 -1.94577e-06 -0.214856 8.03987e-07 1 -5.40155e-06 0.214856 5.10266e-06 0.976646)
    Linear velocity: (1.93416e+13 -7.8135e+09 -1.44365e+15)
    Angular velocity: (-4.69789e+08 -5.84324e+14 3.15625e+09)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.000965759, No Iterations 24
time step continuity errors : sum local = 5.71193e-07, global = -6.43151e-09, cumulative = -1.86003e-09
smoothSolver:  Solving for alpha.water, Initial residual = 5.37648e-07, Final residual = 2.25343e-13, No Iterations 1
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -7.49442e-13  Max(alpha.water) = 1.00001
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.53662  Min(alpha.water) = -2.57565e-09  Max(alpha.water) = 1.00001
GAMG:  Solving for p_rgh, Initial residual = 0.98905, Final residual = 4.95097e-08, No Iterations 252
time step continuity errors : sum local = 9.00697e-09, global = 1.07553e-10, cumulative = -1.75248e-09
smoothSolver:  Solving for omega, Initial residual = 5.76433e-06, Final residual = 2.7095e-08, No Iterations 1
smoothSolver:  Solving for k, Initial residual = 2.84242e-06, Final residual = 4.47391e-09, No Iterations 1
ExecutionTime = 5833.65 s  ClockTime = 5852 s

forces totForce write:
    sum of forces:
        pressure : (2.52119e+37 -4.06381e+35 7.84299e+37)
        viscous  : (-6.91233e+14 1.18334e+13 4.0348e+15)
        porous   : (0 0 0)
    sum of moments:
        pressure : (-4.5691e+34 -2.45591e+37 -3.64627e+35)
        viscous  : (-8.40372e+12 -2.20753e+15 -1.81932e+13)
        porous   : (0 0 0)

Courant Number mean: 0.000253879 max: 210.617
Interface Courant Number mean: 3.62046e-07 max: 0.0692958
deltaT = 3.43069e-24
Time = 0.861924

PIMPLE: Iteration 1
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (1.86617 -4.86051e-07 0.170202)
    Orientation: (0.976646 -1.94577e-06 -0.214856 8.03987e-07 1 -5.40155e-06 0.214856 5.10266e-06 0.976646)
    Linear velocity: (1.94766e+13 -7.86804e+09 -1.45373e+15)
    Angular velocity: (-4.73072e+08 -5.88408e+14 3.17831e+09)
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.00099803, No Iterations 165
time step continuity errors : sum local = 30.0278, global = -0.333752, cumulative = -0.333752
smoothSolver:  Solving for alpha.water, Initial residual = 0.000638797, Final residual = 9.17205e-11, No Iterations 190
Phase-1 volume fraction = 0.231448  Min(alpha.water) = -7.66302e-64  Max(alpha.water) = 1
Applying the previous iteration compression flux
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 673.428  Min(alpha.water) = -2.7388e+10  Max(alpha.water) = 6.97531e+10
GAMG:  Solving for p_rgh, Initial residual = 1, Final residual = 4.73618e-08, No Iterations 150
time step continuity errors : sum local = 0.00172731, global = -0.000342047, cumulative = -0.334094
PIMPLE: Iteration 2
forces forces:
    Not including porosity effects
Rigid-body motion of the hull
    Centre of rotation: (114902 -46.4165 -8.57608e+06)
    Orientation: (0.542239 -4.97447e-06 -0.840224 8.03987e-07 1 -5.40155e-06 0.840224 2.2534e-06 0.542239)
    Linear velocity: (6.69835e+28 -2.70595e+25 -4.99962e+30)
    Angular velocity: (-1.78355e+24 -2.21838e+30 1.19827e+25)
[0] #[57] #00    Foam::error::printStack(Foam::Ostream&)Foam::error::printStack(Foam::Ostream&)--------------------------------------------------------------------------
A process has executed an operation involving a call to the
"fork()" system call to create a child process.  Open MPI is currently
operating in a condition that could result in memory corruption or
other system errors; your job may hang, crash, or produce silent
data corruption.  The use of fork() (or system() or other calls that
create child processes) is strongly discouraged.

The process that invoked fork was:

  Local host:          [[59254,1],0] (PID 211497)

If you are *absolutely sure* that your application will successfully
and correctly survive a call to fork(), you may disable this warning
by setting the mpi_warn_on_fork MCA parameter to 0.
I only upload partial log.interFoam file, because I cannot upload it all here due to its large size, if possible please provide one email-address then I can send it to you.From the error, it can be seen that similar error occurs as the previous one, the force and moment on the hull became so big, also the linear & angular velocity and Courant number also became very big.

The case is simulated at calm water condition Fn=0.26 (U=1.59m/s) with free surface effect without any appendage. Please let me know if you need more information from me.Also, I am not sure about the difference between rigidBodyMotion and sixDoFRigidBodyMotion. Should I try sixDoFRigidBodyMotion for motionsolver? Finally, I need to compute the force on the hull, so I have to use the force function in controlDict.

Please let me know if you have any idea to resolve this error and many thanks for your generous help!

Best regards,
Tony
zyfsoton is offline   Reply With Quote

Old   March 15, 2022, 12:29
Default
  #18
Senior Member
 
Join Date: Sep 2013
Posts: 353
Rep Power: 20
Bloerb will become famous soon enough
You probably do not understand the basics, is my feeling. Sooo a rigid body is a point with inertia...Your centre of Mass is at roughly (1.86 0.00 0.17) in the coordinate system of your rigidBody...you now transform those into the coordinate system of your mesh. With a matrix (you rotate it):

1 0 0
0 1 0
0 0 1

I.e you do nothing to it. No rotation of that point. Now you shift everything once again by (1.86 0.00 0.17). Why? My feeling is that your centre of mass or transform is wrong. Why aren't you using this:

Code:
 transform       (1 0 0 0 1 0 0 0 1) (0 0 0);
Code:
centreOfMass    (1.86617 -4.84309e-07 0.170);
mass            229.5;
inertia         (-0.012 -3.05e-07 0.0056 -0.4300 -3.191e-08 -0.4351);
transform       (1 0 0 0 1 0 0 0 1) (1.8661 -4.8457e-07 0.17331);
Now you solve your rigid body...which is basically m*a=F_damper1+F_damper2+F_fluid
F_fluid is the force on the patch hull...which you now move by the solution of your rigid body....depending on your fvSolution setting once or several times per time step....3 or 4 outer correctors are still not enough...add atleast another 50. You need to converge both solutions. Before you enter the next time step. Once you can solve it you can lower them and increase the maxCo etc. But you need to have a good solution first. In other words look into your logfile...and make sure that within a time step...so within the outerCorrector loops...pimple 1 2 3 4 etc your linear velocity, angular velocity etc is not changing anymore...you need atleast so many outerCorrectors that those values become steady. Simply use grep or awk to search the logfile for those lines, and plot them. Make sure they are fully converged each time step. Simple example....
Code:
awk '{if ($0 ~ "Linear velocity") print}' logfile
awk goes through a file line by line. if the entire line ($0) contains "Linear velocity"...plot the entire line. You could also just plot $3,$4,$5 to only plot the columns containing the velocity. Or add a replacement like gsub(/\(/,""); in front...to remove the brackets.

You also need to keep in mind that you are using accelerationRelaxation 0.4; which means you are using a relaxation of the acceleration of your rigid body. Which is not realistic for transient runs. Hence you need to iterate this out using outer loops to get a realistic solution.

A rigid body has 6 degrees of freedom. You restrict those with a Pz and Ry joint. So your body can only move in z and rotate around y. Hence 2 degrees of freedom are left. rigidBody or sixDoF libaries are both possible. The rigidBody one you are using allows for multiple bodies. But you can't set initial velocites etc or use it as a boundary condition. Hence It does not really matter which one you use.




Since you are solving a rigid body with no internal forces...like a spring and are only adding dampers You ship shouldn't move at all if there is no force from the fluid. You are also not adding gravity. Since there is no line like g (0 0 9.81) in your dynamicMeshDict file...why are you not using gravity in your rigid body solution?....It is therefore highly likely that your fluid force on the patch hull is either wrong, or due to your center of mass being not correctly placed you are spinning your rigid body up.....


So please create a plot of position, velocity and angular velocity from your logfile (by using the awk line as shown) and post those here.
marxioxyz likes this.
Bloerb is offline   Reply With Quote

Old   March 17, 2022, 09:45
Default
  #19
Member
 
Tony Zhang
Join Date: Nov 2019
Location: soton
Posts: 45
Rep Power: 6
zyfsoton is on a distinguished road
Quote:
Originally Posted by Bloerb View Post
You probably do not understand the basics, is my feeling. Sooo a rigid body is a point with inertia...Your centre of Mass is at roughly (1.86 0.00 0.17) in the coordinate system of your rigidBody...you now transform those into the coordinate system of your mesh. With a matrix (you rotate it):

1 0 0
0 1 0
0 0 1

I.e you do nothing to it. No rotation of that point. Now you shift everything once again by (1.86 0.00 0.17). Why? My feeling is that your centre of mass or transform is wrong. Why aren't you using this:

Code:
 transform       (1 0 0 0 1 0 0 0 1) (0 0 0);
Code:
centreOfMass    (1.86617 -4.84309e-07 0.170);
mass            229.5;
inertia         (-0.012 -3.05e-07 0.0056 -0.4300 -3.191e-08 -0.4351);
transform       (1 0 0 0 1 0 0 0 1) (1.8661 -4.8457e-07 0.17331);
Now you solve your rigid body...which is basically m*a=F_damper1+F_damper2+F_fluid
F_fluid is the force on the patch hull...which you now move by the solution of your rigid body....depending on your fvSolution setting once or several times per time step....3 or 4 outer correctors are still not enough...add atleast another 50. You need to converge both solutions. Before you enter the next time step. Once you can solve it you can lower them and increase the maxCo etc. But you need to have a good solution first. In other words look into your logfile...and make sure that within a time step...so within the outerCorrector loops...pimple 1 2 3 4 etc your linear velocity, angular velocity etc is not changing anymore...you need atleast so many outerCorrectors that those values become steady. Simply use grep or awk to search the logfile for those lines, and plot them. Make sure they are fully converged each time step. Simple example....
Code:
awk '{if ($0 ~ "Linear velocity") print}' logfile
awk goes through a file line by line. if the entire line ($0) contains "Linear velocity"...plot the entire line. You could also just plot $3,$4,$5 to only plot the columns containing the velocity. Or add a replacement like gsub(/\(/,""); in front...to remove the brackets.

You also need to keep in mind that you are using accelerationRelaxation 0.4; which means you are using a relaxation of the acceleration of your rigid body. Which is not realistic for transient runs. Hence you need to iterate this out using outer loops to get a realistic solution.

A rigid body has 6 degrees of freedom. You restrict those with a Pz and Ry joint. So your body can only move in z and rotate around y. Hence 2 degrees of freedom are left. rigidBody or sixDoF libaries are both possible. The rigidBody one you are using allows for multiple bodies. But you can't set initial velocites etc or use it as a boundary condition. Hence It does not really matter which one you use.




Since you are solving a rigid body with no internal forces...like a spring and are only adding dampers You ship shouldn't move at all if there is no force from the fluid. You are also not adding gravity. Since there is no line like g (0 0 9.81) in your dynamicMeshDict file...why are you not using gravity in your rigid body solution?....It is therefore highly likely that your fluid force on the patch hull is either wrong, or due to your center of mass being not correctly placed you are spinning your rigid body up.....


So please create a plot of position, velocity and angular velocity from your logfile (by using the awk line as shown) and post those here.
Dear Bloerb,

Many thanks for your detailed explanation and patience with me...

I think now it works since I choose to use 6dof rigid body motion solver instead, as follows:

Code:
/*--------------------------------*- C++ -*----------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     | Website:  https://openfoam.org
    \\  /    A nd           | Version:  7
     \\/     M anipulation  |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      dynamicMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dynamicFvMesh   dynamicMotionSolverFvMesh;

motionSolverLibs    ("libsixDoFRigidBodyMotion.so");

motionSolver    sixDoFRigidBodyMotion;

patches         (hull);
innerDistance   0.19;
outerDistance   0.63;

centreOfMass    (1.86617058305 -4.84309040657e-07 0.170524016331);
mass            229.5;
momentOfInertia (6.20581170707 217.653741789 220.217667153);
rhoInf          1;
report          on;

value           uniform (0 0 0);

accelerationRelaxation 1;

solver
{
    type Newmark;
}

constraints
{
    zAxis
    {
        sixDoFRigidBodyMotionConstraint line;
        direction     (0 0 1);
    }
    yPlane
    {
        sixDoFRigidBodyMotionConstraint axis;
        axis          (0 1 0);
    }
}

restraints
{
    translationDamper
    {
        sixDoFRigidBodyMotionRestraint linearDamper;
        coeff         7217;
    }
    rotationDamper
    {
        sixDoFRigidBodyMotionRestraint sphericalAngularDamper;
        coeff         9727;
    }
}


// ************************************************************************* //
From the logfile, the linear and angular velocities are steady. The simulation still runs on HPC, I attached the first 15s ship resistance coefficients.

Since the hull is free to heave and pitch, I want to plot the heave and pitch RAO vs simulation time, do you have any idea how to do this? Should I set this in controlDict file or can by postprocessing after the simulation is done?

Many thanks for your help as always and best regards,
Tony
Attached Images
File Type: jpg Picture 2.jpg (65.0 KB, 38 views)
zyfsoton is offline   Reply With Quote

Old   March 17, 2022, 10:37
Default
  #20
Senior Member
 
Join Date: Sep 2013
Posts: 353
Rep Power: 20
Bloerb will become famous soon enough
You are using rhoInf=1. This is the density of air. It's been a long time since i used these solvers. It is likely that it isn't read and the density field of interFoam is used....RhoInf used to calculate the forces on the patch. Since many OpenFOAM solvers use a density normalized pressure p/rho. Hence to calculate the Force on the patch F=p*n you'd need to multiply it by a density. F=rhoInf*p*n. Please make sure this isn't the case. Otherwise there is basically 1000 times less force on the patch as realisitc (density of water=1000). From looking into the code it seems that it isn't used. It is basically a forces function object that calculates the forces on the patch which are given to the rigidbody solver. Hence removing that line altogether should work.


For plotting you have to keep in mind that you are only saving your result every so often. Hence if you need something for plotting over time make sure to use functionObjects to output the needed values. I am not familiar with RAOs what quantity do you need to calculate it?




Try this for starters in your controlDict:
Code:
functions 
{
      sixDoFRigidBodyState
    {
        type           sixDoFRigidBodyState;
        libs           ("libsixDoFRigidBodyState.so");
        angleFormat    degrees;
    }
}
Bloerb is offline   Reply With Quote

Reply


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Fatal error with rho in rigidBodyMotion! monnda OpenFOAM Running, Solving & CFD 1 November 3, 2023 09:39
Defining centre of mass/rotation in rigidBodyMotion solver tok3rat0r OpenFOAM Running, Solving & CFD 1 November 26, 2018 09:44
linearAxialAngularSpring restraint in rigidBodyMotion solver jiadongw OpenFOAM Running, Solving & CFD 4 August 13, 2018 21:33
Inertia in 6DOF pabsandoval FLUENT 0 March 13, 2013 18:06
Mass moment of inertia Cluain CFX 4 July 18, 2012 09:36


All times are GMT -4. The time now is 15:50.