basicThermo.C 12.7 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
    \\  /    A nd           | www.openfoam.com
OpenFOAM bot's avatar
OpenFOAM bot committed
6 7
     \\/     M anipulation  |
-------------------------------------------------------------------------------
OpenFOAM bot's avatar
OpenFOAM bot committed
8 9
    Copyright (C) 2011-2016 OpenFOAM Foundation
    Copyright (C) 2017-2019 OpenCFD Ltd.
10 11 12 13
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

14 15 16 17
    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.
18 19 20 21 22 23 24

    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
25
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
26 27 28 29

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

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

39

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

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

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

51 52 53 54 55

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

Foam::wordList Foam::basicThermo::heBoundaryBaseTypes()
{
Andrew Heather's avatar
Andrew Heather committed
56
    const volScalarField::Boundary& tbf = this->T_.boundaryField();
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 86

    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
87
    const volScalarField::Boundary& tbf = this->T_.boundaryField();
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

    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;
        }
    }

    return hbt;
}


123 124
// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

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

135
    isOwner = !ptr;
136

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

        // Transfer ownership of this object to the objectRegistry
153
        ptr->store();
154
    }
155 156

    return *ptr;
157 158 159 160 161 162 163 164
}


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

178 179
    phaseName_(phaseName),

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

182 183
    T_(lookupOrConstruct(mesh, phasePropertyName("T"), TOwner_)),
    TOwner_(lookupOrDefault<Switch>("updateT", TOwner_)),
184 185 186 187 188

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

199
    dpdt_(lookupOrDefault<Switch>("dpdt", true))
200
{}
201 202


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

223 224
    phaseName_(phaseName),

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

227 228
    T_(lookupOrConstruct(mesh, phasePropertyName("T"), TOwner_)),
    TOwner_(lookupOrDefault<Switch>("updateT", TOwner_)),
229 230

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

    dpdt_(lookupOrDefault<Switch>("dpdt", true))
245 246 247
{}


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
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
        (
278
            phasePropertyName("thermo:alpha"),
279 280 281 282 283 284
            mesh.time().timeName(),
            mesh,
            IOobject::READ_IF_PRESENT,
            IOobject::NO_WRITE
        ),
        mesh,
285
        dimensionedScalar(dimensionSet(1, -1, -1, 0, 0), Zero)
286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302
    ),

    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;
    }
}


303 304
// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //

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


315 316
// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //

317
Foam::basicThermo::~basicThermo()
318
{
319
    db().checkOut("p");
320
}
321 322


323 324
// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //

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

    if (thermo)
333
    {
334
        return *thermo;
335 336
    }

337 338 339 340 341 342
    HashTable<const basicThermo*> thermos =
        pf.db().lookupClass<basicThermo>();

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

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


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

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

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

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


450 451 452 453 454 455 456 457
Foam::wordList Foam::basicThermo::splitThermoName
(
    const word& thermoName,
    const int nCmpt
)
{
    wordList cmpts(nCmpt);

458
    string::size_type beg=0, end=0, endb=0, endc=0;
459 460 461 462
    int i = 0;

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

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

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

495 496 497 498 499 500 501
    // 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();
    }

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

    return cmpts;
}


512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529
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_;
}


530 531 532 533 534 535
Foam::volScalarField& Foam::basicThermo::T()
{
    return T_;
}


536 537 538 539 540 541
const Foam::volScalarField& Foam::basicThermo::alpha() const
{
    return alpha_;
}


542 543 544 545 546 547
const Foam::scalarField& Foam::basicThermo::alpha(const label patchi) const
{
    return alpha_.boundaryField()[patchi];
}


548 549 550 551 552 553
bool Foam::basicThermo::read()
{
    return regIOobject::read();
}


554
// ************************************************************************* //