View Issue Details

IDProjectCategoryView StatusLast Update
0001793OpenFOAMBugpublic2015-08-24 18:36
Reporteruser1190Assigned Tohenry  
PrioritynormalSeveritytrivialReproducibilityalways
Status resolvedResolutionfixed 
PlatformGNU/LinuxOSUbuntuOS Version14.04
Summary0001793: activePressureForceBaffleVelocityFvPatchVectorField seems wrong
DescriptionFrom line 299 we have :

if ((mag(valueDiff) > mag(minThresholdValue_) || baffleActivated_))
        {
            openFraction_ =
                max(
                    min(
                        openFraction_
                      + max
                        (
                          this->db().time().deltaT().value()/openingTime_,
                          maxOpenFractionDelta_
                        )*(orientation_),
                        1 - 1e-6
                        ),
                    1e-6
                    );

             baffleActivated_ = true;
        }


I have the feeling it should be :

if ((mag(valueDiff) > mag(minThresholdValue_) || baffleActivated_))
        {
            openFraction_ =
                max(
                    min(
                        openFraction_
                      + min
                        (
                          this->db().time().deltaT().value()/openingTime_,
                          maxOpenFractionDelta_
                        )*(orientation_*sign(valueDiff)),
                        1 - 1e-6
                        ),
                    1e-6
                    );

             baffleActivated_ = true;
        }
TagsNo tags attached.

Activities

wyldckat

2015-08-02 22:30

updater   ~0005169

Having a feeling in science/engineering/mathematics should always be tested, otherwise it's just a feeling.

Can you provide more specific details as to why you believe your proposition is correct?
What is the logic or tests you've done for pointing that consolidated this feeling?

will

2015-08-03 09:22

manager   ~0005170

I think they're saying that, max(/*expression*/, maxValue) is fundamentally faulty logic. It's clipping above the maxValue, not below it. min(/*expression*/, maxValue) is typical. Either that, or the name of the variable is misleading.

I'm not sure where sign(valueDiff) comes from, though. I agree, some proof would be nice, empirical or otherwise.

user1190

2015-08-03 11:31

  ~0005173

The sign(valueDiff) comes from what the bc is supposed to simulate : a check valve.

The valueDiff is the force or pressure difference between the two sides of the patch. Not including it's sign in the calculus of the openfraction would imply the valve to either constantly be opening or closing depending on the chosen value of the orientation.

I've tested a case twice, once with the sign(valueDiff), and once without it.
Without it the valve could only open, and with it the valve was properly oppening and closing.

will

2015-08-03 15:49

manager   ~0005176

Could you upload your case?

user1190

2015-08-04 07:49

  ~0005182

I'm sorry, but i'm not allowed by my company to share any of these cases.

If i find the time, i can make a simplified case though.

alexeym

2015-08-04 08:45

reporter   ~0005184

According to documentation in header file, maxOpenFractionDelta is max open fraction change per timestep. So this

                      + max
                        (
                          this->db().time().deltaT().value()/openingTime_,
                          maxOpenFractionDelta_
                        )*(orientation_),

should be really

                      + min
                        (
                          this->db().time().deltaT().value()/openingTime_,
                          maxOpenFractionDelta_
                        )*(orientation_),

Or in addition to misleading variable name, documentation in header file is wrong.

user21

2015-08-06 21:40

  ~0005203

I agree the code should look like:

        if ((valueDiff > minThresholdValue_) || baffleActivated_)
        {
            openFraction_ =
                max
                (
                    min
                    (
                        openFraction_
                      + min
                        (
                          this->db().time().deltaT().value()/openingTime_,
                          maxOpenFractionDelta_
                        )*(orientation_),
                        1 - 1e-6
                    ),
                    1e-6
                );

             baffleActivated_ = true;
        }

henry

2015-08-18 10:56

manager   ~0005257

With this change is everyone agreed that the implementation is correct and corresponds to the description in the header? If not I propose that this class is removed until we have a specification, implementation and description that are correct and consistent.

alexeym

2015-08-20 20:14

