[lukeweaver]

I'm trying to use an arbitrary "box" to calculate the flux around an actuator disc (reduced-order wind turbine model). However, the calculation for massflow over each surface is not being conserved.

The domain is a rectangular box with axial flow over a plane. The solver is simpleFoam.

The interior box is 80 x 160 x 160 m, and the surfaces are selected using the topoSet with slave cells and faceset as follows:

{
        name    north1SlaveCells;
        type    cellSet;
        action  new;
        source expressionToCell;
        sourceInfo
        {
            expression "(pos().x<=60 && pos().y>=20 && pos().y<=180 && pos().z>=20 && pos().z<=180)";
        }
}
{
        name    north1;
        type    faceSet;
        action  new;
        source  expressionToFace;
        sourceInfo
        {
            expression "(pos().x<=60 && pos().y>=20 && pos().y<=180 && pos().z>=20 && pos().z<=180)";
        }
}

This is done for all six sides (top, bottom,south,east,west). Here's three of the faceSets shown in paraView.

faceSet.jpg

Imagine the box fitting inside these three with the other three faceSets filling in the open space.

In the controlDict, I calculate the massflow:

PhiNorth1
        {
                type swakExpression;
                valueType faceSet;
                variables "direction=vector(0,1,0);dir=direction/mag(direction);";
                setName north1;
                //expression "U & Sf()";
                //expression "phi * flip()";
                expression "(U & direction) * area()";
                autoInterpolate true;
                accumulations (
                        sum
                        min //just for reference
                        max //just for reference
                );
                outputControlMode timeStep;
                outputInterval 10;
                verbose true;
        }

I used the swak tutorial angledDuctImplicitTransient files to build my topoSetDict and controlDict files.

I calculated the mass flow on each faceSet, and summing them should give Phi_In - Phi_Out = 0. After summing up the phi values on all faceSets, the values were Phi_In = 203,860 and sumPhi_Out 206,270, for a total difference of -2410.7 kg/s.

Any help/tips/explanations for this error?

EDIT: For the upstream and downstream faces, I use the U & Sf(). The phi*flip() condition didn't seem to work.

[gschaider]

Originally Posted by lukeweaver

I'm trying to use an arbitrary "box" to calculate the flux around an actuator disc (reduced-order wind turbine model). However, the calculation for massflow over each surface is not being conserved.

The domain is a rectangular box with axial flow over a plane. The solver is simpleFoam.

The interior box is 80 x 160 x 160 m, and the surfaces are selected using the topoSet with slave cells and faceset as follows:

{
        name    north1SlaveCells;
        type    cellSet;
        action  new;
        source expressionToCell;
        sourceInfo
        {
            expression "(pos().x<=60 && pos().y>=20 && pos().y<=180 && pos().z>=20 && pos().z<=180)";
        }
}
{
        name    north1;
        type    faceSet;
        action  new;
        source  expressionToFace;
        sourceInfo
        {
            expression "(pos().x<=60 && pos().y>=20 && pos().y<=180 && pos().z>=20 && pos().z<=180)";
        }
}

This is done for all six sides (top, bottom,south,east,west). Here's three of the faceSets shown in paraView.

Attachment 32645

Imagine the box fitting inside these three with the other three faceSets filling in the open space.

In the controlDict, I calculate the massflow:

PhiNorth1
        {
                type swakExpression;
                valueType faceSet;
                variables "direction=vector(0,1,0);dir=direction/mag(direction);";
                setName north1;
                //expression "U & Sf()";
                //expression "phi * flip()";
                expression "(U & direction) * area()";
                autoInterpolate true;
                accumulations (
                        sum
                        min //just for reference
                        max //just for reference
                );
                outputControlMode timeStep;
                outputInterval 10;
                verbose true;
        }

I used the swak tutorial angledDuctImplicitTransient files to build my topoSetDict and controlDict files.

I calculated the mass flow on each faceSet, and summing them should give Phi_In - Phi_Out = 0. After summing up the phi values on all faceSets, the values were Phi_In = 203,860 and sumPhi_Out 206,270, for a total difference of -2410.7 kg/s.

Any help/tips/explanations for this error?

EDIT: For the upstream and downstream faces, I use the U & Sf(). The phi*flip() condition didn't seem to work.

My guess is that you're counting flows that are not perpendicular to the boundary. If your north1 is the red faceSet in your picture then there is one large patch for which your approach would work but several small ones where for instance a velocity (1 0.1 0.5) would yield a contribution "0.1*area()" although either "1*area()" or "0.5*area()" (depending on the side) would be appropriate

[lukeweaver]

Bernard, thanks for the reply. It turned out that was the issue. I fixed it by using (U & dir) * (Sf() & dir), which gave the correct mass flow from each side (confirmed using Paraview calculators). Now for another question: can you think of a reason that phi * flip() doesn't work?

[gschaider]

Originally Posted by lukeweaver

Bernard, thanks for the reply. It turned out that was the issue. I fixed it by using (U & dir) * (Sf() & dir), which gave the correct mass flow from each side (confirmed using Paraview calculators). Now for another question: can you think of a reason that phi * flip() doesn't work?

No idea.

Check how many faces are flipped by evaluating "flip()" with accumulations "min", "max" and "average"

[lukeweaver]

I'm now running into issues with different mesh configurations. When the cell size is something like (4.0 2.0 2.0), the massflow conservation calculation is fine. However, when the cell size isn't an integer, like (4.1667 2.6667 2.6667) or something similar, I'm getting some weird results in the massflow and area calculations. The attached image shows the faceSets that have misplaced cell faces.

error_faceset.jpg

It seems that this is an issue of the faceSet boundaries that are defined in the topoSetDict not landing evenly on the cell boundaries.

So,
a) is there a way to change the topoSetDict so that the faceSet moves until it can read the faces evenly, and/or

b) is there a way to manipulate the swakExpression so that the cell irregularities don't affect the massflow calculations?

[gschaider]

Originally Posted by lukeweaver

I'm now running into issues with different mesh configurations. When the cell size is something like (4.0 2.0 2.0), the massflow conservation calculation is fine. However, when the cell size isn't an integer, like (4.1667 2.6667 2.6667) or something similar, I'm getting some weird results in the massflow and area calculations. The attached image shows the faceSets that have misplaced cell faces.

Attachment 32778

It seems that this is an issue of the faceSet boundaries that are defined in the topoSetDict not landing evenly on the cell boundaries.

So,
a) is there a way to change the topoSetDict so that the faceSet moves until it can read the faces evenly, and/or

b) is there a way to manipulate the swakExpression so that the cell irregularities don't affect the massflow calculations?

No idea what exactly happens but your "Lego-mesh" looks like a problem with the numeric comparison (if a number "looks" like 60 it may be something slightly different and 60==60 will fail). And as a consequence (this is a speculation) your face-zone and the slave-cell-zone may be inconsistent.

One way to be sure that they are consistent is to first create the cell-zone and instead of duplicating the rather complicated logical expression use it to create the face-zone with the simple expression "set(mySlaveZone)" (set is true if a cell belongs to the cellSet)