[FlyBob91]

Hello to all,
i have two questions about the adjoint shape optimization.

I try to add this equation as boundary condition for the outlet in the file adjointOutletPressureFvPatchScalarField.C but when i compile it gives me an error

the equation is

  operator == ((Up & Uap) + nuw*snGradUan +(phiap/patch().magSf())*(phip/patch().magSf()));

and the error

adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C: In member function ‘virtual void Foam::adjointOutletPressureFvPatchScalarField::updateCoeffs()’:
adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C:107:31: error: ‘nuw’ was not declared in this scope
     operator == ((Up & Uap) + nuw*snGradUan +(phiap/patch().magSf())*(phip/patc
                               ^
adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C:107:35: error: ‘snGradUan’ was not declared in this scope
     operator == ((Up & Uap) + nuw*snGradUan +(phiap/patch().magSf())*(phip/patc
                                   ^

but i think the equation is correct because i took it from a scientific paper.

Second i want to ask the meaning of the terms Up, Uap, phiap and phip because they are not declared in the file.

Thanks for your help and time.

[rhythm]

Hi,

First of all a disclaimer that I am not fully familiar with adjoint optimisation formulation, nonetheless I believe I can give some input as to what is wrong with your implementation.

As I understand you have introduced an additional term to your operator. Unfortunately, as it says in the error message the variables nuw and snGradUan have not been declared in your code. To fix that you simply need to introduce them either in that *.C file or somewhere else in your code.

Second i want to ask the meaning of the terms Up, Uap, phiap and phip because they are not declared in the file.

Again my understanding of adjoint is shaky, so you should consult C. Othmer paper A continuous adjoint formulation for the computation of topological and surface sensitivities of ducted flows, 2008, and see the derivation procedure. (http://onlinelibrary.wiley.com/doi/1....1770/abstract)

Anyways, to my understanding from what is written in the code Up and phip are patch values for U and phi fields.

Similarly Uap and phiap are patch values for adjoint velocity and adjoint pressure variables. If my memory doesn't lie adjoint variables fall out in derivation when Lagrange multiplier method is applied.

All of these parameters are declared in the same updateCoeffs() function just above the operator declaration.


Hope this helps.
Cheers,
Ben

[FlyBob91]

Hi Ben,
First of all thank you for your reply.

I think i started the thread very bad, because i have noticed only now that in the .C file phi is renamed as phip, phia as phiap and so on. That's the reason of my second question and i ask sorry.

about the nuw parameter you have right i haven't defined it yet, but i'm more worried about the snGrad.
Now i want to be more precise.

The equation i want to implement is

q=u*v+un*vn+nu(n*grad)ut

where (v,q) are the lagrangian multipliers and (u,p) the original field variables.
I made this

operator == ((Up & Uap) + nu*snGradUt + (phiap/patch().magSf())*(phip/patch().magSf()))

and i added these lines based on a previous post

    scalarField Un(mag(patch().nf() & Up));         
    scalarField Uan(mag(patch().nf() & Uap));

    vectorField Ut(Up - patch().nf()*Un); 
    vectorField Uat(Uap - patch().nf()*Uan);

So as you said i have to declared nu which is my viscosity whose value is defined in the transportProperties file, but i don't know the right syntax. However i think it's not a very big problem, while i don't understand why it gives me error in snGradUt yet

adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C: In member function ‘virtual void Foam::adjointOutletPressureFvPatchScalarField::updateCoeffs()’:
adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C:116:31: error: ‘nu’ was not declared in this scope
     operator == ((Up & Uap) + nu*snGradUt + (phiap/patch().magSf())*(phip/patch
                               ^
adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C:116:34: error: ‘snGradUt’ was not declared in this scope
     operator == ((Up & Uap) + nu*snGradUt + (phiap/patch().magSf())*(phip/patch
                                  ^

tanks for help

[rhythm]

Hi,

Again, the error message that you are getting is the same as previously. You are asking the code to use variables that haven't been defined.

If you are providing viscosity values in your transportProperties dictionary, then you can use the following code written inside updateCoeffs() function to tell the class where to look for it:

scalar nu = readScalar(db().lookupObject<IOdictionary>("transportProperties").lookup("nu"));

Regarding the other part of your problem, I'm afraid I don't have an answer at the moment.

Couple of comments though:
1) First of all your use of C++ syntax is incorrect, because at the moment it is written as snGradUt which is interpreted as a variable and not a function snGrad() applied on a Ut field. Ideally you would like it to be written in the form

fvc::snGrad(Ut)

or (idea from https://github.com/OpenFOAM/OpenFOAM.../yPlus/yPlus.C)

Ut.boundaryField()[patchi].snGrad()

Both of these examples require Ut to be of a type volVectorField which in your case is not.

2) Your calculation of tangential velocity component is mathematically correct. Unfortunately, snGrad() function requires input of a GeometricField type that has both boundaryField and an internalField, as it needs to know the field values away from the boundary in order to be able to calculate surface normal gradient. One solution might be to use

volVectorField& U = db().lookupObject<volVectorField>("U")

Using it you would need to figure out how to calculate tangential values of the field with respect to the patch and then calculate snGrad() value at the patch. (though, I'm not fully sure if it would a correct approach mathematically)

Anyways here's an example of how someone else extracted surface normal gradients of the velocity field (Taken from http://www.cfd-online.com/Forums/ope...chfield-u.html):

// in the custom boundary condition of U
const volVectorField& Vs = this->db().lookupObject<volVectorField>("Vs");
const surfaceVectorField snGradVs = fvc::snGrad(Vs);
const label& patchIndex(this->patch().index());
const vectorField& patchSnGradVs = snGradVs.boundaryField()[patchIndex];

I know this doesn't solve your problem, so I guess it's more like food for thought rather than a proper answer.

Cheers,
Ben

[FlyBob91]

Hi Ben,
thanks for your reply; your suggestions were really helpfull.

I solved the question about the normal gradient. As you said the command snGrad can't be used and it's necessary to apply a first order discretization

(n*nabla)u_n= (u_nw-u_nw-1)/h

that i implemented in this manner

    // DISTANCE^(-1)
    const scalarField& deltainv = patch().deltaCoeffs();

    // nu / h
    const scalarField& nd = - nuw * deltainv;

    // NORMAL VECTOR
    vectorField n = patch().nf();

    // TANGENT VECTOR
   // vectorField r = 

    // SCALAR NORMAL VELOCITY COMPONENT
    scalarField Un = Up & patch().nf();

    // PATCH-ADIACENT TANGENT ADJOINT VELOCITY (VECTOR)
    vectorField Uac   = Uap.patchInternalField();
    vectorField Uac_n = ( Uac & patch().nf() ) * patch().nf();
    vectorField Uac_t = Uac - Uac_n;

    // TANGENT VELOCITY COMPONENT (VECTOR)
    vectorField U_n = ( Up & patch().nf() ) * patch().nf();
    vectorField U_t = Up - U_n;


    vectorField :: operator = (
                              phiap*patch().Sf()/sqr(patch().magSf())
                              + (1/(Un + nd))
                              * (nd * Uac_t + 2 * (n ^ r) / (pow(mag(r^U_t),3) + SMALL))
                              );

based on this article

http://www.tfd.chalmers.se/~hani/kur...ortAdjoint.pdf


Now i have to implement the equation of the velocity

u_tw = 1/(v_n - nu/h)*(-nu/h*u_(tw-1) + 2*(n x r)/|| r x v_t ||^(3) )

i defined n as

vectorField n = patch().nf();

while r is the vector starting from the center of the outlet face. This term is the only i'm not able to define.

[YUGU]

Hi,

in this " operator == ((Up & Uap) + nuw*snGradUan +(phiap/patch().magSf())*(phip/patch().magSf()));" equation, what's the meaning of "operator ==" ?

thanks