STFS-TUDa / blastAMR

Load-balanced adaptive mesh refinement libraries from blastFoam ported to ESI OpenFOAM
Other
37 stars 12 forks source link

Create Protected Cells - e.g., for wall refinement #4

Closed JanGaertner closed 7 months ago

JanGaertner commented 7 months ago

For wall bounded flows it may be beneficial to design the boundary layer at the wall and prevent mesh changes in this region, e.g., if the wall distance of Y+>20 should be kept that the selected wall model is working correctly.

Within the codedErrorEstimator it is possible to implement this with a cellSet defined in the topoSetDict in the system folder. However, this requires adding #include "cellSet.H" to the errorEstimator.H file.

An example how this can be implemented is shown below:

code            #{

    // =======================================================================
    //                          IMPORTANT NOTES
    // =======================================================================
    // this coded error block does not use the scale/normalize function
    // The error_ field is normalized within this coded block and 
    // protected cells created over the topoSet function in OpenFOAM are 
    // excluded from refinement
    //
    // To create the set of protected cells use the topoSet function of 
    // OpenFOAM to create a set of cells e.g.,
    //
    // system/topoSetDict:
    // ...
    // actions
    // (
    //      {
    //          name    protectedCells;     // Important: DO NOT CHANGE THE NAME
    //          type    cellSet;
    //          action  new;
    //          source      cylinderToCell;
    //          p1         (0 0 0);
    //          p2         (0 0 0.0046);
    //          radius          0.1;
    //      }
    // );

    Info<< "---->! custom error estimator !<----" << endl;
    error_ = 0.0;
    /*
        End result must be (error_ is a volScalarField):
        - error_ == 1  if the cell needs to be refined
        - error_ == 0  if the cell is to be left alone
        - error_ == -1 if the vertices of the cell are to be unrefined
    */

    const auto& T = mesh_.lookupObject<volVectorField>("U");
    error_ == mag(fvc::grad(T)) * dimensionedScalar("one",dimTime,1.0);

    // =======================================================================
    // Enforce that no near wall cells are refined
    // =======================================================================

    // check if it is already in the databank
    if (!mesh_.foundObject<volScalarField>("error:protectedCells"))
    {
        Info << "Create initial error:protectedCells"<<endl;
        cellSet protectedCellsSet
        (
            IOobject
            (
                "protectedCells",
                mesh_.pointsInstance()/polyMesh::meshSubDir/"sets",
                mesh_,
                IOobject::MUST_READ,
                IOobject::NO_WRITE
            )
        );

        volScalarField* protectedCellsPtr = new volScalarField
        (
            IOobject
            (
                "error:protectedCells",
                mesh_.time().timeName(),
                mesh_,
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
            ),
            mesh_,
            dimensionedScalar("temp",dimless,0)
        );
        // Transfer ownership to object registry
        protectedCellsPtr->store();

        volScalarField& protectedCells = *protectedCellsPtr;

        // Set the internal values according to the cellSet
        forAll(mesh_.C(),celli)
        {
            if (protectedCellsSet[celli])
                protectedCells[celli] = 1.0;
            else
                protectedCells[celli] = 0.0;
        }
    }

    // =======================================================================

    /*
        Here is the deal:
        - T Gradient is not recomputed if already available - !!
        - T will form a diffusive interface, 50% and up of the maxGrad are refined
          40%-50% are left alone, and less than 40% are unrefined
        - Not well tested, so play with the values.
    */

    scalar maxGradT = gMax(error_);
    scalar minGradT = gMin(error_);

    lowerRefine_ = minGradT + 0.05*(maxGradT-minGradT); // orig: minGradT + 0.5*(maxGradT-minGradT);
    upperRefine_ =  GREAT;
    lowerUnrefine_ = minGradT + 0.01*(maxGradT-minGradT);// orig: minGradT + 0.4*(maxGradT-minGradT);
    upperUnrefine_ =  GREAT;

    // =======================================================================

    const volScalarField& protectedCells = mesh_.lookupObject<volScalarField>("error:protectedCells");

    // Do not call normalize but do it manually here:
    forAll(error_,celli)
    {
        if (protectedCells[celli] > 0.0)
            error_[celli] = 0.0;
        else if (error_[celli] < lowerUnrefine_)
            error_[celli] = -1.0;
        else if (error_[celli] > lowerRefine_)
            error_[celli] = 1.0;
        else
            error_[celli] = 0.0;
    }
    // error_.correctBoundaryConditions() is called after this code in case you
    // need to manipulate boundary values
    Info<< "---->! end    error estimator !<----" << endl;
#};
JanGaertner commented 7 months ago

Note to restart

When the simulation is restarted the topoSet has to be executed beforehand

FoamScience commented 7 months ago

Yes, we have something like this in the Foam-Extend version already, thanks for reminding me to add it to this library too 😄

FoamScience commented 7 months ago

0b6fe84 takes care of this feature. Instead of relying on cell sets (which would make it possible to operate based on "distance from a patch") I opted for marking boundary cells and their neighbors in layers. It's much cheaper for polyhedral meshes to do this.

This unit test case documents how the new keywords are supposed to be used.

Restarts should pose no problems because this approach does not rely on cell remapping...

The unit tests seem to do fine, so I'm happy with it; closing this issue don't hesitate to re-open if more issues appear