<div dir="ltr"><div dir="ltr">Hello,<div><br></div><div>I am trying to use petsc4py to conduct a quantum mechanics simulation. I've been able to construct all of the relevant matrices, however I am reaching a gigantic bottleneck. </div><div><br></div><div>For the simplest problem I am running I have a few matrices each about 5000x5000. In order to begin time propagation I need to add these matrices together. However, on 6 cores of my local machine it is taking approximately 1-2 seconds per addition. Since I need to do this for each time step in my simulation it is prohibitively slow since there could be upwards of 10K time steps. </div><div><br></div><div>Below is the relevant code:</div><div><br></div><div><div style="color:rgb(166,172,205);background-color:rgb(41,45,62);font-family:Consolas,"Courier New",monospace;font-size:14px;line-height:19px;white-space:pre"><div>structure <span style="color:rgb(199,146,234)">=</span> structure<span style="color:rgb(199,146,234)">=</span>PETSc<span style="color:rgb(137,221,255)">.</span>Mat<span style="color:rgb(137,221,255)">.</span>Structure<span style="color:rgb(137,221,255)">.</span>DIFFERENT_NONZERO_PATTERN</div><div>    <span style="color:rgb(137,221,255);font-style:italic">if</span> test2<span style="color:rgb(137,221,255)">:</span></div><div>        <span style="color:rgb(199,146,234)">def</span> <span style="color:rgb(130,170,255)">makeLeft</span><span style="color:rgb(137,221,255)">(</span><span style="color:rgb(255,83,112)">S</span><span style="color:rgb(137,221,255)">,</span><span style="color:rgb(255,83,112)">MIX</span><span style="color:rgb(137,221,255)">,</span><span style="color:rgb(255,83,112)">ANG</span><span style="color:rgb(137,221,255)">,</span><span style="color:rgb(255,83,112)">ATOM</span><span style="color:rgb(137,221,255)">,</span><span style="color:rgb(255,83,112)">i</span><span style="color:rgb(137,221,255)">):</span></div><div>            </div><br><div>            S<span style="color:rgb(137,221,255)">.</span>axpy<span style="color:rgb(137,221,255)">(</span><span style="color:rgb(199,146,234)">-</span>Field<span style="color:rgb(137,221,255)">.</span>pulse<span style="color:rgb(137,221,255)">[</span>i<span style="color:rgb(137,221,255)">],</span>MIX<span style="color:rgb(137,221,255)">,</span>structure<span style="color:rgb(137,221,255)">)</span></div><div>            S<span style="color:rgb(137,221,255)">.</span>axpy<span style="color:rgb(137,221,255)">(</span><span style="color:rgb(199,146,234)">-</span>Field<span style="color:rgb(137,221,255)">.</span>pulse<span style="color:rgb(137,221,255)">[</span>i<span style="color:rgb(137,221,255)">],</span>ANG<span style="color:rgb(137,221,255)">,</span>structure<span style="color:rgb(137,221,255)">)</span></div><div>            S<span style="color:rgb(137,221,255)">.</span>axpy<span style="color:rgb(137,221,255)">(</span><span style="color:rgb(199,146,234)">-</span><span style="color:rgb(247,140,108)">1</span><span style="color:rgb(137,221,255)">,</span>ATOM<span style="color:rgb(137,221,255)">,</span>structure<span style="color:rgb(137,221,255)">)</span></div><div>            <span style="color:rgb(137,221,255);font-style:italic">return</span> S</div><div>        <span style="color:rgb(199,146,234)">def</span> <span style="color:rgb(130,170,255)">makeRight</span><span style="color:rgb(137,221,255)">(</span><span style="color:rgb(255,83,112)">S</span><span style="color:rgb(137,221,255)">,</span><span style="color:rgb(255,83,112)">MIX</span><span style="color:rgb(137,221,255)">,</span><span style="color:rgb(255,83,112)">ANG</span><span style="color:rgb(137,221,255)">,</span><span style="color:rgb(255,83,112)">ATOM</span><span style="color:rgb(137,221,255)">,</span><span style="color:rgb(255,83,112)">i</span><span style="color:rgb(137,221,255)">):</span></div><div>            </div><div>            </div><div>            </div><div>            S<span style="color:rgb(137,221,255)">.</span>axpy<span style="color:rgb(137,221,255)">(</span>Field<span style="color:rgb(137,221,255)">.</span>pulse<span style="color:rgb(137,221,255)">[</span>i<span style="color:rgb(137,221,255)">],</span>MIX<span style="color:rgb(137,221,255)">,</span>structure<span style="color:rgb(137,221,255)">)</span></div><div>            S<span style="color:rgb(137,221,255)">.</span>axpy<span style="color:rgb(137,221,255)">(</span>Field<span style="color:rgb(137,221,255)">.</span>pulse<span style="color:rgb(137,221,255)">[</span>i<span style="color:rgb(137,221,255)">],</span>ANG<span style="color:rgb(137,221,255)">,</span>structure<span style="color:rgb(137,221,255)">)</span></div><div>            S<span style="color:rgb(137,221,255)">.</span>axpy<span style="color:rgb(137,221,255)">(</span><span style="color:rgb(247,140,108)">1</span><span style="color:rgb(137,221,255)">,</span>ATOM<span style="color:rgb(137,221,255)">,</span>structure<span style="color:rgb(137,221,255)">)</span></div><div>            </div><div>            <span style="color:rgb(137,221,255);font-style:italic">return</span> S</div><br><div>        H_mix <span style="color:rgb(199,146,234)">=</span> Int<span style="color:rgb(137,221,255)">.</span>H_mix</div><div>        H_mix<span style="color:rgb(137,221,255)">.</span>scale<span style="color:rgb(137,221,255)">(</span><span style="color:rgb(247,140,108)">1</span><span style="color:rgb(199,146,234)">j</span> <span style="color:rgb(199,146,234)">*</span> dt <span style="color:rgb(199,146,234)">/</span><span style="color:rgb(247,140,108)">2</span><span style="color:rgb(137,221,255)">)</span></div><div>  </div><div>        H_ang <span style="color:rgb(199,146,234)">=</span> Int<span style="color:rgb(137,221,255)">.</span>H_ang</div><div>        H_ang<span style="color:rgb(137,221,255)">.</span>scale<span style="color:rgb(137,221,255)">(</span><span style="color:rgb(247,140,108)">1</span><span style="color:rgb(199,146,234)">j</span> <span style="color:rgb(199,146,234)">*</span> dt <span style="color:rgb(199,146,234)">/</span><span style="color:rgb(247,140,108)">2</span><span style="color:rgb(137,221,255)">)</span></div><br><div>        H_atom <span style="color:rgb(199,146,234)">=</span> Int<span style="color:rgb(137,221,255)">.</span>H_atom</div><div>        H_atom<span style="color:rgb(137,221,255)">.</span>scale<span style="color:rgb(137,221,255)">(</span><span style="color:rgb(247,140,108)">1</span><span style="color:rgb(199,146,234)">j</span> <span style="color:rgb(199,146,234)">*</span> dt <span style="color:rgb(199,146,234)">/</span><span style="color:rgb(247,140,108)">2</span><span style="color:rgb(137,221,255)">)</span></div><br><div>        S <span style="color:rgb(199,146,234)">=</span> Int<span style="color:rgb(137,221,255)">.</span>S_total</div><br><div>        psi_initial <span style="color:rgb(199,146,234)">=</span> psi<span style="color:rgb(137,221,255)">.</span>psi_initial<span style="color:rgb(137,221,255)">.</span>copy<span style="color:rgb(137,221,255)">()</span></div><div>        ksp <span style="color:rgb(199,146,234)">=</span> PETSc<span style="color:rgb(137,221,255)">.</span>KSP<span style="color:rgb(137,221,255)">().</span>create<span style="color:rgb(137,221,255)">(</span>PETSc<span style="color:rgb(137,221,255)">.</span>COMM_WORLD<span style="color:rgb(137,221,255)">)</span></div><br><div>   </div><div>        <span style="color:rgb(137,221,255);font-style:italic">for</span> i<span style="color:rgb(137,221,255)">,</span>t <span style="color:rgb(137,221,255);font-style:italic">in</span> <span style="color:rgb(130,170,255)">enumerate</span><span style="color:rgb(137,221,255)">(</span>box<span style="color:rgb(137,221,255)">.</span>t<span style="color:rgb(137,221,255)">):</span></div><div>            <span style="color:rgb(130,170,255)">print</span><span style="color:rgb(137,221,255)">(</span>i<span style="color:rgb(137,221,255)">,</span>L<span style="color:rgb(137,221,255)">)</span></div><br><div>            </div><br><div>            O_L <span style="color:rgb(199,146,234)">=</span> makeLeft<span style="color:rgb(137,221,255)">(</span>S<span style="color:rgb(137,221,255)">,</span>H_mix<span style="color:rgb(137,221,255)">,</span>H_ang<span style="color:rgb(137,221,255)">,</span>H_atom<span style="color:rgb(137,221,255)">,</span>i<span style="color:rgb(137,221,255)">)</span></div><div>            O_R <span style="color:rgb(199,146,234)">=</span> makeRight<span style="color:rgb(137,221,255)">(</span>S<span style="color:rgb(137,221,255)">,</span>H_mix<span style="color:rgb(137,221,255)">,</span>H_ang<span style="color:rgb(137,221,255)">,</span>H_atom<span style="color:rgb(137,221,255)">,</span>i<span style="color:rgb(137,221,255)">)</span></div><br><div>            </div><br><div>            <span style="color:rgb(137,221,255);font-style:italic">if</span> i <span style="color:rgb(199,146,234)">==</span> <span style="color:rgb(247,140,108)">0</span><span style="color:rgb(137,221,255)">:</span></div><div>                known <span style="color:rgb(199,146,234)">=</span> O_R<span style="color:rgb(137,221,255)">.</span>getVecRight<span style="color:rgb(137,221,255)">()</span></div><div>                sol <span style="color:rgb(199,146,234)">=</span> O_L<span style="color:rgb(137,221,255)">.</span>getVecRight<span style="color:rgb(137,221,255)">()</span></div><div>        </div><div>            O_R<span style="color:rgb(137,221,255)">.</span>mult<span style="color:rgb(137,221,255)">(</span>psi_initial<span style="color:rgb(137,221,255)">,</span>known<span style="color:rgb(137,221,255)">)</span></div><br><div>            ksp<span style="color:rgb(137,221,255)">.</span>setOperators<span style="color:rgb(137,221,255)">(</span>O_L<span style="color:rgb(137,221,255)">)</span></div><div>        </div><div>            </div><div>            ksp<span style="color:rgb(137,221,255)">.</span>solve<span style="color:rgb(137,221,255)">(</span>known<span style="color:rgb(137,221,255)">,</span>sol<span style="color:rgb(137,221,255)">)</span></div><br><div>            </div><div>            </div><div>    </div><div>            psi_initial<span style="color:rgb(137,221,255)">.</span>copy<span style="color:rgb(137,221,255)">(</span>sol<span style="color:rgb(137,221,255)">)</span></div><div><br></div></div></div><div><br></div><div>I need to clean it up a bit, but the main point is that I need to add the matrices many times for a single time step. I can't preallocate memory very well since some of the matrices aren't the most sparse either. It seems if I cannot speed up the addition it will be difficult to continue so I was wondering if you had any insights.</div><div><br></div><div>Thank you for your time</div></div></div>