[petsc-users] nonlinear eigenvalue problem with structured mesh

[Barry Smith]

> On Jul 25, 2015, at 3:47 PM, Matthew Knepley <knepley at gmail.com> wrote:
> 
> On Sat, Jul 25, 2015 at 2:56 PM, Gideon Simpson <gideon.simpson at gmail.com> wrote:
> Are there any tools for constructing the Jacobian of a function when using DMComposite and DMRedundant?
> 
> No good ones right now, since it would be a sparse Jacobian + a low rank update. We have talked about putting
> this in forever, but no one was screaming for it.

  Tao has some code for low rank update matrices MatLMVM but it is unfortunately a bit to tied to its current usage and eventually needs refactoring for general usability.

   I would start with a MPIAIJ matrix for the full problem and just fill it up and use it. This will allow you to get a code that actually solves the eigenvalue you problem you want to solve. Premature optimization is the root of all evil in HPC.

   Barry

  If you end up needing to solve huge ass problems and this part becomes the bottle neck we can help you improve the performance by treating it in a more sophisticated way but it is crazy to delay your science when you don't even know if  it will be a problem. And like most of these types of optimizations it would not require you rewriting your code totally differently later.

> 
>    Matt
>  
> -gideon
> 
>> On Jul 13, 2015, at 8:25 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>> 
>> 
>>  You can use a DMCOMPOSITE with a DMDA and a DMREDUNDANT see src/snes/examples/tutorials/ex21.c  It is intended for this type of situation.
>> 
>>  Barry
>> 
>> 
>>> On Jul 13, 2015, at 4:05 PM, Gideon Simpson <gideon.simpson at gmail.com> wrote:
>>> 
>>> Thanks for the suggestions.  Suppose though that I wanted to avoid using a self-consistent approach and just attack the full nonlinear problem.  Is there a “smart” way to handle my data where the Vec q, which stores both the function on the mesh, and the eigenvalue parameter.  Suppose the eigenvalue parameter was in the last entry of the q vector, and the N-1 entries represnted the solution.  Is there a way to mask this sub-vector so it can be manipulated with a da?
>>> 
>>> -gideon
>>> 
>>>> On Jul 7, 2015, at 10:38 AM, Matthew Knepley <knepley at gmail.com> wrote:
>>>> 
>>>> On Tue, Jul 7, 2015 at 2:06 AM, Jose E. Roman <jroman at dsic.upv.es> wrote:
>>>> 
>>>> El 07/07/2015, a las 03:34, Matthew Knepley escribió:
>>>> 
>>>>> On Mon, Jul 6, 2015 at 7:43 PM, Gideon Simpson <gideon.simpson at gmail.com> wrote:
>>>>> I have a nonlinear eigenvalue problem for a system of equations, of the form,
>>>>> 
>>>>> -Delta u + f(u) = E * u,
>>>>> 
>>>>> where E is my nonlinear eigenvalue parameter, and u and f are vector valued.
>>>>> 
>>>>> I thus have the following two things to contend with:
>>>>> 
>>>>> 1.  E is a scalar which needs to be distributed across all the processes when the right hand side is formed
>>>>> 
>>>>> 2.  I would like to be able to use a da to manage the spatial and multicomponent nature of u
>>>>> 
>>>>> Obviously the Vec that my nonlinear solver is going to search for has to store both bits of data.  Is there a clever petsc way to handle this, or will I need to do all the indexing and broadcasting by hand?
>>>>> 
>>>>>> My first suggestion would be to investigate SLEPc, which does have some support for nonlinear
>>>>>> eigenvalue problems. Failing that, E normally comes out of the algorithm as the result of some
>>>>>> vector algebra which automatically distributes the constant (like VecDot).
>>>> 
>>>> Nonlinear eigenvalue problems can be generally classified in two groups: those that depend nonlinearly on the eigenvalue parameter T(lambda)*x=0, and those that are nonlinear with respect to the eigenvector H(X)*X=X*Lambda. It seems that Gideon's problem belongs to the seocond type.
>>>> 
>>>> Currently, SLEPc provides solvers for the first type only. We do not foresee to have solvers for the latter type anytime soon.
>>>> 
>>>> I guess what Gideon needs to do is compute E from an initial guess of u, then use SNES to update u, and iterate until a self-consistent state is reached.
>>>> 
>>>> Thanks for the clarification, Jose. I think what you could do is to use SLEPc as the subsolver for each
>>>> eigenproblem with frozen coefficients (so this would be a type of 'Picard' method according to Jed or
>>>> an Inexact Newton Method according to me) and then update using SNES, exactly as Jose says above.
>>>> Even if you do not know the actual Jacobian, you can do the obvious thing in your FormJacobian, which
>>>> is to just reassemble the matrix with the latest 'u' and you will get linear convergence if its contractive. This
>>>> also allows you to try out the nice line searches.
>>>> 
>>>>  Thanks,
>>>> 
>>>>     Matt
>>>> 
>>>> Jose
>>>> 
>>>> 
>>>>>> 
>>>>>>  Thanks,
>>>>>> 
>>>>>>     Matt
>>>>>> 
>>>>> -gideon
>>>>> 
>>>>> 
>>>>> 
>>>>> 
>>>>> --
>>>>> 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
>>> 
>> 
> 
> 
> 
> 
> -- 
> 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