CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM Announcements from Other Sources (http://www.cfd-online.com/Forums/openfoam-news-announcements-other/)
-   -   Equation Reader released (http://www.cfd-online.com/Forums/openfoam-news-announcements-other/78439-equation-reader-released.html)

marupio July 21, 2010 22:54

Equation Reader released
 
Hello Foamers,

I am pleased to announce the release of equationReader, an extension to OpenFOAM that allows you to use equations in a dictionary file.

For instance, constant/transport properties can now be:

Code:

nu    [0 2 -1 0 0 0 0] "3/(1+exp(0.5^2))";
It uses the same syntax as Excel. It works out-of-the-box for all your solvers (although somewhat limited). If you design a solver to use it, it can perform substitutions, and automatic variable updates at every timestep.

For more details, please see the page on the wiki:

http://openfoamwiki.net/index.php/Co...equationReader

Share and enjoy!

[EDIT: Sorry for the duplicate post. I thought this one had failed. Please use this thread for commenting.]

kiran August 8, 2010 11:59

help me for decoding sonicFoam solver
 
hi marupio

i am working sonicFoam solver in OpenFoam 1.6 the problem is iam unable to decode the equations in the solvers. especially p.eqn U.eqn h.eqn.

if u can help me in this regard i will be thank full to you and i did not find any equation decoder from link which you have given.

please help me out.

marupio August 8, 2010 16:43

Hi kiran, thanks for the question. The equationReader extension allows OpenFOAM to read equations from a dictionary, such as transportProperties, or fvSchemes. It has nothing to do with the fvc and fvm namespaces that contain those elegant equations for the solver to work with. I was working on a few articles to explain these equations, but they are still in working form.

The equations are vector operators, so if you can find the vector or tensor forms of the equations for your solver, you'll have a better time translating them. Have a look at:

http://www.openfoam.com/features/creating-solvers.php

It shows some translation.

Good luck!

calebamiles June 26, 2011 08:16

Questions about equationReader
 
Hello Marupio,

I have a few questions about equationReader. I would like to use equationReader's functionality to write a solver with some simple time-dependent source terms in the momentum given by an analytical expression, like cos(a*t) where t is the current time, since the source term is a vector quantity would it be possible to add a source term of the form (A1*vector(1,0,0) + A2*vector(0,1,0) + A3*vector(0,0,1)) where A = (A1, A2, A3) is the vector I would like to add as a source term and A1, A2, A3 are defined in a dictionary. Also what is the best way to create a time dependent expression in the dictionary.

Caleb

marupio June 26, 2011 10:51

Hi Caleb,

(I was just looking through the wiki documentation, and it needs some improving!)

First: adding time as a variable. I believe runTime is derived from dimensionedScalar, but I don't think you should add it directly, because its name changes at every timestep. I'd recommend something like this:

Code:

    eqns.addDataSource
    (
        runTime.value(),
        "time",
        dimTime
    );

Then in your equations, the variable name is "time". Or you could rename the second argument above to "t", or whatever, and that is your variable name.

As for vector construction, you have the right principles... I haven't worked with the vectors much so I don't know if you have the right grammar. A1, A2 and A3 all get values from equationReader, so either they are actively updated, and you are using eqns.update(); or they are passively updated, such as:

Code:

A1 = eqns.evaluate("A1");
A2 = eqns.evaluate("A2");
A3 = eqns.evaluate("A3");

Or use their indices for faster evaluates:


Code:

// Before the solver loop, after equations are read
label a1index = eqns.lookup("A1");
label a2index = eqns.lookup("A2");
label a3index = eqns.lookup("A3");

// Within solver loop
A1 = eqns.evaluate(a1index);
A2 = eqns.evaluate(a2index);
A3 = eqns.evaluate(a3index);

To create a time dependent expression in the dictionary, have a dictionary with the equations defined, such as:

Code:

    A1      "cos(piByTwo_ * time)";
    A2      "sin(pi_ * time)";
    A3      "tan(sqrt(time))";

And add the dictionary as a source, and read the equations in:

Code:

eqns.addDataSource(myDict);
eqns.readEquation(myDict, "A1");
eqns.readEquation(myDict, "A2");
eqns.readEquation(myDict, "A3");

I hope that helps. Any more questions, feel free to ask!

Now, I need to update those instructions...

calebamiles June 26, 2011 20:22

Hello Marupio,

Thank you for the quick reply. Just to be clear the code


Code:

eqns.addDataSource(       
            runTime.value(),       
            "time",       
            dimTime    );

will go directly into my new solver after the code snippet from the "Programming with equationReader" from the wiki

Code:


Code:

#include "IOequationReader.H" // ...   
IOEquationReader eqns
(       
 IOobject
(           
"eqns",           
runTime.
timeName(),           
 runTime,           
IOobject::
READ_IF_PRESENT,           
IOobject::
AUTO_WRITE // Set to NO_WRITE to suppress output       
),        false // set to true to show data sources in output file   
);



I will then add the dictionary containing my equations in the similar manner demonstrated in the demo after which I will read the equations using the syntax described above and then within the solver I will evaluate the equations? Thanks so much for your help.

Caleb

marupio June 27, 2011 08:05

That's right. Let me know if it doesn't work!

calebamiles July 9, 2011 01:30

Position Dependent Equations
 
Hello Marupio,

Thank you so much for all your previous help. I have another simple question, would this be the right way to go about implementing an equation that depends on position

Code:

    eqns.addDataSource(       
        mesh.C().component(0),       
        "x",       
        dimLength   
    );



Thank you so much for all of your help. Also nice changes to the wiki.

Caleb

marupio July 9, 2011 10:51

Hi caleb,

Thanks for the question. I think I should include time/space variables as default sources in future versions. Better yet, I should include vectors and tensors in the equations.

When adding a data source, the one thing you have to ask yourself is: "Will that value always be at that same memory address?" Looking at your code, we have:

Code:

    eqns.addDataSource(       
        mesh.C().component(0),       
        "x",       
        dimLength   
    );

First I checked out mesh.C(), and yes, this is a permanent object, once it is created. (I'd be worried if your solver does any form of remeshing, though.) Then I checked out component(), and this one is a problem. It returns a temporary field.

If your mesh isn't moving, the simple solution is to make it a permanent field, e.g.:

Code:

// Before your main solver loop... perhaps in createFields.H
volScalarField meshX("x", mesh.C().component(0));
eqns.addDataSource(meshX);

If your mesh is moving, you also have to add something like:

Code:

// At the start of an iteration, just after the mesh has moved
meshX = volScalarField("x", mesh.C().component(0));

I haven't tried this. Let me know if something doesn't work!

calebamiles July 9, 2011 12:29

Hello Marupio,

Thanks for the help, you suggestion compiles fine. I have another related question, after using "x" in an equation, for example "cos(x)", defined in a dictionary would I then have to use "eqns.evaluateField" to evaluate the "cos(x)" equation.

Caleb

marupio July 9, 2011 13:05

Hi Caleb,

I really need to fix up the wiki documentation (haven't touched it yet)! Apparently my example for using the eqns.evaluateField() is totally wrong! Here's how you use it:

Code:

// vs is a volScalarField that you want the result in
eqns.evaluateField("nameOfCosXequation", vs);

Apparently it doesn't return a value like the other "evaluate()" functions do... instead, it returns void. The result of the evaluation is returned into the second parameter (vs, above). This is contrary to OpenFOAM conventions (no returning values through parameters). I've put that on the list of things to improve - but don't worry, I'll make it backwards compatible.

calebamiles July 9, 2011 23:59

Hello Marupio,

Thanks for the guidance. Do you know offhand how to define an empty volScalarField in a more elegant way than something like

Code:

volScalarField Null1 = 0*mesh.C().component(0)


Preferably, for readability, I would also like to avoid using an IOobject if at all possible.

Caleb

calebamiles July 10, 2011 03:51

Hello Marupio,

I also have an error that perhaps you can sort out. I have a volScalarField called "bOneX" that I have initialized to all zeros. When I try

Code:

eqns.evaluateField("b1X", bOneX);
I get the error

Quote:

my_icoFoam.C: In function ‘int main(int, char**)’:
my_icoFoam.C:255:32: error: no matching function for call to ‘Foam::equationReader::evaluateField(const char [4], Foam::volScalarField&)’
/opt/openfoam171/src/OpenFOAM/lnInclude/equationReader.H:717:18: note: candidate is: void Foam::equationReader::evaluateField(const Foam::label&, Foam::scalarList&, Foam::dimensionSet&)
any idea what's wrong? Thanks again for all of your help.

Caleb

marupio July 10, 2011 11:13

Hi Caleb,

Sorry, my bad. There's presently no evaluateField("equationName", putResultsHere) function. But there is a evaluateField(equationIndexNumber, putResultsHere) function. You can get the equationIndexNumber (which never changes once the equations are read in) using equationIndexNumber = lookup("equationName")

As for the null volScalarField, you can lookup the constructors in GeometricField.H - they start at line 272. Most of them require IO objects... I think you should get used to seeing the IO object constructor. If you really don't want to see it, you could create a dummy volScalarField with the correct dimensions in createFields.H (using the full IOobject constructor), then use the renaming copy constructor in the solver body.


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