CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM Programming & Development (http://www.cfd-online.com/Forums/openfoam-programming-development/)
-   -   Adding new Transport properties model (http://www.cfd-online.com/Forums/openfoam-programming-development/118841-adding-new-transport-properties-model.html)

dvcauwe June 5, 2013 10:20

Adding new Transport properties model
 
Hello Foamers!

I am currently trying to implement a method for calculating viscosities and thermal conductivities based on the Chapman-Enskog assumption (~Kinetic Theory method in Fluent). Starting from sutherlandTransport, this is very simple but I still have 2 questions:


1) Is it necessary to rebuild psiThermo, psiReactionThermo and psiChemistryModel equivalent classes that can handle this method of calculating transport properties or can I just make this method selectable somewhere? Since obviously at the moment the kinetic theory transport + janaf thermo is not considered a valid combination.

2) In case there's no way to make 1) work, I already went ahead and tried to build "pyThermo", ... as exact copies of the psi-equivalents, while making the necessary changes in thermoPhysicsTypes.H, PyThermos.C, pyReactionThermos.C, ... This new thermo lib compiles without a problem, as well as the pyChemistryModel lib but now comes the strange part: when I try to recompile my custom solver to use this pyChemistryModel instead of psyChemistryModel (everything exact copy, only adjusted so it can use kineticTheoryTransport), I get an undefined reference error:

Code:

Make/linux64GccDPOpt/QSSAFoam.o: In function `main':
QSSAFoam.C:(.text.startup+0x31c): undefined reference to `Foam::pyChemistryModel::New(Foam::fvMesh const&)'
collect2: ld returned 1 exit status


My solver Make/options file looks like this, so I assume the libraries are correctly linked?
Code:

EXE_INC = \
    -I$(WM_PROJECT_USER_DIR)/applications/lib/pyrolysisModel/lnInclude \
    -I$(WM_PROJECT_USER_DIR)/applications/lib/pyThermo/lnInclude \
    -I$(LIB_SRC)/finiteVolume/lnInclude \
    -I$(LIB_SRC)/fvOptions/lnInclude \
    -I$(LIB_SRC)/meshTools/lnInclude \
    -I$(LIB_SRC)/sampling/lnInclude \
    -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)/ODE/lnInclude \

EXE_LIBS = \
    -L$(FOAM_USER_LIBBIN)/libpyrolysisModel.so \
    -L$(FOAM_USER_LIBBIN)/libpyThermo.so \
    -lfiniteVolume \
    -lfvOptions \
    -lmeshTools \
    -lsampling \
    -lcompressibleTurbulenceModel \
    -lcompressibleRASModels \
    -lcompressibleLESModels \
    -lreactionThermophysicalModels \
    -lspecie \
    -lfluidThermophysicalModels \
    -lODE \
    -lcombustionModels

Huge thank you in advance for any help you may be able to offer :o

Best regards,
David

dvcauwe June 10, 2013 03:55

I'm still struggling with this matter...

No one ever tried to implement a new method of calculating specie transport or thermo properties?

Best regards,

David

wyldckat June 10, 2013 08:55

Greetings David,

It's not as much of a matter as to if someone has implemented this, it's a matter of someone who might have done it to see your question And to have time to answer! ;)

... it took me a while to figure out what the problem was, because I didn't read the whole post the first couple of times... OK, so the problem is this:
Code:

-L$(FOAM_USER_LIBBIN)/libpyrolysisModel.so \
-L$(FOAM_USER_LIBBIN)/libpyThermo.so \

Instead, it should be this:
Code:

-L$(FOAM_USER_LIBBIN) \
-lpyrolysisModel \
-lpyThermo \

Best regards,
Bruno

dvcauwe June 11, 2013 12:03

Ah of course, thanks for pointing that out!

Apparently it's also not necessary to create a copy of psiThermo, psiReactionThermo, ... but all it takes is just to add the new thermo combinations to the runTimeSelectionTable. For my case this was achieved by the following code in the ktTransport.H file:

Code:

makeThermo
(
    psiThermo,
    hePsiThermo,
    pureMixture,
    ktTransport,
    sensibleEnthalpy,
    janafThermo,
    perfectGas,
    specie
);

makeReactionThermo
(
    psiThermo,
    psiReactionThermo,
    hePsiThermo,
    reactingMixture,
    ktTransport,
    sensibleEnthalpy,
    janafThermo,
    perfectGas,
    specie
);

typedef ktTransport<species::thermo<janafThermo<perfectGas<specie> >,sensibleEnthalpy> > ktGasThermoPhysics;

