1 #include <../src/ksp/pc/impls/deflation/deflation.h> 2 3 #if defined(HAVE_SLEPC) 4 #include <slepceps.h> 5 #include <slepcbv.h> 6 #endif 7 8 PetscScalar db2[] = {0.7071067811865476,0.7071067811865476}; 9 10 //dec low high 11 PetscScalar db4[] = {-0.12940952255092145,0.22414386804185735,0.836516303737469,0.48296291314469025}; 12 //rec low high 13 //PetscScalar db4[] = {0.48296291314469025,0.836516303737469,0.22414386804185735,-0.12940952255092145}; 14 15 PetscScalar db8[] = {-0.010597401784997278, 16 0.032883011666982945, 17 0.030841381835986965, 18 -0.18703481171888114, 19 -0.02798376941698385, 20 0.6308807679295904, 21 0.7148465705525415, 22 0.23037781330885523}; 23 24 //PetscScalar db8[] = {0.23037781330885523, 25 //0.7148465705525415, 26 //0.6308807679295904, 27 //-0.02798376941698385, 28 //-0.18703481171888114, 29 //0.030841381835986965, 30 //0.032883011666982945, 31 //-0.010597401784997278}; 32 33 PetscScalar db16[] = {-0.00011747678400228192, 34 0.0006754494059985568, 35 -0.0003917403729959771, 36 -0.00487035299301066, 37 0.008746094047015655, 38 0.013981027917015516, 39 -0.04408825393106472, 40 -0.01736930100202211, 41 0.128747426620186, 42 0.00047248457399797254, 43 -0.2840155429624281, 44 -0.015829105256023893, 45 0.5853546836548691, 46 0.6756307362980128, 47 0.3128715909144659, 48 0.05441584224308161}; 49 50 //PetscScalar db16[] = {0.05441584224308161, 51 //0.3128715909144659, 52 //0.6756307362980128, 53 //0.5853546836548691, 54 //-0.015829105256023893, 55 //-0.2840155429624281, 56 //0.00047248457399797254, 57 //0.128747426620186, 58 //-0.01736930100202211, 59 //-0.04408825393106472, 60 //0.013981027917015516, 61 //0.008746094047015655, 62 //-0.00487035299301066, 63 //-0.0003917403729959771, 64 //0.0006754494059985568, 65 //-0.00011747678400228192}; 66 67 PetscScalar biorth22[] = {0.0, 68 -0.1767766952966369, 69 0.3535533905932738, 70 1.0606601717798214, 71 0.3535533905932738, 72 -0.1767766952966369}; 73 74 //PetscScalar biorth22[] = {0.0, 75 //0.3535533905932738, 76 //0.7071067811865476, 77 //0.3535533905932738, 78 //0.0, 79 //0.0}; 80 81 PetscReal meyer[] = {0.0,-1.009999956941423e-12,8.519459636796214e-09,-1.111944952595278e-08,-1.0798819539621958e-08,6.066975741351135e-08,-1.0866516536735883e-07,8.200680650386481e-08,1.1783004497663934e-07,-5.506340565252278e-07,1.1307947017916706e-06,-1.489549216497156e-06,7.367572885903746e-07,3.20544191334478e-06,-1.6312699734552807e-05,6.554305930575149e-05,-0.0006011502343516092,-0.002704672124643725,0.002202534100911002,0.006045814097323304,-0.006387718318497156,-0.011061496392513451,0.015270015130934803,0.017423434103729693,-0.03213079399021176,-0.024348745906078023,0.0637390243228016,0.030655091960824263,-0.13284520043622938,-0.035087555656258346,0.44459300275757724,0.7445855923188063,0.44459300275757724,-0.035087555656258346,-0.13284520043622938,0.030655091960824263,0.0637390243228016,-0.024348745906078023,-0.03213079399021176,0.017423434103729693,0.015270015130934803,-0.011061496392513451,-0.006387718318497156,0.006045814097323304,0.002202534100911002,-0.002704672124643725,-0.0006011502343516092,6.554305930575149e-05,-1.6312699734552807e-05,3.20544191334478e-06,7.367572885903746e-07,-1.489549216497156e-06,1.1307947017916706e-06,-5.506340565252278e-07,1.1783004497663934e-07,8.200680650386481e-08,-1.0866516536735883e-07,6.066975741351135e-08,-1.0798819539621958e-08,-1.111944952595278e-08,8.519459636796214e-09,-1.009999956941423e-12}; 82 83 static PetscErrorCode PCDeflationCreateSpaceJacketHaar(MPI_Comm comm,PetscInt m,PetscInt n,PetscInt M,PetscInt N,PetscBool jacket,Mat *H) 84 { 85 PetscErrorCode ierr; 86 Mat defl; 87 PetscInt i,j,ilo,ihi,alloc=2,*Iidx; 88 PetscReal val,*row; 89 90 PetscFunctionBegin; 91 if (jacket) alloc = 3; 92 ierr = PetscMalloc1(alloc,&row);CHKERRQ(ierr); 93 ierr = PetscMalloc1(alloc,&Iidx);CHKERRQ(ierr); 94 95 val = 1./pow(2,0.5); 96 row[0] = val; 97 row[1] = val; 98 99 /* TODO pass A instead of PC? */ 100 ierr = MatCreate(comm,&defl);CHKERRQ(ierr); 101 ierr = MatSetSizes(defl,m,n,M,N);CHKERRQ(ierr); 102 ierr = MatSetUp(defl);CHKERRQ(ierr); 103 ierr = MatSeqAIJSetPreallocation(defl,alloc,NULL);CHKERRQ(ierr); 104 ierr = MatMPIAIJSetPreallocation(defl,alloc,NULL,alloc,NULL);CHKERRQ(ierr); 105 ierr = MatSetOption(defl,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE);CHKERRQ(ierr); 106 ierr = MatSetOption(defl,MAT_IGNORE_OFF_PROC_ENTRIES,PETSC_TRUE);CHKERRQ(ierr); 107 108 ierr = MatGetOwnershipRange(defl,&ilo,&ihi);CHKERRQ(ierr); 109 for (i=0; i<2; i++) Iidx[i] = i+ilo*2; 110 if (jacket && ihi==M) ihi -=2; 111 if (ihi<ilo) SETERRQ1(comm,PETSC_ERR_ARG_WRONG,"To many cores to assemble Jacket Haar matrix with %d rows",M); 112 for (i=ilo; i<ihi; i++) { 113 ierr = MatSetValues(defl,1,&i,2,Iidx,row,INSERT_VALUES);CHKERRQ(ierr); 114 for (j=0; j<2; j++) Iidx[j] += 2; 115 } 116 if (jacket && ihi == M-2) { 117 for (i=0; i<3; i++) Iidx[i] = i+ilo*2; 118 row[0] = 0.5; row[1] = 0.5; row[2] = val; 119 ierr = MatSetValues(defl,1,&ihi,3,Iidx,row,INSERT_VALUES);CHKERRQ(ierr); 120 ihi += 1; 121 row[2] = -row[2]; 122 ierr = MatSetValues(defl,1,&ihi,3,Iidx,row,INSERT_VALUES);CHKERRQ(ierr); 123 } 124 125 ierr = MatAssemblyBegin(defl,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 126 ierr = MatAssemblyEnd(defl,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 127 128 ierr = PetscFree(row);CHKERRQ(ierr); 129 ierr = PetscFree(Iidx);CHKERRQ(ierr); 130 *H = defl; 131 PetscFunctionReturn(0); 132 } 133 134 static PetscErrorCode PCDeflationCreateSpaceWave(MPI_Comm comm,PetscInt m,PetscInt n,PetscInt M,PetscInt N,PetscInt ncoeffs,PetscScalar *coeffs,PetscBool trunc,Mat *H) 135 { 136 PetscErrorCode ierr; 137 Mat defl; 138 PetscInt i,j,k,ilo,ihi,*Iidx; 139 140 PetscFunctionBegin; 141 ierr = PetscMalloc1(ncoeffs,&Iidx);CHKERRQ(ierr); 142 143 /* TODO pass A instead of PC? */ 144 ierr = MatCreate(comm,&defl);CHKERRQ(ierr); 145 ierr = MatSetSizes(defl,m,n,M,N);CHKERRQ(ierr); 146 ierr = MatSetUp(defl);CHKERRQ(ierr); 147 ierr = MatSeqAIJSetPreallocation(defl,ncoeffs,NULL);CHKERRQ(ierr); 148 ierr = MatMPIAIJSetPreallocation(defl,ncoeffs,NULL,ncoeffs,NULL);CHKERRQ(ierr); 149 ierr = MatSetOption(defl,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE);CHKERRQ(ierr); 150 ierr = MatSetOption(defl,MAT_IGNORE_OFF_PROC_ENTRIES,PETSC_TRUE);CHKERRQ(ierr); 151 152 /* Alg 735 Taswell: fvecmat */ 153 k = ncoeffs -2; 154 if (trunc) k = k/2; 155 156 ierr = MatGetOwnershipRange(defl,&ilo,&ihi);CHKERRQ(ierr); 157 for (i=0; i<ncoeffs; i++) { 158 Iidx[i] = i+ilo*2 -k; 159 if (Iidx[i] >= N) Iidx[i] = PETSC_MIN_INT; 160 } 161 for (i=ilo; i<ihi; i++) { 162 ierr = MatSetValues(defl,1,&i,ncoeffs,Iidx,coeffs,INSERT_VALUES);CHKERRQ(ierr); 163 for (j=0; j<ncoeffs; j++) { 164 Iidx[j] += 2; 165 if (Iidx[j] >= N) Iidx[j] = PETSC_MIN_INT; 166 } 167 } 168 169 ierr = MatAssemblyBegin(defl,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 170 ierr = MatAssemblyEnd(defl,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 171 172 ierr = PetscFree(Iidx);CHKERRQ(ierr); 173 *H = defl; 174 PetscFunctionReturn(0); 175 } 176 177 PetscErrorCode PCDeflationGetSpaceHaar(PC pc,Mat *W,PetscInt size) 178 { 179 PetscErrorCode ierr; 180 Mat A,defl; 181 PetscInt i,j,len,ilo,ihi,*Iidx,m,M; 182 PetscReal *col,val; 183 184 PetscFunctionBegin; 185 /* Haar basis wavelet, level=size */ 186 len = pow(2,size); 187 ierr = PetscMalloc1(len,&col);CHKERRQ(ierr); 188 ierr = PetscMalloc1(len,&Iidx);CHKERRQ(ierr); 189 val = 1./pow(2,size/2.); 190 for (i=0; i<len; i++) col[i] = val; 191 192 /* TODO pass A instead of PC? */ 193 ierr = PCGetOperators(pc,&A,NULL);CHKERRQ(ierr); /* NOTE: Get Pmat instead? */ 194 ierr = MatGetLocalSize(A,&m,NULL);CHKERRQ(ierr); 195 ierr = MatGetSize(A,&M,NULL);CHKERRQ(ierr); 196 ierr = MatCreate(PetscObjectComm((PetscObject)A),&defl);CHKERRQ(ierr); 197 ierr = MatSetSizes(defl,m,PETSC_DECIDE,M,ceil(M/(float)len));CHKERRQ(ierr); 198 ierr = MatSetUp(defl);CHKERRQ(ierr); 199 ierr = MatSeqAIJSetPreallocation(defl,size,NULL);CHKERRQ(ierr); 200 ierr = MatMPIAIJSetPreallocation(defl,size,NULL,size,NULL);CHKERRQ(ierr); 201 ierr = MatSetOption(defl,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE);CHKERRQ(ierr); 202 203 ierr = MatGetOwnershipRangeColumn(defl,&ilo,&ihi);CHKERRQ(ierr); 204 for (i=0; i<len; i++) Iidx[i] = i+ilo*len; 205 if (M%len && ihi == (int)ceil(M/(float)len)) ihi -= 1; 206 for (i=ilo; i<ihi; i++) { 207 ierr = MatSetValues(defl,len,Iidx,1,&i,col,INSERT_VALUES);CHKERRQ(ierr); 208 for (j=0; j<len; j++) Iidx[j] += len; 209 } 210 if (M%len && ihi+1 == ceil(M/(float)len)) { 211 len = M%len; 212 val = 1./pow(pow(2,len),0.5); 213 for (i=0; i<len; i++) col[i] = val; 214 ierr = MatSetValues(defl,len,Iidx,1,&ihi,col,INSERT_VALUES);CHKERRQ(ierr); 215 } 216 217 ierr = MatAssemblyBegin(defl,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 218 ierr = MatAssemblyEnd(defl,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 219 220 ierr = PetscFree(col);CHKERRQ(ierr); 221 ierr = PetscFree(Iidx);CHKERRQ(ierr); 222 *W = defl; 223 PetscFunctionReturn(0); 224 } 225 226 PetscErrorCode PCDeflationGetSpaceJacketHaar(PC pc,Mat *W,PetscInt size) 227 { 228 PetscErrorCode ierr; 229 Mat A,*H,defl; 230 PetscInt i,m,M,Mdefl,Ndefl; 231 PetscBool jh; 232 MPI_Comm comm; 233 234 PetscFunctionBegin; 235 ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr); 236 ierr = PetscMalloc1(size,&H);CHKERRQ(ierr); 237 ierr = PCGetOperators(pc,&A,NULL);CHKERRQ(ierr); /* NOTE: Get Pmat instead? */ 238 ierr = MatGetLocalSize(A,&m,NULL);CHKERRQ(ierr); 239 ierr = MatGetSize(A,&M,NULL);CHKERRQ(ierr); 240 Mdefl = M; 241 Ndefl = M; 242 for (i=0; i<size; i++) { 243 if (Mdefl%2) { 244 jh=PETSC_TRUE; 245 Mdefl = Mdefl/2 +1; 246 } else { 247 jh=PETSC_FALSE; 248 Mdefl = Mdefl/2; 249 } 250 ierr = PCDeflationCreateSpaceJacketHaar(comm,PETSC_DECIDE,m,Mdefl,Ndefl,jh,&H[i]);CHKERRQ(ierr); 251 ierr = MatGetLocalSize(H[i],&m,NULL);CHKERRQ(ierr); 252 Ndefl = Mdefl; 253 } 254 //ierr = MatCreateProd(comm,size,H,&defl);CHKERRQ(ierr); 255 //ierr = MatCreateComposite(comm,size,H,&defl);CHKERRQ(ierr); 256 //ierr = MatCompositeSetType(defl,MAT_COMPOSITE_MULTIPLICATIVE);CHKERRQ(ierr); 257 /* TODO allow implicit */ 258 //ierr = MatCompositeMerge(defl);CHKERRQ(ierr); 259 Mat newmat; 260 defl = H[0]; 261 for (i=0; i<size-1; i++) { 262 ierr = MatMatMult(H[i+1],defl,MAT_INITIAL_MATRIX,PETSC_DECIDE,&newmat);CHKERRQ(ierr); 263 ierr = MatDestroy(&defl);CHKERRQ(ierr); 264 defl = newmat ; 265 } 266 267 ierr = MatTranspose(defl,MAT_INITIAL_MATRIX,W);CHKERRQ(ierr); 268 269 ierr = MatDestroy(&defl);CHKERRQ(ierr); 270 for (i=1; i<size; i++) { 271 ierr = MatDestroy(&H[i]);CHKERRQ(ierr); 272 } 273 ierr = PetscFree(H);CHKERRQ(ierr); 274 PetscFunctionReturn(0); 275 } 276 277 PetscErrorCode PCDeflationGetSpaceWave(PC pc,Mat *W,PetscInt size,PetscInt ncoeffs,PetscScalar *coeffs,PetscBool trunc) 278 { 279 PetscErrorCode ierr; 280 Mat A,*H,defl; 281 PetscInt i,m,M,Mdefl,Ndefl; 282 MPI_Comm comm; 283 284 PetscFunctionBegin; 285 ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr); 286 ierr = PetscMalloc1(size,&H);CHKERRQ(ierr); 287 ierr = PCGetOperators(pc,&A,NULL);CHKERRQ(ierr); /* NOTE: Get Pmat instead? */ 288 ierr = MatGetLocalSize(A,&m,NULL);CHKERRQ(ierr); 289 ierr = MatGetSize(A,&M,NULL);CHKERRQ(ierr); 290 Mdefl = M; 291 Ndefl = M; 292 for (i=0; i<size; i++) { 293 if (Mdefl%2) { 294 if (trunc) { 295 Mdefl = (PetscInt)PetscCeilReal(Mdefl/2.); 296 } else { 297 Mdefl = (PetscInt)PetscFloorReal((ncoeffs+Mdefl-1)/2.); 298 } 299 } else { 300 Mdefl = Mdefl/2; 301 } 302 ierr = PCDeflationCreateSpaceWave(comm,PETSC_DECIDE,m,Mdefl,Ndefl,ncoeffs,coeffs,trunc,&H[i]);CHKERRQ(ierr); 303 ierr = MatGetLocalSize(H[i],&m,NULL);CHKERRQ(ierr); 304 Ndefl = Mdefl; 305 } 306 ierr = MatCreateComposite(comm,size,H,&defl);CHKERRQ(ierr); 307 ierr = MatCompositeSetType(defl,MAT_COMPOSITE_MULTIPLICATIVE);CHKERRQ(ierr); 308 *W = defl; 309 310 for (i=0; i<size; i++) { 311 ierr = MatDestroy(&H[i]);CHKERRQ(ierr); 312 } 313 ierr = PetscFree(H);CHKERRQ(ierr); 314 PetscFunctionReturn(0); 315 } 316 317 PetscErrorCode PCDeflationGetSpaceAggregation(PC pc,Mat *W) 318 { 319 PetscErrorCode ierr; 320 Mat A,defl; 321 PetscInt i,ilo,ihi,*Iidx,m,M; 322 PetscReal *col; 323 MPI_Comm comm; 324 325 PetscFunctionBegin; 326 /* TODO pass A instead of PC? */ 327 ierr = PCGetOperators(pc,&A,NULL);CHKERRQ(ierr); /* NOTE: Get Pmat instead? */ 328 ierr = MatGetOwnershipRangeColumn(A,&ilo,&ihi);CHKERRQ(ierr); 329 ierr = MatGetSize(A,&M,NULL);CHKERRQ(ierr); 330 ierr = PetscObjectGetComm((PetscObject)A,&comm);CHKERRQ(ierr); 331 ierr = MPI_Comm_size(comm,&m);CHKERRQ(ierr); 332 ierr = MatCreate(comm,&defl);CHKERRQ(ierr); 333 ierr = MatSetSizes(defl,ihi-ilo,1,M,m);CHKERRQ(ierr); 334 ierr = MatSetUp(defl);CHKERRQ(ierr); 335 ierr = MatSeqAIJSetPreallocation(defl,1,NULL);CHKERRQ(ierr); 336 ierr = MatMPIAIJSetPreallocation(defl,1,NULL,0,NULL);CHKERRQ(ierr); 337 ierr = MatSetOption(defl,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE);CHKERRQ(ierr); 338 ierr = MatSetOption(defl,MAT_IGNORE_OFF_PROC_ENTRIES,PETSC_TRUE);CHKERRQ(ierr); 339 340 341 ierr = PetscMalloc1(ihi-ilo,&col);CHKERRQ(ierr); 342 ierr = PetscMalloc1(ihi-ilo,&Iidx);CHKERRQ(ierr); 343 for (i=ilo; i<ihi; i++) { 344 Iidx[i-ilo] = i; 345 col[i-ilo] = 1; 346 } 347 ierr = MPI_Comm_rank(comm,&i);CHKERRQ(ierr); 348 ierr = MatSetValues(defl,ihi-ilo,Iidx,1,&i,col,INSERT_VALUES);CHKERRQ(ierr); 349 350 ierr = MatAssemblyBegin(defl,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 351 ierr = MatAssemblyEnd(defl,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 352 353 ierr = PetscFree(col);CHKERRQ(ierr); 354 ierr = PetscFree(Iidx);CHKERRQ(ierr); 355 *W = defl; 356 PetscFunctionReturn(0); 357 } 358 359 PetscErrorCode PCDeflationGetSpaceSLEPc(PC pc,Mat *W,PetscInt size,PetscBool cheapCP) 360 { 361 #if defined(HAVE_SLEPC) 362 PetscErrorCode ierr; 363 PC_Deflation *def = (PC_Deflation*)pc->data; 364 Mat A,defl; 365 Vec vec; 366 EPS eps; 367 PetscScalar *data,*dataScaled,eigval; 368 PetscInt i,nconv,m,M,n=PETSC_DECIDE; 369 PetscBool slepcinit; 370 MPI_Comm comm; 371 372 PetscFunctionBegin; 373 ierr = SlepcInitialized(&slepcinit);CHKERRQ(ierr); 374 if (!slepcinit) { 375 ierr = SlepcInitialize(NULL,NULL,(char*)0,(char*)0);CHKERRQ(ierr); 376 slepcinit = PETSC_TRUE; 377 } 378 ierr = PCGetOperators(pc,&A,NULL);CHKERRQ(ierr); /* NOTE: Get Pmat instead? */ 379 ierr = PetscObjectGetComm((PetscObject)pc,&comm);CHKERRQ(ierr); 380 ierr = EPSCreate(comm,&eps);CHKERRQ(ierr); 381 ierr = EPSSetOperators(eps,A,NULL);CHKERRQ(ierr); 382 ierr = EPSSetProblemType(eps,EPS_HEP);CHKERRQ(ierr); /* Implemented only for def */ 383 ierr = EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);CHKERRQ(ierr); 384 ierr = EPSSetDimensions(eps,size,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr); 385 ierr = EPSSetFromOptions(eps);CHKERRQ(ierr); 386 387 ierr = EPSSolve(eps);CHKERRQ(ierr); 388 ierr = EPSGetConverged(eps,&nconv);CHKERRQ(ierr); 389 if (nconv < size) SETERRQ2(comm,PETSC_ERR_CONV_FAILED,"SLEPc: Number of converged eigenpairs (%d) is less than requested (%d)",nconv,size); 390 ierr = MatCreateVecs(A,NULL,&vec);CHKERRQ(ierr); 391 ierr = MatGetSize(A,&M,NULL);CHKERRQ(ierr); 392 ierr = MatGetLocalSize(A,&m,NULL);CHKERRQ(ierr); 393 ierr = PetscSplitOwnership(comm,&n,&size);CHKERRQ(ierr); 394 ierr = PetscMalloc1(m*size,&data);CHKERRQ(ierr); 395 /* TODO check that eigenvalue is not 0 -> vec is not in Ker A */ 396 for (i=0; i<size; i++) { 397 ierr = VecPlaceArray(vec,&data[i*m]);CHKERRQ(ierr); 398 ierr = EPSGetEigenvector(eps,i,vec,NULL);CHKERRQ(ierr); 399 ierr = VecResetArray(vec);CHKERRQ(ierr); 400 } 401 ierr = MatCreateDense(comm,m,n,M,size,data,&defl);CHKERRQ(ierr); 402 403 if (cheapCP) { 404 ierr = PetscMalloc1(m*size,&dataScaled);CHKERRQ(ierr); 405 for (i=0; i<size; i++) { 406 ierr = VecPlaceArray(vec,&dataScaled[i*m]);CHKERRQ(ierr); 407 ierr = EPSGetEigenpair(eps,i,&eigval,NULL,vec,NULL);CHKERRQ(ierr); 408 ierr = VecScale(vec,eigval);CHKERRQ(ierr); 409 ierr = VecResetArray(vec);CHKERRQ(ierr); 410 } 411 ierr = MatCreateDense(comm,m,n,M,size,dataScaled,&def->AW);CHKERRQ(ierr); 412 } 413 414 //ierr = EPSGetBV(eps,&bv);CHKERRQ(ierr); 415 //ierr = BVCreateMat(bv,&defl);CHKERRQ(ierr); 416 *W = defl; 417 418 ierr = EPSDestroy(&eps);CHKERRQ(ierr); 419 if (slepcinit) ierr = SlepcFinalize();CHKERRQ(ierr); 420 PetscFunctionReturn(0); 421 #else 422 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_CONV_FAILED,"Not compiled with SLEPc support (call make HAVE_SLEPC)"); 423 #endif 424 } 425 426 PetscErrorCode PCDeflationComputeSpace(PC pc) 427 { 428 PetscErrorCode ierr; 429 Mat defl; 430 PetscBool transp=PETSC_TRUE; 431 PC_Deflation *def = (PC_Deflation*)pc->data; 432 433 /* TODO valid header */ 434 PetscFunctionBegin; 435 if (def->spacesize < 1) SETERRQ1(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG,"Wrong PC_DEFLATION Space size specified: %d",def->spacesize); 436 switch (def->spacetype) { 437 case PC_DEFLATION_SPACE_HAAR: 438 transp = PETSC_FALSE; 439 ierr = PCDeflationGetSpaceHaar(pc,&defl,def->spacesize);CHKERRQ(ierr);break; 440 case PC_DEFLATION_SPACE_JACKET_HAAR: 441 transp = PETSC_FALSE; 442 ierr = PCDeflationGetSpaceJacketHaar(pc,&defl,def->spacesize);CHKERRQ(ierr);break; 443 case PC_DEFLATION_SPACE_DB2: 444 ierr = PCDeflationGetSpaceWave(pc,&defl,def->spacesize,2,db2,!def->extendsp);CHKERRQ(ierr);break; 445 case PC_DEFLATION_SPACE_DB4: 446 ierr = PCDeflationGetSpaceWave(pc,&defl,def->spacesize,4,db4,!def->extendsp);CHKERRQ(ierr);break; 447 case PC_DEFLATION_SPACE_DB8: 448 ierr = PCDeflationGetSpaceWave(pc,&defl,def->spacesize,8,db8,!def->extendsp);CHKERRQ(ierr);break; 449 case PC_DEFLATION_SPACE_DB16: 450 ierr = PCDeflationGetSpaceWave(pc,&defl,def->spacesize,16,db16,!def->extendsp);CHKERRQ(ierr);break; 451 case PC_DEFLATION_SPACE_BIORTH22: 452 ierr = PCDeflationGetSpaceWave(pc,&defl,def->spacesize,6,biorth22,!def->extendsp);CHKERRQ(ierr);break; 453 case PC_DEFLATION_SPACE_MEYER: 454 ierr = PCDeflationGetSpaceWave(pc,&defl,def->spacesize,62,meyer,!def->extendsp);CHKERRQ(ierr);break; 455 case PC_DEFLATION_SPACE_AGGREGATION: 456 transp = PETSC_FALSE; 457 ierr = PCDeflationGetSpaceAggregation(pc,&defl);CHKERRQ(ierr);break; 458 case PC_DEFLATION_SPACE_SLEPC: 459 transp = PETSC_FALSE; 460 ierr = PCDeflationGetSpaceSLEPc(pc,&defl,def->spacesize,PETSC_FALSE);CHKERRQ(ierr);break; 461 case PC_DEFLATION_SPACE_SLEPC_CHEAP: 462 transp = PETSC_FALSE; 463 ierr = PCDeflationGetSpaceSLEPc(pc,&defl,def->spacesize,PETSC_TRUE);CHKERRQ(ierr);break; 464 default: SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONG,"Wrong PC_DEFLATION Space Type specified"); 465 } 466 467 ierr = PCDeflationSetSpace(pc,defl,transp);CHKERRQ(ierr); 468 ierr = MatDestroy(&defl);CHKERRQ(ierr); 469 PetscFunctionReturn(0); 470 } 471 472