Matrix operation

[billy at dem.uminho.pt]

Thank you very much for your answers.

I was thinking why doesn't PETSc allow functions inside functions (like LASPack)?

h += (Diag(A) - A) * u <=> AddAsgn_VV (h, Mul_QV (Sub_QQ (Diag_Q (A), A), u));

Billy



Quoting Matthew Knepley <knepley at mcs.anl.gov>:

> Hong Zhang <hzhang at mcs.anl.gov> writes:
> 
>      I would recommend the later
> 
>   MatGetDiagonal(A, diag);
>   VecPointwiseMult(w, diag, u);
>   MatMult(A, u, v);
>   VecAXPY(w, -1.0, v);
> 
>              Matt
> 
> > Billy,
> >
> > You can do
> > MatGetDiagonal(A,diag);
> > MatMult(A,u,&vtmp);
> > VecGetArray(diag,diag_array);
> > VecGetArray(vtmp,vtmp_array);
> > VecGetArray(u,u_arrar);
> > VecGetArray(h,h_array);
> > VecGetLocalSize(h,&n);
> > for (i=0; i<n; i++){
> >   h_array[i] += diag_array[i]*u_array[i] - vtmp_array[i];
> > }
> > VecRestoreArray(h,h_array);
> > VecRestoreArray()
> > ...
> >
> > or you can use VecPointwiseMult() got get
> > Diag(A)*u.
> >
> > Hong
> >
> > On Mon, 20 Feb 2006 billy at dem.uminho.pt wrote:
> >
> >>
> >> Hello,
> >>
> >> I would like to do the following matrix-vector operation in PETSc:
> >>
> >> h += (Diag(A) - A) * u
> >>
> >> where A is a matrix, Diag(A) is a matrix and h, u are vectors. It is the
> >> multiplication of a matrix with diagonals set to zero with a vector.
> >>
> >>
> >> How can I do this in PETSc? because MatGetDiagonal gives a vector not a
> matrix.
> >>
> >>
> >> Billy.
> >>
> >>
> >>
> >
> >
> >
> 
> -- 
> "Failure has a thousand explanations. Success doesn't need one" -- Sir Alec
> Guiness
> 
>