<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:11.0pt;
        font-family:"Calibri",sans-serif;
        mso-ligatures:standardcontextual;}
span.EmailStyle17
        {mso-style-type:personal-compose;
        font-family:"Calibri",sans-serif;
        color:windowtext;}
.MsoChpDefault
        {mso-style-type:export-only;
        font-family:"Calibri",sans-serif;}
@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="#0563C1" vlink="#954F72" style="word-wrap:break-word">
<div class="WordSection1">
<p class="MsoNormal">Hi all,<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">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.<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">I am trying to implement in SLEPc/NEP.  To get started, I have defined a much smaller and simpler system of the form
<o:p></o:p></p>
<p class="MsoNormal">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.<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">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:<o:p></o:p></p>
<p class="MsoNormal"><span style="font-size:10.0pt"><o:p> </o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:#3F7F5F">!.... Create matrix-free operators for A and B</span><span style="font-size:10.0pt;color:black"><o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black">     
</span><span style="font-size:10.0pt;color:#0000C0">PetscCall</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">MatCreateShell</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">PETSC_COMM_SELF</span><span style="font-size:10.0pt;color:black">,</span><span style="font-size:10.0pt;color:#0000C0">n</span><span style="font-size:10.0pt;color:black">,</span><span style="font-size:10.0pt;color:#0000C0">n</span><span style="font-size:10.0pt;color:black">,</span><span style="font-size:10.0pt;color:#0000C0">PETSC_DETERMINE</span><span style="font-size:10.0pt;color:black">,</span><span style="font-size:10.0pt;color:#0000C0">PETSC_DETERMINE</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">PETSC_NULL_INTEGER</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">A</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">ierr</span><span style="font-size:10.0pt;color:black">))<o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black">     
</span><span style="font-size:10.0pt;color:#0000C0">PetscCall</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">MatCreateShell</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">PETSC_COMM_SELF</span><span style="font-size:10.0pt;color:black">,</span><span style="font-size:10.0pt;color:#0000C0">n</span><span style="font-size:10.0pt;color:black">,</span><span style="font-size:10.0pt;color:#0000C0">n</span><span style="font-size:10.0pt;color:black">,</span><span style="font-size:10.0pt;color:#0000C0">PETSC_DETERMINE</span><span style="font-size:10.0pt;color:black">,</span><span style="font-size:10.0pt;color:#0000C0">PETSC_DETERMINE</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">PETSC_NULL_INTEGER</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">B</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">ierr</span><span style="font-size:10.0pt;color:black">))<o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black">     
</span><span style="font-size:10.0pt;color:#0000C0">PetscCall</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">MatShellSetOperation</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">A</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">MATOP_MULT</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">MatMult_A</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">ierr</span><span style="font-size:10.0pt;color:black">))<o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black">     
</span><span style="font-size:10.0pt;color:#0000C0">PetscCall</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">MatShellSetOperation</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">B</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">MATOP_MULT</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">MatMult_B</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">ierr</span><span style="font-size:10.0pt;color:black">))<o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black"><o:p> </o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:#3F7F5F">!.... Create nonlinear eigensolver</span><span style="font-size:10.0pt;color:black"><o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black">     
</span><span style="font-size:10.0pt;color:#0000C0">PetscCall</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">NEPCreate</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">PETSC_COMM_SELF</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">nep</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">ierr</span><span style="font-size:10.0pt;color:black">))<o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black"><o:p> </o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:#3F7F5F">!.... Set the problem type</span><span style="font-size:10.0pt;color:black"><o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black">     
</span><span style="font-size:10.0pt;color:#0000C0">PetscCall</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">NEPSetProblemType</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">nep</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">NEP_GENERAL</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">ierr</span><span style="font-size:10.0pt;color:black">))<o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:#3F7F5F">!</span><span style="font-size:10.0pt;color:black"><o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:#3F7F5F">!.... set the solver type</span><span style="font-size:10.0pt;color:black"><o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black">     
</span><span style="font-size:10.0pt;color:#0000C0">PetscCall</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">NEPSetType</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">nep</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">NEPNLEIGS</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">ierr</span><span style="font-size:10.0pt;color:black">))<o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:#3F7F5F">!</span><span style="font-size:10.0pt;color:black"><o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:#3F7F5F">!.... Set functions and Jacobians for NEP</span><span style="font-size:10.0pt;color:black"><o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black">     
</span><span style="font-size:10.0pt;color:#0000C0">PetscCall</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">NEPSetFunction</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">nep</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">A</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">A</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">MyNEPFunction</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">PETSC_NULL_INTEGER</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">ierr</span><span style="font-size:10.0pt;color:black">))<o:p></o:p></span></p>
<p style="margin:0in;background:white"><span style="font-size:10.0pt;color:black">     
</span><span style="font-size:10.0pt;color:#0000C0">PetscCall</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">NEPSetJacobian</span><span style="font-size:10.0pt;color:black">(</span><span style="font-size:10.0pt;color:#0000C0">nep</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">B</span><span style="font-size:10.0pt;color:black">,   
</span><span style="font-size:10.0pt;color:#0000C0">MyNEPJacobian</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">PETSC_NULL_INTEGER</span><span style="font-size:10.0pt;color:black">,
</span><span style="font-size:10.0pt;color:#0000C0">ierr</span><span style="font-size:10.0pt;color:black">))<o:p></o:p></span></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">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.<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">Any help on getting this problem properly set up would be greatly appreciated.<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<p class="MsoNormal">Kenneth Hall<o:p></o:p></p>
<p class="MsoNormal">ATTACHMENTS: <o:p></o:p></p>
<p class="MsoNormal" style="text-indent:.5in">test_nep.f90<o:p></o:p></p>
<p class="MsoNormal" style="text-indent:.5in">code_output<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
</body>
</html>