Actual source code: lmvmutils.c
petsc-3.10.5 2019-03-28
1: #include <../src/ksp/ksp/utils/lmvm/lmvm.h>
3: /*@
4: MatLMVMUpdate - Adds (X-Xprev) and (F-Fprev) updates to an LMVM matrix.
5: The first time the function is called for an LMVM matrix, no update is
6: applied, but the given X and F vectors are stored for use as Xprev and
7: Fprev in the next update.
8:
9: If the user has provided another LMVM matrix in place of J0, the J0
10: matrix is also updated recursively.
12: Input Parameters:
13: + B - An LMVM-type matrix
14: . X - Solution vector
15: - F - Function vector
17: Level: intermediate
19: .seealso: MatLMVMReset(), MatLMVMAllocate()
20: @*/
21: PetscErrorCode MatLMVMUpdate(Mat B, Vec X, Vec F)
22: {
23: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
24: PetscErrorCode ierr;
25: PetscBool same;
31: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
32: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
33: if (!lmvm->allocated) {
34: MatLMVMAllocate(B, X, F);
35: } else {
36: VecCheckMatCompatible(B, X, 2, F, 3);
37: }
38: if (lmvm->J0) {
39: /* If the user provided an LMVM-type matrix as J0, then trigger its update as well */
40: PetscObjectBaseTypeCompare((PetscObject)lmvm->J0, MATLMVM, &same);
41: if (same) {
42: MatLMVMUpdate(lmvm->J0, X, F);
43: }
44: }
45: (*lmvm->ops->update)(B, X, F);
46: return(0);
47: }
49: /*------------------------------------------------------------*/
51: /*@
52: MatLMVMClearJ0 - Removes all definitions of J0 and reverts to
53: an identity matrix (scale = 1.0).
55: Input Parameters:
56: . B - An LMVM-type matrix
58: Level: advanced
60: .seealso: MatLMVMSetJ0()
61: @*/
62: PetscErrorCode MatLMVMClearJ0(Mat B)
63: {
64: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
65: PetscErrorCode ierr;
66: PetscBool same;
67: MPI_Comm comm = PetscObjectComm((PetscObject)B);
71: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
72: if (!same) SETERRQ(comm, PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
73: lmvm->user_pc = PETSC_FALSE;
74: lmvm->user_ksp = PETSC_FALSE;
75: lmvm->user_scale = PETSC_FALSE;
76: lmvm->J0scalar = 1.0;
77: VecDestroy(&lmvm->J0diag);
78: MatDestroy(&lmvm->J0);
79: PCDestroy(&lmvm->J0pc);
80: return(0);
81: }
83: /*------------------------------------------------------------*/
85: /*@
86: MatLMVMSetJ0Scale - Allows the user to define a scalar value
87: mu such that J0 = mu*I.
89: Input Parameters:
90: + B - An LMVM-type matrix
91: - scale - Scalar value mu that defines the initial Jacobian
93: Level: advanced
95: .seealso: MatLMVMSetDiagScale(), MatLMVMSetJ0()
96: @*/
97: PetscErrorCode MatLMVMSetJ0Scale(Mat B, PetscReal scale)
98: {
99: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
100: PetscErrorCode ierr;
101: PetscBool same;
102: MPI_Comm comm = PetscObjectComm((PetscObject)B);
106: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
107: if (!same) SETERRQ(comm, PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
108: if (!lmvm->square) SETERRQ(comm, PETSC_ERR_SUP, "Scaling is available only for square LMVM matrices");
109: MatLMVMClearJ0(B);
110: lmvm->J0scalar = scale;
111: lmvm->user_scale = PETSC_TRUE;
112: return(0);
113: }
115: /*------------------------------------------------------------*/
117: /*@
118: MatLMVMSetJ0Diag - Allows the user to define a vector
119: V such that J0 = diag(V).
121: Input Parameters:
122: + B - An LMVM-type matrix
123: - V - Vector that defines the diagonal of the initial Jacobian
125: Level: advanced
127: .seealso: MatLMVMSetScale(), MatLMVMSetJ0()
128: @*/
129: PetscErrorCode MatLMVMSetJ0Diag(Mat B, Vec V)
130: {
131: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
132: PetscErrorCode ierr;
133: PetscBool same;
134: MPI_Comm comm = PetscObjectComm((PetscObject)B);
139: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
140: if (!same) SETERRQ(comm, PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
141: if (!lmvm->allocated) SETERRQ(comm, PETSC_ERR_ORDER, "Matrix must be allocated before setting diagonal scaling");
142: if (!lmvm->square) SETERRQ(comm, PETSC_ERR_SUP, "Diagonal scaling is available only for square LMVM matrices");
143: VecCheckSameSize(V, 2, lmvm->Fprev, 3);
144: MatLMVMClearJ0(B);
145: if (!lmvm->J0diag) {
146: VecDuplicate(V, &lmvm->J0diag);
147: }
148: VecCopy(V, lmvm->J0diag);
149: lmvm->user_scale = PETSC_TRUE;
150: return(0);
151: }
153: /*------------------------------------------------------------*/
155: /*@
156: MatLMVMSetJ0 - Allows the user to define the initial
157: Jacobian matrix from which the LMVM approximation is
158: built up. Inverse of this initial Jacobian is applied
159: using an internal KSP solver, which defaults to GMRES.
160: This internal KSP solver has the "mat_lmvm_" option
161: prefix.
162:
163: Note that another LMVM matrix can be used in place of
164: J0, in which case updating the outer LMVM matrix will
165: also trigger the update for the inner LMVM matrix. This
166: is useful in cases where a full-memory diagonal approximation
167: such as MATLMVMDIAGBRDN is used in place of J0.
169: Input Parameters:
170: + B - An LMVM-type matrix
171: - J0 - The initial Jacobian matrix
173: Level: advanced
175: .seealso: MatLMVMSetJ0PC(), MatLMVMSetJ0KSP()
176: @*/
177: PetscErrorCode MatLMVMSetJ0(Mat B, Mat J0)
178: {
179: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
180: PetscErrorCode ierr;
181: PetscBool same;
182: MPI_Comm comm = PetscObjectComm((PetscObject)B);
187: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
188: if (!same) SETERRQ(comm, PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
189: if (!J0->assembled) SETERRQ(comm, PETSC_ERR_ARG_WRONGSTATE, "J0 is not assembled.");
190: if (B->symmetric && (!J0->symmetric)) SETERRQ(comm, PETSC_ERR_ARG_INCOMP, "J0 and J0pre must be symmetric when B is symmetric");
191: if (lmvm->allocated) {
192: MatCheckSameSize(B, 1, J0, 2);
193: }
194: MatLMVMClearJ0(B);
195: MatDestroy(&lmvm->J0);
196: PetscObjectReference((PetscObject)J0);
197: lmvm->J0 = J0;
198: PetscObjectBaseTypeCompare((PetscObject)lmvm->J0, MATLMVM, &same);
199: if (!same && lmvm->square) {
200: KSPSetOperators(lmvm->J0ksp, lmvm->J0, lmvm->J0);
201: }
202: return(0);
203: }
205: /*------------------------------------------------------------*/
207: /*@
208: MatLMVMSetJ0PC - Allows the user to define a PC object that
209: acts as the initial inverse-Jacobian matrix. This PC should
210: already contain all the operators necessary for its application.
211: The LMVM matrix only calls PCApply() without changing any other
212: options.
214: Input Parameters:
215: + B - An LMVM-type matrix
216: - J0pc - PC object where PCApply defines an inverse application for J0
218: Level: advanced
220: .seealso: MatLMVMGetJ0PC()
221: @*/
222: PetscErrorCode MatLMVMSetJ0PC(Mat B, PC J0pc)
223: {
224: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
225: PetscErrorCode ierr;
226: PetscBool same;
227: MPI_Comm comm = PetscObjectComm((PetscObject)B);
228: Mat J0, J0pre;
233: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
234: if (!same) SETERRQ(comm, PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
235: if (!lmvm->square) SETERRQ(comm, PETSC_ERR_SUP, "Inverse J0 can be defined only for square LMVM matrices");
236: PCGetOperators(J0pc, &J0, &J0pre);
237: if (B->symmetric && (!J0->symmetric || !J0pre->symmetric)) SETERRQ(comm, PETSC_ERR_ARG_INCOMP, "J0 and J0pre must be symmetric when B is symmetric");
238: if (lmvm->allocated) {
239: MatCheckSameSize(B, 1, J0, 2);
240: MatCheckSameSize(B, 1, J0pre, 3);
241: }
242: MatDestroy(&J0);
243: MatDestroy(&J0pre);
244: MatLMVMClearJ0(B);
245: PetscObjectReference((PetscObject)J0pc);
246: lmvm->J0pc = J0pc;
247: lmvm->user_pc = PETSC_TRUE;
248: return(0);
249: }
251: /*------------------------------------------------------------*/
253: /*@
254: MatLMVMSetJ0KSP - Allows the user to provide a pre-configured
255: KSP solver for the initial inverse-Jacobian approximation.
256: This KSP solver should already contain all the operators
257: necessary to perform the inversion. The LMVM matrix only
258: calls KSPSolve() without changing any other options.
260: Input Parameters:
261: + B - An LMVM-type matrix
262: - J0ksp - KSP solver for the initial inverse-Jacobian application
264: Level: advanced
266: .seealso: MatLMVMGetJ0KSP()
267: @*/
268: PetscErrorCode MatLMVMSetJ0KSP(Mat B, KSP J0ksp)
269: {
270: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
271: PetscErrorCode ierr;
272: PetscBool same;
273: MPI_Comm comm = PetscObjectComm((PetscObject)B);
274: Mat J0, J0pre;
279: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
280: if (!same) SETERRQ(comm, PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
281: if (!lmvm->square) SETERRQ(comm, PETSC_ERR_SUP, "Inverse J0 can be defined only for square LMVM matrices");
282: KSPGetOperators(J0ksp, &J0, &J0pre);
283: if (B->symmetric && (!J0->symmetric || !J0pre->symmetric)) SETERRQ(comm, PETSC_ERR_ARG_INCOMP, "J0 and J0pre must be symmetric when B is symmetric");
284: if (lmvm->allocated) {
285: MatCheckSameSize(B, 1, J0, 2);
286: MatCheckSameSize(B, 1, J0pre, 3);
287: }
288: MatLMVMClearJ0(B);
289: KSPDestroy(&lmvm->J0ksp);
290: PetscObjectReference((PetscObject)J0ksp);
291: lmvm->J0ksp = J0ksp;
292: lmvm->user_ksp = PETSC_TRUE;
293: return(0);
294: }
296: /*------------------------------------------------------------*/
298: /*@
299: MatLMVMGetJ0 - Returns a pointer to the internal J0 matrix.
301: Input Parameters:
302: . B - An LMVM-type matrix
304: Output Parameter:
305: . J0 - Mat object for defining the initial Jacobian
307: Level: advanced
309: .seealso: MatLMVMSetJ0()
310: @*/
311: PetscErrorCode MatLMVMGetJ0(Mat B, Mat *J0)
312: {
313: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
314: PetscErrorCode ierr;
315: PetscBool same;
319: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
320: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
321: *J0 = lmvm->J0;
322: return(0);
323: }
325: /*------------------------------------------------------------*/
327: /*@
328: MatLMVMGetJ0PC - Returns a pointer to the internal PC object
329: associated with the initial Jacobian.
331: Input Parameters:
332: . B - An LMVM-type matrix
334: Output Parameter:
335: . J0pc - PC object for defining the initial inverse-Jacobian
337: Level: advanced
339: .seealso: MatLMVMSetJ0PC()
340: @*/
341: PetscErrorCode MatLMVMGetJ0PC(Mat B, PC *J0pc)
342: {
343: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
344: PetscErrorCode ierr;
345: PetscBool same;
349: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
350: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
351: if (lmvm->J0pc) {
352: *J0pc = lmvm->J0pc;
353: } else {
354: KSPGetPC(lmvm->J0ksp, J0pc);
355: }
356: return(0);
357: }
359: /*------------------------------------------------------------*/
361: /*@
362: MatLMVMGetJ0KSP - Returns a pointer to the internal KSP solver
363: associated with the initial Jacobian.
365: Input Parameters:
366: . B - An LMVM-type matrix
368: Output Parameter:
369: . J0ksp - KSP solver for defining the initial inverse-Jacobian
371: Level: advanced
373: .seealso: MatLMVMSetJ0KSP()
374: @*/
375: PetscErrorCode MatLMVMGetJ0KSP(Mat B, KSP *J0ksp)
376: {
377: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
378: PetscErrorCode ierr;
379: PetscBool same;
383: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
384: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
385: *J0ksp = lmvm->J0ksp;
386: return(0);
387: }
389: /*------------------------------------------------------------*/
391: /*@
392: MatLMVMApplyJ0Fwd - Applies an approximation of the forward
393: matrix-vector product with the initial Jacobian.
395: Input Parameters:
396: + B - An LMVM-type matrix
397: - X - vector to multiply with J0
399: Output Parameter:
400: . Y - resulting vector for the operation
402: Level: advanced
404: .seealso: MatLMVMSetJ0(), MatLMVMSetJ0Scale(), MatLMVMSetJ0ScaleDiag(),
405: MatLMVMSetJ0PC(), MatLMVMSetJ0KSP(), MatLMVMApplyJ0Inv()
406: @*/
407: PetscErrorCode MatLMVMApplyJ0Fwd(Mat B, Vec X, Vec Y)
408: {
409: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
410: PetscErrorCode ierr;
411: PetscBool same, hasMult;
412: MPI_Comm comm = PetscObjectComm((PetscObject)B);
413: Mat Amat, Pmat;
419: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
420: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
421: if (!lmvm->allocated) SETERRQ(comm, PETSC_ERR_ORDER, "LMVM matrix must be allocated first");
422: VecCheckMatCompatible(B, X, 2, Y, 3);
423: if (lmvm->user_pc || lmvm->user_ksp || lmvm->J0) {
424: /* User may have defined a PC or KSP for J0^{-1} so let's try to use its operators. */
425: if (lmvm->user_pc){
426: PCGetOperators(lmvm->J0pc, &Amat, &Pmat);
427: } else if (lmvm->user_ksp) {
428: KSPGetOperators(lmvm->J0ksp, &Amat, &Pmat);
429: } else {
430: Amat = lmvm->J0;
431: }
432: MatHasOperation(Amat, MATOP_MULT, &hasMult);
433: if (hasMult) {
434: /* product is available, use it */
435: MatMult(Amat, X, Y);
436: } else {
437: /* there's no product, so treat J0 as identity */
438: VecCopy(X, Y);
439: }
440: } else if (lmvm->user_scale) {
441: if (lmvm->J0diag) {
442: /* User has defined a diagonal vector for J0 */
443: VecPointwiseMult(X, lmvm->J0diag, Y);
444: } else {
445: /* User has defined a scalar value for J0 */
446: VecCopy(X, Y);
447: VecScale(Y, lmvm->J0scalar);
448: }
449: } else {
450: /* There is no J0 representation so just apply an identity matrix */
451: VecCopy(X, Y);
452: }
453: return(0);
454: }
456: /*------------------------------------------------------------*/
458: /*@
459: MatLMVMApplyJ0Inv - Applies some estimation of the initial Jacobian
460: inverse to the given vector. The specific form of the application
461: depends on whether the user provided a scaling factor, a J0 matrix,
462: a J0 PC, or a J0 KSP object. If no form of the initial Jacobian is
463: provided, the function simply does an identity matrix application
464: (vector copy).
466: Input Parameters:
467: + B - An LMVM-type matrix
468: - X - vector to "multiply" with J0^{-1}
470: Output Parameter:
471: . Y - resulting vector for the operation
473: Level: advanced
475: .seealso: MatLMVMSetJ0(), MatLMVMSetJ0Scale(), MatLMVMSetJ0ScaleDiag(),
476: MatLMVMSetJ0PC(), MatLMVMSetJ0KSP(), MatLMVMApplyJ0Fwd()
477: @*/
478: PetscErrorCode MatLMVMApplyJ0Inv(Mat B, Vec X, Vec Y)
479: {
480: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
481: PetscErrorCode ierr;
482: PetscBool same, hasSolve;
483: MPI_Comm comm = PetscObjectComm((PetscObject)B);
489: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
490: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
491: if (!lmvm->allocated) SETERRQ(comm, PETSC_ERR_ORDER, "LMVM matrix must be allocated first");
492: VecCheckMatCompatible(B, X, 2, Y, 3);
493: /* Invert the initial Jacobian onto q (or apply scaling) */
494: if (lmvm->user_pc) {
495: /* User has defined a J0 inverse so we can directly apply it as a preconditioner */
496: PCApply(lmvm->J0pc, X, Y);
497: } else if (lmvm->user_ksp) {
498: /* User has defined a J0 or a custom KSP so just perform a solution */
499: KSPSolve(lmvm->J0ksp, X, Y);
500: } else if (lmvm->J0) {
501: MatHasOperation(lmvm->J0, MATOP_SOLVE, &hasSolve);
502: if (hasSolve) {
503: MatSolve(lmvm->J0, X, Y);
504: } else {
505: KSPSolve(lmvm->J0ksp, X, Y);
506: }
507: } else if (lmvm->user_scale) {
508: if (lmvm->J0diag) {
509: VecPointwiseDivide(X, Y, lmvm->J0diag);
510: } else {
511: VecCopy(X, Y);
512: VecScale(Y, 1.0/lmvm->J0scalar);
513: }
514: } else {
515: /* There is no J0 representation so just apply an identity matrix */
516: VecCopy(X, Y);
517: }
518: return(0);
519: }
521: /*------------------------------------------------------------*/
523: /*@
524: MatLMVMIsAllocated - Returns a boolean flag that shows whether
525: the necessary data structures for the underlying matrix is allocated.
527: Input Parameters:
528: . B - An LMVM-type matrix
530: Output Parameter:
531: . flg - PETSC_TRUE if allocated, PETSC_FALSE otherwise
533: Level: intermediate
535: .seealso: MatLMVMAllocate(), MatLMVMReset()
536: @*/
538: PetscErrorCode MatLMVMIsAllocated(Mat B, PetscBool *flg)
539: {
540: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
541: PetscErrorCode ierr;
542: PetscBool same;
543:
546: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
547: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
548: *flg = PETSC_FALSE;
549: if (lmvm->allocated && B->preallocated && B->assembled) *flg = PETSC_TRUE;
550: return(0);
551: }
553: /*------------------------------------------------------------*/
555: /*@
556: MatLMVMAllocate - Produces all necessary common memory for
557: LMVM approximations based on the solution and function vectors
558: provided. If MatSetSizes() and MatSetUp() have not been called
559: before MatLMVMAllocate(), the allocation will read sizes from
560: the provided vectors and update the matrix.
562: Input Parameters:
563: + B - An LMVM-type matrix
564: . X - Solution vector
565: - F - Function vector
567: Level: intermediate
569: .seealso: MatLMVMReset(), MatLMVMUpdate()
570: @*/
571: PetscErrorCode MatLMVMAllocate(Mat B, Vec X, Vec F)
572: {
573: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
574: PetscErrorCode ierr;
575: PetscBool same;
581: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
582: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
583: (*lmvm->ops->allocate)(B, X, F);
584: if (lmvm->J0) {
585: PetscObjectBaseTypeCompare((PetscObject)lmvm->J0, MATLMVM, &same);
586: if (same) {
587: MatLMVMAllocate(lmvm->J0, X, F);
588: }
589: }
590: return(0);
591: }
593: /*------------------------------------------------------------*/
595: /*@
596: MatLMVMResetShift - Zero the shift factor.
598: Input Parameters:
599: . B - An LMVM-type matrix
601: Level: intermediate
603: .seealso: MatLMVMAllocate(), MatLMVMUpdate()
604: @*/
605: PetscErrorCode MatLMVMResetShift(Mat B)
606: {
607: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
608: PetscErrorCode ierr;
609: PetscBool same;
613: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
614: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
615: lmvm->shift = 0.0;
616: return(0);
617: }
619: /*------------------------------------------------------------*/
621: /*@
622: MatLMVMReset - Flushes all of the accumulated updates out of
623: the LMVM approximation. In practice, this will not actually
624: destroy the data associated with the updates. It simply resets
625: counters, which leads to existing data being overwritten, and
626: MatSolve() being applied as if there are no updates. A boolean
627: flag is available to force destruction of the update vectors.
628:
629: If the user has provided another LMVM matrix as J0, the J0
630: matrix is also reset in this function.
632: Input Parameters:
633: + B - An LMVM-type matrix
634: - destructive - flag for enabling destruction of data structures
636: Level: intermediate
638: .seealso: MatLMVMAllocate(), MatLMVMUpdate()
639: @*/
640: PetscErrorCode MatLMVMReset(Mat B, PetscBool destructive)
641: {
642: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
643: PetscErrorCode ierr;
644: PetscBool same;
648: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
649: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
650: (*lmvm->ops->reset)(B, destructive);
651: if (lmvm->J0) {
652: PetscObjectBaseTypeCompare((PetscObject)lmvm->J0, MATLMVM, &same);
653: if (same) {
654: MatLMVMReset(lmvm->J0, destructive);
655: }
656: }
657: return(0);
658: }
660: /*------------------------------------------------------------*/
662: /*@
663: MatLMVMGetUpdateCount - Returns the number of accepted updates.
664: This number may be greater than the total number of update vectors
665: stored in the matrix. The counters are reset when MatLMVMReset()
666: is called.
668: Input Parameters:
669: . B - An LMVM-type matrix
671: Output Parameter:
672: . nupdates - number of accepted updates
674: Level: intermediate
676: .seealso: MatLMVMGetRejectCount(), MatLMVMReset()
677: @*/
678: PetscErrorCode MatLMVMGetUpdateCount(Mat B, PetscInt *nupdates)
679: {
680: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
681: PetscErrorCode ierr;
682: PetscBool same;
683:
686: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
687: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
688: *nupdates = lmvm->nupdates;
689: return(0);
690: }
692: /*------------------------------------------------------------*/
694: /*@
695: MatLMVMGetRejectCount - Returns the number of rejected updates.
696: The counters are reset when MatLMVMReset() is called.
698: Input Parameters:
699: . B - An LMVM-type matrix
701: Output Parameter:
702: . nrejects - number of rejected updates
704: Level: intermediate
706: .seealso: MatLMVMGetRejectCount(), MatLMVMReset()
707: @*/
708: PetscErrorCode MatLMVMGetRejectCount(Mat B, PetscInt *nrejects)
709: {
710: Mat_LMVM *lmvm = (Mat_LMVM*)B->data;
711: PetscErrorCode ierr;
712: PetscBool same;
713:
716: PetscObjectBaseTypeCompare((PetscObject)B, MATLMVM, &same);
717: if (!same) SETERRQ(PetscObjectComm((PetscObject)B), PETSC_ERR_ARG_WRONG, "Matrix must be an LMVM-type.");
718: *nrejects = lmvm->nrejects;
719: return(0);
720: }