
[Sponsors] 
July 13, 2010, 06:55 
Transport equation in chtMultiRegionSimpleFoam solver

#1 
Senior Member
Aurelien Thinat
Join Date: Jul 2010
Posts: 152
Rep Power: 6 
Hello everybody,
First of all, I am a new user of OpenFoam. I need to develop a scalar transport equation (which would be mass fraction) in the chtMultiRegionSimpleFoam solver in order to study a multispecies case. I'm using the reactingFoam solver as an example. Here are my problems : 1) I added the YEqn.H file in the "fluids" folder and I adapted the equation for a steady state solver and nonreactive species (no source term). I have a problem with the constructor : fv::convectionScheme<scalar>:New ( fluidRegions[i], // initally "mesh", fields, phi, mesh.divScheme("div(phi,Yi_h)") ) What is "fields" related to ? Should I modify "mesh.divScheme" into fluidRegions[i].divScheme ? 2) How can I initialize my scalars ? In the createFluidFields.H ? I didn't find the reader of the dictionnary "Ydefaut" in the reactingFoam source code. Thank you very much for your help. Aurélien 

July 15, 2010, 12:03 

#2 
Senior Member
Aurelien Thinat
Join Date: Jul 2010
Posts: 152
Rep Power: 6 
I have added 3 transport equations in the chtMultiRegionSimpleFoam solver.
I defined scalar[i]Fluid, i = 1, 2 or 3 as a ptrList<volScalarField> Then I define scalar[i]Fluid[j] , j= number of FluidRegion and I call solve { fvm::div(phi,scalar1Fluidi) } But then, when I launch the computation I got this message : "Time = 0 Solving for fluid region fluid DILUPBiCG: Solving for Ux, Initial residual = 1, Final residual = 0.00068726968, No Iterations 1 DILUPBiCG: Solving for Uy, Initial residual = 1, Final residual = 0.00057863562, No Iterations 1 DILUPBiCG: Solving for Uz, Initial residual = 1, Final residual = 0.00097519241, No Iterations 1 #0 Foam::error:rintStack(Foam::Ostream&) in "/home/trap/OpenFOAM/OpenFOAM1.7.0/lib/linux64GccDPOpt/libOpenFOAM.so" #1 Foam::sigFpe::sigFpeHandler(int) in "/home/trap/OpenFOAM/OpenFOAM1.7.0/lib/linux64GccDPOpt/libOpenFOAM.so" #2 ?? in "/lib64/libc.so.6" #3 Foam:ILUPreconditioner::calcReciprocalD(Foam::Fi eld<double>&, Foam::lduMatrix const&) in "/home/trap/OpenFOAM/OpenFOAM1.7.0/lib/linux64GccDPOpt/libOpenFOAM.so" #4 Foam:ILUPreconditioner:ILUPreconditioner(Foam: :lduMatrix::solver const&, Foam::dictionary const&) in "/home/trap/OpenFOAM/OpenFOAM1.7.0/lib/linux64GccDPOpt/libOpenFOAM.so" #5 Foam::lduMatrix:reconditioner::addasymMatrixCons tructorToTable<Foam:ILUPreconditioner>::New(Foam ::lduMatrix::solver const&, Foam::dictionary const&) in "/home/trap/OpenFOAM/OpenFOAM1.7.0/lib/linux64GccDPOpt/libOpenFOAM.so" #6 Foam::lduMatrix:reconditioner::New(Foam::lduMatr ix::solver const&, Foam::dictionary const&) in "/home/trap/OpenFOAM/OpenFOAM1.7.0/lib/linux64GccDPOpt/libOpenFOAM.so" #7 Foam::PBiCG::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const in "/home/trap/OpenFOAM/OpenFOAM1.7.0/lib/linux64GccDPOpt/libOpenFOAM.so" #8 Foam::fvMatrix<double>::solve(Foam::dictionary const&) in "/home/trap/OpenFOAM/OpenFOAM1.7.0/lib/linux64GccDPOpt/libfiniteVolume.so" #9 Foam::lduMatrix::solverPerformance Foam::solve<double>(Foam::tmp<Foam::fvMatrix<doubl e> > const&) in "/home/trap/OpenFOAM/OpenFOAM1.7.0/applications/bin/linux64GccDPOpt/mychtMultiRegionSimpleFoam" #10 main in "/home/trap/OpenFOAM/OpenFOAM1.7.0/applications/bin/linux64GccDPOpt/mychtMultiRegionSimpleFoam" #11 __libc_start_main in "/lib64/libc.so.6" #12 _start at /usr/src/packages/BUILD/glibc2.10.1/csu/../sysdeps/x86_64/elf/start.S:116 Floating point exception" Do you have any idea ? EDIT : I solve UEqn first, then YEqn1, YEqn2, YEqn3, hEqn... 

