[chery1986]

I am using porousInterFoam and porousWaveFoam in OF221, and want to output the drag force inside the porous media but cannot figure out the method. I define the porous media in blockMeshDict with a block
hex (1 2 6 5 13 14 18 17) porosity (73 93 1) simpleGrading (1 1 1)
and in porosityZone file define the property of it

1
(
    porosity
    {
        coordinateSystem
        {
            e1  (1 0 0);
            e2  (0 0 1);
        }

        resistanceFormulation vanGent1995;

        porosity 0.49;
        KC KC [0 0 0 0 0 0 0] 128;
        gammaAddedMass 0.34;
        
        d50 d50 [0 1 0 0 0 0 0] 0.0159;
        alpha alpha [0 0 0 0 0 0 0] 500;
        beta beta [0 0 0 0 0 0 0] 2.0;
    }
)

I found that if I use libforces.so,

    porousforces
    {
        type forces;
        functionObjectLibs ("libforces.so");
        outputControl timeStep;
        outputInterval     1;

        // Patches to sample
        patches (lowerWall);
        // Name of fields
        pName p;
        Uname U;
        // Density
        rhoName rho;
//        nuName    nuInf;
        rhoInf  1000;
//        nuInf 1e-06;
        // Dump to file
        log true;
        // Centre of rotation
        CofR (0 0 0);
    }

Then it will output the the force history with the first line # Time forces(pressure, viscous, porous) moment(pressure, viscous, porous), but the porous force is always (0, 0, 0).

Anyone know how to output the porous force in the porous media zone.

[RygeltheXVI]

I too am looking to work out porousity drag in openfoam and was wondering if you ever found a solution?

It doesn't make sense that a porous region is a cellzone to work out the forces on a patch (unless the zone is 1 cell wide), but even then you can't have the face bounding your porous region be patches else the flow wont go though them...

[styleo86]

Do you have found any solution to the problem?

I have the same request but I don't find a workaround.

[wyldckat]

Greetings to all!

I had this thread on my to-do list for a while now and finally managed to look into this.

I have to admit I was going to say at first that this wasn't possible, but after some searching and looking into the source code, here's what I've found:

1. Although the feature of calculating the porous part of the forces was already available in version 2.2.1, it was only announced as a feature in 2.3.0: http://www.openfoam.org/version2.3.0...processing.php
   

   forces and forceCoeffs include porosity contribution, and updated writing of bin data to single file(s);

2. If you look at this file: https://github.com/OpenFOAM/OpenFOAM.../forces.C#L681 - in line 681, you'll find that the missing optional option is named "porosity".

Therefore, since I don't have a test case myself, I can only state that the following line is missing in the function object you use in "controlDict":

porosity true;

Using chery1986's example:

    porousforces
    {
        type forces;
        functionObjectLibs ("libforces.so");
        outputControl timeStep;
        outputInterval     1;

        // Patches to sample
        patches (lowerWall);
        // Name of fields
        pName p;
        Uname U;
        // Density
        rhoName rho;
//        nuName    nuInf;
        rhoInf  1000;
//        nuInf 1e-06;
        // Dump to file
        log true;
        // Centre of rotation
        CofR (0 0 0);
        // Include porosity effects
        porosity   true;
    }

Best regards,
Bruno

[peppino]

I would like to revive this thread since I have not yet found an answer to the main question in the forum.
I attach a simple case to work with. The domain is a simple box with a porous zone in it and the flow solution seems reasonable.
The outputted forces though seems enormously high.
This may not be a surprise considering that i specified a cellZone (named framework) instead of a patch in the forces dictionary (see below).
Can I tell the functionObject that framework is a porous cellZone and not a patch? Or is there something else that is wrong?

forces_1
{type forces;
functionObjectLibs ("libforces.so"); 
porosity true;
patches (framework );    //not a patch!!
rhoName rhoInf;
rhoInf 1.2;
CofR (0 0 0);
outputControl   timeStep;
outputInterval  1;}

[peppino]

I seem to have found the reason for the extremely high forces, simply i forgot to scale the mesh before running...

Also, I noticed an error message when running the solver that complained about framework not being a patch. So I tested to remove all entries from the patch list in the dictionary (see below) and the error message disappeared. the porous forces are outputted anyway.

I presume this means that the functionObject will print out forces exherted on all porous zones in the domain.

  forces_1
{type forces;
functionObjectLibs ("libforces.so"); 
porosity true;
patches ( );    
rhoName rhoInf;
rhoInf 1.2;
CofR (0 0 0);
outputControl   timeStep;
outputInterval  1;}

[peppino]

next question is then: if there are several porous zones in the domain, is there a way to output the forces on a single porous zone rather than the sum of the forces on all zones? any suggestions are appreciated!

[wyldckat]

Originally Posted by peppino

next question is then: if there are several porous zones in the domain, is there a way to output the forces on a single porous zone rather than the sum of the forces on all zones? any suggestions are appreciated!

Quick question: Which OpenFOAM version are you using?