reporter   ~0005278

In fact "if not" branch is quite interesting, since in the report comments everybody was appealing to documentation in the header file, and there were no comments from people, who are actually using the feature.

henry

2015-08-20 20:33

manager   ~0005279

@oim: could you provide some feedback on the correctness/usefulness of this class?

Do the changes described here resolve the issues you were having?

Do you agree with the current documentation for the class?

If not can you propose corrections?

If none of the above I will remove the class.

wyldckat

2015-08-22 22:04

updater  

activePressureForceBaffleVelocityFvPatchVectorField.C (11,976 bytes)   
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | Copyright (C) 2011 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 "activePressureForceBaffleVelocityFvPatchVectorField.H"
#include "addToRunTimeSelectionTable.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "cyclicFvPatch.H"

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

Foam::activePressureForceBaffleVelocityFvPatchVectorField::
activePressureForceBaffleVelocityFvPatchVectorField
(
    const fvPatch& p,
    const DimensionedField<vector, volMesh>& iF
)
:
    fixedValueFvPatchVectorField(p, iF),
    pName_("p"),
    cyclicPatchName_(),
    cyclicPatchLabel_(-1),
    orientation_(1),
    initWallSf_(0),
    initCyclicSf_(0),
    nbrCyclicSf_(0),
    openFraction_(0),
    openingTime_(0),
    maxOpenFractionDelta_(0),
    curTimeIndex_(-1),
    minThresholdValue_(0),
    fBased_(1),
    baffleActivated_(0)
{}


Foam::activePressureForceBaffleVelocityFvPatchVectorField::
activePressureForceBaffleVelocityFvPatchVectorField
(
    const activePressureForceBaffleVelocityFvPatchVectorField& ptf,
    const fvPatch& p,
    const DimensionedField<vector, volMesh>& iF,
    const fvPatchFieldMapper& mapper
)
:
    fixedValueFvPatchVectorField(ptf, p, iF, mapper),
    pName_(ptf.pName_),
    cyclicPatchName_(ptf.cyclicPatchName_),
    cyclicPatchLabel_(ptf.cyclicPatchLabel_),
    orientation_(ptf.orientation_),
    initWallSf_(ptf.initWallSf_),
    initCyclicSf_(ptf.initCyclicSf_),
    nbrCyclicSf_(ptf.nbrCyclicSf_),
    openFraction_(ptf.openFraction_),
    openingTime_(ptf.openingTime_),
    maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
    curTimeIndex_(-1),
    minThresholdValue_(ptf.minThresholdValue_),
    fBased_(ptf.fBased_),
    baffleActivated_(ptf.baffleActivated_)
{}


Foam::activePressureForceBaffleVelocityFvPatchVectorField::
activePressureForceBaffleVelocityFvPatchVectorField
(
    const fvPatch& p,
    const DimensionedField<vector, volMesh>& iF,
    const dictionary& dict
)
:
    fixedValueFvPatchVectorField(p, iF),
    pName_(dict.lookupOrDefault<word>("p", "p")),
    cyclicPatchName_(dict.lookup("cyclicPatch")),
    cyclicPatchLabel_(p.patch().boundaryMesh().findPatchID(cyclicPatchName_)),
    orientation_(readLabel(dict.lookup("orientation"))),
    initWallSf_(0),
    initCyclicSf_(0),
    nbrCyclicSf_(0),
    openFraction_(readScalar(dict.lookup("openFraction"))),
    openingTime_(readScalar(dict.lookup("openingTime"))),
    maxOpenFractionDelta_(readScalar(dict.lookup("maxOpenFractionDelta"))),
    curTimeIndex_(-1),
    minThresholdValue_(readScalar(dict.lookup("minThresholdValue"))),
    fBased_(readBool(dict.lookup("forceBased"))),
    baffleActivated_(0)
{
    fvPatchVectorField::operator=(vector::zero);

    if (p.size() > 0)
    {
        initWallSf_ = p.Sf();
        initCyclicSf_ = p.boundaryMesh()[cyclicPatchLabel_].Sf();
        nbrCyclicSf_ =  refCast<const cyclicFvPatch>
        (
            p.boundaryMesh()[cyclicPatchLabel_]
        ).neighbFvPatch().Sf();
    }

    if (dict.found("p"))
    {
        dict.lookup("p") >> pName_;
    }
}


