code: def _initialize_b_x(self):
          r
          # Get vector(s) compatible with the matrix,
          # i.e., with the same parallel layout.
          self.petsc_x, self.petsc_b = self.petsc_A.getVecs()
  
          #  Set the solution vector to zeros.
          self.petsc_x.set(0)
  
          # Define the petsc rhs vector from the numpy one.
          # If the rhs is defined by blocks, use this:
          PETSc.Vec.setValuesBlocked(self.petsc_b, [sp.arange(self.m)], self.b)
docs: r"""
          Initialize the solution vector (self.petsc_x), which is a dense
          matrix (1D vector) and defines the rhs vector (self.petsc_b) from
          the existing data.