[qkrtjdrns]

I am trying to make a PID model to control valve using an immersed boundary method of foam-extend. The inlet velocity is fixed and the valve is controlled to adjust P_goal through PID algorithm. To do so, I need to calculate the pressure difference before and after the valve.
I calculated the pressure and the result is written /postProcessing/cuttingPlane/0/faceSourcePressure.dat, which is as below

#   Source : sampledSurface averagingPlane
#   Faces  : 1543
#   Area   : 2.825372e-03
#        Time     areaAverage(p)
2.000000e-03      -2.124464e+00
4.000000e-03      -3.515062e-01
6.000000e-03      2.320372e-01
8.000000e-03      2.114872e-01
1.000000e-02      2.317919e-01 ...

Here, I would like to read the pressure data at every time step to the source code, where PID control algorithm will be added.

#include "valveControl.H"
#include "addToRunTimeSelectionTable.H"
#include "mathematicalConstants.H"

using namespace Foam::mathematicalConstant;

// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

namespace Foam
{
namespace solidBodyMotionFunctions
{   
    defineTypeNameAndDebug(valvePIDControl, 0);
    addToRunTimeSelectionTable
    (   
        solidBodyMotionFunction,
        valvePIDControl,
        dictionary
    );
};
};


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

Foam::vector
Foam::solidBodyMotionFunctions::valvePIDControl::calcPosition
(
    const scalar t
) const
{
    return kP_*sin(2*pi*t/pGoal_);
}


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

Foam::solidBodyMotionFunctions::valvePIDControl::valvePIDControl
(
    const dictionary& SBMFCoeffs,
    const Time& runTime
)
:
    solidBodyMotionFunction(SBMFCoeffs, runTime),
    pGoal_(readScalar(SBMFCoeffs_.lookup("pGoal"))),
    kP_(SBMFCoeffs_.lookup("kP")),
    kI_(SBMFCoeffs_.lookup("kI")),
    kD_(SBMFCoeffs_.lookup("kD"))
{}


// * * * * * * * * * * * * * * * * Destructors * * * * * * * * * * * * * * * //

Foam::solidBodyMotionFunctions::valvePIDControl::~valvePIDControl()
{}


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

Foam::septernion
Foam::solidBodyMotionFunctions::valvePIDControl::transformation() const
{
    scalar t = time_.value();

    septernion TR(calcPosition(t), quaternion::I);

    Info<< "solidBodyMotionFunctions::valvePIDControl::transformation(): "
        << "Time = " << t << " transformation: " << TR << endl;

    return TR;
}


Foam::septernion
Foam::solidBodyMotionFunctions::valvePIDControl::velocity() const
{
    scalar t = time_.value();
    scalar dt = time_.deltaT().value();

    septernion TV
    (
        (calcPosition(t) - calcPosition(t - dt))/dt,
        quaternion::zero
    );

    return TV;
}


bool Foam::solidBodyMotionFunctions::valvePIDControl::read
(
    const dictionary& SBMFCoeffs
)
{
    solidBodyMotionFunction::read(SBMFCoeffs);

    SBMFCoeffs_.lookup("pGoal") >> pGoal_;
    SBMFCoeffs_.lookup("kP") >> kP_;
    SBMFCoeffs_.lookup("kI") >> kI_;
    SBMFCoeffs_.lookup("kD") >> kD_;

    return true;
}

In this case, how can I transfer pressure data at every time step written in /postProcessing/cuttingPlane/0/faceSourcePressure.dat to the above source code.

Sincerely
Goon

[Shadab_Ilahi]

Hi Goon,
Instead of reading the average pressure data from the dat file written by function objects, you can calculate it directly inside the code by performing the area average of pressure by looping over the all face of inlet patch/ or the face of any face zone.

Or

Incude the functionality of reading the file by ifstream class by giving the address of that file and read always the latest and last line of file.

I hope it solves your problem.

Cheers,
Shadab

[Shadab_Ilahi]

Here is the demo code you can use it something like this.

//Average pressure calculation at 2D location
        // Retrieve a reference to the pressure field


        const volScalarField& p = obr_.lookupObject<volScalarField>(pName_);

        // itnerpolate onto the faces
        surfaceScalarField pface = fvc::interpolate(p);

        scalar twoDPressure(0.0);
        scalar twoDArea(0.0);

        // create a twoDfaceId_ by directly reading the patch faces and then perform this calculations

        forAll(twoDfaceId_, faceI){

        twoDPressure += pface[twoDfaceId_[faceI]]* mesh_.magSf()[twoDfaceId_[faceI]];
        twoDArea += mesh_.magSf()[twoDfaceId_[faceI]];

        }


        // reduce for parallel running
        reduce(twoDPressure, sumOp<scalar>());


        reduce(twoDArea, sumOp<scalar>());


        twoDPressure =twoDPressure/twoDArea;


        Info << "2D average pressure on the "
             << returnReduce(twoDfaceId_.size(), sumOp<label>()) << " faces is [m2/s2]: " << twoDPressure << nl << endl;

I am assuming you are implementing your code as a function objet.

I hope it will solve your problem.


Regards,
Shadab

[qkrtjdrns]

I appreciate for your reply. I have three further questions.

1. Right. I think that I am implementing my code as a function object. I compile valve.c (the code that I attached) and make library file at $FOAM_USER_LIBBIN. And then add the library in the controlDict. Is it right what you mean?

2. If so, where can I add the codes that you wrote for me to my code valve.C?

3. In the part [// create a twoDfaceId_ by directly reading the patch faces], should I add additional dummy patches only to calculate the pressure? Because I need to calculate the averaged pressure at the two places (before and after valve), which are not inlet and outlet. When I use faceSource library to calculate the averaged pressure, I separately generated the patches at constant/triSurface/averagingPlane1.stl & constant/triSurface/averagingPlane2.stl. Can I still use the two patches there?

I appreciate for your reply again.

Thanks