Foam::activePressureForceBaffleVelocityFvPatchVectorField::
activePressureForceBaffleVelocityFvPatchVectorField
(
    const activePressureForceBaffleVelocityFvPatchVectorField& ptf
)
:
    fixedValueFvPatchVectorField(ptf),
    pName_(ptf.pName_),
    cyclicPatchName_(ptf.cyclicPatchName_),
    cyclicPatchLabel_(ptf.cyclicPatchLabel_),
    orientation_(ptf.orientation_),
    initWallSf_(ptf.initWallSf_),
    initCyclicSf_(ptf.initCyclicSf_),
    nbrCyclicSf_(ptf.nbrCyclicSf_),
    openFraction_(ptf.openFraction_),
    openingTime_(ptf.openingTime_),
    maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
    curTimeIndex_(-1),
    minThresholdValue_(ptf.minThresholdValue_),
    fBased_(ptf.fBased_),
    baffleActivated_(ptf.baffleActivated_)
{}


Foam::activePressureForceBaffleVelocityFvPatchVectorField::
activePressureForceBaffleVelocityFvPatchVectorField
(
    const activePressureForceBaffleVelocityFvPatchVectorField& ptf,
    const DimensionedField<vector, volMesh>& iF
)
:
    fixedValueFvPatchVectorField(ptf, iF),
    pName_(ptf.pName_),
    cyclicPatchName_(ptf.cyclicPatchName_),
    cyclicPatchLabel_(ptf.cyclicPatchLabel_),
    orientation_(ptf.orientation_),
    initWallSf_(ptf.initWallSf_),
    initCyclicSf_(ptf.initCyclicSf_),
    nbrCyclicSf_(ptf.nbrCyclicSf_),
    openFraction_(ptf.openFraction_),
    openingTime_(ptf.openingTime_),
    maxOpenFractionDelta_(ptf.maxOpenFractionDelta_),
    curTimeIndex_(-1),
    minThresholdValue_(ptf.minThresholdValue_),
    fBased_(ptf.fBased_),
    baffleActivated_(ptf.baffleActivated_)
{}


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

void Foam::activePressureForceBaffleVelocityFvPatchVectorField::autoMap
(
    const fvPatchFieldMapper& m
)
{
    fixedValueFvPatchVectorField::autoMap(m);

    //- Note: cannot map field from cyclic patch anyway so just recalculate
    //  Areas should be consistent when doing autoMap except in case of
    //  topo changes.
    //- Note: we don't want to use Sf here since triggers rebuilding of
    //  fvMesh::S() which will give problems when mapped (since already
    //  on new mesh)
    forAll (patch().boundaryMesh().mesh().faceAreas(), i)
    {
        if (mag(patch().boundaryMesh().mesh().faceAreas()[i]) == 0)
        {
            Info << "faceArea[active] "<< i << endl;
        }
    }
    if (patch().size() > 0)
    {
        const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
        initWallSf_ = patch().patchSlice(areas);
        initCyclicSf_ = patch().boundaryMesh()
        [
            cyclicPatchLabel_
        ].patchSlice(areas);
        nbrCyclicSf_ = refCast<const cyclicFvPatch>
        (
            patch().boundaryMesh()
            [
                cyclicPatchLabel_
            ]
        ).neighbFvPatch().patch().patchSlice(areas);
    }
}

