[CedricVH]

It is not possible to calculate the mass flow rate through a faceSet in standard OpenFOAM. This is because the normal vectors of the faces have a random position, giving a result close to 0. Even after converting the set to cyclic patches using createBaffles, this problems persists as the normal vectors of the created patches do not have the same direction, unlike normal boundary patches.

However, I have read that swak4foam can help by using the flip() function:

    flow_pipe                    
    {                
        type    swakExpression;            
        valueType    faceSet;            
        setName    pipe_plane_faces;            
        expression    "phi*flip()";            
        accumulations    (            
                        sum
        );                
        verbose    true;    
        autoInterpolate true;            
    }

The flip() function however needs slaveCells, and this gives some problems.

My workflow is:

• Create a mesh using snappyHexMesh
• Create a cutting plane using a graphical package, and export as stl
• Start setSet
• Create a cellSet using surfaceToCell with that stl and keep the cells that are cut by the stl
• Convert cellSet to faceSet, keep all faces
• The stl plane has a thickness, so it cuts through multiple cells, leaving me with n layers of faces. If I would perform sum(phi) on these faces, the result would be n times to high. That's why I can only keep one layer of connected faces.
• Invert the cellSet and take a subSet of the faceSet using that cellSet leaves me with only the boundary faces
• Create a second cellSet using surfaceToCell with the stl and keep the cells inside
• Take a subSet of the faceSet using that cellSet leaves me with only one layer of faces that is connected to that cellSet

The result is now a faceSet consisting of one layer of faces. Now I have to create slaveCells. As I onderstand correctly, this should be a set of cells that enclose the faces at both directions. The simplest way would be creating a cellSet from this faceSet as every face would then connect to two cells.

cellSet pipe_plane_facesSlaveCells new faceToCell pipe_plane_faces all

However, when running my case, I get this error

--> FOAM FATAL ERROR: 
One of owner or neighbour of internal face 108712 should be in cellSet pipe_plane_facesSlaveCells to be able to determine orientation.
Face:108712 own:38634 OwnInCellSet:1 nei:39797 NeiInCellSet:1

    From function scalarField *FaceSetValueExpressionDriver::makeFaceFlipField()
    in file FaceSetValueExpressionDriver.C at line 233.

It seems that the slaveCells are not correctly defined.

Is my workflow correct or is there a much better workflow to correctly measure the mass flow rate through a faceSet/faceZone? If my workflow is correct, how can I manage to make the case run and get correct results?

[gschaider]

Originally Posted by CedricVH

It is not possible to calculate the mass flow rate through a faceSet in standard OpenFOAM. This is because the normal vectors of the faces have a random position, giving a result close to 0. Even after converting the set to cyclic patches using createBaffles, this problems persists as the normal vectors of the created patches do not have the same direction, unlike normal boundary patches.

However, I have read that swak4foam can help by using the flip() function:

    flow_pipe                    
    {                
        type    swakExpression;            
        valueType    faceSet;            
        setName    pipe_plane_faces;            
        expression    "phi*flip()";            
        accumulations    (            
                        sum
        );                
        verbose    true;    
        autoInterpolate true;            
    }

The flip() function however needs slaveCells, and this gives some problems.

My workflow is:

• Create a mesh using snappyHexMesh
• Create a cutting plane using a graphical package, and export as stl
• Start setSet
• Create a cellSet using surfaceToCell with that stl and keep the cells that are cut by the stl
• Convert cellSet to faceSet, keep all faces
• The stl plane has a thickness, so it cuts through multiple cells, leaving me with n layers of faces. If I would perform sum(phi) on these faces, the result would be n times to high. That's why I can only keep one layer of connected faces.
• Invert the cellSet and take a subSet of the faceSet using that cellSet leaves me with only the boundary faces
• Create a second cellSet using surfaceToCell with the stl and keep the cells inside
• Take a subSet of the faceSet using that cellSet leaves me with only one layer of faces that is connected to that cellSet

The result is now a faceSet consisting of one layer of faces. Now I have to create slaveCells. As I onderstand correctly, this should be a set of cells that enclose the faces at both directions. The simplest way would be creating a cellSet from this faceSet as every face would then connect to two cells.

