12a6744ebSBarry Smith 22a6744ebSBarry Smith /* 32a6744ebSBarry Smith This provides a matrix that applies a VecScatter to a vector. 42a6744ebSBarry Smith */ 52a6744ebSBarry Smith 6af0996ceSBarry Smith #include <petsc/private/matimpl.h> /*I "petscmat.h" I*/ 7af0996ceSBarry Smith #include <petsc/private/vecimpl.h> 82a6744ebSBarry Smith 92a6744ebSBarry Smith typedef struct { 102a6744ebSBarry Smith VecScatter scatter; 112a6744ebSBarry Smith } Mat_Scatter; 122a6744ebSBarry Smith 132a6744ebSBarry Smith /*@ 142a6744ebSBarry Smith MatScatterGetVecScatter - Returns the user-provided scatter set with MatScatterSetVecScatter() 152a6744ebSBarry Smith 162a6744ebSBarry Smith Not Collective, but not cannot use scatter if not used collectively on Mat 172a6744ebSBarry Smith 182a6744ebSBarry Smith Input Parameter: 192a6744ebSBarry Smith . mat - the matrix, should have been created with MatCreateScatter() or have type MATSCATTER 202a6744ebSBarry Smith 212a6744ebSBarry Smith Output Parameter: 222a6744ebSBarry Smith . scatter - the scatter context 232a6744ebSBarry Smith 242a6744ebSBarry Smith Level: intermediate 252a6744ebSBarry Smith 262a6744ebSBarry Smith .seealso: MatCreateScatter(), MatScatterSetVecScatter(), MATSCATTER 272a6744ebSBarry Smith @*/ 287087cfbeSBarry Smith PetscErrorCode MatScatterGetVecScatter(Mat mat,VecScatter *scatter) 292a6744ebSBarry Smith { 302a6744ebSBarry Smith Mat_Scatter *mscatter; 312a6744ebSBarry Smith 322a6744ebSBarry Smith PetscFunctionBegin; 330700a824SBarry Smith PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 342a6744ebSBarry Smith PetscValidPointer(scatter,2); 352a6744ebSBarry Smith mscatter = (Mat_Scatter*)mat->data; 362a6744ebSBarry Smith *scatter = mscatter->scatter; 372a6744ebSBarry Smith PetscFunctionReturn(0); 382a6744ebSBarry Smith } 392a6744ebSBarry Smith 402a6744ebSBarry Smith PetscErrorCode MatDestroy_Scatter(Mat mat) 412a6744ebSBarry Smith { 422a6744ebSBarry Smith PetscErrorCode ierr; 432a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter*)mat->data; 442a6744ebSBarry Smith 452a6744ebSBarry Smith PetscFunctionBegin; 466bf464f9SBarry Smith ierr = VecScatterDestroy(&scatter->scatter);CHKERRQ(ierr); 476bf464f9SBarry Smith ierr = PetscFree(mat->data);CHKERRQ(ierr); 482a6744ebSBarry Smith PetscFunctionReturn(0); 492a6744ebSBarry Smith } 502a6744ebSBarry Smith 512a6744ebSBarry Smith PetscErrorCode MatMult_Scatter(Mat A,Vec x,Vec y) 522a6744ebSBarry Smith { 532a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter*)A->data; 542a6744ebSBarry Smith PetscErrorCode ierr; 552a6744ebSBarry Smith 562a6744ebSBarry Smith PetscFunctionBegin; 572c71b3e2SJacob Faibussowitsch PetscCheckFalse(!scatter->scatter,PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_WRONGSTATE,"Need to first call MatScatterSetScatter()"); 5874d95942SJed Brown ierr = VecZeroEntries(y);CHKERRQ(ierr); 5974d95942SJed Brown ierr = VecScatterBegin(scatter->scatter,x,y,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 6074d95942SJed Brown ierr = VecScatterEnd(scatter->scatter,x,y,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 612a6744ebSBarry Smith PetscFunctionReturn(0); 622a6744ebSBarry Smith } 632a6744ebSBarry Smith 642a6744ebSBarry Smith PetscErrorCode MatMultAdd_Scatter(Mat A,Vec x,Vec y,Vec z) 652a6744ebSBarry Smith { 662a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter*)A->data; 672a6744ebSBarry Smith PetscErrorCode ierr; 682a6744ebSBarry Smith 692a6744ebSBarry Smith PetscFunctionBegin; 702c71b3e2SJacob Faibussowitsch PetscCheckFalse(!scatter->scatter,PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_WRONGSTATE,"Need to first call MatScatterSetScatter()"); 712a6744ebSBarry Smith if (z != y) {ierr = VecCopy(y,z);CHKERRQ(ierr);} 72ca9f406cSSatish Balay ierr = VecScatterBegin(scatter->scatter,x,z,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 73ca9f406cSSatish Balay ierr = VecScatterEnd(scatter->scatter,x,z,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 742a6744ebSBarry Smith PetscFunctionReturn(0); 752a6744ebSBarry Smith } 762a6744ebSBarry Smith 772a6744ebSBarry Smith PetscErrorCode MatMultTranspose_Scatter(Mat A,Vec x,Vec y) 782a6744ebSBarry Smith { 792a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter*)A->data; 802a6744ebSBarry Smith PetscErrorCode ierr; 812a6744ebSBarry Smith 822a6744ebSBarry Smith PetscFunctionBegin; 832c71b3e2SJacob Faibussowitsch PetscCheckFalse(!scatter->scatter,PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_WRONGSTATE,"Need to first call MatScatterSetScatter()"); 8474d95942SJed Brown ierr = VecZeroEntries(y);CHKERRQ(ierr); 8574d95942SJed Brown ierr = VecScatterBegin(scatter->scatter,x,y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 8674d95942SJed Brown ierr = VecScatterEnd(scatter->scatter,x,y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 872a6744ebSBarry Smith PetscFunctionReturn(0); 882a6744ebSBarry Smith } 892a6744ebSBarry Smith 902a6744ebSBarry Smith PetscErrorCode MatMultTransposeAdd_Scatter(Mat A,Vec x,Vec y,Vec z) 912a6744ebSBarry Smith { 922a6744ebSBarry Smith Mat_Scatter *scatter = (Mat_Scatter*)A->data; 932a6744ebSBarry Smith PetscErrorCode ierr; 942a6744ebSBarry Smith 952a6744ebSBarry Smith PetscFunctionBegin; 962c71b3e2SJacob Faibussowitsch PetscCheckFalse(!scatter->scatter,PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_WRONGSTATE,"Need to first call MatScatterSetScatter()"); 972a6744ebSBarry Smith if (z != y) {ierr = VecCopy(y,z);CHKERRQ(ierr);} 98ca9f406cSSatish Balay ierr = VecScatterBegin(scatter->scatter,x,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 99ca9f406cSSatish Balay ierr = VecScatterEnd(scatter->scatter,x,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1002a6744ebSBarry Smith PetscFunctionReturn(0); 1012a6744ebSBarry Smith } 1022a6744ebSBarry Smith 103f4259b30SLisandro Dalcin static struct _MatOps MatOps_Values = {NULL, 104f4259b30SLisandro Dalcin NULL, 105f4259b30SLisandro Dalcin NULL, 1062a6744ebSBarry Smith MatMult_Scatter, 1072a6744ebSBarry Smith /* 4*/ MatMultAdd_Scatter, 1082a6744ebSBarry Smith MatMultTranspose_Scatter, 1092a6744ebSBarry Smith MatMultTransposeAdd_Scatter, 110f4259b30SLisandro Dalcin NULL, 111f4259b30SLisandro Dalcin NULL, 112f4259b30SLisandro Dalcin NULL, 113f4259b30SLisandro Dalcin /* 10*/ NULL, 114f4259b30SLisandro Dalcin NULL, 115f4259b30SLisandro Dalcin NULL, 116f4259b30SLisandro Dalcin NULL, 117f4259b30SLisandro Dalcin NULL, 118f4259b30SLisandro Dalcin /* 15*/ NULL, 119f4259b30SLisandro Dalcin NULL, 120f4259b30SLisandro Dalcin NULL, 121f4259b30SLisandro Dalcin NULL, 122f4259b30SLisandro Dalcin NULL, 123f4259b30SLisandro Dalcin /* 20*/ NULL, 124f4259b30SLisandro Dalcin NULL, 125f4259b30SLisandro Dalcin NULL, 126f4259b30SLisandro Dalcin NULL, 127f4259b30SLisandro Dalcin /* 24*/ NULL, 128f4259b30SLisandro Dalcin NULL, 129f4259b30SLisandro Dalcin NULL, 130f4259b30SLisandro Dalcin NULL, 131f4259b30SLisandro Dalcin NULL, 132f4259b30SLisandro Dalcin /* 29*/ NULL, 133f4259b30SLisandro Dalcin NULL, 134f4259b30SLisandro Dalcin NULL, 135f4259b30SLisandro Dalcin NULL, 136f4259b30SLisandro Dalcin NULL, 137f4259b30SLisandro Dalcin /* 34*/ NULL, 138f4259b30SLisandro Dalcin NULL, 139f4259b30SLisandro Dalcin NULL, 140f4259b30SLisandro Dalcin NULL, 141f4259b30SLisandro Dalcin NULL, 142f4259b30SLisandro Dalcin /* 39*/ NULL, 143f4259b30SLisandro Dalcin NULL, 144f4259b30SLisandro Dalcin NULL, 145f4259b30SLisandro Dalcin NULL, 146f4259b30SLisandro Dalcin NULL, 147f4259b30SLisandro Dalcin /* 44*/ NULL, 148f4259b30SLisandro Dalcin NULL, 1497d68702bSBarry Smith MatShift_Basic, 150f4259b30SLisandro Dalcin NULL, 151f4259b30SLisandro Dalcin NULL, 152f4259b30SLisandro Dalcin /* 49*/ NULL, 153f4259b30SLisandro Dalcin NULL, 154f4259b30SLisandro Dalcin NULL, 155f4259b30SLisandro Dalcin NULL, 156f4259b30SLisandro Dalcin NULL, 157f4259b30SLisandro Dalcin /* 54*/ NULL, 158f4259b30SLisandro Dalcin NULL, 159f4259b30SLisandro Dalcin NULL, 160f4259b30SLisandro Dalcin NULL, 161f4259b30SLisandro Dalcin NULL, 162f4259b30SLisandro Dalcin /* 59*/ NULL, 1632a6744ebSBarry Smith MatDestroy_Scatter, 164f4259b30SLisandro Dalcin NULL, 165f4259b30SLisandro Dalcin NULL, 166f4259b30SLisandro Dalcin NULL, 167f4259b30SLisandro Dalcin /* 64*/ NULL, 168f4259b30SLisandro Dalcin NULL, 169f4259b30SLisandro Dalcin NULL, 170f4259b30SLisandro Dalcin NULL, 171f4259b30SLisandro Dalcin NULL, 172f4259b30SLisandro Dalcin /* 69*/ NULL, 173f4259b30SLisandro Dalcin NULL, 174f4259b30SLisandro Dalcin NULL, 175f4259b30SLisandro Dalcin NULL, 176f4259b30SLisandro Dalcin NULL, 177f4259b30SLisandro Dalcin /* 74*/ NULL, 178f4259b30SLisandro Dalcin NULL, 179f4259b30SLisandro Dalcin NULL, 180f4259b30SLisandro Dalcin NULL, 181f4259b30SLisandro Dalcin NULL, 182f4259b30SLisandro Dalcin /* 79*/ NULL, 183f4259b30SLisandro Dalcin NULL, 184f4259b30SLisandro Dalcin NULL, 185f4259b30SLisandro Dalcin NULL, 186f4259b30SLisandro Dalcin NULL, 187f4259b30SLisandro Dalcin /* 84*/ NULL, 188f4259b30SLisandro Dalcin NULL, 189f4259b30SLisandro Dalcin NULL, 190f4259b30SLisandro Dalcin NULL, 191f4259b30SLisandro Dalcin NULL, 192f4259b30SLisandro Dalcin /* 89*/ NULL, 193f4259b30SLisandro Dalcin NULL, 194f4259b30SLisandro Dalcin NULL, 195f4259b30SLisandro Dalcin NULL, 196f4259b30SLisandro Dalcin NULL, 197f4259b30SLisandro Dalcin /* 94*/ NULL, 198f4259b30SLisandro Dalcin NULL, 199f4259b30SLisandro Dalcin NULL, 200f4259b30SLisandro Dalcin NULL, 201f4259b30SLisandro Dalcin NULL, 202f4259b30SLisandro Dalcin /*99*/ NULL, 203f4259b30SLisandro Dalcin NULL, 204f4259b30SLisandro Dalcin NULL, 205f4259b30SLisandro Dalcin NULL, 206f4259b30SLisandro Dalcin NULL, 207f4259b30SLisandro Dalcin /*104*/ NULL, 208f4259b30SLisandro Dalcin NULL, 209f4259b30SLisandro Dalcin NULL, 210f4259b30SLisandro Dalcin NULL, 211f4259b30SLisandro Dalcin NULL, 212f4259b30SLisandro Dalcin /*109*/ NULL, 213f4259b30SLisandro Dalcin NULL, 214f4259b30SLisandro Dalcin NULL, 215f4259b30SLisandro Dalcin NULL, 216f4259b30SLisandro Dalcin NULL, 217f4259b30SLisandro Dalcin /*114*/ NULL, 218f4259b30SLisandro Dalcin NULL, 219f4259b30SLisandro Dalcin NULL, 220f4259b30SLisandro Dalcin NULL, 221f4259b30SLisandro Dalcin NULL, 222f4259b30SLisandro Dalcin /*119*/ NULL, 223f4259b30SLisandro Dalcin NULL, 224f4259b30SLisandro Dalcin NULL, 225f4259b30SLisandro Dalcin NULL, 226f4259b30SLisandro Dalcin NULL, 227f4259b30SLisandro Dalcin /*124*/ NULL, 228f4259b30SLisandro Dalcin NULL, 229f4259b30SLisandro Dalcin NULL, 230f4259b30SLisandro Dalcin NULL, 231f4259b30SLisandro Dalcin NULL, 232f4259b30SLisandro Dalcin /*129*/ NULL, 233f4259b30SLisandro Dalcin NULL, 234f4259b30SLisandro Dalcin NULL, 235f4259b30SLisandro Dalcin NULL, 236f4259b30SLisandro Dalcin NULL, 237f4259b30SLisandro Dalcin /*134*/ NULL, 238f4259b30SLisandro Dalcin NULL, 239f4259b30SLisandro Dalcin NULL, 240f4259b30SLisandro Dalcin NULL, 241f4259b30SLisandro Dalcin NULL, 242f4259b30SLisandro Dalcin /*139*/NULL, 243*d70f29a3SPierre Jolivet NULL, 244*d70f29a3SPierre Jolivet NULL, 245*d70f29a3SPierre Jolivet NULL, 246*d70f29a3SPierre Jolivet NULL, 247*d70f29a3SPierre Jolivet /*144*/NULL, 248*d70f29a3SPierre Jolivet NULL, 249f4259b30SLisandro Dalcin NULL, 250f4259b30SLisandro Dalcin NULL 2513964eb88SJed Brown }; 2522a6744ebSBarry Smith 2532a6744ebSBarry Smith /*MC 2544a440b9bSBarry Smith MATSCATTER - MATSCATTER = "scatter" - A matrix type that simply applies a VecScatterBegin/End() 2552a6744ebSBarry Smith 2562a6744ebSBarry Smith Level: advanced 2572a6744ebSBarry Smith 2582a6744ebSBarry Smith .seealso: MatCreateScatter(), MatScatterSetVecScatter(), MatScatterGetVecScatter() 2592a6744ebSBarry Smith 2602a6744ebSBarry Smith M*/ 2612a6744ebSBarry Smith 2628cc058d9SJed Brown PETSC_EXTERN PetscErrorCode MatCreate_Scatter(Mat A) 2632a6744ebSBarry Smith { 2642a6744ebSBarry Smith Mat_Scatter *b; 2652a6744ebSBarry Smith PetscErrorCode ierr; 2662a6744ebSBarry Smith 2672a6744ebSBarry Smith PetscFunctionBegin; 2682a6744ebSBarry Smith ierr = PetscMemcpy(A->ops,&MatOps_Values,sizeof(struct _MatOps));CHKERRQ(ierr); 269b00a9115SJed Brown ierr = PetscNewLog(A,&b);CHKERRQ(ierr); 2702a6744ebSBarry Smith 2712a6744ebSBarry Smith A->data = (void*)b; 2722a6744ebSBarry Smith 27326283091SBarry Smith ierr = PetscLayoutSetUp(A->rmap);CHKERRQ(ierr); 27426283091SBarry Smith ierr = PetscLayoutSetUp(A->cmap);CHKERRQ(ierr); 2752a6744ebSBarry Smith 2762a6744ebSBarry Smith A->assembled = PETSC_TRUE; 2772a6744ebSBarry Smith A->preallocated = PETSC_FALSE; 27874d95942SJed Brown 27974d95942SJed Brown ierr = PetscObjectChangeTypeName((PetscObject)A,MATSCATTER);CHKERRQ(ierr); 2802a6744ebSBarry Smith PetscFunctionReturn(0); 2812a6744ebSBarry Smith } 2822a6744ebSBarry Smith 28397929ea7SJunchao Zhang #include <petsc/private/sfimpl.h> 2842a6744ebSBarry Smith /*@C 2852a6744ebSBarry Smith MatCreateScatter - Creates a new matrix based on a VecScatter 2862a6744ebSBarry Smith 287d083f849SBarry Smith Collective 2882a6744ebSBarry Smith 2892a6744ebSBarry Smith Input Parameters: 2902a6744ebSBarry Smith + comm - MPI communicator 2912a6744ebSBarry Smith - scatter - a VecScatterContext 2922a6744ebSBarry Smith 2932a6744ebSBarry Smith Output Parameter: 2942a6744ebSBarry Smith . A - the matrix 2952a6744ebSBarry Smith 2962a6744ebSBarry Smith Level: intermediate 2972a6744ebSBarry Smith 2982a6744ebSBarry Smith PETSc requires that matrices and vectors being used for certain 2992a6744ebSBarry Smith operations are partitioned accordingly. For example, when 3002a6744ebSBarry Smith creating a scatter matrix, A, that supports parallel matrix-vector 3012a6744ebSBarry Smith products using MatMult(A,x,y) the user should set the number 3022a6744ebSBarry Smith of local matrix rows to be the number of local elements of the 3032a6744ebSBarry Smith corresponding result vector, y. Note that this is information is 3042a6744ebSBarry Smith required for use of the matrix interface routines, even though 3052a6744ebSBarry Smith the scatter matrix may not actually be physically partitioned. 3062a6744ebSBarry Smith 3076eb45d04SBarry Smith Developer Notes: This directly accesses information inside the VecScatter associated with the matrix-vector product 3086eb45d04SBarry Smith for this matrix. This is not desirable.. 3096eb45d04SBarry Smith 3102a6744ebSBarry Smith .seealso: MatScatterSetVecScatter(), MatScatterGetVecScatter(), MATSCATTER 3112a6744ebSBarry Smith @*/ 3127087cfbeSBarry Smith PetscErrorCode MatCreateScatter(MPI_Comm comm,VecScatter scatter,Mat *A) 3132a6744ebSBarry Smith { 3142a6744ebSBarry Smith PetscErrorCode ierr; 3152a6744ebSBarry Smith 3162a6744ebSBarry Smith PetscFunctionBegin; 3172a6744ebSBarry Smith ierr = MatCreate(comm,A);CHKERRQ(ierr); 31897929ea7SJunchao Zhang ierr = MatSetSizes(*A,scatter->vscat.to_n,scatter->vscat.from_n,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr); 3192a6744ebSBarry Smith ierr = MatSetType(*A,MATSCATTER);CHKERRQ(ierr); 3202a6744ebSBarry Smith ierr = MatScatterSetVecScatter(*A,scatter);CHKERRQ(ierr); 32137e1895aSJed Brown ierr = MatSetUp(*A);CHKERRQ(ierr); 3222a6744ebSBarry Smith PetscFunctionReturn(0); 3232a6744ebSBarry Smith } 3242a6744ebSBarry Smith 3252a6744ebSBarry Smith /*@ 3262a6744ebSBarry Smith MatScatterSetVecScatter - sets that scatter that the matrix is to apply as its linear operator 3272a6744ebSBarry Smith 3282a6744ebSBarry Smith Collective on Mat 3292a6744ebSBarry Smith 3302a6744ebSBarry Smith Input Parameters: 3312a6744ebSBarry Smith + mat - the scatter matrix 3329448b7f1SJunchao Zhang - scatter - the scatter context create with VecScatterCreate() 3332a6744ebSBarry Smith 3342a6744ebSBarry Smith Level: advanced 3352a6744ebSBarry Smith 3362a6744ebSBarry Smith .seealso: MatCreateScatter(), MATSCATTER 3372a6744ebSBarry Smith @*/ 3387087cfbeSBarry Smith PetscErrorCode MatScatterSetVecScatter(Mat mat,VecScatter scatter) 3392a6744ebSBarry Smith { 3402a6744ebSBarry Smith Mat_Scatter *mscatter = (Mat_Scatter*)mat->data; 3412a6744ebSBarry Smith PetscErrorCode ierr; 3422a6744ebSBarry Smith 3432a6744ebSBarry Smith PetscFunctionBegin; 3440700a824SBarry Smith PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 34597929ea7SJunchao Zhang PetscValidHeaderSpecific(scatter,PETSCSF_CLASSID,2); 346064a246eSJacob Faibussowitsch PetscCheckSameComm((PetscObject)scatter,2,(PetscObject)mat,1); 3472c71b3e2SJacob Faibussowitsch PetscCheckFalse(mat->rmap->n != scatter->vscat.to_n,PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Number of local rows in matrix %" PetscInt_FMT " not equal local scatter size %" PetscInt_FMT,mat->rmap->n,scatter->vscat.to_n); 3482c71b3e2SJacob Faibussowitsch PetscCheckFalse(mat->cmap->n != scatter->vscat.from_n,PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Number of local columns in matrix %" PetscInt_FMT " not equal local scatter size %" PetscInt_FMT,mat->cmap->n,scatter->vscat.from_n); 3492a6744ebSBarry Smith 350c3122656SLisandro Dalcin ierr = PetscObjectReference((PetscObject)scatter);CHKERRQ(ierr); 3516bf464f9SBarry Smith ierr = VecScatterDestroy(&mscatter->scatter);CHKERRQ(ierr); 35226fbe8dcSKarl Rupp 3532a6744ebSBarry Smith mscatter->scatter = scatter; 3542a6744ebSBarry Smith PetscFunctionReturn(0); 3552a6744ebSBarry Smith } 3562a6744ebSBarry Smith 357