July 16, 2010, 12:09 

#3 
Senior Member
Aurelien Thinat
Join Date: Jul 2010
Posts: 152
Rep Power: 6 
I'm still on it. I have a problem to declare the diffusivity of the scalar.
In a multiregion solver, the createFluidFields.H file is a loop on each fluid region and the variables are defined as ptrList with a component for each region. I managed to declare correctly the scalar and to read the BCs' files. But when I'm trying to do quite the same for the diffusivity scalar, ie :  defining a ptrList<dimensionedScalar> Diffusivity(fluidRegions.size());  reading the transportProperties dictionnary, there is a problem of "redeclaration" of the diffusivity. Does anyone have an idea or any example ? 

July 17, 2010, 16:30 

#4 
Senior Member
David Gaden
Join Date: Apr 2009
Location: Winnipeg, Canada
Posts: 396
Rep Power: 11 
So your problem is now with defining a ptrList < dimensionedScalar > ?
Not sure why it's complaining about redeclaration... are you using the ptrList object correctly? I wrote a solver that did something similar. My scalars are not volume fractions, though... and I called it "gamma" instead of "diffusivity". Heres the code: Code:
// scalars is a subDictionary containing admScalar data dictionary admScalarsDict(admDict.subDict("scalars")); wordList admScalarKeywords(admScalarsDict.toc()); label nScalars(admScalarKeywords.size()); if (admDebug > 0) { Info << "Found " << nScalars << " scalars listed under keywords:\n" << endl; for (label i=0; i<(nScalars1); i++) { Info << "[" << admScalarKeywords[i] << "], "; } Info << "[" << admScalarKeywords[nScalars1] << "]\n" << endl; } PtrList<volScalarField> admVar(nScalars); PtrList<dimensionedScalar> gamma(nScalars); admScalar tempScalarID; forAll(admScalarKeywords, i) { tempScalarID = admScalarsDict.lookup(admScalarKeywords[i]); Info << "Reading field named " << tempScalarID.name(); Info << " from keyword " << admScalarKeywords[i] << " in admDict. \n" << endl; gamma.set ( i, new dimensionedScalar ( "gamma(" + tempScalarID.name() + ")", dimensionSet(0, 2, 1, 0, 0), tempScalarID.gamma() ) ); admVar.set ( i, new volScalarField ( IOobject ( tempScalarID.name(), runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE ), mesh ) ); } 

July 17, 2010, 20:55 

#5 
Member
Paulo Bufacchi
Join Date: Jun 2010
Posts: 50
Rep Power: 6 
Aurélien,
Are you working on diffusion of species from one region to another? I'm not sure I understood your model correctly. Paulo 

July 17, 2010, 22:12 

#6 
Senior Member
David Gaden
Join Date: Apr 2009
Location: Winnipeg, Canada
Posts: 396
Rep Power: 11 
This section of code is only in the initial setup (ie createFields). Since you were running into duplicate definition problems, I suspected you were using the ptrList wrong. The code is an example for that. As for the model, these are scalars that are used in a parallel biochemical model... so they are being transported by convection and diffusion, but they are not volume fractions because the biochemical model isn't mass preserving.
In other words, as far as OpenFOAM is concerned, these are just scalars. 

July 19, 2010, 03:40 

