[petsc-users] Providing matrix vector product yourself for matrix free in SNES

[E.T.A.vanderWeide at utwente.nl]

Hi Matt,

Thanks, this is working. The only problem I still have is that I have to store much more than just a vector in context of MATSHELL in order to carry out the matrix vector multiplication. I basically need all the contents of an object of a certain class. I tried to store the this pointer of the object in the context of MATSHELL, but then I get an error. So it looks like it is not allowed to pass the this pointer to MatShellSetContext, is that correct?

My current workaround is to use a global variable. That works, but I don't like to use global variables. Is there a better way to make available a pointer to an object in MATSHELL?

Thanks again,

Edwin


-----Original Message-----
From: petsc-users-bounces at mcs.anl.gov on behalf of Matthew Knepley
Sent: Tue 12/20/2011 3:46 PM
To: PETSc users list
Subject: Re: [petsc-users] Providing matrix vector product yourself for matrix free in SNES
 
On Tue, Dec 20, 2011 at 7:24 AM, <E.T.A.vanderWeide at utwente.nl> wrote:

> **
>
> Hi,
>
> I was wondering if it is possible to provide PETSc the matrix vector
> product yourself for a matrix free computation inside the SNES environment.
> The reason for asking is that the truncation error in the standard finite
> difference approach may lead to some noisy behavior for our application. We
> have found a nice way to get rid off the truncation error completely by
> using dual numbers for the creation of the Jacobian matrices. The same
> technique can be used for the matrix vector product inside SNES, but this
> means that the user has to provide this product to PETSc.
>
> To be more specific. I would like to solve the nonlinear problem F(U) = 0
> via the SNES environment using the matrix free option. So similarly to the
> function SNESSetFunction, which sets the function evaluation routine, I'm
> seeking for a function that sets the function to compute the matrix vector
> product dF/dU times V, where U and V are inputs. Is such a function
> available or is there a different way to provide the matrix vector product
> yourself?
>
No need to use the matrix-free option. Just provide a MATSHELL (
http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateShell.html
)
as the Jacobian using SNESSetJacobian(). Then your FormJacobian() function
will get a vector 'x' as input. Update the context for the MATSHELL with
this vector in that function, which is all it would really do. Then when
the MatMult() is called on your shell matrix, compute the action using that
'x'
which is stored in the context.

   Matt


> Thanks,
>
> Edwin
> ---------------------------------------------------
> Edwin van der Weide
> Department of Mechanical Engineering
> University of Twente
> Enschede, the Netherlands
>



-- 
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/20111221/137d2de4/attachment.htm>