basicThermo.C 12.8 KB
Newer Older
1
2
3
4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
OpenFOAM bot's avatar
OpenFOAM bot committed
5
6
7
8
    \\  /    A nd           | Copyright (C) 2017-2019 OpenCFD Ltd.
     \\/     M anipulation  |
-------------------------------------------------------------------------------
                            | Copyright (C) 2011-2016 OpenFOAM Foundation
9
10
11
12
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

13
14
15
16
    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.
17
18
19
20
21
22
23

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
24
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
25
26
27
28

\*---------------------------------------------------------------------------*/

#include "basicThermo.H"
29
30
31
32
33
34
35
36
#include "zeroGradientFvPatchFields.H"
#include "fixedEnergyFvPatchScalarField.H"
#include "gradientEnergyFvPatchScalarField.H"
#include "mixedEnergyFvPatchScalarField.H"
#include "fixedJumpFvPatchFields.H"
#include "fixedJumpAMIFvPatchFields.H"
#include "energyJumpFvPatchScalarField.H"
#include "energyJumpAMIFvPatchScalarField.H"
37

38

39
/* * * * * * * * * * * * * * * private static data * * * * * * * * * * * * * */
40
41
42

namespace Foam
{
43
    defineTypeNameAndDebug(basicThermo, 0);
44
    defineRunTimeSelectionTable(basicThermo, fvMesh);
45
    defineRunTimeSelectionTable(basicThermo, fvMeshDictPhase);
46
}
47

48
const Foam::word Foam::basicThermo::dictName("thermophysicalProperties");
49

50
51
52
53
54

// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //

Foam::wordList Foam::basicThermo::heBoundaryBaseTypes()
{
Andrew Heather's avatar
Andrew Heather committed
55
    const volScalarField::Boundary& tbf = this->T_.boundaryField();
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85

    wordList hbt(tbf.size(), word::null);

    forAll(tbf, patchi)
    {
        if (isA<fixedJumpFvPatchScalarField>(tbf[patchi]))
        {
            const fixedJumpFvPatchScalarField& pf =
                dynamic_cast<const fixedJumpFvPatchScalarField&>(tbf[patchi]);

            hbt[patchi] = pf.interfaceFieldType();
        }
        else if (isA<fixedJumpAMIFvPatchScalarField>(tbf[patchi]))
        {
            const fixedJumpAMIFvPatchScalarField& pf =
                dynamic_cast<const fixedJumpAMIFvPatchScalarField&>
                (
                    tbf[patchi]
                );

            hbt[patchi] = pf.interfaceFieldType();
        }
    }

    return hbt;
}


Foam::wordList Foam::basicThermo::heBoundaryTypes()
{
Andrew Heather's avatar
Andrew Heather committed
86
    const volScalarField::Boundary& tbf = this->T_.boundaryField();
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125

    wordList hbt = tbf.types();

    forAll(tbf, patchi)
    {
        if (isA<fixedValueFvPatchScalarField>(tbf[patchi]))
        {
            hbt[patchi] = fixedEnergyFvPatchScalarField::typeName;
        }
        else if
        (
            isA<zeroGradientFvPatchScalarField>(tbf[patchi])
         || isA<fixedGradientFvPatchScalarField>(tbf[patchi])
        )
        {
            hbt[patchi] = gradientEnergyFvPatchScalarField::typeName;
        }
        else if (isA<mixedFvPatchScalarField>(tbf[patchi]))
        {
            hbt[patchi] = mixedEnergyFvPatchScalarField::typeName;
        }
        else if (isA<fixedJumpFvPatchScalarField>(tbf[patchi]))
        {
            hbt[patchi] = energyJumpFvPatchScalarField::typeName;
        }
        else if (isA<fixedJumpAMIFvPatchScalarField>(tbf[patchi]))
        {
            hbt[patchi] = energyJumpAMIFvPatchScalarField::typeName;
        }
        else if (tbf[patchi].type() == "energyRegionCoupledFvPatchScalarField")
        {
            hbt[patchi] = "energyRegionCoupledFvPatchScalarField";
        }
    }

    return hbt;
}


126
127
// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

128
129
130
Foam::volScalarField& Foam::basicThermo::lookupOrConstruct
(
    const fvMesh& mesh,
131
132
133
    const word& name,
    bool& isOwner
)
134
{
135
136
    volScalarField* ptr =
        mesh.objectRegistry::getObjectPtr<volScalarField>(name);
137

138
    isOwner = !ptr;
139

140
    if (!ptr)
141
    {
142
        ptr = new volScalarField
143
        (
144
            IOobject
145
            (
146
147
148
149
150
151
152
                name,
                mesh.time().timeName(),
                mesh,
                IOobject::MUST_READ,
                IOobject::AUTO_WRITE
            ),
            mesh
153
154
155
        );

        // Transfer ownership of this object to the objectRegistry
156
        ptr->store();
157
    }
158
159

    return *ptr;
160
161
162
163
164
165
166
167
}


