financialFoam.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  financialFoam
26 
27 Description
28  Solves the Black-Scholes equation to price commodities.
29 
30 \*---------------------------------------------------------------------------*/
31 
32 #include "argList.H"
33 #include "timeSelector.H"
34 #include "setWriter.H"
35 #include "writeFile.H"
36 
37 #include "fvcGrad.H"
38 
39 #include "fvmDdt.H"
40 #include "fvmDiv.H"
41 #include "fvmLaplacian.H"
42 
43 using namespace Foam;
44 
45 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
46 
47 int main(int argc, char *argv[])
48 {
49  #define NO_CONTROL
50  #include "postProcess.H"
51 
52  #include "setRootCase.H"
53  #include "createTime.H"
54  #include "createMesh.H"
55  #include "createFields.H"
56 
57  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
58 
59  Info<< nl << "Calculating value(price of commodities)" << endl;
60 
61  surfaceScalarField phi("phi", (sigmaSqr - r)*(Pf & mesh.Sf()));
62 
63  volScalarField DV("DV", 0.5*sigmaSqr*sqr(P.component(Foam::vector::X)));
64 
65  Info<< "Starting time loop\n" << endl;
66 
67  while (runTime.loop())
68  {
70 
71  solve
72  (
73  fvm::ddt(V)
74  + fvm::div(phi, V)
75  - fvm::Sp(fvc::div(phi), V)
76  - fvm::laplacian(DV, V)
77  ==
78  - fvm::Sp(r, V)
79  );
80 
81  runTime.write();
82 
83  if (runTime.writeTime())
84  {
85  setWriter::New(runTime.controlDict().lookup("graphFormat"))->write
86  (
87  runTime.globalPath()
89  /args.executable()
90  /runTime.name(),
91 
92  args.executable(),
93 
94  coordSet(true, word::null, mesh.C().primitiveField(), "x"),
95 
96  "V", V.primitiveField(),
97  "delta", delta.primitiveField()
98  );
99  }
100 
101  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
102  << " ClockTime = " << runTime.elapsedClockTime() << " s"
103  << nl << endl;
104  }
105 
106  Info<< "End\n" << endl;
107 
108  return 0;
109 }
110 
111 
112 // ************************************************************************* //
scalar delta
Generic GeometricField class.
const word & executable() const
Name of executable without the path.
Definition: argListI.H:36
Holds list of sampling positions.
Definition: coordSet.H:51
static const word outputPrefix
Directory prefix.
Definition: writeFile.H:72
static autoPtr< setWriter > New(const word &writeType, const IOstream::streamFormat writeFormat=IOstream::ASCII, const IOstream::compressionType writeCompression=IOstream::UNCOMPRESSED)
Select given write options.
Definition: setWriter.C:286
static const word null
An empty word.
Definition: word.H:77
int main(int argc, char *argv[])
Definition: financialFoam.C:44
Calculate the gradient of the given field.
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for the laplacian of the field.
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
tmp< VolField< Type > > div(const SurfaceField< Type > &ssf)
Definition: fvcDiv.C:47
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const VolField< Type > &vf, const word &name)
Definition: fvmDiv.C:48
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
Definition: fvmDdt.C:46
Namespace for OpenFOAM.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
messageStream Info
static const char nl
Definition: Ostream.H:266
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
Execute application functionObjects to post-process existing results.
Foam::argList args(argc, argv)