<html dir="ltr">
<head>
<meta http-equiv="Content-Type" content="text/html; charset=iso-8859-1">
<style id="owaParaStyle" type="text/css">P {margin-top:0;margin-bottom:0;}</style>
</head>
<body ocsi="0" fpstyle="1">
<div style="direction: ltr;font-family: Tahoma;color: #000000;font-size: 10pt;">Hello,<br>
<br>
I am trying to figure out a faster way of doing in petsc4py and thought I would ask to see how this should be done in petsc. If I can figure that out then I think I can figure out how to translate it but I am having issues with the basic algorithm I should
be using.<br>
<br>
I am solving a linear system of ODEs of the form <br>
<br>
B du/dt = A u<br>
<div><br>
<div style="font-family:Tahoma; font-size:13px">
<div style="font-family:Tahoma; font-size:13px">using Crank-Nicholson. So the update equation can be written as<br>
<br>
(B - dt/2*A)*u^(n+1) = (B + dt/2*A)*u^n<br>
<br>
or, if you don't mind using the backslash operator,<br>
<br>
u^(n+1) = (B - dt/2*A)\(B + dt/2*A)*u^n<br>
<br>
The integration is faster of course after I build the matrix <br>
<br>
C = (B - dt/2*A)\(B + dt/2*A)<br>
<br>
however building this is very slow and I am not sure how to build it efficiently. At the moment the majority of my code is just doing that.<br>
<br>
In petsc4py I do the following: <br>
<br>
(you can find the full code at https://github.com/francispoulin/qg_1L_ts/blob/master/qg_1L_ts.py)<br>
<br>
<div class="line js-file-line" id="LC122"><span class="n"> AeT</span> <span class="o">
=</span> <span class="n">PETSc</span><span class="o">.</span><span class="n">Mat</span><span class="p">()</span><span class="o">.</span><span class="n">createAIJ</span><span class="p">([</span><span class="n">Ny</span><span class="o">-</span><span class="mi">1</span><span class="p">,</span><span class="n">Ny</span><span class="o">-</span><span class="mi">1</span><span class="p">])</span></div>
<div class="line js-file-line" id="LC123"> <span class="n">AeT</span><span class="o">.</span><span class="n">setUp</span><span class="p">();</span>
<span class="n">AeT</span><span class="o">.</span><span class="n">assemble</span><span class="p">()</span></div>
<div class="line js-file-line" id="LC124"> <span class="n">ksp</span> <span class="o">
=</span> <span class="n">PETSc</span><span class="o">.</span><span class="n">KSP</span><span class="p">()</span><span class="o">.</span><span class="n">create</span><span class="p">(</span><span class="n">PETSc</span><span class="o">.</span><span class="n">COMM_WORLD</span><span class="p">)</span></div>
<div class="line js-file-line" id="LC125"> <span class="n">ksp</span><span class="o">.</span><span class="n">setOperators</span><span class="p">(</span><span class="n">B</span><span class="o">-</span><span class="n">dto2</span><span class="o">*</span><span class="n">A</span><span class="p">)</span></div>
<div class="line js-file-line" id="LC126"> <span class="n">ksp</span><span class="o">.</span><span class="n">setTolerances</span><span class="p">(</span><span class="mf">1e-16</span><span class="p">)</span></div>
<div class="line js-file-line" id="LC127"> <span class="n">pc</span> <span class="o">
=</span> <span class="n">ksp</span><span class="o">.</span><span class="n">getPC</span><span class="p">()</span></div>
<div class="line js-file-line" id="LC128"> <span class="n">pc</span><span class="o">.</span><span class="n">setType</span><span class="p">(</span><span class="s">'none'</span><span class="p">)</span></div>
<div class="line js-file-line" id="LC129"> <span class="n">ksp</span><span class="o">.</span><span class="n">setFromOptions</span><span class="p">()</span></div>
<div class="line js-file-line" id="LC130"><br>
</div>
<div class="line js-file-line" id="LC131"> <span class="n">btemp</span> <span class="o">
=</span> <span class="n">B</span><span class="o">+</span><span class="n">dto2</span><span class="o">*</span><span class="n">A</span></div>
<div class="line js-file-line" id="LC132"> <span class="n">btemp</span><span class="o">.</span><span class="n">assemble</span><span class="p">()</span></div>
<div class="line js-file-line" id="LC133"><br>
</div>
<div class="line js-file-line" id="LC134"> <span class="n">bcol</span> <span class="o">
=</span> <span class="n">PETSc</span><span class="o">.</span><span class="n">Vec</span><span class="p">()</span><span class="o">.</span><span class="n">createMPI</span><span class="p">(</span><span class="n">Ny</span><span class="o">-</span><span class="mi">1</span><span class="p">)</span></div>
<div class="line js-file-line" id="LC135"> <span class="n">xcol</span> <span class="o">
=</span> <span class="n">PETSc</span><span class="o">.</span><span class="n">Vec</span><span class="p">()</span><span class="o">.</span><span class="n">createMPI</span><span class="p">(</span><span class="n">Ny</span><span class="o">-</span><span class="mi">1</span><span class="p">)</span></div>
<div class="line js-file-line" id="LC136"> <span class="n">bcol</span><span class="o">.</span><span class="n">setUp</span><span class="p">();</span>
<span class="n">xcol</span><span class="o">.</span><span class="n">setUp</span><span class="p">()</span></div>
<div class="line js-file-line" id="LC137"><br>
</div>
<div class="line js-file-line" id="LC138"> <span class="n">sc</span><span class="p">,</span><span class="n">ec</span>
<span class="o">=</span> <span class="n">bcol</span><span class="o">.</span><span class="n">getOwnershipRange</span><span class="p">()</span></div>
<div class="line js-file-line" id="LC139"> <span class="k">for</span> <span class="n">
i</span> <span class="ow">in</span> <span class="nb">range</span><span class="p">(</span><span class="mi">0</span><span class="p">,</span><span class="n">Ny</span><span class="o">-</span><span class="mi">1</span><span class="p">):</span></div>
<div class="line js-file-line" id="LC140"> <span class="n">bcol</span><span class="p">[</span><span class="n">sc</span><span class="p">:</span><span class="n">ec</span><span class="p">]</span>
<span class="o">=</span> <span class="n">btemp</span><span class="p">[</span><span class="n">sc</span><span class="p">:</span><span class="n">ec</span><span class="p">,</span><span class="n">i</span><span class="p">]</span></div>
<div class="line js-file-line" id="LC141"> <span class="n">bcol</span><span class="o">.</span><span class="n">assemble</span><span class="p">();</span>
<span class="n">xcol</span><span class="o">.</span><span class="n">assemble</span><span class="p">()</span></div>
<div class="line js-file-line" id="LC142"> <span class="n">ksp</span><span class="o">.</span><span class="n">solve</span><span class="p">(</span><span class="n">bcol</span><span class="p">,</span><span class="n">xcol</span><span class="p">)</span></div>
<div class="line js-file-line" id="LC143"> <span class="n">AeT</span><span class="p">[</span><span class="n">sc</span><span class="p">:</span><span class="n">ec</span><span class="p">,</span><span class="n">i</span><span class="p">]</span>
<span class="o">=</span> <span class="n">xcol</span><span class="p">[</span><span class="n">sc</span><span class="p">:</span><span class="n">ec</span><span class="p">]</span></div>
<div class="line js-file-line" id="LC144"><br>
</div>
<div class="line js-file-line" id="LC145"> <span class="n">AeT</span><span class="o">.</span><span class="n">assemble</span><span class="p">()</span></div>
<br>
We are building the matrix one row at a time and, as accurate as that maybe, it is very slow. Is there perhaps a fast way of doing this kind of operation in petsc?<br>
<br>
Sorry for the bother but any advice would be greatly appreciated.<br>
<br>
Cheers, Francis<br>
<br>
------------------<br>
Francis Poulin <br>
Associate Professor<br>
Associate Chair, Undergraduate Studies <br>
Department of Applied Mathematics<br>
University of Waterloo<br>
<br>
email: fpoulin@uwaterloo.ca<br>
Web: https://uwaterloo.ca/poulin-research-group/<br>
Telephone: +1 519 888 4567 x32637<br>
<br>
</div>
</div>
</div>
</div>
</body>
</html>