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