Multi species mass transport library [update]
 

I suggest also putting a link on the wiki: http://openfoamwiki.net/index.php/Main_ContribOther

Dear Anton,

I will follow your suggestion soon.

Many thanks

Valerio

Library conversion to be compatible with OF-12.1.1
 

Hi All,
that are really good news, but that will be also nice getting able to run the library under OF-2.1.1.

Many will be happy to read the programmer's suggestions regarding this.

BR

Im in the process of updating your library to 2.1.x and I have a small question.

For instance, the correct routine in Fick.C show that there is a slight chance of
making an easy mistake, (without knowing probably)

the inert species mass fraction is grabbed from thermo like this.
// Calculate inert species
volScalarField& yInert = thermo_.composition().Y()[inertIndex_];

which is set to 1 - sum of the first N-1 species, like this
yInert = 1 - yt;

since yt is just calculated in the loop
forAll(this->D_, i)

where D is
D_.setSize(species().size()-1);

and it gets the name
"D_" + species()[i],

so D is 1 less than the number of species

Now...if the inert specie is not placed in the last position in the species list,
the inertIndex will point to the wrong specie and it will use the inertspecie to calculate yt...
so basically everything will be a mess.

why dont you do it like in the YEqn.H where you check the name of the species to see if its the same as
the inert specie name
if (Y[i].name() != inertSpecie)
or
if (i != inertIndex_)

to prevent people from making this mistake and allow the inert specie to be placed anywhere in the species list.

just a thought,
N

Dear Niklas,

you are completely right! Actually the library works only if the inert species occupies the last position. But (and I'm sorry for this) I forgot to mention it.

This implementation is a kind of "misprint" of the firsts code versions . This bug fix was on to my to-do list, but I completely forgot it. I will add it again in my new to-do list!

Many thanks to have discovered (and shared) this problem. Surely it will be useful for all people that want to experiment the library.

Thanks again

Valerio

Dear Niklas, I canot open the pdf file after downloaded , can you check the file or send me a copy, many thanks!
my email is suningchen@gmail.com

Sorry, but I cannot read the pdf either. It says its damaged.

Dear Niklas,
Did you manage to update Novaresio's library into OF 2.1.x? I'd be really interested on it. Does it use cantera for the transport properties?
Regards
Jose

Dear Valerio,
Could you kindly re-upload the pdf file you posted on the 1st post to this thread? The original seems to be corrupt.
Thank you!
Regards
Jose

Greetings to all!

Well, this was funny... guess what: the 1st PDF file is in fact a DVI file!!!

I opened it easily with Ubuntu's default "Document Viewer" from the Evince authors and was finding it strange that everyone wasn't able to open it. So I then tried with Gimp and it did complain about the bad PDF file.
Tried with pdfedit to open it, same problem. pdf2ps, same deal.

So I then did what I usually do: open the file with a text editor (or hex editor, if necessary) and saw several bytes then something about "TeX" and more bytes... so I then tried dvipdf et voilà, there you have it, the real PDF :D :cool:

Have fun!
Bruno

LOL... funny

Thanks Bruno! Looks much better now!

Hi Again everyone!
I managed to compile Novaresio's library and the modifiedReactingFoam solver.
Now I want to set up a new simulation with Smooke's chemistry scheme but the library needs addtional info in the transportDictionary, namely the collision diameter and the characteristic Lennard-Jones energy of each involving species, to compute the binary diffusion coefficient.

I found in Novaresio's paper the some of these data is tabulated in "Hirshfelder et al., Molecular Theory of Gases and Liquids". However, neither the values he has in example case (supplied with the library) correspond to the ones in the former reference (despeite units), nor I am able to find data for all the species, namely, OH, H2O2, H, O, HCO, etc.

Does anybody have any idea where can I find this data?
Novaresio, could you enlighten me on this? Where did you get your data?

Regards everyone!

I've been monitoring this thread and am quite interested in experimenting with it. Jose, other than Niklas' comment about the inert species, were there any other fixes required to get it running in 2.1.x? From your earlier comments it appears you were gunning for that.

Thanks!

Well, I actually managed to compile it under OF-1.6-ext and not OF-2.1.x. My bad!

If you have any news on this please report them. I am very much interested on it. Anyways, if I have the time and patience, I will try to do it myself.

Regards

No worries. I'm sure Niklas would have mentioned something else had he come across it. I've done work with both versions and have some experience adapting between the two, so maybe I'll also tackle it if I find the time (been working on polyTopoChangers in 21x and that's been taking a LOT of time).

