CFDEMcoupling  2.4
 All Classes
cfdemCloud.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  cloud class managing DEM data for CFD-DEM coupling
32 
33 Class
34  Foam::cfdemCloud
35 
36 SourceFiles
37  cfdemCloud.C
38  cfdemCloudIO.C
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #ifndef cfdemCloud_H
43 #define cfdemCloud_H
44 
45 // choose version
46 #include "OFversion.H"
47 #include <vector>
48 
49 #include "fvCFD.H"
50 #include "IFstream.H"
51 
52 #if defined(version21) || defined(version16ext)
53  #include "turbulenceModel.H"
54 #elif defined(version15)
55  #include "RASModel.H"
56 #endif
57 
58 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
59 
60 namespace Foam
61 {
62 
63 // forward declarations
64 class forceModel;
65 class locateModel;
66 class voidFractionModel;
67 class dataExchangeModel;
68 class IOModel;
69 class probeModel;
70 class averagingModel;
71 class clockModel;
72 class smoothingModel;
73 class momCoupleModel;
74 class meshMotionModel;
75 class liggghtsCommandModel;
76 
77 /*---------------------------------------------------------------------------*\
78  Class cfdemCloud Declaration
79 \*---------------------------------------------------------------------------*/
80 
82 {
83 
84 // protected data
85 protected:
86  const fvMesh& mesh_;
87 
88  IOdictionary couplingProperties_;
89 
90  IOdictionary liggghtsCommandDict_;
91 
92  Switch solveFlow_;
93 
94  bool verbose_;
95 
96  bool ignore_;
97 
98  const word modelType_;
99 
100  mutable double **positions_;
101 
102  mutable double **velocities_;
103 
104  mutable double **fluidVel_;
105 
106  mutable double **fAcc_;
107 
108  mutable double **impForces_;
109 
110  mutable double **expForces_;
111 
112  mutable double **DEMForces_;
113 
114  mutable double **Cds_;
115 
116  mutable double **radii_;
117 
118  mutable double **voidfractions_;
119 
120  mutable double **cellIDs_;
121 
122  mutable double **particleWeights_;
123 
124  mutable double **particleVolumes_;
125 
126  mutable double **particleV_;
127 
128  int numberOfParticles_;
129 
130  scalar d32_;
131 
132  bool numberOfParticlesChanged_;
133 
134  mutable bool arraysReallocated_;
135 
136  const wordList forceModels_;
137 
138  const wordList momCoupleModels_;
139 
140  const wordList liggghtsCommandModelList_;
141 
142  const word turbulenceModelType_;
143 
144  mutable scalar cg_;
145 
146  bool cgOK_;
147 
148  bool impDEMdrag_;
149 
150  bool impDEMdragAcc_;
151 
152  mutable scalar imExSplitFactor_;
153 
154  mutable bool treatVoidCellsAsExplicitForce_; //will treat the coupling force in cells with no Us data explicitly
155 
156  bool useDDTvoidfraction_;
157 
158  mutable volScalarField ddtVoidfraction_;
159 
160  #if defined(version21) || defined(version16ext)
161  #ifdef compre
162  const compressible::turbulenceModel& turbulence_;
163  #else
164  const incompressible::turbulenceModel& turbulence_;
165  #endif
166  #elif defined(version15)
167  const incompressible::RASModel& turbulence_;
168  #endif
169 
170  autoPtr<forceModel>* forceModel_;
171 
172  autoPtr<locateModel> locateModel_;
173 
174  autoPtr<momCoupleModel>* momCoupleModel_;
175 
176  autoPtr<dataExchangeModel> dataExchangeModel_;
177 
178  autoPtr<IOModel> IOModel_;
179 
180  autoPtr<probeModel> probeModel_;
181 
182  autoPtr<voidFractionModel> voidFractionModel_;
183 
184  autoPtr<averagingModel> averagingModel_;
185 
186  autoPtr<clockModel> clockModel_;
187 
188  autoPtr<smoothingModel> smoothingModel_;
189 
190  autoPtr<meshMotionModel> meshMotionModel_;
191 
192  autoPtr<liggghtsCommandModel>* liggghtsCommand_;
193 
194 // Protected member functions
195  virtual void getDEMdata();
196 
197  virtual void giveDEMdata();
198 
199  virtual void setNumberOfParticles(int);
200 
201  virtual void findCells();
202 
203  virtual void setForces();
204 
205  virtual void setParticleForceField();
206 
207  virtual void setVectorAverages();
208 
209 public:
210 
211  friend class dataExchangeModel;
212  friend class voidFractionModel;
213  friend class forceModel;
214  friend class forceSubModel;
215 
216 // Constructors
217 
218  //- Construct from mesh and a list of particles
219  cfdemCloud
220  (
221  const fvMesh& mesh
222  );
223 
224  //- Destructor
225  virtual ~cfdemCloud();
226 
227 // public Member Functions
228 
229  // Access
230  void checkCG(bool);
231 
232  void setPos(double **&);
233 
234  word modelType(){ return modelType_; };
235 
236  label particleCell(int);
237 
238  vector position(int);
239 
240  vector velocity(int);
241 
242  vector fluidVel(int);
243 
244  virtual const forceModel& forceM(int);
245 
246  virtual int nrForceModels();
247 
248  scalar voidfraction(int);
249 
250  label liggghtsCommandModelIndex(word);
251 
252  inline void setCG(double) const;
253 
254  inline const scalar& cg() const;
255 
256  inline const bool& impDEMdrag() const;
257 
258  inline const bool& impDEMdragAcc() const;
259 
260  inline const scalar& imExSplitFactor() const;
261 
262  inline const bool& treatVoidCellsAsExplicitForce() const;
263 
264  inline const bool& ignore() const;
265 
266  inline const fvMesh& mesh() const;
267 
268  inline bool solveFlow() const;
269 
270  inline bool verbose() const;
271 
272  inline const IOdictionary& couplingProperties() const;
273 
274  inline double ** positions() const;
275 
276  inline double ** velocities() const;
277 
278  inline double ** fluidVels() const;
279 
280  inline double ** fAccs() const;
281 
282  inline double ** impForces() const;
283 
284  inline double ** expForces() const;
285 
286  inline double ** DEMForces() const;
287 
288  inline double ** Cds() const;
289 
290  inline double ** radii() const;
291 
292  inline double ** voidfractions() const;
293 
294  inline void get_radii(double**&) const;
295 
296  inline double ** cellIDs() const;
297 
298  inline void get_cellIDs(double**&) const;
299 
300  inline double ** particleWeights() const;
301 
302  virtual inline label body(int);
303 
304  virtual inline double particleVolume(int);
305 
306  inline scalar radius(int);
307 
308  virtual inline double d(int);
309 
310  inline scalar d32(bool recalc=true);
311  virtual inline double dMin() {return -1;}
312  virtual inline double dMax() {return -1;}
313  virtual inline int minType() {return -1;}
314  virtual inline int maxType() {return -1;}
315  virtual inline bool multipleTypesDMax() {return false;}
316  virtual inline bool multipleTypesDMin() {return false;}
317  virtual inline double ** particleDensity() const {return NULL;};
318  virtual inline int ** particleTypes() const {return NULL;};
319  virtual label particleType(label index) const {return -1;};
320 
321  //access to the particle's rotation and torque data
322  virtual inline double ** DEMTorques() const {return NULL;};
323  virtual inline double ** omegaArray() const {return NULL;};
324  virtual vector omega(int) const {return Foam::vector(0,0,0);};
325 
326  //access to the particles' orientation information
327  virtual inline double ** exArray() const {return NULL;};
328  virtual vector ex(int) const {return Foam::vector(0,0,0);};
329 
330  //Detector if SRF module is enable or not
331  virtual inline bool SRFOn(){return false;}
332 
333  inline int numberOfParticles() const;
334 
335  inline bool numberOfParticlesChanged() const;
336 
337  inline int numberOfClumps() const;
338 
339  inline bool arraysReallocated() const;
340 
341  inline const wordList& forceModels();
342 
343  inline const voidFractionModel& voidFractionM() const;
344 
345  inline const locateModel& locateM() const;
346 
347  inline const momCoupleModel& momCoupleM(int) const;
348 
349  inline const dataExchangeModel& dataExchangeM() const;
350 
351  inline const IOModel& IOM() const;
352 
353  inline const probeModel& probeM() const;
354 
355  inline const averagingModel& averagingM() const;
356 
357  inline const clockModel& clockM() const;
358 
359  inline const smoothingModel& smoothingM() const;
360 
361  inline const meshMotionModel& meshMotionM() const;
362 
363  inline const wordList& liggghtsCommandModelList() const;
364 
365  inline autoPtr<liggghtsCommandModel>* liggghtsCommand() const;
366 
367  #if defined(version21) || defined(version16ext)
368  #ifdef compre
369  inline const compressible::turbulenceModel& turbulence() const;
370  #else
371  inline const incompressible::turbulenceModel& turbulence() const;
372  #endif
373  #elif defined(version15)
374  inline const incompressible::RASModel& turbulence() const;
375  #endif
376 
377  // Write
378 
379  // write cfdemCloud internal data
380  virtual bool evolve(volScalarField&,volVectorField&,volVectorField&);
381 
382  virtual bool reAllocArrays() const;
383 
384  virtual bool reAllocArrays(int nP, bool forceRealloc) const; //force number of particles during reallocation
385 
386 
387  // IO
388  void writeScalarFieldToTerminal(double**&);
389 
390  void writeVectorFieldToTerminal(double**&);
391 
392  // functions
393  tmp<fvVectorMatrix> divVoidfractionTau(volVectorField& ,volScalarField&) const;
394 
395  tmp<volScalarField> ddtVoidfraction() const;
396 
397  void calcDdtVoidfraction(volScalarField& voidfraction) const;
398 
399  //tmp<fvVectorMatrix> ddtVoidfractionU(volVectorField& ,volScalarField&) const;
400 
401  tmp<volScalarField> voidfractionNuEff(volScalarField&) const;
402 
403  void resetArray(double**&,int,int,double resetVal=0.);
404 
405  std::vector< std::vector<double*> >* getVprobe();
406 
407  std::vector< std::vector<double> >* getSprobe();
408 
409 };
410 
411 
412 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
413 
414 } // End namespace Foam
415 
416 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
417 
418 #include "cfdemCloudI.H"
419 
420 #endif
421 
422 // ************************************************************************* //
Definition: averagingModel.H:54
Definition: meshMotionModel.H:53
Definition: dataExchangeModel.H:53
Definition: IOModel.H:54
Definition: momCoupleModel.H:54
Definition: clockModel.H:58
Definition: forceSubModel.H:54
Definition: smoothingModel.H:54
Definition: cfdemCloud.H:81
Definition: voidFractionModel.H:56
Definition: forceModel.H:55
Definition: probeModel.H:54
Definition: locateModel.H:53