Hi,

As said in the title, I am trying to compute `X=X+a*Y`

where `X`

and `Y`

are `Mat`

PETSc matrices and `a`

is a scalar. Is there a way to do this ? Here a stupid example :

```
load "PETSc"
mesh Th = square(10,10);
fespace Vh(Th,P1);
varf vX(u,v) = int2d(Th)(u*v) + on(1,u=0);
matrix X = vX(Vh,Vh);
matrix Y = vX(Vh,Vh);
Mat Xp(X);
Mat Yp(Y);
Mat Zp;
real a=2.0;
Yp = a*Xp + Yp;
```

I see something called `MatAXPY`

in PETSc doc but I don’t know if it is available in FreeFem ? (I see it in `PETSc-code.hpp`

)

I can do it with `matrix`

format and then import the results in `Mat`

but I loose the advantage of distributed matrices in large problem ?

Thanks for your help