Cheers

Hi Again

I am getting a strange error with modifiedReactingFoam solver:

If I try to use a Chemkin file as input for chemical scheme I get this error:

I wonder what this is...
You can replicate this (i guess) by simply replacing the chemkin file (chem.inp) used in the dieselFoam tutorial distributed with OF by the chemistry used in the modifiedReactingFoam test case. Everything else is untouched.

Cheers

I think you also need the therm.dat file to get the thermophysical data for C7H16. chem.inp and therm.dat are tied closely together. Check to see if C7H16 is missing from your therm.dat.

Nop, C7H16 is there.
That's why it's really strange. I have no idea why that table has those species loaded. H2, CH4, OH, H2O2 are not even in the reaction file. So it makes no sense that they are in the table. Did you try to replicate the set up?

I tried to use reactingFoam instead of modifiedReactingFoam and I didnt got the error. So I guess the problem is with the multispecies mass diffusion library. Maybe this table is hardcoded with this species. I will try to investigate the code.

cheers

Hi!
I found a data file with all this info. For those interested, you can find it here:

http://diesel.me.berkeley.edu/~gri_m...12/text12.html

It has all the data necessary for methane combustion, using the grimech chemical mechanism. It is used for the sandia flame code.

Cheers
Jose

Ok.. I was able to track down this bug and came up with a simple fix.

So, I found out the code was crashing as it ran lines 65 and 66 in file Wilke.C:

I guess the problem is with function "molecularWeights" and so I replaced those lines with the following to avoid the fatal error I was getting so far.

const scalar& W1 = readScalar(dic.subDict("molarWeight").lookup(name1 ));
This will look for the molecular weights of each species in a subdictionary in the constant/transportProperties file.

This worked for me. If you're using other Binary diffusivity model than Wilke, you should do the same procedure in the respective .C file (i. e. ChapmanEnskog.C or Fuller.C).

I haven't look inside the molecularWeights.C file, where I think the prob is, to find out what is wrong with it and fix it. However, I will had this task to my todo list and, whenever I´ll have the time, I will try to have a look on it.

Hope I could help.

Cheers

Jose

