waveTransmissiveFvPatchField.C 4.74 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
6
     \\/     M anipulation  |
OpenFOAM bot's avatar
OpenFOAM bot committed
7
-------------------------------------------------------------------------------
OpenFOAM bot's avatar
OpenFOAM bot committed
8
    Copyright (C) 2011-2012 OpenFOAM Foundation
9
    Copyright (C) 2020 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 30 31 32 33

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

#include "waveTransmissiveFvPatchField.H"
#include "addToRunTimeSelectionTable.H"
#include "fvPatchFieldMapper.H"
#include "volFields.H"
#include "EulerDdtScheme.H"
34
#include "CrankNicolsonDdtScheme.H"
35 36 37 38 39
#include "backwardDdtScheme.H"

// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

template<class Type>
40
Foam::waveTransmissiveFvPatchField<Type>::waveTransmissiveFvPatchField
41 42 43 44 45 46
(
    const fvPatch& p,
    const DimensionedField<Type, volMesh>& iF
)
:
    advectiveFvPatchField<Type>(p, iF),
47
    psiName_("thermo:psi"),
48 49 50 51 52
    gamma_(0.0)
{}


template<class Type>
53
Foam::waveTransmissiveFvPatchField<Type>::waveTransmissiveFvPatchField
54 55 56 57 58 59 60 61 62 63 64 65 66 67
(
    const waveTransmissiveFvPatchField& ptf,
    const fvPatch& p,
    const DimensionedField<Type, volMesh>& iF,
    const fvPatchFieldMapper& mapper
)
:
    advectiveFvPatchField<Type>(ptf, p, iF, mapper),
    psiName_(ptf.psiName_),
    gamma_(ptf.gamma_)
{}


template<class Type>
68
Foam::waveTransmissiveFvPatchField<Type>::waveTransmissiveFvPatchField
69 70 71 72 73 74 75
(
    const fvPatch& p,
    const DimensionedField<Type, volMesh>& iF,
    const dictionary& dict
)
:
    advectiveFvPatchField<Type>(p, iF, dict),
76
    psiName_(dict.getOrDefault<word>("psi", "thermo:psi")),
77
    gamma_(dict.get<scalar>("gamma"))
78 79 80 81
{}


template<class Type>
82
Foam::waveTransmissiveFvPatchField<Type>::waveTransmissiveFvPatchField
83 84 85 86 87 88 89 90 91 92 93
(
    const waveTransmissiveFvPatchField& ptpsf
)
:
    advectiveFvPatchField<Type>(ptpsf),
    psiName_(ptpsf.psiName_),
    gamma_(ptpsf.gamma_)
{}


template<class Type>
94
Foam::waveTransmissiveFvPatchField<Type>::waveTransmissiveFvPatchField
95 96 97 98 99 100 101 102 103 104 105 106 107 108
(
    const waveTransmissiveFvPatchField& ptpsf,
    const DimensionedField<Type, volMesh>& iF
)
:
    advectiveFvPatchField<Type>(ptpsf, iF),
    psiName_(ptpsf.psiName_),
    gamma_(ptpsf.gamma_)
{}


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

template<class Type>
109 110
Foam::tmp<Foam::scalarField>
Foam::waveTransmissiveFvPatchField<Type>::advectionSpeed() const
111 112
{
    // Lookup the velocity and compressibility of the patch
113
    const fvPatchField<scalar>& psip =
114 115
        this->patch().template
            lookupPatchField<volScalarField, scalar>(psiName_);
116

117
    const surfaceScalarField& phi =
118
        this->db().template lookupObject<surfaceScalarField>(this->phiName_);
119

120
    fvsPatchField<scalar> phip =
121 122
        this->patch().template
            lookupPatchField<surfaceScalarField, scalar>(this->phiName_);
123 124 125

    if (phi.dimensions() == dimDensity*dimVelocity*dimArea)
    {
126
        const fvPatchScalarField& rhop =
127 128
            this->patch().template
                lookupPatchField<volScalarField, scalar>(this->rhoName_);
129 130 131 132 133 134 135 136 137 138 139 140

        phip /= rhop;
    }

    // Calculate the speed of the field wave w
    // by summing the component of the velocity normal to the boundary
    // and the speed of sound (sqrt(gamma_/psi)).
    return phip/this->patch().magSf() + sqrt(gamma_/psip);
}


template<class Type>
141
void Foam::waveTransmissiveFvPatchField<Type>::write(Ostream& os) const
142 143 144
{
    fvPatchField<Type>::write(os);

145 146 147
    os.writeEntryIfDifferent<word>("phi", "phi", this->phiName_);
    os.writeEntryIfDifferent<word>("rho", "rho", this->rhoName_);
    os.writeEntryIfDifferent<word>("psi", "thermo:psi", psiName_);
148

149
    os.writeEntry("gamma", gamma_);
150 151 152

    if (this->lInf_ > SMALL)
    {
153 154
        os.writeEntry("fieldInf", this->fieldInf_);
        os.writeEntry("lInf", this->lInf_);
155 156 157 158 159 160 161
    }

    this->writeEntry("value", os);
}


// ************************************************************************* //