financialFoam.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 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 "fvCFD.H"
33 #include "writeCellGraph.H"
34 #include "OSspecific.H"
35 
36 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37 
38 int main(int argc, char *argv[])
39 {
40  #include "setRootCase.H"
41 
42  #include "createTime.H"
43  #include "createMesh.H"
44  #include "createFields.H"
45 
46  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47 
48  Info<< nl << "Calculating value(price of comodities)" << endl;
49 
50  surfaceScalarField phi("phi", (sigmaSqr - r)*(Pf & mesh.Sf()));
51 
52  volScalarField DV("DV", 0.5*sigmaSqr*sqr(P.component(Foam::vector::X)));
53 
54  Info<< "Starting time loop\n" << endl;
55 
56  while (runTime.loop())
57  {
59 
60  solve
61  (
62  fvm::ddt(V)
63  + fvm::div(phi, V)
64  - fvm::Sp(fvc::div(phi), V)
65  - fvm::laplacian(DV, V)
66  ==
67  - fvm::Sp(r, V)
68  );
69 
70  runTime.write();
71 
72  if (runTime.outputTime())
73  {
74  writeCellGraph(V, runTime.graphFormat());
75  writeCellGraph(delta, runTime.graphFormat());
76  }
77 
78  Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
79  << " ClockTime = " << runTime.elapsedClockTime() << " s"
80  << nl << endl;
81  }
82 
83  Info<< "End\n" << endl;
84 
85  return 0;
86 }
87 
88 
89 // ************************************************************************* //
rhoEqn solve()
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
tmp< GeometricField< Type, fvPatchField, volMesh > > Sp(const volScalarField &sp, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:67
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
messageStream Info
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
dynamicFvMesh & mesh
void writeCellGraph(const volScalarField &vsf, const word &graphFormat)
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
scalar delta
surfaceScalarField & phi
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
int main(int argc, char *argv[])
Definition: postCalc.C:54
dimensionedSymmTensor sqr(const dimensionedVector &dv)