> On Feb 8, 2015, at 7:20 PM, Ronal Celaya <[email protected]> wrote: > > I know that. I want to have all the vector x replicated in all processes and > update it in each iteration
You will have to do communication when you "update it each iteration", no less communication then is done for the MatMult() so you will save no communication. Believe me, people have been doing parallel Krylov methods for 25 years, there is no savings on communication that can be had; it can only be shifted around to different points in the algorithm. Barry > , so I don't need to communicate the vector x each time MatMult() is called. > I'm not sure I'm making myself clear, sorry > > On Sun, Feb 8, 2015 at 7:52 PM, Barry Smith <[email protected]> wrote: > > > On Feb 8, 2015, at 6:14 PM, Ronal Celaya <[email protected]> wrote: > > > > Thank you Barry. > > Is there a way to reuse the vector x? I don't want to gather the vector in > > each iteration, I'd rather replicate the vector x in each process. > > I don't understand. With each new matrix vector product there are new > values in x (which are updated from other parts of the CG algorithm), these > new values need to be communicated in the MatMult() to where they are needed; > you can't just reuse the old values in x. > > Barry > > > > > Thanks in advance. > > > > On Sun, Feb 8, 2015 at 7:17 PM, Barry Smith <[email protected]> wrote: > > > > > On Feb 8, 2015, at 5:41 PM, Ronal Celaya <[email protected]> > > > wrote: > > > > > > Hello > > > If I have a MatMult operation inside a for loop (e. g. CG algorithm), and > > > the matrix A is MPIAIJ, vector x is gathered to local process in every > > > loop? > > > > Yes, internal to MatMult() it calls MatMult_MPIAIJ() which is in > > src/mat/impls/aij/mpi/mpiaij,c which has the following code: > > > > PetscErrorCode MatMult_MPIAIJ(Mat A,Vec xx,Vec yy) > > { > > Mat_MPIAIJ *a = (Mat_MPIAIJ*)A->data; > > PetscErrorCode ierr; > > PetscInt nt; > > > > PetscFunctionBegin; > > ierr = VecGetLocalSize(xx,&nt);CHKERRQ(ierr); > > if (nt != A->cmap->n) > > SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Incompatible partition of A > > (%D) and xx (%D)",A->cmap->n,nt); > > ierr = > > VecScatterBegin(a->Mvctx,xx,a->lvec,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); > > ierr = (*a->A->ops->mult)(a->A,xx,yy);CHKERRQ(ierr); > > ierr = > > VecScatterEnd(a->Mvctx,xx,a->lvec,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); > > ierr = (*a->B->ops->multadd)(a->B,a->lvec,yy,yy);CHKERRQ(ierr); > > PetscFunctionReturn(0); > > > > The needed values of x are communicated in the VecScatterBegin() to > > VecScatterEnd(). Note only exactly those values needed by each process are > > communicated in the scatter so not all values are communicated to all > > processes. Since the matrix is very sparse (normally) only a small > > percentage of the values need to be communicated. > > > > Barry > > > > > > > > I'm sorry for my English. > > > > > > Regards, > > > > > > -- > > > Ronal Celaya > > > > > > > > > > -- > > Ronal Celaya > > > > > -- > Ronal Celaya
