
[Sponsors] 
interFoam with scalar transport in single phase 

LinkBack  Thread Tools  Search this Thread  Display Modes 
March 3, 2021, 09:50 
interFoam with scalar transport in single phase

#1 
Senior Member
Reviewer #2
Join Date: Jul 2015
Location: Gainesville,FL
Posts: 134
Rep Power: 9 
Hi everyone,
I am interested in simulating a scalar transport within the primary phase (i.e., water) in the interFoam. I follow the approach of phaseScalarTransport (https://github.com/OpenFOAM/OpenFOAM...larTransport.H) and successfully added a scalar transport to the primary phase with MULES compressed flux (i.e., alphaPhiUn). Code:
fvScalarMatrix YEqn ( fvm::ddt(Y) + fvm::div(alphaPhiUn, Y)  fvm::laplacian(Deff, Y) == alpha1*fvOptions(Y) ); YEqn.relax(); fvOptions.constrain(YEqn); YEqn.solve(); The code runs well and as expected the scalar is transported only within the water phase. When the water surface arises, the scalar will arise too. However, when the water surface falls, the part of the scalar will "stuck" above the surface (see attachment). I would like to have scalar fall as the water falls. Could you please shed some light on this issue? Thanks, Rdf 

March 5, 2021, 15:52 

#2 
Senior Member
Reviewer #2
Join Date: Jul 2015
Location: Gainesville,FL
Posts: 134
Rep Power: 9 
Okay. I think I figured it out. I will post my findings here for future reference. There are several post concerns about this issue and provide valuable information.
InterFoam: Add an equation to only one phase interFoam with transport: add a transport equation to only one phase Numerical problem in species Transport in InterFoam at the interface interFoam + a scalar transport equation New phaseScalarTransport function object for a single phase in interFoam In summary, the simplest way (at this moment) to simulate a scalar transport in a single phase (e.g., water) with interFoam is using the "phaseScalarTransport" utility from the OpenFOAM 7 or 8 (org version). The field alphaS.water = alpha.water * S.water should be used as the results, not the scalar itself (S.water). This is because the continuity equation is calculated for alphaS.water, as shown in the following code. If one looks at the S.water field, the scalar itself will be dispersed into the other phase. This is not a bug. I use the dambreak as a testing case and found the mass of alphaS.water is generally conserved with an error of less than 1% while the mass of the scalar itself (S.water) will have an error range from 100% to 400%. Code:
fvScalarMatrix fieldEqn ( fvm::ddt(alpha, s_) + fvm::div(alphaPhi, s_, divScheme)  fvm::laplacian ( fvc::interpolate(alpha)*fvc::interpolate(D), s_, laplacianScheme ) == fvOptions_(alpha, s_)  fvm::ddt(residualAlpha_, s_) + fvc::ddt(residualAlpha_, s_) ); fieldEqn.relax(relaxCoeff); fvOptions_.constrain(fieldEqn); fieldEqn.solve(schemesField_); In the ESIOpenFOAM, although there is a phase option in the "scalarTransport" utility and they even provide tutorials (e.g, angledDuct, waterChannel), it does not seems that phase constrained transport is implemented at all. In these tutorials, the scalars only diffusive but not convective, as observed in this post (Foam::error::printStack(Foam::Ostream&) with interFoam parallel). This is likely because the default option for phasePhiCompressed, alphaPhiUn is zero. Although alphaPhiUn is created in the interFoam, the interFoam solver does not pass a value to alphaPhiUn. This seems quite odd to me and please correct me if I am wrong. There is a method suggested for convecting the scalar in the "scalarTransport" utility of ESI OpenFOAM. Post (Foam::error::printStack(Foam::Ostream&) with interFoam parallel) suggest that defined the "alphaPhi0.water" for the "phasePhiCompressed". However, after some tests, this method will let some scalar stuck in the air. The "scalarTransport" utility of ESI OpenFOAM tries to conserve the mass for s.water instead of alphaS.water, as shown in the following code. Code:
tmp<surfaceScalarField> tTPhiUD; for (label i = 0; i <= nCorr_; i++) { fvScalarMatrix sEqn ( fvm::ddt(s) + fvm::div(limitedPhiAlpha, s, divScheme)  fvm::laplacian(D, s, laplacianScheme) == alpha*fvOptions_(s) ); sEqn.relax(relaxCoeff); fvOptions_.constrain(sEqn); sEqn.solve(mesh_.solverDict(schemesField_)); tTPhiUD = sEqn.flux(); } Rdf 

May 26, 2021, 08:34 

#3  
New Member
Sudarshan
Join Date: Sep 2020
Posts: 3
Rep Power: 4 
Hi,
I'm also trying to simulate a case (similar to sand in water + air) however I do not understand how to implement the phaseTransportFoam (in particular could you please explain what you meant by interfoam + phaseTransportFoam). I tried to implement this in the dambreak case. In the controlDict, file I included the following lines before the end. phaseScalarTransport1 { type phaseScalarTransport; libs ("libsolverFunctionObjects.so"); field s.water; p p_rgh; } phaseScalarTransport2 { type phaseScalarTransport; libs ("libsolverFunctionObjects.so"); field alpha.water * s.water; p p_rgh; } Then I run the command "interfoam" which ends well but when I enter the command "phaseTransportFoam" it is not found. Could you please tell the steps clearly on how to implement phaseTransportFoam. Thanks. Sudarshan. Quote:


July 13, 2021, 09:04 

#4 
New Member
Join Date: Feb 2016
Posts: 18
Rep Power: 8 
It looks like they are working on a fix for this in the ESIversion:
https://develop.openfoam.com/Develop...//issues/2053 

July 13, 2021, 09:33 

#5 
Senior Member
Reviewer #2
Join Date: Jul 2015
Location: Gainesville,FL
Posts: 134
Rep Power: 9 
Hi,
In the case of fox may encounter this issue. Here is one of my juggling by integrating the "phaseScalarTransport" (org version) with ESI. https://github.com/Rdfing/interAdsFoam Thanks, Rdf 

September 23, 2021, 04:12 

#6 
New Member
Shengjie Lu
Join Date: Sep 2020
Location: Nanjing,China
Posts: 12
Rep Power: 4 
Hi, randolph.
Your comment is of great help in imposing phaseScalarTransport.However, I'm still confused about its usage. For instance, what do you mean by 'The field alphaS.water = alpha.water * S.water should be used'? Should I specify a new variable named alphaS.water? I post part of my functions in controlDict file below(interFoam), could you help checking its correctness? Thank you in advance! functions { phaseScalarTransport1 { type phaseScalarTransport; libs ("libsolverFunctionObjects.so"); alpha alpha.water; field alpha.water * s.water; p p_rgh; writeAlphaField true; nCorr 1; resetOnStartUp false; writeControl outputTime; // write scalar field results writeInterval 20; // D 0.001; //difussion coefficient // nut nut; //name of field to use as diffusivity, default = 'none' // alphaDt 2; // schemesField U; //name of field to use when looking up schemes from fvSchemes bounded01 true; 

October 4, 2021, 21:52 

#7  
Senior Member
Reviewer #2
Join Date: Jul 2015
Location: Gainesville,FL
Posts: 134
Rep Power: 9 
Quote:
You can just use phaseScalarTransport as explained in the header file (https://cpp.openfoam.org/v7/classFoa...Transport.html). alphaS.water will be the field for postprocessing. Thanks, Rdf 

November 2, 2021, 23:25 

#8 
Member
Join Date: Aug 2018
Posts: 47
Rep Power: 6 
Thanks Randolph for explaining phaseScalarTransport.
I'm new to this function object. It looks your conclusion is alphaS.water is the scalar retained in water phase, and should be used for processing. I wonder why S.water is present? Is it involved in the scalar transport equation? What's it's physical meaning and under what situation it should be used in processing? Looking forward to your reply! 

November 17, 2021, 19:51 

#9  
Senior Member
Reviewer #2
Join Date: Jul 2015
Location: Gainesville,FL
Posts: 134
Rep Power: 9 
Quote:
This is not the exact reason why alphaS.water should be used. The alphaS.water should be used is because (i) "the continuity equation is calculated for alphaS.water" in the phaseScalarTransport as illustrated in my previous post and (ii) adhoc flux correction method implemented in phaseScalarTransport (OpenFOAM 7 or 8 [org version]). In the ESIOpenFOAM (e.g., v1912) formulation (i.e., phase option in the "scalarTransport" utility), the mass of the scalar will not be conserved even if alphaS.water is used for the processing. The way that phaseScalarTransport formulates the phaseconstrained scalar transport for interFoam is very similar to how the scalar transport handled in a typical dispersed multiphase flow simulation. The only added effort is on flux correction routine, which is also very similar to the (momentum) flux correction in a typical NS solver (if my memory serves). Thanks, Rdf 

January 10, 2022, 01:50 

#10 
Member
Thomas Sprich
Join Date: Mar 2015
Posts: 75
Rep Power: 9 
Hi Randolph,
I have read your posts with interest because I am attempting to model the scalar transport between phases, whereas you are interested in the scalar transport limited to a single phase. I thought I would share my method with you and am interested in your comments. I am sure my approach would work in your case as well. In my approach, I search through the mesh and compare the values of alpha. If alpha is between 0.05 and 0.95, then I set the diffusivity to a high value and if outside these bounds I set it to 0. This basically meant I have a diffusivity field for which I can solve the scalar transport equation. I actually took the scalar equation directly from scalarTransportFoam. To make this work for your case, all that would be required would be to set a high diffusivity for the phase that you want and set diffusivity everywhere else to zero. I created a test case based on the dam break case. I was particularly looking for the problem you identified in your first post, i.e. that the scalar was trapped in the one phase. For my approach, I do not seem to have this problem. I was then interested in your evaluation of the continuity of the scalar. How did you do this so I can verify for myself with my approach? Anyway, seeing as for interFoam, both phases share the velocity field and the scalar equation comes straight from scalarTransportFoam, I can't imagine there would be continuity issues. Maybe I'm being niave though... Let me know your thoughts. I've put the code as it is so far on github as well as a test case for you to evaluate. I compiled in on the 2106 ESI version. Here is the link to it: It is a work in progress, so it may not be as well documented as it could be. Let me know if you have any questions. My intention is to extend the code and make the diffusivity at the interface dependent on the turbulence as well as alpha and finally to have a maximum limit for the scalar in the phases. In this case, even if there is a difference in value between two cells, but one cell has reached its limit, no more scalar will be diffused across the faces. This is analogous to oxygen transfer of air to water. Air can have 23% O2, whereas the maximum concentration of O2 in water is about 8.5% I'd appreciate your (or anyone else's) thoughts and comments on my approach. Thanks, Thomas 

May 24, 2022, 15:29 

#11 
New Member
Join Date: May 2022
Posts: 2
Rep Power: 0 
Hi everyone! I'm on the brink of a nervous breakdown and about to torch my pc.
Hopefully some of you can give me a hand. I'm trying to determine the mixing time of a fermenter (stirred vessel with water and bubbles) To do so I was hoping to use phaseScalarTransport on twoPhaseEulerFoam. Inject some known quatity of scalar S (topoSet) and probe some points in paraview to do the post processing. The simulation runs normally but no s.water file is created in each of the timesteps. When looking at the log file I see that no calculations are done for the scalar. Here what I consider the relevant parts of my code: Regarding boundary conditions the scalar has zerogrdient everywhere since it's being generated by the injectionRateSuSp In controlDict I put the following (based on https://cpp.openfoam.org/v7/classFoa...ransport.html) Code:
functions { s { type phaseScalarTransport; libs ("libsolverFunctionObjects.so"); field s.water; alphaPhi alphaRhoPhi.water; rho thermo:rho.water; resetOnStartUp no; fvOptions { s { type scalarSemiImplicitSource; timeStart 0; duration 1; active true; selectionMode cellSet; cellSet tracer; volumeMode absolute; injectionRateSuSp { s.water (10 0); // kg/s } } } }; }; Additionally I added my fvScheme and my fvSolution Code:
fvSchemes div(phi,s.water) Gauss linearUpwind grad(s.water); fvSolution "s.*" { solver PBiCGStab; preconditioner DILU; tolerance 1e08; relTol 0; minIter 1; } 

June 20, 2022, 07:35 
phaseScalarTransport

#12  
Member
sadra mahmoudi
Join Date: Feb 2021
Location: Austria
Posts: 34
Rep Power: 3 
Hello dear JVNR,
I am also simulating a single bubble moving in a solution of sugar and water. In order to solve the transport of sugar in one phase, I am trying to use phaseScalarTransport. My problem is: I dont know how exactly should I implement phaseScalarTransport in the interFoam. I would be more than happy if you give me some tips for implementation. Best regards, Sadra Quote:


June 21, 2022, 12:26 

#13 
New Member
Join Date: Dec 2021
Posts: 19
Rep Power: 3 
Hey JVNR,
you probably already solved your problem, but I'll try to help you anyway. Make sure you have a s.water file in your /0. I'm pretty sure your fvSchemes is not correct and that you must add a div scheme for alphaPhi.water, s.water, e.g.: Code:
fvSchemes div(alphaPhi.water, s.water) Gauss linear When you start the simulation, abort it immediately and check the beginning of the log in your terminal. For me, the rest of the simulation is running well, too, even if phaseScalarTransport is crashing due to missing entries or something similar. However, the beginning of the log always said something about what's wrong with phaseScalarTransport.  I got a question of my own as well: Did anybody manage to simulate a constant injection of a concentration in the water phase through a boundary condition? I'm injecting a concetration into an already steady state from a boundary that is fully submerged in alpha.water=1. I've tried InletOutlet, but the simulation crashes after a few time steps the concentration explodes all over the domain. The mesh is rather coarse, so I can quickly test it on my laptop, but I'm not sure if that plays an important role here. Any help or tips would be appreciated. Best regards Finn Last edited by finn_amann; June 21, 2022 at 12:30. Reason: Added some details 

June 22, 2022, 05:12 
PhaseScalarTransport

#14  
Member
sadra mahmoudi
Join Date: Feb 2021
Location: Austria
Posts: 34
Rep Power: 3 
Dear Finn,
I would be thankful if you give me a hand. I am working on single bubble behavior in a solution of sugar and water with InterFoam. In order not to solve the concentration equation inside the bubble, I decided to use PhaseScalarTransport. Since a month ago, I am really struggling to implement it. I would be appreciative if you have a look at my case to see whether something is wrong. Best regards, Sadra Quote:


June 22, 2022, 11:26 

#15 
New Member
Join Date: Dec 2021
Posts: 19
Rep Power: 3 
Hey Sadra,
unfortunately I really don't know much about neither OpenFOAM nor phaseScalarTransport. But looking at your controlDict, I think you need something like Code:
functions { s { type phaseScalarTransport; functionObjectLibs ("libsolverFunctionObjects.so"); enabled true; writeControl outputTime; p p_rgh; field s.water; // name of your field, e.g. s.water // nCorr 1; D 0.001; schemeField U; log yes; } }; I hope that helps you. All my simulations are not working anyway and explode after a few time steps, so I'm not really much of a help to you, I guess. Best regards Finn 

June 24, 2022, 05:52 
phaseScalarTransport

#16  
Member
sadra mahmoudi
Join Date: Feb 2021
Location: Austria
Posts: 34
Rep Power: 3 
Hi Finn,
Thanks a lot for your guidance. I really appreciate it. The problem was resolved and I can simulate it properly. Although I still have some other problems. As you know, when someone use phaseScalarTransport, the final result (distribution of concentration in the primary phase) will be printed in a variable called "alphaS.water" in the paraview. In my base solver (InterFoam), I would like to make the density dependant on the concentration (alphaS.water). In fact, I would need something like: rho=alpha2*mu2 + alpha1*(0.0008*exp(0.0486*alphaS.water)) in which alpha2 is air and alpha1 is water. Now, the problem is, in spite of defining alphaS.water in the creatFiled of solver, the solver can not be compiled because it does not recongnise alphaS.water. I would be thankful if anybody give me a hand. Best regards, Sadra Quote:


July 12, 2022, 12:26 

#17 
New Member
Join Date: Dec 2021
Posts: 19
Rep Power: 3 
Hello people!
I wanted to ask what kind of schemes and solution settings you guys have been using for your simulations in this context so far? Have you been using rectangular or triangular cells and how succesful were youg with either one? 

July 13, 2022, 06:53 

#18 
Member
sadra mahmoudi
Join Date: Feb 2021
Location: Austria
Posts: 34
Rep Power: 3 
Hi finn_amann,
This is my fvScheme: ddtSchemes { default Euler; } gradSchemes { default Gauss linear; } divSchemes { //div(rhoPhi,U) Gauss limitedLinear 1; div(rhoPhi,U) Gauss linear 1; div(phi,alpha) Gauss linear 1; div(phirb,alpha) Gauss interfaceCompression; //Gauss linear; div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; div(phi,C) Gauss limitedLinear 1; div((phi*interpolate(alpha.air)),C) Gauss upwind 1; div(rhoPhi,C) Gauss upwind 1; div((rhoPhi*interpolate(alpha.air)),C) Gauss upwind; div(((100*dc)*grad(alpha.air))) Gauss linear 1; div((interpolate((1alpha.air))*phi),C) Gauss linear 1; div(alphaPhi,C) Gauss linear 1; div(alphaPhi,C_) Gauss linear 1; div(alphaPhi.water,s.water) Gauss linear 1; div(alphaPhi.water,alpha.water) Gauss linear 1; } laplacianSchemes { default Gauss linear corrected; } interpolationSchemes { default linear; } snGradSchemes { default corrected; } and this one is my fv solution: /** C++ **\  =========    \\ / F ield  OpenFOAM: The Open Source CFD Toolbox   \\ / O peration  Version: v2006   \\ / A nd  Website: www.openfoam.com   \\/ M anipulation   \**/ FoamFile { version 2.0; format ascii; class dictionary; location "system"; object fvSolution; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // solvers { "alpha.*.*" { interfaceMethod "isoAdvector"; isoFaceTol 1e6; surfCellTol 1e6; snapTol 1e12; nAlphaBounds 3; clip true; reconstructionScheme plicRDF; nAlphaCorr 2; nAlphaSubCycles 1; cAlpha 1; MULESCorr yes; nLimiterIter 3; solver smoothSolver; smoother symGaussSeidel; tolerance 1e8; relTol 0; } /////////////////// "s.*" { solver PBiCGStab; preconditioner DILU; tolerance 1e08; relTol 0; minIter 1; } /////////////////// "pcorr.*" { solver PCG; preconditioner { preconditioner GAMG; tolerance 1e05; relTol 0; smoother DICGaussSeidel; } tolerance 1e05; relTol 0; maxIter 100; } p_rgh { solver GAMG; tolerance 1e08; relTol 0.001; smoother DIC; } p_rghFinal { solver PCG; preconditioner { preconditioner GAMG; tolerance 1e07; relTol 0; nVcycles 2; smoother DICGaussSeidel; nPreSweeps 2; } tolerance 1e07; relTol 0; maxIter 20; } "(UkC)" { solver smoothSolver; smoother GaussSeidel; tolerance 1e06; relTol 0.1; nSweeps 1; } "(UkC)Final" { solver smoothSolver; smoother symGaussSeidel; tolerance 1e08; relTol 0; } } PIMPLE { momentumPredictor no; nCorrectors 2; nOuterCorrectors 1; nNonOrthogonalCorrectors 0; pRefCell 0; pRefValue 0; } // ************************************************** *********************** // 

Tags 
interfoam, scalar transport 
Thread Tools  Search this Thread 
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Temperature calculation from total enthalpy in OpenFOAM (XiFOAM)  sharifi  OpenFOAM Running, Solving & CFD  1  October 8, 2020 10:16 
Division by zero exception  loop over scalarField  Pat84  OpenFOAM Programming & Development  6  February 18, 2017 06:57 
problems concerning mass conservativity in bubbleFoam with custom scalar transport  cutter  OpenFOAM Programming & Development  3  February 10, 2015 05:25 
dieselFoam problem!! trying to introduce a new heat transfer model  vivek070176  OpenFOAM Programming & Development  10  December 24, 2014 00:48 
Variable diffusion coefficient for scalar transport in interFoam  amwitt  OpenFOAM Programming & Development  5  September 10, 2013 20:12 