CFDEMcoupling  2.4
 All Classes
forceModel.H
1 /*---------------------------------------------------------------------------*\
2  CFDEMcoupling - Open Source CFD-DEM coupling
3 
4  CFDEMcoupling is part of the CFDEMproject
5  www.cfdem.com
6  Christoph Goniva, christoph.goniva@cfdem.com
7  Copyright 2009-2012 JKU Linz
8  Copyright 2012- DCS Computing GmbH, Linz
9 -------------------------------------------------------------------------------
10 License
11  This file is part of CFDEMcoupling.
12 
13  CFDEMcoupling is free software; you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by the
15  Free Software Foundation; either version 3 of the License, or (at your
16  option) any later version.
17 
18  CFDEMcoupling is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with CFDEMcoupling; if not, write to the Free Software Foundation,
25  Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
26 
27 Description
28  This code is designed to realize coupled CFD-DEM simulations using LIGGGHTS
29  and OpenFOAM(R). Note: this code is not part of OpenFOAM(R) (see DISCLAIMER).
30 
31 Class
32  forceModel
33 
34 SourceFiles
35  forceModel.C
36 
37 \*---------------------------------------------------------------------------*/
38 
39 #ifndef forceModel_H
40 #define forceModel_H
41 
42 #include "fvCFD.H"
43 #include "cfdemCloud.H"
44 #include "probeModel.H"
45 #include "forceSubModel.H"
46 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47 
48 namespace Foam
49 {
50 
51 /*---------------------------------------------------------------------------*\
52  Class forceModel Declaration
53 \*---------------------------------------------------------------------------*/
54 
56 {
57 
58 protected:
59 
60  // Protected data
61  const dictionary& dict_;
62 
63  cfdemCloud& particleCloud_;
64 
65  //Switch treatExplicit_; // marker to treat force in implicit way (otherwise explicit)
66  //Switch treatDEM_; // marker to use the force only on DEM side
67  //Switch implDEM_; // marker to use the implicit force on DEM side
68 
69  mutable volVectorField impParticleForces_; // sum of implicit particle forces [N]
70 
71  mutable volVectorField expParticleForces_; // sum of explicit particle forces [N]
72 
73  bool coupleForce_;
74 
75  const word modelType_;
76 
77  bool probeIt_;
78 
79  bool requiresEx_; //requires a orientation vector
80 
81  wordList forceSubModels_;
82 
83  autoPtr<forceSubModel>* forceSubModel_;
84 
85 public:
86 
87  //- Runtime type information
88  TypeName("forceModel");
89 
90  // Declare runtime constructor selection table
91 
92  declareRunTimeSelectionTable
93  (
94  autoPtr,
95  forceModel,
96  dictionary,
97  (
98  const dictionary& dict,
99  cfdemCloud& sm
100  ),
101  (dict,sm)
102  );
103 
104 
105  // Constructors
106 
107  //- Construct from components
108  forceModel
109  (
110  const dictionary& dict,
111  cfdemCloud& sm
112  );
113 
114 
115  // Destructor
116 
117  virtual ~forceModel();
118 
119 
120  // Selector
121 
122  static autoPtr<forceModel> New
123  (
124  const dictionary& dict,
125  cfdemCloud& sm,
126  word forceType
127  );
128 
129 
130  // Member Functions
131  virtual void setForce() const = 0;
132 
133  //tmp<volScalarField> provideScalarField();
134 
135  virtual void manipulateScalarField(volScalarField&) const;
136 
137  // Access
138  word modelType(){ return modelType_; };
139 
140  inline volVectorField& impParticleForces() const { return impParticleForces_;};
141 
142  inline volVectorField& expParticleForces() const { return expParticleForces_;};
143 
144  inline double ** impForces() const { return particleCloud_.impForces_;};
145 
146  inline double ** expForces() const { return particleCloud_.expForces_;};
147 
148  inline double ** DEMForces() const { return particleCloud_.DEMForces_;};
149 
150  inline double ** Cds() const { return particleCloud_.Cds_;};
151 
152  inline double ** fluidVel() const { return particleCloud_.fluidVel_;};
153 
154  inline const bool& coupleForce() const { return coupleForce_;};
155 
156  virtual inline bool& requiresEx() { return requiresEx_;};
157 
158  void repartitionImExForces() const; //Repartition Implixit/Explicit forces
159 
160  void treatVoidCells() const;
161 
162  inline const wordList& forceSubModels(){ return forceSubModels_; };
163 
164  inline const forceSubModel& forceSubM(int i) const { return forceSubModel_[i]; };
165 
166  inline int nrForceSubModels(){ return forceSubModels_.size(); };
167 
168  void setForceSubModels(dictionary& dict);
169 };
170 
171 
172 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
173 
174 } // End namespace Foam
175 
176 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
177 
178 #endif
179 
180 // ************************************************************************* //
Definition: forceSubModel.H:54
Definition: cfdemCloud.H:81
Definition: forceModel.H:55