[petsc-users] Is it possible to use -snes_ksp_ew in the framework of KSP?

[Song Gao]

Dear Petsc users,

We are working on a matrix-free solver for NS equations. We first developed
the solver in the SNES framework. But because we want to take control of
the Newton iterations and reuse some of our existing code, we stop using
SNES and begin developing the matrix-free solver in KSP framework.

But we have an issue when using MFFD and VecSetValuesBlocked together.
My code looks like

call MatMFFDSetFunction(pet_mat_mf, flowsol_ng_mf,  ctx, ierpetsc)
call MatMFFDSetBase(pet_mat_mf, pet_solu_snes, PETSC_NULL_OBJECT, ierpetsc)

The vec pet_solu_snes are set to blocksize 5.

flowsol_ng_mf looks like:
subroutine flowsol_ng_mf ( ctxx, pet_solu_snesa, pet_rhs_snesa,  ierrpet )
...............
    call VecGetBlockSize(pet_rhs_snesa,i,ierpetsc)
    write(*,*) "pet_rhs_snesa",i
.............
    call VecSetValuesBlocked ( pet_rhs_snesa, ndperl, pirow2, rhsloc,
ADD_VALUES, ierpetsc )
..........
end

I'm worried about the block size of pet_rhs_snesa is not properly
initialized. So I checked the blocksize in the function evaluation.

The output with single CPU looks like:

  0 KSP preconditioned resid norm 7.818306841541e+00 true resid norm
9.619278462343e-03 ||r(i)||/||b|| 1.000000000000e+00
 pet_rhs_snesa           1
 pet_rhs_snesa           1
 pet_rhs_snesa           5
  1 KSP preconditioned resid norm 1.739193723080e+00 true resid norm
2.564460053330e-03 ||r(i)||/||b|| 2.665958848545e-01
 pet_rhs_snesa           1
 pet_rhs_snesa           5
  2 KSP preconditioned resid norm 3.816590845074e-01 true resid norm
8.222646506071e-04 ||r(i)||/||b|| 8.548090730777e-02
 pet_rhs_snesa           1
 pet_rhs_snesa           5

I think the block size is not initialized correctly for the first call,
which leads to the misuse of  VecSetValuesBlocked.

One way to circumvent is to create a temporary vec with correct blocksize,
assemble it with VecSetValuesBlocked and then copy back. But Is there other
solutions? Thank you very much.



Song



On Sat, Mar 8, 2014 at 8:51 PM, Matthew Knepley <knepley at gmail.com> wrote:

> On Sat, Mar 8, 2014 at 7:45 PM, Song Gao <song.gao2 at mail.mcgill.ca> wrote:
>
>> Thank you all.
>>
>> Most of the difficulties come from the current structure of the program.
>> Its size and "age" make any radical modification a challenge. The large
>> uses of global variables and deprecated "goto" statements call for a
>> revision, which however is unlikely to occur, against our better
>> judgement...
>>
>> That being said, the current tools provided byPETSc are sufficient,
>> however not necessarily convenient for our purposes. As a wishful thinking
>> we can say that the implementation of PETSc SNES features into existing
>> codes would be easier if the outer Newton-like loop could be managed by the
>> original code, out of the SNES context. On the other hand, we realize that
>> this might be in contrast with the requirements of PETSc itself.
>>
>
> I guess I wanted to know specifically why you cannot just dump all of the
> other things in the legacy code into SNESSetUpdate()
> or SNESSetConvergenceTest(). Was there something about the way we
> structured the hooks that made this impossible?
>
>   Thanks,
>
>      Matt
>
>
>> The application of the Eistenstat-Walker method together with a
>> Matrix-Free approach are of key importance. Thus, as kindly suggested by
>> Jed Brown, implementing our own version of EW scheme might turn out to be
>> the only way. We would be happy to provide more details if you think that
>> they might be helpful for the future developments of PETSc.
>>
>>
>> On Fri, Mar 7, 2014 at 2:39 PM, Matthew Knepley <knepley at gmail.com>wrote:
>>
>>> On Fri, Mar 7, 2014 at 1:21 PM, Song Gao <song.gao2 at mail.mcgill.ca>wrote:
>>>
>>>> Hello,
>>>>
>>>> We are working on a legacy codes which solves the NS equations. The
>>>> codes take control of newton iterations itself and use KSP to solve the
>>>> linear system.
>>>>
>>>> We modified the codes, used SNES to control the newton iteration and
>>>> changed the code to matrix free fashion. But the legacy codes did a lot of
>>>> other things between two newton iterations (such as output solution, update
>>>> variables....). I know we could use linesearchpostcheck but it is
>>>> difficulty to do that correctly. Therefore, we decide to go back to the KSP
>>>> framework but still use matrix free.
>>>>
>>>
>>> What makes it difficult?
>>>
>>>   Thanks,
>>>
>>>      Matt
>>>
>>>
>>>> When using SNES, we always use the runtime option -snes_ksp_ew, we
>>>> observe that for some test cases, the residual stalls without -snes_ksp_ew,
>>>> but converges with-snes_ksp_ew. So I'm thinking if it is possible to
>>>> use -snes_ksp_ew in KSP? Thanks in advance.
>>>>
>>>>
>>>>
>>>
>>>
>>> --
>>> What most experimenters take for granted before they begin their
>>> experiments is infinitely more interesting than any results to which their
>>> experiments lead.
>>> -- Norbert Wiener
>>>
>>
>>
>
>
> --
> What most experimenters take for granted before they begin their
> experiments is infinitely more interesting than any results to which their
> experiments lead.
> -- Norbert Wiener
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140429/a6c63dc1/attachment.html>