1 /* TODOLIST 2 3 Solvers 4 - Add support for cholesky for coarse solver (similar to local solvers) 5 - Propagate ksp prefixes for solvers to mat objects? 6 7 User interface 8 - ** DM attached to pc? 9 10 Debugging output 11 - * Better management of verbosity levels of debugging output 12 13 Extra 14 - *** Is it possible to work with PCBDDCGraph on boundary indices only (less memory consumed)? 15 - BDDC with MG framework? 16 17 FETIDP 18 - Move FETIDP code to its own classes 19 20 MATIS related operations contained in BDDC code 21 - Provide general case for subassembling 22 23 */ 24 25 #include <../src/ksp/pc/impls/bddc/bddc.h> /*I "petscpc.h" I*/ /* includes for fortran wrappers */ 26 #include <../src/ksp/pc/impls/bddc/bddcprivate.h> 27 #include <petscblaslapack.h> 28 29 /* temporarily declare it */ 30 PetscErrorCode PCApply_BDDC(PC,Vec,Vec); 31 32 /* -------------------------------------------------------------------------- */ 33 #undef __FUNCT__ 34 #define __FUNCT__ "PCSetFromOptions_BDDC" 35 PetscErrorCode PCSetFromOptions_BDDC(PetscOptions *PetscOptionsObject,PC pc) 36 { 37 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 38 PetscErrorCode ierr; 39 40 PetscFunctionBegin; 41 ierr = PetscOptionsHead(PetscOptionsObject,"BDDC options");CHKERRQ(ierr); 42 /* Verbose debugging */ 43 ierr = PetscOptionsInt("-pc_bddc_check_level","Verbose output for PCBDDC (intended for debug)","none",pcbddc->dbg_flag,&pcbddc->dbg_flag,NULL);CHKERRQ(ierr); 44 /* Primal space cumstomization */ 45 ierr = PetscOptionsBool("-pc_bddc_use_local_mat_graph","Use or not adjacency graph of local mat for interface analysis","none",pcbddc->use_local_adj,&pcbddc->use_local_adj,NULL);CHKERRQ(ierr); 46 ierr = PetscOptionsBool("-pc_bddc_use_vertices","Use or not corner dofs in coarse space","none",pcbddc->use_vertices,&pcbddc->use_vertices,NULL);CHKERRQ(ierr); 47 ierr = PetscOptionsBool("-pc_bddc_use_edges","Use or not edge constraints in coarse space","none",pcbddc->use_edges,&pcbddc->use_edges,NULL);CHKERRQ(ierr); 48 ierr = PetscOptionsBool("-pc_bddc_use_faces","Use or not face constraints in coarse space","none",pcbddc->use_faces,&pcbddc->use_faces,NULL);CHKERRQ(ierr); 49 ierr = PetscOptionsBool("-pc_bddc_use_true_nnsp","Use near null space attached to the matrix without modifications","none",pcbddc->use_nnsp_true,&pcbddc->use_nnsp_true,NULL);CHKERRQ(ierr); 50 ierr = PetscOptionsBool("-pc_bddc_use_qr_single","Use QR factorization for single constraints on cc (QR is used when multiple constraints are present)","none",pcbddc->use_qr_single,&pcbddc->use_qr_single,NULL);CHKERRQ(ierr); 51 /* Change of basis */ 52 ierr = PetscOptionsBool("-pc_bddc_use_change_of_basis","Use or not internal change of basis on local edge nodes","none",pcbddc->use_change_of_basis,&pcbddc->use_change_of_basis,NULL);CHKERRQ(ierr); 53 ierr = PetscOptionsBool("-pc_bddc_use_change_on_faces","Use or not internal change of basis on local face nodes","none",pcbddc->use_change_on_faces,&pcbddc->use_change_on_faces,NULL);CHKERRQ(ierr); 54 if (!pcbddc->use_change_of_basis) { 55 pcbddc->use_change_on_faces = PETSC_FALSE; 56 } 57 /* Switch between M_2 (default) and M_3 preconditioners (as defined by C. Dohrmann in the ref. article) */ 58 ierr = PetscOptionsBool("-pc_bddc_switch_static","Switch on static condensation ops around the interface preconditioner","none",pcbddc->switch_static,&pcbddc->switch_static,NULL);CHKERRQ(ierr); 59 ierr = PetscOptionsInt("-pc_bddc_coarse_redistribute","Number of procs where to redistribute coarse problem","none",pcbddc->redistribute_coarse,&pcbddc->redistribute_coarse,NULL);CHKERRQ(ierr); 60 ierr = PetscOptionsInt("-pc_bddc_coarsening_ratio","Set coarsening ratio used in multilevel coarsening","none",pcbddc->coarsening_ratio,&pcbddc->coarsening_ratio,NULL);CHKERRQ(ierr); 61 ierr = PetscOptionsInt("-pc_bddc_levels","Set maximum number of levels for multilevel","none",pcbddc->max_levels,&pcbddc->max_levels,NULL);CHKERRQ(ierr); 62 ierr = PetscOptionsBool("-pc_bddc_use_coarse_estimates","Use estimated eigenvalues for coarse problem","none",pcbddc->use_coarse_estimates,&pcbddc->use_coarse_estimates,NULL);CHKERRQ(ierr); 63 ierr = PetscOptionsBool("-pc_bddc_use_deluxe_scaling","Use deluxe scaling for BDDC","none",pcbddc->use_deluxe_scaling,&pcbddc->use_deluxe_scaling,NULL);CHKERRQ(ierr); 64 ierr = PetscOptionsBool("-pc_bddc_schur_rebuild","Whether or not the interface graph for Schur principal minors has to be rebuilt (i.e. define the interface without any adjacency)","none",pcbddc->sub_schurs_rebuild,&pcbddc->sub_schurs_rebuild,NULL);CHKERRQ(ierr); 65 ierr = PetscOptionsInt("-pc_bddc_schur_layers","Number of dofs' layers for the computation of principal minors (i.e. -1 uses all dofs)","none",pcbddc->sub_schurs_layers,&pcbddc->sub_schurs_layers,NULL);CHKERRQ(ierr); 66 ierr = PetscOptionsBool("-pc_bddc_schur_use_useradj","Whether or not the CSR graph specified by the user should be used for computing successive layers (default is to use adj of local mat)","none",pcbddc->sub_schurs_use_useradj,&pcbddc->sub_schurs_use_useradj,NULL);CHKERRQ(ierr); 67 ierr = PetscOptionsBool("-pc_bddc_schur_exact","Whether or not to use the exact Schur complement instead of the reduced one (which excludes size 1 cc)","none",pcbddc->sub_schurs_exact_schur,&pcbddc->sub_schurs_exact_schur,NULL);CHKERRQ(ierr); 68 ierr = PetscOptionsBool("-pc_bddc_deluxe_zerorows","Zero rows and columns of deluxe operators associated with primal dofs","none",pcbddc->deluxe_zerorows,&pcbddc->deluxe_zerorows,NULL);CHKERRQ(ierr); 69 ierr = PetscOptionsBool("-pc_bddc_deluxe_faster","Faster application of deluxe scaling (requires extra work during setup)","none",pcbddc->faster_deluxe,&pcbddc->faster_deluxe,NULL);CHKERRQ(ierr); 70 ierr = PetscOptionsBool("-pc_bddc_adaptive_userdefined","Use user-defined constraints (should be attached via MatSetNearNullSpace to pmat) in addition to those adaptively generated","none",pcbddc->adaptive_userdefined,&pcbddc->adaptive_userdefined,NULL);CHKERRQ(ierr); 71 ierr = PetscOptionsReal("-pc_bddc_adaptive_threshold","Threshold to be used for adaptive selection of constraints","none",pcbddc->adaptive_threshold,&pcbddc->adaptive_threshold,NULL);CHKERRQ(ierr); 72 ierr = PetscOptionsInt("-pc_bddc_adaptive_nmin","Minimum number of constraints per connected components","none",pcbddc->adaptive_nmin,&pcbddc->adaptive_nmin,NULL);CHKERRQ(ierr); 73 ierr = PetscOptionsInt("-pc_bddc_adaptive_nmax","Maximum number of constraints per connected components","none",pcbddc->adaptive_nmax,&pcbddc->adaptive_nmax,NULL);CHKERRQ(ierr); 74 ierr = PetscOptionsBool("-pc_bddc_symmetric","Symmetric computation of primal basis functions","none",pcbddc->symmetric_primal,&pcbddc->symmetric_primal,NULL);CHKERRQ(ierr); 75 ierr = PetscOptionsInt("-pc_bddc_coarse_adj","Number of processors where to map the coarse adjacency list","none",pcbddc->coarse_adj_red,&pcbddc->coarse_adj_red,NULL);CHKERRQ(ierr); 76 ierr = PetscOptionsBool("-pc_bddc_benign_trick","Apply the benign subspace trick to a class of saddle point problems","none",pcbddc->benign_saddle_point,&pcbddc->benign_saddle_point,NULL);CHKERRQ(ierr); 77 ierr = PetscOptionsBool("-pc_bddc_detect_disconnected","Detects disconnected subdomains","none",pcbddc->detect_disconnected,&pcbddc->detect_disconnected,NULL);CHKERRQ(ierr); 78 ierr = PetscOptionsTail();CHKERRQ(ierr); 79 PetscFunctionReturn(0); 80 } 81 /* -------------------------------------------------------------------------- */ 82 #undef __FUNCT__ 83 #define __FUNCT__ "PCBDDCSetChangeOfBasisMat_BDDC" 84 static PetscErrorCode PCBDDCSetChangeOfBasisMat_BDDC(PC pc, Mat change) 85 { 86 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 87 PetscErrorCode ierr; 88 89 PetscFunctionBegin; 90 ierr = PetscObjectReference((PetscObject)change);CHKERRQ(ierr); 91 ierr = MatDestroy(&pcbddc->user_ChangeOfBasisMatrix);CHKERRQ(ierr); 92 pcbddc->user_ChangeOfBasisMatrix = change; 93 PetscFunctionReturn(0); 94 } 95 #undef __FUNCT__ 96 #define __FUNCT__ "PCBDDCSetChangeOfBasisMat" 97 /*@ 98 PCBDDCSetChangeOfBasisMat - Set user defined change of basis for dofs 99 100 Collective on PC 101 102 Input Parameters: 103 + pc - the preconditioning context 104 - change - the change of basis matrix 105 106 Level: intermediate 107 108 Notes: 109 110 .seealso: PCBDDC 111 @*/ 112 PetscErrorCode PCBDDCSetChangeOfBasisMat(PC pc, Mat change) 113 { 114 PetscErrorCode ierr; 115 116 PetscFunctionBegin; 117 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 118 PetscValidHeaderSpecific(change,MAT_CLASSID,2); 119 PetscCheckSameComm(pc,1,change,2); 120 if (pc->mat) { 121 PetscInt rows_c,cols_c,rows,cols; 122 ierr = MatGetSize(pc->mat,&rows,&cols);CHKERRQ(ierr); 123 ierr = MatGetSize(change,&rows_c,&cols_c);CHKERRQ(ierr); 124 if (rows_c != rows) { 125 SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Invalid number of rows for change of basis matrix! %d != %d",rows_c,rows); 126 } 127 if (cols_c != cols) { 128 SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Invalid number of columns for change of basis matrix! %d != %d",cols_c,cols); 129 } 130 ierr = MatGetLocalSize(pc->mat,&rows,&cols);CHKERRQ(ierr); 131 ierr = MatGetLocalSize(change,&rows_c,&cols_c);CHKERRQ(ierr); 132 if (rows_c != rows) { 133 SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Invalid number of local rows for change of basis matrix! %d != %d",rows_c,rows); 134 } 135 if (cols_c != cols) { 136 SETERRQ2(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Invalid number of local columns for change of basis matrix! %d != %d",cols_c,cols); 137 } 138 } 139 ierr = PetscTryMethod(pc,"PCBDDCSetChangeOfBasisMat_C",(PC,Mat),(pc,change));CHKERRQ(ierr); 140 PetscFunctionReturn(0); 141 } 142 /* -------------------------------------------------------------------------- */ 143 #undef __FUNCT__ 144 #define __FUNCT__ "PCBDDCSetPrimalVerticesIS_BDDC" 145 static PetscErrorCode PCBDDCSetPrimalVerticesIS_BDDC(PC pc, IS PrimalVertices) 146 { 147 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 148 PetscBool isequal = PETSC_FALSE; 149 PetscErrorCode ierr; 150 151 PetscFunctionBegin; 152 ierr = PetscObjectReference((PetscObject)PrimalVertices);CHKERRQ(ierr); 153 if (pcbddc->user_primal_vertices) { 154 ierr = ISEqual(PrimalVertices,pcbddc->user_primal_vertices,&isequal);CHKERRQ(ierr); 155 } 156 ierr = ISDestroy(&pcbddc->user_primal_vertices);CHKERRQ(ierr); 157 ierr = ISDestroy(&pcbddc->user_primal_vertices_local);CHKERRQ(ierr); 158 pcbddc->user_primal_vertices = PrimalVertices; 159 if (!isequal) pcbddc->recompute_topography = PETSC_TRUE; 160 PetscFunctionReturn(0); 161 } 162 #undef __FUNCT__ 163 #define __FUNCT__ "PCBDDCSetPrimalVerticesIS" 164 /*@ 165 PCBDDCSetPrimalVerticesIS - Set additional user defined primal vertices in PCBDDC 166 167 Collective 168 169 Input Parameters: 170 + pc - the preconditioning context 171 - PrimalVertices - index set of primal vertices in global numbering (can be empty) 172 173 Level: intermediate 174 175 Notes: 176 Any process can list any global node 177 178 .seealso: PCBDDC 179 @*/ 180 PetscErrorCode PCBDDCSetPrimalVerticesIS(PC pc, IS PrimalVertices) 181 { 182 PetscErrorCode ierr; 183 184 PetscFunctionBegin; 185 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 186 PetscValidHeaderSpecific(PrimalVertices,IS_CLASSID,2); 187 PetscCheckSameComm(pc,1,PrimalVertices,2); 188 ierr = PetscTryMethod(pc,"PCBDDCSetPrimalVerticesIS_C",(PC,IS),(pc,PrimalVertices));CHKERRQ(ierr); 189 PetscFunctionReturn(0); 190 } 191 /* -------------------------------------------------------------------------- */ 192 #undef __FUNCT__ 193 #define __FUNCT__ "PCBDDCSetPrimalVerticesLocalIS_BDDC" 194 static PetscErrorCode PCBDDCSetPrimalVerticesLocalIS_BDDC(PC pc, IS PrimalVertices) 195 { 196 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 197 PetscBool isequal = PETSC_FALSE; 198 PetscErrorCode ierr; 199 200 PetscFunctionBegin; 201 ierr = PetscObjectReference((PetscObject)PrimalVertices);CHKERRQ(ierr); 202 if (pcbddc->user_primal_vertices_local) { 203 ierr = ISEqual(PrimalVertices,pcbddc->user_primal_vertices_local,&isequal);CHKERRQ(ierr); 204 } 205 ierr = ISDestroy(&pcbddc->user_primal_vertices);CHKERRQ(ierr); 206 ierr = ISDestroy(&pcbddc->user_primal_vertices_local);CHKERRQ(ierr); 207 pcbddc->user_primal_vertices_local = PrimalVertices; 208 if (!isequal) pcbddc->recompute_topography = PETSC_TRUE; 209 PetscFunctionReturn(0); 210 } 211 #undef __FUNCT__ 212 #define __FUNCT__ "PCBDDCSetPrimalVerticesLocalIS" 213 /*@ 214 PCBDDCSetPrimalVerticesLocalIS - Set additional user defined primal vertices in PCBDDC 215 216 Collective 217 218 Input Parameters: 219 + pc - the preconditioning context 220 - PrimalVertices - index set of primal vertices in local numbering (can be empty) 221 222 Level: intermediate 223 224 Notes: 225 226 .seealso: PCBDDC 227 @*/ 228 PetscErrorCode PCBDDCSetPrimalVerticesLocalIS(PC pc, IS PrimalVertices) 229 { 230 PetscErrorCode ierr; 231 232 PetscFunctionBegin; 233 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 234 PetscValidHeaderSpecific(PrimalVertices,IS_CLASSID,2); 235 PetscCheckSameComm(pc,1,PrimalVertices,2); 236 ierr = PetscTryMethod(pc,"PCBDDCSetPrimalVerticesLocalIS_C",(PC,IS),(pc,PrimalVertices));CHKERRQ(ierr); 237 PetscFunctionReturn(0); 238 } 239 /* -------------------------------------------------------------------------- */ 240 #undef __FUNCT__ 241 #define __FUNCT__ "PCBDDCSetCoarseningRatio_BDDC" 242 static PetscErrorCode PCBDDCSetCoarseningRatio_BDDC(PC pc,PetscInt k) 243 { 244 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 245 246 PetscFunctionBegin; 247 pcbddc->coarsening_ratio = k; 248 PetscFunctionReturn(0); 249 } 250 251 #undef __FUNCT__ 252 #define __FUNCT__ "PCBDDCSetCoarseningRatio" 253 /*@ 254 PCBDDCSetCoarseningRatio - Set coarsening ratio used in multilevel 255 256 Logically collective on PC 257 258 Input Parameters: 259 + pc - the preconditioning context 260 - k - coarsening ratio (H/h at the coarser level) 261 262 Options Database Keys: 263 . -pc_bddc_coarsening_ratio 264 265 Level: intermediate 266 267 Notes: 268 Approximatively k subdomains at the finer level will be aggregated into a single subdomain at the coarser level 269 270 .seealso: PCBDDC, PCBDDCSetLevels() 271 @*/ 272 PetscErrorCode PCBDDCSetCoarseningRatio(PC pc,PetscInt k) 273 { 274 PetscErrorCode ierr; 275 276 PetscFunctionBegin; 277 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 278 PetscValidLogicalCollectiveInt(pc,k,2); 279 ierr = PetscTryMethod(pc,"PCBDDCSetCoarseningRatio_C",(PC,PetscInt),(pc,k));CHKERRQ(ierr); 280 PetscFunctionReturn(0); 281 } 282 283 /* The following functions (PCBDDCSetUseExactDirichlet PCBDDCSetLevel) are not public */ 284 #undef __FUNCT__ 285 #define __FUNCT__ "PCBDDCSetUseExactDirichlet_BDDC" 286 static PetscErrorCode PCBDDCSetUseExactDirichlet_BDDC(PC pc,PetscBool flg) 287 { 288 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 289 290 PetscFunctionBegin; 291 pcbddc->use_exact_dirichlet_trick = flg; 292 PetscFunctionReturn(0); 293 } 294 295 #undef __FUNCT__ 296 #define __FUNCT__ "PCBDDCSetUseExactDirichlet" 297 PetscErrorCode PCBDDCSetUseExactDirichlet(PC pc,PetscBool flg) 298 { 299 PetscErrorCode ierr; 300 301 PetscFunctionBegin; 302 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 303 PetscValidLogicalCollectiveBool(pc,flg,2); 304 ierr = PetscTryMethod(pc,"PCBDDCSetUseExactDirichlet_C",(PC,PetscBool),(pc,flg));CHKERRQ(ierr); 305 PetscFunctionReturn(0); 306 } 307 308 #undef __FUNCT__ 309 #define __FUNCT__ "PCBDDCSetLevel_BDDC" 310 static PetscErrorCode PCBDDCSetLevel_BDDC(PC pc,PetscInt level) 311 { 312 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 313 314 PetscFunctionBegin; 315 pcbddc->current_level = level; 316 PetscFunctionReturn(0); 317 } 318 319 #undef __FUNCT__ 320 #define __FUNCT__ "PCBDDCSetLevel" 321 PetscErrorCode PCBDDCSetLevel(PC pc,PetscInt level) 322 { 323 PetscErrorCode ierr; 324 325 PetscFunctionBegin; 326 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 327 PetscValidLogicalCollectiveInt(pc,level,2); 328 ierr = PetscTryMethod(pc,"PCBDDCSetLevel_C",(PC,PetscInt),(pc,level));CHKERRQ(ierr); 329 PetscFunctionReturn(0); 330 } 331 332 #undef __FUNCT__ 333 #define __FUNCT__ "PCBDDCSetLevels_BDDC" 334 static PetscErrorCode PCBDDCSetLevels_BDDC(PC pc,PetscInt levels) 335 { 336 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 337 338 PetscFunctionBegin; 339 pcbddc->max_levels = levels; 340 PetscFunctionReturn(0); 341 } 342 343 #undef __FUNCT__ 344 #define __FUNCT__ "PCBDDCSetLevels" 345 /*@ 346 PCBDDCSetLevels - Sets the maximum number of levels for multilevel 347 348 Logically collective on PC 349 350 Input Parameters: 351 + pc - the preconditioning context 352 - levels - the maximum number of levels (max 9) 353 354 Options Database Keys: 355 . -pc_bddc_levels 356 357 Level: intermediate 358 359 Notes: 360 Default value is 0, i.e. traditional one-level BDDC 361 362 .seealso: PCBDDC, PCBDDCSetCoarseningRatio() 363 @*/ 364 PetscErrorCode PCBDDCSetLevels(PC pc,PetscInt levels) 365 { 366 PetscErrorCode ierr; 367 368 PetscFunctionBegin; 369 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 370 PetscValidLogicalCollectiveInt(pc,levels,2); 371 if (levels > 99) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Maximum number of levels for bddc is 99\n"); 372 ierr = PetscTryMethod(pc,"PCBDDCSetLevels_C",(PC,PetscInt),(pc,levels));CHKERRQ(ierr); 373 PetscFunctionReturn(0); 374 } 375 /* -------------------------------------------------------------------------- */ 376 377 #undef __FUNCT__ 378 #define __FUNCT__ "PCBDDCSetNullSpace_BDDC" 379 static PetscErrorCode PCBDDCSetNullSpace_BDDC(PC pc,MatNullSpace NullSpace) 380 { 381 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 382 PetscErrorCode ierr; 383 384 PetscFunctionBegin; 385 ierr = PetscObjectReference((PetscObject)NullSpace);CHKERRQ(ierr); 386 ierr = MatNullSpaceDestroy(&pcbddc->NullSpace);CHKERRQ(ierr); 387 pcbddc->NullSpace = NullSpace; 388 PetscFunctionReturn(0); 389 } 390 391 #undef __FUNCT__ 392 #define __FUNCT__ "PCBDDCSetNullSpace" 393 /*@ 394 PCBDDCSetNullSpace - Set nullspace for BDDC operator 395 396 Logically collective on PC and MatNullSpace 397 398 Input Parameters: 399 + pc - the preconditioning context 400 - NullSpace - Null space of the linear operator to be preconditioned (Pmat) 401 402 Level: intermediate 403 404 Notes: 405 406 .seealso: PCBDDC 407 @*/ 408 PetscErrorCode PCBDDCSetNullSpace(PC pc,MatNullSpace NullSpace) 409 { 410 PetscErrorCode ierr; 411 412 PetscFunctionBegin; 413 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 414 PetscValidHeaderSpecific(NullSpace,MAT_NULLSPACE_CLASSID,2); 415 PetscCheckSameComm(pc,1,NullSpace,2); 416 ierr = PetscTryMethod(pc,"PCBDDCSetNullSpace_C",(PC,MatNullSpace),(pc,NullSpace));CHKERRQ(ierr); 417 PetscFunctionReturn(0); 418 } 419 /* -------------------------------------------------------------------------- */ 420 421 #undef __FUNCT__ 422 #define __FUNCT__ "PCBDDCSetDirichletBoundaries_BDDC" 423 static PetscErrorCode PCBDDCSetDirichletBoundaries_BDDC(PC pc,IS DirichletBoundaries) 424 { 425 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 426 PetscBool isequal = PETSC_FALSE; 427 PetscErrorCode ierr; 428 429 PetscFunctionBegin; 430 ierr = PetscObjectReference((PetscObject)DirichletBoundaries);CHKERRQ(ierr); 431 if (pcbddc->DirichletBoundaries) { 432 ierr = ISEqual(DirichletBoundaries,pcbddc->DirichletBoundaries,&isequal);CHKERRQ(ierr); 433 } 434 /* last user setting takes precendence -> destroy any other customization */ 435 ierr = ISDestroy(&pcbddc->DirichletBoundariesLocal);CHKERRQ(ierr); 436 ierr = ISDestroy(&pcbddc->DirichletBoundaries);CHKERRQ(ierr); 437 pcbddc->DirichletBoundaries = DirichletBoundaries; 438 if (!isequal) pcbddc->recompute_topography = PETSC_TRUE; 439 PetscFunctionReturn(0); 440 } 441 442 #undef __FUNCT__ 443 #define __FUNCT__ "PCBDDCSetDirichletBoundaries" 444 /*@ 445 PCBDDCSetDirichletBoundaries - Set IS defining Dirichlet boundaries for the global problem. 446 447 Collective 448 449 Input Parameters: 450 + pc - the preconditioning context 451 - DirichletBoundaries - parallel IS defining the Dirichlet boundaries 452 453 Level: intermediate 454 455 Notes: 456 Provide the information if you used MatZeroRows/Columns routines. Any process can list any global node 457 458 .seealso: PCBDDC, PCBDDCSetDirichletBoundariesLocal(), MatZeroRows(), MatZeroRowsColumns() 459 @*/ 460 PetscErrorCode PCBDDCSetDirichletBoundaries(PC pc,IS DirichletBoundaries) 461 { 462 PetscErrorCode ierr; 463 464 PetscFunctionBegin; 465 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 466 PetscValidHeaderSpecific(DirichletBoundaries,IS_CLASSID,2); 467 PetscCheckSameComm(pc,1,DirichletBoundaries,2); 468 ierr = PetscTryMethod(pc,"PCBDDCSetDirichletBoundaries_C",(PC,IS),(pc,DirichletBoundaries));CHKERRQ(ierr); 469 PetscFunctionReturn(0); 470 } 471 /* -------------------------------------------------------------------------- */ 472 473 #undef __FUNCT__ 474 #define __FUNCT__ "PCBDDCSetDirichletBoundariesLocal_BDDC" 475 static PetscErrorCode PCBDDCSetDirichletBoundariesLocal_BDDC(PC pc,IS DirichletBoundaries) 476 { 477 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 478 PetscBool isequal = PETSC_FALSE; 479 PetscErrorCode ierr; 480 481 PetscFunctionBegin; 482 ierr = PetscObjectReference((PetscObject)DirichletBoundaries);CHKERRQ(ierr); 483 if (pcbddc->DirichletBoundariesLocal) { 484 ierr = ISEqual(DirichletBoundaries,pcbddc->DirichletBoundariesLocal,&isequal);CHKERRQ(ierr); 485 } 486 /* last user setting takes precendence -> destroy any other customization */ 487 ierr = ISDestroy(&pcbddc->DirichletBoundariesLocal);CHKERRQ(ierr); 488 ierr = ISDestroy(&pcbddc->DirichletBoundaries);CHKERRQ(ierr); 489 pcbddc->DirichletBoundariesLocal = DirichletBoundaries; 490 if (!isequal) pcbddc->recompute_topography = PETSC_TRUE; 491 PetscFunctionReturn(0); 492 } 493 494 #undef __FUNCT__ 495 #define __FUNCT__ "PCBDDCSetDirichletBoundariesLocal" 496 /*@ 497 PCBDDCSetDirichletBoundariesLocal - Set IS defining Dirichlet boundaries for the global problem in local ordering. 498 499 Collective 500 501 Input Parameters: 502 + pc - the preconditioning context 503 - DirichletBoundaries - parallel IS defining the Dirichlet boundaries (in local ordering) 504 505 Level: intermediate 506 507 Notes: 508 509 .seealso: PCBDDC, PCBDDCSetDirichletBoundaries(), MatZeroRows(), MatZeroRowsColumns() 510 @*/ 511 PetscErrorCode PCBDDCSetDirichletBoundariesLocal(PC pc,IS DirichletBoundaries) 512 { 513 PetscErrorCode ierr; 514 515 PetscFunctionBegin; 516 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 517 PetscValidHeaderSpecific(DirichletBoundaries,IS_CLASSID,2); 518 PetscCheckSameComm(pc,1,DirichletBoundaries,2); 519 ierr = PetscTryMethod(pc,"PCBDDCSetDirichletBoundariesLocal_C",(PC,IS),(pc,DirichletBoundaries));CHKERRQ(ierr); 520 PetscFunctionReturn(0); 521 } 522 /* -------------------------------------------------------------------------- */ 523 524 #undef __FUNCT__ 525 #define __FUNCT__ "PCBDDCSetNeumannBoundaries_BDDC" 526 static PetscErrorCode PCBDDCSetNeumannBoundaries_BDDC(PC pc,IS NeumannBoundaries) 527 { 528 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 529 PetscBool isequal = PETSC_FALSE; 530 PetscErrorCode ierr; 531 532 PetscFunctionBegin; 533 ierr = PetscObjectReference((PetscObject)NeumannBoundaries);CHKERRQ(ierr); 534 if (pcbddc->NeumannBoundaries) { 535 ierr = ISEqual(NeumannBoundaries,pcbddc->NeumannBoundaries,&isequal);CHKERRQ(ierr); 536 } 537 /* last user setting takes precendence -> destroy any other customization */ 538 ierr = ISDestroy(&pcbddc->NeumannBoundariesLocal);CHKERRQ(ierr); 539 ierr = ISDestroy(&pcbddc->NeumannBoundaries);CHKERRQ(ierr); 540 pcbddc->NeumannBoundaries = NeumannBoundaries; 541 if (!isequal) pcbddc->recompute_topography = PETSC_TRUE; 542 PetscFunctionReturn(0); 543 } 544 545 #undef __FUNCT__ 546 #define __FUNCT__ "PCBDDCSetNeumannBoundaries" 547 /*@ 548 PCBDDCSetNeumannBoundaries - Set IS defining Neumann boundaries for the global problem. 549 550 Collective 551 552 Input Parameters: 553 + pc - the preconditioning context 554 - NeumannBoundaries - parallel IS defining the Neumann boundaries 555 556 Level: intermediate 557 558 Notes: 559 Any process can list any global node 560 561 .seealso: PCBDDC, PCBDDCSetNeumannBoundariesLocal() 562 @*/ 563 PetscErrorCode PCBDDCSetNeumannBoundaries(PC pc,IS NeumannBoundaries) 564 { 565 PetscErrorCode ierr; 566 567 PetscFunctionBegin; 568 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 569 PetscValidHeaderSpecific(NeumannBoundaries,IS_CLASSID,2); 570 PetscCheckSameComm(pc,1,NeumannBoundaries,2); 571 ierr = PetscTryMethod(pc,"PCBDDCSetNeumannBoundaries_C",(PC,IS),(pc,NeumannBoundaries));CHKERRQ(ierr); 572 PetscFunctionReturn(0); 573 } 574 /* -------------------------------------------------------------------------- */ 575 576 #undef __FUNCT__ 577 #define __FUNCT__ "PCBDDCSetNeumannBoundariesLocal_BDDC" 578 static PetscErrorCode PCBDDCSetNeumannBoundariesLocal_BDDC(PC pc,IS NeumannBoundaries) 579 { 580 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 581 PetscBool isequal = PETSC_FALSE; 582 PetscErrorCode ierr; 583 584 PetscFunctionBegin; 585 ierr = PetscObjectReference((PetscObject)NeumannBoundaries);CHKERRQ(ierr); 586 if (pcbddc->NeumannBoundariesLocal) { 587 ierr = ISEqual(NeumannBoundaries,pcbddc->NeumannBoundariesLocal,&isequal);CHKERRQ(ierr); 588 } 589 /* last user setting takes precendence -> destroy any other customization */ 590 ierr = ISDestroy(&pcbddc->NeumannBoundariesLocal);CHKERRQ(ierr); 591 ierr = ISDestroy(&pcbddc->NeumannBoundaries);CHKERRQ(ierr); 592 pcbddc->NeumannBoundariesLocal = NeumannBoundaries; 593 if (!isequal) pcbddc->recompute_topography = PETSC_TRUE; 594 PetscFunctionReturn(0); 595 } 596 597 #undef __FUNCT__ 598 #define __FUNCT__ "PCBDDCSetNeumannBoundariesLocal" 599 /*@ 600 PCBDDCSetNeumannBoundariesLocal - Set IS defining Neumann boundaries for the global problem in local ordering. 601 602 Collective 603 604 Input Parameters: 605 + pc - the preconditioning context 606 - NeumannBoundaries - parallel IS defining the subdomain part of Neumann boundaries (in local ordering) 607 608 Level: intermediate 609 610 Notes: 611 612 .seealso: PCBDDC, PCBDDCSetNeumannBoundaries() 613 @*/ 614 PetscErrorCode PCBDDCSetNeumannBoundariesLocal(PC pc,IS NeumannBoundaries) 615 { 616 PetscErrorCode ierr; 617 618 PetscFunctionBegin; 619 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 620 PetscValidHeaderSpecific(NeumannBoundaries,IS_CLASSID,2); 621 PetscCheckSameComm(pc,1,NeumannBoundaries,2); 622 ierr = PetscTryMethod(pc,"PCBDDCSetNeumannBoundariesLocal_C",(PC,IS),(pc,NeumannBoundaries));CHKERRQ(ierr); 623 PetscFunctionReturn(0); 624 } 625 /* -------------------------------------------------------------------------- */ 626 627 #undef __FUNCT__ 628 #define __FUNCT__ "PCBDDCGetDirichletBoundaries_BDDC" 629 static PetscErrorCode PCBDDCGetDirichletBoundaries_BDDC(PC pc,IS *DirichletBoundaries) 630 { 631 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 632 633 PetscFunctionBegin; 634 *DirichletBoundaries = pcbddc->DirichletBoundaries; 635 PetscFunctionReturn(0); 636 } 637 638 #undef __FUNCT__ 639 #define __FUNCT__ "PCBDDCGetDirichletBoundaries" 640 /*@ 641 PCBDDCGetDirichletBoundaries - Get parallel IS for Dirichlet boundaries 642 643 Collective 644 645 Input Parameters: 646 . pc - the preconditioning context 647 648 Output Parameters: 649 . DirichletBoundaries - index set defining the Dirichlet boundaries 650 651 Level: intermediate 652 653 Notes: 654 The IS returned (if any) is the same passed in earlier by the user with PCBDDCSetDirichletBoundaries 655 656 .seealso: PCBDDC 657 @*/ 658 PetscErrorCode PCBDDCGetDirichletBoundaries(PC pc,IS *DirichletBoundaries) 659 { 660 PetscErrorCode ierr; 661 662 PetscFunctionBegin; 663 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 664 ierr = PetscUseMethod(pc,"PCBDDCGetDirichletBoundaries_C",(PC,IS*),(pc,DirichletBoundaries));CHKERRQ(ierr); 665 PetscFunctionReturn(0); 666 } 667 /* -------------------------------------------------------------------------- */ 668 669 #undef __FUNCT__ 670 #define __FUNCT__ "PCBDDCGetDirichletBoundariesLocal_BDDC" 671 static PetscErrorCode PCBDDCGetDirichletBoundariesLocal_BDDC(PC pc,IS *DirichletBoundaries) 672 { 673 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 674 675 PetscFunctionBegin; 676 *DirichletBoundaries = pcbddc->DirichletBoundariesLocal; 677 PetscFunctionReturn(0); 678 } 679 680 #undef __FUNCT__ 681 #define __FUNCT__ "PCBDDCGetDirichletBoundariesLocal" 682 /*@ 683 PCBDDCGetDirichletBoundariesLocal - Get parallel IS for Dirichlet boundaries (in local ordering) 684 685 Collective 686 687 Input Parameters: 688 . pc - the preconditioning context 689 690 Output Parameters: 691 . DirichletBoundaries - index set defining the subdomain part of Dirichlet boundaries 692 693 Level: intermediate 694 695 Notes: 696 The IS returned could be the same passed in earlier by the user (if provided with PCBDDCSetDirichletBoundariesLocal) or a global-to-local map of the global IS (if provided with PCBDDCSetDirichletBoundaries). 697 In the latter case, the IS will be available after PCSetUp. 698 699 .seealso: PCBDDC 700 @*/ 701 PetscErrorCode PCBDDCGetDirichletBoundariesLocal(PC pc,IS *DirichletBoundaries) 702 { 703 PetscErrorCode ierr; 704 705 PetscFunctionBegin; 706 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 707 ierr = PetscUseMethod(pc,"PCBDDCGetDirichletBoundariesLocal_C",(PC,IS*),(pc,DirichletBoundaries));CHKERRQ(ierr); 708 PetscFunctionReturn(0); 709 } 710 /* -------------------------------------------------------------------------- */ 711 712 #undef __FUNCT__ 713 #define __FUNCT__ "PCBDDCGetNeumannBoundaries_BDDC" 714 static PetscErrorCode PCBDDCGetNeumannBoundaries_BDDC(PC pc,IS *NeumannBoundaries) 715 { 716 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 717 718 PetscFunctionBegin; 719 *NeumannBoundaries = pcbddc->NeumannBoundaries; 720 PetscFunctionReturn(0); 721 } 722 723 #undef __FUNCT__ 724 #define __FUNCT__ "PCBDDCGetNeumannBoundaries" 725 /*@ 726 PCBDDCGetNeumannBoundaries - Get parallel IS for Neumann boundaries 727 728 Collective 729 730 Input Parameters: 731 . pc - the preconditioning context 732 733 Output Parameters: 734 . NeumannBoundaries - index set defining the Neumann boundaries 735 736 Level: intermediate 737 738 Notes: 739 The IS returned (if any) is the same passed in earlier by the user with PCBDDCSetNeumannBoundaries 740 741 .seealso: PCBDDC 742 @*/ 743 PetscErrorCode PCBDDCGetNeumannBoundaries(PC pc,IS *NeumannBoundaries) 744 { 745 PetscErrorCode ierr; 746 747 PetscFunctionBegin; 748 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 749 ierr = PetscUseMethod(pc,"PCBDDCGetNeumannBoundaries_C",(PC,IS*),(pc,NeumannBoundaries));CHKERRQ(ierr); 750 PetscFunctionReturn(0); 751 } 752 /* -------------------------------------------------------------------------- */ 753 754 #undef __FUNCT__ 755 #define __FUNCT__ "PCBDDCGetNeumannBoundariesLocal_BDDC" 756 static PetscErrorCode PCBDDCGetNeumannBoundariesLocal_BDDC(PC pc,IS *NeumannBoundaries) 757 { 758 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 759 760 PetscFunctionBegin; 761 *NeumannBoundaries = pcbddc->NeumannBoundariesLocal; 762 PetscFunctionReturn(0); 763 } 764 765 #undef __FUNCT__ 766 #define __FUNCT__ "PCBDDCGetNeumannBoundariesLocal" 767 /*@ 768 PCBDDCGetNeumannBoundariesLocal - Get parallel IS for Neumann boundaries (in local ordering) 769 770 Collective 771 772 Input Parameters: 773 . pc - the preconditioning context 774 775 Output Parameters: 776 . NeumannBoundaries - index set defining the subdomain part of Neumann boundaries 777 778 Level: intermediate 779 780 Notes: 781 The IS returned could be the same passed in earlier by the user (if provided with PCBDDCSetNeumannBoundariesLocal) or a global-to-local map of the global IS (if provided with PCBDDCSetNeumannBoundaries). 782 In the latter case, the IS will be available after PCSetUp. 783 784 .seealso: PCBDDC 785 @*/ 786 PetscErrorCode PCBDDCGetNeumannBoundariesLocal(PC pc,IS *NeumannBoundaries) 787 { 788 PetscErrorCode ierr; 789 790 PetscFunctionBegin; 791 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 792 ierr = PetscUseMethod(pc,"PCBDDCGetNeumannBoundariesLocal_C",(PC,IS*),(pc,NeumannBoundaries));CHKERRQ(ierr); 793 PetscFunctionReturn(0); 794 } 795 /* -------------------------------------------------------------------------- */ 796 797 #undef __FUNCT__ 798 #define __FUNCT__ "PCBDDCSetLocalAdjacencyGraph_BDDC" 799 static PetscErrorCode PCBDDCSetLocalAdjacencyGraph_BDDC(PC pc, PetscInt nvtxs,const PetscInt xadj[],const PetscInt adjncy[], PetscCopyMode copymode) 800 { 801 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 802 PCBDDCGraph mat_graph = pcbddc->mat_graph; 803 PetscBool same_data = PETSC_FALSE; 804 PetscErrorCode ierr; 805 806 PetscFunctionBegin; 807 if (mat_graph->nvtxs_csr == nvtxs) { 808 if (mat_graph->xadj == xadj && mat_graph->adjncy == adjncy) same_data = PETSC_TRUE; 809 if (!same_data && mat_graph->xadj[nvtxs] == xadj[nvtxs]) { 810 ierr = PetscMemcmp(adjncy,mat_graph->adjncy,nvtxs*sizeof(PetscInt),&same_data);CHKERRQ(ierr); 811 } 812 } 813 if (!same_data) { 814 /* free old CSR */ 815 ierr = PCBDDCGraphResetCSR(mat_graph);CHKERRQ(ierr); 816 /* get CSR into graph structure */ 817 if (copymode == PETSC_COPY_VALUES) { 818 ierr = PetscMalloc1(nvtxs+1,&mat_graph->xadj);CHKERRQ(ierr); 819 ierr = PetscMalloc1(xadj[nvtxs],&mat_graph->adjncy);CHKERRQ(ierr); 820 ierr = PetscMemcpy(mat_graph->xadj,xadj,(nvtxs+1)*sizeof(PetscInt));CHKERRQ(ierr); 821 ierr = PetscMemcpy(mat_graph->adjncy,adjncy,xadj[nvtxs]*sizeof(PetscInt));CHKERRQ(ierr); 822 } else if (copymode == PETSC_OWN_POINTER) { 823 mat_graph->xadj = (PetscInt*)xadj; 824 mat_graph->adjncy = (PetscInt*)adjncy; 825 } 826 mat_graph->nvtxs_csr = nvtxs; 827 pcbddc->recompute_topography = PETSC_TRUE; 828 } 829 PetscFunctionReturn(0); 830 } 831 832 #undef __FUNCT__ 833 #define __FUNCT__ "PCBDDCSetLocalAdjacencyGraph" 834 /*@ 835 PCBDDCSetLocalAdjacencyGraph - Set adjacency structure (CSR graph) of the local matrix 836 837 Not collective 838 839 Input Parameters: 840 + pc - the preconditioning context 841 . nvtxs - number of local vertices of the graph (i.e., the size of the local problem) 842 . xadj, adjncy - the CSR graph 843 - copymode - either PETSC_COPY_VALUES or PETSC_OWN_POINTER. 844 845 Level: intermediate 846 847 Notes: 848 849 .seealso: PCBDDC,PetscCopyMode 850 @*/ 851 PetscErrorCode PCBDDCSetLocalAdjacencyGraph(PC pc,PetscInt nvtxs,const PetscInt xadj[],const PetscInt adjncy[], PetscCopyMode copymode) 852 { 853 void (*f)(void) = 0; 854 PetscErrorCode ierr; 855 856 PetscFunctionBegin; 857 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 858 PetscValidIntPointer(xadj,3); 859 PetscValidIntPointer(adjncy,4); 860 if (copymode != PETSC_COPY_VALUES && copymode != PETSC_OWN_POINTER) { 861 SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_SUP,"Unsupported copy mode %d in %s\n",copymode,__FUNCT__); 862 } 863 ierr = PetscTryMethod(pc,"PCBDDCSetLocalAdjacencyGraph_C",(PC,PetscInt,const PetscInt[],const PetscInt[],PetscCopyMode),(pc,nvtxs,xadj,adjncy,copymode));CHKERRQ(ierr); 864 /* free arrays if PCBDDC is not the PC type */ 865 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCBDDCSetLocalAdjacencyGraph_C",&f);CHKERRQ(ierr); 866 if (!f && copymode == PETSC_OWN_POINTER) { 867 ierr = PetscFree(xadj);CHKERRQ(ierr); 868 ierr = PetscFree(adjncy);CHKERRQ(ierr); 869 } 870 PetscFunctionReturn(0); 871 } 872 /* -------------------------------------------------------------------------- */ 873 874 #undef __FUNCT__ 875 #define __FUNCT__ "PCBDDCSetDofsSplittingLocal_BDDC" 876 static PetscErrorCode PCBDDCSetDofsSplittingLocal_BDDC(PC pc,PetscInt n_is, IS ISForDofs[]) 877 { 878 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 879 PetscInt i; 880 PetscBool isequal = PETSC_FALSE; 881 PetscErrorCode ierr; 882 883 PetscFunctionBegin; 884 if (pcbddc->n_ISForDofsLocal == n_is) { 885 for (i=0;i<n_is;i++) { 886 PetscBool isequalt; 887 ierr = ISEqual(ISForDofs[i],pcbddc->ISForDofsLocal[i],&isequalt);CHKERRQ(ierr); 888 if (!isequalt) break; 889 } 890 if (i == n_is) isequal = PETSC_TRUE; 891 } 892 for (i=0;i<n_is;i++) { 893 ierr = PetscObjectReference((PetscObject)ISForDofs[i]);CHKERRQ(ierr); 894 } 895 /* Destroy ISes if they were already set */ 896 for (i=0;i<pcbddc->n_ISForDofsLocal;i++) { 897 ierr = ISDestroy(&pcbddc->ISForDofsLocal[i]);CHKERRQ(ierr); 898 } 899 ierr = PetscFree(pcbddc->ISForDofsLocal);CHKERRQ(ierr); 900 /* last user setting takes precendence -> destroy any other customization */ 901 for (i=0;i<pcbddc->n_ISForDofs;i++) { 902 ierr = ISDestroy(&pcbddc->ISForDofs[i]);CHKERRQ(ierr); 903 } 904 ierr = PetscFree(pcbddc->ISForDofs);CHKERRQ(ierr); 905 pcbddc->n_ISForDofs = 0; 906 /* allocate space then set */ 907 if (n_is) { 908 ierr = PetscMalloc1(n_is,&pcbddc->ISForDofsLocal);CHKERRQ(ierr); 909 } 910 for (i=0;i<n_is;i++) { 911 pcbddc->ISForDofsLocal[i] = ISForDofs[i]; 912 } 913 pcbddc->n_ISForDofsLocal = n_is; 914 if (n_is) pcbddc->user_provided_isfordofs = PETSC_TRUE; 915 if (!isequal) pcbddc->recompute_topography = PETSC_TRUE; 916 PetscFunctionReturn(0); 917 } 918 919 #undef __FUNCT__ 920 #define __FUNCT__ "PCBDDCSetDofsSplittingLocal" 921 /*@ 922 PCBDDCSetDofsSplittingLocal - Set index sets defining fields of the local subdomain matrix 923 924 Collective 925 926 Input Parameters: 927 + pc - the preconditioning context 928 . n_is - number of index sets defining the fields 929 - ISForDofs - array of IS describing the fields in local ordering 930 931 Level: intermediate 932 933 Notes: 934 n_is should be the same among processes. Not all nodes need to be listed: unlisted nodes will belong to the complement field. 935 936 .seealso: PCBDDC 937 @*/ 938 PetscErrorCode PCBDDCSetDofsSplittingLocal(PC pc,PetscInt n_is, IS ISForDofs[]) 939 { 940 PetscInt i; 941 PetscErrorCode ierr; 942 943 PetscFunctionBegin; 944 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 945 PetscValidLogicalCollectiveInt(pc,n_is,2); 946 for (i=0;i<n_is;i++) { 947 PetscCheckSameComm(pc,1,ISForDofs[i],3); 948 PetscValidHeaderSpecific(ISForDofs[i],IS_CLASSID,3); 949 } 950 ierr = PetscTryMethod(pc,"PCBDDCSetDofsSplittingLocal_C",(PC,PetscInt,IS[]),(pc,n_is,ISForDofs));CHKERRQ(ierr); 951 PetscFunctionReturn(0); 952 } 953 /* -------------------------------------------------------------------------- */ 954 955 #undef __FUNCT__ 956 #define __FUNCT__ "PCBDDCSetDofsSplitting_BDDC" 957 static PetscErrorCode PCBDDCSetDofsSplitting_BDDC(PC pc,PetscInt n_is, IS ISForDofs[]) 958 { 959 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 960 PetscInt i; 961 PetscBool isequal = PETSC_FALSE; 962 PetscErrorCode ierr; 963 964 PetscFunctionBegin; 965 if (pcbddc->n_ISForDofs == n_is) { 966 for (i=0;i<n_is;i++) { 967 PetscBool isequalt; 968 ierr = ISEqual(ISForDofs[i],pcbddc->ISForDofs[i],&isequalt);CHKERRQ(ierr); 969 if (!isequalt) break; 970 } 971 if (i == n_is) isequal = PETSC_TRUE; 972 } 973 for (i=0;i<n_is;i++) { 974 ierr = PetscObjectReference((PetscObject)ISForDofs[i]);CHKERRQ(ierr); 975 } 976 /* Destroy ISes if they were already set */ 977 for (i=0;i<pcbddc->n_ISForDofs;i++) { 978 ierr = ISDestroy(&pcbddc->ISForDofs[i]);CHKERRQ(ierr); 979 } 980 ierr = PetscFree(pcbddc->ISForDofs);CHKERRQ(ierr); 981 /* last user setting takes precendence -> destroy any other customization */ 982 for (i=0;i<pcbddc->n_ISForDofsLocal;i++) { 983 ierr = ISDestroy(&pcbddc->ISForDofsLocal[i]);CHKERRQ(ierr); 984 } 985 ierr = PetscFree(pcbddc->ISForDofsLocal);CHKERRQ(ierr); 986 pcbddc->n_ISForDofsLocal = 0; 987 /* allocate space then set */ 988 if (n_is) { 989 ierr = PetscMalloc1(n_is,&pcbddc->ISForDofs);CHKERRQ(ierr); 990 } 991 for (i=0;i<n_is;i++) { 992 pcbddc->ISForDofs[i] = ISForDofs[i]; 993 } 994 pcbddc->n_ISForDofs = n_is; 995 if (n_is) pcbddc->user_provided_isfordofs = PETSC_TRUE; 996 if (!isequal) pcbddc->recompute_topography = PETSC_TRUE; 997 PetscFunctionReturn(0); 998 } 999 1000 #undef __FUNCT__ 1001 #define __FUNCT__ "PCBDDCSetDofsSplitting" 1002 /*@ 1003 PCBDDCSetDofsSplitting - Set index sets defining fields of the global matrix 1004 1005 Collective 1006 1007 Input Parameters: 1008 + pc - the preconditioning context 1009 . n_is - number of index sets defining the fields 1010 - ISForDofs - array of IS describing the fields in global ordering 1011 1012 Level: intermediate 1013 1014 Notes: 1015 Any process can list any global node. Not all nodes need to be listed: unlisted nodes will belong to the complement field. 1016 1017 .seealso: PCBDDC 1018 @*/ 1019 PetscErrorCode PCBDDCSetDofsSplitting(PC pc,PetscInt n_is, IS ISForDofs[]) 1020 { 1021 PetscInt i; 1022 PetscErrorCode ierr; 1023 1024 PetscFunctionBegin; 1025 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 1026 PetscValidLogicalCollectiveInt(pc,n_is,2); 1027 for (i=0;i<n_is;i++) { 1028 PetscCheckSameComm(pc,1,ISForDofs[i],3); 1029 PetscValidHeaderSpecific(ISForDofs[i],IS_CLASSID,3); 1030 } 1031 ierr = PetscTryMethod(pc,"PCBDDCSetDofsSplitting_C",(PC,PetscInt,IS[]),(pc,n_is,ISForDofs));CHKERRQ(ierr); 1032 PetscFunctionReturn(0); 1033 } 1034 1035 /* -------------------------------------------------------------------------- */ 1036 #undef __FUNCT__ 1037 #define __FUNCT__ "PCPreSolve_BDDC" 1038 /* -------------------------------------------------------------------------- */ 1039 /* 1040 PCPreSolve_BDDC - Changes the right hand side and (if necessary) the initial 1041 guess if a transformation of basis approach has been selected. 1042 1043 Input Parameter: 1044 + pc - the preconditioner contex 1045 1046 Application Interface Routine: PCPreSolve() 1047 1048 Notes: 1049 The interface routine PCPreSolve() is not usually called directly by 1050 the user, but instead is called by KSPSolve(). 1051 */ 1052 static PetscErrorCode PCPreSolve_BDDC(PC pc, KSP ksp, Vec rhs, Vec x) 1053 { 1054 PetscErrorCode ierr; 1055 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 1056 PC_IS *pcis = (PC_IS*)(pc->data); 1057 Vec used_vec; 1058 PetscBool copy_rhs = PETSC_TRUE; 1059 PetscBool benign_correction_is_zero = PETSC_FALSE; 1060 PetscBool iscg = PETSC_FALSE; 1061 1062 PetscFunctionBegin; 1063 /* if we are working with cg, one dirichlet solve can be avoided during Krylov iterations */ 1064 if (ksp) { 1065 ierr = PetscObjectTypeCompare((PetscObject)ksp,KSPCG,&iscg);CHKERRQ(ierr); 1066 if (!iscg) { 1067 ierr = PCBDDCSetUseExactDirichlet(pc,PETSC_FALSE);CHKERRQ(ierr); 1068 } 1069 } 1070 /* Creates parallel work vectors used in presolve */ 1071 if (!pcbddc->original_rhs) { 1072 ierr = VecDuplicate(pcis->vec1_global,&pcbddc->original_rhs);CHKERRQ(ierr); 1073 } 1074 if (!pcbddc->temp_solution) { 1075 ierr = VecDuplicate(pcis->vec1_global,&pcbddc->temp_solution);CHKERRQ(ierr); 1076 } 1077 1078 if (x) { 1079 ierr = PetscObjectReference((PetscObject)x);CHKERRQ(ierr); 1080 used_vec = x; 1081 } else { /* it can only happen when calling PCBDDCMatFETIDPGetRHS */ 1082 ierr = PetscObjectReference((PetscObject)pcbddc->temp_solution);CHKERRQ(ierr); 1083 used_vec = pcbddc->temp_solution; 1084 ierr = VecSet(used_vec,0.0);CHKERRQ(ierr); 1085 } 1086 1087 /* hack into ksp data structure since PCPreSolve comes earlier than setting to zero the guess in src/ksp/ksp/interface/itfunc.c */ 1088 if (ksp) { 1089 /* store the flag for the initial guess since it will be restored back during PCPostSolve_BDDC */ 1090 ierr = KSPGetInitialGuessNonzero(ksp,&pcbddc->ksp_guess_nonzero);CHKERRQ(ierr); 1091 if (!pcbddc->ksp_guess_nonzero) { 1092 ierr = VecSet(used_vec,0.0);CHKERRQ(ierr); 1093 } 1094 } 1095 1096 pcbddc->rhs_change = PETSC_FALSE; 1097 /* Take into account zeroed rows -> change rhs and store solution removed */ 1098 if (rhs) { 1099 IS dirIS = NULL; 1100 1101 /* DirichletBoundariesLocal may not be consistent among neighbours; gets a dirichlet dofs IS from graph (may be cached) */ 1102 ierr = PCBDDCGraphGetDirichletDofs(pcbddc->mat_graph,&dirIS);CHKERRQ(ierr); 1103 if (dirIS) { 1104 Mat_IS *matis = (Mat_IS*)pc->pmat->data; 1105 PetscInt dirsize,i,*is_indices; 1106 PetscScalar *array_x; 1107 const PetscScalar *array_diagonal; 1108 1109 ierr = MatGetDiagonal(pc->pmat,pcis->vec1_global);CHKERRQ(ierr); 1110 ierr = VecPointwiseDivide(pcis->vec1_global,rhs,pcis->vec1_global);CHKERRQ(ierr); 1111 ierr = VecScatterBegin(matis->rctx,pcis->vec1_global,pcis->vec2_N,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1112 ierr = VecScatterEnd(matis->rctx,pcis->vec1_global,pcis->vec2_N,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1113 ierr = VecScatterBegin(matis->rctx,used_vec,pcis->vec1_N,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1114 ierr = VecScatterEnd(matis->rctx,used_vec,pcis->vec1_N,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1115 ierr = ISGetLocalSize(dirIS,&dirsize);CHKERRQ(ierr); 1116 ierr = VecGetArray(pcis->vec1_N,&array_x);CHKERRQ(ierr); 1117 ierr = VecGetArrayRead(pcis->vec2_N,&array_diagonal);CHKERRQ(ierr); 1118 ierr = ISGetIndices(dirIS,(const PetscInt**)&is_indices);CHKERRQ(ierr); 1119 for (i=0; i<dirsize; i++) array_x[is_indices[i]] = array_diagonal[is_indices[i]]; 1120 ierr = ISRestoreIndices(dirIS,(const PetscInt**)&is_indices);CHKERRQ(ierr); 1121 ierr = VecRestoreArrayRead(pcis->vec2_N,&array_diagonal);CHKERRQ(ierr); 1122 ierr = VecRestoreArray(pcis->vec1_N,&array_x);CHKERRQ(ierr); 1123 ierr = VecScatterBegin(matis->rctx,pcis->vec1_N,used_vec,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1124 ierr = VecScatterEnd(matis->rctx,pcis->vec1_N,used_vec,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1125 pcbddc->rhs_change = PETSC_TRUE; 1126 ierr = ISDestroy(&dirIS);CHKERRQ(ierr); 1127 } 1128 } 1129 1130 /* remove the computed solution or the initial guess from the rhs */ 1131 if (pcbddc->rhs_change || (ksp && pcbddc->ksp_guess_nonzero) ) { 1132 /* store the original rhs */ 1133 if (copy_rhs) { 1134 ierr = VecCopy(rhs,pcbddc->original_rhs);CHKERRQ(ierr); 1135 copy_rhs = PETSC_FALSE; 1136 } 1137 pcbddc->rhs_change = PETSC_TRUE; 1138 ierr = VecScale(used_vec,-1.0);CHKERRQ(ierr); 1139 ierr = MatMultAdd(pc->mat,used_vec,rhs,rhs);CHKERRQ(ierr); 1140 ierr = VecScale(used_vec,-1.0);CHKERRQ(ierr); 1141 ierr = VecCopy(used_vec,pcbddc->temp_solution);CHKERRQ(ierr); 1142 if (ksp) { 1143 ierr = KSPSetInitialGuessNonzero(ksp,PETSC_FALSE);CHKERRQ(ierr); 1144 } 1145 } 1146 ierr = VecDestroy(&used_vec);CHKERRQ(ierr); 1147 1148 /* When using the benign trick: (TODO: what about FETI-DP?) 1149 - change rhs on pressures (iteration matrix is surely of type MATIS) 1150 - compute initial vector in benign space 1151 */ 1152 if (rhs && pcbddc->benign_saddle_point) { 1153 Mat_IS *matis = (Mat_IS*)(pc->mat->data); 1154 PetscInt i; 1155 1156 /* store the original rhs */ 1157 if (copy_rhs) { 1158 ierr = VecCopy(rhs,pcbddc->original_rhs);CHKERRQ(ierr); 1159 copy_rhs = PETSC_FALSE; 1160 } 1161 1162 /* now change (locally) the basis */ 1163 ierr = VecScatterBegin(matis->rctx,rhs,pcis->vec1_N,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1164 ierr = VecScatterEnd(matis->rctx,rhs,pcis->vec1_N,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1165 if (pcbddc->benign_change) { 1166 ierr = MatMultTranspose(pcbddc->benign_change,pcis->vec1_N,pcis->vec2_N);CHKERRQ(ierr); 1167 ierr = VecScatterBegin(matis->rctx,pcis->vec2_N,rhs,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1168 ierr = VecScatterEnd(matis->rctx,pcis->vec2_N,rhs,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1169 /* swap local iteration matrix (with the benign trick, amat == pmat) */ 1170 if (pcbddc->benign_n && !pcbddc->benign_change_explicit) { 1171 /* benign shell mat mult for matis->A and pcis->A_IB*/ 1172 ierr = PCBDDCBenignShellMat(pc,PETSC_FALSE);CHKERRQ(ierr); 1173 } else { 1174 ierr = PCBDDCBenignPopOrPushB0(pc,PETSC_FALSE);CHKERRQ(ierr); 1175 pcbddc->benign_original_mat = matis->A; 1176 matis->A = pcbddc->local_mat; 1177 } 1178 } else { 1179 ierr = VecScatterBegin(matis->rctx,pcis->vec1_N,rhs,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1180 ierr = VecScatterEnd(matis->rctx,pcis->vec1_N,rhs,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1181 if (pcbddc->ChangeOfBasisMatrix) { 1182 pcbddc->benign_original_mat = matis->A; 1183 matis->A = pcbddc->local_mat; 1184 } 1185 } 1186 pcbddc->rhs_change = PETSC_TRUE; 1187 1188 /* compute u^*_h as in Xuemin Tu's thesis (see Section 4.8.1) */ 1189 /* TODO: what about Stokes? */ 1190 if (!pcbddc->benign_vec) { 1191 ierr = VecDuplicate(rhs,&pcbddc->benign_vec);CHKERRQ(ierr); 1192 } 1193 ierr = PetscMemzero(pcbddc->benign_p0,pcbddc->benign_n*sizeof(PetscScalar));CHKERRQ(ierr); 1194 if (pcbddc->benign_n) { 1195 const PetscScalar *array; 1196 1197 ierr = VecGetArrayRead(pcis->vec2_N,&array);CHKERRQ(ierr); 1198 for (i=0;i<pcbddc->benign_n;i++) pcbddc->benign_p0[i] = array[pcbddc->benign_p0_lidx[i]]; 1199 ierr = VecRestoreArrayRead(pcis->vec2_N,&array);CHKERRQ(ierr); 1200 } 1201 if (pcbddc->benign_null && iscg) { /* this is a workaround, need to understand more */ 1202 PetscBool iszero_l = PETSC_TRUE; 1203 for (i=0;i<pcbddc->benign_n;i++) { 1204 iszero_l = (iszero_l && (PetscAbsScalar(pcbddc->benign_p0[i]) < PETSC_SMALL ? PETSC_TRUE : PETSC_FALSE)); 1205 } 1206 ierr = MPI_Allreduce(&iszero_l,&benign_correction_is_zero,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)pc));CHKERRQ(ierr); 1207 } 1208 if (!benign_correction_is_zero) { /* use the coarse solver only */ 1209 ierr = PCBDDCBenignGetOrSetP0(pc,pcis->vec1_global,PETSC_FALSE);CHKERRQ(ierr); 1210 pcbddc->benign_apply_coarse_only = PETSC_TRUE; 1211 ierr = PCApply_BDDC(pc,pcis->vec1_global,pcbddc->benign_vec);CHKERRQ(ierr); 1212 pcbddc->benign_apply_coarse_only = PETSC_FALSE; 1213 ierr = PetscMemzero(pcbddc->benign_p0,pcbddc->benign_n*sizeof(PetscScalar));CHKERRQ(ierr); 1214 ierr = PCBDDCBenignGetOrSetP0(pc,pcbddc->benign_vec,PETSC_FALSE);CHKERRQ(ierr); 1215 } 1216 } 1217 1218 /* change rhs and iteration matrix if using the change of basis */ 1219 if (pcbddc->ChangeOfBasisMatrix) { 1220 PCBDDCChange_ctx change_ctx; 1221 1222 /* get change ctx */ 1223 ierr = MatShellGetContext(pcbddc->new_global_mat,&change_ctx);CHKERRQ(ierr); 1224 1225 /* set current iteration matrix inside change context (change of basis has been already set into the ctx during PCSetUp) */ 1226 if (!pcbddc->benign_original_mat) { 1227 ierr = MatDestroy(&change_ctx->original_mat);CHKERRQ(ierr); 1228 ierr = PetscObjectReference((PetscObject)pc->mat);CHKERRQ(ierr); 1229 change_ctx->original_mat = pc->mat; 1230 1231 /* change iteration matrix */ 1232 ierr = MatDestroy(&pc->mat);CHKERRQ(ierr); 1233 ierr = PetscObjectReference((PetscObject)pcbddc->new_global_mat);CHKERRQ(ierr); 1234 pc->mat = pcbddc->new_global_mat; 1235 } 1236 /* store the original rhs */ 1237 if (copy_rhs) { 1238 ierr = VecCopy(rhs,pcbddc->original_rhs);CHKERRQ(ierr); 1239 copy_rhs = PETSC_FALSE; 1240 } 1241 1242 /* change rhs */ 1243 ierr = MatMultTranspose(change_ctx->global_change,rhs,pcis->vec1_global);CHKERRQ(ierr); 1244 ierr = VecCopy(pcis->vec1_global,rhs);CHKERRQ(ierr); 1245 pcbddc->rhs_change = PETSC_TRUE; 1246 } 1247 1248 /* remove non-benign solution from the rhs */ 1249 if (pcbddc->benign_saddle_point) { 1250 /* store the original rhs */ 1251 if (copy_rhs) { 1252 ierr = VecCopy(rhs,pcbddc->original_rhs);CHKERRQ(ierr); 1253 copy_rhs = PETSC_FALSE; 1254 } 1255 if (benign_correction_is_zero) { /* still need to understand why it works great */ 1256 ierr = PCBDDCSetUseExactDirichlet(pc,PETSC_FALSE);CHKERRQ(ierr); 1257 ierr = PCApply_BDDC(pc,rhs,pcbddc->benign_vec);CHKERRQ(ierr); 1258 } 1259 ierr = VecScale(pcbddc->benign_vec,-1.0);CHKERRQ(ierr); 1260 ierr = MatMultAdd(pc->mat,pcbddc->benign_vec,rhs,rhs);CHKERRQ(ierr); 1261 pcbddc->rhs_change = PETSC_TRUE; 1262 } 1263 1264 /* set initial guess if using PCG */ 1265 if (x && pcbddc->use_exact_dirichlet_trick) { 1266 ierr = VecSet(x,0.0);CHKERRQ(ierr); 1267 ierr = VecScatterBegin(pcis->global_to_D,rhs,pcis->vec1_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1268 ierr = VecScatterEnd(pcis->global_to_D,rhs,pcis->vec1_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1269 ierr = KSPSolve(pcbddc->ksp_D,pcis->vec1_D,pcis->vec2_D);CHKERRQ(ierr); 1270 ierr = VecScatterBegin(pcis->global_to_D,pcis->vec2_D,x,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1271 ierr = VecScatterEnd(pcis->global_to_D,pcis->vec2_D,x,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1272 if (ksp) { 1273 ierr = KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);CHKERRQ(ierr); 1274 } 1275 } 1276 if (pcbddc->benign_null) { 1277 MatNullSpace null_space; 1278 Vec nullv; 1279 PetscBool isnull; 1280 PetscInt i; 1281 1282 for (i=0;i<pcbddc->benign_n;i++) pcbddc->benign_p0[i] = 1.; 1283 ierr = VecDuplicate(pcis->vec1_global,&nullv);CHKERRQ(ierr); 1284 ierr = VecSet(nullv,0.);CHKERRQ(ierr); 1285 ierr = PCBDDCBenignGetOrSetP0(pc,nullv,PETSC_FALSE);CHKERRQ(ierr); 1286 ierr = VecNormalize(nullv,NULL);CHKERRQ(ierr); 1287 ierr = MatNullSpaceCreate(PetscObjectComm((PetscObject)pc),PETSC_FALSE,1,&nullv,&null_space);CHKERRQ(ierr); 1288 ierr = MatNullSpaceTest(null_space,pc->mat,&isnull);CHKERRQ(ierr); 1289 if (isnull) { 1290 ierr = MatSetNullSpace(pc->mat,null_space);CHKERRQ(ierr); 1291 } 1292 ierr = MatNullSpaceDestroy(&null_space);CHKERRQ(ierr); 1293 ierr = VecDestroy(&nullv);CHKERRQ(ierr); 1294 } 1295 1296 /* remove nullspace if present */ 1297 if (ksp && x && pcbddc->NullSpace) { 1298 ierr = MatNullSpaceRemove(pcbddc->NullSpace,x);CHKERRQ(ierr); 1299 /* store the original rhs */ 1300 if (copy_rhs) { 1301 ierr = VecCopy(rhs,pcbddc->original_rhs);CHKERRQ(ierr); 1302 copy_rhs = PETSC_FALSE; 1303 } 1304 pcbddc->rhs_change = PETSC_TRUE; 1305 ierr = MatNullSpaceRemove(pcbddc->NullSpace,rhs);CHKERRQ(ierr); 1306 } 1307 PetscFunctionReturn(0); 1308 } 1309 1310 /* -------------------------------------------------------------------------- */ 1311 #undef __FUNCT__ 1312 #define __FUNCT__ "PCPostSolve_BDDC" 1313 /* -------------------------------------------------------------------------- */ 1314 /* 1315 PCPostSolve_BDDC - Changes the computed solution if a transformation of basis 1316 approach has been selected. Also, restores rhs to its original state. 1317 1318 Input Parameter: 1319 + pc - the preconditioner contex 1320 1321 Application Interface Routine: PCPostSolve() 1322 1323 Notes: 1324 The interface routine PCPostSolve() is not usually called directly by 1325 the user, but instead is called by KSPSolve(). 1326 */ 1327 static PetscErrorCode PCPostSolve_BDDC(PC pc, KSP ksp, Vec rhs, Vec x) 1328 { 1329 PetscErrorCode ierr; 1330 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 1331 1332 PetscFunctionBegin; 1333 if (pcbddc->ChangeOfBasisMatrix && !pcbddc->benign_original_mat) { 1334 PCBDDCChange_ctx change_ctx; 1335 1336 /* get change ctx */ 1337 ierr = MatShellGetContext(pcbddc->new_global_mat,&change_ctx);CHKERRQ(ierr); 1338 1339 /* restore iteration matrix */ 1340 ierr = MatDestroy(&pc->mat);CHKERRQ(ierr); 1341 ierr = PetscObjectReference((PetscObject)change_ctx->original_mat);CHKERRQ(ierr); 1342 pc->mat = change_ctx->original_mat; 1343 } 1344 1345 /* need to restore the local matrices */ 1346 if (pcbddc->benign_original_mat) { 1347 Mat_IS *matis = (Mat_IS*)pc->mat->data; 1348 1349 if (pcbddc->benign_n && !pcbddc->benign_change_explicit) { 1350 ierr = PCBDDCBenignShellMat(pc,PETSC_TRUE);CHKERRQ(ierr); 1351 } else { 1352 pcbddc->local_mat = matis->A; 1353 matis->A = pcbddc->benign_original_mat; 1354 pcbddc->benign_original_mat = NULL; 1355 } 1356 } 1357 1358 /* get solution in original basis */ 1359 if (x) { 1360 PC_IS *pcis = (PC_IS*)(pc->data); 1361 1362 /* restore solution on pressures */ 1363 if (pcbddc->benign_saddle_point) { 1364 Mat_IS *matis = (Mat_IS*)pc->mat->data; 1365 1366 /* add non-benign solution */ 1367 ierr = VecAXPY(x,-1.0,pcbddc->benign_vec);CHKERRQ(ierr); 1368 1369 /* change basis on pressures for x */ 1370 ierr = VecScatterBegin(matis->rctx,x,pcis->vec1_N,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1371 ierr = VecScatterEnd(matis->rctx,x,pcis->vec1_N,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1372 if (pcbddc->benign_change) { 1373 ierr = MatMult(pcbddc->benign_change,pcis->vec1_N,pcis->vec2_N);CHKERRQ(ierr); 1374 ierr = VecScatterBegin(matis->rctx,pcis->vec2_N,x,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1375 ierr = VecScatterEnd(matis->rctx,pcis->vec2_N,x,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1376 } else { 1377 ierr = VecScatterBegin(matis->rctx,pcis->vec1_N,x,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1378 ierr = VecScatterEnd(matis->rctx,pcis->vec1_N,x,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1379 } 1380 } 1381 1382 /* change basis on x */ 1383 if (pcbddc->ChangeOfBasisMatrix) { 1384 PCBDDCChange_ctx change_ctx; 1385 1386 ierr = MatShellGetContext(pcbddc->new_global_mat,&change_ctx);CHKERRQ(ierr); 1387 ierr = MatMult(change_ctx->global_change,x,pcis->vec1_global);CHKERRQ(ierr); 1388 ierr = VecCopy(pcis->vec1_global,x);CHKERRQ(ierr); 1389 } 1390 } 1391 1392 /* add solution removed in presolve */ 1393 if (x && pcbddc->rhs_change) { 1394 ierr = VecAXPY(x,1.0,pcbddc->temp_solution);CHKERRQ(ierr); 1395 } 1396 1397 /* restore rhs to its original state */ 1398 if (rhs && pcbddc->rhs_change) { 1399 ierr = VecCopy(pcbddc->original_rhs,rhs);CHKERRQ(ierr); 1400 } 1401 pcbddc->rhs_change = PETSC_FALSE; 1402 1403 /* restore ksp guess state */ 1404 if (ksp) { 1405 ierr = KSPSetInitialGuessNonzero(ksp,pcbddc->ksp_guess_nonzero);CHKERRQ(ierr); 1406 } 1407 PetscFunctionReturn(0); 1408 } 1409 /* -------------------------------------------------------------------------- */ 1410 #undef __FUNCT__ 1411 #define __FUNCT__ "PCSetUp_BDDC" 1412 /* -------------------------------------------------------------------------- */ 1413 /* 1414 PCSetUp_BDDC - Prepares for the use of the BDDC preconditioner 1415 by setting data structures and options. 1416 1417 Input Parameter: 1418 + pc - the preconditioner context 1419 1420 Application Interface Routine: PCSetUp() 1421 1422 Notes: 1423 The interface routine PCSetUp() is not usually called directly by 1424 the user, but instead is called by PCApply() if necessary. 1425 */ 1426 PetscErrorCode PCSetUp_BDDC(PC pc) 1427 { 1428 PetscErrorCode ierr; 1429 PC_BDDC* pcbddc = (PC_BDDC*)pc->data; 1430 Mat_IS* matis; 1431 MatNullSpace nearnullspace; 1432 IS zerodiag = NULL; 1433 PetscInt nrows,ncols; 1434 PetscBool computetopography,computesolvers,computesubschurs; 1435 PetscBool computeconstraintsmatrix; 1436 PetscBool new_nearnullspace_provided,ismatis; 1437 1438 PetscFunctionBegin; 1439 ierr = PetscObjectTypeCompare((PetscObject)pc->pmat,MATIS,&ismatis);CHKERRQ(ierr); 1440 if (!ismatis) { 1441 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG,"PCBDDC preconditioner requires matrix of type MATIS"); 1442 } 1443 ierr = MatGetSize(pc->pmat,&nrows,&ncols);CHKERRQ(ierr); 1444 if (nrows != ncols) { 1445 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"PCBDDC preconditioner requires a square preconditioning matrix"); 1446 } 1447 matis = (Mat_IS*)pc->pmat->data; 1448 /* the following lines of code should be replaced by a better logic between PCIS, PCNN, PCBDDC and other future nonoverlapping preconditioners */ 1449 /* For BDDC we need to define a local "Neumann" problem different to that defined in PCISSetup 1450 Also, BDDC directly build the Dirichlet problem */ 1451 /* split work */ 1452 if (pc->setupcalled) { 1453 if (pc->flag == SAME_NONZERO_PATTERN) { 1454 computetopography = PETSC_FALSE; 1455 computesolvers = PETSC_TRUE; 1456 } else { /* DIFFERENT_NONZERO_PATTERN */ 1457 computetopography = PETSC_TRUE; 1458 computesolvers = PETSC_TRUE; 1459 } 1460 } else { 1461 computetopography = PETSC_TRUE; 1462 computesolvers = PETSC_TRUE; 1463 } 1464 if (pcbddc->recompute_topography) { 1465 computetopography = PETSC_TRUE; 1466 } 1467 pcbddc->recompute_topography = computetopography; 1468 computeconstraintsmatrix = PETSC_FALSE; 1469 1470 /* check parameters' compatibility */ 1471 if (pcbddc->adaptive_threshold > 0.0 && !pcbddc->use_deluxe_scaling) { 1472 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Cannot compute adaptive constraints without deluxe scaling. Rerun with -pc_bddc_use_deluxe_scaling"); 1473 } 1474 pcbddc->adaptive_selection = (PetscBool)(pcbddc->adaptive_threshold > 0.0 && pcbddc->use_deluxe_scaling); 1475 pcbddc->adaptive_userdefined = (PetscBool)(pcbddc->adaptive_selection && pcbddc->adaptive_userdefined); 1476 if (pcbddc->adaptive_selection) pcbddc->use_faces = PETSC_TRUE; 1477 1478 computesubschurs = (PetscBool)(pcbddc->adaptive_selection || pcbddc->use_deluxe_scaling); 1479 if (pcbddc->faster_deluxe && pcbddc->adaptive_selection && pcbddc->use_change_of_basis) { 1480 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Cannot compute faster deluxe if adaptivity and change of basis are both requested. Rerun with -pc_bddc_deluxe_faster false"); 1481 } 1482 1483 /* check if the iteration matrix is of type MATIS in case the benign trick has been requested */ 1484 ierr = PetscObjectTypeCompare((PetscObject)pc->mat,MATIS,&ismatis);CHKERRQ(ierr); 1485 if (pcbddc->benign_saddle_point && !ismatis) { 1486 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG,"PCBDDC preconditioner with benign subspace trick requires the iteration matrix to be of type MATIS"); 1487 } 1488 if (pcbddc->benign_saddle_point && pc->mat != pc->pmat) { 1489 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG,"PCBDDC preconditioner with benign subspace trick requires Amat == Pmat"); 1490 } 1491 1492 /* Get stdout for dbg */ 1493 if (pcbddc->dbg_flag) { 1494 if (!pcbddc->dbg_viewer) { 1495 pcbddc->dbg_viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pc)); 1496 ierr = PetscViewerASCIIPushSynchronized(pcbddc->dbg_viewer);CHKERRQ(ierr); 1497 } 1498 ierr = PetscViewerASCIIAddTab(pcbddc->dbg_viewer,2*pcbddc->current_level);CHKERRQ(ierr); 1499 } 1500 1501 if (pcbddc->user_ChangeOfBasisMatrix) { 1502 /* use_change_of_basis flag is used to automatically compute a change of basis from constraints */ 1503 pcbddc->use_change_of_basis = PETSC_FALSE; 1504 ierr = PCBDDCComputeLocalMatrix(pc,pcbddc->user_ChangeOfBasisMatrix);CHKERRQ(ierr); 1505 } else { 1506 ierr = MatDestroy(&pcbddc->local_mat);CHKERRQ(ierr); 1507 ierr = PetscObjectReference((PetscObject)matis->A);CHKERRQ(ierr); 1508 pcbddc->local_mat = matis->A; 1509 } 1510 1511 /* detect local disconnected subdomains if requested and not done before */ 1512 if (pcbddc->detect_disconnected && !pcbddc->n_local_subs) { 1513 ierr = MatDetectDisconnectedComponents(pcbddc->local_mat,PETSC_FALSE,&pcbddc->n_local_subs,&pcbddc->local_subs);CHKERRQ(ierr); 1514 } 1515 1516 /* 1517 change basis on local pressures (aka zerodiag dofs) 1518 This should come earlier then PCISSetUp for extracting the correct subdomain matrices 1519 */ 1520 if (pcbddc->benign_saddle_point) { 1521 PC_IS* pcis = (PC_IS*)pc->data; 1522 1523 if (pcbddc->user_ChangeOfBasisMatrix || pcbddc->use_change_of_basis) pcbddc->benign_change_explicit = PETSC_TRUE; 1524 /* detect local saddle point and change the basis in pcbddc->local_mat */ 1525 ierr = PCBDDCBenignDetectSaddlePoint(pc,&zerodiag);CHKERRQ(ierr); 1526 /* pop B0 mat from local mat */ 1527 ierr = PCBDDCBenignPopOrPushB0(pc,PETSC_TRUE);CHKERRQ(ierr); 1528 /* set flag in pcis to not reuse submatrices during PCISCreate */ 1529 if (pc->flag == SAME_NONZERO_PATTERN && pcbddc->benign_n && !pcbddc->benign_change_explicit && !pcbddc->dbg_flag) { 1530 pcis->reusesubmatrices = PETSC_TRUE; 1531 } else { 1532 pcis->reusesubmatrices = PETSC_FALSE; 1533 } 1534 } 1535 /* propagate relevant information */ 1536 #if !defined(PETSC_USE_COMPLEX) /* workaround for reals */ 1537 if (matis->A->symmetric_set) { 1538 ierr = MatSetOption(pcbddc->local_mat,MAT_HERMITIAN,matis->A->symmetric);CHKERRQ(ierr); 1539 } 1540 #endif 1541 if (matis->A->symmetric_set) { 1542 ierr = MatSetOption(pcbddc->local_mat,MAT_SYMMETRIC,matis->A->symmetric);CHKERRQ(ierr); 1543 } 1544 if (matis->A->spd_set) { 1545 ierr = MatSetOption(pcbddc->local_mat,MAT_SPD,matis->A->spd);CHKERRQ(ierr); 1546 } 1547 1548 /* Set up all the "iterative substructuring" common block without computing solvers */ 1549 { 1550 Mat temp_mat; 1551 1552 temp_mat = matis->A; 1553 matis->A = pcbddc->local_mat; 1554 ierr = PCISSetUp(pc,PETSC_FALSE);CHKERRQ(ierr); 1555 pcbddc->local_mat = matis->A; 1556 matis->A = temp_mat; 1557 } 1558 1559 /* Analyze interface */ 1560 if (computetopography) { 1561 ierr = PCBDDCAnalyzeInterface(pc);CHKERRQ(ierr); 1562 computeconstraintsmatrix = PETSC_TRUE; 1563 } 1564 1565 /* check existence of a divergence free extension, i.e. 1566 b(v_I,p_0) = 0 for all v_I (raise error if not). 1567 Also, check that PCBDDCBenignGetOrSetP0 works */ 1568 #if defined(PETSC_USE_DEBUG) 1569 if (pcbddc->benign_saddle_point) { 1570 ierr = PCBDDCBenignCheck(pc,zerodiag);CHKERRQ(ierr); 1571 } 1572 #endif 1573 ierr = ISDestroy(&zerodiag);CHKERRQ(ierr); 1574 1575 /* Setup local dirichlet solver ksp_D and sub_schurs solvers */ 1576 if (computesolvers) { 1577 PCBDDCSubSchurs sub_schurs=pcbddc->sub_schurs; 1578 1579 if (computesubschurs && computetopography) { 1580 ierr = PCBDDCInitSubSchurs(pc);CHKERRQ(ierr); 1581 } 1582 if (sub_schurs->schur_explicit) { 1583 if (computesubschurs) { 1584 ierr = PCBDDCSetUpSubSchurs(pc);CHKERRQ(ierr); 1585 } 1586 ierr = PCBDDCSetUpLocalSolvers(pc,PETSC_TRUE,PETSC_FALSE);CHKERRQ(ierr); 1587 } else { 1588 ierr = PCBDDCSetUpLocalSolvers(pc,PETSC_TRUE,PETSC_FALSE);CHKERRQ(ierr); 1589 if (computesubschurs) { 1590 ierr = PCBDDCSetUpSubSchurs(pc);CHKERRQ(ierr); 1591 } 1592 } 1593 if (pcbddc->adaptive_selection) { 1594 ierr = PCBDDCAdaptiveSelection(pc);CHKERRQ(ierr); 1595 computeconstraintsmatrix = PETSC_TRUE; 1596 } 1597 } 1598 1599 /* infer if NullSpace object attached to Mat via MatSetNearNullSpace has changed */ 1600 new_nearnullspace_provided = PETSC_FALSE; 1601 ierr = MatGetNearNullSpace(pc->pmat,&nearnullspace);CHKERRQ(ierr); 1602 if (pcbddc->onearnullspace) { /* already used nearnullspace */ 1603 if (!nearnullspace) { /* near null space attached to mat has been destroyed */ 1604 new_nearnullspace_provided = PETSC_TRUE; 1605 } else { 1606 /* determine if the two nullspaces are different (should be lightweight) */ 1607 if (nearnullspace != pcbddc->onearnullspace) { 1608 new_nearnullspace_provided = PETSC_TRUE; 1609 } else { /* maybe the user has changed the content of the nearnullspace so check vectors ObjectStateId */ 1610 PetscInt i; 1611 const Vec *nearnullvecs; 1612 PetscObjectState state; 1613 PetscInt nnsp_size; 1614 ierr = MatNullSpaceGetVecs(nearnullspace,NULL,&nnsp_size,&nearnullvecs);CHKERRQ(ierr); 1615 for (i=0;i<nnsp_size;i++) { 1616 ierr = PetscObjectStateGet((PetscObject)nearnullvecs[i],&state);CHKERRQ(ierr); 1617 if (pcbddc->onearnullvecs_state[i] != state) { 1618 new_nearnullspace_provided = PETSC_TRUE; 1619 break; 1620 } 1621 } 1622 } 1623 } 1624 } else { 1625 if (!nearnullspace) { /* both nearnullspaces are null */ 1626 new_nearnullspace_provided = PETSC_FALSE; 1627 } else { /* nearnullspace attached later */ 1628 new_nearnullspace_provided = PETSC_TRUE; 1629 } 1630 } 1631 1632 /* Setup constraints and related work vectors */ 1633 /* reset primal space flags */ 1634 pcbddc->new_primal_space = PETSC_FALSE; 1635 pcbddc->new_primal_space_local = PETSC_FALSE; 1636 if (computeconstraintsmatrix || new_nearnullspace_provided) { 1637 /* It also sets the primal space flags */ 1638 ierr = PCBDDCConstraintsSetUp(pc);CHKERRQ(ierr); 1639 /* Allocate needed local vectors (which depends on quantities defined during ConstraintsSetUp) */ 1640 ierr = PCBDDCSetUpLocalWorkVectors(pc);CHKERRQ(ierr); 1641 } 1642 1643 if (computesolvers || pcbddc->new_primal_space) { 1644 if (pcbddc->use_change_of_basis) { 1645 PC_IS *pcis = (PC_IS*)(pc->data); 1646 Mat temp_mat = NULL; 1647 1648 if (pcbddc->benign_change) { 1649 /* insert B0 in pcbddc->local_mat */ 1650 ierr = PCBDDCBenignPopOrPushB0(pc,PETSC_FALSE);CHKERRQ(ierr); 1651 /* swap local matrices */ 1652 ierr = MatISGetLocalMat(pc->pmat,&temp_mat);CHKERRQ(ierr); 1653 ierr = PetscObjectReference((PetscObject)temp_mat);CHKERRQ(ierr); 1654 ierr = MatISSetLocalMat(pc->pmat,pcbddc->local_mat);CHKERRQ(ierr); 1655 ierr = MatDestroy(&pcbddc->local_mat);CHKERRQ(ierr); 1656 } 1657 ierr = PCBDDCComputeLocalMatrix(pc,pcbddc->ChangeOfBasisMatrix);CHKERRQ(ierr); 1658 if (pcbddc->benign_change) { 1659 /* restore original matrix */ 1660 ierr = MatISSetLocalMat(pc->pmat,temp_mat);CHKERRQ(ierr); 1661 ierr = PetscObjectDereference((PetscObject)temp_mat);CHKERRQ(ierr); 1662 /* pop B0 from pcbddc->local_mat */ 1663 ierr = PCBDDCBenignPopOrPushB0(pc,PETSC_TRUE);CHKERRQ(ierr); 1664 } 1665 /* get submatrices */ 1666 ierr = MatDestroy(&pcis->A_IB);CHKERRQ(ierr); 1667 ierr = MatDestroy(&pcis->A_BI);CHKERRQ(ierr); 1668 ierr = MatDestroy(&pcis->A_BB);CHKERRQ(ierr); 1669 ierr = MatGetSubMatrix(pcbddc->local_mat,pcis->is_B_local,pcis->is_B_local,MAT_INITIAL_MATRIX,&pcis->A_BB);CHKERRQ(ierr); 1670 ierr = MatGetSubMatrix(pcbddc->local_mat,pcis->is_I_local,pcis->is_B_local,MAT_INITIAL_MATRIX,&pcis->A_IB);CHKERRQ(ierr); 1671 ierr = MatGetSubMatrix(pcbddc->local_mat,pcis->is_B_local,pcis->is_I_local,MAT_INITIAL_MATRIX,&pcis->A_BI);CHKERRQ(ierr); 1672 /* set flag in pcis to not reuse submatrices during PCISCreate */ 1673 pcis->reusesubmatrices = PETSC_FALSE; 1674 } else if (!pcbddc->user_ChangeOfBasisMatrix && !pcbddc->benign_change) { 1675 ierr = MatDestroy(&pcbddc->local_mat);CHKERRQ(ierr); 1676 ierr = PetscObjectReference((PetscObject)matis->A);CHKERRQ(ierr); 1677 pcbddc->local_mat = matis->A; 1678 } 1679 /* SetUp coarse and local Neumann solvers */ 1680 ierr = PCBDDCSetUpSolvers(pc);CHKERRQ(ierr); 1681 /* SetUp Scaling operator */ 1682 ierr = PCBDDCScalingSetUp(pc);CHKERRQ(ierr); 1683 } 1684 pcbddc->recompute_topography = PETSC_FALSE; 1685 1686 if (pcbddc->dbg_flag) { 1687 ierr = PetscViewerASCIISubtractTab(pcbddc->dbg_viewer,2*pcbddc->current_level);CHKERRQ(ierr); 1688 } 1689 PetscFunctionReturn(0); 1690 } 1691 1692 /* -------------------------------------------------------------------------- */ 1693 /* 1694 PCApply_BDDC - Applies the BDDC operator to a vector. 1695 1696 Input Parameters: 1697 + pc - the preconditioner context 1698 - r - input vector (global) 1699 1700 Output Parameter: 1701 . z - output vector (global) 1702 1703 Application Interface Routine: PCApply() 1704 */ 1705 #undef __FUNCT__ 1706 #define __FUNCT__ "PCApply_BDDC" 1707 PetscErrorCode PCApply_BDDC(PC pc,Vec r,Vec z) 1708 { 1709 PC_IS *pcis = (PC_IS*)(pc->data); 1710 PC_BDDC *pcbddc = (PC_BDDC*)(pc->data); 1711 PetscInt n_B = pcis->n_B, n_D = pcis->n - n_B; 1712 PetscErrorCode ierr; 1713 const PetscScalar one = 1.0; 1714 const PetscScalar m_one = -1.0; 1715 const PetscScalar zero = 0.0; 1716 1717 /* This code is similar to that provided in nn.c for PCNN 1718 NN interface preconditioner changed to BDDC 1719 Added support for M_3 preconditioner in the reference article (code is active if pcbddc->switch_static == PETSC_TRUE) */ 1720 1721 PetscFunctionBegin; 1722 if (pcbddc->benign_saddle_point) { /* get p0 from r */ 1723 ierr = PCBDDCBenignGetOrSetP0(pc,r,PETSC_TRUE);CHKERRQ(ierr); 1724 } 1725 if (!pcbddc->use_exact_dirichlet_trick) { 1726 ierr = VecCopy(r,z);CHKERRQ(ierr); 1727 /* First Dirichlet solve */ 1728 ierr = VecScatterBegin(pcis->global_to_D,r,pcis->vec1_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1729 ierr = VecScatterEnd(pcis->global_to_D,r,pcis->vec1_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1730 /* 1731 Assembling right hand side for BDDC operator 1732 - pcis->vec1_D for the Dirichlet part (if needed, i.e. pcbddc->switch_static == PETSC_TRUE) 1733 - pcis->vec1_B the interface part of the global vector z 1734 */ 1735 if (n_D) { 1736 ierr = KSPSolve(pcbddc->ksp_D,pcis->vec1_D,pcis->vec2_D);CHKERRQ(ierr); 1737 ierr = VecScale(pcis->vec2_D,m_one);CHKERRQ(ierr); 1738 if (pcbddc->switch_static) { ierr = MatMultAdd(pcis->A_II,pcis->vec2_D,pcis->vec1_D,pcis->vec1_D);CHKERRQ(ierr); } 1739 ierr = MatMult(pcis->A_BI,pcis->vec2_D,pcis->vec1_B);CHKERRQ(ierr); 1740 } else { 1741 ierr = VecSet(pcis->vec1_B,zero);CHKERRQ(ierr); 1742 } 1743 ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_B,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1744 ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_B,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1745 ierr = PCBDDCScalingRestriction(pc,z,pcis->vec1_B);CHKERRQ(ierr); 1746 } else { 1747 if (!pcbddc->benign_apply_coarse_only) { 1748 if (pcbddc->switch_static) { 1749 ierr = VecSet(pcis->vec1_D,zero);CHKERRQ(ierr); 1750 } 1751 ierr = PCBDDCScalingRestriction(pc,r,pcis->vec1_B);CHKERRQ(ierr); 1752 } 1753 } 1754 1755 /* Apply interface preconditioner 1756 input/output vecs: pcis->vec1_B and pcis->vec1_D */ 1757 ierr = PCBDDCApplyInterfacePreconditioner(pc,PETSC_FALSE);CHKERRQ(ierr); 1758 1759 /* Apply transpose of partition of unity operator */ 1760 ierr = PCBDDCScalingExtension(pc,pcis->vec1_B,z);CHKERRQ(ierr); 1761 1762 /* Second Dirichlet solve and assembling of output */ 1763 ierr = VecScatterBegin(pcis->global_to_B,z,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1764 ierr = VecScatterEnd(pcis->global_to_B,z,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1765 if (n_B) { 1766 ierr = MatMult(pcis->A_IB,pcis->vec1_B,pcis->vec3_D);CHKERRQ(ierr); 1767 if (pcbddc->switch_static) { ierr = MatMultAdd(pcis->A_II,pcis->vec1_D,pcis->vec3_D,pcis->vec3_D);CHKERRQ(ierr); } 1768 } else if (pcbddc->switch_static) { 1769 ierr = MatMult(pcis->A_II,pcis->vec1_D,pcis->vec3_D);CHKERRQ(ierr); 1770 } 1771 ierr = KSPSolve(pcbddc->ksp_D,pcis->vec3_D,pcis->vec4_D);CHKERRQ(ierr); 1772 1773 if (!pcbddc->use_exact_dirichlet_trick) { 1774 if (pcbddc->switch_static) { 1775 ierr = VecAXPBYPCZ(pcis->vec2_D,m_one,one,m_one,pcis->vec4_D,pcis->vec1_D);CHKERRQ(ierr); 1776 } else { 1777 ierr = VecAXPBY(pcis->vec2_D,m_one,m_one,pcis->vec4_D);CHKERRQ(ierr); 1778 } 1779 ierr = VecScatterBegin(pcis->global_to_D,pcis->vec2_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1780 ierr = VecScatterEnd(pcis->global_to_D,pcis->vec2_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1781 } else { 1782 if (pcbddc->switch_static) { 1783 ierr = VecAXPBY(pcis->vec4_D,one,m_one,pcis->vec1_D);CHKERRQ(ierr); 1784 } else { 1785 ierr = VecScale(pcis->vec4_D,m_one);CHKERRQ(ierr); 1786 } 1787 ierr = VecScatterBegin(pcis->global_to_D,pcis->vec4_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1788 ierr = VecScatterEnd(pcis->global_to_D,pcis->vec4_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1789 } 1790 1791 if (pcbddc->benign_saddle_point) { /* set p0 (computed in PCBDDCApplyInterface) */ 1792 ierr = PCBDDCBenignGetOrSetP0(pc,z,PETSC_FALSE);CHKERRQ(ierr); 1793 } 1794 PetscFunctionReturn(0); 1795 } 1796 1797 /* -------------------------------------------------------------------------- */ 1798 /* 1799 PCApplyTranspose_BDDC - Applies the transpose of the BDDC operator to a vector. 1800 1801 Input Parameters: 1802 + pc - the preconditioner context 1803 - r - input vector (global) 1804 1805 Output Parameter: 1806 . z - output vector (global) 1807 1808 Application Interface Routine: PCApplyTranspose() 1809 */ 1810 #undef __FUNCT__ 1811 #define __FUNCT__ "PCApplyTranspose_BDDC" 1812 PetscErrorCode PCApplyTranspose_BDDC(PC pc,Vec r,Vec z) 1813 { 1814 PC_IS *pcis = (PC_IS*)(pc->data); 1815 PC_BDDC *pcbddc = (PC_BDDC*)(pc->data); 1816 PetscInt n_B = pcis->n_B, n_D = pcis->n - n_B; 1817 PetscErrorCode ierr; 1818 const PetscScalar one = 1.0; 1819 const PetscScalar m_one = -1.0; 1820 const PetscScalar zero = 0.0; 1821 1822 PetscFunctionBegin; 1823 if (pcbddc->benign_saddle_point) { /* get p0 from r */ 1824 ierr = PCBDDCBenignGetOrSetP0(pc,r,PETSC_TRUE);CHKERRQ(ierr); 1825 } 1826 if (!pcbddc->use_exact_dirichlet_trick) { 1827 ierr = VecCopy(r,z);CHKERRQ(ierr); 1828 /* First Dirichlet solve */ 1829 ierr = VecScatterBegin(pcis->global_to_D,r,pcis->vec1_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1830 ierr = VecScatterEnd(pcis->global_to_D,r,pcis->vec1_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1831 /* 1832 Assembling right hand side for BDDC operator 1833 - pcis->vec1_D for the Dirichlet part (if needed, i.e. pcbddc->switch_static == PETSC_TRUE) 1834 - pcis->vec1_B the interface part of the global vector z 1835 */ 1836 if (n_D) { 1837 ierr = KSPSolveTranspose(pcbddc->ksp_D,pcis->vec1_D,pcis->vec2_D);CHKERRQ(ierr); 1838 ierr = VecScale(pcis->vec2_D,m_one);CHKERRQ(ierr); 1839 if (pcbddc->switch_static) { ierr = MatMultTransposeAdd(pcis->A_II,pcis->vec2_D,pcis->vec1_D,pcis->vec1_D);CHKERRQ(ierr); } 1840 ierr = MatMultTranspose(pcis->A_IB,pcis->vec2_D,pcis->vec1_B);CHKERRQ(ierr); 1841 } else { 1842 ierr = VecSet(pcis->vec1_B,zero);CHKERRQ(ierr); 1843 } 1844 ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_B,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1845 ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_B,z,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1846 ierr = PCBDDCScalingRestriction(pc,z,pcis->vec1_B);CHKERRQ(ierr); 1847 } else { 1848 if (pcbddc->switch_static) { 1849 ierr = VecSet(pcis->vec1_D,zero);CHKERRQ(ierr); 1850 } 1851 ierr = PCBDDCScalingRestriction(pc,r,pcis->vec1_B);CHKERRQ(ierr); 1852 } 1853 1854 /* Apply interface preconditioner 1855 input/output vecs: pcis->vec1_B and pcis->vec1_D */ 1856 ierr = PCBDDCApplyInterfacePreconditioner(pc,PETSC_TRUE);CHKERRQ(ierr); 1857 1858 /* Apply transpose of partition of unity operator */ 1859 ierr = PCBDDCScalingExtension(pc,pcis->vec1_B,z);CHKERRQ(ierr); 1860 1861 /* Second Dirichlet solve and assembling of output */ 1862 ierr = VecScatterBegin(pcis->global_to_B,z,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1863 ierr = VecScatterEnd(pcis->global_to_B,z,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1864 if (n_B) { 1865 ierr = MatMultTranspose(pcis->A_BI,pcis->vec1_B,pcis->vec3_D);CHKERRQ(ierr); 1866 if (pcbddc->switch_static) { ierr = MatMultTransposeAdd(pcis->A_II,pcis->vec1_D,pcis->vec3_D,pcis->vec3_D);CHKERRQ(ierr); } 1867 } else if (pcbddc->switch_static) { 1868 ierr = MatMultTranspose(pcis->A_II,pcis->vec1_D,pcis->vec3_D);CHKERRQ(ierr); 1869 } 1870 ierr = KSPSolveTranspose(pcbddc->ksp_D,pcis->vec3_D,pcis->vec4_D);CHKERRQ(ierr); 1871 if (!pcbddc->use_exact_dirichlet_trick) { 1872 if (pcbddc->switch_static) { 1873 ierr = VecAXPBYPCZ(pcis->vec2_D,m_one,one,m_one,pcis->vec4_D,pcis->vec1_D);CHKERRQ(ierr); 1874 } else { 1875 ierr = VecAXPBY(pcis->vec2_D,m_one,m_one,pcis->vec4_D);CHKERRQ(ierr); 1876 } 1877 ierr = VecScatterBegin(pcis->global_to_D,pcis->vec2_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1878 ierr = VecScatterEnd(pcis->global_to_D,pcis->vec2_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1879 } else { 1880 if (pcbddc->switch_static) { 1881 ierr = VecAXPBY(pcis->vec4_D,one,m_one,pcis->vec1_D);CHKERRQ(ierr); 1882 } else { 1883 ierr = VecScale(pcis->vec4_D,m_one);CHKERRQ(ierr); 1884 } 1885 ierr = VecScatterBegin(pcis->global_to_D,pcis->vec4_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1886 ierr = VecScatterEnd(pcis->global_to_D,pcis->vec4_D,z,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 1887 } 1888 if (pcbddc->benign_saddle_point) { /* set p0 (computed in PCBDDCApplyInterface) */ 1889 ierr = PCBDDCBenignGetOrSetP0(pc,z,PETSC_FALSE);CHKERRQ(ierr); 1890 } 1891 PetscFunctionReturn(0); 1892 } 1893 /* -------------------------------------------------------------------------- */ 1894 1895 #undef __FUNCT__ 1896 #define __FUNCT__ "PCDestroy_BDDC" 1897 PetscErrorCode PCDestroy_BDDC(PC pc) 1898 { 1899 PC_BDDC *pcbddc = (PC_BDDC*)pc->data; 1900 PetscErrorCode ierr; 1901 1902 PetscFunctionBegin; 1903 /* free data created by PCIS */ 1904 ierr = PCISDestroy(pc);CHKERRQ(ierr); 1905 /* free BDDC custom data */ 1906 ierr = PCBDDCResetCustomization(pc);CHKERRQ(ierr); 1907 /* destroy objects related to topography */ 1908 ierr = PCBDDCResetTopography(pc);CHKERRQ(ierr); 1909 /* free allocated graph structure */ 1910 ierr = PetscFree(pcbddc->mat_graph);CHKERRQ(ierr); 1911 /* free allocated sub schurs structure */ 1912 ierr = PetscFree(pcbddc->sub_schurs);CHKERRQ(ierr); 1913 /* destroy objects for scaling operator */ 1914 ierr = PCBDDCScalingDestroy(pc);CHKERRQ(ierr); 1915 ierr = PetscFree(pcbddc->deluxe_ctx);CHKERRQ(ierr); 1916 /* free solvers stuff */ 1917 ierr = PCBDDCResetSolvers(pc);CHKERRQ(ierr); 1918 /* free global vectors needed in presolve */ 1919 ierr = VecDestroy(&pcbddc->temp_solution);CHKERRQ(ierr); 1920 ierr = VecDestroy(&pcbddc->original_rhs);CHKERRQ(ierr); 1921 /* free stuff for change of basis hooks */ 1922 if (pcbddc->new_global_mat) { 1923 PCBDDCChange_ctx change_ctx; 1924 ierr = MatShellGetContext(pcbddc->new_global_mat,&change_ctx);CHKERRQ(ierr); 1925 ierr = MatDestroy(&change_ctx->original_mat);CHKERRQ(ierr); 1926 ierr = MatDestroy(&change_ctx->global_change);CHKERRQ(ierr); 1927 ierr = VecDestroyVecs(2,&change_ctx->work);CHKERRQ(ierr); 1928 ierr = PetscFree(change_ctx);CHKERRQ(ierr); 1929 } 1930 ierr = MatDestroy(&pcbddc->new_global_mat);CHKERRQ(ierr); 1931 /* remove functions */ 1932 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetChangeOfBasisMat_C",NULL);CHKERRQ(ierr); 1933 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetPrimalVerticesLocalIS_C",NULL);CHKERRQ(ierr); 1934 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetPrimalVerticesIS_C",NULL);CHKERRQ(ierr); 1935 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetCoarseningRatio_C",NULL);CHKERRQ(ierr); 1936 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetLevel_C",NULL);CHKERRQ(ierr); 1937 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetUseExactDirichlet_C",NULL);CHKERRQ(ierr); 1938 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetLevels_C",NULL);CHKERRQ(ierr); 1939 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetNullSpace_C",NULL);CHKERRQ(ierr); 1940 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetDirichletBoundaries_C",NULL);CHKERRQ(ierr); 1941 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetDirichletBoundariesLocal_C",NULL);CHKERRQ(ierr); 1942 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetNeumannBoundaries_C",NULL);CHKERRQ(ierr); 1943 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetNeumannBoundariesLocal_C",NULL);CHKERRQ(ierr); 1944 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCGetDirichletBoundaries_C",NULL);CHKERRQ(ierr); 1945 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCGetDirichletBoundariesLocal_C",NULL);CHKERRQ(ierr); 1946 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCGetNeumannBoundaries_C",NULL);CHKERRQ(ierr); 1947 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCGetNeumannBoundariesLocal_C",NULL);CHKERRQ(ierr); 1948 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetDofsSplitting_C",NULL);CHKERRQ(ierr); 1949 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetDofsSplittingLocal_C",NULL);CHKERRQ(ierr); 1950 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetLocalAdjacencyGraph_C",NULL);CHKERRQ(ierr); 1951 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCCreateFETIDPOperators_C",NULL);CHKERRQ(ierr); 1952 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCMatFETIDPGetRHS_C",NULL);CHKERRQ(ierr); 1953 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCMatFETIDPGetSolution_C",NULL);CHKERRQ(ierr); 1954 /* Free the private data structure */ 1955 ierr = PetscFree(pc->data);CHKERRQ(ierr); 1956 PetscFunctionReturn(0); 1957 } 1958 /* -------------------------------------------------------------------------- */ 1959 1960 #undef __FUNCT__ 1961 #define __FUNCT__ "PCBDDCMatFETIDPGetRHS_BDDC" 1962 static PetscErrorCode PCBDDCMatFETIDPGetRHS_BDDC(Mat fetidp_mat, Vec standard_rhs, Vec fetidp_flux_rhs) 1963 { 1964 FETIDPMat_ctx mat_ctx; 1965 Vec copy_standard_rhs; 1966 PC_IS* pcis; 1967 PC_BDDC* pcbddc; 1968 PetscErrorCode ierr; 1969 1970 PetscFunctionBegin; 1971 ierr = MatShellGetContext(fetidp_mat,&mat_ctx);CHKERRQ(ierr); 1972 pcis = (PC_IS*)mat_ctx->pc->data; 1973 pcbddc = (PC_BDDC*)mat_ctx->pc->data; 1974 1975 /* 1976 change of basis for physical rhs if needed 1977 It also changes the rhs in case of dirichlet boundaries 1978 TODO: better management when FETIDP will have its own class 1979 */ 1980 ierr = VecDuplicate(standard_rhs,©_standard_rhs);CHKERRQ(ierr); 1981 ierr = VecCopy(standard_rhs,copy_standard_rhs);CHKERRQ(ierr); 1982 ierr = PCPreSolve_BDDC(mat_ctx->pc,NULL,copy_standard_rhs,NULL);CHKERRQ(ierr); 1983 /* store vectors for computation of fetidp final solution */ 1984 ierr = VecScatterBegin(pcis->global_to_D,copy_standard_rhs,mat_ctx->temp_solution_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1985 ierr = VecScatterEnd(pcis->global_to_D,copy_standard_rhs,mat_ctx->temp_solution_D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1986 /* scale rhs since it should be unassembled */ 1987 /* TODO use counter scaling? (also below) */ 1988 ierr = VecScatterBegin(pcis->global_to_B,copy_standard_rhs,mat_ctx->temp_solution_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1989 ierr = VecScatterEnd(pcis->global_to_B,copy_standard_rhs,mat_ctx->temp_solution_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 1990 /* Apply partition of unity */ 1991 ierr = VecPointwiseMult(mat_ctx->temp_solution_B,pcis->D,mat_ctx->temp_solution_B);CHKERRQ(ierr); 1992 /* ierr = PCBDDCScalingRestriction(mat_ctx->pc,copy_standard_rhs,mat_ctx->temp_solution_B);CHKERRQ(ierr); */ 1993 if (!pcbddc->switch_static) { 1994 /* compute partially subassembled Schur complement right-hand side */ 1995 ierr = KSPSolve(pcbddc->ksp_D,mat_ctx->temp_solution_D,pcis->vec1_D);CHKERRQ(ierr); 1996 ierr = MatMult(pcis->A_BI,pcis->vec1_D,pcis->vec1_B);CHKERRQ(ierr); 1997 ierr = VecAXPY(mat_ctx->temp_solution_B,-1.0,pcis->vec1_B);CHKERRQ(ierr); 1998 ierr = VecSet(copy_standard_rhs,0.0);CHKERRQ(ierr); 1999 ierr = VecScatterBegin(pcis->global_to_B,mat_ctx->temp_solution_B,copy_standard_rhs,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 2000 ierr = VecScatterEnd(pcis->global_to_B,mat_ctx->temp_solution_B,copy_standard_rhs,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 2001 /* ierr = PCBDDCScalingRestriction(mat_ctx->pc,copy_standard_rhs,mat_ctx->temp_solution_B);CHKERRQ(ierr); */ 2002 ierr = VecScatterBegin(pcis->global_to_B,copy_standard_rhs,mat_ctx->temp_solution_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 2003 ierr = VecScatterEnd(pcis->global_to_B,copy_standard_rhs,mat_ctx->temp_solution_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 2004 ierr = VecPointwiseMult(mat_ctx->temp_solution_B,pcis->D,mat_ctx->temp_solution_B);CHKERRQ(ierr); 2005 } 2006 ierr = VecDestroy(©_standard_rhs);CHKERRQ(ierr); 2007 /* BDDC rhs */ 2008 ierr = VecCopy(mat_ctx->temp_solution_B,pcis->vec1_B);CHKERRQ(ierr); 2009 if (pcbddc->switch_static) { 2010 ierr = VecCopy(mat_ctx->temp_solution_D,pcis->vec1_D);CHKERRQ(ierr); 2011 } 2012 /* apply BDDC */ 2013 ierr = PCBDDCApplyInterfacePreconditioner(mat_ctx->pc,PETSC_FALSE);CHKERRQ(ierr); 2014 /* Application of B_delta and assembling of rhs for fetidp fluxes */ 2015 ierr = VecSet(fetidp_flux_rhs,0.0);CHKERRQ(ierr); 2016 ierr = MatMult(mat_ctx->B_delta,pcis->vec1_B,mat_ctx->lambda_local);CHKERRQ(ierr); 2017 ierr = VecScatterBegin(mat_ctx->l2g_lambda,mat_ctx->lambda_local,fetidp_flux_rhs,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 2018 ierr = VecScatterEnd(mat_ctx->l2g_lambda,mat_ctx->lambda_local,fetidp_flux_rhs,ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 2019 PetscFunctionReturn(0); 2020 } 2021 2022 #undef __FUNCT__ 2023 #define __FUNCT__ "PCBDDCMatFETIDPGetRHS" 2024 /*@ 2025 PCBDDCMatFETIDPGetRHS - Compute the right-hand side for FETI-DP linear system using the physical right-hand side 2026 2027 Collective 2028 2029 Input Parameters: 2030 + fetidp_mat - the FETI-DP matrix object obtained by a call to PCBDDCCreateFETIDPOperators 2031 - standard_rhs - the right-hand side of the original linear system 2032 2033 Output Parameters: 2034 . fetidp_flux_rhs - the right-hand side for the FETI-DP linear system 2035 2036 Level: developer 2037 2038 Notes: 2039 2040 .seealso: PCBDDC, PCBDDCCreateFETIDPOperators, PCBDDCMatFETIDPGetSolution 2041 @*/ 2042 PetscErrorCode PCBDDCMatFETIDPGetRHS(Mat fetidp_mat, Vec standard_rhs, Vec fetidp_flux_rhs) 2043 { 2044 FETIDPMat_ctx mat_ctx; 2045 PetscErrorCode ierr; 2046 2047 PetscFunctionBegin; 2048 ierr = MatShellGetContext(fetidp_mat,&mat_ctx);CHKERRQ(ierr); 2049 ierr = PetscTryMethod(mat_ctx->pc,"PCBDDCMatFETIDPGetRHS_C",(Mat,Vec,Vec),(fetidp_mat,standard_rhs,fetidp_flux_rhs));CHKERRQ(ierr); 2050 PetscFunctionReturn(0); 2051 } 2052 /* -------------------------------------------------------------------------- */ 2053 2054 #undef __FUNCT__ 2055 #define __FUNCT__ "PCBDDCMatFETIDPGetSolution_BDDC" 2056 static PetscErrorCode PCBDDCMatFETIDPGetSolution_BDDC(Mat fetidp_mat, Vec fetidp_flux_sol, Vec standard_sol) 2057 { 2058 FETIDPMat_ctx mat_ctx; 2059 PC_IS* pcis; 2060 PC_BDDC* pcbddc; 2061 PetscErrorCode ierr; 2062 2063 PetscFunctionBegin; 2064 ierr = MatShellGetContext(fetidp_mat,&mat_ctx);CHKERRQ(ierr); 2065 pcis = (PC_IS*)mat_ctx->pc->data; 2066 pcbddc = (PC_BDDC*)mat_ctx->pc->data; 2067 2068 /* apply B_delta^T */ 2069 ierr = VecScatterBegin(mat_ctx->l2g_lambda,fetidp_flux_sol,mat_ctx->lambda_local,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 2070 ierr = VecScatterEnd (mat_ctx->l2g_lambda,fetidp_flux_sol,mat_ctx->lambda_local,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 2071 ierr = MatMultTranspose(mat_ctx->B_delta,mat_ctx->lambda_local,pcis->vec1_B);CHKERRQ(ierr); 2072 /* compute rhs for BDDC application */ 2073 ierr = VecAYPX(pcis->vec1_B,-1.0,mat_ctx->temp_solution_B);CHKERRQ(ierr); 2074 if (pcbddc->switch_static) { 2075 ierr = VecCopy(mat_ctx->temp_solution_D,pcis->vec1_D);CHKERRQ(ierr); 2076 } 2077 /* apply BDDC */ 2078 ierr = PCBDDCApplyInterfacePreconditioner(mat_ctx->pc,PETSC_FALSE);CHKERRQ(ierr); 2079 /* put values into standard global vector */ 2080 ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_B,standard_sol,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 2081 ierr = VecScatterEnd (pcis->global_to_B,pcis->vec1_B,standard_sol,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 2082 if (!pcbddc->switch_static) { 2083 /* compute values into the interior if solved for the partially subassembled Schur complement */ 2084 ierr = MatMult(pcis->A_IB,pcis->vec1_B,pcis->vec1_D);CHKERRQ(ierr); 2085 ierr = VecAXPY(mat_ctx->temp_solution_D,-1.0,pcis->vec1_D);CHKERRQ(ierr); 2086 ierr = KSPSolve(pcbddc->ksp_D,mat_ctx->temp_solution_D,pcis->vec1_D);CHKERRQ(ierr); 2087 } 2088 ierr = VecScatterBegin(pcis->global_to_D,pcis->vec1_D,standard_sol,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 2089 ierr = VecScatterEnd (pcis->global_to_D,pcis->vec1_D,standard_sol,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 2090 /* final change of basis if needed 2091 Is also sums the dirichlet part removed during RHS assembling */ 2092 ierr = PCPostSolve_BDDC(mat_ctx->pc,NULL,NULL,standard_sol);CHKERRQ(ierr); 2093 PetscFunctionReturn(0); 2094 } 2095 2096 #undef __FUNCT__ 2097 #define __FUNCT__ "PCBDDCMatFETIDPGetSolution" 2098 /*@ 2099 PCBDDCMatFETIDPGetSolution - Compute the physical solution using the solution of the FETI-DP linear system 2100 2101 Collective 2102 2103 Input Parameters: 2104 + fetidp_mat - the FETI-DP matrix obtained by a call to PCBDDCCreateFETIDPOperators 2105 - fetidp_flux_sol - the solution of the FETI-DP linear system 2106 2107 Output Parameters: 2108 . standard_sol - the solution defined on the physical domain 2109 2110 Level: developer 2111 2112 Notes: 2113 2114 .seealso: PCBDDC, PCBDDCCreateFETIDPOperators, PCBDDCMatFETIDPGetRHS 2115 @*/ 2116 PetscErrorCode PCBDDCMatFETIDPGetSolution(Mat fetidp_mat, Vec fetidp_flux_sol, Vec standard_sol) 2117 { 2118 FETIDPMat_ctx mat_ctx; 2119 PetscErrorCode ierr; 2120 2121 PetscFunctionBegin; 2122 ierr = MatShellGetContext(fetidp_mat,&mat_ctx);CHKERRQ(ierr); 2123 ierr = PetscTryMethod(mat_ctx->pc,"PCBDDCMatFETIDPGetSolution_C",(Mat,Vec,Vec),(fetidp_mat,fetidp_flux_sol,standard_sol));CHKERRQ(ierr); 2124 PetscFunctionReturn(0); 2125 } 2126 /* -------------------------------------------------------------------------- */ 2127 2128 extern PetscErrorCode FETIDPMatMult(Mat,Vec,Vec); 2129 extern PetscErrorCode FETIDPMatMultTranspose(Mat,Vec,Vec); 2130 extern PetscErrorCode PCBDDCDestroyFETIDPMat(Mat); 2131 extern PetscErrorCode FETIDPPCApply(PC,Vec,Vec); 2132 extern PetscErrorCode FETIDPPCApplyTranspose(PC,Vec,Vec); 2133 extern PetscErrorCode PCBDDCDestroyFETIDPPC(PC); 2134 2135 #undef __FUNCT__ 2136 #define __FUNCT__ "PCBDDCCreateFETIDPOperators_BDDC" 2137 static PetscErrorCode PCBDDCCreateFETIDPOperators_BDDC(PC pc, Mat *fetidp_mat, PC *fetidp_pc) 2138 { 2139 2140 FETIDPMat_ctx fetidpmat_ctx; 2141 Mat newmat; 2142 FETIDPPC_ctx fetidppc_ctx; 2143 PC newpc; 2144 MPI_Comm comm; 2145 PetscErrorCode ierr; 2146 2147 PetscFunctionBegin; 2148 ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr); 2149 /* FETIDP linear matrix */ 2150 ierr = PCBDDCCreateFETIDPMatContext(pc,&fetidpmat_ctx);CHKERRQ(ierr); 2151 ierr = PCBDDCSetupFETIDPMatContext(fetidpmat_ctx);CHKERRQ(ierr); 2152 ierr = MatCreateShell(comm,PETSC_DECIDE,PETSC_DECIDE,fetidpmat_ctx->n_lambda,fetidpmat_ctx->n_lambda,fetidpmat_ctx,&newmat);CHKERRQ(ierr); 2153 ierr = MatShellSetOperation(newmat,MATOP_MULT,(void (*)(void))FETIDPMatMult);CHKERRQ(ierr); 2154 ierr = MatShellSetOperation(newmat,MATOP_MULT_TRANSPOSE,(void (*)(void))FETIDPMatMultTranspose);CHKERRQ(ierr); 2155 ierr = MatShellSetOperation(newmat,MATOP_DESTROY,(void (*)(void))PCBDDCDestroyFETIDPMat);CHKERRQ(ierr); 2156 ierr = MatSetUp(newmat);CHKERRQ(ierr); 2157 /* FETIDP preconditioner */ 2158 ierr = PCBDDCCreateFETIDPPCContext(pc,&fetidppc_ctx);CHKERRQ(ierr); 2159 ierr = PCBDDCSetupFETIDPPCContext(newmat,fetidppc_ctx);CHKERRQ(ierr); 2160 ierr = PCCreate(comm,&newpc);CHKERRQ(ierr); 2161 ierr = PCSetType(newpc,PCSHELL);CHKERRQ(ierr); 2162 ierr = PCShellSetContext(newpc,fetidppc_ctx);CHKERRQ(ierr); 2163 ierr = PCShellSetApply(newpc,FETIDPPCApply);CHKERRQ(ierr); 2164 ierr = PCShellSetApplyTranspose(newpc,FETIDPPCApplyTranspose);CHKERRQ(ierr); 2165 ierr = PCShellSetDestroy(newpc,PCBDDCDestroyFETIDPPC);CHKERRQ(ierr); 2166 ierr = PCSetOperators(newpc,newmat,newmat);CHKERRQ(ierr); 2167 ierr = PCSetUp(newpc);CHKERRQ(ierr); 2168 /* return pointers for objects created */ 2169 *fetidp_mat=newmat; 2170 *fetidp_pc=newpc; 2171 PetscFunctionReturn(0); 2172 } 2173 2174 #undef __FUNCT__ 2175 #define __FUNCT__ "PCBDDCCreateFETIDPOperators" 2176 /*@ 2177 PCBDDCCreateFETIDPOperators - Create FETI-DP operators 2178 2179 Collective 2180 2181 Input Parameters: 2182 . pc - the BDDC preconditioning context (setup should have been called before) 2183 2184 Output Parameters: 2185 + fetidp_mat - shell FETI-DP matrix object 2186 - fetidp_pc - shell Dirichlet preconditioner for FETI-DP matrix 2187 2188 Options Database Keys: 2189 . -fetidp_fullyredundant <false> - use or not a fully redundant set of Lagrange multipliers 2190 2191 Level: developer 2192 2193 Notes: 2194 Currently the only operations provided for FETI-DP matrix are MatMult and MatMultTranspose 2195 2196 .seealso: PCBDDC, PCBDDCMatFETIDPGetRHS, PCBDDCMatFETIDPGetSolution 2197 @*/ 2198 PetscErrorCode PCBDDCCreateFETIDPOperators(PC pc, Mat *fetidp_mat, PC *fetidp_pc) 2199 { 2200 PetscErrorCode ierr; 2201 2202 PetscFunctionBegin; 2203 PetscValidHeaderSpecific(pc,PC_CLASSID,1); 2204 if (pc->setupcalled) { 2205 ierr = PetscUseMethod(pc,"PCBDDCCreateFETIDPOperators_C",(PC,Mat*,PC*),(pc,fetidp_mat,fetidp_pc));CHKERRQ(ierr); 2206 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"You must call PCSetup_BDDC() first \n"); 2207 PetscFunctionReturn(0); 2208 } 2209 /* -------------------------------------------------------------------------- */ 2210 /*MC 2211 PCBDDC - Balancing Domain Decomposition by Constraints. 2212 2213 An implementation of the BDDC preconditioner based on 2214 2215 .vb 2216 [1] C. R. Dohrmann. "An approximate BDDC preconditioner", Numerical Linear Algebra with Applications Volume 14, Issue 2, pages 149-168, March 2007 2217 [2] A. Klawonn and O. B. Widlund. "Dual-Primal FETI Methods for Linear Elasticity", http://cs.nyu.edu/csweb/Research/TechReports/TR2004-855/TR2004-855.pdf 2218 [3] J. Mandel, B. Sousedik, C. R. Dohrmann. "Multispace and Multilevel BDDC", http://arxiv.org/abs/0712.3977 2219 [4] C. Pechstein and C. R. Dohrmann. "Modern domain decomposition methods BDDC, deluxe scaling, and an algebraic approach", Seminar talk, Linz, December 2013, http://people.ricam.oeaw.ac.at/c.pechstein/pechstein-bddc2013.pdf 2220 .ve 2221 2222 The matrix to be preconditioned (Pmat) must be of type MATIS. 2223 2224 Currently works with MATIS matrices with local matrices of type MATSEQAIJ, MATSEQBAIJ or MATSEQSBAIJ, either with real or complex numbers. 2225 2226 It also works with unsymmetric and indefinite problems. 2227 2228 Unlike 'conventional' interface preconditioners, PCBDDC iterates over all degrees of freedom, not just those on the interface. This allows the use of approximate solvers on the subdomains. 2229 2230 Approximate local solvers are automatically adapted for singular linear problems (see [1]) if the user has provided the nullspace using PCBDDCSetNullSpace() 2231 2232 Boundary nodes are split in vertices, edges and faces classes using information from the local to global mapping of dofs and the local connectivity graph of nodes. The latter can be customized by using PCBDDCSetLocalAdjacencyGraph() 2233 Additional information on dofs can be provided by using PCBDDCSetDofsSplitting(), PCBDDCSetDirichletBoundaries(), PCBDDCSetNeumannBoundaries(), and PCBDDCSetPrimalVerticesIS() and their local counterparts. 2234 2235 Constraints can be customized by attaching a MatNullSpace object to the MATIS matrix via MatSetNearNullSpace(). Non-singular modes are retained via SVD. 2236 2237 Change of basis is performed similarly to [2] when requested. When more than one constraint is present on a single connected component (i.e. an edge or a face), a robust method based on local QR factorizations is used. 2238 User defined change of basis can be passed to PCBDDC by using PCBDDCSetChangeOfBasisMat() 2239 2240 The PETSc implementation also supports multilevel BDDC [3]. Coarse grids are partitioned using a MatPartitioning object. 2241 2242 Adaptive selection of primal constraints [4] is supported for SPD systems with high-contrast in the coefficients if MUMPS or MKL_PARDISO are present. Future versions of the code will also consider using PASTIX. 2243 2244 An experimental interface to the FETI-DP method is available. FETI-DP operators could be created using PCBDDCCreateFETIDPOperators(). A stand-alone class for the FETI-DP method will be provided in the next releases. 2245 Deluxe scaling is not supported yet for FETI-DP. 2246 2247 Options Database Keys (some of them, run with -h for a complete list): 2248 2249 . -pc_bddc_use_vertices <true> - use or not vertices in primal space 2250 . -pc_bddc_use_edges <true> - use or not edges in primal space 2251 . -pc_bddc_use_faces <false> - use or not faces in primal space 2252 . -pc_bddc_symmetric <true> - symmetric computation of primal basis functions. Specify false for unsymmetric problems 2253 . -pc_bddc_use_change_of_basis <false> - use change of basis approach (on edges only) 2254 . -pc_bddc_use_change_on_faces <false> - use change of basis approach on faces if change of basis has been requested 2255 . -pc_bddc_switch_static <false> - switches from M_2 (default) to M_3 operator (see reference article [1]) 2256 . -pc_bddc_levels <0> - maximum number of levels for multilevel 2257 . -pc_bddc_coarsening_ratio <8> - number of subdomains which will be aggregated together at the coarser level (e.g. H/h ratio at the coarser level, significative only in the multilevel case) 2258 . -pc_bddc_redistribute <0> - size of a subset of processors where the coarse problem will be remapped (the value is ignored if not at the coarsest level) 2259 . -pc_bddc_use_deluxe_scaling <false> - use deluxe scaling 2260 . -pc_bddc_schur_layers <-1> - select the economic version of deluxe scaling by specifying the number of layers (-1 corresponds to the original deluxe scaling) 2261 . -pc_bddc_adaptive_threshold <0.0> - when a value greater than one is specified, adaptive selection of constraints is performed on edges and faces (requires deluxe scaling and MUMPS or MKL_PARDISO installed) 2262 - -pc_bddc_check_level <0> - set verbosity level of debugging output 2263 2264 Options for Dirichlet, Neumann or coarse solver can be set with 2265 .vb 2266 -pc_bddc_dirichlet_ 2267 -pc_bddc_neumann_ 2268 -pc_bddc_coarse_ 2269 .ve 2270 e.g -pc_bddc_dirichlet_ksp_type richardson -pc_bddc_dirichlet_pc_type gamg. PCBDDC uses by default KPSPREONLY and PCLU. 2271 2272 When using a multilevel approach, solvers' options at the N-th level (N > 1) can be specified as 2273 .vb 2274 -pc_bddc_dirichlet_lN_ 2275 -pc_bddc_neumann_lN_ 2276 -pc_bddc_coarse_lN_ 2277 .ve 2278 Note that level number ranges from the finest (0) to the coarsest (N). 2279 In order to specify options for the BDDC operators at the coarser levels (and not for the solvers), prepend -pc_bddc_coarse_ or -pc_bddc_coarse_l to the option, e.g. 2280 .vb 2281 -pc_bddc_coarse_pc_bddc_adaptive_threshold 5 -pc_bddc_coarse_l1_pc_bddc_redistribute 3 2282 .ve 2283 will use a threshold of 5 for constraints' selection at the first coarse level and will redistribute the coarse problem of the first coarse level on 3 processors 2284 2285 Level: intermediate 2286 2287 Developer notes: 2288 2289 Contributed by Stefano Zampini 2290 2291 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, MATIS 2292 M*/ 2293 2294 #undef __FUNCT__ 2295 #define __FUNCT__ "PCCreate_BDDC" 2296 PETSC_EXTERN PetscErrorCode PCCreate_BDDC(PC pc) 2297 { 2298 PetscErrorCode ierr; 2299 PC_BDDC *pcbddc; 2300 2301 PetscFunctionBegin; 2302 /* Creates the private data structure for this preconditioner and attach it to the PC object. */ 2303 ierr = PetscNewLog(pc,&pcbddc);CHKERRQ(ierr); 2304 pc->data = (void*)pcbddc; 2305 2306 /* create PCIS data structure */ 2307 ierr = PCISCreate(pc);CHKERRQ(ierr); 2308 2309 /* BDDC customization */ 2310 pcbddc->use_local_adj = PETSC_TRUE; 2311 pcbddc->use_vertices = PETSC_TRUE; 2312 pcbddc->use_edges = PETSC_TRUE; 2313 pcbddc->use_faces = PETSC_FALSE; 2314 pcbddc->use_change_of_basis = PETSC_FALSE; 2315 pcbddc->use_change_on_faces = PETSC_FALSE; 2316 pcbddc->switch_static = PETSC_FALSE; 2317 pcbddc->use_nnsp_true = PETSC_FALSE; 2318 pcbddc->use_qr_single = PETSC_FALSE; 2319 pcbddc->symmetric_primal = PETSC_TRUE; 2320 pcbddc->benign_saddle_point = PETSC_FALSE; 2321 pcbddc->dbg_flag = 0; 2322 /* private */ 2323 pcbddc->local_primal_size = 0; 2324 pcbddc->local_primal_size_cc = 0; 2325 pcbddc->local_primal_ref_node = 0; 2326 pcbddc->local_primal_ref_mult = 0; 2327 pcbddc->n_vertices = 0; 2328 pcbddc->primal_indices_local_idxs = 0; 2329 pcbddc->recompute_topography = PETSC_FALSE; 2330 pcbddc->coarse_size = -1; 2331 pcbddc->new_primal_space = PETSC_FALSE; 2332 pcbddc->new_primal_space_local = PETSC_FALSE; 2333 pcbddc->global_primal_indices = 0; 2334 pcbddc->onearnullspace = 0; 2335 pcbddc->onearnullvecs_state = 0; 2336 pcbddc->user_primal_vertices = 0; 2337 pcbddc->user_primal_vertices_local = 0; 2338 pcbddc->NullSpace = 0; 2339 pcbddc->temp_solution = 0; 2340 pcbddc->original_rhs = 0; 2341 pcbddc->local_mat = 0; 2342 pcbddc->ChangeOfBasisMatrix = 0; 2343 pcbddc->user_ChangeOfBasisMatrix = 0; 2344 pcbddc->new_global_mat = 0; 2345 pcbddc->coarse_vec = 0; 2346 pcbddc->coarse_ksp = 0; 2347 pcbddc->coarse_phi_B = 0; 2348 pcbddc->coarse_phi_D = 0; 2349 pcbddc->coarse_psi_B = 0; 2350 pcbddc->coarse_psi_D = 0; 2351 pcbddc->vec1_P = 0; 2352 pcbddc->vec1_R = 0; 2353 pcbddc->vec2_R = 0; 2354 pcbddc->local_auxmat1 = 0; 2355 pcbddc->local_auxmat2 = 0; 2356 pcbddc->R_to_B = 0; 2357 pcbddc->R_to_D = 0; 2358 pcbddc->ksp_D = 0; 2359 pcbddc->ksp_R = 0; 2360 pcbddc->NeumannBoundaries = 0; 2361 pcbddc->NeumannBoundariesLocal = 0; 2362 pcbddc->DirichletBoundaries = 0; 2363 pcbddc->DirichletBoundariesLocal = 0; 2364 pcbddc->user_provided_isfordofs = PETSC_FALSE; 2365 pcbddc->n_ISForDofs = 0; 2366 pcbddc->n_ISForDofsLocal = 0; 2367 pcbddc->ISForDofs = 0; 2368 pcbddc->ISForDofsLocal = 0; 2369 pcbddc->ConstraintMatrix = 0; 2370 pcbddc->use_exact_dirichlet_trick = PETSC_TRUE; 2371 pcbddc->coarse_loc_to_glob = 0; 2372 pcbddc->coarsening_ratio = 8; 2373 pcbddc->coarse_adj_red = 0; 2374 pcbddc->current_level = 0; 2375 pcbddc->max_levels = 0; 2376 pcbddc->use_coarse_estimates = PETSC_FALSE; 2377 pcbddc->redistribute_coarse = 0; 2378 pcbddc->coarse_subassembling = 0; 2379 pcbddc->coarse_subassembling_init = 0; 2380 pcbddc->detect_disconnected = PETSC_FALSE; 2381 pcbddc->n_local_subs = 0; 2382 pcbddc->local_subs = NULL; 2383 2384 /* benign subspace trick */ 2385 pcbddc->benign_change = 0; 2386 pcbddc->benign_vec = 0; 2387 pcbddc->benign_original_mat = 0; 2388 pcbddc->benign_sf = 0; 2389 pcbddc->benign_B0 = 0; 2390 pcbddc->benign_n = 0; 2391 pcbddc->benign_p0 = NULL; 2392 pcbddc->benign_p0_lidx = NULL; 2393 pcbddc->benign_p0_gidx = NULL; 2394 pcbddc->benign_null = PETSC_FALSE; 2395 2396 /* create local graph structure */ 2397 ierr = PCBDDCGraphCreate(&pcbddc->mat_graph);CHKERRQ(ierr); 2398 2399 /* scaling */ 2400 pcbddc->work_scaling = 0; 2401 pcbddc->use_deluxe_scaling = PETSC_FALSE; 2402 pcbddc->faster_deluxe = PETSC_FALSE; 2403 2404 /* create sub schurs structure */ 2405 ierr = PCBDDCSubSchursCreate(&pcbddc->sub_schurs);CHKERRQ(ierr); 2406 pcbddc->sub_schurs_rebuild = PETSC_FALSE; 2407 pcbddc->sub_schurs_layers = -1; 2408 pcbddc->sub_schurs_use_useradj = PETSC_FALSE; 2409 2410 pcbddc->computed_rowadj = PETSC_FALSE; 2411 2412 /* adaptivity */ 2413 pcbddc->adaptive_threshold = 0.0; 2414 pcbddc->adaptive_nmax = 0; 2415 pcbddc->adaptive_nmin = 0; 2416 2417 /* function pointers */ 2418 pc->ops->apply = PCApply_BDDC; 2419 pc->ops->applytranspose = PCApplyTranspose_BDDC; 2420 pc->ops->setup = PCSetUp_BDDC; 2421 pc->ops->destroy = PCDestroy_BDDC; 2422 pc->ops->setfromoptions = PCSetFromOptions_BDDC; 2423 pc->ops->view = 0; 2424 pc->ops->applyrichardson = 0; 2425 pc->ops->applysymmetricleft = 0; 2426 pc->ops->applysymmetricright = 0; 2427 pc->ops->presolve = PCPreSolve_BDDC; 2428 pc->ops->postsolve = PCPostSolve_BDDC; 2429 2430 /* composing function */ 2431 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetChangeOfBasisMat_C",PCBDDCSetChangeOfBasisMat_BDDC);CHKERRQ(ierr); 2432 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetPrimalVerticesLocalIS_C",PCBDDCSetPrimalVerticesLocalIS_BDDC);CHKERRQ(ierr); 2433 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetPrimalVerticesIS_C",PCBDDCSetPrimalVerticesIS_BDDC);CHKERRQ(ierr); 2434 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetCoarseningRatio_C",PCBDDCSetCoarseningRatio_BDDC);CHKERRQ(ierr); 2435 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetLevel_C",PCBDDCSetLevel_BDDC);CHKERRQ(ierr); 2436 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetUseExactDirichlet_C",PCBDDCSetUseExactDirichlet_BDDC);CHKERRQ(ierr); 2437 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetLevels_C",PCBDDCSetLevels_BDDC);CHKERRQ(ierr); 2438 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetNullSpace_C",PCBDDCSetNullSpace_BDDC);CHKERRQ(ierr); 2439 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetDirichletBoundaries_C",PCBDDCSetDirichletBoundaries_BDDC);CHKERRQ(ierr); 2440 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetDirichletBoundariesLocal_C",PCBDDCSetDirichletBoundariesLocal_BDDC);CHKERRQ(ierr); 2441 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetNeumannBoundaries_C",PCBDDCSetNeumannBoundaries_BDDC);CHKERRQ(ierr); 2442 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetNeumannBoundariesLocal_C",PCBDDCSetNeumannBoundariesLocal_BDDC);CHKERRQ(ierr); 2443 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCGetDirichletBoundaries_C",PCBDDCGetDirichletBoundaries_BDDC);CHKERRQ(ierr); 2444 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCGetDirichletBoundariesLocal_C",PCBDDCGetDirichletBoundariesLocal_BDDC);CHKERRQ(ierr); 2445 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCGetNeumannBoundaries_C",PCBDDCGetNeumannBoundaries_BDDC);CHKERRQ(ierr); 2446 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCGetNeumannBoundariesLocal_C",PCBDDCGetNeumannBoundariesLocal_BDDC);CHKERRQ(ierr); 2447 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetDofsSplitting_C",PCBDDCSetDofsSplitting_BDDC);CHKERRQ(ierr); 2448 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetDofsSplittingLocal_C",PCBDDCSetDofsSplittingLocal_BDDC);CHKERRQ(ierr); 2449 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCSetLocalAdjacencyGraph_C",PCBDDCSetLocalAdjacencyGraph_BDDC);CHKERRQ(ierr); 2450 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCCreateFETIDPOperators_C",PCBDDCCreateFETIDPOperators_BDDC);CHKERRQ(ierr); 2451 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCMatFETIDPGetRHS_C",PCBDDCMatFETIDPGetRHS_BDDC);CHKERRQ(ierr); 2452 ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCMatFETIDPGetSolution_C",PCBDDCMatFETIDPGetSolution_BDDC);CHKERRQ(ierr); 2453 PetscFunctionReturn(0); 2454 } 2455 2456