createFields.H
Go to the documentation of this file.
1  Info<< "Reading field psi\n" << endl;
3  (
4  IOobject
5  (
6  "psi",
7  runTime.name(),
8  mesh,
9  IOobject::MUST_READ,
10  IOobject::AUTO_WRITE
11  ),
12  mesh
13  );
14 
15  Info<< "Reading physicalProperties\n" << endl;
16 
17  IOdictionary physicalProperties
18  (
19  IOobject
20  (
21  "physicalProperties",
22  runTime.constant(),
23  mesh,
24  IOobject::MUST_READ_IF_MODIFIED,
25  IOobject::NO_WRITE
26  )
27  );
28 
29  List<magnet> magnets(physicalProperties.lookup("magnets"));
30 
32  (
33  IOobject
34  (
35  "murf",
36  runTime.name(),
37  mesh
38  ),
39  mesh,
40  1
41  );
42 
44  (
45  IOobject
46  (
47  "Mrf",
48  runTime.name(),
49  mesh
50  ),
51  mesh,
52  dimensionedScalar(dimensionSet(0, 1, 0, 0, 0, 1, 0), 0)
53  );
54 
56  {
57  label magnetZonei = mesh.faceZones().findIndex(magnets[i].name());
58 
59  if (magnetZonei == -1)
60  {
61  FatalIOErrorIn(args.executable().c_str(), physicalProperties)
62  << "Cannot find faceZone for magnet " << magnets[i].name()
63  << exit(FatalIOError);
64  }
65 
66  const labelList& faces = mesh.faceZones()[magnetZonei];
67 
68  const scalar muri = magnets[i].mur();
69  const scalar Mri = magnets[i].Mr().value();
70  const vector& orientationi = magnets[i].orientation();
71 
72  const surfaceVectorField& Sf = mesh.Sf();
73 
74  forAll(faces, i)
75  {
76  label facei = faces[i];
77  murf[facei] = muri;
78  Mrf[facei] = Mri*(orientationi & Sf[facei]);
79  }
80  }
const word & executable() const
Name of executable without the path.
Definition: argListI.H:36
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:340
const volScalarField & psi
surfaceScalarField Mrf(IOobject("Mrf", runTime.name(), mesh), mesh, dimensionedScalar(dimensionSet(0, 1, 0, 0, 0, 1, 0), 0))
Info<< "Reading field psi\n"<< endl;volScalarField psi(IOobject("psi", runTime.name(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading physicalProperties\n"<< endl;IOdictionary physicalProperties(IOobject("physicalProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));List< magnet > magnets(physicalProperties.lookup("magnets"))
surfaceScalarField murf(IOobject("murf", runTime.name(), mesh), mesh, 1)
forAll(magnets, i)
Definition: createFields.H:55
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
List< label > labelList
A List of labels.
Definition: labelList.H:56
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
SurfaceField< scalar > surfaceScalarField
messageStream Info
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
IOerror FatalIOError
SurfaceField< vector > surfaceVectorField
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Foam::argList args(argc, argv)