[juggler]

Hello Everyone!

My current Project (using OpenFOAM 4.0) demands (for the Sake of Stability) the spatial Filtering of a scalar Field to remove features smaller than a specified length scale. The Values of this Field should be blurred by averaging them with the Values of Cells within a certain radius. My first attempt at this (as specified by the Paper I was reading: http://www.topopt.dtu.dk/files/TopOpt88.pdf) was doing a weighted Average over all Cells within a certain radius. My Implementation of involved a for-loop over all Cells with max() and gsum() functions. While this yielded some Results (see below), it proved to be unusably slow. I therefore wanted to replace this by something faster, preferably leveraging the PDE-Solving powers of OpenFOAM.

The Paper proposed an alternative Approach for filtering by solving a Reaction-Diffusion-Type Equation, which I wrote as follows:

RFilter = pow(rFilter/(2*::sqrt(3.0)),2.0);
filteredObjectSensitivity = ObjectSensitivity;

Info<< "rFilter = " << rFilter << nl
    << "RFilter = " << RFilter << nl << endl;
    
fvScalarMatrix PDEFilter
        (
            laplacianDimKiller*RFilter*fvm::laplacian(filteredObjectSensitivity)
            -filteredObjectSensitivity
        ==
            -ObjectSensitivity
        );
        
PDEFilter.solve();

filteredObjectSensitivity = filteredObjectSensitivity/(max(filteredObjectSensitivity).value()*0.000000001);

with

rFilter as the original Filter length scale (set by user, type scalar)
ObjectSensitivity as the original (unfiltered) volScalarField
filteredObjectSensitivity as the filtered volScalarField
laplacianDimKiller as a scalar of unity value to make FOAM stop complaining about dimensions

The last Line is to normalize values to a certain Range.

The Paper specifies BCs as Gradients normal to the Boundary to be zero, so I set zeroGradient BCs for the filteredObjectSensitivity field.

Entries in fvSolution and fvSchemes were:

    filteredObjectSensitivity
    {
        solver          GAMG;
        tolerance       1e-06;
        relTol          0.001;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 20;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

and

    laplacian(filteredObjectSensitivity) Gauss linear corrected;

To test this, I disabled the usual Input to the filter and made it load a test Field created with setFields, then ran the Code for one Iteration on my 2D-Testcase. The test Field looks like this:

unfiltered.png

Filtering this with rFilter set to 0.5 yields following output on the CLI:

rFilter = 0.5 [0 0 0 0 0 0 0] 0.5
RFilter = pow((0.5|3.4641),2) [0 0 0 0 0 0 0] 0.02083333333

GAMG:  Solving for filteredObjectSensitivity, Initial residual = 1, Final residual = 0.0008298633416, No Iterations 6

Resulting filteredObjectSensitivity looks like following:

filtered05.png

For Comparison, Result of old filter Approach:

oldfilter05.png

Filtering with rFilter set to 0.05 gives this:

rFilter = 0.05 [0 0 0 0 0 0 0] 0.05
RFilter = pow((0.05|3.4641),2) [0 0 0 0 0 0 0] 0.0002083333333

GAMG:  Solving for filteredObjectSensitivity, Initial residual = 1, Final residual = 0.0008298633416, No Iterations 6

Again, Result of PDE Filter:

filtered005.png

And old Filter using the same Parameters:

oldfilter005.png

As can be seen, the PDE filter Approach does not really care for the Values of rFilter, as the old Filter does. What can I do to make the Output of the PDE-Filter at least similar to the Output of the old Filter?
What I tried:
- Flipping signs in the Equation -> Does not do anything exept inverting the Solution
- Tried very large/small values for rFilter -> Same Results
- Setting one BC as fixedValue, value 0 -> Distorts solution, still no Impact of rFilter to filteredObjectSensitivity field

Now I am pretty much out of Ideas. Maybe I did not quite understand something?

Best Regards,

Juggler