magneticFoam.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  magneticFoam
26 
27 Description
28  Solver for the magnetic field generated by permanent magnets.
29 
30  A Poisson's equation for the magnetic scalar potential psi is solved
31  from which the magnetic field intensity H and magnetic flux density B
32  are obtained. The paramagnetic particle force field (H dot grad(H))
33  is optionally available.
34 
35 \*---------------------------------------------------------------------------*/
36 
37 #include "argList.H"
38 #include "magnet.H"
40 #include "simpleControl.H"
41 
42 #include "fvcGrad.H"
43 #include "fvcSnGrad.H"
44 #include "fvcReconstruct.H"
45 
46 #include "fvmLaplacian.H"
47 
48 using namespace Foam;
49 
50 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51 
52 int main(int argc, char *argv[])
53 {
55  (
56  "noH",
57  "do not write the magnetic field intensity field"
58  );
59 
61  (
62  "noB",
63  "do not write the magnetic flux density field"
64  );
65 
67  (
68  "HdotGradH",
69  "write the paramagnetic particle force field"
70  );
71 
72  #include "setRootCase.H"
73  #include "createTime.H"
74  #include "createMesh.H"
75 
76  simpleControl simple(mesh);
77 
78  #include "createFields.H"
79 
80  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
81 
82  Info<< "Calculating the magnetic field potential" << endl;
83 
84  runTime++;
85 
87  {
89  }
90 
91  psi.write();
92 
93  if (!args.optionFound("noH") || args.optionFound("HdotGradH"))
94  {
96  (
97  IOobject
98  (
99  "H",
100  runTime.name(),
101  mesh
102  ),
103  fvc::reconstruct(fvc::snGrad(psi)*mesh.magSf())
104  );
105 
106  if (!args.optionFound("noH"))
107  {
108  Info<< nl
109  << "Creating field H for time "
110  << runTime.name() << endl;
111 
112  H.write();
113  }
114 
115  if (args.optionFound("HdotGradH"))
116  {
117  Info<< nl
118  << "Creating field HdotGradH for time "
119  << runTime.name() << endl;
120 
121  volVectorField HdotGradH
122  (
123  IOobject
124  (
125  "HdotGradH",
126  runTime.name(),
127  mesh
128  ),
129  H & fvc::grad(H)
130  );
131 
132  HdotGradH.write();
133  }
134  }
135 
136  if (!args.optionFound("noB"))
137  {
138  Info<< nl
139  << "Creating field B for time "
140  << runTime.name() << endl;
141 
143  (
144  IOobject
145  (
146  "B",
147  runTime.name(),
148  mesh
149  ),
151  *fvc::reconstruct(murf*fvc::snGrad(psi)*mesh.magSf() + murf*Mrf)
152  );
153 
154  B.write();
155  }
156 
157  Info<< "\nEnd\n" << endl;
158 
159  return 0;
160 }
161 
162 
163 // ************************************************************************* //
static const Foam::dimensionedScalar B("B", Foam::dimless, 18.678)
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
Definition: argList.C:118
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:114
bool correctNonOrthogonal(const bool finalIter=false)
Non-orthogonal corrector loop.
virtual bool write(const bool write=true) const
Write using setting from DB.
Simple control class. Provides time-loop control methods which exit the simulation once convergence c...
Definition: simpleControl.H:67
simpleControl simple(mesh)
Calculate the gradient of the given field.
Reconstruct volField from a face flux field.
Calculate the snGrad of the given volField.
Calculate the matrix for the laplacian of the field.
const volScalarField & psi
surfaceScalarField Mrf(IOobject("Mrf", runTime.name(), mesh), mesh, dimensionedScalar(dimensionSet(0, 1, 0, 0, 0, 1, 0), 0))
surfaceScalarField murf(IOobject("murf", runTime.name(), mesh), mesh, 1)
int main(int argc, char *argv[])
Definition: magneticFoam.C:49
const dimensionedScalar mu0
Magnetic constant/permeability of free space: default SI units: [H/m].
tmp< VolField< typename outerProduct< vector, Type >::type > > reconstruct(const SurfaceField< Type > &ssf)
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< SurfaceField< Type > > snGrad(const VolField< Type > &vf, const word &name)
Definition: fvcSnGrad.C:45
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
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
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
Foam::argList args(argc, argv)