# [petsc-dev] [petsc-maint] SNESMF

panourg at mech.upatras.gr panourg at mech.upatras.gr
Thu Feb 13 22:26:34 CST 2014

```Mr.Matthew, yes,

It is something like you said. Basically, I augmented the diagonal
elements of the operator.

Also, Yes, I want to find the solution of n+1 time using the steady state
solution of non-linear system. That is to reply to the other guy.

Thanks a lot for your help

Kostas

>
> On Feb 13, 2014, at 9:44 PM, <panourg at mech.upatras.gr>
> <panourg at mech.upatras.gr> wrote:
>
>> Barry, I am trying to make a dual time implicit scheme where essentially
>> I
>> have a pseudo-time step which appears only in the Jacobian.
>
>    Hmm, so if this is a reasonable thing (Jed and Matt?) shouldn’t we add
> support for it for TSPSEUDO? This is also much like a trust region in
> that you are adding some scaling of the identity (or some other matrix)
> to the Jacobian, maybe.
>
>    Barry
>
> Request-assigned: who,   add support to TSPSEUDO for pseudo-time step
> which appears only in the Jacobian
>
>>
>> For example: So far I had the classical implicit scheme where I had to
>> solve at each time step the following non-linear scheme after
>> linearizing
>> it:
>>
>> [M/Dt + dR/dt ] * DU = -F   where [M/Dt + dR/dt ] = dF/dU.
>>
>>
>> Now I have resulted into a new linearized system using pseudo-time step:
>>
>> [M/Dtau +  M/Dt + dR/dt ] * DU = -F  where the F is the same as above.
>>
>> As you see if I differentiate the F function at RHS I will not take the
>> operator due to the fact that have an additional term (M/Dtau).
>>
>> Of course this problem is coming up due to the numerical
>> differentiation.
>> If I could do it analytically by hand I will have the same function for
>> both RHS and operator.
>>
>> Ok. From what you said to me petsc works as I need. Therefore I am ok
>> for
>> the time being.
>>
>> Kostas.
>>
>>
>>>   Just curious but why for your problem do you wish to use a different
>>> function for differencing for the Jacobian than that you want to find
>>> the zero of?
>>>
>>>   Barry
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>>
>>> On Feb 13, 2014, at 8:58 PM, <panourg at mech.upatras.gr>
>>> <panourg at mech.upatras.gr> wrote:
>>>
>>>> Thank you a lot. But, As far as I saw I cant do what I need with that
>>>> because I noticed after implementing it that the first function
>>>> (Non_linear_function_1 which has been set by SNESsetfunction) is
>>>> called
>>>> when petsc wants just to show me the non-linear residual. To be
>>>> understandable take a look at the following print:
>>>>
>>>> As you can see the Non_linear_function_1 is called only at the
>>>> beginning
>>>> of each newton step and not at any intermediate gmres stage.
>>>
>>>   Intermediate steps of GMRES never need to recompute the RHS b since
>>> it
>>> is stored and stays the same for the entire GMRES solve
>>>
>>>> And here I am
>>>> wondering whether when petsc goes to compute the RHS (b vector of
>>>> Ax=b)
>>>> of
>>>> the linear system at the k Newton iteration uses the
>>>> Non_linear_function_1
>>>
>>>  It is using 1
>>>
>>>> or the Non_linear_function_2.
>>>>
>>>> If petsc uses the Non_linear_function_1 whenever it needs a value of
>>>> (b)
>>>> then I am ok. Otherwise I cant use it.
>>>>
>>>> Basically, with that I want to solve a dual time implicit method for
>>>> compressible flows where the operator is computed by slighlty
>>>> different
>>>> function than the RHS. Just to know what I try to do.
>>>>
>>>>
>>>> Non_linear_function_1
>>>> Non_linear_function_1
>>>> Non_linear_function_1
>>>> Non_linear_function_1
>>>
>>>   Above 1 is computing the rhs
>>>
>>>> 0 SNES Function norm 1.942206120770e-02
>>>>   0 KSP unpreconditioned resid norm 1.942206120770e-02 true resid norm
>>>> 1.942206120770e-02 ||r(i)||/||b|| 1.000000000000e+00
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>>   1 KSP unpreconditioned resid norm 3.959771148164e-04 true resid norm
>>>> 3.959771148164e-04 ||r(i)||/||b|| 2.038800674047e-02
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>> Non_linear_function_2
>>>>   2 KSP unpreconditioned resid norm 8.488671743678e-06 true resid norm
>>>> 8.488680869120e-06 ||r(i)||/||b|| 4.370638511711e-04
>>>
>>> Above 2 is being used to do matrix vector products
>>>
>>>> Non_linear_function_1
>>>> Non_linear_function_1
>>>> Non_linear_function_1
>>>> Non_linear_function_1
>>>
>>>   Above 1 is being used for the line search.
>>>
>>>> 1 SNES Function norm 1.938487707733e-02
>>>
>>>   Based on the output above it is doing exactly what it is suppose to
>>> be
>>> doing.
>>>
>>>   Just curious but why for your problem do you wish to use a different
>>> function for differencing for the Jacobian than that you want to find
>>> the zero of?
>>>
>>>   Barry
>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>> Thanks and I am sorry for this confusion that maybe I have caused to
>>>> you.
>>>>
>>>> Kostas
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>>
>>>>>  Ok, it is now fixed and ready for your use. You will need to obtain
>>>>> PETSc with
>>>>>
>>>>> git clone https://bitbucket.org/petsc/petsc.git petsc
>>>>> cd petsc
>>>>> git branch barry/fix-setfunction-for-snesmf
>>>>>
>>>>> (if you are not currently getting PETSc from git)
>>>>>
>>>>> then configure and make it in the usual manner.
>>>>>
>>>>> You can then pass your ctx through the call to MatMFFDSetFunction()
>>>>> and
>>>>> not snes.
>>>>>
>>>>> See src/snes/examples/tests/ex7.c for a sample use of
>>>>> MatMFFDSetFunction()
>>>>> though your code should be ok.
>>>>>
>>>>> Any problems please let us know
>>>>>
>>>>> Barry
>>>>>
>>>>>
>>>>> On Feb 13, 2014, at 9:26 AM, <panourg at mech.upatras.gr>
>>>>> <panourg at mech.upatras.gr> wrote:
>>>>>
>>>>>> Ok I corrected the snes arguments but now how can I pass into the
>>>>>> Non_linear_function_2 some structs of mine which have some
>>>>>> information
>>>>>> that I need?
>>>>>>
>>>>>> when I use the SNESSetFunction(SNES snes,Vec r,PetscErrorCode
>>>>>> (*SNESFunction)(SNES,Vec,Vec,void*),void *ctx) ,  I pass my
>>>>>> information
>>>>>> through *ctx.
>>>>>>
>>>>>>
>>>>>> Thanks
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>>
>>>>>>> On Feb 13, 2014, at 12:23 AM, Jed Brown <jed at jedbrown.org> wrote:
>>>>>>>
>>>>>>>> panourg at mech.upatras.gr writes:
>>>>>>>>
>>>>>>>>> Hi.
>>>>>>>>> As you can see at the previous e-mail that I sent you and you
>>>>>>>>> replied,
>>>>>>>>> I
>>>>>>>>> need to define different nonlinear function for the computation
>>>>>>>>> of
>>>>>>>>> dF/du
>>>>>>>>> than I have set for the computation of F (RHS of eq (1)).
>>>>>>>>> Concisely, I use (apart from other petsc functions) the
>>>>>>>>> following:
>>>>>>>>>
>>>>>>>>> Mat Jmf;
>>>>>>>>>
>>>>>>>>> SNESSetFunction(snes_self,Solution1_self,Non_Linear_Residual,run);
>>>>>>>>> MatCreateSNESMF(snes,&Jmf);
>>>>>>>>> MatMFFDSetFunction(Jmf,Non_Linear_Residual_2,run);
>>>>>>>>> SNESSetJacobian(snes,Jmf,P,preconditioner,run)
>>>>>>>>
>>>>>>>
>>>>>>> Currently the final argument to MatMFFDSetFunction() and the first
>>>>>>> argument to your function Non_Linear_Residual_2 has to be the snes
>>>>>>> object. Sorry I did not explain this.
>>>>>>>
>>>>>>>  Please try again and if it doesn’t work send us the entire error
>>>>>>> message.
>>>>>>>
>>>>>>> Barry
>>>>>>>
>>>>>>>
>>>>>>>> Does your function "preconditioner" assemble only into P, then
>>>>>>>> call
>>>>>>>> MatMFFDSetBase(Jmf,X,NULL)?
>>>>>>>
>>>>>>> Jed,
>>>>>>> With  MatCreateSNESMF(), (as opposed to MatCreateMFFD()) the user
>>>>>>> need
>>>>>>> not call MatMFFDSetBase(), they only need to call
>>>>>>> MatAssemblyBegin/End()
>>>>>>>
>>>>>>>>
>>>>>>>>> and Neither do I put the option -snes_mf_operator nor -snes_mf.
>>>>>>>>> In MatMFFDSetFunction, I have set different function as you can
>>>>>>>>> see.
>>>>>>>>>
>>>>>>>>> But here I take an error
>>>>>>>>
>>>>>>>> Paraphrasing errors kills kittens!
>>>>>>>>
>>>>>>>> Always, always, always send the ENTIRE error message.
>>>>>>>
>>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>
>>>
>>>
>>
>
>

```