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/)
-   -   How to interface a Fortran thermodynamic tool with OpenFOAM? (http://www.cfd-online.com/Forums/openfoam-programming-development/81534-how-interface-fortran-thermodynamic-tool-openfoam.html)

Cyp October 29, 2010 09:18

How to interface a Fortran thermodynamic tool with OpenFOAM?
 
Hi everybody!

I am looking to interface a specific thermodynamic tool with OpenFOAM. This program is written in Fortran. From a Pressure, a temperature and a composition it can provide a density : rho=f(P,T,x)

I wonder how to bound these different programs. I would like to keep the OpenFOAM structure and evaluate the density at each time step with a command as : thermo.rho();

Is anybody already performed such an interfacing ? Can you give some hints to start ? from which thermophysicalModels can I be inspired?

Regards,
Cyp

Chris Lucas November 1, 2010 06:13

Hi,

compile the fortran programm as a dynamic library. Then, use that library in OpenFoam.

Is the programm you want the use in OpenFoam Refprop (NIST)?

Regards,
Christian

Cyp November 1, 2010 11:44

Hi Chris!

Thank you for your answer!! Indeed it is exactly the NIST REFPROP I want to use with OpenFOAM! Did you already experienced to call REFPROP from OpenFOAM ?


For the moment, I compiled the Fortran subroutine as a shared object (.so). How can I call it from OF ??

Best regards,
Cyp

Chris Lucas November 1, 2010 11:53

Hi,

yes, I have implemented Refprop in OpenFoam. However, I noticed that Refprop is far too slow to be used in OpenFOAM. The simulation time will increase in order of many magnitudes.

Regards,
Christian

kdneroorkar November 1, 2010 12:03

Hi
another way to do this (we use this with REFPROP) is to create a lookup table using REFPROP with rho as a function of (P,T,x). Then write a class called, say, thermoprops. when you call thermo.rho(P,T,x) (where thermo is of type thermoprops) , it would interpolate in this table and return the required rho.
Hope this helps
Kshitij

Cyp November 2, 2010 05:18

Hi!

Quote:

Hi,

yes, I have implemented Refprop in OpenFoam. However, I noticed that Refprop is far too slow to be used in OpenFOAM. The simulation time will increase in order of many magnitudes.

Regards,
Christian
It is not really a good news....

Quote:

Originally Posted by kdneroorkar (Post 281681)
Hi
another way to do this (we use this with REFPROP) is to create a lookup table using REFPROP with rho as a function of (P,T,x). Then write a class called, say, thermoprops. when you call thermo.rho(P,T,x) (where thermo is of type thermoprops) , it would interpolate in this table and return the required rho.
Hope this helps
Kshitij

Do you mean I have to generate a table from REFPROP and then pick up the right value of rho in this table from P, T and x ?

I think my main problem is either I chose a direct call of REFPROP or an interpolation from a table, is the creation of the thermoprops class... From which existed class can I start ??

Regards,
Cyp

kdneroorkar November 2, 2010 09:47

yes. so you would have something like
#P T x rho
1000 273 0.01 900
1200 273 0.01 890
.....
and you would interpolate to get rho=f(P,T,x). another idea you could try ( I dont know if this is feasible in your case), instead of using P, T and x, you could use P,H. This way the interpolation would be simplified.

I dont know of any existing class that does this interpolation. We had to write our own.

Cyp November 2, 2010 10:37

Thank you for your answers! For the moment, I would like to test to implement a direct call from REFPROP instead of a reading in a table.

As a first test, and to become familiar with the thermodynamic class programmation in OpenFOAM, I try to create my own class from basicRhoThermo : rho=f(T,P). However, when I read the .C and .H files, I can't where rho is calculated...

Do you have any idea ?

Regards,
Cyp

kdneroorkar November 2, 2010 11:38

basicRhoThermo is a base class. since rho() is a virtual function, it would be calculated in some class that is derived from basicRhoThermo. If you look at some class that is derived from basicRhoThermo, you might get an idea.
I am still using 1.5-dev so I dont have basicRhoThermo. Sorry that I cant be of more help

Chris Lucas November 2, 2010 16:55

Hi,

might I ask, what exactly you want to simulated?

Ok, here is what you have to do to implement refprop in OpenFoam. This is not the best way to do it, but it works.

Start from the file hPsiThermo (or hRhoThermo). In this file, you have to call the Refprop functions. In the constuctor, you have to create the Refprop object. You also need a second "rho" field. This is needed because rho not directly connected to psi anymore (as it is for perfect gases). In this class, all the functions have to be changed so that the Refprop functions are called.

Now comes the tricky part. Firstly, you have to change the enthalpy BC's. The "normal" BC are programmed for perfect gases, so we have h(T) instead of h(T,p). You also need to change the functions in basicThermo called by the enthalpy BC.



As I said, Refprop is very slow and you have to solve the fundamental equations for each grid point (without changing refprop even more than once).

About using tables, it is a possibility. But one problem is that you get an error and I don't mean the "normal" interpolation error. The problem is that thermo physical properties are connected to each other e.g. c_p=dh/dT. So, the interpolated value of c_p might not correspond correctly with you enthalpy.


Regards,
Christian

Cyp November 3, 2010 06:12

Hi Christian!

Thank you for your extended answer! In my investigation fields, the final step will be to simulate a two-phase multicomponent flow. I already have developed models for the mass fractions in each phase and I need to obtain the density of each phase from REFPROP. For each phase, density depends on Temperature, Pressure and Composition. It is my final objectives. I will consider a rough mesh. So I hope that the use of REFPROP will not have a consequent impact on my CPU time...

For the moment, to become familiar with OpenFoam/REFPROP interfacing, I tried the implementation on a simpler model : just one phase is considered. There is no question of enthalpy, Cp and other thermodynamic properties. Just the density (I call the TPRHO() Refprop subroutine).

As you advised me, I started from hPsiThermo. I saw enthalpy, Cp and Cv variables. Since I not consider such variables, I suppose I can rid of them ? And so I can rid of the part where I should change the enthalpy BC, can't I ?

I am not sure to well understand the part with the second density. Do you mean I have to create a rho2 variable in hPsiThermo.H file ?

Moreover, if I understand well, it is in the hPsiThermo.C file that I call the REFPROP function ?

Thank you very much for your help.

Chris Lucas November 3, 2010 06:48

Hi,

I'm also interested in multiphase flow. So you have phase change in your application? How did you get your pressure equation?

But back to your question:

"As you advised me, I started from hPsiThermo. I saw enthalpy, Cp and Cv variable. Since I not consider such variables, I suppose I can rid of them ?"

For now, yes. But for a real compressible flow, no, you need this functions like
"template<class MixtureType>
Foam::tmp<Foam::scalarField> Foam::hPsiThermo<MixtureType>::Cp"
even so the variables given to the function have to be changed. Within, this function, you must call the refprop routines.

"I am not sure to well understand the part with the second density. Do you mean I have to create a rho2 variable in hPsiThermo.H file?"

About the second density field. If you look at the rho function in basicPsiThermo.C
( virtual tmp<volScalarField> rho() const
{
return p_*psi();
})

You see that rho is connected to psi for perfect gases, so that no rho field has to be defined to save the results for rho. For real fluids, there is no relationship as the one above between rho and psi, so you have the save the results for rho separately.

"And so I can rid of the part where I should change the enthalpy BC, can't I ?"

If you don't have a energy equation, yes, you bypass the BC problem. However this is not correct for compressible flows
This BC functions are in different classes (within a subsubfolder of the basic folder). This classes call functions in basicThermo, which are later overwriten by the functions in your new "hPsiThermo".

"(I call the TPRHO() Refprop subroutine)".

Unless you rewrite the energy equations (I would not use energy equations based on temperature for real fluids), you need a function like rho(p,h) and rho(p,T) (for BC's).

"Moreover, if I understand well, it is in the hPsiThermo.C file that I call the REFPROP function?"

Yes, you can.
Regards,
Christian

Cyp November 3, 2010 10:38

Yes I have phase change in my application. For the pressure equation, since I work in a porous medium, I developed an IMPES method (Implicit Pressure, Explicit Saturation). Pressures within each phase are related by a capillary pressure which depends on the saturation.

I was thinking about the creation of a second rho to rid of the perfect gas law. I tried to define a new base class which is a copy of basicPsiThermo (I called it testThermo1). In the equivalent of the makeBasicPsiThermo.H file, I could precised only two function : the Cthermo and the Mixture. So I could have only one rho variable. What do you think about that ?

In a second step, I defined a new class deriving from the previous class and MixtureType. (exactly what have been done for the hPsiThermo class, I called it testThermo2). Then I modify the equivalent of the hPsiThermos.C file in order to define only Cthermo and pureMixture.... but it doesn't compile.....

I get :

Code:

cyp@cyp-laptop:~/OpenFOAM/cyp-1.7.1/applications/lib$ wmake libso
Making dependency list for source file testThermo2/testThermos2.C
SOURCE=testThermo2/testThermos2.C  ;  g++ -m64 -Dlinux64 -DWM_DP -Wall -Wextra -Wno-unused-parameter  -Wold-style-cast -Wnon-virtual-dtor -O3  -DNoRepository  -ftemplate-depth-40 -I/opt/openfoam171/src/finiteVolume/lnInclude    -I/opt/openfoam171/src/thermophysicalModels/specie/lnInclude    -I/opt/openfoam171/src/thermophysicalModels/basic/lnInclude -IlnInclude -I. -I/opt/openfoam171/src/OpenFOAM/lnInclude -I/opt/openfoam171/src/OSspecific/POSIX/lnInclude  -fPIC -c $SOURCE -o Make/linux64GccDPOpt/testThermos2.o
testThermo2/testThermos2.C:49: error: type/value mismatch at argument 1  in template parameter list for ‘template<class MixtureType> class  Foam::testThermo2’
testThermo2/testThermos2.C:49: error:  expected a type, got ‘pureMixture’
testThermo2/testThermos2.C:49: error: invalid type in declaration before ‘;’ token
testThermo2/testThermos2.C:49: error: ‘testThermo2pureMixture’ is not a class or namespace
testThermo2/testThermos2.C:49: error: ‘testThermo2pureMixture’ is not a class or namespace
testThermo2/testThermos2.C:49: warning: unused variable ‘Foam::debug’
testThermo2/testThermos2.C: In constructor ‘Foam::testThermo1::addfvMeshConstructorToTable<testThermo1Type>::addfvMeshConstructorToTable(const Foam::word&) [with testThermo1Type = int]’:
testThermo2/testThermos2.C:56: error: ‘typeName’ is not a member of ‘int’
In file included from lnInclude/makeTestThermo1.H:35,
                from testThermo2/testThermos2.C:28:
lnInclude/testThermo1.H: In static member function ‘static Foam::autoPtr<Foam::testThermo1> Foam::testThermo1::addfvMeshConstructorToTable<testThermo1Type>::New(const Foam::fvMesh&) [with testThermo1Type = int]’:
lnInclude/testThermo1.H:78:  instantiated from ‘Foam::testThermo1::addfvMeshConstructorToTable<testThermo1Type>::addfvMeshConstructorToTable(const Foam::word&) [with testThermo1Type = int]’
testThermo2/testThermos2.C:56:  instantiated from here
lnInclude/testThermo1.H:71: error: cannot convert ‘const Foam::fvMesh’ to ‘int’ in initialization

I do not understand why testThermo2pureMixture is written a single word...

I think it should come from my makeTestThermo.H file :
Quote:

#ifndef makeTestThermo1_H
#define makeTestThermo1_H

#include "testThermo1.H"
#include "addToRunTimeSelectionTable.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#define makeTestThermo1(Cthermo,Mixture) \
\
typedef Cthermo<Mixture> \
Cthermo##Mixture; \
\
defineTemplateTypeNameAndDebugWithName \
( \
Cthermo##Mixture, \
#Cthermo \
"<"#Mixture">", \
0 \
); \
\
addToRunTimeSelectionTable \
( \
testThermo1, \
Cthermo##Mixture, \
fvMesh \
)


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif
I can't see my mistake...

Best Regards,
Cyp

nimasam November 16, 2010 21:40

how can i make dynamic file from a fortan code and use it in openFoam ?

karthik1414 January 20, 2011 04:49

Basicthermo.H
 
hello everyone,
I am new to the OpenFoam software.
Can anyone tell me how I can make use of Basicthermo.H while using Interfoam as solver??

nimasam January 20, 2011 05:23

may be its good to open new thread for ur issue, however you should look at user guide chapter 3 but a brief description can be like that:
in interfoam.C u should include ur new class:
# include "Basicthermo.H"
then in case directory use "wmake" to compile it.

karthik1414 January 24, 2011 09:43

hey nimasam,
Thank you for your reply, i tried that.
Any idea how to vary density with pressure change??

nimasam January 24, 2011 09:56

look in compressible solver for example "compressibleInterFoam" and compare it with interFoam

pbe_cfd March 7, 2011 15:15

Hi,

I have couple of fortran.90 subroutines which I would like to call from OpenFOAM1.6. What I know direct call of fortran subroutine from c++ is not good. Thus, one has to have C interface to routines and some more technical difficulties as, fortran stores rowwise C stores columnwise :confused: Whatever, these are general problems of C++ and Fortran mixing. My questions are more OpenFOAM specific:

1) What should be the compiled fortran source codes, dynamic library, shared objects (.so), objects (.o) ... ? If more than one is good, which one performs better?

2) How should be this "executables" (compiled files) linked to OpenFOAM, let's say scalarFoam? How should the <file> and <options> in Make folder be modified?

Please no suggestion like "why don't you rewrite your subroutines in C++" ;-)

All the best,

@ Cyp: You see, I am chasing you ;-)

Cyp March 9, 2011 06:15

Hi!

To link a library with OpenFOAM you can follow this 2 steps :

- First ,you create a *.so file using "wmake libso"
- Then you link this *.so file to OpenFOAM by including the following command in your controlDict file :
Code:

libs
(
    "*.so "
);

Best,
Cyp


All times are GMT -4. The time now is 05:02.