[Andreas M.]

Hi OpenFoam Community,

i want to implement an inlet boundary control for a given ducted flow fluid system.
For more context:
I use the given adjointShapeOptimizationFoam solver. I adapted the primal and
adjoint pde, so it fits my conditions.
I got/solved the Design equation and want to run a steepest descent algorithm like the adjointShapeOptimizationFoam solver does.

The issue:
Where do i implement this steepest descent algorithm/equation?
Is there a way to do it in the first place?
Can anyone give my some advice/papertips/books/whatever that treat such cases?

Thanks in Advance

[Andreas M.]

Since no one replied yet, i go a bit further in Detail:

I have to calculate something like:

alpha = -n*(U & Ua) - Ua*(U & n) - nu*(grad(U)n) + qn (<= Design equation)

Where n is the surface normal (vector) of the inlet.
U is the primal velocity-vectorfield,
Ua is the adjoint velocity-vectorfield,
( . & . ) is scalar/inner product,
nu is the kinematic viscocity
q is the adjoint pressure scalar field

alpha itself is the velocity vector of the inlet. The values of alpha are copied in the 0/U file before the simulation is started.

At the end of the simulation, the objective function value is calculated and given to a third party optimizer (Dakota) for further processing.

My issue:
I have to calculate a gradient of alpha, and return it to dakota as well.
In order to do so, i have to find a way, how i can calculate the above mentioned alpha.
Alpha itself has no influence what so ever on the current simulation.

At best OpenFoam would calculate alpha (or update alpha) in each iteration.

It would be sufficient, if i would be able to calculate the alpha at the end of the simulation (post processing).

Anyway, since the calculation has to be done for the inlet only and the solver im using is adapted to my concerns, im still clueless how to do it.

Still hoping someone can give me a hint.

[Andreas M.]

Hi again,

i digged a bit into it. Got some solutions and some (no clue, why they appear in the first place) problems:

Basically i added a new IOField into the adapted adjointShapeOptimizationFoamTest.C in order to calculate my needed equation.

I embedded the following line

int main(int argc, char *argv[])
{
label PatchID = mesh.boundaryMesh().findPatchID("inlet");

and in the simple.loop right after runtime.write the following code:

runTime.write();
    if (runTime.write())
    {
        IOField<vector> cfOut
        (
        IOobject
        (
            "cf",
            mesh.time().timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
        -mesh.Cf().boundaryField()[PatchID]*(U.boundaryField()[PatchID] & Ua.boundaryField()[PatchID]) - Ua.boundaryField()[PatchID]*(U.boundaryField()[PatchID] & mesh.Cf().boundaryField()[PatchID]) + pa.boundaryField()[PatchID] * mesh.Cf().boundaryField()[PatchID]
        );
        cfOut.write();
}

in adjointShapeOptimizationFoamTest.C

After this changes (and wclean&wmake) the solver calculates the cf vector (or better my needed and above mentioned alpha in code called cf) without any further problems (besides the unpleasant way of implementation).

In the next step i tried to implement the approximation of my derived objective function (which is an indicator function of the velocity field).
I changed the scalar field alpha into a velocity field in createFields.H and added several variables to get the dimensional units in place. Then I deleted every alpha in adjointShapeOptimizationFoamTest.C and added my line for the update of alpha:

In createFields.H:

volVectorField alpha
    (
        IOobject
        (
            "alpha",
            runTime.timeName(),
            mesh,
            IOobject::READ_IF_PRESENT,
            IOobject::AUTO_WRITE
        ),
        -100000.0*U*lambda2*1.0/(1.0*alphaMax + 100000.0*sqr((U & U) - lambda))//lambda*max(Ua & U, zeroSensitivity)
    );

in adjointShapeOptimizationTest.C, instead of the "alpha += ":

alpha = -100000.0*U*lambda2*1.0/(1.0*alphaMax + 100000.0*sqr((U & U) - lambda));

My errors so far:
if i keep it like that, my cf explodes in no time (values get so huge that it cant get past to 70-80 iterations, while without any changes to alpha it looks like a possible solution)

if i comment every alpha existing in the solver (CreateFields and adjointShapeOptimizationTest.C), my cf explodes as well and it calculates the time directories in the wrong order (seems to skip alot)

I added the working and not working files. Those in the folder are my working backup.

IF i get the solver working, i would try to add the alpha via fvm::Sp(1.0, alpha ) into Uaqn <fvVectorMatrix> since U and Ua are not compatible vector fields, but since my objective function has a velocity part inside the tank, it seems impossible to get the adjoint velocity in here!?

I would be glad, if someone could tell me, why the solver stops working when i change the alpha as mentioned, since it has nothing to do (codewise) with my cf or with the solution of the primal and adjoint equations (yet).

[Andreas M.]

Im using OpenFoam 2.4 for those who are interested.