void Foam::activePressureForceBaffleVelocityFvPatchVectorField::rmap
(
    const fvPatchVectorField& ptf,
    const labelList& addr
)
{
    fixedValueFvPatchVectorField::rmap(ptf, addr);

    // See autoMap.
    const vectorField& areas = patch().boundaryMesh().mesh().faceAreas();
    initWallSf_ = patch().patchSlice(areas);
    initCyclicSf_ = patch().boundaryMesh()
    [
        cyclicPatchLabel_
    ].patchSlice(areas);
    nbrCyclicSf_ = refCast<const cyclicFvPatch>
    (
        patch().boundaryMesh()
        [
            cyclicPatchLabel_
        ]
    ).neighbFvPatch().patch().patchSlice(areas);
}


void Foam::activePressureForceBaffleVelocityFvPatchVectorField::updateCoeffs()
{
    if (updated())
    {
        return;
    }
    // Execute the change to the openFraction only once per time-step
    if (curTimeIndex_ != this->db().time().timeIndex())
    {
        const volScalarField& p = db().lookupObject<volScalarField>
        (
            pName_
        );

        const fvPatch& cyclicPatch = patch().boundaryMesh()[cyclicPatchLabel_];
        const labelList& cyclicFaceCells = cyclicPatch.patch().faceCells();
        const fvPatch& nbrPatch = refCast<const cyclicFvPatch>
        (
            cyclicPatch
        ).neighbFvPatch();

        const labelList& nbrFaceCells = nbrPatch.patch().faceCells();

        scalar valueDiff = 0;

        if (fBased_)
        {
             // Add this side
            forAll(cyclicFaceCells, facei)
            {
                valueDiff +=p[cyclicFaceCells[facei]]*mag(initCyclicSf_[facei]);
            }

            // Remove other side
            forAll(nbrFaceCells, facei)
            {
                valueDiff -=p[nbrFaceCells[facei]]*mag(initCyclicSf_[facei]);
            }

            Info<< "Force difference = " << valueDiff << endl;
        }
        else //pressure based
        {
            forAll(cyclicFaceCells, facei)
            {
                valueDiff += p[cyclicFaceCells[facei]];
            }

            forAll(nbrFaceCells, facei)
            {
                valueDiff -= p[nbrFaceCells[facei]];
            }

            Info<< "Pressure difference = " << valueDiff << endl;
        }

        if ((mag(valueDiff) > mag(minThresholdValue_)) || baffleActivated_)
        {
            openFraction_ =
                max(
                    min(
                        openFraction_
                      + min
                        (
                          this->db().time().deltaT().value()/openingTime_,
                          maxOpenFractionDelta_
                        )*(orientation_),
                        1 - 1e-6
                        ),
                    1e-6
                    );

             baffleActivated_ = true;
        }
        else
        {
            openFraction_ = max(min(1 - 1e-6, openFraction_), 1e-6);
        }

        Info<< "Open fraction = " << openFraction_ << endl;

        vectorField::subField Sfw = patch().patch().faceAreas();
        vectorField newSfw((1 - openFraction_)*initWallSf_);
        forAll(Sfw, facei)
        {
            Sfw[facei] = newSfw[facei];
        }
        const_cast<scalarField&>(patch().magSf()) = mag(patch().Sf());

        // Update owner side of cyclic
        const_cast<vectorField&>(cyclicPatch.Sf()) =
            openFraction_*initCyclicSf_;

        const_cast<scalarField&>(cyclicPatch.magSf()) =
            mag(cyclicPatch.Sf());

        // Update neighbour side of cyclic
        const_cast<vectorField&>(nbrPatch.Sf()) =
            openFraction_*nbrCyclicSf_;

        const_cast<scalarField&>(nbrPatch.magSf()) =
            mag(nbrPatch.Sf());

        curTimeIndex_ = this->db().time().timeIndex();
    }

    fixedValueFvPatchVectorField::updateCoeffs();
}