Foam::basicThermo::basicThermo
(
    const fvMesh& mesh,
    const word& phaseName
)
168
:
169
170
171
172
    IOdictionary
    (
        IOobject
        (
173
            phasePropertyName(dictName, phaseName),
174
175
176
177
178
179
180
            mesh.time().constant(),
            mesh,
            IOobject::MUST_READ_IF_MODIFIED,
            IOobject::NO_WRITE
        )
    ),

181
182
    phaseName_(phaseName),

183
    p_(lookupOrConstruct(mesh, "p", pOwner_)),
184

185
186
    T_(lookupOrConstruct(mesh, phasePropertyName("T"), TOwner_)),
    TOwner_(lookupOrDefault<Switch>("updateT", TOwner_)),
187
188
189
190
191

    alpha_
    (
        IOobject
        (
192
            phasePropertyName("thermo:alpha"),
193
194
            mesh.time().timeName(),
            mesh,
195
            IOobject::READ_IF_PRESENT,
196
197
198
            IOobject::NO_WRITE
        ),
        mesh,
199
        dimensionedScalar(dimensionSet(1, -1, -1, 0, 0), Zero)
200
    ),
201

202
    dpdt_(lookupOrDefault<Switch>("dpdt", true))
203
{}
204
205


206
207
208
Foam::basicThermo::basicThermo
(
    const fvMesh& mesh,
209
210
    const dictionary& dict,
    const word& phaseName
211
)
212
:
213
214
215
216
    IOdictionary
    (
        IOobject
        (
217
            phasePropertyName(dictName, phaseName),
218
219
220
221
222
223
224
225
            mesh.time().constant(),
            mesh,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
        dict
    ),

226
227
    phaseName_(phaseName),

228
    p_(lookupOrConstruct(mesh, "p", pOwner_)),
229

230
231
    T_(lookupOrConstruct(mesh, phasePropertyName("T"), TOwner_)),
    TOwner_(lookupOrDefault<Switch>("updateT", TOwner_)),
232
233

    alpha_
234
235
236
    (
        IOobject
        (
237
            phasePropertyName("thermo:alpha"),
238
239
            mesh.time().timeName(),
            mesh,
240
241
            IOobject::NO_READ,
            IOobject::NO_WRITE
242
        ),
243
        mesh,
244
        dimensionedScalar(dimensionSet(1, -1, -1, 0, 0), Zero)
245
246
247
    ),

    dpdt_(lookupOrDefault<Switch>("dpdt", true))
248
249
250
{}


251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
Foam::basicThermo::basicThermo
(
    const fvMesh& mesh,
    const word& phaseName,
    const word& dictionaryName
)
:
    IOdictionary
    (
        IOobject
        (
            dictionaryName,
            mesh.time().constant(),
            mesh,
            IOobject::MUST_READ_IF_MODIFIED,
            IOobject::NO_WRITE
        )
    ),

    phaseName_(phaseName),

    p_(lookupOrConstruct(mesh, "p", pOwner_)),

    T_(lookupOrConstruct(mesh, "T", TOwner_)),
    TOwner_(lookupOrDefault<Switch>("updateT", TOwner_)),

    alpha_
    (
        IOobject
        (
            "thermo:alpha",
            mesh.time().timeName(),
            mesh,
            IOobject::READ_IF_PRESENT,
            IOobject::NO_WRITE
        ),
        mesh,
288
        dimensionedScalar(dimensionSet(1, -1, -1, 0, 0), Zero)
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
    ),

    dpdt_(lookupOrDefault<Switch>("dpdt", true))
{
    if (debug)
    {
        Pout<< "Constructed shared thermo : mesh:" << mesh.name()
            << " phase:" << phaseName
            << " dictionary:" << dictionaryName
            << " T:" << T_.name()
            << " updateT:" << TOwner_
            << " alphaName:" << alpha_.name()
            << endl;
    }
}


306
307
// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //

308
309
Foam::autoPtr<Foam::basicThermo> Foam::basicThermo::New
(
310
311
    const fvMesh& mesh,
    const word& phaseName
312
313
)
{
314
    return New<basicThermo>(mesh, phaseName);
315
316
317
}


318
319
// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //

320
Foam::basicThermo::~basicThermo()
321
{
322
    db().checkOut("p");
323
}
324
325


326
327
// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //

328
329
330
331
332
const Foam::basicThermo& Foam::basicThermo::lookupThermo
(
    const fvPatchScalarField& pf
)
{
333
334
335
    const basicThermo* thermo = pf.db().findObject<basicThermo>(dictName);

    if (thermo)
336
    {
337
        return *thermo;
338
339
    }

340
341
342
343
344
345
    HashTable<const basicThermo*> thermos =
        pf.db().lookupClass<basicThermo>();

    forAllConstIters(thermos, iter)
    {
        if
346
        (
347
348
            &(iter()->he().internalField())
         == &(pf.internalField())
349
350
        )
        {
351
            return *iter();
352
353
354
355
356
357
358
        }
    }

    return pf.db().lookupObject<basicThermo>(dictName);
}


359
360
void Foam::basicThermo::validate
(
361
    const string& app,
362
363
364
    const word& a
) const
{
365
    if (!(he().name() == phasePropertyName(a)))
366
    {
367
        FatalErrorInFunction
368
            << "Supported energy type is " << phasePropertyName(a)
369
370
371
372
373
374
375
            << ", thermodynamics package provides " << he().name()
            << exit(FatalError);
    }
}

void Foam::basicThermo::validate
(
376
    const string& app,
377
378
379
380
    const word& a,
    const word& b
) const
{
381
382
383
384
385
386
387
    if
    (
       !(
            he().name() == phasePropertyName(a)
         || he().name() == phasePropertyName(b)
        )
    )
388
    {
389
        FatalErrorInFunction
390
391
            << "Supported energy types are " << phasePropertyName(a)
            << " and " << phasePropertyName(b)
392
393
394
395
396
397
398
            << ", thermodynamics package provides " << he().name()
            << exit(FatalError);
    }
}

void Foam::basicThermo::validate
(
399
    const string& app,
400
401
402
403
404
405
406
407
    const word& a,
    const word& b,
    const word& c
) const
{
    if
    (
       !(
408
409
410
            he().name() == phasePropertyName(a)
         || he().name() == phasePropertyName(b)
         || he().name() == phasePropertyName(c)
411
412
413
        )
    )
    {
414
        FatalErrorInFunction
415
416
417
            << "Supported energy types are " << phasePropertyName(a)
            << ", " << phasePropertyName(b)
            << " and " << phasePropertyName(c)
418
419
420
421
422
423
424
            << ", thermodynamics package provides " << he().name()
            << exit(FatalError);
    }
}

void Foam::basicThermo::validate
(
425
    const string& app,
426
427
428
429
430
431
432
433
434
    const word& a,
    const word& b,
    const word& c,
    const word& d
) const
{
    if
    (
       !(
435
436
437
438
            he().name() == phasePropertyName(a)
         || he().name() == phasePropertyName(b)
         || he().name() == phasePropertyName(c)
         || he().name() == phasePropertyName(d)
439
440
441
        )
    )
    {
442
        FatalErrorInFunction
443
444
445
446
            << "Supported energy types are " << phasePropertyName(a)
            << ", " << phasePropertyName(b)
            << ", " << phasePropertyName(c)
            << " and " << phasePropertyName(d)
447
448
449
450
451
452
            << ", thermodynamics package provides " << he().name()
            << exit(FatalError);
    }
}


453
454
455
456
457
458
459
460
Foam::wordList Foam::basicThermo::splitThermoName
(
    const word& thermoName,
    const int nCmpt
)
{
    wordList cmpts(nCmpt);

461
    string::size_type beg=0, end=0, endb=0, endc=0;
462
463
464
465
    int i = 0;

    while
    (
466
467
        (endb = thermoName.find('<', beg)) != string::npos
     || (endc = thermoName.find(',', beg)) != string::npos
468
469
    )
    {
470
471
472
473
474
475
        if (endb == string::npos)
        {
            end = endc;
        }
        else if ((endc = thermoName.find(',', beg)) != string::npos)
        {
476
            end = std::min(endb, endc);
477
478
479
480
481
482
        }
        else
        {
            end = endb;
        }

483
484
485
486
        if (beg < end)
        {
            cmpts[i] = thermoName.substr(beg, end-beg);
            cmpts[i++].replaceAll(">","");
487
488
489
490
491
492
493

            // If the number of number of components in the name
            // is greater than nCmpt return an empty list
            if (i == nCmpt)
            {
                return wordList::null();
            }
494
495
496
497
        }
        beg = end + 1;
    }

498
499
500
501
502
503
504
    // If the number of number of components in the name is not equal to nCmpt
    // return an empty list
    if (i + 1 != nCmpt)
    {
        return wordList::null();
    }

505
506
507
    if (beg < thermoName.size())
    {
        cmpts[i] = thermoName.substr(beg, string::npos);
508
        cmpts[i].replaceAll(">","");
509
510
511
512
513
514
    }

    return cmpts;
}


515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
Foam::volScalarField& Foam::basicThermo::p()
{
    return p_;
}


const Foam::volScalarField& Foam::basicThermo::p() const
{
    return p_;
}


const Foam::volScalarField& Foam::basicThermo::T() const
{
    return T_;
}


533
534
535
536
537
538
Foam::volScalarField& Foam::basicThermo::T()
{
    return T_;
}


539
540
541
542
543
544
const Foam::volScalarField& Foam::basicThermo::alpha() const
{
    return alpha_;
}


545
546
547
548
549
550
const Foam::scalarField& Foam::basicThermo::alpha(const label patchi) const
{
    return alpha_.boundaryField()[patchi];
}


551
552
553
554
555
556
bool Foam::basicThermo::read()
{
    return regIOobject::read();
}


557
// ************************************************************************* //