MoDeNa  1.0
Software framework facilitating sequential multi-scale modelling
PUFoam.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2013 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 Application
25  MODENAFoam
26 
27 Description
28  Solver for 2 compressible, non-isothermal immiscible fluids using a VOF
29  (volume of fluid) phase-fraction based interface capturing approach.
30 
31  The momentum and other fluid properties are of the "mixture" and a single
32  momentum equation is solved.
33 
34  Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
35 
36  The solver has been adapted for modeling of Polyurethane foams (PU). In that,
37  it includes the Quadrature Method of Moments (QMOM) to solve a 'Population Balance Equation'
38  (PBE) determining the bubble size distribution inside PU foams.
39 
40  Moreover, the kinetics of the reactions including gelling, blowing and evaporation
41  of the physical blowing agent are incorporated into the solver.
42 
43  Finally, the kinetics and PBE has been coupled to describe the time evolution
44  of foaming process.
45 
46 \*---------------------------------------------------------------------------*/
57 #include "fvCFD.H"
59 #include "MULES.H"
60 #include "subCycle.H"
61 #include "rhoThermo.H"
62 #include "interfaceProperties.H"
63 #include "twoPhaseMixture.H"
64 #include "twoPhaseMixtureThermo.H"
65 #include "turbulenceModel.H"
66 #include "turbulentFluidThermoModel.H"
67 #include "pimpleControl.H"
68 #include "fixedFluxPressureFvPatchScalarField.H"
69 
70 extern "C"{void dsteqr_(char &, int *, double *, double *, double *, int *, double *, int *); }
71 // MoDeNa
72 #include "modena.h"
73 #include "modenaData.H"
74 // Kinetics headers
75 #include "KineticsFunctions.H"
76 // Moments headers
77 #include "PDA.H"
78 #include "growthRate.H"
79 #include "growthSource.H"
80 #include "coalescenceKernel.H"
81 #include "coalescenceSource.H"
82 
83 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
84 
85 int main(int argc, char *argv[])
86 {
87  #include "modenaCalls.H"
88  #include "setRootCase.H"
89  #include "createTime.H"
90  #include "createMesh.H"
91  #include "createTimeControls.H"
92 
93  #include "readGravitationalAcceleration.H"
94 
95  pimpleControl pimple(mesh);
96 
97  #include "readTimeControls.H"
98  #include "initContinuityErrs.H"
99  #include "createFields.H"
100  #include "CourantNo.H"
101  #include "setInitialDeltaT.H"
102  bool gellingPoint = false;
103  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
104 
105  Info<< "\nStarting time loop\n" << endl;
106  Info<< "initialFoamMass : " << initialFoamMass << endl;
107  while (runTime.run())
108  {
109  #include "readTimeControls.H"
110  #include "CourantNo.H"
111  #include "setDeltaT.H"
112 
113  runTime++;
114 
115  Info<< "Time = " << runTime.timeName() << nl << endl;
116 
117  // --- Pressure-velocity PIMPLE corrector loop
118  while (pimple.loop())
119  {
120  if (simulationTarget == "mold-filling")
121  {
122  if
123  (
124  alpha1.weightedAverage(mesh.V()).value() > 0.01
125  && gellingPoint == false
126  )
127  {
128  #include "alphaEqnsSubCycle.H"
129  // correct interface on first PIMPLE corrector
130  if (pimple.corr() == 1)
131  {
132  interface.correct();
133  }
134  #include "alphaCorrection.H"
135  #include "checkGellingPoint.H"
136  solve(fvm::ddt(rho) + fvc::div(rhoPhi));
137 
138 
139  #include "conversionSources.H"
140  #include "conversionEqns.H"
141  // #include "conversionCheck.H"
142 
143  #include "rheologyModel.H"
144 
145  #include "UEqn.H"
146 
147  #include "TSSource.H"
148  #include "TSEqn.H"
149 
150  #include "TCheck.H"
151  #include "densityEqns.H"
152 
153  #include "MomSources.H"
154  #include "MomEqns.H"
155  #include "MomConvert.H"
156 
157  // --- Pressure corrector loop
158  while (pimple.correct())
159  {
160  #include "pEqn.H"
161  }
162  #include "BASources.H"
163  #include "BAEqns.H"
164  #include "BACheck.H"
165 
166  if (pimple.turbCorr())
167  {
168  turbulence->correct();
169  }
170  }
171  }
172  else
173  {
174  if (gellingPoint == false)
175  {
176  #include "alphaEqnsSubCycle.H"
177 
178  // correct interface on first PIMPLE corrector
179  if (pimple.corr() == 1)
180  {
181  interface.correct();
182  }
183  #include "alphaCorrection.H"
184 
185 
186  #include "checkGellingPoint.H"
187  solve(fvm::ddt(rho) + fvc::div(rhoPhi));
188  #include "conversionSources.H"
189  #include "conversionEqns.H"
190  // #include "conversionCheck.H"
191 
192  #include "rheologyModel.H"
193 
194  #include "UEqn.H"
195 
196  #include "TSSource.H"
197  #include "TSEqn.H"
198 
199  #include "TCheck.H"
200  #include "densityEqns.H"
201 
202  #include "MomSources.H"
203  #include "MomEqns.H"
204  #include "MomConvert.H"
205 
206  // --- Pressure corrector loop
207  while (pimple.correct())
208  {
209  #include "pEqn.H"
210  }
211  #include "BASources.H"
212  #include "BAEqns.H"
213  #include "BACheck.H"
214 
215  if (pimple.turbCorr())
216  {
217  turbulence->correct();
218  }
219  }
220  }
221 
222  #include "continuityErrors.H"
223 
224  Info<< "\nMass of foam: "
225  << fvc::domainIntegrate(alpha2*rho_foam)
226  << endl;
227  }
228 
229  runTime.write();
230 
231  Info<< "ExecutionTime = "
232  << runTime.elapsedCpuTime()
233  << " s\n\n" << endl;
234  }
235  Info<< "finalFoamMass : "
236  << fvc::domainIntegrate(rho_foam*alpha2)
237  << endl;
238  Info<< "End\n" << endl;
239 
240  return 0;
241 }
242 // ************************************************************************* //
243 ///@endcond
builds the RHS of PDEs for the moments.
allocate memory for the surrogate models
momentum equation
converts moments based on the unit volume of foam
three different functions for growth rate
builds and solves moments equations
builds and solves conversion equations based on the kinetic model chosen.
functions for kinetics calculations
checks the boundedness of field variables related to the physical and chemical blowing agents...
source term due to bubble growth
checks the boundedness of alpha (phase volume fraction)
builds the RHS of conversion equations based on the kinetics model chosen.
source term due to bubbles coalescence
builds and solves the temperature equation
implementatino of three rheology models including constant, Newtonian, and non-Newtonian ...
creates all the required field variables.
calculates and reports the continuity error
builds the RHS of PDEs for the blowing agents.
calculates the PU foam density
int main(int argc, char *argv[])
Reads parameters. Creates struts and walls. Saves foam morphology to a file.
Definition: foams.cc:23
constant coalescence kernel
builds the RHS of temperature equation
checks if the gelling point reaches.
checks the boundedness of temperature values
instantiates the surrogate models
builds the PDEs for the presence of different blowing agents. It also considers the user choice for t...
builds and solves pressure equation