[joel.lehikoinen]

Hi,

I need to create a solver based on buoyantPimpleFoam which solves an arbitrary (around 10) scalar transport equations. I have written a solver which solves one passive scalar transport equation. I could just copy-paste that code enough many times, but I would like to know if a more elegant solution is possible: Namely, how easy it would be to implement a for-loop that solves N scalar transport equations, where N is specified by the user? The problem is, I don't know C++ and I'm not too familiar with the classes present in OpenFOAM.

As far as I know, the vectors in OpenFOAM are always 3-dimensional. But could I use a tensor of rank (1,N) or (N,1) to store the scalars? Or is using a tensor a bad idea, if I want to have different BCs for the scalars? Tensor would be nice because then I wouldn't have to worry about how to get OF to read/write files with indices in their names (scalar1, scalar2, etc.).

This is not a critical problem, as I said I can just copy-paste the code snippet that solves one scalar transport equation N times, but if I want to change the number of scalars later it becomes cumbersome to recompile the solver every time. I was just wondering if anyone with more knowledge of C++ and OF source code knows these things.

Regards,
Joel

[l_r_mcglashan]

Sure, you can use a PtrList. Here's how you initialise them:

//  Create a list of pointers to the mass fraction field of each species.
    PtrList<volScalarField> chemicalSpecies(numberOfSpecies);
    for (label i=0; i<numberOfSpecies; ++i)
    {
        Info << "Creating Species " << namesOfSpecies[i] << endl;
        chemicalSpecies.set
        (
            i,
            new volScalarField
            (
                IOobject
                (
                    namesOfSpecies[i],
                    runTime.timeName(),
                    mesh,
                    IOobject::NO_READ,
                    IOobject::AUTO_WRITE
                ),
                mesh,
                dimensionedScalar("zero", dimless, 0.0),
                mixtureFraction.boundaryField().types()
            )
        );
    }

Then use a for loop to solve the equations, as done in, for example, YEqn.H in solvers/combustion/chemFoam

[joel.lehikoinen]

Thank you very much for your help, I managed to create the solver I will need.

[dinksy]

I am also trying to solve a similar problem of solving transport equation of n scalars, say c, and would like to define them as an array, such as c[n] instead of defining n separate scalarfields. The above thread was of some help but if you can tell me which solver it is a part of then I could look into it to clarify my doubts. The main problem I face is with the definition of c as an array in the initial time folder [0].
Regards

[marc.immer]

Hi,
I'm asking myself the same question (how to initialize the array in the 0 folder).
Did you make any progress?

Marc

[olivierG]

hello,

Why you don't use reactingFoam or rhoReacting(Buoyant)Foam ?
At least take a look at this solver.

regards,
olivier

[marc.immer]

Thanks, reactingFoam is helpful. But I guess there is no easier way than duplicating/renaming the files in the 0 Folder for each species. But this at least can be automated with a script.

Cheers
Marc

[olivierG]

yes and no

If you take a look at reactingFoam tutorial, for product of reaction, you have a "Ydefault" file, in which you specify the defaut boundary conditions for the species. So you may try to use this.

regards,
olivier

[marc.immer]

I took a look at the reactingFoam solver. Since I'm only interested in a passive scalar transport, having a whole thermodynamic model might be a bit overkill... I'm not sure how I could only extract the classes I need, since it all seems to dependend on the classes basicMultiComponentMixture and chemistryReader.

Regards
Marc

[marc.immer]

Hi,

I managed to get it to work. I used parts of the reactingFoam for the loop structure,

tmp<fv::convectionScheme<scalar> > mvConvection
	(
	    fv::convectionScheme<scalar>::New
	    (
		mesh,
		fields,
		phi,
		mesh.divScheme("div(phi,si_h)")
	    )
	);
	volScalarField kappaEff
	(
	    "kappaEff",
	    turbulence->nu()/Pr + turbulence->nut()/Prt
	);

	forAll(s, i)
    	{
	    volScalarField& si = s[i];

	    tmp<fvScalarMatrix> siEqn
	    (
		    fvm::ddt(si)
		  + mvConvection->fvmDiv(phi, si)
		  - fvm::laplacian(kappaEff, si)
	    );

	    sources.constrain(siEqn());
            solve(siEqn() == sources(si),mesh.solver("si"));
        }

and copy pasted parts of the chemistryReader to get the initial files.
To define the species, I added a list into the transportProperties dict. (makes it easier to read in the createFields.H

wordList names(transportProperties.lookup("scalars"));
    PtrList<volScalarField> s(names.size());

    forAll(s, i)
    {
        IOobject header
        (
            names[i],
            mesh.time().timeName(),
            mesh,
            IOobject::NO_READ
        );

        // check if field exists and can be read
        if (header.headerOk())
        {
            s.set
            (
                i,
                new volScalarField
                (
                    IOobject
                    (
                        names[i],
                        mesh.time().timeName(),
                        mesh,
                        IOobject::MUST_READ,
                        IOobject::AUTO_WRITE
                    ),
                    mesh
                )
            );
        }
        else
        {
            volScalarField sdefault
            (
                IOobject
                (
                    "sdefault",
                    mesh.time().timeName(),
                    mesh,
                    IOobject::MUST_READ,
                    IOobject::NO_WRITE
                ),
                mesh
            );

            s.set
            (
                i,
                new volScalarField
                (
                    IOobject
                    (
                        names[i],
                        mesh.time().timeName(),
                        mesh,
                        IOobject::NO_READ,
                        IOobject::AUTO_WRITE
                    ),
                    sdefault
                )
            );
        }
    }
    multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;

    forAll(s, i)
	{
	    fields.add(s[i]);
	}