federicabi |
September 14, 2016 08:44 |
Customized restraint
Hi!
I'm pretty new with openFOAM and with C++, and I'm having troubles with a customized restraint that I modified from linearSpring.
These are my .C and .H files
Code:
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "variableLinearSpring.H"
#include "addToRunTimeSelectionTable.H"
#include "sixDoFRigidBodyMotion.H"
#include "forces.H"
#include "IOobject.H"
#include "dictionary.H"
#include "Time.H"
#include "IOmanip.H"
#include "Pstream.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
namespace sixDoFRigidBodyMotionRestraints
{
defineTypeNameAndDebug(variableLinearSpring, 0);
addToRunTimeSelectionTable
(
sixDoFRigidBodyMotionRestraint,
variableLinearSpring,
dictionary
);
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::sixDoFRigidBodyMotionRestraints::variableLinearSpring::variableLinearSpring
(
const word& name,
const dictionary& sDoFRBMRDict
)
:
sixDoFRigidBodyMotionRestraint(name, sDoFRBMRDict),
refAttachmentPt_(),
alfaini_(),
forzaxini_()
{
read(sDoFRBMRDict);
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::sixDoFRigidBodyMotionRestraints::variableLinearSpring::~variableLinearSpring()
{}
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
void Foam::sixDoFRigidBodyMotionRestraints::variableLinearSpring::restrain
(
const sixDoFRigidBodyMotion& motion,
vector& restraintPosition,
vector& restraintForce,
vector& restraintMoment
) const
{
static bool firstIt_ = true;
scalar alfa;
scalar forzax;
scalar forzaxPrevIter;
scalar alfaPrevIter;
if (firstIt_)
{
alfa = (alfaini_ * 3.1415 / 180);
forzax = forzaxini_;
}
else
{
FILE *forze;
scalar timestep;
scalar fxp;
scalar fyp;
scalar fzp;
scalar fxv;
scalar fyv;
scalar fzv;
scalar fxpo;
scalar fypo;
scalar fzpo;
scalar mxp;
scalar myp;
scalar mzp;
scalar mxv;
scalar myv;
scalar mzv;
scalar mxpo;
scalar mypo;
scalar mzpo;
forze = fopen("postProcessing/forces/0.3/forces.dat", "r");
scanf("%*[^\n]"); /* Skip to the End of the Line */
for (int i = 1; i <= 3; i++)
{
scanf("%*1[\n]"); /* Skip Three Newline */
}
fscanf(forze, "%le ((%le %le %le) (%le %le %le) (%le %le %le)) ((%le %le %le) (%le %le %le) (%le %le %le))", ×tep, &fxp, &fyp, &fzp, &fxv, &fyv, &fzv, &fxpo, &fypo, &fzpo, &mxp, &myp, &mzp, &mxv, &myv, &mzv, &mxpo, &mypo, &mzpo);
fclose(forze);
forzaxPrevIter = fxp + fxv + fxpo;
alfaPrevIter = (alfaini_ * 3.1415 / 180) + asin(motion.orientation()[3]);
alfa = alfaPrevIter;
forzax = forzaxPrevIter;
}
restraintPosition = motion.transform(refAttachmentPt_);
vector v = motion.velocity(restraintPosition);
restraintForce = vector(forzax,0,forzax*tan(alfa));
restraintMoment = vector::zero;
if (motion.report())
{
Info<< " attachmentPt" << restraintPosition
<< " force " << restraintForce
<< endl;
}
firstIt_ = false;
}
bool Foam::sixDoFRigidBodyMotionRestraints::variableLinearSpring::read
(
const dictionary& sDoFRBMRDict
)
{
sixDoFRigidBodyMotionRestraint::read(sDoFRBMRDict);
sDoFRBMRCoeffs_.lookup("refAttachmentPt") >> refAttachmentPt_;
sDoFRBMRCoeffs_.lookup("alfaini") >> alfaini_;
sDoFRBMRCoeffs_.lookup("forzaxini") >> forzaxini_;
return true;
}
void Foam::sixDoFRigidBodyMotionRestraints::variableLinearSpring::write
(
Ostream& os
) const
{
os.writeKeyword("refAttachmentPt")
<< refAttachmentPt_ << token::END_STATEMENT << nl;
os.writeKeyword("alfaini")
<< alfaini_ << token::END_STATEMENT << nl;
os.writeKeyword("forzaxini")
<< forzaxini_ << token::END_STATEMENT << nl;
}
// ************************************************************************* //
Code:
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::sixDoFRigidBodyMotionRestraints::variableLinearSpring
Description
sixDoFRigidBodyMotionRestraints model. Variable Linear spring.
SourceFiles
variableLinearSpring.C
\*---------------------------------------------------------------------------*/
#ifndef variableLinearSpring_H
#define variableLinearSpring_H
#include "sixDoFRigidBodyMotionRestraint.H"
#include "point.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace sixDoFRigidBodyMotionRestraints
{
/*---------------------------------------------------------------------------*\
Class variableLinearSpring Declaration
\*---------------------------------------------------------------------------*/
class variableLinearSpring
:
public sixDoFRigidBodyMotionRestraint
{
// Private data
//- Punto di applicazione del traino (punto elica)
point refAttachmentPt_;
//- Angolo di inclinazione della linea d'assi
scalar alfaini_;
//- Forza di resistenza iniziale
scalar forzaxini_;
public:
//- Runtime type information
TypeName("variableLinearSpring");
// Constructors
//- Construct from components
variableLinearSpring
(
const word& name,
const dictionary& sDoFRBMRDict
);
//- Construct and return a clone
virtual autoPtr<sixDoFRigidBodyMotionRestraint> clone() const
{
return autoPtr<sixDoFRigidBodyMotionRestraint>
(
new variableLinearSpring(*this)
);
}
//- Destructor
virtual ~variableLinearSpring();
// Member Functions
//- Calculate the restraint position, force and moment.
// Global reference frame vectors.
virtual void restrain
(
const sixDoFRigidBodyMotion& motion,
vector& restraintPosition,
vector& restraintForce,
vector& restraintMoment
) const;
//- Update properties from given dictionary
virtual bool read(const dictionary& sDoFRBMRCoeff);
//- Write
virtual void write(Ostream&) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace solidBodyMotionFunctions
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
I compiled it without errors and I applied it to a simulation with interDyMFoam solver.
I run the code and it runs for the first iteration and then it stops without any error warning, just it remains like below for ever and I don't know why.
Code:
/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 3.0.1 |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : 3.0.1-d8a290b55d28
Exec : interDyMFoam
Date : Sep 14 2016
Time : 14:17:37
Host : "rd-cfd"
PID : 45770
Case : /home/jolly/OpenFOAM/Stage_Federica/AZ95rph/AZ95rph_snappyHexMesh/AZ95rph_snappyHexMesh_2DoFsenzaForze_26kn
nProcs : 1
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster
allowSystemOperations : Allowing user-supplied system call operations
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time
Create mesh for time = 0.3
Selecting dynamicFvMesh dynamicMotionSolverFvMesh
Selecting motion solver: sixDoFRigidBodyMotion
Selecting sixDoFSolver Newmark
Translational constraint tensor (0 0 0 0 0 0 0 0 1)
Rotational constraint tensor (0 0 0 0 1 0 0 0 0)
PIMPLE: Operating solver in PISO mode
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
kOmegaSSTCoeffs
{
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.000767548, No Iterations 9
GAMG: Solving for pcorr, Initial residual = 0.000969538, Final residual = 8.08389e-07, No Iterations 5
time step continuity errors : sum local = 4.04061e-11, global = -3.40519e-12, cumulative = -3.40519e-12
Reading/calculating face velocity Uf
Courant Number mean: 0.00117436 max: 2.90517
Starting time loop
forces forces:
Not including porosity effects
Courant Number mean: 0.00117436 max: 2.90517
Interface Courant Number mean: 1.75985e-06 max: 1.32479
deltaT = 0.000178571
Time = 0.300179
PIMPLE: iteration 1
Restraint externalForces: attachmentPt(1.58696 1.05 -1.73213) force (-55645 0 -2902.51)
6-DoF rigid body motion
Centre of rotation: (9.3 0 3.18594)
Centre of mass: (9.3 0 3.18594)
Orientation: (1 0 -0.000808245 0 1 0 0.000808245 0 1)
Linear velocity: (0 0 0.152943)
Angular velocity: (0 0.0157576 0)
Execution time for mesh.update() = 9.18 s
GAMG: Solving for pcorr, Initial residual = 1, Final residual = 0.000991738, No Iterations 11
GAMG: Solving for pcorr, Initial residual = 0.000365091, Final residual = 7.69852e-07, No Iterations 4
time step continuity errors : sum local = 9.41351e-12, global = -6.12812e-13, cumulative = -4.01801e-12
smoothSolver: Solving for alpha.water, Initial residual = 4.79871e-08, Final residual = 1.22406e-16, No Iterations 3
Phase-1 volume fraction = 0.545033 Min(alpha.water) = -2.91408e-11 Max(alpha.water) = 1.00003
MULES: Correcting alpha.water
MULES: Correcting alpha.water
MULES: Correcting alpha.water
Phase-1 volume fraction = 0.545033 Min(alpha.water) = -3.92276e-10 Max(alpha.water) = 1.00003
GAMG: Solving for p_rgh, Initial residual = 0.00530249, Final residual = 4.33916e-06, No Iterations 10
GAMG: Solving for p_rgh, Initial residual = 0.000367444, Final residual = 3.22317e-07, No Iterations 5
time step continuity errors : sum local = 4.18558e-12, global = -2.79785e-13, cumulative = -4.29779e-12
GAMG: Solving for p_rgh, Initial residual = 9.06052e-05, Final residual = 5.46971e-08, No Iterations 6
GAMG: Solving for p_rgh, Initial residual = 4.83781e-05, Final residual = 7.13157e-09, No Iterations 7
time step continuity errors : sum local = 9.60132e-14, global = -4.5178e-16, cumulative = -4.29824e-12
smoothSolver: Solving for omega, Initial residual = 1.11214e-05, Final residual = 5.23818e-12, No Iterations 1
smoothSolver: Solving for k, Initial residual = 6.08195e-08, Final residual = 5.60729e-12, No Iterations 1
ExecutionTime = 219.07 s ClockTime = 222 s
forces forces output:
sum of forces:
pressure : (-1.49126e+06 -1.20981e+07 2.89098e+07)
viscous : (-16593.9 2035.3 -26.9414)
porous : (0 0 0)
sum of moments:
pressure : (-9.04111e+06 3.01946e+07 89945)
viscous : (5158.34 53168.4 42987.5)
porous : (0 0 0)
fieldMinMax minmaxdomain output:
min(p) = -279303 at location (0.139124 0.399286 -1.14294)
max(p) = 787024 at location (5.64798 0.333083 -1.02371)
min(U) = (-41.2283 -18.1915 -41.8936) at location (18.8947 2.47892 1.3215)
max(U) = (5.22655 7.86842 8.85684) at location (21.4319 0.100175 0.228842)
Courant Number mean: 0.000873884 max: 2.21732
Interface Courant Number mean: 1.31434e-06 max: 0.781705
deltaT = 0.000200437
Time = 0.300379
PIMPLE: iteration 1
Can anybody help me to solve my problem?
Thanks in advance
UPDATE: I tried to remove the if cycle and the firstIt dependancy and the simulation get stuck at the first iteration with no errors. I assume that the problem is about the reading of the forces.dat file but I have no idea about how to solve this problem.
Federica
|