It seems a shame to have to specify the molecular weights in a dictionary when we should already have the data available through the thermo model. Can we (in general) call , where specI is the label of the specie in questions (we'll have to do a search for it every time). I've tried to compile this but I'm told that composition (which is a basic) has no member W, when it does:

http://foam.sourceforge.net/docs/cpp...87805d01d0cc28

I'm sure I'm missing something about pointers and casting, but I'm not seeing it. Anyone have ideas?

Greetings to all!
From the link, it says it's a pure virtual method. This means that if the method is not reimplemented, then it's as good as non-existent!

Now, I haven't checked/looked at the code, but I'm assuming you can't simply use the basic thermo dictionary and will have to use the complex one, namely one of those long ones with the various characteristics for the fluid/gas...

Best regards,
Bruno

Ah, so I need to go deeper in the templatization. Following the inheritance diagram, multiComponentMixture<thermoType> implements the W() function. In createFields for the modified reactingFoam, which I think would be the general purpose type. The modifiedReactingFoam test case uses the thermo model.

Is there any way to explicitly to call the version of the function we want, maybe calling it straight from the psiChemistry object?

Perhaps something similar to:
Not very pretty, but it's the most direct thing I can think of, although this might not work exactly this way...

Additionally, the instantiated class might not like the dictionary file if said dictionary is using any other terms from the ones above :(

Bummer. So there isn't a way to call the implemented version of a pure virtual function of a template class without knowing what the implementation will be?

It shouldn't be necessary to bother with the details in the solver itself. Other solvers should show that it's as simple as using the correct thermo dictionary options for what you need to solve. What might happen is that some thermo options disregard completely molecular weight and such :rolleyes:

For example, buoyantSimpleFoam:
Seems simple enough here....

Tutorial "heatTransfer/buoyantSimpleFoam/hotRoom":
Weight should be easily accessible in this case. But it needs one big template hierarchy in there for it to work out of the box.

Thank you for the explanation. It just seems a bit off to me that a class named basicMultiComponentMixture would have an access function for the molecular weight and not be able to use it.

Doing a quick comparison between versions (as I'm switching between 1.6-ext and 2.1.x until the polyTopoChange/polyTopoChanger mechanisms are fixed), it seems that the function is available in 2.1.x (to accomodate the SLGThermo class I'll bet) but not in 1.6.ext.

So if you are using this with 2.1.x, you can get the molecular weights using

. I've compiled the code and it seems to work as desired. as previously mentioned. In 1.6-ext it appears we are stuck with the dictionary.

Hi, I'm interested in using this model for OpenFoam2.1.1. Did you guys happen to make it compatible with this version?
I tried to compile it using OF2.1.1 but got this error:
diffusivityModel/diffusivityModel.C:104:69: error: ‘const class Foam::porousZone’ has no member named ‘zoneId’
It might be because the zoneId is named something else in Of2.1.1! Does anybody have any clue?

Thanks

After looking at the file porousZone.H, I guess I should change "zoneId" to "zoneIds" which I did and now I get the below error:

diffusivityModel/diffusivityModel.C:104:76: error: no match for ‘operator[]’ in ‘(&(& T)->Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>::<anonymous>.Foam::DimensionedField <Type, GeoMesh>::mesh [with Type = double, GeoMesh = Foam::volMesh, typename GeoMesh::Mesh = Foam::fvMesh]())->Foam::fvMesh::<anonymous>.Foam::polyMesh::cellZon es()[(& pZones)->Foam::PorousZones<Foam::porousZone>::<anonymous>. Foam::IOPtrList<Foam::porousZone>::<anonymous>.Foa m::PtrList<T>::operator[] [with T = Foam::porousZone, Foam::label = int](zoneI).Foam::porousZone::zoneIds]’


The error refers to the part of the code which looks like this:

forAll(pZones, zoneI)
{
const scalar& porosity = pZones[zoneI].porosity();
const labelList& cells = T.mesh().cellZones()[pZones[zoneI].zoneIds()];
forAll(cells, cellI)
{
eps_[cells[cellI]] = porosity;
}
}

Can anybody help how I can resolve this error?
Thanks

Hi mhsn,

you have do add a loop over all porosity Zones.

forAll(pZones, zoneI)
{
const labelList& cellZoneIds = pZones[zoneI].zoneIds();
const scalar& porosity = pZones[zoneI].porosity();
forAll(cellZoneIds, zoneJ)
{
const labelList& cells = T.mesh().cellZones()[cellZoneIds[zoneJ]];
forAll(cells, cellI)
{
eps_[cells[cellI]] = porosity;
}
}
}

BR

Thanks, it compiled :) Now I should check what error I will see when compiling the other library "multiSpeciesTransportModels".

Well, now I'm trying to compile "multiSpeciesTransportModels" and after resolving some errors regarding universal gas constant "R", I see some errors when it tries to read submodels.
Here is what I get:
lnInclude/Fick.C:149:19: error: ‘const class Foam::fvMesh’ has no member named ‘relax’
lnInclude/Fick.C:167:19: error: ‘const class Foam::fvMesh’ has no member named ‘relax’
lnInclude/Fick.C:169:28: error: ‘const class Foam::fvMesh’ has no member named ‘relaxationFactor’

and here is what the code looks like at those lines:

iline149:
f (mesh_.relax("Yi"))
{
yi.storePrevIter();
}


line167:
if (mesh_.relax("Yi"))
{
yi.relax(mesh_.relaxationFactor("Yi"));
}


I guess from of1.6 to of2.1.1 that relax may be defined differently! IS that right? anyone knows how this can be resolved?
Tetraeder, can you help?

Hi mhsn,

with


if (mesh_.relaxField("Yi"))
{
yi.storePrevIter();
}

.....

if (mesh_.relaxField("Yi"))
{
yi.relax(mesh_.fieldRelaxationFactor("Yi"));
}

it should compile in OF 2.1.x

BR

Yes, it compiled. Thanks a lot

Just to learn more for future, how did you find the equivalent member functions in of2.1? I looked into everything and couldn't find them and it was because I didn't know where and what to look for!

Another problem, now that I have both libraries compiled, I went to to $FOAM_RUN/../applications and typed wmake modifiedReactingFoam and got this error:
modifiedReactingFoam.C:41:40: fatal error: multiSpeciesTransportModel.H: No such file or directory

Why can't it find the .H file although I have the my Make/options looks like it addresses that. Here is my Make/options:
EXE_INC = \
-I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
-I$(LIB_SRC)/ODE/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(FOAM_USER_LIBBIN)/../diffusivityModels/lnInclude \
-I$(FOAM_USER_LIBBIN)/../multiSpeciesTransportModels/lnInclude


EXE_LIBS = \
-lcompressibleTurbulenceModel \
-lcompressibleRASModels \
-lcompressibleLESModels \
-lreactionThermophysicalModels \
-lbasicThermophysicalModels \
-lchemistryModel \
-lODE \
-lfiniteVolume \
-llduSolvers \
-L$(FOAM_USER_LIBBIN) \
-lmultiSpeciesTransportModels

Will be appreciated if you give me a hint on that too! :)
Thanks again