void Foam::activePressureForceBaffleVelocityFvPatchVectorField::
write(Ostream& os) const
{
    fvPatchVectorField::write(os);
    writeEntryIfDifferent<word>(os, "p", "p", pName_);
    os.writeKeyword("cyclicPatch")
        << cyclicPatchName_ << token::END_STATEMENT << nl;
    os.writeKeyword("orientation")
        << orientation_ << token::END_STATEMENT << nl;
    os.writeKeyword("openingTime")
        << openingTime_ << token::END_STATEMENT << nl;
    os.writeKeyword("maxOpenFractionDelta")
        << maxOpenFractionDelta_ << token::END_STATEMENT << nl;
    os.writeKeyword("openFraction")
        << openFraction_ << token::END_STATEMENT << nl;
    os.writeKeyword("minThresholdValue")
        << minThresholdValue_ << token::END_STATEMENT << nl;
    os.writeKeyword("forceBased")
        << fBased_ << token::END_STATEMENT << nl;
    writeEntry("value", os);
}


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

namespace Foam
{
    makePatchTypeField
    (
        fvPatchVectorField,
        activePressureForceBaffleVelocityFvPatchVectorField
    );
}

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

wyldckat

2015-08-22 22:04

updater  

activePressureForceBaffleVelocityFvPatchVectorField.H (8,504 bytes)   
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | Copyright (C) 2011-2012 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/>.

Class
    Foam::activePressureForceBaffleVelocityFvPatchVectorField

Group
    grpCoupledBoundaryConditions

Description
    This boundary condition is applied to the flow velocity, to simulate the
    opening or closure of a baffle due to local pressure or force changes,
    by merging the behaviours of wall and cyclic conditions.

    The baffle joins two mesh regions, where the open fraction determines
    the interpolation weights applied to each cyclic- and neighbour-patch
    contribution. This means that this is boundary condition is meant to be
    used in an extra wall beyond an existing cyclic patch pair. See PDRMesh
    for more details.

    Once the threshold is crossed, this condition activated and continues to
    open or close at a fixed rate using

        \f[
            x = x_{old} + s \times \frac{dt}{DT}
        \f]

    where

    \vartable
        x       | baffle open fraction [0-1]
        x_{old} | baffle open fraction on previous evaluation
        s       | sign for orientation: 1 to open or -1 to close
        dt      | simulation time step
        DT      | time taken to open the baffle
    \endvartable

    The open fraction is then applied to scale the patch areas.

    \heading Patch usage

    \table
        Property     | Description             | Required    | Default value
        p            | pressure field name     | no          | p
        cyclicPatch  | cyclic patch name       | yes         |
        orientation  | 1 to open or -1 to close | yes|
        openFraction | current open fraction [0-1] | yes |
        openingTime  | time taken to open or close the baffle | yes |
        maxOpenFractionDelta | max fraction change per timestep | yes |
        minThresholdValue | minimum absolute pressure or force difference for activation | yes |
        forceBased   | force (true) or pressure-based (false) activation | yes |
    \endtable

    Example of the boundary condition specification:
    \verbatim
    myPatch
    {
        type            activePressureForceBaffleVelocity;
        p               p;
        cyclicPatch     cyclic1;
        orientation     1;
        openFraction    0.2;
        openingTime     5.0;
        maxOpenFractionDelta 0.1;
        minThresholdValue 0.01;
        forceBased      false;
    }
    \endverbatim
    
SourceFiles
    activePressureForceBaffleVelocityFvPatchVectorField.C

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

#ifndef activePressureForceBaffleVelocityFvPatchVectorField_H
#define activePressureForceBaffleVelocityFvPatchVectorField_H

#include "fvPatchFields.H"
#include "fixedValueFvPatchFields.H"

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

