# Porosity and the energy equation

 User Name Remember Me Password
 Register Blogs Members List Search Today's Posts Mark Forums Read

 January 20, 2020, 10:39 Porosity and the energy equation #1 Senior Member   Join Date: Dec 2019 Posts: 215 Rep Power: 7 Hi, I am trying to model a porous media using the darcy forchheimer model. I am using Openfoam 4 and the solver sonicFoam. For this purpose I have added a fvOptions dict: Code: ```porosity_wall { type explicitPorositySource; active true; explicitPorositySourceCoeffs { type DarcyForchheimer; selectionMode cellZone; cellZone porosity; // Specify the name of the cellZone DarcyForchheimerCoeffs { // Negative coeffs are multiplied by largest positive coeff, // taking the magnitude, e.g. for -1000, coeff = |1e7*-1000| = 1e10 d [0 -2 0 0 0 0 0] (1e10 1e10 1e10); f [0 -1 0 0 0 0 0] (0 0 0); coordinateSystem // Cartesian coordinates for the cellZone { x (1 0 0); y (0 1 0); #includeEtc "caseDicts/general/coordinateSystem/cartesianXY" } } } }``` It works fine concerning the pressure drop. The darcy forchheimer modell is an additional source term in the momentum equation. But the problem is, that the temperature seems not to be affected by the porosity. I have looked up EEqn.H for sonicFoam: Code: ```fvScalarMatrix EEqn ( fvm::ddt(rho, e) + fvm::div(phi, e) + fvc::ddt(rho, K) + fvc::div(phi, K) + fvc::div(fvc::absolute(phi/fvc::interpolate(rho), U), p, "div(phiv,p)") - fvm::laplacian(turbulence->alphaEff(), e) == fvOptions(rho, e) );``` It seems the solver considers fvOptions in the energy equation too. Unfortunately I cannot find how or whether there is an additional source term in the energy equation aswell. Does anyone know whether there is a source term in the energy equation? Kind regards, shock77 alainislas likes this.

 January 20, 2020, 17:07 #2 Senior Member   Lasse Brams Vinther Join Date: Oct 2015 Posts: 112 Rep Power: 10 Hello shock77, The fvOptions in EEqn, regards the potential sources, constraints or corrections you could add to your energy equation, however, the thermal part of the porous media is not described in openFoam. You can do this your self by developing the fvOptions library required, or implement it directly in the solver (this is what I have done previously). Lastly you could do a coded source that does it all, however it requires you develop the entire sink term yourself from a empirical formula or similarly. I have done something alike the last part as well, just a source term that depended on the geometrical location of the cell in the domain, volume of the cell, velocity in the cell and so on. I have attached an example of the code I have used, where a1 through a5 are just numerical constants from a fit. Code: ```energySource { type scalarCodedSource; //scalarSemiImplicitSource active true; name sourceTime; scalarCodedSourceCoeffs //scalarSemiImplicitSourceCoeffs S(x) = Su + Sp*x // q in [W]; or in [W/m³] if you use specific mode { selectionMode cellZone; cellZone porousity1; fields (h); fieldNames (h); name sourceTime; codeInclude #{ #}; codeCorrect #{ // Pout<< "**codeCorrect**" << endl; #}; codeAddSup #{ // const Time& time = mesh().time(); const scalarField& V = mesh_.V(); const vectorField& C = mesh_.C(); const volVectorField& U = mesh().lookupObject("U"); const volScalarField z = U.mesh().C() & vector(0,0,1); scalarField& hSource = eqn.source(); forAll(C, i) { hSource[i] -= a1*((a2- a3*exp(-a4*mag(a5-z[i])))*V[i]); } // Pout << "***codeAddSup***" << endl; #}; codeSetValue #{ // Pout<< "**codeSetValue**" << endl; #}; // Dummy entry. Make dependent on above to trigger recompilation code #{ \$codeInclude \$codeCorrect \$codeAddSup \$codeSetValue #}; } sourceTimeCoeffs { // Dummy entry } }``` Best regards, Lasse

 January 24, 2020, 09:53 #3 Member   Join Date: Jun 2017 Posts: 73 Rep Power: 8 Hi, could you explain how you derived your change in enthalpie due to porosity?

 January 24, 2020, 10:23 #4 Senior Member   Lasse Brams Vinther Join Date: Oct 2015 Posts: 112 Rep Power: 10 My approach is somewhat based upon the following source: https://www.slideshare.net/Cypiii/op...training-v51en With the heat transfer coefficient being based upon empirical formulas. I believe I found applicable formulas in Chemical reactor analysis and design by Froment, Bischoff. And then only executing the catalyst code if in the porous zone. Let me know if you have any additional questions.

 January 31, 2020, 09:00 #5 Senior Member   Join Date: Dec 2019 Posts: 215 Rep Power: 7 Dear Swagga5aur, thank you very much for your reply and sry for my very late response. So what you do is treating the heat transfer within the porous media when the gas and the solid have different temperatures right? What I am wondering about is, that I cant find any treatment of the temperature change due to the increased friction inside a porous material. It seems to be always negleted. Maybe I am mistaken and its not important. What I also cant find is how turbulence models are affected by porosity in OpenFOAM. Do you know how it works? Kind regards, shock77

 January 31, 2020, 15:15 #6 Senior Member   Lasse Brams Vinther Join Date: Oct 2015 Posts: 112 Rep Power: 10 Yes, heat transfer requires a difference in temperatures to exist. Do you suggest that due to the lower cross sectional flow area of the fluid a frictional heat generation occurs? To my understanding frictional heat generation are usually neglect-able compared to the convective heat transfer. I have mainly dealt with low velocities, however, perhaps some the of the following links may help. http://scientiairanica.sharif.edu/ar...f69f2ffd37.pdf https://www.cfd-online.com/Wiki/V2-f_models In general I would consider the resistance of the porous media outweighs the influence of the turbulent flow. https://vbn.aau.dk/ws/portalfiles/po..._equations.pdf Check page 242 for a description of the critical Reynolds number. Hope this helps. Best regards Lasse