[peppino]

Originally Posted by wyldckat

Quick question: Which OpenFOAM version are you using?

Hi Bruno,
I am using OF 2.3.1.
By the way, thank you for your previous post where you discussed the keyword "porosity"

[wyldckat]

Quick answers:

Originally Posted by peppino

I am using OF 2.3.1.
By the way, thank you for your previous post where you discussed the keyword "porosity"

That makes it a bit clearer and you're welcome!

Originally Posted by peppino

next question is then: if there are several porous zones in the domain, is there a way to output the forces on a single porous zone rather than the sum of the forces on all zones? any suggestions are appreciated!

I've checked the code and the only way would be to create a custom function object. I don't have the time necessary to do it myself, but you can find a very similar example here: https://github.com/wyldckat/forceDirCoeffs/tree/of22x

edit: Sorry, forgot to mention that the piece of code you're looking for is here: https://github.com/OpenFOAM/OpenFOAM.../forces.C#L923 - starts in line 923, "fPTot" and "Md" is what you're looking for.

[peppino]

Thanks for the suggestion and the links. I doubt that I will be able to modify the code but I will give it a try later on. If i can get it to work i will post the result.

[Tempest]

Originally Posted by wyldckat

Quick answers:

That makes it a bit clearer and you're welcome!


I've checked the code and the only way would be to create a custom function object. I don't have the time necessary to do it myself, but you can find a very similar example here: https://github.com/wyldckat/forceDirCoeffs/tree/of22x

edit: Sorry, forgot to mention that the piece of code you're looking for is here: https://github.com/OpenFOAM/OpenFOAM.../forces.C#L923 - starts in line 923, "fPTot" and "Md" is what you're looking for.

Is it possible to calculate the porous Forces in OF without using the forceCoeffs macro in controlDict but to add a certain piece of code in the main solver?

[mechy]

Hi All
I have the following question about porosity part of forces
How OF calculated the forces applied on the porous zones ?
what formula is used for calculating the forces

thanks for any help

[wyldckat]

@mechy: Quick answer... Given my old answer:

Originally Posted by wyldckat

https://github.com/OpenFOAM/OpenFOAM.../forces.C#L923 - starts in line 923, "fPTot" and "Md" is what you're looking for.

1. https://cpp.openfoam.org/v6/
2. Search "porosityModel" --> https://cpp.openfoam.org/v6/classFoa...sityModel.html
3. Example, clicked on "DarcyForchheimer" --> https://cpp.openfoam.org/v6/DarcyFor...8C_source.html
4. Found the calcForce method that is used there: https://cpp.openfoam.org/v6/classFoa...fa6632b88fd61b
5. Clicked on the line number where it's defined:
   

   Definition at line 172 of file DarcyForchheimer.C.

6. And there you go: https://cpp.openfoam.org/v6/DarcyFor...ce.html#l00172

The calculation depends on the model.

[mechy]

Dear Bruno
thanks so much for your reply.
I have found that the force is calculated by integrating of Darcy-Forchimmer term over the porous volumes.

I need to calculate the force on the porous-fluid interface by integrating the pressure and viscous stress on the interface. however, the force library only gives the patches as its input. and this library can not give the interface for calculating the force.

I will be so thankful if you can help me to calculate the force on the fluid-porous interface

Best Regards

[edo2822]

Originally Posted by chery1986

I am using porousInterFoam and porousWaveFoam in OF221, and want to output the drag force inside the porous media but cannot figure out the method. I define the porous media in blockMeshDict with a block
hex (1 2 6 5 13 14 18 17) porosity (73 93 1) simpleGrading (1 1 1)
and in porosityZone file define the property of it

1
(
    porosity
    {
        coordinateSystem
        {
            e1  (1 0 0);
            e2  (0 0 1);
        }

        resistanceFormulation vanGent1995;

        porosity 0.49;
        KC KC [0 0 0 0 0 0 0] 128;
        gammaAddedMass 0.34;
        
        d50 d50 [0 1 0 0 0 0 0] 0.0159;
        alpha alpha [0 0 0 0 0 0 0] 500;
        beta beta [0 0 0 0 0 0 0] 2.0;
    }
)

I found that if I use libforces.so,

    porousforces
    {
        type forces;
        functionObjectLibs ("libforces.so");
        outputControl timeStep;
        outputInterval     1;

        // Patches to sample
        patches (lowerWall);
        // Name of fields
        pName p;
        Uname U;
        // Density
        rhoName rho;
//        nuName    nuInf;
        rhoInf  1000;
//        nuInf 1e-06;
        // Dump to file
        log true;
        // Centre of rotation
        CofR (0 0 0);
    }

Then it will output the the force history with the first line # Time forces(pressure, viscous, porous) moment(pressure, viscous, porous), but the porous force is always (0, 0, 0).

Anyone know how to output the porous force in the porous media zone.

Hey! I don't have the solution for this problem, but can I ask tyou which boundary condition you imposed for the boundaries of your porous zones?