[Nek5000-users] Clarification on Adjoint solver
nek5000-users at lists.mcs.anl.gov
nek5000-users at lists.mcs.anl.gov
Wed Feb 17 14:37:21 CST 2016
Hi Philipp ,
Many thanks for your replies. I think what we have in mind is very
similar to what has been done in the paper you provided for the JCF,
with one major difference. In our problem an optimal disturbance with a
finite amplitude is sought after, which results in a time-dependent U_b.
>From previous comments looks like it can't be done automatically in Nek
(and the user needs to implement this manually). Does this manual way
mean the user has to load velocity at each time step (say in a loop that
executes the Nek solver)? Is there any advice on how can one implement
the classical check-pointing approach?
My second question is boundary condition. From Nek literature, it seems
line "inhomogeneous" boundary conditions can't be implemented (which is
usually the case for the outlet). Do you confirm that? In JCF paper you
provided I see $P_{a}n+1/Re \nabla u_{a}.n=(U_{b}.n)u_{a}$ with u_{a}
and p_{a} as adjoint variables is considered as inhomogeneous.
Neglecting the pressure term, this BC more looks like a Robin type to
me. So what's the jargon here for inhomogeneous BC?
Last question, is equation (7) of that paper the kind of equation solved
by Nek? Or one should change the terms slightly in the source code?
Many thanks,
Saleh
On 2/10/2016 4:10 PM, nek5000-users at lists.mcs.anl.gov wrote:
> Dear all,
>
> Let me add some answers below.
>
> On 2016-02-10 18:58, nek5000-users at lists.mcs.anl.gov wrote:
>> Adam and Oana,
>> Thanks again for your responses.
>> Let me summarize what I understood.
>> 1) We may have either evolving (unsteady) or non-evolving (steady) base
>> flow. As for the latter, we may use restart option in rea file.
>
> The repo version is for linear adjoint, i.e. for a non-evolving base
> flows.
>
>> 2) For non-evolving case (which is what I have in mind in the end) Nek
>> "is able to" handle and to do so parameter #31 in rea file needs to be
>> modified. Can you please clarify how. Suppose that there is velocity and
>> a passive scalar or even temperature, so there are three adjoint fields
>> (velocity, pressure and the scalar).
>
> The repo version only has the adjoint equations implemented for the
> velocity. The temperature adjoint has first been implemented at ANL,
> and we have further validated it. We hope that we can commit it back
> to the repo within a few weeks.
> Parameter 31 gives the number of disturbance fields, negative for a
> fixed base flow. For typical adjoint runs it should be -1 (one adjoint
> variable for fixed base flow).
>
>> Also, is what Oana mentioned, "one has to do the change of the
>> velocities himself" related to remark number 2 above?
>
> Exactly, for a changing base flow, one had to do that manually.
>
>> 3) There's also a confusion in userchk: First Adam mentioned that we
>> need to include 'ADJOINT' plus if(ISTEP.eq.0) ADJOINT=.TRUE. for adjoint
>> solver in order to i) changing the time and pressure sign and 2)
>> modification of advection terms. But in the second example, while he
>> still uses include 'ADJOINT' in userchk but then there is this statement
>> IFADJ = .FALSE. Can you please clarify what's the difference here?
>
> It was a typo, the relevant variable is IFADJ and nothing else. This
> one had to be true or false depending on the equation to be solved.
>
>> I'm really grateful that you take the time to answer these questions.
>> Please let me know if we can take this conversation off the mailing
>> list.
>>
>> Brief summary of the problem we have in mind: we'd like to optimize the
>> flow with respect to some parameters (either BCs or ICs) in terms of a
>> defined cost function. So this can be regarded as either design or
>> optimal control problem. Calculus of variations and Lagrange multipliers
>> lead to adjoint equations which are very similar to NS equations. We
>> already solved the problem using NEK and we'd like to use this again for
>> adjoint equations. All the results will be in an outer loop for
>> optimization (based on, say, steepest gradient). In general, the problem
>> is unsteady, but we can start off with steady and non-evolving cases at
>> early stages.
>
> We have calculated optimal initial conditions etc. using Nek, see e.g.
> the results of the jet in crossflow,
> http://www.sciencedirect.com/science/article/pii/S0997754614001009. Is
> your intended optimisation similar to that?
>
> Best regards,
> Philipp
>
>
>
>>
>> Best,
>> Saleh
>>
>> On 2/8/2016 4:26 AM, nek5000-users at lists.mcs.anl.gov wrote:
>>> Hi Saleh
>>>
>>> If I understand you correctly, you ask about a base flow, i.e. the
>>> velocity field you do your linearisation about (the one you use to
>>> advect prturbation). If this is the case, than it is possible to do
>>> with nek. Otherwise you have to wait for the code Oana mentioned. What
>>> we do in our stability calculations, we take the solution of
>>> non-linear NS (steady state reached with e.g. SFD) and next we load it
>>> into nek as non-linear filed. One can do it by setting restart options
>>> in ###.rea file. It should look like:
>>> 3.33333 3.33333 -0.833333 -1.16667
>>> XFAC,YFAC,XZERO,YZERO
>>> **MESH DATA** 1st line is X of corner 1,2,3,4. 2nd line is Y.
>>> -1472 2 1472 NELT,NDIM,NELV
>>> 1 PRESOLVE/RESTART OPTIONS *****
>>> bfext_cyl?.f00000
>>> 7 INITIAL CONDITIONS *****
>>> C Default
>>> C Default
>>>
>>>
>>> Where bfext_cyl0.f00000 is our restart file (we usually save this
>>> velocity field in double precision). Next you have to set linear
>>> solver with non-evolving base flow (nek is able to evolve both so be
>>> careful). You freeze baseflow by setting negative value of parameter
>>> 31 in ###.rea:
>>> 0.00000 p030 > 0 ==> properties set in uservp()
>>> -1.00000 p031 NPERT: #perturbation modes
>>> 0.00000 p032 #BCs in re2 file, if > 0
>>>
>>> Nek is able of evolving number of perturbation fields set by
>>> abs(param(31)) (in ###.rea) and lpert in SIZE. I have to mention that
>>> nek by itself outputs the non-linear fields only, so to see
>>> perturbation you should add to your userchk something like this
>>>
>>>
>>> subroutine userchk
>>>
>>> include 'SIZE' ! NX1, NY1, NZ1, NELV, NIO
>>> include 'TSTEP' ! ISTEP, IOSTEP
>>> include 'SOLN' ! V[XYZ]P, TP
>>> include 'INPUT' ! IF3D, PARAM
>>> include 'ADJOINT' ! IFADJ
>>>
>>> ! direct or adjoint run
>>> IFADJ = .FALSE.
>>>
>>> ! write perturbation field
>>> if (mod(ISTEP,int(PARAM(15))).eq.0) then
>>> call outpost2(VXP,VYP,VZP,PRP,TP,0,'prt')
>>> endif
>>>
>>> I hope it is an answer to your question.
>>> regards
>>> adam
>>>
>>> On 2016-02-06 19:15, nek5000-users at lists.mcs.anl.gov wrote:
>>>> Hi Saleh,
>>>>
>>>> There is no easily accessible implementation for this, one has to do
>>>> the change of the velocities himself. In fact we have ongoing work on
>>>> this and once we are done I guess we can share it.
>>>> What kind of case are you looking at, long time horizons are quite
>>>> tricky to handle. We can take this conversation off the mailing list
>>>> if needed.
>>>>
>>>> Oana
>>>>
>>>>
>>>>
>>>> ________________________________________
>>>> From: nek5000-users-bounces at lists.mcs.anl.gov
>>>> [nek5000-users-bounces at lists.mcs.anl.gov] on behalf of
>>>> nek5000-users at lists.mcs.anl.gov [nek5000-users at lists.mcs.anl.gov]
>>>> Sent: Friday, February 05, 2016 2:15 PM
>>>> To: nek5000-users at lists.mcs.anl.gov
>>>> Subject: Re: [Nek5000-users] Clarification on Adjoint solver
>>>>
>>>> Adams,
>>>>
>>>> I really appreciate your detailed answer. I definitely clarified
>>>> the Nek
>>>> functionality for adjoint solvers. One question remains is how to
>>>> include the velocity of the forward solver into adjoint equations. I
>>>> assume this is also the case when you do linear stability analysis
>>>> with
>>>> Nek, since the velocity in convective terms is not completely part of
>>>> the solution. Mathematically, my question is:
>>>>
>>>> in $V. \nabla u_a$ term, as the advection, how can we substitute
>>>> $V$ (as
>>>> the forward velocity which is already solved for in a regular NS
>>>> solution) in adjoint equations. suppose that the velocity, along with
>>>> other quantities are saved in afield out put such as fld.f0000X for X
>>>> time-step. So we have access to that velocity but how can we
>>>> include it
>>>> in adjoint equations?
>>>>
>>>> Many thanks again,
>>>> Saleh
>>>>
>>>> On 2/2/2016 3:09 AM, nek5000-users at lists.mcs.anl.gov wrote:
>>>>> Hi
>>>>>
>>>>> The adjoint implementation in nek assumes change of the direction of
>>>>> time evolution. The goal is to use the standard nek solver with
>>>>> minimal modifications (only advection terms have to be changed). You
>>>>> get these equations from you standard one by changing the time and
>>>>> pressure sign (t -> -t; p-> -p). This of course changes your time
>>>>> integration limits, but it shouldn't be a problem. So when you have
>>>>> your linear direct simulation, to turn it into dual one it is enough
>>>>> to include in your userchk include files ADJOINT and TSTEP, and for
>>>>> ISTEP.eq.0 set IFADJOINT to true. Something like:
>>>>> subroutine userchk
>>>>>
>>>>> include 'SIZE'
>>>>> include 'TSTEP' ! ISTEP
>>>>> include 'ADJOINT' ! IFADJ
>>>>>
>>>>> if(ISTEP.eq.0) ADJOINT=.TRUE.
>>>>>
>>>>> There is no need to modify ###.rea or ###.map files. The only problem
>>>>> are boundary conditions for open flows, as they would require
>>>>> non-homogeneous bc. In this case to avoid problems we usually use
>>>>> sufficiently large domain together with zero Dirichlet bc. The last
>>>>> issue is sufficient resolution for both direct and adjoin runs. I
>>>>> hope
>>>>> this gives you some idea about nek implementation.
>>>>> Regards
>>>>> Adam
>>>>>
>>>>> On 2016-02-01 22:17, nek5000-users at lists.mcs.anl.gov wrote:
>>>>>> Hi Neks,
>>>>>>
>>>>>> To my knowledge, Nek5000 is able to solve adjoint equations of the
>>>>>> form:
>>>>>> $\partial u_a/\partial t + V.\nabla u_a - u_a.(\nabla V)^T + \nabla
>>>>>> \p_a + 1/Re \nabla^2 u_a=0$ and $\nabla.u_a=0$
>>>>>> where u_a and p_a are adjoint velocity and pressure. (Slightly
>>>>>> different formulation may be seen in literature).
>>>>>>
>>>>>> Let's assume the forward problem is already solved so that $V$, i.e.
>>>>>> velocity, (and possibly $p$ as the pressure) is (are) already know.
>>>>>> Thus, the initial and boundary conditions for adjoint NS are also
>>>>>> known. How can we now solve the dual/adjoint NS problem with
>>>>>> Nek5000?
>>>>>> Is there any example on how to modify .rea files? Specailly, the
>>>>>> "convection term of $V.\nabla u_a - u_a.(\nabla V)^T$ needs to be
>>>>>> modified compared to forward NS but am not sure how.
>>>>>> My search of mailing list entails some modifications in "perturb.f"
>>>>>> subroutine, but a little bit more clarification on that would be
>>>>>> really appreciated.
>>>>>>
>>>>>> Yours,
>>>>>> Saleh
>>>>>> _______________________________________________
>>>>>> Nek5000-users mailing list
>>>>>> Nek5000-users at lists.mcs.anl.gov
>>>>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>>>> _______________________________________________
>>>>> Nek5000-users mailing list
>>>>> Nek5000-users at lists.mcs.anl.gov
>>>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>>> _______________________________________________
>>>> Nek5000-users mailing list
>>>> Nek5000-users at lists.mcs.anl.gov
>>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>>> _______________________________________________
>>>> Nek5000-users mailing list
>>>> Nek5000-users at lists.mcs.anl.gov
>>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>>
>>> _______________________________________________
>>> Nek5000-users mailing list
>>> Nek5000-users at lists.mcs.anl.gov
>>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
>>
>> _______________________________________________
>> Nek5000-users mailing list
>> Nek5000-users at lists.mcs.anl.gov
>> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
> _______________________________________________
> Nek5000-users mailing list
> Nek5000-users at lists.mcs.anl.gov
> https://lists.mcs.anl.gov/mailman/listinfo/nek5000-users
More information about the Nek5000-users
mailing list