Calculation of integral boundary values

 Register Blogs Members List Search Today's Posts Mark Forums Read

 January 12, 2005, 04:01 How can I calculate the - #1 Joern Beilke (Beilke) Guest   Posts: n/a How can I calculate the - mass flux averaged pressure for a boundary region - mass flow rate for a region (sum, inflow, outflow) during the calculation?

 January 12, 2005, 06:54 Have a look at this code examp #2 Eugene de Villiers (Eugene) Guest   Posts: n/a Have a look at this code example: /********************************************/ scalar volumeFluxWeightedPressureInlet = 0.0; scalar volumeFluxInlet = 0.0; word inletPatchName = "inlet"; //access boundary elements const fvPatchList& patches = mesh.boundary(); //loop over boundaries forAll(patches, patchI) {
const fvPatch cPatch = patches[patchI]; //check boundary name if(cPatch.name() == inletPatchName) {
//Access boundary face area vectors vectorField inletAreaVectors = cPatch.faceAreas(); //loop over all faces of selected //boundary forAll(cPatch, faceI) {
//calculate boundary face flux scalar cFaceFlux = U.boundaryField()[patchI][faceI] & inletAreaVectors[faceI]; //sum face fluxes volumeFluxInlet += cFaceFlux; //sum flux weighted face pressures volumeFluxWeightedPressureInlet += p.boundaryField()[patchI][faceI] * cFaceFlux;
}
volumeFluxWeightedPressureInlet \= volumeFluxInlet;
}
} /**********************************************/ Here "mesh" is the fvMesh grid, U is the velocity field (vector), p is the pressure field and "&" is the dot product operator. The same general procedure can be used to access, modify and otherwise manipulate all types of boundary field values. For example, it can also be used to specify inlet velocities through face centre locations and or time dependent functions. See "src/OpenFOAM/lnInclude/fvPatch.H" and its base classes in the OpenFOAM directory or the Doxygen documentation for more relevant access functions. Eugene iLyrics likes this.

 January 12, 2005, 15:32 Hi Eugene, I tried to get #3 Joern Beilke (Beilke) Guest   Posts: n/a Hi Eugene, I tried to get it running but I still get one error message: error: 'const class Foam::fvPatch' has no member named 'faceAreas' I think that I got an idea what you are doing in the example, but I'm a total c++ beginner ;-) I will also send you the example via email. Thanks very much for your help Jörn

 January 12, 2005, 19:22 Try the function Sf() instead #4 Eugene de Villiers (Eugene) Guest   Posts: n/a Try the function Sf() instead (faceAreas is a polyPatch member function, which I mistakenly thought would be inherited by fvPatch). For a C++ beginner, the Doxygen documentation can be extremely valuable; providing an easy to understand summary of class member functions and inheritance. It also has a built-in search engine, which saved me loads of time when I first started to use Foam and was still struggling with the basics of C++. In response to your e-mail, yes "rho.boundaryField()[patchI][faceI]" will give you the density at a specified face on the specified boundary patch. Eugene

 January 13, 2005, 04:50 I don't agree with Eugene's e #5 Henry Weller (Henry) Guest   Posts: n/a I don't agree with Eugene's example as a method of calculating the patch mass flow rate. It is the flux field phi which contains the mass flux for compressible codes or the volumetic flux for incompressible codes. Also it is not usually neccessary to loop over faces; there is a large set of field functions included in FOAM which sould be used in preference to low-level looping where ever possible. Thus the inlet mass flux can simply be obtained from: scalar massFlux = sum(phi.boundaryField()[inletPatchi]); where the inletPatchi is the index of the inlet patch which can be obtined by searching as in Eugene's example or from label inletPatchi = mesh.boundaryMesh().findPatchID(nameOfInlet); where nameOfInlet is the name of the inlet patch which could be provided from a dictionary etc. sharonyue likes this.

 January 13, 2005, 08:33 Thanks Eugene and Henry, a #6 Joern Beilke (Beilke) Guest   Posts: n/a Thanks Eugene and Henry, after removing some typos and changing to the Sf() function the compilation and calculation work perfectly well. I already had a look into the Doxygen documentation but for a c++ beginner it looks a bit like "chinese backwards" ;-) For those who want to use it I will include the running version here: /********************************************/ scalar massFluxWeightedPressureInlet = 0.0 ; scalar massFluxInlet = 0.0 ; word inletPatchName = "INLET"; //access boundary elements const fvPatchList& patches = mesh.boundary(); //loop over boundaries forAll(patches, patchI) { const fvPatch& cPatch = patches[patchI]; //check boundary name if(cPatch.name() == inletPatchName) { //Access boundary face area vectors const vectorField& inletAreaVectors = cPatch.Sf(); //loop over all faces of selected //boundary forAll(cPatch, faceI) { //calculate boundary face mass flux scalar cFaceFlux = U.boundaryField()[patchI][faceI] & inletAreaVectors[faceI] * rho.boundaryField()[patchI][faceI]; //sum face mass fluxes massFluxInlet += cFaceFlux; //sum mass flux weighted face pressures massFluxWeightedPressureInlet += p.boundaryField()[patchI][faceI] * cFaceFlux; } massFluxWeightedPressureInlet /= massFluxInlet; } } Info<< "MassFlowRate: " << massFluxInlet << " kg/s" << endl; Info<< "Pressure : " << massFluxWeightedPressureInlet << " Pa " << endl; /**********************************************/

 January 13, 2005, 08:46 I still think this is an abso #7 Henry Weller (Henry) Guest   Posts: n/a I still think this is an absolutely terrible way of doing things in FOAM. I posted a much more elegant and appropriate way of calculating the patch-flux and I must warn you again that for outlets in particular the flux obtained from the velocity in this way does not obey continuity.

 January 13, 2005, 09:28 Ok Henry, I replaced the #8 Joern Beilke (Beilke) Guest   Posts: n/a Ok Henry, I replaced the //calculate boundary face mass flux scalar cFaceFlux = U.boundaryField()[patchI][faceI] & inletAreaVectors[faceI] * rho.boundaryField()[patchI][faceI]; by scalar cFaceFlux = phi.boundaryField()[patchI][faceI]; Your suggestion works for calculating the flux without looping the faces. But is there a way to calculate the mass flow averaged pressure without this loop? PS: the pressure difference in now below 2% with ke and ud.

 January 13, 2005, 10:15 You can do nearly anything in #9 Henry Weller (Henry) Guest   Posts: n/a You can do nearly anything in FOAM without looping. Imagine you have a field psi and want a flux-weighted average of it over patch patchi: sum(phi.boundaryField()[patchi]*psi.boundaryField()[patchi])/sum(phi.boundaryField()[patchi])

 January 13, 2005, 11:08 I guess that "...nearly anyth #10 Nils Smeds (Smeds) Guest   Posts: n/a I guess that "...nearly anything in FOAM without looping" is supposed to be read as "...nearly anything in FOAM without explicit looping"

 March 3, 2005, 23:46 I am using the following code #11 Jarrod Sinclair (Sinclair) Guest   Posts: n/a I am using the following code to summarise the flux at each boundary patch, along with the net, after each timestep. { // loop through all boundary patches int patchItr; label patchIndex; scalar flux, netFlux = 0.0; forAll (mesh.boundaryMesh(), patchItr) { // get the patch index patchIndex = mesh.boundaryMesh()[patchItr].index(); // compute the mass flux and net flux flux = sum(phi.boundaryField()[patchIndex]); netFlux += flux; // compute and display the mass flux with patch name Info<< "Mass flux at " << mesh.boundaryMesh()[patchIndex].name() << " = " << flux << endl; } Info<< "Net mass flux = " << netFlux << endl << endl; } This works fine, but in parallel it can't resolve some of the boundary names that are not available on the master processor (giving a name of "procBoundary0to1" for example). How can the master processor get access to these names via the mesh.boundaryMesh()[patchIndex].name() statement? (I couldn't find a solution in the Doxygen notes). However, the flux appears to be calculated fine.

 March 4, 2005, 03:39 "procBoundary0to1" is the pro #12 Henry Weller (Henry) Guest   Posts: n/a "procBoundary0to1" is the processor-processor "boundary" patch which isn't a true boundary, it's the internal patch created by the decomposition and probably not important to your flux calculation. In order to sum over the "real" boundary patches in parallel you will need to test that the patches you are interested are on the processor and after the loop do a reduce operation to sum up the contributions from the processors to the boundary regions. Have a look through src/cfdTools/adjustPhi/adjustPhi.C which uses the logic you need to sum inlet and outlet fluxes to correct mass-balance for steady-state cases.

 March 7, 2005, 03:11 Thanks for that Henry, I have #13 Jarrod Sinclair (Sinclair) Guest   Posts: n/a Thanks for that Henry, I have made those modifications and am using the reduce command to sum the contributions from each processor. Is there a way to get a list of all boundaries on the original mesh while running in parallel? Or is this not possible due to running on decomposed grids in seperate processes? At the moment I am defining a patch list explicity, but I could potentially parse the original constant/polyMesh/boundary file to automate this (assuming the files are available during runtime).

 March 7, 2005, 04:48 You can send the information #14 Henry Weller (Henry) Guest   Posts: n/a You can send the information about the boundary patches to the master, sort the lists and send them back. However, the patch numbering on the processors will be different because each only hold the list of it's own patches so I am not sure how this complete list will help you.

 March 8, 2005, 22:02 I think I've got a solution th #15 New Member   Jarrod Sinclair Join Date: Mar 2009 Posts: 6 Rep Power: 9 I think I've got a solution that works in both serial and parallel, and uses the communcation and sorting you suggested. Firstly, you need to #include "SortableList.H" as one of the top level includes in a solver application. Then #include "buildGlobalBoundaryList.H" once "createMesh.H" and other initialisations have been run in the main function. (This only needs to be run once.) The file is attached below. buildGlobalBoundaryList.H The result of this is a wordList variable called "globalPatchNames" that is available on each processor, and contains a sorted list of all real patches in the model. To compute mass fluxes after each timestep, #include "computeMassFlux.H" (see attachment). computeMassFlux.H I've tested this in both serial and parallel, so hopefully it works in other people's cases. mm.abdollahzadeh likes this.

 March 9, 2005, 04:51 Thanks for these pieces of sam #16 Senior Member   Join Date: Mar 2009 Posts: 854 Rep Power: 14 Thanks for these pieces of sample code, it all seems fine except for one point, did you need to include Pstream.H? If so it should be included outside the "main" function at the top of the code although it is usually already included via some other includes.

 March 16, 2005, 20:07 You are right about Pstream.H #17 New Member   Jarrod Sinclair Join Date: Mar 2009 Posts: 6 Rep Power: 9 You are right about Pstream.H Henry, I've made some changes to the above code. buildGlobalBoundaryList.H is now a function that returns a wordList of patch names rather than inline code, so to avoid any variable redefinition problems. This also makes the headers a bit cleaner. However, when I perform a sort on a SortableList with OpenFOAM 1.1, I get a segmentation fault (this worked in 1.0.2). Do you know why this is so?

 March 17, 2005, 04:31 The difference between 1.0.2 a #18 Senior Member   Join Date: Mar 2009 Posts: 854 Rep Power: 14 The difference between 1.0.2 and 1.1 with respect to sortable list is that I have changed all calls to the C-library qsort and replaced them the equivalent STL sort calls which is certainly a more elegant method but we have had problems with it in the past. Could you please send us the code that produces the failure and we will look into it.

 March 17, 2005, 09:08 Sure, the new version of "buil #19 New Member   Jarrod Sinclair Join Date: Mar 2009 Posts: 6 Rep Power: 9 Sure, the new version of "buildGlobalBoundaryList.H" is attached below: buildGlobalBoundaryList.H This should be included outside the main function after top level includes like "fvCFD.H". It can then be called inside the main function like so, wordList globalBoundaryList; globalBoundaryList = buildGlobalBoundaryList(mesh); Then, mass fluxes can be calculated via the following include (similar to the previous example). computeMassFlux.H I've commented out the sort() operation in "buildGlobalBoundaryList.H" for now.

 March 17, 2005, 09:58 Hello Jarrod, we located fo #20 Super Moderator   Mattijs Janssens Join Date: Mar 2009 Posts: 1,416 Rep Power: 17 Hello Jarrod, we located found the problem in SortableList. It is because setSize when it grows the list does not adapt the list of indices. Until we fix it you can work around it by creating your list big enough and shrinking it (i.e. setSize with something smaller than the current size) Mattijs

 Thread Tools Display Modes Linear Mode

 Posting Rules You may not post new threads You may not post replies You may not post attachments You may not edit your posts BB code is On Smilies are On [IMG] code is On HTML code is OffTrackbacks are On Pingbacks are On Refbacks are On Forum Rules

 Similar Threads Thread Thread Starter Forum Replies Last Post Bloshchitsyn Vladimir CFX 0 November 18, 2007 23:19 aunola Main CFD Forum 0 May 8, 2006 07:54 CFDtoy Main CFD Forum 0 October 5, 2004 13:08 StarCD CD-adapco 2 May 24, 2004 14:33 R.Sripriya FLUENT 4 June 12, 2003 12:49

All times are GMT -4. The time now is 12:35.