namespace Foam
{

/*---------------------------------------------------------------------------*\
    Class activePressureForceBaffleVelocityFvPatchVectorField Declaration
\*---------------------------------------------------------------------------*/

class activePressureForceBaffleVelocityFvPatchVectorField
:
    public fixedValueFvPatchVectorField
{
    // Private data

        //- Name of the pressure field used to calculate the force
        //  on the active baffle
        word pName_;

        //- Name of the cyclic patch used when the active baffle is open
        word cyclicPatchName_;

        //- Index of the cyclic patch used when the active baffle is open
        label cyclicPatchLabel_;

        //- Orientation (1 or -1) of the active baffle patch.
        //  Used to change the direction of opening without the need for
        //  reordering the patch faces
        label orientation_;

        //- Initial wall patch areas
        vectorField initWallSf_;

        //- Initial cyclic patch areas
        vectorField initCyclicSf_;

        //- Initial neighbour-side cyclic patch areas
        vectorField nbrCyclicSf_;

        //- Current fraction of the active baffle which is open
        scalar openFraction_;

        //- Time taken for the active baffle to open
        scalar openingTime_;

        //- Maximum fractional change to the active baffle openness
        //  per time-step
        scalar maxOpenFractionDelta_;

        label curTimeIndex_;

        //- Minimum value for the active baffle to start opening
        scalar minThresholdValue_;

        //- Force based active baffle
        bool fBased_;

        //- Baffle is activated
        bool baffleActivated_;


public:

    //- Runtime type information
    TypeName("activePressureForceBaffleVelocity");


    // Constructors

        //- Construct from patch and internal field
        activePressureForceBaffleVelocityFvPatchVectorField
        (
            const fvPatch&,
            const DimensionedField<vector, volMesh>&
        );

        //- Construct from patch, internal field and dictionary
        activePressureForceBaffleVelocityFvPatchVectorField
        (
            const fvPatch&,
            const DimensionedField<vector, volMesh>&,
            const dictionary&
        );

        //- Construct by mapping
        activePressureForceBaffleVelocityFvPatchVectorField
        (
            const activePressureForceBaffleVelocityFvPatchVectorField&,
            const fvPatch&,
            const DimensionedField<vector, volMesh>&,
            const fvPatchFieldMapper&
        );

        //- Construct as copy
        activePressureForceBaffleVelocityFvPatchVectorField
        (
            const activePressureForceBaffleVelocityFvPatchVectorField&
        );

        //- Construct and return a clone
        virtual tmp<fvPatchVectorField> clone() const
        {
            return tmp<fvPatchVectorField>
            (
                new activePressureForceBaffleVelocityFvPatchVectorField(*this)
            );
        }

        //- Construct as copy setting internal field reference
        activePressureForceBaffleVelocityFvPatchVectorField
        (
            const activePressureForceBaffleVelocityFvPatchVectorField&,
            const DimensionedField<vector, volMesh>&
        );

        //- Construct and return a clone setting internal field reference
        virtual tmp<fvPatchVectorField> clone
        (
            const DimensionedField<vector, volMesh>& iF
        ) const
        {
            return tmp<fvPatchVectorField>
            (
                new activePressureForceBaffleVelocityFvPatchVectorField
                (
                    *this,
                    iF
                )
            );
        }


    // Member functions

        // Mapping functions

            //- Map (and resize as needed) from self given a mapping object
            virtual void autoMap
            (
                const fvPatchFieldMapper&
            );

            //- Reverse map the given fvPatchField onto this fvPatchField
            virtual void rmap
            (
                const fvPatchVectorField&,
                const labelList&
            );


        //- Update the coefficients associated with the patch field
        virtual void updateCoeffs();

        //- Write
        virtual void write(Ostream&) const;
};


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

} // End namespace Foam

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

#endif

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

wyldckat

2015-08-22 22:10

updater   ~0005281

Last edited: 2015-08-22 22:14

I was curious about this and I've taken a better look into this. Specially because I'm concerned that this kind of boundary condition (wall+cyclic hybrid) could be forgotten about if removed.