typedef Reaction<ktGasThermoPhysics> ktGasReaction;
makeReactions(ktGasThermoPhysics, ktGasReaction);

makeChemistryModel(chemistryModel, psiChemistryModel, ktGasThermoPhysics);

makeChemistrySolverType(ode, psiChemistryModel, ktGasThermoPhysics);

makeChemistryReader(ktGasThermoPhysics);
makeChemistryReaderType(foamChemistryReader, ktGasThermoPhysics);

Hope this might be useful for someone else some day, thanks again for your help!

Kind regards,
David

dvcauwe August 12, 2013 11:58

On a related note, I'd like to implement the Wilke equation for calculating the viscosity of a gas mixture (cf. Fluent implementation: http://www.sharcnet.ca/Software/Flue...ug/node296.htm). I see that others have done this before me and although apparently this should have become easier in 2.2.x, I can't quite figure out where to make the changes.

In my case I'm using psiReactionThermo, so I suppose the essential part is in hePsiThermo, correct?
Code:

template<class BasicPsiThermo, class MixtureType>
void Foam::hePsiThermo<BasicPsiThermo, MixtureType>::calculate()
{
    const scalarField& hCells = this->he_.internalField();
    const scalarField& pCells = this->p_.internalField();

    scalarField& TCells = this->T_.internalField();
    scalarField& psiCells = this->psi_.internalField();
    scalarField& muCells = this->mu_.internalField();
    scalarField& alphaCells = this->alpha_.internalField();

    forAll(TCells, celli)
    {
        const typename MixtureType::thermoType& mixture_ =
            this->cellMixture(celli);

        TCells[celli] = mixture_.THE
        (
            hCells[celli],
            pCells[celli],
            TCells[celli]
        );

        psiCells[celli] = mixture_.psi(pCells[celli], TCells[celli]);

        muCells[celli] = mixture_.mu(pCells[celli], TCells[celli]);
        alphaCells[celli] = mixture_.alphah(pCells[celli], TCells[celli]);

    }

At this point, can I somehow access the species values of mu and alphah and simply write a different way of calculating these mixture properties or do I need to go deeper into multiComponentMixture where "cellMixture" is defined?

Thank you in advance for any advice (and/or working code :o) to help me solve this problem!

Best regards,
David

ChrisA August 12, 2013 19:40

There's plenty of mixing discussion on the forums in how the mixture quantities are handled... A thread of my own that helped me understand whats going on greatly is here (for 2.1.1, I'm not sure what 2.2.x has done with this but first glances show its very similar)

That's a good starting point, as far as I understand it it's not quite as simple as you stated above, I'm not familiar with the wilke model, but I'm guessing you can't change the += operator to fit with the model? I've been secretly hoping someone would change around the mixture methods so we could select a mixture method at run time rather than be forced into the molar weighted coefficients but alas.

Also, your typedefs inside the .h file for the new viscosity function is absolute genius and I'm not sure why I didn't think of it before, also not sure why I haven't found it documented elsewhere (I've been appending the thermo selection code like you had been trying to avoids)... Don't mind me, just expressing my shameless gratitude in saving me hours of time when I convert my 2.1.1 viscosity libraries over to 2.2.x

dvcauwe August 13, 2013 03:52

Thanks for pointing my attention to that thread again, I must've missed some parts the first time I read it. ;) If I understand correctly, the coefficients are mixed into a new "average specie" and not the actual properties. In that case, I suppose the changes need to be made at an earlier stage than hePsiThermo, probably in multiComponentMixture...
Code:

template<class ThermoType>
const ThermoType& Foam::multiComponentMixture<ThermoType>::cellMixture
(
    const label celli
) const
{
    mixture_ = Y_[0][celli]/speciesData_[0].W()*speciesData_[0];
    for (label n=1; n<Y_.size(); n++)
    {
        mixture_ += Y_[n][celli]/speciesData_[n].W()*speciesData_[n];
    }
    return mixture_;
}

So I guess the solution would be to edit this part to differentiate between the thermo and the transport coefficients. With the additional problem that for Wilke's formula the weighting coefficients depend on the specie viscosity relative to that of all other species. This is also the reason that I don't think it's possible to include it in the += operator...

Best regards,
David

ChrisA August 13, 2013 13:35

Quote:

Originally Posted by dvcauwe (Post 445351)
If I understand correctly, the coefficients are mixed into a new "average specie" and not the actual properties.

