[ahcai007]

Hi,

I am using the tutorial dieselFoam with its global n-heptane mechanism to simulate diesel combustion. I wanted to output the reaction rates for CO2 as the fuel combusts. I checked that in basicChemistryModel.H and found that RR() refers to Reaction rates. However, when I type "chemistry.RR()" it does not work.

If I write it this way where I include the chemistry.RR() in a for loop.

PHP Code:

forAll (Y , celli)
{
       scalarField rr=chemistry.RR(celli);
       Info<< "Reaction rate = " << rr << endl << endl;
{ 


it will produce a results as seen in log.txt.

If I were to write it as such:

PHP Code:

forAll (Y , celli)
{
       scalarField rr=chemistry.RR(celli);
       Info<< "Reaction rate = " << rr[3] << endl << endl;
{ 


it will produce an output as seen in log2.txt.

My guess is that "rr" will output the reaction rates of ALL the species at ALL the cells while "rr[3]" will only output the reaction rates of ALL the species at cell[3].

What I had in mind is to extract only the RR for one of the species (for my case, CO2) over all the cells.

Thanks.


Ahcai

[mturcios777]

The RR in basicChemistryModel is a virtual function and is implemented later, in ODEChemistry (this isn't immediately obvious, but by looking at the thermyphysicalType in the thermphysicalProperties file you can see it). It also returns a scalar field, which is why you are able to initialize rr in both cases. If you look at ODEChemistryModelI.H at line 71, what RR(const label i) actually does is return the reation rate for species i.

In both cases, rr contains the entire RR field for whatever species happens to be in index cellI. Since Y is a PtrList for scalarfields, what you are actually doing is stepping through all 5 species. Have a look at YEqn.H to see what I mean. In the first case, you are returning the entire RR field for each species (for each iteration, you should see 5 outputs of "Reaction Rate". In the second, you are obtaining the reation rate for each species at cell index 3 (again, you have 5 outputs corresponding to the 5 species).

What you need to find is the index of CO2 in your thermo setup, then call chemistry.RR(CO2 index), which will give what you are looking for.

[ahcai007]

Originally Posted by mturcios777

The RR in basicChemistryModel is a virtual function and is implemented later, in ODEChemistry (this isn't immediately obvious, but by looking at the thermyphysicalType in the thermphysicalProperties file you can see it). It also returns a scalar field, which is why you are able to initialize rr in both cases. If you look at ODEChemistryModelI.H at line 71, what RR(const label i) actually does is return the reation rate for species i.

In both cases, rr contains the entire RR field for whatever species happens to be in index cellI. Since Y is a PtrList for scalarfields, what you are actually doing is stepping through all 5 species. Have a look at YEqn.H to see what I mean. In the first case, you are returning the entire RR field for each species (for each iteration, you should see 5 outputs of "Reaction Rate". In the second, you are obtaining the reation rate for each species at cell index 3 (again, you have 5 outputs corresponding to the 5 species).

What you need to find is the index of CO2 in your thermo setup, then call chemistry.RR(CO2 index), which will give what you are looking for.

Hi Marco,

Thanks for the explanation. It helped a lot and I manage to make it work =)

Cheers~

[lx882211]

Originally Posted by ahcai007

Hi Marco,

Thanks for the explanation. It helped a lot and I manage to make it work =)

Cheers~

I meet a problem in Openfoam 2.3.0, reactingFoam, i need the data of CH4 reaction rate in contour(figure), but i failed to autowrite it like U or Yi.

First i added this code in the CreatField.H,

volScalarField Rrate
(
IOobject
(
"Rrate",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedScalar("Rrate", dimMass/dimVolume/dimTime, 0.0)
);

Then, i added this code in YEqn,

forAll(Y, i)
{
if (Y[i].name() != "CH4")
RR = reaction->R(Yi);
}

After wmake, it shows

YEqn.H:26:14: error: no match for ‘operator=’ (operand types are ‘Foam::volScalarField {aka Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>}’ and ‘Foam::tmp<Foam::fvMatrix<double> >’)
Rrate = reaction->R(Yi);

Can you tell me how to do that? Thank you very much.