[petsc-users] A matrix-free Jacobian for SNES using MatShell

[Barry Smith]

On Aug 31, 2010, at 3:16 PM, Jeremy Roberts wrote:

> Ah, I understand now what you mean. I have u in a module used by all subroutines that need u, including the main program where it is an argument of SNESSolve and the subroutine jac_shell.  I assumed having all routines use the same variable makes it the same everywhere at all times...maybe I'm just lucky with this particular problem.  Putting it as you do, I see why the user might not want to assume what goes into solve is consistent with (points to?) what enters FOO at any given call.

   It does not neccessarily match, you should not depend on it.

  Barry

> 
> Thanks,
> Jeremy
> 
> P.S.  I've attached the working example, where "x" is the unknown vector.  
> 
> On Tue, Aug 31, 2010 at 2:27 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
> 
> On Aug 31, 2010, at 2:17 PM, Jeremy Roberts wrote:
> 
>> To implement the safe method---i.e. copying the unknown vector u within FOO---what is proper generic form of FOO?  
>> 
>> I guess I'm not clear exactly how FOO and jac_shell (the J*x routine) are linked.  It appears they aren't, given I can use an empty subroutine, but how then is a local copy of u within FOO used in for J(u) in J(U)*x of jac_shell?
> 
>   Don't ask me. You wrote the code, I have no idea how you are accessing u.  But you need to compute J(u)*x so you need to use u.
> 
>   subroutine yourfoo(snes,u,A,B,flg,ierr)
>   common /mycommonblock/ v 
> 
>   call veccopy(u,v,ierr)
>   return
> 
>    and create the v upfront in your main program. Of course, you don't need to use common blocks but that is one way of maintaining this global variables.
> 
>    Barry
> 
>    
> 
>> 
>> 
>> On Tue, Aug 31, 2010 at 1:14 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>> 
>> On Aug 31, 2010, at 1:03 PM, Jeremy Roberts wrote:
>> 
>> >
>> > Hi Barry,
>> >
>> > Thanks for the quick response.
>> >
>> > I have the unknown vector u globally accessible and so forgo use of the user-defined context (I tried to get it to work---can it be used in Fortran?).
>> >
>> > I used
>> >
>> > call SNESSetJacobian( snes, Jshell, Jshell, FOO, PETSC_NULL_OBJECT, ierr )
>> >
>> > where FOO is an empty function having the correct argument types.  It seems to work now, and I get the exact same answer using an explicit Jacobian matrix and the shell version.  Is there a way to get around needing FOO at all?
>> >
>> 
>>   If it works great. But I how you know what the "unknown vector u " that is "globally accessible"  in SNES is? How do you know it is the same each time a new Jacobian is needed? The only save way is to copy the vector passed into the FOO funcion.
>> 
>>   Barry
>> 
>> 
>> 
>> > Jeremy
>> >
>> >
>> 
>> 
> 
> 
> <exJac.F><exJac_mod.F><exJac_mod2.F>

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20100831/9dcdaf30/attachment.htm>