[dmoroian]

Hello everyone,
I have this issue that I could not find an answer for it on the forum: how to set a value in a field variable, but only in cells belonging to a certain cellSet?

A simple example is that I want to patch a scalar field with value 1 in a subdomain of my computational domain. This is of course done easy with setFields, but I need a snippet, and until know I have this:

    //generate the scalar field, initialized with 0.0
    volScalarField Tij(
        IOobject
        (
                "Tij",
                runTime.timeName(),
                mesh,
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
        ),
        (mesh.C() & mesh.C())*0.0
    );
 
    //define the limits of the subdomain where I want to patch with 1.0
    scalar xmin = -0.015;
    scalar xmax = 0.035;
    scalar ymin = -0.015;
    scalar ymax = 0.035;
    scalar zmin = -0.0;
    scalar zmax = 0.02;
 
    //patch the cells that have the centroid inside the desired domain
    forAll(Tij,cellI)
    {
        vector centroid = Tij.mesh().C()[cellI];
        if(centroid[0] >= xmin &&
                centroid[0] <= xmax &&
                centroid[1] >= ymin &&
                centroid[1] <= ymax &&
                centroid[2] >= zmin &&
                centroid[2] <= zmax)
 
                Tij[cellI] = 1.0;
    }

The above code works but I would like to use a cellSet instead. Therfore I tried the following:

    //generate the scalar field, initialized with 0.0
    volScalarField Tij(
        IOobject
        (
                "Tij",
                runTime.timeName(),
                mesh,
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
        ),
        (mesh.C() & mesh.C())*0.0
    );
 
    //read a previously defined cellSet (the set was defined using the 
    //utility setSet)
    cellSet cSet(mesh,"box");
 
    //patch all the cells contained within cSet with 1.0
    forAll(cSet,cellI)
    {
        Tij[cSet[cellI]] = 1.0;
    }

...without succes , meaning the cells are not patched with 1, they remain with 0. The cSet[cellI] does not return the cell index.
Any hints?

Dragos

[tomislav_maric]

A cellSet is a HashSet is a HashTable is a HashTableCore. If you follow the inheritance and the typedef used for the HashSet, you'll see that the cellSet relies on the uniqueness of keys in the HashSet to prevent duplication of the stored cell IDs. You can access the cell labels by accessing the HashTable:toc() table of contents.

Here is the working example:

#include "fvCFD.H"

#include "cellSet.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Main program:

int main(int argc, char *argv[])
{
    #include "setRootCase.H"
    #include "createTime.H"
    #include "createMesh.H"

    // The cellSet requires a simulation case, since all constructors take an
    // IOobject as the argument.
    cellSet testCells(
        IOobject
        (
            "testCells", 
            "constant/polyMesh/sets",
            mesh, 
            IOobject::NO_READ, 
            IOobject::AUTO_WRITE
        )
    ); 

    // HashTable<-HashSet<-cellSet, so there is no need to first create the 
    // cell set using the 'setSet' app, you can do it for testing purposes manually.
    testCells.insert(0);
    testCells.insert(1);
    testCells.insert(3);
    testCells.insert(4);
    testCells.insert(5);
    testCells.insert(6);
    testCells.insert(7); // Note the repeating ID.
    testCells.insert(7); // Note the repeating ID.
    testCells.insert(7); // Note the repeating ID.
    testCells.insert(34);

    // HashSet uses unique keys, so the cellSet is based on it, to prevent repeating
    // cell IDs. 

    // Write the cellSet using the regIOobject inheritance path. 
    testCells.write(); 

    // Since a cellSet is a HashTable is a HashTableCore, use toc() to access keys.
    const labelList& cells = testCells.toc(); 

    // Here are your cell IDs. 
    forAll(cells, I)
    {
        Info << cells[I] << endl;
    }

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

    Info<< "\nEnd\n" << endl;
    return 0;
}


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

[dmoroian]

It works great!

Thanks!

[tomislav_maric]

@Dragos Np, I'm glad it helped!