Unable to perform a nonlinear solve with SNES when using a MatShell for J (the Jacobian matrix)

[GoogleCodeExporter]

The problem is that it is quite necessary to sometimes perform a nonlinear
solve of F(x) = b, when only the Matvec (matrix-vector product) of J*v is
known. This is quite different than the "Matrix-Free" approach.  Even
though they sound similar, MF involves a finite difference approximation of
 J. Anyhow, I have tried for a considerably long time to get this to work
and I have no idea how.  The standard approach of solving F(x) = b does
indeed work for me, but is not what I need.

What steps will reproduce the problem?
Create a class MyShell that inherits all the methods of MatShell, but
includes a new definition of mult(x,y).  Then create an instance, J, of
MyShell. Use this J to be the J that provides the desired matvec without
having to actually explicitly have a jacobian matrix. 

What is the expected output? What do you see instead?
The expected output is a solved problem. However, in the call to the
jacobian function needed by SNES, this J shell somehow gets reassigned to
be a Mat and loses the mult function that we gave it.  The code is attached
so that you can have a look.

What version of the product are you using? On what operating system?
petsc4py 0.7.5, PETSc 2.3.3, Python 2.5, Suse linux 10.2

Please provide any additional information below.
Please, there may be errors or stupid approaches to doing things in the
attached file, but I'm new to this game, so go easy on me ;)

Thanks in advance!

Original issue reported on code.google.com by frank.la...@gmail.com on 12 Jan 2008 at 5:50

Attachments:

• SNES_Shell_Test.py

[GoogleCodeExporter]

Here you have example script for solving x**2=4 using a shell matrix. Run it 
with:

$ python SNES_MatShell.py -snes_monitor -snes_view

Basically, your previous code had two problems:

1 - DO NOT inherit from MatShell class. This is counter-intuitive, but I've not 
found a easy way to support this. Instead, you have to write your class, next 
create 
an instance of this class, and then set that instance inside your PETSc shell 
matrix 
(or create a shell matrix and pass the instance).

2 - Your matvec computation was deadly wrong. If your nonlinear problem is 
x**2=4, 
then J*v = 2*x_i*v, where 'x_i' is the current nonlinear iterate, but 'v' is a 
given 
vector (it generally changes during the inner linear solver loop). So you have 
to 
save some way the current nonlinear iterate in your Jacobian routine in order 
to 
properly define your matvec operation.

Original comment by dalcinl on 14 Jan 2008 at 7:13

• Changed state: Accepted

Attachments:

• SNES_MatShell.py

[GoogleCodeExporter]

Terrific! Thank you so much!

Yes, I see that the jacobian matvec was wrong. I was more focused on getting the
programming to work out at the time, I guess :)

For anyone else who may care: 
You may need to change the '...' expressions seen in dalcinl's example code 
with the
':' to run.  I needed this, no problem though.

Thanks again!
Frank

Original comment by frank.la...@gmail.com on 14 Jan 2008 at 7:23

[GoogleCodeExporter]

One last thing -- also should be snes.solve(x,b) not (b,x). 

Original comment by frank.la...@gmail.com on 14 Jan 2008 at 9:31

[GoogleCodeExporter]

A note about last comment:

- In version 0.7.5 and below, use snes.solve(x,b)

- In current SVN and future versions, use snes.solve(b,x)

Original comment by dalcinl on 14 Jan 2008 at 9:48

• Changed state: Fixed

[GoogleCodeExporter]

Hi!,

In the current version of petsc4py, what would be the changes to make it work?

Original comment by romveg...@gmail.com on 20 Apr 2014 at 9:21

[GoogleCodeExporter]

See attached, tested with current petsc4py 
(https://bitbucket.org/petsc/petsc4py).

Original comment by five...@gmail.com on 21 Apr 2014 at 12:52

Attachments:

• SNES_MatShellv2.py