displacementLinearMotionMotionSolver.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2018 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33  defineTypeNameAndDebug(displacementLinearMotionMotionSolver, 0);
34 
36  (
37  motionSolver,
38  displacementLinearMotionMotionSolver,
39  dictionary
40  );
41 }
42 
43 
44 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45 
46 Foam::displacementLinearMotionMotionSolver::
47 displacementLinearMotionMotionSolver
48 (
49  const polyMesh& mesh,
50  const IOdictionary& dict
51 )
52 :
53  points0MotionSolver(mesh, dict, typeName),
54  axis_(normalised(vector(coeffDict().lookup("axis")))),
55  xFixed_(readScalar(coeffDict().lookup("xFixed"))),
56  xMoving_(readScalar(coeffDict().lookup("xMoving"))),
57  displacement_(Function1<scalar>::New("displacement", coeffDict()))
58 {}
59 
60 
61 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
62 
65 {}
66 
67 
68 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
69 
72 {
73  tmp<pointField> tcurPoints(new pointField(points0()));
74  pointField& curPoints = tcurPoints.ref();
75 
76  const scalar t = time().value();
77 
78  const scalar displacement = displacement_->value(t);
79 
80  forAll(curPoints, i)
81  {
82  const scalar lambda =
83  (xFixed_ - (axis_ & curPoints[i]))/(xFixed_ - xMoving_);
84 
85  if (lambda > 1)
86  {
87  curPoints[i] += axis_*displacement;
88  }
89  else if (lambda > 0)
90  {
91  curPoints[i] += axis_*lambda*displacement;
92  }
93  }
94 
95  return tcurPoints;
96 }
97 
98 
99 // ************************************************************************* //
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
Definition: VectorSpaceI.H:413
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
Macros for easy insertion into run-time selection tables.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
stressControl lookup("compactNormalStress") >> compactNormalStress
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
Virtual base class for displacement motion solvers.
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.