laplacianFoam.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) 2011-2023 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  laplacianFoam
26 
27 Description
28  Solves a simple Laplace equation, e.g. for thermal diffusion in a solid.
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #include "argList.H"
33 #include "fvModels.H"
34 #include "fvConstraints.H"
35 #include "simpleControl.H"
36 
37 #include "fvmLaplacian.H"
38 
39 using namespace Foam;
40 
41 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42 
43 int main(int argc, char *argv[])
44 {
45  #include "setRootCase.H"
46 
47  #include "createTime.H"
48  #include "createMesh.H"
49 
50  simpleControl simple(mesh);
51 
52  #include "createFields.H"
53 
54  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55 
56  Info<< "\nCalculating temperature distribution\n" << endl;
57 
58  while (simple.loop(runTime))
59  {
60  Info<< "Time = " << runTime.userTimeName() << nl << endl;
61 
62  fvModels.correct();
63 
65  {
66  fvScalarMatrix TEqn
67  (
68  fvm::ddt(T) - fvm::laplacian(DT, T)
69  ==
71  );
72 
74  TEqn.solve();
76  }
77 
78  #include "write.H"
79 
80  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
81  << " ClockTime = " << runTime.elapsedClockTime() << " s"
82  << nl << endl;
83  }
84 
85  Info<< "End\n" << endl;
86 
87  return 0;
88 }
89 
90 
91 // ************************************************************************* //
Finite volume constraints.
Definition: fvConstraints.H:67
bool constrain(fvMatrix< Type > &eqn) const
Apply constraints to an equation.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
Finite volume models.
Definition: fvModels.H:65
tmp< fvMatrix< Type > > source(const VolField< Type > &field) const
Return source for an equation.
virtual void correct()
Correct the fvModels.
Definition: fvModels.C:380
bool correctNonOrthogonal(const bool finalIter=false)
Non-orthogonal corrector loop.
Simple control class. Provides time-loop control methods which exit the simulation once convergence c...
Definition: simpleControl.H:67
bool loop(Time &time)
Time loop loop.
Definition: simpleControl.C:80
simpleControl simple(mesh)
Calculate the matrix for the laplacian of the field.
int main(int argc, char *argv[])
Definition: laplacianFoam.C:40
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
Definition: fvmDdt.C:46
Namespace for OpenFOAM.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
messageStream Info
static const char nl
Definition: Ostream.H:266
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)