I had to change two the below lines:
-I$(FOAM_USER_LIBBIN)/../diffusivityModels/lnInclude \
-I$(FOAM_USER_LIBBIN)/../multiSpeciesTransportModels/lnInclude

to:
-I$(FOAM_RUN)l/../lib/diffusivityModels/lnInclude \
-I$(FOAM_RUN)/../lib/multiSpeciesTransportModels/lnInclude

I had to add some more lines as well, so the error resolved.
Now I face this error which I don't understand what it is about:

/usr/bin/ld: cannot find -llduSolvers

Do you know what it relates to? let me know if you needed more information.

Thanks

Hey guys, thanks to your help, I could make modifiedReactingFoam work with openFoam2.1.1.
you may have to change some parts to be adapted to your openFoam especially in modifiedReactingFoam/Make/options and modifiedReactingFoam/Make/files

let me know if you needed that. I couldn't upload them regarding the maximum limit for the file size.

fixedFlux boundary
 

Dear novyno, dear community,


Great work and thanks for this library!

In your publication in Computer Physics Communications, you mention a generic boundary condition for species mass flux, called fixedFlux.
Could you please post it here? I would like to use it for part of my bachelor thesis on silicium wafer processing.

I am trying to implement a simple surface reaction boundary with a sticking coefficient 0..1 (only valid for Ficks law, similar to the boundary fixedGradient):

// Get Species name
const word specie_(dimensionedInternalField().name());

// Get p, T, rho, moleFrac, massFrac, Alpha,diffAlpha, molarWeight
const fvPatchField<scalar>& p = patch().lookupPatchField<volScalarField, scalar>("p");
const fvPatchField<scalar>& T = patch().lookupPatchField<volScalarField, scalar>("T");
const fvPatchField<scalar>& rho = patch().lookupPatchField<volScalarField, scalar>("rho");
// const fvPatchField<scalar>& massFrac = patch().lookupPatchField<volScalarField, scalar>(specie_); // Mass fraction
const fvPatchField<scalar>& moleFrac = patch().lookupPatchField<volScalarField, scalar>("x_"+specie_); // Mole fraction
const fvPatchField<scalar>& diffAlpha = patch().lookupPatchField<volScalarField, scalar>("D_"+specie_); // Diffusion constant
const scalar molarWeight(molecularWeights[specie_]); // kg/mol

// Constants
const scalar R(8.3144621); // J/(K*mol) = kg^2*m^2*s^-2*kgmol^-1
const scalar twoPi(6.283185307179586);

gradient() = -stickingCoeff_*p*moleFrac / (sqrt(twoPi*molarWeight*R*T) * (rho*diffAlpha));

Unfortunately, in every second iteration, this formula produces a gradient only containing zeros "gradient() 25{-0}" which probably means, moleFrac is exactly zero (unphysical I think).
Can somebody help me out please?

Thanks a lot,
OpenFoamStudent