Adding rows of Matrix

Linear Algebra
1

I want to add two rows of an assembled matrix, and reinsert the result in the original Matrix in place of the first of the added rows, and set the second row to the identity matrix, see the MWE below.
Unfortunately I get an error when trying to set the new values in the Matrix, see below, everything else seems to work.
Am I using the set() function from the Matrix class wrong, or is there another way to do this?

from dolfin import *
import numpy as np

mesh = UnitSquareMesh(5,5)

V = VectorFunctionSpace(mesh, 'CG', 2)

u = TrialFunction(V)
v = TestFunction(V)

A = assemble( inner(u, v) * dx )

row_block_0 = np.ndarray( ( 1, V.dim() ), dtype=np.float_)
row_block_1 = np.ndarray( ( 1, V.dim() ), dtype=np.float_)

A.get( row_block_0, [0], list(range(V.dim())) )
A.get( row_block_1, [1], list(range(V.dim())) )

A.set( row_block_0+row_block_1, [1], list(range(V.dim())) )

A.ident([1])

Traceback (most recent call last):
File “960_Test.py”, line 19, in
A.set( row_block_0+row_block_1, [1], list(range(V.dim())) )
RuntimeError:

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at

*** fenics-support@googlegroups.com

*** Remember to include the error message listed below and, if possible,
*** include a minimal running example to reproduce the error.

*** -------------------------------------------------------------------------
*** Error: Unable to successfully call PETSc function ‘MatSetValues’.
*** Reason: PETSc error code is: 63 (Argument out of range).
*** Where: This error was encountered inside /build/dolfin-ltzKbC/dolfin-2019.1.0/dolfin/la/PETScMatrix.cpp.
*** Process: 0

*** DOLFIN version: 2019.1.0
*** Git changeset: unknown
*** ---------------------------------------------

2

Hi,

Have you tried using the petsc4py API as mentioned here? Something like

from dolfin import *
import numpy as np
from petsc4py import PETSc
addMode = PETSc.InsertMode.ADD

mesh = UnitSquareMesh(5,5)

V = VectorFunctionSpace(mesh, 'CG', 2)

u = TrialFunction(V)
v = TestFunction(V)
A = assemble( inner(u, v) * dx )

row_block_0 = np.ndarray( ( 1, V.dim() ), dtype=np.float_)
row_block_1 = np.ndarray( ( 1, V.dim() ), dtype=np.float_)
identityRow = np.identity(V.dim())[1] #row to set to identity

A.get( row_block_0, [0], list(range(V.dim())) )
A.get( row_block_1, [1], list(range(V.dim())) )
print(A.array())  # original


Am = as_backend_type(A).mat()
Am.setOption(PETSc.Mat.Option.NEW_NONZERO_ALLOCATION_ERR, False) 

# setting the first row to the sum of first two rows
Am.setValues([0], list(range(V.dim())), row_block_0 + row_block_1)

# setting the row-1 to identity: assuming this is what A.ident([1]) does
Am.setValues([1], list(range(V.dim())), identityRow)

Am.assemble()

print(A.array())  # modified
3

Thank You, @bhaveshshrimali, this works and solves this part of my problem.

However I would still like to know how the set() function for matrices is supposed to be used, so if someone knows I would appreciate an explanation.