This is correct. (I should note I'm still using 2.1.1 so this might be slightly different in 2.2.x?) You're correct in that you are essentially going to have to make your own version of multicomponentmixture (Again, I think, there keeps being talk of easier implementation of different mixture models but I haven't seen it yet). I don't quite 100% understand the way speciesData_ is organized so I'm not going to be much more help on that front :(.

If you figure it out I'd be very interested to hear how you do it, my 2.1.1 code calculates and mixes the viscosity in the solver, which is a non ideal solution if I ever want to compare viscosity models or change it up. I might take another look at the whole thing as well when I finish with the current paper I'm on but that might be awhile.

ChrisA August 13, 2013 19:49

Quote:

Originally Posted by dvcauwe (Post 433420)
Ah of course, thanks for pointing that out!

Apparently it's also not necessary to create a copy of psiThermo, psiReactionThermo, ... but all it takes is just to add the new thermo combinations to the runTimeSelectionTable. For my case this was achieved by the following code in the ktTransport.H file

I gave this a shot and there's a couple questions I'd like to ask you... firstly, did you have to recompile Basic/ReactionThermo/ChemistryModel as well? or is it enough to add the new library/source to the solver and the compiler sorts all that out when the solver is compiled?

Secondly, I'm noticing a long list of "duplicate entry in runtime selection table" errors/warnings at runtime when I do this:

Code:

Duplicate entry irreversibleLangmuirHinshelwoodReaction in runtime selection table Reaction
#0    /opt/openfoam221/platforms/linux64GccDPOpt/lib/libOpenFOAM.so(_ZN4Foam5error14safePrintStackERSo+0x25) [0x7fe3b19693b5]
#1    /home/carisman/OpenFOAM/carisman-2.2.1/platforms/linux64GccDPOpt/lib/libRCBFV2.0Specie.so(+0x2cac28) [0x7fe3b3628c28]
#2    /lib64/ld-linux-x86-64.so.2(+0xf306) [0x7fe3b6483306]
#3    /lib64/ld-linux-x86-64.so.2(+0xf3df) [0x7fe3b64833df]
#4    /lib64/ld-linux-x86-64.so.2(+0x16ea) [0x7fe3b64756ea]

Did you encounter these as well? Or is there something horribly wrong with what I've done...


Never mind the above, there's some other screwy things going on causing that I believe... Back on topic, I cannot get the code to compile by adding the appropriate selection combinations to the ktTransport.H file, did you include anything else in your ktTransport to get your method to work? It seems like I'm missing a vital step. The error I get is:

Code:

lnInclude/hollisTransport.H:52:1: error: expected constructor, destructor, or type conversion before ( token

dvcauwe August 14, 2013 04:09

1 Attachment(s)
Quote:

Originally Posted by ChrisA (Post 445558)
I gave this a shot and there's a couple questions I'd like to ask you... firstly, did you have to recompile Basic/ReactionThermo/ChemistryModel as well? or is it enough to add the new library/source to the solver and the compiler sorts all that out when the solver is compiled?

Secondly, I'm noticing a long list of "duplicate entry in runtime selection table" errors/warnings at runtime when I do this:

I didn't need to touch ChemistryModel but now that I look back on it, I did need a different version of makeReactions.H and I'm pretty sure that was because of the same duplicate entries.

Maybe it's easiest if I just send you the files. Basically it's a full copy of sutherland with some alterations to ktTransport.H and ktTransportI.H and a bunch of includes to allow these typedefs. Let me know if it works for you!

Best regards,
David

ChrisA August 14, 2013 19:09

Ah I see, I thought it might have been some missing includes... what got me was that I didn't quite walk through where those functions are coming from...

I see one issue with this (forgive any miss wording C++ is a side thing for me), you're including files higher up in the chain, things that rely on specie to be there in the first place... this kind of cyclic referencing can lead to nastyness pretty easily (I think, feel free to correct me oh great C++ gurus).

The compiler seems to be having a terrible time compiling my attempt to replicate your code, it would appear that the cyclic references are throwing it for a loop... Thanks for your help though, I think I'll just pawn editing the Basic/ReactionThermo/ChemistryModel files off on a summer student and avoid the cyclic-ness :p.

As a side note, if you run into any issues with your mixture implementation let me know, I'll be glad to help as much as I can.

wyldckat August 16, 2013 13:27

Greetings to all!

I honestly have to say that I'm a bit lost on this thread. I'm not sure how I can give some help here, but perhaps the following threads/posts might assist you in getting a few more ideas:
Both of those threads address similar topics, when it comes to the header inclusion and usage of macros, so perhaps it can help you to get some more ideas on how OpenFOAM handles these headers and macros.

Best regards,
Bruno


All times are GMT -4. The time now is 10:18.