February 4, 2020, 07:10
#7
Senior Member

Join Date: Dec 2019
Posts: 215
Rep Power: 7

Quote:
 Do you suggest that due to the lower cross sectional flow area of the fluid a frictional heat generation occurs?
Yes excactly. I have a high velocity flow and I am not sure whether the friction plays a role or not. I have read lots of papers and every author neglects the additional friction. Maybe you are right and its not very important.

Quote:
 In general I would consider the resistance of the porous media outweighs the influence of the turbulent flow.
I have read both. Either you neglect the turbulence, whenever you have a very fine porosity or you have to implement it into your turbulence model.

It is a difficult topic, at least for me. I am lacking experience here.

But thank you very much for your kind help! It helped me a lot.

Kind regards,
shock77

 March 7, 2020, 08:00 #8 New Member   Kunqing Jiang Join Date: Mar 2020 Posts: 2 Rep Power: 0 Dear Swagga5aur and shock77 I have met the same problem,and I benefited a lot from your discussion.But it's hard for me to do this "And then only executing the catalyst code if in the porous zone".Because I haven't learned openfoam long. Could you please tall me how to do that ,thanks.

 March 10, 2020, 09:15 #9 Senior Member   Lasse Brams Vinther Join Date: Oct 2015 Posts: 112 Rep Power: 10 Hello Kunqing, I currently don't have the time to make this for you, however, just a brief description of what I would do. I would add cells from a porosity zone as the following. Code: ``` const labelList& cells = mesh.cellZones() [mesh.cellZones().findZoneID("porosity")];``` And then add a for loop or while loop for the cells where the normal y equation is solved and then everywhere else the y equation is solved with removed reactions such that no reactions occur, but the mixing of the specie still does. Hope it helps.

 March 14, 2020, 09:33 #10 New Member   Kunqing Jiang Join Date: Mar 2020 Posts: 2 Rep Power: 0 Dear Swagga5aur, Thanks for your reply,it helps me a lot. I will have a try myself. Thank you very much indeed.

 March 28, 2020, 15:59 #11 Senior Member   Lasse Brams Vinther Join Date: Oct 2015 Posts: 112 Rep Power: 10 To anyone of interest here is an example of how to make a zone specific Y solver where you can remove for example the reaction term. Its in no way improved for speed just briefly tested. Code: ```tmp > mvConvection ( fv::convectionScheme::New ( mesh, fields, phi, mesh.divScheme("div(phi,Yi)") ) ); { reaction.correct(); Qdot = reaction.Qdot(); volScalarField Yt(0.0*Y[0]); const labelList& cellZ = mesh.cellZones() [mesh.cellZones().findZoneID("catalyst")]; const label zoneID = mesh.cellZones().findZoneID("catalyst"); if (zoneID != -1) { forAll(cellZ, i) { zF=1; } forAll(Y, i) { if (Y[i].name() != inertSpecie) { volScalarField& Yi = Y[i]; fvScalarMatrix YiEqn ( mvConvection->fvmDiv(phi, Yi) - fvm::laplacian(turb.muEff(), Yi) == reaction.R(Yi) + fvOptions(rho, Yi) ); YiEqn.relax(); YiEqn.relax(mesh.equationRelaxationFactor("Yi")); fvOptions.constrain(YiEqn); YiEqn.solve(mesh.solver("Yi")); fvOptions.correct(Yi); Yi.max(0.0); Yt += Yi; } } } forAll(Y, i) { if (Y[i].name() != inertSpecie) { volScalarField& Yi = Y[i]; fvScalarMatrix YiEqn ( mvConvection->fvmDiv(phi, Yi) - fvm::laplacian(turb.muEff(), Yi) == fvOptions(rho, Yi) ); YiEqn.relax(); YiEqn.relax(mesh.equationRelaxationFactor("Yi")); fvOptions.constrain(YiEqn); YiEqn.solve(mesh.solver("Yi")); fvOptions.correct(Yi); Yi.max(0.0); Yt += Yi; } } Y[inertIndex] = scalar(1) - Yt; Y[inertIndex].max(0.0); }``` Just remove the following code to make it compile able. Code: ``` forAll(cellZ, i) { zF=1; }``` If interested for a simple indicator of the functionality, leave this code and add the following to createFields directory to make a field zF being 1 if its the specified zone or 0 if outside it. Note this is made for a custom multiRegion solver but the main idea is the same,. Code: `PtrList zFReactor(reactorRegions.size());` Code: ``` Info<< " Adding to zFReactor\n" << endl; zFReactor.set ( i, new volScalarField ( IOobject ( "zF", runTime.timeName(), reactorRegions[i], IOobject::NO_READ, IOobject::AUTO_WRITE ), reactorRegions[i], dimensionedScalar("zF", dimless, scalar(0.0)) ) );``` Let me know if there is any questions or issues. Best regards, Lasse

 May 9, 2020, 15:52 #12 Senior Member   Lasse Brams Vinther Join Date: Oct 2015 Posts: 112 Rep Power: 10 Hello all, Turned out this doesn't work, however, the following does but is less elegant regarding skipping unnecessary steps. Added this to setRegionFields: Code: ``` word cellSetName = "catalyst"; label zoneID = mesh.cellZones().findZoneID(cellSetName); if (zoneID == -1) { FatalErrorIn("yourFunctionName") << "Cannot find cellZone " << cellSetName << endl << "Valid cellZones are " << mesh.cellZones().names() << exit(FatalError); } const labelList& cells = mesh.cellZones()[zoneID]; Info << "Cells in cellzone " << cellSetName << ":" << endl; forAll(cells, i) { const label cell = cells[i]; zF[cell]=1; }``` The following to createFields: Code: ``` Info<< " Adding to zF\n" << endl; zF.set ( i, new volScalarField ( IOobject ( "zF", runTime.timeName(), reactorRegions[i], IOobject::NO_READ, IOobject::AUTO_WRITE ), Regions[i], dimensionedScalar("zF", dimless, scalar(0.0)) ) );``` and then just uses as a multiplier to disable fx. reaction.R(Yi) call in YEqn to remove reactions in cells outside of the catalyst cellzone. Best regards, Lasse

 Thread Tools Search this Thread Search this Thread: Advanced Search 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 Off Pingbacks are On Refbacks are On Forum Rules

All times are GMT -4. The time now is 19:56.

 Contact Us - CFD Online - Privacy Statement - Top