1 2 /* 3 This provides a matrix that applies a VecScatter to a vector. 4 */ 5 6 #include <petsc/private/matimpl.h> /*I "petscmat.h" I*/ 7 #include <petsc/private/vecimpl.h> 8 9 typedef struct { 10 VecScatter scatter; 11 } Mat_Scatter; 12 13 /*@ 14 MatScatterGetVecScatter - Returns the user-provided scatter set with MatScatterSetVecScatter() 15 16 Not Collective, but not cannot use scatter if not used collectively on Mat 17 18 Input Parameter: 19 . mat - the matrix, should have been created with MatCreateScatter() or have type MATSCATTER 20 21 Output Parameter: 22 . scatter - the scatter context 23 24 Level: intermediate 25 26 .seealso: MatCreateScatter(), MatScatterSetVecScatter(), MATSCATTER 27 @*/ 28 PetscErrorCode MatScatterGetVecScatter(Mat mat,VecScatter *scatter) 29 { 30 Mat_Scatter *mscatter; 31 32 PetscFunctionBegin; 33 PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 34 PetscValidPointer(scatter,2); 35 mscatter = (Mat_Scatter*)mat->data; 36 *scatter = mscatter->scatter; 37 PetscFunctionReturn(0); 38 } 39 40 PetscErrorCode MatDestroy_Scatter(Mat mat) 41 { 42 PetscErrorCode ierr; 43 Mat_Scatter *scatter = (Mat_Scatter*)mat->data; 44 45 PetscFunctionBegin; 46 ierr = VecScatterDestroy(&scatter->scatter);CHKERRQ(ierr); 47 ierr = PetscFree(mat->data);CHKERRQ(ierr); 48 PetscFunctionReturn(0); 49 } 50 51 PetscErrorCode MatMult_Scatter(Mat A,Vec x,Vec y) 52 { 53 Mat_Scatter *scatter = (Mat_Scatter*)A->data; 54 PetscErrorCode ierr; 55 56 PetscFunctionBegin; 57 if (!scatter->scatter) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_WRONGSTATE,"Need to first call MatScatterSetScatter()"); 58 ierr = VecZeroEntries(y);CHKERRQ(ierr); 59 ierr = VecScatterBegin(scatter->scatter,x,y,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 60 ierr = VecScatterEnd(scatter->scatter,x,y,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 61 PetscFunctionReturn(0); 62 } 63 64 PetscErrorCode MatMultAdd_Scatter(Mat A,Vec x,Vec y,Vec z) 65 { 66 Mat_Scatter *scatter = (Mat_Scatter*)A->data; 67 PetscErrorCode ierr; 68 69 PetscFunctionBegin; 70 if (!scatter->scatter) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_WRONGSTATE,"Need to first call MatScatterSetScatter()"); 71 if (z != y) {ierr = VecCopy(y,z);CHKERRQ(ierr);} 72 ierr = VecScatterBegin(scatter->scatter,x,z,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 73 ierr = VecScatterEnd(scatter->scatter,x,z,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 74 PetscFunctionReturn(0); 75 } 76 77 PetscErrorCode MatMultTranspose_Scatter(Mat A,Vec x,Vec y) 78 { 79 Mat_Scatter *scatter = (Mat_Scatter*)A->data; 80 PetscErrorCode ierr; 81 82 PetscFunctionBegin; 83 if (!scatter->scatter) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_WRONGSTATE,"Need to first call MatScatterSetScatter()"); 84 ierr = VecZeroEntries(y);CHKERRQ(ierr); 85 ierr = VecScatterBegin(scatter->scatter,x,y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 86 ierr = VecScatterEnd(scatter->scatter,x,y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 87 PetscFunctionReturn(0); 88 } 89 90 PetscErrorCode MatMultTransposeAdd_Scatter(Mat A,Vec x,Vec y,Vec z) 91 { 92 Mat_Scatter *scatter = (Mat_Scatter*)A->data; 93 PetscErrorCode ierr; 94 95 PetscFunctionBegin; 96 if (!scatter->scatter) SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_WRONGSTATE,"Need to first call MatScatterSetScatter()"); 97 if (z != y) {ierr = VecCopy(y,z);CHKERRQ(ierr);} 98 ierr = VecScatterBegin(scatter->scatter,x,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 99 ierr = VecScatterEnd(scatter->scatter,x,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 100 PetscFunctionReturn(0); 101 } 102 103 static struct _MatOps MatOps_Values = {0, 104 0, 105 0, 106 MatMult_Scatter, 107 /* 4*/ MatMultAdd_Scatter, 108 MatMultTranspose_Scatter, 109 MatMultTransposeAdd_Scatter, 110 0, 111 0, 112 0, 113 /* 10*/ 0, 114 0, 115 0, 116 0, 117 0, 118 /* 15*/ 0, 119 0, 120 0, 121 0, 122 0, 123 /* 20*/ 0, 124 0, 125 0, 126 0, 127 /* 24*/ 0, 128 0, 129 0, 130 0, 131 0, 132 /* 29*/ 0, 133 0, 134 0, 135 0, 136 0, 137 /* 34*/ 0, 138 0, 139 0, 140 0, 141 0, 142 /* 39*/ 0, 143 0, 144 0, 145 0, 146 0, 147 /* 44*/ 0, 148 0, 149 MatShift_Basic, 150 0, 151 0, 152 /* 49*/ 0, 153 0, 154 0, 155 0, 156 0, 157 /* 54*/ 0, 158 0, 159 0, 160 0, 161 0, 162 /* 59*/ 0, 163 MatDestroy_Scatter, 164 0, 165 0, 166 0, 167 /* 64*/ 0, 168 0, 169 0, 170 0, 171 0, 172 /* 69*/ 0, 173 0, 174 0, 175 0, 176 0, 177 /* 74*/ 0, 178 0, 179 0, 180 0, 181 0, 182 /* 79*/ 0, 183 0, 184 0, 185 0, 186 0, 187 /* 84*/ 0, 188 0, 189 0, 190 0, 191 0, 192 /* 89*/ 0, 193 0, 194 0, 195 0, 196 0, 197 /* 94*/ 0, 198 0, 199 0, 200 0, 201 0, 202 /*99*/ 0, 203 0, 204 0, 205 0, 206 0, 207 /*104*/ 0, 208 0, 209 0, 210 0, 211 0, 212 /*109*/ 0, 213 0, 214 0, 215 0, 216 0, 217 /*114*/ 0, 218 0, 219 0, 220 0, 221 0, 222 /*119*/ 0, 223 0, 224 0, 225 0, 226 0, 227 /*124*/ 0, 228 0, 229 0, 230 0, 231 0, 232 /*129*/ 0, 233 0, 234 0, 235 0, 236 0, 237 /*134*/ 0, 238 0, 239 0, 240 0, 241 0, 242 /*139*/ 0, 243 0, 244 0 245 }; 246 247 /*MC 248 MATSCATTER - MATSCATTER = "scatter" - A matrix type that simply applies a VecScatterBegin/End() 249 250 Level: advanced 251 252 .seealso: MatCreateScatter(), MatScatterSetVecScatter(), MatScatterGetVecScatter() 253 254 M*/ 255 256 PETSC_EXTERN PetscErrorCode MatCreate_Scatter(Mat A) 257 { 258 Mat_Scatter *b; 259 PetscErrorCode ierr; 260 261 PetscFunctionBegin; 262 ierr = PetscMemcpy(A->ops,&MatOps_Values,sizeof(struct _MatOps));CHKERRQ(ierr); 263 ierr = PetscNewLog(A,&b);CHKERRQ(ierr); 264 265 A->data = (void*)b; 266 267 ierr = PetscLayoutSetUp(A->rmap);CHKERRQ(ierr); 268 ierr = PetscLayoutSetUp(A->cmap);CHKERRQ(ierr); 269 270 A->assembled = PETSC_TRUE; 271 A->preallocated = PETSC_FALSE; 272 273 ierr = PetscObjectChangeTypeName((PetscObject)A,MATSCATTER);CHKERRQ(ierr); 274 PetscFunctionReturn(0); 275 } 276 277 #include <petsc/private/vecscatterimpl.h> 278 /*@C 279 MatCreateScatter - Creates a new matrix based on a VecScatter 280 281 Collective on MPI_Comm 282 283 Input Parameters: 284 + comm - MPI communicator 285 - scatter - a VecScatterContext 286 287 Output Parameter: 288 . A - the matrix 289 290 Level: intermediate 291 292 PETSc requires that matrices and vectors being used for certain 293 operations are partitioned accordingly. For example, when 294 creating a scatter matrix, A, that supports parallel matrix-vector 295 products using MatMult(A,x,y) the user should set the number 296 of local matrix rows to be the number of local elements of the 297 corresponding result vector, y. Note that this is information is 298 required for use of the matrix interface routines, even though 299 the scatter matrix may not actually be physically partitioned. 300 301 Developer Notes: This directly accesses information inside the VecScatter associated with the matrix-vector product 302 for this matrix. This is not desirable.. 303 304 305 .seealso: MatScatterSetVecScatter(), MatScatterGetVecScatter(), MATSCATTER 306 @*/ 307 PetscErrorCode MatCreateScatter(MPI_Comm comm,VecScatter scatter,Mat *A) 308 { 309 PetscErrorCode ierr; 310 311 PetscFunctionBegin; 312 ierr = MatCreate(comm,A);CHKERRQ(ierr); 313 ierr = MatSetSizes(*A,scatter->to_n,scatter->from_n,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr); 314 ierr = MatSetType(*A,MATSCATTER);CHKERRQ(ierr); 315 ierr = MatScatterSetVecScatter(*A,scatter);CHKERRQ(ierr); 316 ierr = MatSetUp(*A);CHKERRQ(ierr); 317 PetscFunctionReturn(0); 318 } 319 320 /*@ 321 MatScatterSetVecScatter - sets that scatter that the matrix is to apply as its linear operator 322 323 Collective on Mat 324 325 Input Parameters: 326 + mat - the scatter matrix 327 - scatter - the scatter context create with VecScatterCreate() 328 329 Level: advanced 330 331 332 .seealso: MatCreateScatter(), MATSCATTER 333 @*/ 334 PetscErrorCode MatScatterSetVecScatter(Mat mat,VecScatter scatter) 335 { 336 Mat_Scatter *mscatter = (Mat_Scatter*)mat->data; 337 PetscErrorCode ierr; 338 339 PetscFunctionBegin; 340 PetscValidHeaderSpecific(mat,MAT_CLASSID,1); 341 PetscValidHeaderSpecific(scatter,VEC_SCATTER_CLASSID,2); 342 PetscCheckSameComm((PetscObject)scatter,1,(PetscObject)mat,2); 343 if (mat->rmap->n != scatter->to_n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Number of local rows in matrix %D not equal local scatter size %D",mat->rmap->n,scatter->to_n); 344 if (mat->cmap->n != scatter->from_n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Number of local columns in matrix %D not equal local scatter size %D",mat->cmap->n,scatter->from_n); 345 346 ierr = PetscObjectReference((PetscObject)scatter);CHKERRQ(ierr); 347 ierr = VecScatterDestroy(&mscatter->scatter);CHKERRQ(ierr); 348 349 mscatter->scatter = scatter; 350 PetscFunctionReturn(0); 351 } 352 353 354