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

Barry Smith bsmith at mcs.anl.gov
Thu Feb 13 22:07:38 CST 2014

```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.
>>>>>>
>>>>>>
>>>>>
>>>>
>>>>
>>>
>>
>>
>

```