[petsc-users] SLEPc/NEP for shell matrice T(lambda) and T'(lambda)

[Jose E. Roman]

I am getting an error in a different place than you. I started to debug, but don't have much time at the moment.
Can you try something? Comparing to ex21.c, I see that a difference that may be relevant is the MATOP_DUPLICATE operation. Can you try defining it for your A matrix?

Note: If you plan to use the NLEIGS solver, there is no need to define the derivative T' so you can skip the call to NEPSetJacobian().

Jose


> El 6 oct 2023, a las 0:37, Kenneth C Hall <kenneth.c.hall at duke.edu> escribió:
> 
> Hi all,
>  
> I have a very large eigenvalue problem of the form T(\lambda).x = 0. The eigenvalues appear in a complicated way, and I must use a matrix-free approach to compute the products T.x and T’.x.
>  
> I am trying to implement in SLEPc/NEP.  To get started, I have defined a much smaller and simpler system of the form
> A.x - \lambda x = 0 where A is a 10x10 matrix. This is of course a simple standard eigenvalue problem, but I am using it as a surrogate to understand how to use NEP.
>  
> I have set the problem up using shell matrices (as that is my ultimate goal).  The full code is attached, but here is a smaller snippet of code:
>  
> !.... Create matrix-free operators for A and B
>       PetscCall(MatCreateShell(PETSC_COMM_SELF,n,n,PETSC_DETERMINE,PETSC_DETERMINE, PETSC_NULL_INTEGER, A, ierr))
>       PetscCall(MatCreateShell(PETSC_COMM_SELF,n,n,PETSC_DETERMINE,PETSC_DETERMINE, PETSC_NULL_INTEGER, B, ierr))
>       PetscCall(MatShellSetOperation(A, MATOP_MULT, MatMult_A, ierr))
>       PetscCall(MatShellSetOperation(B, MATOP_MULT, MatMult_B, ierr))
>  
> !.... Create nonlinear eigensolver
>       PetscCall(NEPCreate(PETSC_COMM_SELF, nep, ierr))
>  
> !.... Set the problem type
>       PetscCall(NEPSetProblemType(nep, NEP_GENERAL, ierr))
> !
> !.... set the solver type
>       PetscCall(NEPSetType(nep, NEPNLEIGS, ierr))
> !
> !.... Set functions and Jacobians for NEP
>       PetscCall(NEPSetFunction(nep, A, A, MyNEPFunction, PETSC_NULL_INTEGER, ierr))
>       PetscCall(NEPSetJacobian(nep, B,    MyNEPJacobian, PETSC_NULL_INTEGER, ierr))
>  
> The code runs, calls MyNEPFunction and MatMult_A multiple times, sweeping over the prescribed RG range, but crashes before it ever calls MyNEPJacobian or MatMult_B.  The NEP viewer and error messages are attached.
>  
> Any help on getting this problem properly set up would be greatly appreciated.
>  
> Kenneth Hall
> ATTACHMENTS: 
> test_nep.f90
> code_output
>  
> <code_output><test_nep.f90>