<html xmlns:o="urn:schemas-microsoft-com:office:office" xmlns:w="urn:schemas-microsoft-com:office:word" xmlns:m="http://schemas.microsoft.com/office/2004/12/omml" xmlns="http://www.w3.org/TR/REC-html40">
<head>
<meta http-equiv="Content-Type" content="text/html; charset=Windows-1252">
<meta name="Generator" content="Microsoft Word 15 (filtered medium)">
<style><!--
/* Font Definitions */
@font-face
{font-family:"Cambria Math";
panose-1:2 4 5 3 5 4 6 3 2 4;}
@font-face
{font-family:Calibri;
panose-1:2 15 5 2 2 2 4 3 2 4;}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
{margin:0in;
font-size:10.0pt;
font-family:"Calibri",sans-serif;}
a:link, span.MsoHyperlink
{mso-style-priority:99;
color:blue;
text-decoration:underline;}
span.EmailStyle19
{mso-style-type:personal-reply;
font-family:"Calibri",sans-serif;
color:windowtext;}
.MsoChpDefault
{mso-style-type:export-only;
font-size:10.0pt;
mso-ligatures:none;}
@page WordSection1
{size:8.5in 11.0in;
margin:1.0in 1.0in 1.0in 1.0in;}
div.WordSection1
{page:WordSection1;}
--></style>
</head>
<body lang="EN-US" link="blue" vlink="purple" style="word-wrap:break-word">
<div class="WordSection1">
<p class="MsoNormal"><span style="font-size:11.0pt">Jose,<o:p></o:p></span></p>
<p class="MsoNormal"><span style="font-size:11.0pt"><o:p> </o:p></span></p>
<p class="MsoNormal"><span style="font-size:11.0pt">Thanks very much for your help with this. Greatly appreciated. I will look at the MR. Please let me know if you do get the Fortran example working.<o:p></o:p></span></p>
<p class="MsoNormal"><span style="font-size:11.0pt"><o:p> </o:p></span></p>
<p class="MsoNormal"><span style="font-size:11.0pt">Thanks, and best regards,<o:p></o:p></span></p>
<p class="MsoNormal"><span style="font-size:11.0pt">Kenneth<o:p></o:p></span></p>
<p class="MsoNormal"><span style="font-size:11.0pt"><o:p> </o:p></span></p>
<p class="MsoNormal"><span style="font-size:11.0pt"><o:p> </o:p></span></p>
<div id="mail-editor-reference-message-container">
<div>
<div style="border:none;border-top:solid #B5C4DF 1.0pt;padding:3.0pt 0in 0in 0in">
<p class="MsoNormal" style="margin-bottom:12.0pt"><b><span style="font-size:12.0pt;color:black">From:
</span></b><span style="font-size:12.0pt;color:black">Jose E. Roman <jroman@dsic.upv.es><br>
<b>Date: </b>Wednesday, October 11, 2023 at 2:41 AM<br>
<b>To: </b>Kenneth C Hall <kenneth.c.hall@duke.edu><br>
<b>Cc: </b>petsc-users@mcs.anl.gov <petsc-users@mcs.anl.gov><br>
<b>Subject: </b>Re: [petsc-users] SLEPc/NEP for shell matrice T(lambda) and T'(lambda)<o:p></o:p></span></p>
</div>
<div>
<p class="MsoNormal" style="margin-bottom:12.0pt"><span style="font-size:11.0pt">Kenneth,<br>
<br>
The MatDuplicate issue should be fixed in the following MR <a href="https://urldefense.com/v3/__https:/gitlab.com/petsc/petsc/-/merge_requests/6912__;!!OToaGQ!p1tu1lzpyqM4wU-3WRzXN9bH3sFnXjyJvwQZh4PQBG5GNgB472qfxKOASyjxsg23AUQGusU-HpzI855ViaFfRCI$">
https://urldefense.com/v3/__https://gitlab.com/petsc/petsc/-/merge_requests/6912__;!!OToaGQ!p1tu1lzpyqM4wU-3WRzXN9bH3sFnXjyJvwQZh4PQBG5GNgB472qfxKOASyjxsg23AUQGusU-HpzI855ViaFfRCI$</a>
<br>
<br>
Note that the NLEIGS solver internally uses MatDuplicate for creating multiple copies of the shell matrix, each one with its own value of lambda. Hence your implementation of the shell matrix is not appropriate, since you have a single global lambda within
the module. I have attempted to write a Fortran example that duplicates the lambda correctly (see the MR), but does not work yet.<br>
<br>
Jose<br>
<br>
<br>
> El 6 oct 2023, a las 22:28, Kenneth C Hall <kenneth.c.hall@duke.edu> escribió:<br>
> <br>
> Jose,<br>
> <br>
> Unfortunately, I was unable to implement the MATOP_DUPLICATE operation in fortran (and I do not know enough c to work in c). Here is the error message I get:<br>
> <br>
> [0]PETSC ERROR: #1 MatShellSetOperation_Fortran() at /Users/hall/Documents/Fortran_Codes/Packages/petsc/src/mat/impls/shell/ftn-custom/zshellf.c:283<br>
> [0]PETSC ERROR: #2 src/test_nep.f90:62<br>
> <br>
> When I look at zshellf.c, MATOP_DUPLICATE is not one of the supported operations. See below.<br>
> <br>
> Kenneth<br>
> <br>
> <br>
> /**
<br>
> * Subset of MatOperation that is supported by the Fortran wrappers.
<br>
> */<br>
> enum FortranMatOperation {<br>
> FORTRAN_MATOP_MULT = 0,<br>
> FORTRAN_MATOP_MULT_ADD = 1,<br>
> FORTRAN_MATOP_MULT_TRANSPOSE = 2,<br>
> FORTRAN_MATOP_MULT_TRANSPOSE_ADD = 3,<br>
> FORTRAN_MATOP_SOR = 4,<br>
> FORTRAN_MATOP_TRANSPOSE = 5,<br>
> FORTRAN_MATOP_GET_DIAGONAL = 6,<br>
> FORTRAN_MATOP_DIAGONAL_SCALE = 7,<br>
> FORTRAN_MATOP_ZERO_ENTRIES = 8,<br>
> FORTRAN_MATOP_AXPY = 9,<br>
> FORTRAN_MATOP_SHIFT = 10,<br>
> FORTRAN_MATOP_DIAGONAL_SET = 11,<br>
> FORTRAN_MATOP_DESTROY = 12,<br>
> FORTRAN_MATOP_VIEW = 13,<br>
> FORTRAN_MATOP_CREATE_VECS = 14,<br>
> FORTRAN_MATOP_GET_DIAGONAL_BLOCK = 15,<br>
> FORTRAN_MATOP_COPY = 16,<br>
> FORTRAN_MATOP_SCALE = 17,<br>
> FORTRAN_MATOP_SET_RANDOM = 18,<br>
> FORTRAN_MATOP_ASSEMBLY_BEGIN = 19,<br>
> FORTRAN_MATOP_ASSEMBLY_END = 20,<br>
> FORTRAN_MATOP_SIZE = 21<br>
> };<br>
> <br>
> <br>
> From: Jose E. Roman <jroman@dsic.upv.es><br>
> Date: Friday, October 6, 2023 at 7:01 AM<br>
> To: Kenneth C Hall <kenneth.c.hall@duke.edu><br>
> Cc: petsc-users@mcs.anl.gov <petsc-users@mcs.anl.gov><br>
> Subject: Re: [petsc-users] SLEPc/NEP for shell matrice T(lambda) and T'(lambda)<br>
> <br>
> I am getting an error in a different place than you. I started to debug, but don't have much time at the moment.<br>
> 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?<br>
> <br>
> 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().<br>
> <br>
> Jose<br>
> <br>
> <br>
> > El 6 oct 2023, a las 0:37, Kenneth C Hall <kenneth.c.hall@duke.edu> escribió:<br>
> > <br>
> > Hi all,<br>
> > <br>
> > 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.<br>
> > <br>
> > I am trying to implement in SLEPc/NEP. To get started, I have defined a much smaller and simpler system of the form<br>
> > 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.<br>
> > <br>
> > 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:<br>
> > <br>
> > !.... Create matrix-free operators for A and B<br>
> > PetscCall(MatCreateShell(PETSC_COMM_SELF,n,n,PETSC_DETERMINE,PETSC_DETERMINE, PETSC_NULL_INTEGER, A, ierr))<br>
> > PetscCall(MatCreateShell(PETSC_COMM_SELF,n,n,PETSC_DETERMINE,PETSC_DETERMINE, PETSC_NULL_INTEGER, B, ierr))<br>
> > PetscCall(MatShellSetOperation(A, MATOP_MULT, MatMult_A, ierr))<br>
> > PetscCall(MatShellSetOperation(B, MATOP_MULT, MatMult_B, ierr))<br>
> > <br>
> > !.... Create nonlinear eigensolver<br>
> > PetscCall(NEPCreate(PETSC_COMM_SELF, nep, ierr))<br>
> > <br>
> > !.... Set the problem type<br>
> > PetscCall(NEPSetProblemType(nep, NEP_GENERAL, ierr))<br>
> > !<br>
> > !.... set the solver type<br>
> > PetscCall(NEPSetType(nep, NEPNLEIGS, ierr))<br>
> > !<br>
> > !.... Set functions and Jacobians for NEP<br>
> > PetscCall(NEPSetFunction(nep, A, A, MyNEPFunction, PETSC_NULL_INTEGER, ierr))<br>
> > PetscCall(NEPSetJacobian(nep, B, MyNEPJacobian, PETSC_NULL_INTEGER, ierr))<br>
> > <br>
> > 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.<br>
> > <br>
> > Any help on getting this problem properly set up would be greatly appreciated.<br>
> > <br>
> > Kenneth Hall<br>
> > ATTACHMENTS: <br>
> > test_nep.f90<br>
> > code_output<br>
> > <br>
> > <code_output><test_nep.f90><br>
> <o:p></o:p></span></p>
</div>
</div>
</div>
</div>
</body>
</html>