Add Matrix to Assembled Matrix

2

If I need to edit values of an assembled matrix, I usually go through the petsc4py API. See the following example:

from dolfin import *

mesh = UnitIntervalMesh(2)
V = FunctionSpace(mesh,"Lagrange",1)
A = assemble(TrialFunction(V)*TestFunction(V)*dx)

# If by "+ B" you mean adding numbers directly to specific values of A, I
# usually do this through petsc4py:
import petsc4py, sys
petsc4py.init(sys.argv)
from petsc4py import PETSc

# To add to existing values:
ADD_MODE = PETSc.InsertMode.ADD

# To replace values instead:
#ADD_MODE = PETSc.InsertMode.INSERT

print(A.array())  # Original A

# (Note: Printed output is somewhat unclear in parallel)

# What value to add/insert at matrix entry i,j (i.e., the definition of B):
i = 0
j = 1
value = 1.0

Am = as_backend_type(A).mat()

# If the value you want to modify is not allocated as a nonzero, you need to
# set this option (with some cost to performance).  Ideally, you would
# set up a matrix with the necessary space allocated, assemble into that one,
# and then edit values without this.
Am.setOption(PETSc.Mat.Option.NEW_NONZERO_ALLOCATION_ERR, False)

# In parallel, you want to make sure a given value is added only once, ideally
# on the process that owns that row of the matrix.  (Can skip the check for
# ownership range in serial.)
Istart, Iend = Am.getOwnershipRange()
if(i<Iend and i>=Istart):
    Am.setValue(i,j,value,addv=ADD_MODE)
Am.assemble()

print(A.array())  # A with 1.0 added to the i,j entry

A complete listing of petsc4py functions and arguments can be found here:

https://www.mcs.anl.gov/petsc/petsc4py-current/docs/apiref/index.html

(You can usually figure out what all the functions/arguments mean by looking up the regular PETSc documentation, since the function names are similar.)

1 Like