Hi,<br><br>I noticed that PETSc doesn't have a routine to compute<br>the outer product of two Vec objects.<br><br>Specifically, I want to compute the outer product for two<br>Vec objects that are large enough to require storage on<br>
multiple processors.<br><br>If I write my own routine in such a way that it performs<br>local computations of the entries of the outer product<br>matrix (called OP from hereon), say using <br>VecGetOwnershipRange to retrieve local entries of the<br>
two vectors (Vec a and Vec b) and then computing local<br>blocks of OP, then there will be parts of OP that will not<br>be computed.<br><br>Say we have 3 processors, and<br>Vec a is (2,4,8,16,32), Vec b is (1, 1/2, 1/4). Then the<br>
above approach will produce OP as<br><br>2 - -<br>4 - -<br>- 4 -<br>- 8 -<br>- - 8<br><br>instead of <br><br>2 1 0.5<br>4 2 1<br>8 4 2<br>16 8 4<br>32 16 8<br><br>So, I chose to implement my routine in this way instead >><br>
I would convert the "row" vector {vector b in a*transpose(b) }<br>into an array of PetscScalar values, using PetscMalloc.<br>Then, I run a for-loop traversing this array-version of vector b<br>where in each iteration of the loop, the processor-specific<br>
part of vector a is scaled with the element of b that is being<br>accessed in the loop. Then MatSetValues is called to set<br>the scaled part-of-vector-a in the right locations in OP.<br><br>However, this approach is also producing the same result as<br>
above. My guess is that the array-version of Vec b created<br>using PetscMalloc is not a globally visible array that is being accessed<br>by each processor -- instead, each processor is creating its own<br>version of array-vector-b and then again computing only parts of<br>
OP. Further, I think this probably has to do with the fact that<br>the routine PetscMalloc is not collective.<br><br>How do I get an array of PetscScalar values that is not<br>processor-specific but is visible to all processors? If this<br>
cannot be done, how do I go about computing the outer<br>product of two vectors?<br><br>I am appending my code (routine to compute outer product and<br>the "main" function) and sample output with this mail.<br><br>
Thanks very much,<br>Harsh<br><br>---- code and sample output ----<br><br><span style="font-family: courier new,monospace;">static char helpMsg[] = "\nComputes outer product of two vectors.\n";</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;">#include "petscmat.h"</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"><br>// function to compute outer-product matrix of two vectors</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;">PetscErrorCode MyMPIVecOuterProd(Mat OP, Vec a, Vec b, InsertMode addV)</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;">{</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /*</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> for vectors a and b, computes a.transpose(b) and <br>
adds/stores the resulting</span> <span style="font-family: courier new,monospace;">outer-product matrix<br> to/in the matrix OP</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"></span><span style="font-family: courier new,monospace;"> */ </span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"></span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInt nRows = 0; /* number of rows of OP matrix */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscInt nCols = 0; /* number of columns of OP matrix */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInt nA = 0; /* length of vector a */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscInt nB = 0; /* length of vector b */</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscScalar *locAVals; /* array to hold local vector a values */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscScalar *locSAVals; /* array to hold scaled local vector a values */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscScalar *locBVals; /* array to hold local vector b values */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscScalar *bArr; /* array to hold the entire vector b */</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInt aLow,aHigh,bLow,bHigh; /* local index-range limits */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscInt ia,ib; /* for-loop index variables for a and b vectors */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInt * locRowIdxOP; /* locally-set OP column indices for MatSetValues */</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscScalar curBVal; /* value of vector b's component with which to scale vector a */</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* get the dimension of vector a */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecGetSize(a,&nA);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* get the dimension of vector b */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecGetSize(b,&nB);</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">/* get the dimensions of outer-product matrix */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> MatGetSize(OP,&nRows,&nCols);</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* check for dimensional compatibility */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> if ((nRows != nA) || (nCols != nB)) {</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> SETERRQ(1,"Error: MyMPIVecOuterProd: Dimensional Incompatibility!");</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> return(1);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> }</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* --------------------------------------------- */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* first, convert vector b into array of scalars */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* --------------------------------------------- */</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* allocate memory for array-representation of vector b */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscMalloc((nB)*sizeof(PetscScalar),&bArr);</span><span style="font-family: courier new,monospace;"><br><br style="font-family: courier new,monospace;"></span><span style="font-family: courier new,monospace;"> /* do local assignment from vector b values to array-representation */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* first, obtain local range of vector b */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecGetOwnershipRange(b,&bLow,&bHigh); // bHigh is one more than highest local index</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">/* then, obtain pointer to local elements of vector b */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> VecGetArray(b,&locBVals);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">/* then, assign local values of vector b to corresponding locations in bArr */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> for (ib = bLow; ib < bHigh; ib++) {</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> *(bArr + ib) = *(locBVals + ib - bLow);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> } // end of b for loop</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* finally, restore local elements of vector b */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> VecRestoreArray(b,&locBVals);</span><br><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* ------------------------------------------------- */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* next, scale local values of vector a and add them */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* to corresponding locations in the OP matrix */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* ------------------------------------------------- */<br><br style="font-family: courier new,monospace;"></span><span style="font-family: courier new,monospace;"> /* first, obtain local range of vector a */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> VecGetOwnershipRange(a,&aLow,&aHigh); // aHigh is one more than highest local index</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">/* then, obtain pointer to local elements of vector a */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> VecGetArray(a,&locAVals);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* allocate memory for local array-of-scaled-vector-a-values */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscMalloc((aHigh-aLow)*sizeof(PetscScalar),&locSAVals);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* allocate memory for locally-set OP row indices */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscMalloc((aHigh-aLow)*sizeof(PetscInt),&locRowIdxOP);</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* set locally-set OP row indices */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> for (ia = 0; ia < aHigh-aLow; ia++) {</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> *(locRowIdxOP + ia) = ia + aLow;</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> } // end of for : set locally-set OP row indices</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* next, for each component of vector b (bArr), scale local vector a values */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* and set them up in the corresponding locations in the OP matrix */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> for (ib = 0; ib < nB; ib++) {</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* get component of vector b to scale with */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> curBVal = *(bArr + ib);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* scale vector a local values */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> for (ia = 0; ia < aHigh-aLow; ia++) {</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> *(locSAVals + ia) = (*(locAVals + ia)) * curBVal;</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> } // end of for: scale local vector a values</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* set scaled values in appropriate locations in OP */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> MatSetValues(OP,aHigh-aLow,locRowIdxOP,1,&ib,locSAVals,addV)</span><span style="font-family: courier new,monospace;"></span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> } // end of for: set scaled values in OP matrix</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* next, restore local elements of vector a */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> VecRestoreArray(a,&locAVals);</span><span style="font-family: courier new,monospace;"></span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* free memory for local row indices of OP */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscFree(locRowIdxOP);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* free memory for local scaled vector a values */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscFree(locSAVals);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* free memory for array representation of vector b */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscFree(bArr);</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* ------------------------------- */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> /* finally, assemble the OP matrix */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* ------------------------------- */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> MatAssemblyBegin(OP,MAT_FINAL_ASSEMBLY);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"></span><span style="font-family: courier new,monospace;"> MatAssemblyEnd(OP,MAT_FINAL_ASSEMBLY);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"></span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> return(0);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;">}</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;">int main (int argc, char **argv) {</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"></span><span style="font-family: courier new,monospace;"> PetscInt n1=5; /* first vector's dim */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInt n2=3; /* second vector's dim */</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> Vec v1; /* n1x1 supervector */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> Vec v2; /* n2x1 reduced-dimension representation of y */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> Mat OP1; /* n1xn2 outer-product matrix */</span><span style="font-family: courier new,monospace;"></span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscInt ii; /* for-loop index variables */</span><br><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInt * locRowIdx1; /* locally-set row indices for VecSetValues */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscScalar * locRowVals1; /* array to hold local vector values */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInt low1, high1; /* variables to get local-indices' range */</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInt * locRowIdx2; /* locally-set row indices for VecSetValues */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscScalar * locRowVals2; /* array to hold local vector values */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInt low2, high2; /* variables to get local-indices' range */</span><br>
<br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscInitialize(&argc,&argv,(char*)0,helpMsg);</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecCreateMPI(PETSC_COMM_WORLD,PETSC_DECIDE,n1,&v1);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">VecCreateMPI(PETSC_COMM_WORLD,PETSC_DECIDE,n2,&v2);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">VecZeroEntries(v1);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">VecZeroEntries(v2);</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">MatCreateMPIDense(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,n1,n2,PETSC_NULL,&OP1);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">MatZeroEntries(OP1);</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecGetOwnershipRange(v1,&low1,&high1);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">PetscMalloc((high1-low1)*sizeof(PetscInt),&locRowIdx1);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">PetscMalloc((high1-low1)*sizeof(PetscScalar),&locRowVals1);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">/* set locally-set indices and values */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> for (ii = 0; ii < high1-low1; ii++) {</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> *(locRowIdx1 + ii) = ii + low1;</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> *(locRowVals1 + ii) = pow(2,ii+low1+1);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> } // end of for : set locally-set indices and values</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecSetValues(v1,high1-low1,locRowIdx1,locRowVals1,INSERT_VALUES);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">PetscFree(locRowIdx1); errCode += PetscFree(locRowVals1);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">/* now, assemble the vector v1 */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecAssemblyBegin(v1);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">VecAssemblyEnd(v1);</span><br><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">VecGetOwnershipRange(v2,&low2,&high2);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">PetscMalloc((high2-low2)*sizeof(PetscInt),&locRowIdx2);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">PetscMalloc((high2-low2)*sizeof(PetscScalar),&locRowVals2);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">/* set locally-set indices and values */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> for (ii = 0; ii < high2-low2; ii++) {</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> *(locRowIdx2 + ii) = ii + low2;</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> *(locRowVals2 + ii) = 1.0 / pow(2,ii+low2);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> } // end of for : set locally-set indices and values</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecSetValues(v2,high2-low2,locRowIdx2,locRowVals2,INSERT_VALUES);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">PetscFree(locRowIdx2); errCode += PetscFree(locRowVals2);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">/* now, assemble the vector v2 */</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecAssemblyBegin(v2);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">VecAssemblyEnd(v2);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"></span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> MyMPIVecOuterProd(OP1, v1, v2, INSERT_VALUES);</span><span style="font-family: courier new,monospace;"></span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscPrintf(PETSC_COMM_WORLD,"---- vector v1 ----\n");</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> VecView(v1,PETSC_VIEWER_STDOUT_WORLD);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"></span><span style="font-family: courier new,monospace;"> PetscPrintf(PETSC_COMM_WORLD,"---- vector v1 ----\n");</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> PetscPrintf(PETSC_COMM_WORLD,"---- vector v2 ----\n");</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> VecView(v2,PETSC_VIEWER_STDOUT_WORLD);</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">PetscPrintf(PETSC_COMM_WORLD,"---- vector v2 ----\n");</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscPrintf(PETSC_COMM_WORLD,"---- matrix OP1 ----\n");</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> MatView(OP1,PETSC_VIEWER_STDOUT_WORLD);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> </span><span style="font-family: courier new,monospace;">PetscPrintf(PETSC_COMM_WORLD,"---- matrix OP1 ----\n");</span><span style="font-family: courier new,monospace;"></span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* destroy OP1 */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> MatDestroy(OP1);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"></span><span style="font-family: courier new,monospace;"> /* destroy v1 */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> VecDestroy(v1);</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> /* destroy v2 */</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> VecDestroy(v2);</span><br style="font-family: courier new,monospace;"><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;"> PetscFinalize();</span><br style="font-family: courier new,monospace;">
<span style="font-family: courier new,monospace;"> return 0;</span><br style="font-family: courier new,monospace;"><span style="font-family: courier new,monospace;">}</span><br><br>[hsharma@ifp-32]$ petscmpiexec -np 3 ./OuterProductCheckOutput<br>
<br><span style="font-family: verdana,sans-serif;">---- vector v1 ----</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">Process [0]</span><br style="font-family: verdana,sans-serif;">
<span style="font-family: verdana,sans-serif;">2</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">4</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">Process [1]</span><br style="font-family: verdana,sans-serif;">
<span style="font-family: verdana,sans-serif;">8</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">16</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">Process [2]</span><br style="font-family: verdana,sans-serif;">
<span style="font-family: verdana,sans-serif;">32</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">---- vector v1 ----</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">---- vector v2 ----</span><br style="font-family: verdana,sans-serif;">
<span style="font-family: verdana,sans-serif;">Process [0]</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">1</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">Process [1]</span><br style="font-family: verdana,sans-serif;">
<span style="font-family: verdana,sans-serif;">0.5</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">Process [2]</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">0.25</span><br style="font-family: verdana,sans-serif;">
<span style="font-family: verdana,sans-serif;">---- vector v2 ----</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">---- matrix OP1 ----</span><br style="font-family: verdana,sans-serif;">
<span style="font-family: verdana,sans-serif;">2.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">4.0000000000000000e+00 0.0000000000000000e+00 0.0000000000000000e+00</span><br style="font-family: verdana,sans-serif;">
<span style="font-family: verdana,sans-serif;">0.0000000000000000e+00 4.0000000000000000e+00 0.0000000000000000e+00</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">0.0000000000000000e+00 8.0000000000000000e+00 0.0000000000000000e+00</span><br style="font-family: verdana,sans-serif;">
<span style="font-family: verdana,sans-serif;">0.0000000000000000e+00 0.0000000000000000e+00 8.0000000000000000e+00</span><br style="font-family: verdana,sans-serif;"><span style="font-family: verdana,sans-serif;">---- matrix OP1 ----</span><br style="font-family: courier new,monospace;">
<br style="font-family: courier new,monospace;">---- code and sample output ----<br><br>