<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=utf-8">
<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:DengXian;
panose-1:2 1 6 0 3 1 1 1 1 1;}
@font-face
{font-family:Calibri;
panose-1:2 15 5 2 2 2 4 3 2 4;}
@font-face
{font-family:"\@DengXian";
panose-1:2 1 6 0 3 1 1 1 1 1;}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
{margin:0in;
margin-bottom:.0001pt;
font-size:11.0pt;
font-family:"Calibri",sans-serif;}
a:link, span.MsoHyperlink
{mso-style-priority:99;
color:blue;
text-decoration:underline;}
span.apple-converted-space
{mso-style-name:apple-converted-space;}
span.EmailStyle19
{mso-style-type:personal-reply;
font-family:"Calibri",sans-serif;
color:windowtext;}
.MsoChpDefault
{mso-style-type:export-only;
font-size:10.0pt;}
@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">
<div class="WordSection1">
<p class="MsoNormal">Thanks a lot for the example!<o:p></o:p></p>
<p class="MsoNormal"><o:p> </o:p></p>
<div style="border:none;border-top:solid #B5C4DF 1.0pt;padding:3.0pt 0in 0in 0in">
<p class="MsoNormal"><b><span style="font-size:12.0pt;color:black">From: </span></b><span style="font-size:12.0pt;color:black">"Zhang, Hong" <hongzhang@anl.gov><br>
<b>Date: </b>Tuesday, February 18, 2020 at 11:32 PM<br>
<b>To: </b>Yuyun Yang <yyang85@stanford.edu><br>
<b>Cc: </b>Matthew Knepley <knepley@gmail.com>, "petsc-users@mcs.anl.gov" <petsc-users@mcs.anl.gov><br>
<b>Subject: </b>Re: [petsc-users] Matrix-free method in PETSc<o:p></o:p></span></p>
</div>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<div>
<p class="MsoNormal">Here is an TS example using DMDA and matrix-free Jacobians. Though the matrix-free part is faked, it demonstrates the workflow.<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<p class="MsoNormal"><a href="https://gitlab.com/petsc/petsc/-/blob/hongzh/ts-matshell-example/src/ts/examples/tutorials/advection-diffusion-reaction/ex5_mf.c">https://gitlab.com/petsc/petsc/-/blob/hongzh/ts-matshell-example/src/ts/examples/tutorials/advection-diffusion-reaction/ex5_mf.c</a>
<o:p></o:p></p>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
<div>
<p class="MsoNormal">Hong (Mr.)<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal"><o:p> </o:p></p>
<div>
<p class="MsoNormal"><br>
<br>
<o:p></o:p></p>
<blockquote style="margin-top:5.0pt;margin-bottom:5.0pt">
<div>
<p class="MsoNormal">On Feb 18, 2020, at 8:26 AM, Yuyun Yang <<a href="mailto:yyang85@stanford.edu">yyang85@stanford.edu</a>> wrote:<o:p></o:p></p>
</div>
<p class="MsoNormal"><o:p> </o:p></p>
<div>
<div>
<p class="MsoNormal">Thanks. Also, when using KSP, would the syntax be KSPSetOperators(ksp,A,A)? Since you mentioned preconditioners are not generally used for matrix-free operators, I wasn’t sure whether I should still put “A” in the Pmat field.<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
<div>
<p class="MsoNormal">Is it still possible to use TS in conjunction with the matrix-free operator? I’d like to create a simple test case that solves the 1d heat equation implicitly with variable coefficients, but didn’t know how the time stepping can be set
up.<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
<div>
<p class="MsoNormal">Thanks,<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal">Yuyun<o:p></o:p></p>
</div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
<div style="border:none;border-top:solid #B5C4DF 1.0pt;padding:3.0pt 0in 0in 0in">
<div>
<p class="MsoNormal"><b><span style="font-size:12.0pt">From:<span class="apple-converted-space"> </span></span></b><span style="font-size:12.0pt">Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>><br>
<b>Date:<span class="apple-converted-space"> </span></b>Tuesday, February 18, 2020 at 9:23 PM<br>
<b>To:<span class="apple-converted-space"> </span></b>Yuyun Yang <<a href="mailto:yyang85@stanford.edu">yyang85@stanford.edu</a>><br>
<b>Cc:<span class="apple-converted-space"> </span></b>"Smith, Barry F." <<a href="mailto:bsmith@mcs.anl.gov">bsmith@mcs.anl.gov</a>>, "<a href="mailto:petsc-users@mcs.anl.gov">petsc-users@mcs.anl.gov</a>" <<a href="mailto:petsc-users@mcs.anl.gov">petsc-users@mcs.anl.gov</a>><br>
<b>Subject:<span class="apple-converted-space"> </span></b>Re: [petsc-users] Matrix-free method in PETSc</span><o:p></o:p></p>
</div>
</div>
<div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
</div>
<div>
<div>
<div>
<p class="MsoNormal">On Tue, Feb 18, 2020 at 8:20 AM Yuyun Yang <<a href="mailto:yyang85@stanford.edu">yyang85@stanford.edu</a>> wrote:<o:p></o:p></p>
</div>
</div>
<div>
<blockquote style="border:none;border-left:solid #CCCCCC 1.0pt;padding:0in 0in 0in 6.0pt;margin-left:4.8pt;margin-top:5.0pt;margin-right:0in;margin-bottom:5.0pt">
<div>
<p class="MsoNormal">Thanks for the clarification.<br>
<br>
Got one more question: if I have variable coefficients, my stencil will be updated at every time step, so will the coefficients in myMatMult. In that case, is it necessary to destroy the shell matrix and create it all over again, or can I use it as it is, only
calling the stencil update function, assuming the result will be passed into the matrix operation automatically?<o:p></o:p></p>
</div>
</blockquote>
<div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
</div>
<div>
<div>
<p class="MsoNormal">You update the information in the context associated with the shell matrix. No need to destroy it.<o:p></o:p></p>
</div>
</div>
<div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
</div>
<div>
<div>
<p class="MsoNormal"> Thanks,<o:p></o:p></p>
</div>
</div>
<div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
</div>
<div>
<div>
<p class="MsoNormal"> Matt<o:p></o:p></p>
</div>
</div>
<div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
</div>
<blockquote style="border:none;border-left:solid #CCCCCC 1.0pt;padding:0in 0in 0in 6.0pt;margin-left:4.8pt;margin-top:5.0pt;margin-right:0in;margin-bottom:5.0pt">
<p class="MsoNormal" style="margin-bottom:12.0pt">Thanks,<br>
Yuyun<br>
<br>
On 2/18/20, 7:34 AM, "Smith, Barry F." <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>> wrote:<br>
<br>
<br>
<br>
> On Feb 17, 2020, at 7:56 AM, Yuyun Yang <<a href="mailto:yyang85@stanford.edu" target="_blank">yyang85@stanford.edu</a>> wrote:<br>
><span class="apple-converted-space"> </span><br>
> Hello,<br>
><span class="apple-converted-space"> </span><br>
> I actually have a question about the usage of DMDA since I'm quite new to this. I wonder if the DMDA suite of functions can be directly called on vectors created from VecCreate?<br>
<br>
Yes, but you have to make sure the ones you create have the same sizes and parallel layouts. Generally best to get them from the DMDA or VecDuplicate() than the hassle of figuring out sizes.<br>
<br>
> Or the vectors have to be formed by DMDACreateGlobalVector? I'm also not sure about what the dof and stencil width arguments do.<br>
><span class="apple-converted-space"> </span><br>
> I'm still unsure about the usage of MatCreateShell and MatShellSetOperation, since it seems that MyMatMult should still have 3 inputs just like MatMult (the matrix and two vectors). Since I'm not forming the matrix, does that mean the matrix input is
meaningless but still needs to exist for the sake of this format?<br>
<br>
Well the matrix input is your shell matrix so it likely has information you need to do your multiply routine. MatShellGetContext() (No you do not want to put your information about the matrix stencil inside global variables!)<br>
<br>
<br>
><span class="apple-converted-space"> </span><br>
> After I create such a shell matrix, can I use it like a regular matrix in KSP and utilize preconditioners?<br>
><span class="apple-converted-space"> </span><br>
> Thanks!<br>
> Yuyun<br>
> From: petsc-users <<a href="mailto:petsc-users-bounces@mcs.anl.gov" target="_blank">petsc-users-bounces@mcs.anl.gov</a>> on behalf of Yuyun Yang <<a href="mailto:yyang85@stanford.edu" target="_blank">yyang85@stanford.edu</a>><br>
> Sent: Sunday, February 16, 2020 3:12 AM<br>
> To: Smith, Barry F. <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>><br>
> Cc:<span class="apple-converted-space"> </span><a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a><span class="apple-converted-space"> </span><<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a>><br>
> Subject: Re: [petsc-users] Matrix-free method in PETSc<br>
> <span class="apple-converted-space"> </span><br>
> Thank you, that is very helpful information indeed! I will try it and send you my code when it works.<br>
><span class="apple-converted-space"> </span><br>
> Best regards,<br>
> Yuyun<br>
> From: Smith, Barry F. <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>><br>
> Sent: Saturday, February 15, 2020 10:02 PM<br>
> To: Yuyun Yang <<a href="mailto:yyang85@stanford.edu" target="_blank">yyang85@stanford.edu</a>><br>
> Cc:<span class="apple-converted-space"> </span><a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a><span class="apple-converted-space"> </span><<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a>><br>
> Subject: Re: [petsc-users] Matrix-free method in PETSc<br>
> <span class="apple-converted-space"> </span><br>
> Yuyun,<br>
><span class="apple-converted-space"> </span><br>
> If you are speaking about using a finite difference stencil on a structured grid where you provide the Jacobian vector products yourself by looping over the grid doing the stencil operation we unfortunately do not have exactly that kind of example.<span class="apple-converted-space"> </span><br>
><span class="apple-converted-space"> </span><br>
> But it is actually not difficult. I suggest starting with src/ts/examples/tests/ex22.c It computes the sparse matrix explicitly with FormIJacobian()<span class="apple-converted-space"> </span><br>
><span class="apple-converted-space"> </span><br>
> What you need to do is instead in main() use MatCreateShell() and MatShellSetOperation(,MATOP_MULT,(void (*)(void))MyMatMult) then provide the routine MyMatMult() to do your stencil based matrix free product; note that you can create this new routine
by taking the structure of IFunction() and reorganizing it to do the Jacobian product instead. You will need to get the information about the shell matrix size on each process by calling DMDAGetCorners().<span class="apple-converted-space"> </span><br>
><span class="apple-converted-space"> </span><br>
> You will then remove the explicit computation of the Jacobian, and also remove the Event stuff since you don't need it.<br>
><span class="apple-converted-space"> </span><br>
> Extending to 2 and 3d is straight forward.<span class="apple-converted-space"> </span><br>
><span class="apple-converted-space"> </span><br>
> Any questions let us know.<br>
><span class="apple-converted-space"> </span><br>
> Barry<br>
><span class="apple-converted-space"> </span><br>
> If you like this would make a great merge request with your code to improve our examples.<br>
><span class="apple-converted-space"> </span><br>
><span class="apple-converted-space"> </span><br>
> > On Feb 15, 2020, at 9:42 PM, Yuyun Yang <<a href="mailto:yyang85@stanford.edu" target="_blank">yyang85@stanford.edu</a>> wrote:<br>
> ><span class="apple-converted-space"> </span><br>
> > Hello team,<br>
> ><span class="apple-converted-space"> </span><br>
> > I wanted to apply the Krylov subspace method to a matrix-free implementation of a stencil, such that the iterative method acts on the operation without ever constructing the matrix explicitly (for example, when doing backward Euler).<br>
> ><span class="apple-converted-space"> </span><br>
> > I'm not sure whether there is already an example for that somewhere. If so, could you point me to a relevant example?<br>
> ><span class="apple-converted-space"> </span><br>
> > Thank you!<br>
> ><span class="apple-converted-space"> </span><br>
> > Best regards,<br>
> > Yuyun<br>
<br>
<br>
<br>
<o:p></o:p></p>
</blockquote>
</div>
<div>
<p class="MsoNormal"><br clear="all">
<o:p></o:p></p>
</div>
<div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
</div>
<div>
<p class="MsoNormal">--<span class="apple-converted-space"> </span><o:p></o:p></p>
</div>
<div>
<div>
<div>
<div>
<div>
<div>
<div>
<div>
<p class="MsoNormal">What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener<o:p></o:p></p>
</div>
</div>
<div>
<div>
<p class="MsoNormal"> <o:p></o:p></p>
</div>
</div>
<div>
<div>
<p class="MsoNormal"><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><o:p></o:p></p>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</div>
</blockquote>
</div>
<p class="MsoNormal"><o:p> </o:p></p>
</div>
</div>
</body>
</html>