CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM > OpenFOAM Programming & Development

Adjoint Shape Optimization

Register Blogs Members List Search Today's Posts Mark Forums Read

Like Tree2Likes
  • 1 Post By rhythm
  • 1 Post By rhythm

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   November 3, 2016, 06:13
Arrow Adjoint Shape Optimization
  #1
Senior Member
 
FlyBob91's Avatar
 
Join Date: Mar 2016
Location: Bergamo
Posts: 157
Rep Power: 10
FlyBob91 is on a distinguished road
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

Code:
  operator == ((Up & Uap) + nuw*snGradUan +(phiap/patch().magSf())*(phip/patch().magSf()));
and the error

Code:
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.
Attached Files
File Type: c adjointOutletPressureFvPatchScalarField.C (3.7 KB, 11 views)
FlyBob91 is offline   Reply With Quote

Old   November 3, 2016, 08:35
Default
  #2
Member
 
Ben Jankauskas
Join Date: Jun 2013
Location: Exeter
Posts: 34
Rep Power: 12
rhythm is on a distinguished road
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.

Quote:
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 likes this.
rhythm is offline   Reply With Quote

Old   November 3, 2016, 10:06
Default
  #3
Senior Member
 
FlyBob91's Avatar
 
Join Date: Mar 2016
Location: Bergamo
Posts: 157
Rep Power: 10
FlyBob91 is on a distinguished road
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

Code:
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

Code:
operator == ((Up & Uap) + nu*snGradUt + (phiap/patch().magSf())*(phip/patch().magSf()))
and i added these lines based on a previous post

Code:
    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

Code:
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
FlyBob91 is offline   Reply With Quote

Old   November 3, 2016, 17:03
Default
  #4
Member
 
Ben Jankauskas
Join Date: Jun 2013
Location: Exeter
Posts: 34
Rep Power: 12
rhythm is on a distinguished road
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:

Code:
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
Code:
fvc::snGrad(Ut)
or (idea from https://github.com/OpenFOAM/OpenFOAM.../yPlus/yPlus.C)
Code:
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
Code:
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):
Code:
// 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 likes this.
rhythm is offline   Reply With Quote

Old   November 6, 2016, 15:34
Default
  #5
Senior Member
 
FlyBob91's Avatar
 
Join Date: Mar 2016
Location: Bergamo
Posts: 157
Rep Power: 10
FlyBob91 is on a distinguished road
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

Code:
(n*nabla)u_n= (u_nw-u_nw-1)/h
that i implemented in this manner

Code:
    // 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

Code:
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

Code:
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.

Last edited by FlyBob91; November 14, 2016 at 10:19.
FlyBob91 is offline   Reply With Quote

Old   May 24, 2019, 16:28
Default
  #6
New Member
 
liliu
Join Date: Nov 2018
Posts: 12
Rep Power: 7
YUGU is on a distinguished road
Hi,

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

thanks
YUGU is offline   Reply With Quote

Reply

Thread Tools Search this Thread
Search this Thread:

Advanced Search
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
adjoint Shape Optimization mechy OpenFOAM Running, Solving & CFD 11 June 19, 2017 09:05
SU2 Shape optimization questions AdriC SU2 Shape Design 6 January 27, 2016 04:25
hard to convergent with SU2 4.0 for Euler adjoint Xianguu SU2 1 July 14, 2015 17:09
Adjoint method for shape optimization Edris Main CFD Forum 6 November 1, 2012 08:10
SQP method for adjoint based aerodynamic shape optimization cfdbooks Main CFD Forum 2 May 26, 2009 10:40


All times are GMT -4. The time now is 12:42.