No. Basically the slave-zone is "one side" of the faceSet.

An example: assume that your faceSet is the Bruxelles city-limit and you want to calculate the flow of phi out of Bruxelles. For the reason you gave above just summing up won't do it. Then for instance you could define whole Bruxelles as the slaveCells (the name is EXTREMLY unfortunate here). It is only important that each face in the set has ONE cell that is in the slaveCells. That way it knows which is the direction "to town" and calculates the flip() accordingly. Not every slaveCell has to be at the city-limits

Originally Posted by CedricVH

cellSet pipe_plane_facesSlaveCells new faceToCell pipe_plane_faces all

However, when running my case, I get this error

--> FOAM FATAL ERROR: 
One of owner or neighbour of internal face 108712 should be in cellSet pipe_plane_facesSlaveCells to be able to determine orientation.
Face:108712 own:38634 OwnInCellSet:1 nei:39797 NeiInCellSet:1

    From function scalarField *FaceSetValueExpressionDriver::makeFaceFlipField()
    in file FaceSetValueExpressionDriver.C at line 233.

It seems that the slaveCells are not correctly defined.

Is my workflow correct or is there a much better workflow to correctly measure the mass flow rate through a faceSet/faceZone? If my workflow is correct, how can I manage to make the case run and get correct results?

Not completely sure (never worked with the surfaceToCell-topoSource so I don't know what exactly it calculates)

Just one possibility (this might be computationally more expensive and it involves interpolation so you might have a small error): the sampled surfaces in OpenFOAM have one type sampledTriSurfaceMesh which is defined by an STL (which you have). And swak can calculate on sampledSurfaces. Of course you can't use phi (surfaceFields can't be interpolated) but you'll have to calculate the equivalent from U and the area (possibly rho if you're compressible)

[CedricVH]

No. Basically the slave-zone is "one side" of the faceSet.

If I understand correctly, the second cellSet that I have generated using surfaceToCell has to be the slaveCell set as this set includes the one layer faceSet but only on one side (the inside in my case) of the mesh. I will try this tomorrow!

[CedricVH]

I did not manage to get a correct mass flow rate through the faceSet using the flip() function. I'm going to visualise my problem using images.

The first attachment (faceSet_pipe.png) shows the system. At the right side there is an inlet that is connected through a pipe to a faceSet. This faceSet includes one layer of faces, thus sum(phi) inlet should be equal to sum(phi) faceSet. However, this is not the case:

Expression flow_inlet :  sum=-2.94504858156e-07
Expression flow_faceSet_noFlip :  sum=-2.07518040304e-07

This difference is due to the fact that the normals of the faceSet are not in the same direction as the normals of the patch. You can see in the image that the Y normals of the patch are all colored blue while the Y normals of the faceSet are colored blue or red.

The next step is to use the flip() function of swak4foam to correct this problem. In order calculate the flipMap, a slaveCellSet should be defined. This cellSet should be on one side of the faceSet en should include these faces. In the second image (faceSet.jpg), you see a detailed view of the the faceSet, in the third image (faceSetSlaveCellSet.jpg) you see the faceSet and the slaveCellSet. The total system is now given in the fourth image (faceSetSlaveCellSet.png).

If I now calculate the mass flow through this faceSet using the flip() function, I get:

Expression flow_inlet :  sum=-2.94504858156e-07
Expression flow_faceSet_noFlip :  sum=-2.07518040304e-07
Expression flow_faceSet_flip :  sum=-5.78237927163e-08

It can be seen that the difference between the flow though the inlet and the flow through the faceSet has become larger! Also, the sum of the noFlip and flip massflow sum(phi + phi*flip()) is not equal to the flow through the inlet.

Is there something that I missed?

[gschaider]

Originally Posted by CedricVH

I did not manage to get a correct mass flow rate through the faceSet using the flip() function. I'm going to visualise my problem using images.

The first attachment (faceSet_pipe.png) shows the system. At the right side there is an inlet that is connected through a pipe to a faceSet. This faceSet includes one layer of faces, thus sum(phi) inlet should be equal to sum(phi) faceSet. However, this is not the case:

Expression flow_inlet :  sum=-2.94504858156e-07
Expression flow_faceSet_noFlip :  sum=-2.07518040304e-07

This difference is due to the fact that the normals of the faceSet are not in the same direction as the normals of the patch. You can see in the image that the Y normals of the patch are all colored blue while the Y normals of the faceSet are colored blue or red.

The next step is to use the flip() function of swak4foam to correct this problem. In order calculate the flipMap, a slaveCellSet should be defined. This cellSet should be on one side of the faceSet en should include these faces. In the second image (faceSet.jpg), you see a detailed view of the the faceSet, in the third image (faceSetSlaveCellSet.jpg) you see the faceSet and the slaveCellSet. The total system is now given in the fourth image (faceSetSlaveCellSet.png).

If I now calculate the mass flow through this faceSet using the flip() function, I get:

Expression flow_inlet :  sum=-2.94504858156e-07
Expression flow_faceSet_noFlip :  sum=-2.07518040304e-07
Expression flow_faceSet_flip :  sum=-5.78237927163e-08

It can be seen that the difference between the flow though the inlet and the flow through the faceSet has become larger! Also, the sum of the noFlip and flip massflow sum(phi + phi*flip()) is not equal to the flow through the inlet.

Is there something that I missed?

No. It LOOKS alright (but it's a bit hard to tell without being able to manipulate the geometry)

BTW: sum(phi + phi*flip()) can't be equal to the inflow as you're summing up f*phi with f either being 0 or 2 (depending on the flip value)

Currently I have no idea what could be the problem. What I usually recommend is testing your method with a very simple blockMesh-channel because there it is easier to see what the matter is

BTW2: I assume that you're doing a steady simulation. Because for a transient simulation the mass-flows at a point in time do not necessarily have to be equal (in average they are). But even that would not explain the huge diffeence you have

[gschaider]

Originally Posted by CedricVH

I did not manage to get a correct mass flow rate through the faceSet using the flip() function. I'm going to visualise my problem using images.

The first attachment (faceSet_pipe.png) shows the system. At the right side there is an inlet that is connected through a pipe to a faceSet. This faceSet includes one layer of faces, thus sum(phi) inlet should be equal to sum(phi) faceSet. However, this is not the case:

Expression flow_inlet :  sum=-2.94504858156e-07
Expression flow_faceSet_noFlip :  sum=-2.07518040304e-07

This difference is due to the fact that the normals of the faceSet are not in the same direction as the normals of the patch. You can see in the image that the Y normals of the patch are all colored blue while the Y normals of the faceSet are colored blue or red.

The next step is to use the flip() function of swak4foam to correct this problem. In order calculate the flipMap, a slaveCellSet should be defined. This cellSet should be on one side of the faceSet en should include these faces. In the second image (faceSet.jpg), you see a detailed view of the the faceSet, in the third image (faceSetSlaveCellSet.jpg) you see the faceSet and the slaveCellSet. The total system is now given in the fourth image (faceSetSlaveCellSet.png).

If I now calculate the mass flow through this faceSet using the flip() function, I get:

Expression flow_inlet :  sum=-2.94504858156e-07
Expression flow_faceSet_noFlip :  sum=-2.07518040304e-07
Expression flow_faceSet_flip :  sum=-5.78237927163e-08

It can be seen that the difference between the flow though the inlet and the flow through the faceSet has become larger! Also, the sum of the noFlip and flip massflow sum(phi + phi*flip()) is not equal to the flow through the inlet.

Is there something that I missed?

No. It seems that this was a bug (Usage of a wrong data structure).

In the File Libraries/swak4FoamParsers/FaceSetValueExpressionDriver.C look for the line

SortableList<label> faceLabels(faceSet_->toc());

and replace it with

List<label> faceLabels(faceSet_->toc());

The fix has also been pushed to the mercurial-repository.

BTW: took me some time to find this posting again. That's why I appreciate it very much if people report such problems on the Mantis. Because that way I'm sure those interested get notified once I fixed it