#7 
Senior Member
Aurelien Thinat
Join Date: Jul 2010
Posts: 152
Rep Power: 6 
Thank you Marupio. It was a problem with the prtList. It seems it's working now, I just have a problem of dimension conflict.
Edit : I might need your help again. The equation I'm solving is based on the scalarTransport solver ie : { fvm::div(phi,scalar1Fluid)  fvm::laplacian(Dscalar1Fluid,scalar1Fluid) } but the equation I need to solve is : div(U*c)  laplacian(Dc * c) = 0 with Dc = Dscalar1Fluid, c = scalar1Fluid, U = velocity. I saw in the user guide that phi = rho * U, The error I get is : incompatible dimensions for operation [scalar1Fluid[1 3 1 0 0 0 0] ]  [scalar1Fluid[0 0 1 0 0 0 0] ] It confirms my problem with phi. So I tried to solve "fvm::div(U,scalar1Fluid)" which is not working of course : "no matching function for call to 'div(Foam::GeometricField<Foam::Vector<double>,Foa m::fvPatchField, Foam::volMesh>&,...)" Aurélien Last edited by Aurelien Thinat; July 19, 2010 at 05:25. 

July 19, 2010, 10:19 

#8 
Senior Member
David Gaden
Join Date: Apr 2009
Location: Winnipeg, Canada
Posts: 396
Rep Power: 11 
phi is actually velocity flux interpolated to the volume faces, and appears frequently in the discretization process.
All the incompressible solvers I've looked at use transport equations that are divided by density... this changes the units you'd expect to see. What dimensions are you giving Dscalar1Fluid? In my solver, it turns out to be [0 2 1 0 0 0 0]. 

July 19, 2010, 10:28 

#9 
Senior Member
Aurelien Thinat
Join Date: Jul 2010
Posts: 152
Rep Power: 6 
I had a look to the createFluidFields and phi is defined as :
phiFluid.set ( i, new surfaceScalarField ( "phi", runTime.timeName(), fluidRegions[i], IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE ), linearInterpolate(rhoFluid[i]*UFluid[i]) & fluidRegions[i].Sf() ) ); Well, I'm adding rho in my transport equation. Last edited by Aurelien Thinat; July 19, 2010 at 11:01. 

July 19, 2010, 11:07 