Here's what I've been able to diagnose:

 - The initial problem here is that this boundary condition was introduced in 2.0.0, without any clear explanation as for what it's really for. This missing information is critical for properly fixing this boundary condition.

 - The closest details available as to why it was introduced, seems to be because it's part of PDRFoam, PDRMesh and the tutorial case "combustion/PDRFoam/flamePropagationWithObstacles".

 - Studying the tutorial case and the description on the header file, it seems that this boundary condition was originally designed to act as a baffle that is "blown away" by a pressure difference above a threshold (i.e. explosion), regardless of being a positive or negative difference.

 - Based on this analysis, this means that this boundary condition was originally _not_ designed to act as a check valve. If anything, it seems to be something like a sheet of metal that gets torn up by the blast, at least within the tutorial.

 - The "orientation" parameter as it is implemented, is used for defining which way the fraction is "opened". In other words, if positive, the "openFraction" increases; if negative, it decreases. Problem here might be that if it's meant to "close" when a high pressure (or force) difference occurs, then this would mean that the baffle could act as a "blast door" that closes when the pressure difference is above the threshold.

 - The "maxOpenFractionDelta" is meant to act as a pseudo-relaxation factor. In other words, if for some reason the simulation deltaT is too large, this acts as a limiter that will restrain the rate at which the "open fraction" is increased/decreased.

 - The downside of simply removing this boundary condition is that the tutorial "combustion/PDRFoam/flamePropagationWithObstacles" needs to be readjusted to not use nor create the baffle, therefore loosing one of the apparent two reasons as to why PDRMesh exists.


Based on this assessment, I've attached the updated 2 files ("activePressureForceBaffleVelocityFvPatchVectorField.[CH]") that aim to fix what seems to be the original intention for this boundary condition, for placing in the folder "src/finiteVolume/fields/fvPatchFields/derived/activePressureForceBaffleVelocity". It uses Alexey's deduction (comment #0005184) and I've updated the description... oops, I only updated the description, and forgot to update the comments for the variables in the "private" section of the class structure... I'll upload the second updated header file in a few minutes.

I also vote that the name of the class and boundary condition is also changed to something more to the point of its existence, something like: triggeredPressureForceBaffleVelocity

wyldckat

2015-08-22 22:15

updater  

wyldckat

2015-08-22 22:20

updater   ~0005282

Attached "activePressureForceBaffleVelocityFvPatchVectorField.H.v2", which has the comment for "Orientation" updated in the private section of the class structure.

I forgot to mention in the previous comment that I've tested this with the tutorial "combustion/PDRFoam/flamePropagationWithObstacles", to double-check if it's working as intended, namely with "orientation = 1". I have not yet tested with "orientation = -1".

In addition, it's still missing a security check for whether the orientation value is 1 or -1... not sure if it's worth implementing such a check.

henry

2015-08-24 18:36

manager   ~0005293

@Bruno: Thanks for collating all the information and updating the files. I am happy to change the name if you are now convinced of the purpose of the BC and that a new name would be more consistent with this purpose.

I am not entirely convinced that the implementation is correct yet but I do not have a specification for this BC otherwise I would probably start again.

Resolved by commit 8685344a0b1ff2c91652121400fd1daaeeeea5ff

Issue History

Date Modified Username Field Change
2015-07-22 14:01 user1190 New Issue
2015-08-02 22:30 wyldckat Note Added: 0005169
2015-08-03 09:22 will Note Added: 0005170
2015-08-03 11:31 user1190 Note Added: 0005173
2015-08-03 15:49 will Note Added: 0005176
2015-08-04 07:49 user1190 Note Added: 0005182
2015-08-04 08:45 alexeym Note Added: 0005184
2015-08-06 21:40 user21 Note Added: 0005203
2015-08-18 10:56 henry Note Added: 0005257
2015-08-20 20:14 alexeym Note Added: 0005278
2015-08-20 20:33 henry Note Added: 0005279
2015-08-22 22:04 wyldckat File Added: activePressureForceBaffleVelocityFvPatchVectorField.C
2015-08-22 22:04 wyldckat File Added: activePressureForceBaffleVelocityFvPatchVectorField.H
2015-08-22 22:10 wyldckat Note Added: 0005281
2015-08-22 22:14 wyldckat Note Edited: 0005281
2015-08-22 22:15 wyldckat File Added: activePressureForceBaffleVelocityFvPatchVectorField.H.v2
2015-08-22 22:20 wyldckat Note Added: 0005282
2015-08-24 18:36 henry Note Added: 0005293
2015-08-24 18:36 henry Status new => resolved
2015-08-24 18:36 henry Resolution open => fixed
2015-08-24 18:36 henry Assigned To => henry