View Issue Details
ID | Project | Category | View Status | Date Submitted | Last Update |
---|---|---|---|---|---|
0003473 | OpenFOAM | Bug | public | 2020-04-01 12:15 | 2020-04-01 17:05 |
Reporter | jkim2000 | Assigned To | henry | ||
Priority | low | Severity | minor | Reproducibility | always |
Status | resolved | Resolution | fixed | ||
Platform | Intel PC | OS | ubuntu | OS Version | 18.04 |
Fixed in Version | dev | ||||
Summary | 0003473: wallHeatFlux and P1 radiation | ||||
Description | During postProcessing with wallHeatFlux, it was found that P1 radiation heat flux is not added, and only the convective heat flux is calculated and saved in the field of "wallHeatFlux" It was expected that wallHeatFlux in the case with P1 radiation must be larger than WallHeatFlux in the case without radiation. But, the wallHeatFlux function gives same heat flux in the cases without radiation and with P1 radiation The problem is from that in P1.C code read-option of qr is NO_READ. On the contrary in FvDOM.C read-option of qr is READ_IF_PRESENT, which means that an object of the P1 radiation model created by wallHeatFlux function does not read the qr field. When the read-option of qr in P1.C code is changed to READ_IF_PRESENT, wallHeatFlux function gives the sum of convective and radiative heat flux | ||||
Tags | No tags attached. | ||||
|
P1.C (6,606 bytes)
/*---------------------------------------------------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | Website: https://openfoam.org \\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License This file is part of OpenFOAM. 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. 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 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. \*---------------------------------------------------------------------------*/ #include "P1.H" #include "fvmLaplacian.H" #include "fvmSup.H" #include "absorptionEmissionModel.H" #include "scatterModel.H" #include "constants.H" #include "addToRunTimeSelectionTable.H" using namespace Foam::constant; // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // namespace Foam { namespace radiationModels { defineTypeNameAndDebug(P1, 0); addToRadiationRunTimeSelectionTables(P1); } } // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // Foam::radiationModels::P1::P1(const volScalarField& T) : radiationModel(typeName, T), G_ ( IOobject ( "G", mesh_.time().timeName(), mesh_, IOobject::MUST_READ, IOobject::AUTO_WRITE ), mesh_ ), qr_ ( IOobject ( "qr", mesh_.time().timeName(), mesh_, IOobject::READ_IF_PRESENT, //NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar(dimMass/pow3(dimTime), 0) ), a_ ( IOobject ( "a", mesh_.time().timeName(), mesh_, IOobject::NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar(dimless/dimLength, 0) ), e_ ( IOobject ( "e", mesh_.time().timeName(), mesh_, IOobject::NO_READ, IOobject::NO_WRITE ), mesh_, dimensionedScalar(dimless/dimLength, 0) ), E_ ( IOobject ( "E", mesh_.time().timeName(), mesh_, IOobject::NO_READ, IOobject::NO_WRITE ), mesh_, dimensionedScalar(dimMass/dimLength/pow3(dimTime), 0) ) {} Foam::radiationModels::P1::P1(const dictionary& dict, const volScalarField& T) : radiationModel(typeName, dict, T), G_ ( IOobject ( "G", mesh_.time().timeName(), mesh_, IOobject::MUST_READ, IOobject::AUTO_WRITE ), mesh_ ), qr_ ( IOobject ( "qr", mesh_.time().timeName(), mesh_, IOobject::READ_IF_PRESENT, //NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar(dimMass/pow3(dimTime), 0) ), a_ ( IOobject ( "a", mesh_.time().timeName(), mesh_, IOobject::NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar(dimless/dimLength, 0) ), e_ ( IOobject ( "e", mesh_.time().timeName(), mesh_, IOobject::NO_READ, IOobject::NO_WRITE ), mesh_, dimensionedScalar(dimless/dimLength, 0) ), E_ ( IOobject ( "E", mesh_.time().timeName(), mesh_, IOobject::NO_READ, IOobject::NO_WRITE ), mesh_, dimensionedScalar(dimMass/dimLength/pow3(dimTime), 0) ) {} // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // Foam::radiationModels::P1::~P1() {} // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // bool Foam::radiationModels::P1::read() { if (radiationModel::read()) { // nothing to read return true; } else { return false; } } void Foam::radiationModels::P1::calculate() { a_ = absorptionEmission_->a(); e_ = absorptionEmission_->e(); E_ = absorptionEmission_->E(); const volScalarField sigmaEff(scatter_->sigmaEff()); const dimensionedScalar a0 ("a0", a_.dimensions(), rootVSmall); // Construct diffusion const volScalarField gamma ( IOobject ( "gammaRad", G_.mesh().time().timeName(), G_.mesh(), IOobject::NO_READ, IOobject::NO_WRITE ), 1.0/(3.0*a_ + sigmaEff + a0) ); // Solve G transport equation solve ( fvm::laplacian(gamma, G_) - fvm::Sp(a_, G_) == - 4.0*(e_*physicoChemical::sigma*pow4(T_) ) - E_ ); volScalarField::Boundary& qrBf = qr_.boundaryFieldRef(); // Calculate radiative heat flux on boundaries. forAll(mesh_.boundaryMesh(), patchi) { if (!G_.boundaryField()[patchi].coupled()) { qrBf[patchi] = -gamma.boundaryField()[patchi] *G_.boundaryField()[patchi].snGrad(); } } } Foam::tmp<Foam::volScalarField> Foam::radiationModels::P1::Rp() const { return volScalarField::New ( "Rp", 4.0*absorptionEmission_->eCont()*physicoChemical::sigma ); } Foam::tmp<Foam::DimensionedField<Foam::scalar, Foam::volMesh>> Foam::radiationModels::P1::Ru() const { const volScalarField::Internal& G = G_(); const volScalarField::Internal E = absorptionEmission_->ECont()()(); const volScalarField::Internal a = absorptionEmission_->aCont()()(); return a*G - E; } // ************************************************************************* // |
|
Resolved by commit 4867477252f4da773ea4bb5bd03d2740d88f2d00 |
Date Modified | Username | Field | Change |
---|---|---|---|
2020-04-01 12:15 | jkim2000 | New Issue | |
2020-04-01 12:15 | jkim2000 | File Added: P1.C | |
2020-04-01 17:05 | henry | Assigned To | => henry |
2020-04-01 17:05 | henry | Status | new => resolved |
2020-04-01 17:05 | henry | Resolution | open => fixed |
2020-04-01 17:05 | henry | Fixed in Version | => dev |
2020-04-01 17:05 | henry | Note Added: 0011279 |