#10 
Senior Member
David Gaden
Join Date: Apr 2009
Location: Winnipeg, Canada
Posts: 396
Rep Power: 11 
Ah, so phi actually matches convention for once. I'm used to working with incompressible. You don't want to define another phi. I think the problem is simply a units issue with your diffusion coefficient. The code I gave you (above) was for incompressible... so maybe your diffusion coefficient needs to be multiplied by density.
It's something simple like that, I'm sure. Take a moment to work out the units on paper  make sure you are defining the diffusion coefficient properly. (Boussinesq doesn't use density) 

July 19, 2010, 11:21 

#11 
Senior Member
Aurelien Thinat
Join Date: Jul 2010
Posts: 152
Rep Power: 6 
I have added the density inside all my term of the equation. It still computes and the scalar is still unbounded.
The Y1Eqn.H files : { for (label i=0; i>scalar1Fluid.size(),i++) { solve ( fvm::div(phi,scalar1Fluid[i])  fvm::laplacian(rhoFluid[i]*Dscalar1Fluid[i], scalar1Fluid[i]) ); runTime.write(); } } Do I have to add a test to detect the max between scalar1Fluid and 0/1 ? 

July 19, 2010, 11:31 

#12 
Senior Member
David Gaden
Join Date: Apr 2009
Location: Winnipeg, Canada
Posts: 396
Rep Power: 11 
You could try it. I'm not entirely sure. If you want to bound the solution, try using:
Code:
{ for (label i=0; i>scalar1Fluid.size(),i++) { solve ( fvm::div(phi,scalar1Fluid[i])  fvm::laplacian(rhoFluid[i]*Dscalar1Fluid[i], scalar1Fluid[i]) ); scalar1Fluid.max(0.0); scalar1Fluid.min(1.0); runTime.write(); } } 

July 20, 2010, 06:57 

#13 
Senior Member
Aurelien Thinat
Join Date: Jul 2010
Posts: 152
Rep Power: 6 
Hi,
My scalars transport equations are now working well. I didn't add the max and min function. Now I will have a look at the rho/Cp/lambda/mu state equations. Thank you again. Aurélien 

July 21, 2010, 06:53 

#14 
Senior Member
Aurelien Thinat
Join Date: Jul 2010
Posts: 152
Rep Power: 6 
Hello everybody,
I need some help to clear my mind. My purpose is to adapt the chtMultiRegionSimpleFoam solver for multispecies case.  So I started by adding 3 scalar transport equations, which will be use for the mass fraction of each species. The next steps are :  Adapt the EoS for the density and thermal diffusivity,  Add the condition : Sum(Yi) = 1. Am I right ? 

December 4, 2010, 13:09 

#15 
New Member
Join Date: Nov 2010
Posts: 10
Rep Power: 6 
dear Aurélien
i have a problem and i think you can help me. my problem is mass transfer from solid to liquid with electrochemical method and i want simulate this problem with openfoam and i cant find solver for it. how resolvent your problem? use a solver or write a code for it. if you write a code : can i have your code? very very thank in advance 

April 10, 2012, 14:51 

#16 
Senior Member
mahdi abdollahzadeh
Join Date: Mar 2011
Location: Covilha,Portugal
Posts: 149
Rep Power: 5 
Hi all,
I am also want to add a scaler equation to Chtmultiregionfoam. I have two question: 1 How is it possible to add a dimenedScaler in creatsfeild which be readed from a dictionary. ( for example the diffissity ) 2is the boundary condition for conjucate heat transfer is general enough to be used for other scalar equations with out any modification? Best Mahdi 

April 11, 2012, 02:17 

#17  
Senior Member

Quote:
how ever i guess you want to make a variable diffusivity! thats why you want to use setField to assign it values! so first create a volScalarField variable with the name diffusivity in creatFields.H in your solver, then you should compile your modified solver with wmake command (read user guide i think chapter 4 how to compile a solver) then make a diffusivity file in folder 0 , here in this file you can assign uniform or non uniform value ( with setField) for diffusivity Quote:


April 11, 2012, 04:51 

#18 
Senior Member
mahdi abdollahzadeh
Join Date: Mar 2011
Location: Covilha,Portugal
Posts: 149
Rep Power: 5 
Dear Nima
Thanks for your kind reply. I am a little bit familer with icoFoam. in icoFoam I can define a dimeionscaler and read it from a dictionary. the thing is that I want to solve another equations like electric field in both solid and fluid zones. so i want to difine electric diffusitivity in a dictionary and read it seprately for both solid and fluid zone. my second question is related to this fact that my electricpotential feild should be solved in both solid and fluid zones. I was thinking that can i use the turbulentTemperatureCoupledBaffleMixed boundary condition for another electricfield??? there is two option in that boundary ? K and KName which i dont know if can i change them to make this boundary applicable for electricfield? Best Mehdi 

April 11, 2012, 06:36 

#19 
Senior Member

so you should create your own boundary condition, im not familiar with electricpotential feild, but for first step to implement your bounday you may want to look
compressible::turbulentTemperatureCoupledBaffleMix ed to see how it is implemented 

April 11, 2012, 07:35 

#20 
Senior Member
mahdi abdollahzadeh
Join Date: Mar 2011
Location: Covilha,Portugal
Posts: 149
Rep Power: 5 
yes, I an on it, but a guy from chalmers where doing the same with privious versions and claim that its enough general!


Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Species transport equation  Bunmi  FLUENT  1  December 3, 2013 16:39 
UDF for concentration transport equation  heri  FLUENT  1  August 15, 2011 12:23 
Solver Equations  Carola  CFX  9  August 12, 2003 09:27 
compressible two phase flow in CFX4.4  youngan  CFX  0  July 2, 2003 00:32 