Lines Matching +full:- +full:m

13   Mat         M;  member
16 /* sampler for Newton-Schultz */
31 /* Newton-Schultz customization */
54 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCH2OpusInferCoordinates_Private()
59 if (pch2opus->sdim) PetscFunctionReturn(PETSC_SUCCESS); in PCH2OpusInferCoordinates_Private()
61 if (!dm) PetscCall(MatGetDM(pc->useAmat ? pc->mat : pc->pmat, &dm)); in PCH2OpusInferCoordinates_Private()
80 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCReset_H2OPUS()
83 pch2opus->sdim = 0; in PCReset_H2OPUS()
84 pch2opus->nlocc = 0; in PCReset_H2OPUS()
85 PetscCall(PetscFree(pch2opus->coords)); in PCReset_H2OPUS()
91 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCSetCoordinates_H2OPUS()
95 if (pch2opus->sdim && sdim == pch2opus->sdim && nlocc == pch2opus->nlocc) { in PCSetCoordinates_H2OPUS()
96 PetscCall(PetscArraycmp(pch2opus->coords, coords, sdim * nlocc, &reset)); in PCSetCoordinates_H2OPUS()
102 PetscCall(PetscMalloc1(sdim * nlocc, &pch2opus->coords)); in PCSetCoordinates_H2OPUS()
103 PetscCall(PetscArraycpy(pch2opus->coords, coords, sdim * nlocc)); in PCSetCoordinates_H2OPUS()
104 pch2opus->sdim = sdim; in PCSetCoordinates_H2OPUS()
105 pch2opus->nlocc = nlocc; in PCSetCoordinates_H2OPUS()
112 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCDestroy_H2OPUS()
116 PetscCall(MatDestroy(&pch2opus->A)); in PCDestroy_H2OPUS()
117 PetscCall(MatDestroy(&pch2opus->M)); in PCDestroy_H2OPUS()
118 PetscCall(MatDestroy(&pch2opus->T)); in PCDestroy_H2OPUS()
119 PetscCall(VecDestroy(&pch2opus->w)); in PCDestroy_H2OPUS()
120 PetscCall(MatDestroy(&pch2opus->S)); in PCDestroy_H2OPUS()
121 PetscCall(VecDestroy(&pch2opus->wns[0])); in PCDestroy_H2OPUS()
122 PetscCall(VecDestroy(&pch2opus->wns[1])); in PCDestroy_H2OPUS()
123 PetscCall(VecDestroy(&pch2opus->wns[2])); in PCDestroy_H2OPUS()
124 PetscCall(VecDestroy(&pch2opus->wns[3])); in PCDestroy_H2OPUS()
125 PetscCall(MatDestroy(&pch2opus->wnsmat[0])); in PCDestroy_H2OPUS()
126 PetscCall(MatDestroy(&pch2opus->wnsmat[1])); in PCDestroy_H2OPUS()
127 PetscCall(MatDestroy(&pch2opus->wnsmat[2])); in PCDestroy_H2OPUS()
128 PetscCall(MatDestroy(&pch2opus->wnsmat[3])); in PCDestroy_H2OPUS()
130 PetscCall(PetscFree(pc->data)); in PCDestroy_H2OPUS()
136 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCSetFromOptions_H2OPUS()
140 …scCall(PetscOptionsInt("-pc_h2opus_maxits", "Maximum number of iterations for Newton-Schultz", NUL… in PCSetFromOptions_H2OPUS()
141 …PetscCall(PetscOptionsBool("-pc_h2opus_monitor", "Monitor Newton-Schultz convergence", NULL, pch2o… in PCSetFromOptions_H2OPUS()
142 …PetscCall(PetscOptionsReal("-pc_h2opus_atol", "Absolute tolerance", NULL, pch2opus->atol, &pch2opu… in PCSetFromOptions_H2OPUS()
143 …PetscCall(PetscOptionsReal("-pc_h2opus_rtol", "Relative tolerance", NULL, pch2opus->rtol, &pch2opu… in PCSetFromOptions_H2OPUS()
144 …onsEnum("-pc_h2opus_norm_type", "Norm type for convergence monitoring", NULL, NormTypes, (PetscEnu… in PCSetFromOptions_H2OPUS()
145 …PetscCall(PetscOptionsInt("-pc_h2opus_hyperorder", "Hyper power order of sampling", NULL, pch2opus in PCSetFromOptions_H2OPUS()
146 …etscCall(PetscOptionsInt("-pc_h2opus_leafsize", "Leaf size when constructed from kernel", NULL, pc… in PCSetFromOptions_H2OPUS()
147 …PetscCall(PetscOptionsReal("-pc_h2opus_eta", "Admissibility condition tolerance", NULL, pch2opus->… in PCSetFromOptions_H2OPUS()
148 …tscCall(PetscOptionsInt("-pc_h2opus_maxrank", "Maximum rank when constructed from matvecs", NULL, … in PCSetFromOptions_H2OPUS()
149 …OptionsInt("-pc_h2opus_samples", "Number of samples to be taken concurrently when constructing fro… in PCSetFromOptions_H2OPUS()
150 …cCall(PetscOptionsReal("-pc_h2opus_mrtol", "Relative tolerance for construction from sampling", NU… in PCSetFromOptions_H2OPUS()
151 …PetscCall(PetscOptionsBool("-pc_h2opus_forcecpu", "Force construction on CPU", NULL, pch2opus->for… in PCSetFromOptions_H2OPUS()
158 Mat M; member
168 /* PetscCall(MatMultTranspose(aat->M,x,aat->w)); */ in MatMult_AAt()
169 PetscCall(MatMultTranspose(aat->A, x, aat->w)); in MatMult_AAt()
170 PetscCall(MatMult(aat->A, aat->w, y)); in MatMult_AAt()
176 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCH2OpusSetUpInit()
177 Mat A = pc->useAmat ? pc->mat : pc->pmat, AAt; in PCH2OpusSetUpInit()
178 PetscInt M, m; in PCH2OpusSetUpInit() local
185 aat.M = pch2opus->M; /* unused so far */ in PCH2OpusSetUpInit()
187 PetscCall(MatGetSize(A, &M, NULL)); in PCH2OpusSetUpInit()
188 PetscCall(MatGetLocalSize(A, &m, NULL)); in PCH2OpusSetUpInit()
189 PetscCall(MatCreateShell(PetscObjectComm((PetscObject)A), m, m, M, M, &aat, &AAt)); in PCH2OpusSetUpInit()
190 PetscCall(MatBindToCPU(AAt, pch2opus->boundtocpu)); in PCH2OpusSetUpInit()
197 pch2opus->s0 = 1. / n; in PCH2OpusSetUpInit()
205 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCApplyKernel_H2OPUS()
208 if (t) PetscCall(MatMultTranspose(pch2opus->M, x, y)); in PCApplyKernel_H2OPUS()
209 else PetscCall(MatMult(pch2opus->M, x, y)); in PCApplyKernel_H2OPUS()
215 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCApplyMatKernel_H2OPUS()
218 if (t) PetscCall(MatTransposeMatMult(pch2opus->M, X, MAT_REUSE_MATRIX, PETSC_DETERMINE, &Y)); in PCApplyMatKernel_H2OPUS()
219 else PetscCall(MatMatMult(pch2opus->M, X, MAT_REUSE_MATRIX, PETSC_CURRENT, &Y)); in PCApplyMatKernel_H2OPUS()
251 /* used to test the norm of (M^-1 A - I) */
252 static PetscErrorCode MatMultKernel_MAmI(Mat M, Vec x, Vec y, PetscBool t) in MatMultKernel_MAmI() argument
260 PetscCall(MatShellGetContext(M, &pc)); in MatMultKernel_MAmI()
261 pch2opus = (PC_H2OPUS *)pc->data; in MatMultKernel_MAmI()
262 if (!pch2opus->w) PetscCall(MatCreateVecs(pch2opus->M, &pch2opus->w, NULL)); in MatMultKernel_MAmI()
263 A = pch2opus->A; in MatMultKernel_MAmI()
264 PetscCall(VecBindToCPU(pch2opus->w, pch2opus->boundtocpu)); in MatMultKernel_MAmI()
267 PetscCall(PCApplyTranspose_H2OPUS(pc, x, pch2opus->w)); in MatMultKernel_MAmI()
268 PetscCall(MatMultTranspose(A, pch2opus->w, y)); in MatMultKernel_MAmI()
270 PetscCall(MatMultTranspose(A, x, pch2opus->w)); in MatMultKernel_MAmI()
271 PetscCall(PCApplyTranspose_H2OPUS(pc, pch2opus->w, y)); in MatMultKernel_MAmI()
275 PetscCall(MatMult(A, x, pch2opus->w)); in MatMultKernel_MAmI()
276 PetscCall(PCApply_H2OPUS(pc, pch2opus->w, y)); in MatMultKernel_MAmI()
278 PetscCall(PCApply_H2OPUS(pc, x, pch2opus->w)); in MatMultKernel_MAmI()
279 PetscCall(MatMult(A, pch2opus->w, y)); in MatMultKernel_MAmI()
282 PetscCall(VecAXPY(y, -1.0, x)); in MatMultKernel_MAmI()
302 for i = 1 . . . l - 1 do
303 R = (I - A * Xk) * R
307 static PetscErrorCode MatMultKernel_Hyper(Mat M, Vec x, Vec y, PetscBool t) in MatMultKernel_Hyper() argument
314 PetscCall(MatShellGetContext(M, &pc)); in MatMultKernel_Hyper()
315 pch2opus = (PC_H2OPUS *)pc->data; in MatMultKernel_Hyper()
316 A = pch2opus->A; in MatMultKernel_Hyper()
317 …PetscCall(MatCreateVecs(pch2opus->M, pch2opus->wns[0] ? NULL : &pch2opus->wns[0], pch2opus->wns[1]… in MatMultKernel_Hyper()
318 …PetscCall(MatCreateVecs(pch2opus->M, pch2opus->wns[2] ? NULL : &pch2opus->wns[2], pch2opus->wns[3]… in MatMultKernel_Hyper()
319 PetscCall(VecBindToCPU(pch2opus->wns[0], pch2opus->boundtocpu)); in MatMultKernel_Hyper()
320 PetscCall(VecBindToCPU(pch2opus->wns[1], pch2opus->boundtocpu)); in MatMultKernel_Hyper()
321 PetscCall(VecBindToCPU(pch2opus->wns[2], pch2opus->boundtocpu)); in MatMultKernel_Hyper()
322 PetscCall(VecBindToCPU(pch2opus->wns[3], pch2opus->boundtocpu)); in MatMultKernel_Hyper()
323 PetscCall(VecCopy(x, pch2opus->wns[0])); in MatMultKernel_Hyper()
324 PetscCall(VecCopy(x, pch2opus->wns[3])); in MatMultKernel_Hyper()
326 for (PetscInt i = 0; i < pch2opus->hyperorder - 1; i++) { in MatMultKernel_Hyper()
327 PetscCall(MatMultTranspose(A, pch2opus->wns[0], pch2opus->wns[1])); in MatMultKernel_Hyper()
328 PetscCall(PCApplyTranspose_H2OPUS(pc, pch2opus->wns[1], pch2opus->wns[2])); in MatMultKernel_Hyper()
329 PetscCall(VecAXPY(pch2opus->wns[0], -1., pch2opus->wns[2])); in MatMultKernel_Hyper()
330 PetscCall(VecAXPY(pch2opus->wns[3], 1., pch2opus->wns[0])); in MatMultKernel_Hyper()
332 PetscCall(PCApplyTranspose_H2OPUS(pc, pch2opus->wns[3], y)); in MatMultKernel_Hyper()
334 for (PetscInt i = 0; i < pch2opus->hyperorder - 1; i++) { in MatMultKernel_Hyper()
335 PetscCall(PCApply_H2OPUS(pc, pch2opus->wns[0], pch2opus->wns[1])); in MatMultKernel_Hyper()
336 PetscCall(MatMult(A, pch2opus->wns[1], pch2opus->wns[2])); in MatMultKernel_Hyper()
337 PetscCall(VecAXPY(pch2opus->wns[0], -1., pch2opus->wns[2])); in MatMultKernel_Hyper()
338 PetscCall(VecAXPY(pch2opus->wns[3], 1., pch2opus->wns[0])); in MatMultKernel_Hyper()
340 PetscCall(PCApply_H2OPUS(pc, pch2opus->wns[3], y)); in MatMultKernel_Hyper()
345 static PetscErrorCode MatMult_Hyper(Mat M, Vec x, Vec y) in MatMult_Hyper() argument
348 PetscCall(MatMultKernel_Hyper(M, x, y, PETSC_FALSE)); in MatMult_Hyper()
352 static PetscErrorCode MatMultTranspose_Hyper(Mat M, Vec x, Vec y) in MatMultTranspose_Hyper() argument
355 PetscCall(MatMultKernel_Hyper(M, x, y, PETSC_TRUE)); in MatMultTranspose_Hyper()
360 static PetscErrorCode MatMatMultKernel_Hyper(Mat M, Mat X, Mat Y, PetscBool t) in MatMatMultKernel_Hyper() argument
368 PetscCall(MatShellGetContext(M, &pc)); in MatMatMultKernel_Hyper()
369 pch2opus = (PC_H2OPUS *)pc->data; in MatMatMultKernel_Hyper()
370 A = pch2opus->A; in MatMatMultKernel_Hyper()
371 if (pch2opus->wnsmat[0] && pch2opus->wnsmat[0]->cmap->N != X->cmap->N) { in MatMatMultKernel_Hyper()
372 PetscCall(MatDestroy(&pch2opus->wnsmat[0])); in MatMatMultKernel_Hyper()
373 PetscCall(MatDestroy(&pch2opus->wnsmat[1])); in MatMatMultKernel_Hyper()
375 if (!pch2opus->wnsmat[0]) { in MatMatMultKernel_Hyper()
376 PetscCall(MatDuplicate(X, MAT_SHARE_NONZERO_PATTERN, &pch2opus->wnsmat[0])); in MatMatMultKernel_Hyper()
377 PetscCall(MatDuplicate(Y, MAT_SHARE_NONZERO_PATTERN, &pch2opus->wnsmat[1])); in MatMatMultKernel_Hyper()
379 if (pch2opus->wnsmat[2] && pch2opus->wnsmat[2]->cmap->N != X->cmap->N) { in MatMatMultKernel_Hyper()
380 PetscCall(MatDestroy(&pch2opus->wnsmat[2])); in MatMatMultKernel_Hyper()
381 PetscCall(MatDestroy(&pch2opus->wnsmat[3])); in MatMatMultKernel_Hyper()
383 if (!pch2opus->wnsmat[2]) { in MatMatMultKernel_Hyper()
384 PetscCall(MatDuplicate(X, MAT_SHARE_NONZERO_PATTERN, &pch2opus->wnsmat[2])); in MatMatMultKernel_Hyper()
385 PetscCall(MatDuplicate(Y, MAT_SHARE_NONZERO_PATTERN, &pch2opus->wnsmat[3])); in MatMatMultKernel_Hyper()
387 PetscCall(MatCopy(X, pch2opus->wnsmat[0], SAME_NONZERO_PATTERN)); in MatMatMultKernel_Hyper()
388 PetscCall(MatCopy(X, pch2opus->wnsmat[3], SAME_NONZERO_PATTERN)); in MatMatMultKernel_Hyper()
390 for (i = 0; i < pch2opus->hyperorder - 1; i++) { in MatMatMultKernel_Hyper()
391 …PetscCall(MatTransposeMatMult(A, pch2opus->wnsmat[0], MAT_REUSE_MATRIX, PETSC_CURRENT, &pch2opus->… in MatMatMultKernel_Hyper()
392 PetscCall(PCApplyTransposeMat_H2OPUS(pc, pch2opus->wnsmat[1], pch2opus->wnsmat[2])); in MatMatMultKernel_Hyper()
393 PetscCall(MatAXPY(pch2opus->wnsmat[0], -1., pch2opus->wnsmat[2], SAME_NONZERO_PATTERN)); in MatMatMultKernel_Hyper()
394 PetscCall(MatAXPY(pch2opus->wnsmat[3], 1., pch2opus->wnsmat[0], SAME_NONZERO_PATTERN)); in MatMatMultKernel_Hyper()
396 PetscCall(PCApplyTransposeMat_H2OPUS(pc, pch2opus->wnsmat[3], Y)); in MatMatMultKernel_Hyper()
398 for (i = 0; i < pch2opus->hyperorder - 1; i++) { in MatMatMultKernel_Hyper()
399 PetscCall(PCApplyMat_H2OPUS(pc, pch2opus->wnsmat[0], pch2opus->wnsmat[1])); in MatMatMultKernel_Hyper()
400 …PetscCall(MatMatMult(A, pch2opus->wnsmat[1], MAT_REUSE_MATRIX, PETSC_CURRENT, &pch2opus->wnsmat[2]… in MatMatMultKernel_Hyper()
401 PetscCall(MatAXPY(pch2opus->wnsmat[0], -1., pch2opus->wnsmat[2], SAME_NONZERO_PATTERN)); in MatMatMultKernel_Hyper()
402 PetscCall(MatAXPY(pch2opus->wnsmat[3], 1., pch2opus->wnsmat[0], SAME_NONZERO_PATTERN)); in MatMatMultKernel_Hyper()
404 PetscCall(PCApplyMat_H2OPUS(pc, pch2opus->wnsmat[3], Y)); in MatMatMultKernel_Hyper()
409 static PetscErrorCode MatMatMultNumeric_Hyper(Mat M, Mat X, Mat Y, PetscCtx ctx) in MatMatMultNumeric_Hyper() argument
412 PetscCall(MatMatMultKernel_Hyper(M, X, Y, PETSC_FALSE)); in MatMatMultNumeric_Hyper()
416 /* Basic Newton-Schultz sampler: (2 * I - M * A)*M */
417 static PetscErrorCode MatMultKernel_NS(Mat M, Vec x, Vec y, PetscBool t) in MatMultKernel_NS() argument
424 PetscCall(MatShellGetContext(M, &pc)); in MatMultKernel_NS()
425 pch2opus = (PC_H2OPUS *)pc->data; in MatMultKernel_NS()
426 A = pch2opus->A; in MatMultKernel_NS()
427 …PetscCall(MatCreateVecs(pch2opus->M, pch2opus->wns[0] ? NULL : &pch2opus->wns[0], pch2opus->wns[1]… in MatMultKernel_NS()
428 PetscCall(VecBindToCPU(pch2opus->wns[0], pch2opus->boundtocpu)); in MatMultKernel_NS()
429 PetscCall(VecBindToCPU(pch2opus->wns[1], pch2opus->boundtocpu)); in MatMultKernel_NS()
432 PetscCall(MatMultTranspose(A, y, pch2opus->wns[1])); in MatMultKernel_NS()
433 PetscCall(PCApplyTranspose_H2OPUS(pc, pch2opus->wns[1], pch2opus->wns[0])); in MatMultKernel_NS()
434 PetscCall(VecAXPBY(y, -1., 2., pch2opus->wns[0])); in MatMultKernel_NS()
437 PetscCall(MatMult(A, y, pch2opus->wns[0])); in MatMultKernel_NS()
438 PetscCall(PCApply_H2OPUS(pc, pch2opus->wns[0], pch2opus->wns[1])); in MatMultKernel_NS()
439 PetscCall(VecAXPBY(y, -1., 2., pch2opus->wns[1])); in MatMultKernel_NS()
444 static PetscErrorCode MatMult_NS(Mat M, Vec x, Vec y) in MatMult_NS() argument
447 PetscCall(MatMultKernel_NS(M, x, y, PETSC_FALSE)); in MatMult_NS()
451 static PetscErrorCode MatMultTranspose_NS(Mat M, Vec x, Vec y) in MatMultTranspose_NS() argument
454 PetscCall(MatMultKernel_NS(M, x, y, PETSC_TRUE)); in MatMultTranspose_NS()
458 /* Basic Newton-Schultz sampler: (2 * I - M * A)*M, MatMat version */
459 static PetscErrorCode MatMatMultKernel_NS(Mat M, Mat X, Mat Y, PetscBool t) in MatMatMultKernel_NS() argument
466 PetscCall(MatShellGetContext(M, &pc)); in MatMatMultKernel_NS()
467 pch2opus = (PC_H2OPUS *)pc->data; in MatMatMultKernel_NS()
468 A = pch2opus->A; in MatMatMultKernel_NS()
469 if (pch2opus->wnsmat[0] && pch2opus->wnsmat[0]->cmap->N != X->cmap->N) { in MatMatMultKernel_NS()
470 PetscCall(MatDestroy(&pch2opus->wnsmat[0])); in MatMatMultKernel_NS()
471 PetscCall(MatDestroy(&pch2opus->wnsmat[1])); in MatMatMultKernel_NS()
473 if (!pch2opus->wnsmat[0]) { in MatMatMultKernel_NS()
474 PetscCall(MatDuplicate(X, MAT_SHARE_NONZERO_PATTERN, &pch2opus->wnsmat[0])); in MatMatMultKernel_NS()
475 PetscCall(MatDuplicate(Y, MAT_SHARE_NONZERO_PATTERN, &pch2opus->wnsmat[1])); in MatMatMultKernel_NS()
479 PetscCall(MatTransposeMatMult(A, Y, MAT_REUSE_MATRIX, PETSC_CURRENT, &pch2opus->wnsmat[1])); in MatMatMultKernel_NS()
480 PetscCall(PCApplyTransposeMat_H2OPUS(pc, pch2opus->wnsmat[1], pch2opus->wnsmat[0])); in MatMatMultKernel_NS()
482 PetscCall(MatAXPY(Y, -1., pch2opus->wnsmat[0], SAME_NONZERO_PATTERN)); in MatMatMultKernel_NS()
485 PetscCall(MatMatMult(A, Y, MAT_REUSE_MATRIX, PETSC_CURRENT, &pch2opus->wnsmat[0])); in MatMatMultKernel_NS()
486 PetscCall(PCApplyMat_H2OPUS(pc, pch2opus->wnsmat[0], pch2opus->wnsmat[1])); in MatMatMultKernel_NS()
488 PetscCall(MatAXPY(Y, -1., pch2opus->wnsmat[1], SAME_NONZERO_PATTERN)); in MatMatMultKernel_NS()
493 static PetscErrorCode MatMatMultNumeric_NS(Mat M, Mat X, Mat Y, PetscCtx ctx) in MatMatMultNumeric_NS() argument
496 PetscCall(MatMatMultKernel_NS(M, X, Y, PETSC_FALSE)); in MatMatMultNumeric_NS()
502 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCH2OpusSetUpSampler_Private()
503 Mat A = pc->useAmat ? pc->mat : pc->pmat; in PCH2OpusSetUpSampler_Private()
506 if (!pch2opus->S) { in PCH2OpusSetUpSampler_Private()
507 PetscInt M, N, m, n; in PCH2OpusSetUpSampler_Private() local
509 PetscCall(MatGetSize(A, &M, &N)); in PCH2OpusSetUpSampler_Private()
510 PetscCall(MatGetLocalSize(A, &m, &n)); in PCH2OpusSetUpSampler_Private()
511 PetscCall(MatCreateShell(PetscObjectComm((PetscObject)A), m, n, M, N, pc, &pch2opus->S)); in PCH2OpusSetUpSampler_Private()
512 PetscCall(MatSetBlockSizesFromMats(pch2opus->S, A, A)); in PCH2OpusSetUpSampler_Private()
514 PetscCall(MatShellSetVecType(pch2opus->S, VECCUDA)); in PCH2OpusSetUpSampler_Private()
517 if (pch2opus->hyperorder >= 2) { in PCH2OpusSetUpSampler_Private()
518 PetscCall(MatShellSetOperation(pch2opus->S, MATOP_MULT, (PetscErrorCodeFn *)MatMult_Hyper)); in PCH2OpusSetUpSampler_Private()
519 …PetscCall(MatShellSetOperation(pch2opus->S, MATOP_MULT_TRANSPOSE, (PetscErrorCodeFn *)MatMultTrans… in PCH2OpusSetUpSampler_Private()
520 …PetscCall(MatShellSetMatProductOperation(pch2opus->S, MATPRODUCT_AB, NULL, MatMatMultNumeric_Hyper… in PCH2OpusSetUpSampler_Private()
521 …PetscCall(MatShellSetMatProductOperation(pch2opus->S, MATPRODUCT_AB, NULL, MatMatMultNumeric_Hyper… in PCH2OpusSetUpSampler_Private()
523 PetscCall(MatShellSetOperation(pch2opus->S, MATOP_MULT, (PetscErrorCodeFn *)MatMult_NS)); in PCH2OpusSetUpSampler_Private()
524 …PetscCall(MatShellSetOperation(pch2opus->S, MATOP_MULT_TRANSPOSE, (PetscErrorCodeFn *)MatMultTrans… in PCH2OpusSetUpSampler_Private()
525 …PetscCall(MatShellSetMatProductOperation(pch2opus->S, MATPRODUCT_AB, NULL, MatMatMultNumeric_NS, N… in PCH2OpusSetUpSampler_Private()
526 …PetscCall(MatShellSetMatProductOperation(pch2opus->S, MATPRODUCT_AB, NULL, MatMatMultNumeric_NS, N… in PCH2OpusSetUpSampler_Private()
528 PetscCall(MatPropagateSymmetryOptions(A, pch2opus->S)); in PCH2OpusSetUpSampler_Private()
529 PetscCall(MatBindToCPU(pch2opus->S, pch2opus->boundtocpu)); in PCH2OpusSetUpSampler_Private()
531 PetscCall(MatSetOption(pch2opus->S, MAT_SYMMETRIC, PETSC_TRUE)); in PCH2OpusSetUpSampler_Private()
537 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCSetUp_H2OPUS()
538 Mat A = pc->useAmat ? pc->mat : pc->pmat; in PCSetUp_H2OPUS()
539 NormType norm = pch2opus->normtype; in PCSetUp_H2OPUS()
544 if (!pch2opus->T) { in PCSetUp_H2OPUS()
545 PetscInt M, N, m, n; in PCSetUp_H2OPUS() local
549 PetscCall(MatGetSize(A, &M, &N)); in PCSetUp_H2OPUS()
550 PetscCall(MatGetLocalSize(A, &m, &n)); in PCSetUp_H2OPUS()
551 PetscCall(MatCreateShell(PetscObjectComm((PetscObject)pc->pmat), m, n, M, N, pc, &pch2opus->T)); in PCSetUp_H2OPUS()
552 PetscCall(MatSetBlockSizesFromMats(pch2opus->T, A, A)); in PCSetUp_H2OPUS()
553 PetscCall(MatShellSetOperation(pch2opus->T, MATOP_MULT, (PetscErrorCodeFn *)MatMult_MAmI)); in PCSetUp_H2OPUS()
554 …PetscCall(MatShellSetOperation(pch2opus->T, MATOP_MULT_TRANSPOSE, (PetscErrorCodeFn *)MatMultTrans… in PCSetUp_H2OPUS()
555 PetscCall(MatShellSetOperation(pch2opus->T, MATOP_NORM, (PetscErrorCodeFn *)MatNorm_H2OPUS)); in PCSetUp_H2OPUS()
557 PetscCall(MatShellSetVecType(pch2opus->T, VECCUDA)); in PCSetUp_H2OPUS()
559 PetscCall(MatSetOptionsPrefix(pch2opus->T, prefix)); in PCSetUp_H2OPUS()
560 PetscCall(MatAppendOptionsPrefix(pch2opus->T, "pc_h2opus_est_")); in PCSetUp_H2OPUS()
562 PetscCall(MatDestroy(&pch2opus->A)); in PCSetUp_H2OPUS()
566 pch2opus->A = A; in PCSetUp_H2OPUS()
571 if (!pch2opus->sdim) PetscCall(PCH2OpusInferCoordinates_Private(pc)); in PCSetUp_H2OPUS()
572 … pch2opus->sdim, pch2opus->coords, PETSC_FALSE, pch2opus->eta, pch2opus->leafsize, pch2opus->max_r… in PCSetUp_H2OPUS()
574 PetscCall(MatSetOption(pch2opus->A, MAT_SYMMETRIC, PETSC_TRUE)); in PCSetUp_H2OPUS()
576 PetscCall(MatSetOptionsPrefix(pch2opus->A, prefix)); in PCSetUp_H2OPUS()
577 PetscCall(MatAppendOptionsPrefix(pch2opus->A, "pc_h2opus_init_")); in PCSetUp_H2OPUS()
578 PetscCall(MatSetFromOptions(pch2opus->A)); in PCSetUp_H2OPUS()
579 PetscCall(MatAssemblyBegin(pch2opus->A, MAT_FINAL_ASSEMBLY)); in PCSetUp_H2OPUS()
580 PetscCall(MatAssemblyEnd(pch2opus->A, MAT_FINAL_ASSEMBLY)); in PCSetUp_H2OPUS()
582 PetscCall(MatSetOption(pch2opus->A, MAT_SYMMETRIC, PETSC_TRUE)); in PCSetUp_H2OPUS()
585 PetscCall(MatBindToCPU(pch2opus->A, pch2opus->forcecpu)); in PCSetUp_H2OPUS()
588 pch2opus->boundtocpu = pch2opus->forcecpu ? PETSC_TRUE : pch2opus->A->boundtocpu; in PCSetUp_H2OPUS()
590 PetscCall(MatBindToCPU(pch2opus->T, pch2opus->boundtocpu)); in PCSetUp_H2OPUS()
591 if (pch2opus->M) { /* see if we can reuse M as initial guess */ in PCSetUp_H2OPUS()
594 PetscCall(MatBindToCPU(pch2opus->M, pch2opus->boundtocpu)); in PCSetUp_H2OPUS()
595 PetscCall(MatNorm(pch2opus->T, norm, &reuse)); in PCSetUp_H2OPUS()
596 if (reuse >= 1.0) PetscCall(MatDestroy(&pch2opus->M)); in PCSetUp_H2OPUS()
598 if (!pch2opus->M) { in PCSetUp_H2OPUS()
600 PetscCall(MatDuplicate(pch2opus->A, MAT_COPY_VALUES, &pch2opus->M)); in PCSetUp_H2OPUS()
602 PetscCall(MatSetOptionsPrefix(pch2opus->M, prefix)); in PCSetUp_H2OPUS()
603 PetscCall(MatAppendOptionsPrefix(pch2opus->M, "pc_h2opus_inv_")); in PCSetUp_H2OPUS()
604 PetscCall(MatSetFromOptions(pch2opus->M)); in PCSetUp_H2OPUS()
606 PetscCall(MatScale(pch2opus->M, pch2opus->s0)); in PCSetUp_H2OPUS()
608 /* A and M have the same h2 matrix nonzero structure, save on reordering routines */ in PCSetUp_H2OPUS()
609 PetscCall(MatH2OpusSetNativeMult(pch2opus->A, PETSC_TRUE)); in PCSetUp_H2OPUS()
610 PetscCall(MatH2OpusSetNativeMult(pch2opus->M, PETSC_TRUE)); in PCSetUp_H2OPUS()
611 …if (norm == NORM_1 || norm == NORM_2 || norm == NORM_INFINITY) PetscCall(MatNorm(pch2opus->T, norm… in PCSetUp_H2OPUS()
612 if (PetscIsInfOrNanReal(initerr)) pc->failedreason = PC_SETUP_ERROR; in PCSetUp_H2OPUS()
614 …if (pch2opus->monitor) PetscCall(PetscPrintf(PetscObjectComm((PetscObject)pc), "%" PetscInt_FMT ":… in PCSetUp_H2OPUS()
615 if (initerr > pch2opus->atol && !pc->failedreason) { in PCSetUp_H2OPUS()
619 for (i = 0; i < pch2opus->maxits; i++) { in PCSetUp_H2OPUS()
620 Mat M; in PCSetUp_H2OPUS() local
623 PetscCall(MatDuplicate(pch2opus->M, MAT_SHARE_NONZERO_PATTERN, &M)); in PCSetUp_H2OPUS()
624 PetscCall(MatGetOptionsPrefix(pch2opus->M, &prefix)); in PCSetUp_H2OPUS()
625 PetscCall(MatSetOptionsPrefix(M, prefix)); in PCSetUp_H2OPUS()
626 PetscCall(MatH2OpusSetSamplingMat(M, pch2opus->S, PETSC_DECIDE, PETSC_DECIDE)); in PCSetUp_H2OPUS()
627 PetscCall(MatSetFromOptions(M)); in PCSetUp_H2OPUS()
628 PetscCall(MatH2OpusSetNativeMult(M, PETSC_TRUE)); in PCSetUp_H2OPUS()
629 PetscCall(MatAssemblyBegin(M, MAT_FINAL_ASSEMBLY)); in PCSetUp_H2OPUS()
630 PetscCall(MatAssemblyEnd(M, MAT_FINAL_ASSEMBLY)); in PCSetUp_H2OPUS()
631 PetscCall(MatDestroy(&pch2opus->M)); in PCSetUp_H2OPUS()
632 pch2opus->M = M; in PCSetUp_H2OPUS()
633 …if (norm == NORM_1 || norm == NORM_2 || norm == NORM_INFINITY) PetscCall(MatNorm(pch2opus->T, norm… in PCSetUp_H2OPUS()
634 …if (pch2opus->monitor) PetscCall(PetscPrintf(PetscObjectComm((PetscObject)pc), "%" PetscInt_FMT ":… in PCSetUp_H2OPUS()
635 if (PetscIsInfOrNanReal(err)) pc->failedreason = PC_SETUP_ERROR; in PCSetUp_H2OPUS()
636 if (err < pch2opus->atol || err < pch2opus->rtol * initerr || pc->failedreason) break; in PCSetUp_H2OPUS()
640 PetscCall(MatH2OpusSetNativeMult(pch2opus->A, PETSC_FALSE)); in PCSetUp_H2OPUS()
641 PetscCall(MatH2OpusSetNativeMult(pch2opus->M, PETSC_FALSE)); in PCSetUp_H2OPUS()
642 PetscCall(VecDestroy(&pch2opus->wns[0])); in PCSetUp_H2OPUS()
643 PetscCall(VecDestroy(&pch2opus->wns[1])); in PCSetUp_H2OPUS()
644 PetscCall(VecDestroy(&pch2opus->wns[2])); in PCSetUp_H2OPUS()
645 PetscCall(VecDestroy(&pch2opus->wns[3])); in PCSetUp_H2OPUS()
646 PetscCall(MatDestroy(&pch2opus->wnsmat[0])); in PCSetUp_H2OPUS()
647 PetscCall(MatDestroy(&pch2opus->wnsmat[1])); in PCSetUp_H2OPUS()
648 PetscCall(MatDestroy(&pch2opus->wnsmat[2])); in PCSetUp_H2OPUS()
649 PetscCall(MatDestroy(&pch2opus->wnsmat[3])); in PCSetUp_H2OPUS()
655 PC_H2OPUS *pch2opus = (PC_H2OPUS *)pc->data; in PCView_H2OPUS()
661 if (pch2opus->A && pch2opus->A != pc->mat && pch2opus->A != pc->pmat) { in PCView_H2OPUS()
665 PetscCall(MatView(pch2opus->A, viewer)); in PCView_H2OPUS()
669 if (pch2opus->M) { in PCView_H2OPUS()
673 PetscCall(MatView(pch2opus->M, viewer)); in PCView_H2OPUS()
677 PetscCall(PetscViewerASCIIPrintf(viewer, "Initial scaling: %g\n", pch2opus->s0)); in PCView_H2OPUS()
683 …PCH2OPUS = "h2opus" - A preconditioner type for, `MATH2OPUS`, hierarchical matrices using the H2Op…
686 + -pc_type h2opus - pc type to "h2opus" during a call to `PCSetFromOptions()`
687 . -pc_h2opus_maxits - maximum number of iterations for Newton-Schultz
688 . -pc_h2opus_monitor - monitor Newton-Schultz convergence
689 . -pc_h2opus_atol - absolute tolerance
690 . -pc_h2opus_rtol - relative tolerance
691 . -pc_h2opus_norm_type - normtype
692 . -pc_h2opus_hyperorder - Hyper power order of sampling
693 . -pc_h2opus_leafsize - leaf size when constructed from kernel
694 . -pc_h2opus_eta - admissibility condition tolerance
695 . -pc_h2opus_maxrank - maximum rank when constructed from matvecs
696 . -pc_h2opus_samples - number of samples to be taken concurrently when constructing from matvecs
697 . -pc_h2opus_mrtol - relative tolerance for construction from sampling
698 - -pc_h2opus_forcecpu - force construction of preconditioner on CPU
703 M*/
711 pc->data = (void *)pch2opus; in PCCreate_H2OPUS()
713 pch2opus->atol = 1.e-2; in PCCreate_H2OPUS()
714 pch2opus->rtol = 1.e-6; in PCCreate_H2OPUS()
715 pch2opus->maxits = 50; in PCCreate_H2OPUS()
716 pch2opus->hyperorder = 1; /* defaults to basic Newton-Schultz */ in PCCreate_H2OPUS()
717 pch2opus->normtype = NORM_2; in PCCreate_H2OPUS()
720 pch2opus->eta = PETSC_DECIDE; in PCCreate_H2OPUS()
721 pch2opus->leafsize = PETSC_DECIDE; in PCCreate_H2OPUS()
722 pch2opus->max_rank = PETSC_DECIDE; in PCCreate_H2OPUS()
723 pch2opus->bs = PETSC_DECIDE; in PCCreate_H2OPUS()
724 pch2opus->mrtol = PETSC_DECIDE; in PCCreate_H2OPUS()
725 pch2opus->boundtocpu = PetscDefined(H2OPUS_USE_GPU) ? PETSC_FALSE : PETSC_TRUE; in PCCreate_H2OPUS()
726 pc->ops->destroy = PCDestroy_H2OPUS; in PCCreate_H2OPUS()
727 pc->ops->setup = PCSetUp_H2OPUS; in PCCreate_H2OPUS()
728 pc->ops->apply = PCApply_H2OPUS; in PCCreate_H2OPUS()
729 pc->ops->matapply = PCApplyMat_H2OPUS; in PCCreate_H2OPUS()
730 pc->ops->applytranspose = PCApplyTranspose_H2OPUS; in PCCreate_H2OPUS()
731 pc->ops->reset = PCReset_H2OPUS; in PCCreate_H2OPUS()
732 pc->ops->setfromoptions = PCSetFromOptions_H2OPUS; in PCCreate_H2OPUS()
733 pc->ops->view = PCView_H2OPUS; in PCCreate_H2OPUS()