
<html><head></head><body><div class="yahoo-style-wrap" style="font-family:Helvetica Neue, Helvetica, Arial, sans-serif;font-size:13px;"><div dir="ltr" data-setdir="false">Hi All,</div><div dir="ltr" data-setdir="false"><br></div><div dir="ltr" data-setdir="false">I need to fill non-zero values of a Petsc matrix via petsc4py for the domain defined by <span>A.getOwnershipRange() using three Numpy arrays: (1) array containing row indices of non-zero value, (2) array containing column indices of non-zero values and (3) array containing the non-zero matrix values. How can one perform this type of filling operation in <span>petsc4py</span>? The method A.setValues does not appear to allow this since it only works on an individual matrix element or a block of matrix elements.<br></span></div><div dir="ltr" data-setdir="false"><span><br></span></div><div dir="ltr" data-setdir="false"><span>I am using Numpy arrays since they can be computed in loops optimized using Numba on each processor. <span><span>I also cannot pass the Petsc matrix to a Numba compiled function since type information cannot be inferred.</span></span> I absolutely need to avoid looping in standard Python to define Petsc matrix elements due to performance issues. I also need to use a standard petscy4py method and avoid writing new C or Fortran wrappers to minimize language complexity.<br></span></div><div dir="ltr" data-setdir="false"><span><br></span></div><div dir="ltr" data-setdir="false"><span>Example Algorithm Building on <span>Lisandro Dalcin's 2D Poisson Example:</span></span></div><div dir="ltr" data-setdir="false"><span><span>----------------------------------------------------------------------------------------------<br></span></span></div><div dir="ltr" data-setdir="false"><span><span>comm = PETSc.COMM_WORLD</span><br></span></div><div dir="ltr" data-setdir="false"><span><span>rank = comm.getRank()</span><br></span></div><div dir="ltr" data-setdir="false"><span><span><br></span></span></div><div dir="ltr" data-setdir="false"><span><span>dx = 1.0/(xnodes + 1)</span> # xnodes is the number of nodes in the x and y-directions of the grid<br></span></div><div dir="ltr" data-setdir="false"><span><span><span>nnz_max</span></span> = 5 # max number of non-zero values per row<br></span></div><div dir="ltr" data-setdir="false"><span><span><br></span></span></div><div dir="ltr" data-setdir="false"><span><span>A = PETSc.Mat()</span><br></span></div><div dir="ltr" data-setdir="false"><div>A.create(comm=PETSc.COMM_WORLD)<br>A.setSizes((xnodes*ynodes, xnodes*ynodes))<br><div>A.setType(PETSc.Mat.Type.AIJ)</div><div dir="ltr" data-setdir="false"><span>A.setPreallocationNNZ(<span><span><span><span>nnz_max</span></span></span></span>)</span></div><div dir="ltr" data-setdir="false"><span><span><br></span></span></div><div dir="ltr" data-setdir="false"><span><span>rstart, rend = A.getOwnershipRange()</span></span></div><div dir="ltr" data-setdir="false"><span><div><div><br></div><div dir="ltr" data-setdir="false"># Here Anz, Arow and Acol are vectors with size equal to the number of non-zero values <br></div><div>Anz, Arow, Acol = build_nonzero_numpy_arrays_using_numba(rstart, rend, nnz_max, dx, xnodes, ynodes)</div><div><br></div><div>A.setValues(Arow, Acol, Anz) # <--- This does not work.<br></div></div></span><div><br></div><div dir="ltr" data-setdir="false"><div>A.assemblyBegin()<br><div>A.assemblyEnd()</div><div dir="ltr" data-setdir="false"><div>ksp = PETSc.KSP()<br>ksp.create(comm=A.getComm())</div><div dir="ltr" data-setdir="false"><div>ksp.setType(PETSc.KSP.Type.CG)<br>ksp.getPC().setType(PETSc.PC.Type.GAMG)</div><div dir="ltr" data-setdir="false"><div>ksp.setOperators(A)<br>ksp.setFromOptions()</div><div dir="ltr" data-setdir="false"><div>x, b = A.createVecs()<br><div>b.set(1.0)</div><div><br></div></div><div dir="ltr" data-setdir="false"><span>ksp.solve(b, x)</span></div></div></div></div><div><br></div><div dir="ltr" data-setdir="false">Regards,</div><div dir="ltr" data-setdir="false">Erik Kneller, Ph.D.<br></div></div></div><div><br></div></div><span></span></div></div></div></div></body></html>