[petsc-users] Galerkin projection using petsc4py

Thanasis Boutsikakis thanasis.boutsikakis at corintis.com
Thu Oct 5 06:09:11 CDT 2023


Hi everyone,

I am trying a Galerkin projection (see MFE below) and I cannot get the Phi.transposeMatMult(A, A1) work. The error is

    Phi.transposeMatMult(A, A1)
  File "petsc4py/PETSc/Mat.pyx", line 1514, in petsc4py.PETSc.Mat.transposeMatMult
petsc4py.PETSc.Error: error code 56
[0] MatTransposeMatMult() at /Users/boutsitron/firedrake/src/petsc/src/mat/interface/matrix.c:10135
[0] MatProduct_Private() at /Users/boutsitron/firedrake/src/petsc/src/mat/interface/matrix.c:9989
[0] No support for this operation for this object type
[0] Call MatProductCreate() first

Do you know if these exposed to petsc4py or maybe there is another way? I cannot get the MFE to work (neither in sequential nor in parallel)

"""Experimenting with PETSc mat-mat multiplication"""

import time

import numpy as np
from colorama import Fore
from firedrake import COMM_SELF, COMM_WORLD
from firedrake.petsc import PETSc
from mpi4py import MPI
from numpy.testing import assert_array_almost_equal

from utilities import (
    Print,
    create_petsc_matrix,
)

nproc = COMM_WORLD.size
rank = COMM_WORLD.rank

# --------------------------------------------
# EXP: Galerkin projection of an mpi PETSc matrix A with an mpi PETSc matrix Phi
#  A' = Phi.T * A * Phi
# [k x k] <- [k x m] x [m x m] x [m x k]
# --------------------------------------------

m, k = 11, 7
# Generate the random numpy matrices
np.random.seed(0)  # sets the seed to 0
A_np = np.random.randint(low=0, high=6, size=(m, m))
Phi_np = np.random.randint(low=0, high=6, size=(m, k))

# Create A as an mpi matrix distributed on each process
A = create_petsc_matrix(A_np)

# Create Phi as an mpi matrix distributed on each process
Phi = create_petsc_matrix(Phi_np)

A1 = create_petsc_matrix(np.zeros((k, m)))

# Now A1 contains the result of Phi^T * A
Phi.transposeMatMult(A, A1)

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20231005/feb6dc6e/attachment-0001.html>


More information about the petsc-users mailing list