Actual source code: dsops.c
1: /*
2: DS operations: DSSolve(), DSVectors(), etc
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2012, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
9:
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include <slepc-private/dsimpl.h> /*I "slepcds.h" I*/
28: /*@
29: DSGetLeadingDimension - Returns the leading dimension of the allocated
30: matrices.
32: Not Collective
34: Input Parameter:
35: . ds - the direct solver context
37: Output Parameter:
38: . ld - leading dimension (maximum allowed dimension for the matrices)
40: Level: advanced
42: .seealso: DSAllocate(), DSSetDimensions()
43: @*/
44: PetscErrorCode DSGetLeadingDimension(DS ds,PetscInt *ld)
45: {
48: if (ld) *ld = ds->ld;
49: return(0);
50: }
54: /*@
55: DSSetState - Change the state of the DS object.
57: Logically Collective on DS
59: Input Parameters:
60: + ds - the direct solver context
61: - state - the new state
63: Notes:
64: The state indicates that the dense system is in an initial state (raw),
65: in an intermediate state (such as tridiagonal, Hessenberg or
66: Hessenberg-triangular), in a condensed state (such as diagonal, Schur or
67: generalized Schur), or in a truncated state.
69: This function is normally used to return to the raw state when the
70: condensed structure is destroyed.
72: Level: advanced
74: .seealso: DSGetState()
75: @*/
76: PetscErrorCode DSSetState(DS ds,DSStateType state)
77: {
83: switch (state) {
84: case DS_STATE_RAW:
85: case DS_STATE_INTERMEDIATE:
86: case DS_STATE_CONDENSED:
87: case DS_STATE_TRUNCATED:
88: if (ds->state<state) { PetscInfo(ds,"DS state has been increased\n"); }
89: ds->state = state;
90: break;
91: default:
92: SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_WRONG,"Wrong state");
93: }
94: return(0);
95: }
99: /*@
100: DSGetState - Returns the current state.
102: Not Collective
104: Input Parameter:
105: . ds - the direct solver context
107: Output Parameter:
108: . state - current state
110: Level: advanced
112: .seealso: DSSetState()
113: @*/
114: PetscErrorCode DSGetState(DS ds,DSStateType *state)
115: {
118: if (state) *state = ds->state;
119: return(0);
120: }
124: /*@
125: DSSetDimensions - Resize the matrices in the DS object.
127: Logically Collective on DS
129: Input Parameters:
130: + ds - the direct solver context
131: . n - the new size
132: . m - the new column size (only for SVD)
133: . l - number of locked (inactive) leading columns
134: - k - intermediate dimension (e.g., position of arrow)
136: Notes:
137: The internal arrays are not reallocated.
139: The value m is not used except in the case of DSSVD, pass PETSC_IGNORE
140: otherwise.
142: Level: intermediate
144: .seealso: DSGetDimensions(), DSAllocate()
145: @*/
146: PetscErrorCode DSSetDimensions(DS ds,PetscInt n,PetscInt m,PetscInt l,PetscInt k)
147: {
154: if (!ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must call DSAllocate() first");
155: if (n!=PETSC_IGNORE) {
156: if (n==PETSC_DECIDE || n==PETSC_DEFAULT) {
157: ds->n = ds->ld;
158: } else {
159: if (n<1 || n>ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of n. Must be between 1 and ld");
160: if (ds->extrarow && n+1>ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_OUTOFRANGE,"A value of n equal to ld leaves no room for extra row");
161: ds->n = n;
162: }
163: }
164: if (m!=PETSC_IGNORE) {
165: if (m==PETSC_DECIDE || m==PETSC_DEFAULT) {
166: ds->m = ds->ld;
167: } else {
168: if (m<1 || m>ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of m. Must be between 1 and ld");
169: ds->m = m;
170: }
171: }
172: if (l==PETSC_DECIDE || l==PETSC_DEFAULT) {
173: ds->l = 0;
174: } else {
175: if (l<0 || l>ds->n) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of l. Must be between 0 and n");
176: ds->l = l;
177: }
178: if (k==PETSC_DECIDE || k==PETSC_DEFAULT) {
179: ds->k = ds->n/2;
180: } else {
181: if (k<0 || k>ds->n) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of k. Must be between 0 and n");
182: ds->k = k;
183: }
184: return(0);
185: }
189: /*@
190: DSGetDimensions - Returns the current dimensions.
192: Not Collective
194: Input Parameter:
195: . ds - the direct solver context
197: Output Parameter:
198: + n - the current size
199: . m - the current column size (only for SVD)
200: . l - number of locked (inactive) leading columns
201: - k - intermediate dimension (e.g., position of arrow)
203: Level: intermediate
205: .seealso: DSSetDimensions()
206: @*/
207: PetscErrorCode DSGetDimensions(DS ds,PetscInt *n,PetscInt *m,PetscInt *l,PetscInt *k)
208: {
211: if (n) *n = ds->n;
212: if (m) *m = ds->m;
213: if (l) *l = ds->l;
214: if (k) *k = ds->k;
215: return(0);
216: }
220: /*@
221: DSTruncate - Truncates the system represented in the DS object.
223: Logically Collective on DS
225: Input Parameters:
226: + ds - the direct solver context
227: - n - the new size
229: Note:
230: The new size is set to n. In cases where the extra row is meaningful,
231: the first n elements are kept as the extra row for the new system.
233: Level: advanced
235: .seealso: DSSetDimensions(), DSSetExtraRow()
236: @*/
237: PetscErrorCode DSTruncate(DS ds,PetscInt n)
238: {
244: if (!ds->ops->truncate) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_SUP,"DS type %s",((PetscObject)ds)->type_name);
245: if (ds->state<DS_STATE_CONDENSED) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must call DSSolve() first");
246: if (n<ds->l || n>ds->n) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of n. Must be between l and n");
247: PetscLogEventBegin(DS_Other,ds,0,0,0);
248: PetscFPTrapPush(PETSC_FP_TRAP_OFF);
249: (*ds->ops->truncate)(ds,n);
250: PetscFPTrapPop();
251: PetscLogEventEnd(DS_Other,ds,0,0,0);
252: ds->state = DS_STATE_TRUNCATED;
253: PetscObjectStateIncrease((PetscObject)ds);
254: return(0);
255: }
259: /*@C
260: DSGetArray - Returns a pointer to one of the internal arrays used to
261: represent matrices. You MUST call DSRestoreArray() when you no longer
262: need to access the array.
264: Not Collective
266: Input Parameters:
267: + ds - the direct solver context
268: - m - the requested matrix
270: Output Parameter:
271: . a - pointer to the values
273: Level: advanced
275: .seealso: DSRestoreArray(), DSGetArrayReal()
276: @*/
277: PetscErrorCode DSGetArray(DS ds,DSMatType m,PetscScalar *a[])
278: {
282: if (m<0 || m>=DS_NUM_MAT) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_WRONG,"Invalid matrix");
283: if (!ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must call DSAllocate() first");
284: if (!ds->mat[m]) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_WRONGSTATE,"Requested matrix was not created in this DS");
285: *a = ds->mat[m];
286: CHKMEMQ;
287: return(0);
288: }
292: /*@C
293: DSRestoreArray - Restores the matrix after DSGetArray() was called.
295: Not Collective
297: Input Parameters:
298: + ds - the direct solver context
299: . m - the requested matrix
300: - a - pointer to the values
302: Level: advanced
304: .seealso: DSGetArray(), DSGetArrayReal()
305: @*/
306: PetscErrorCode DSRestoreArray(DS ds,DSMatType m,PetscScalar *a[])
307: {
313: if (m<0 || m>=DS_NUM_MAT) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_WRONG,"Invalid matrix");
314: CHKMEMQ;
315: *a = 0;
316: PetscObjectStateIncrease((PetscObject)ds);
317: return(0);
318: }
322: /*@C
323: DSGetArrayReal - Returns a pointer to one of the internal arrays used to
324: represent real matrices. You MUST call DSRestoreArray() when you no longer
325: need to access the array.
327: Not Collective
329: Input Parameters:
330: + ds - the direct solver context
331: - m - the requested matrix
333: Output Parameter:
334: . a - pointer to the values
336: Level: advanced
338: .seealso: DSRestoreArrayReal(), DSGetArray()
339: @*/
340: PetscErrorCode DSGetArrayReal(DS ds,DSMatType m,PetscReal *a[])
341: {
345: if (m<0 || m>=DS_NUM_MAT) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_WRONG,"Invalid matrix");
346: if (!ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must call DSAllocate() first");
347: if (!ds->rmat[m]) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_WRONGSTATE,"Requested matrix was not created in this DS");
348: *a = ds->rmat[m];
349: CHKMEMQ;
350: return(0);
351: }
355: /*@C
356: DSRestoreArrayReal - Restores the matrix after DSGetArrayReal() was called.
358: Not Collective
360: Input Parameters:
361: + ds - the direct solver context
362: . m - the requested matrix
363: - a - pointer to the values
365: Level: advanced
367: .seealso: DSGetArrayReal(), DSGetArray()
368: @*/
369: PetscErrorCode DSRestoreArrayReal(DS ds,DSMatType m,PetscReal *a[])
370: {
376: if (m<0 || m>=DS_NUM_MAT) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_WRONG,"Invalid matrix");
377: CHKMEMQ;
378: *a = 0;
379: PetscObjectStateIncrease((PetscObject)ds);
380: return(0);
381: }
385: /*@
386: DSSolve - Solves the problem.
388: Logically Collective on DS
390: Input Parameters:
391: + ds - the direct solver context
392: . eigr - array to store the computed eigenvalues (real part)
393: - eigi - array to store the computed eigenvalues (imaginary part)
395: Note:
396: This call brings the dense system to condensed form. No ordering
397: of the eigenvalues is enforced (for this, call DSSort() afterwards).
399: Level: intermediate
401: .seealso: DSSort()
402: @*/
403: PetscErrorCode DSSolve(DS ds,PetscScalar *eigr,PetscScalar *eigi)
404: {
410: if (!ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must call DSAllocate() first");
411: if (!ds->ops->solve) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_SUP,"DS type %s",((PetscObject)ds)->type_name);
412: if (ds->state>=DS_STATE_CONDENSED) return(0);
413: PetscLogEventBegin(DS_Solve,ds,0,0,0);
414: PetscFPTrapPush(PETSC_FP_TRAP_OFF);
415: if (!ds->ops->solve[ds->method]) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_OUTOFRANGE,"The specified method number does not exist for this DS");
416: (*ds->ops->solve[ds->method])(ds,eigr,eigi);
417: PetscFPTrapPop();
418: PetscLogEventEnd(DS_Solve,ds,0,0,0);
419: ds->state = DS_STATE_CONDENSED;
420: PetscObjectStateIncrease((PetscObject)ds);
421: return(0);
422: }
426: /*@
427: DSSort - Sorts the result of DSSolve() according to a given sorting
428: criterion.
430: Logically Collective on DS
432: Input Parameters:
433: + ds - the direct solver context
434: . eigr - array containing the computed eigenvalues (real part)
435: . eigi - array containing the computed eigenvalues (imaginary part)
436: . rr - (optional) array containing auxiliary values (real part)
437: - ri - (optional) array containing auxiliary values (imaginary part)
439: Input/Output Parameter:
440: . k - (optional) number of elements in the leading group
442: Notes:
443: This routine sorts the arrays provided in eigr and eigi, and also
444: sorts the dense system stored inside ds (assumed to be in condensed form).
445: The sorting criterion is specified with DSSetEigenvalueComparison().
447: If arrays rr and ri are provided, then a (partial) reordering based on these
448: values rather than on the eigenvalues is performed. In symmetric problems
449: a total order is obtained (parameter k is ignored), but otherwise the result
450: is sorted only partially. In this latter case, it is only guaranteed that
451: all the first k elements satisfy the comparison with any of the last n-k
452: elements. The output value of parameter k is the final number of elements in
453: the first set.
455: Level: intermediate
457: .seealso: DSSolve(), DSSetEigenvalueComparison()
458: @*/
459: PetscErrorCode DSSort(DS ds,PetscScalar *eigr,PetscScalar *eigi,PetscScalar *rr,PetscScalar *ri,PetscInt *k)
460: {
467: if (ds->state<DS_STATE_CONDENSED) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must call DSSolve() first");
468: if (ds->state==DS_STATE_TRUNCATED) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Cannot sort a truncated DS");
469: if (!ds->ops->sort) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_SUP,"DS type %s",((PetscObject)ds)->type_name);
470: if (!ds->comp_fun) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must provide a sorting criterion with DSSetEigenvalueComparison() first");
471: if (k && !rr) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_WRONG,"Argument k can only be used together with rr");
472: PetscLogEventBegin(DS_Other,ds,0,0,0);
473: PetscFPTrapPush(PETSC_FP_TRAP_OFF);
474: (*ds->ops->sort)(ds,eigr,eigi,rr,ri,k);
475: PetscFPTrapPop();
476: PetscLogEventEnd(DS_Other,ds,0,0,0);
477: PetscObjectStateIncrease((PetscObject)ds);
478: return(0);
479: }
483: /*@
484: DSVectors - Compute vectors associated to the dense system such
485: as eigenvectors.
487: Logically Collective on DS
489: Input Parameters:
490: + ds - the direct solver context
491: - mat - the matrix, used to indicate which vectors are required
493: Input/Output Parameter:
494: - j - (optional) index of vector to be computed
496: Output Parameter:
497: . rnorm - (optional) computed residual norm
499: Notes:
500: Allowed values for mat are DS_MAT_X, DS_MAT_Y, DS_MAT_U and DS_MAT_VT, to
501: compute right or left eigenvectors, or left or right singular vectors,
502: respectively.
504: If PETSC_NULL is passed in argument j then all vectors are computed,
505: otherwise j indicates which vector must be computed. In real non-symmetric
506: problems, on exit the index j will be incremented when a complex conjugate
507: pair is found.
509: This function can be invoked after the dense problem has been solved,
510: to get the residual norm estimate of the associated Ritz pair. In that
511: case, the relevant information is returned in rnorm.
513: For computing eigenvectors, LAPACK's _trevc is used so the matrix must
514: be in (quasi-)triangular form, or call DSSolve() first.
516: Level: intermediate
518: .seealso: DSSolve()
519: @*/
520: PetscErrorCode DSVectors(DS ds,DSMatType mat,PetscInt *j,PetscReal *rnorm)
521: {
527: if (!ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must call DSAllocate() first");
528: if (!ds->ops->vectors) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_SUP,"DS type %s",((PetscObject)ds)->type_name);
529: if (rnorm && !j) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must give a value of j");
530: if (!ds->mat[mat]) { DSAllocateMat_Private(ds,mat); }
531: PetscLogEventBegin(DS_Vectors,ds,0,0,0);
532: PetscFPTrapPush(PETSC_FP_TRAP_OFF);
533: (*ds->ops->vectors)(ds,mat,j,rnorm);
534: PetscFPTrapPop();
535: PetscLogEventEnd(DS_Vectors,ds,0,0,0);
536: PetscObjectStateIncrease((PetscObject)ds);
537: return(0);
538: }
543: /*@
544: DSNormalize - Normalize a column or all the columns of a matrix. Considers
545: the case when the columns represent the real and the imaginary part of a vector.
547: Logically Collective on DS
549: Input Parameter:
550: + ds - the direct solver context
551: . mat - the matrix to be modified
552: - col - the column to normalize or -1 to normalize all of them
554: Notes:
555: The columns are normalized with respect to the 2-norm.
557: If col and col+1 (or col-1 and col) represent the real and the imaginary
558: part of a vector, both columns are scaled.
560: Level: advanced
561: @*/
562: PetscErrorCode DSNormalize(DS ds,DSMatType mat,PetscInt col)
563: {
570: if (ds->state<DS_STATE_CONDENSED) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must call DSSolve() first");
571: if (!ds->ops->normalize) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_SUP,"DS type %s",((PetscObject)ds)->type_name);
572: if (col<-1) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_OUTOFRANGE,"col should be at least minus one");
573: PetscLogEventBegin(DS_Other,ds,0,0,0);
574: PetscFPTrapPush(PETSC_FP_TRAP_OFF);
575: (*ds->ops->normalize)(ds,mat,col);
576: PetscFPTrapPop();
577: PetscLogEventEnd(DS_Other,ds,0,0,0);
578: return(0);
579: }
583: /*@C
584: DSUpdateExtraRow - Performs all necessary operations so that the extra
585: row gets up-to-date after a call to DSSolve().
587: Not Collective
589: Input Parameters:
590: . ds - the direct solver context
592: Level: advanced
594: .seealso: DSSolve(), DSSetExtraRow()
595: @*/
596: PetscErrorCode DSUpdateExtraRow(DS ds)
597: {
602: if (!ds->ops->update) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_SUP,"DS type %s",((PetscObject)ds)->type_name);
603: if (!ds->extrarow) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ARG_WRONGSTATE,"Should have called DSSetExtraRow");
604: if (!ds->ld) SETERRQ(((PetscObject)ds)->comm,PETSC_ERR_ORDER,"Must call DSAllocate() first");
605: PetscLogEventBegin(DS_Other,ds,0,0,0);
606: PetscFPTrapPush(PETSC_FP_TRAP_OFF);
607: (*ds->ops->update)(ds);
608: PetscFPTrapPop();
609: PetscLogEventEnd(DS_Other,ds,0,0,0);
610: return(0);
611: }
615: /*@C
616: DSCond - Compute the inf-norm condition number of the first matrix
617: as cond(A) = norm(A)*norm(inv(A)).
619: Not Collective
621: Input Parameters:
622: + ds - the direct solver context
623: - cond - the computed condition number
625: Level: advanced
627: .seealso: DSSolve()
628: @*/
629: PetscErrorCode DSCond(DS ds,PetscReal *cond)
630: {
636: if (!ds->ops->cond) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_SUP,"DS type %s",((PetscObject)ds)->type_name);
637: PetscLogEventBegin(DS_Other,ds,0,0,0);
638: PetscFPTrapPush(PETSC_FP_TRAP_OFF);
639: (*ds->ops->cond)(ds,cond);
640: PetscFPTrapPop();
641: PetscLogEventEnd(DS_Other,ds,0,0,0);
642: return(0);
643: }
647: /*@C
648: DSTranslateHarmonic - Computes a translation of the dense system.
650: Logically Collective on DS
652: Input Parameters:
653: + ds - the direct solver context
654: . tau - the translation amount
655: . beta - last component of vector b
656: - recover - boolean flag to indicate whether to recover or not
658: Output Parameters:
659: + g - the computed vector (optional)
660: - gamma - scale factor (optional)
662: Notes:
663: This function is intended for use in the context of Krylov methods only.
664: It computes a translation of a Krylov decomposition in order to extract
665: eigenpair approximations by harmonic Rayleigh-Ritz.
666: The matrix is updated as A + g*b' where g = (A-tau*eye(n))'\b and
667: vector b is assumed to be beta*e_n^T.
669: The gamma factor is defined as sqrt(1+g'*g) and can be interpreted as
670: the factor by which the residual of the Krylov decomposition is scaled.
672: If the recover flag is activated, the computed translation undoes the
673: translation done previously. In that case, parameter tau is ignored.
675: Level: developer
676: @*/
677: PetscErrorCode DSTranslateHarmonic(DS ds,PetscScalar tau,PetscReal beta,PetscBool recover,PetscScalar *g,PetscReal *gamma)
678: {
683: if (!ds->ops->transharm) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_SUP,"DS type %s",((PetscObject)ds)->type_name);
684: PetscLogEventBegin(DS_Other,ds,0,0,0);
685: PetscFPTrapPush(PETSC_FP_TRAP_OFF);
686: (*ds->ops->transharm)(ds,tau,beta,recover,g,gamma);
687: PetscFPTrapPop();
688: PetscLogEventEnd(DS_Other,ds,0,0,0);
689: ds->state = DS_STATE_RAW;
690: PetscObjectStateIncrease((PetscObject)ds);
691: return(0);
692: }
696: /*@C
697: DSTranslateRKS - Computes a modification of the dense system corresponding
698: to an update of the shift in a rational Krylov method.
700: Logically Collective on DS
702: Input Parameters:
703: + ds - the direct solver context
704: - alpha - the translation amount
706: Notes:
707: This function is intended for use in the context of Krylov methods only.
708: It takes the leading (k+1,k) submatrix of A, containing the truncated
709: Rayleigh quotient of a Krylov-Schur relation computed from a shift
710: sigma1 and transforms it to obtain a Krylov relation as if computed
711: from a different shift sigma2. The new matrix is computed as
712: 1.0/alpha*(eye(k)-Q*inv(R)), where [Q,R]=qr(eye(k)-alpha*A) and
713: alpha = sigma1-sigma2.
715: Matrix Q is placed in DS_MAT_Q so that it can be used to update the
716: Krylov basis.
718: Level: developer
719: @*/
720: PetscErrorCode DSTranslateRKS(DS ds,PetscScalar alpha)
721: {
726: if (!ds->ops->transrks) SETERRQ1(((PetscObject)ds)->comm,PETSC_ERR_SUP,"DS type %s",((PetscObject)ds)->type_name);
727: PetscLogEventBegin(DS_Other,ds,0,0,0);
728: PetscFPTrapPush(PETSC_FP_TRAP_OFF);
729: (*ds->ops->transrks)(ds,alpha);
730: PetscFPTrapPop();
731: PetscLogEventEnd(DS_Other,ds,0,0,0);
732: ds->state = DS_STATE_RAW;
733: ds->compact = PETSC_FALSE;
734: PetscObjectStateIncrease((PetscObject)ds);
735: return(0);
736: }