You should be able to use the PETSc.Mat and PETSc.Vec functions in petsc4py to compute these matrix-matrix-matrix problems.
However, as you seem to be using np.linalg.svd to get your singular value decomposition, the easiest way might be to convert the petsc-matrix A to a numpy matrix: Transform dolfinx assemble matrix to numpy array - #2 by dokken
and compute the matrix-matrix-matrix product using numpy.
1 Like