How to perform simple operations in PETSc like matrix addition, multiplication, constant times a matrix? I’ve tried MatMatMult and MatScale both did not work.
MatMatMult
should work if you have a recent enough FreeFEM installation, cf. https://github.com/FreeFem/FreeFem-sources/blob/develop/plugin/mpi/PETSc-code.hpp#L4460. MatScale
is simply A *= alpha
. Addition is A += alpha * B
where both A
and B
are Mat
. You can also do, assuming your Mat
follow the same distributions:
matrix a = alpha * b + gamma * c; // sequential FreeFEM matrix
Mat B(A, a); // B and A follow the same distribution
This is usually more efficient since the sequential matrix operations do not require any communication.
Is there a direct way to obtain the number of rows and columns like in sequential? ex: B.n or B.m
.n
will give you the number of local rows of a Mat
. .m
does not exist, but since you need to Mat
do define a rectangular one, i.e., Mat C(A, B)
, instead of doing C.m
(which does not exist), you can do B.n
.
Ah I see, thank you for your reply. I would like to add something about my first comment regarding matrix operations, whenever I try to do for instance A+=alpha I get the following error:
any idea on how to proceed?
[0]PETSC ERROR: ------------------------------------------------------------------------
[0]PETSC ERROR: Caught signal number 11 SEGV: Segmentation Violation, probably memory access out of range
[0]PETSC ERROR: Try option -start_in_debugger or -on_error_attach_debugger
[0]PETSC ERROR: or see https://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind
[0]PETSC ERROR: or try http://valgrind.org on GNU/linux and Apple Mac OS X to find memory corruption errors
[0]PETSC ERROR: configure using --with-debugging=yes, recompile, link, and run
[0]PETSC ERROR: to get more information on the crash.
[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
[0]PETSC ERROR: Signal received
[0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.
[0]PETSC ERROR: Petsc Release Version 3.15.0, Mar 30, 2021
[0]PETSC ERROR: /softs/local/freefem/4.9/bin/FreeFem+±mpi on a named node001.cluster by mnoun Mon Sep 6 10:39:22 2021
[0]PETSC ERROR: Configure options --prefix=/softs/local/freefem/4.9/ff-petsc/r MAKEFLAGS= --with-debugging=0 COPTFLAGS="-O3 -mtune=generic" CXXOPTFLAGS="-O3 -mtune=generic" FOPTFLAGS="-O3 -mtune=generic" --with-cxx-dialect=C++11 --with-ssl=0 --with-x=0 --with-fortran-bindings=0 --with-cudac=0 --with-scalar-type=real --with-blaslapack-include= --with-blaslapack-lib="-llapack -lblas" --download-metis --download-ptscotch --download-hypre --download-parmetis --download-mmg --download-parmmg --download-superlu --download-suitesparse --download-tetgen --download-slepc --download-hpddm --with-cc=gcc --with-cxx=g++ --download-mpich --with-fc=gfortran --download-scalapack --download-mumps --download-slepc-configure-arguments=–download-arpack=https://github.com/prj-/arpack-ng/archive/b64dccb.tar.gz PETSC_ARCH=fr
[0]PETSC ERROR: #1 User provided function() at unknown file:0
[0]PETSC ERROR: Run with -malloc_debug to check if memory corruption is causing the crash.
That’s a very generic error. It’s nearly impossible to help without a MWE to reproduce the error.
This is a part of the code : testpetsc.edp (3.0 KB) I defined three matrices K, C and M and I am trying to combine them in the matrix a1.
Neither a1
nor a2
are defined, so it doesn’t make sense to do a1 += something
with an undefined a1
. Just put everything in a matrix linearComb
, and then do Mat a1(K, linearComb)
.
I defined a1 and a2 on line 81 as Mat a1, a2. Could you elaborate a little bit more your method please?
Edit: K, C and M are Mat type also
Sorry, I meant that none of the matrices are initialized, so you cannot do unintializedMat += something
, it’s undefined behavior (and results in a segfault).
Instead of
K = stiffness(Vh, Vh);
C = damping(Vh,Vh);
M = mass(Vh,Vh);
Do
matrix locK = stiffness(Vh, Vh);
K = locK;
matrix locC = damping(Vh,Vh);
C = locC;
matrix locM = damping(Vh,Vh);
M = locM;
matrix locLinearComb = alpha * locK + beta * locC + gamma * locM;
Mat linearComb(M, locLinearComb);