Actual source code: rvector.c
petsc-master 2020-01-26
2: /*
3: Provides the interface functions for vector operations that have PetscScalar/PetscReal in the signature
4: These are the vector functions the user calls.
5: */
6: #include <petsc/private/vecimpl.h>
7: #if defined(PETSC_HAVE_CUDA)
8: #include <../src/vec/vec/impls/dvecimpl.h>
9: #include <../src/vec/vec/impls/seq/seqcuda/cudavecimpl.h>
10: #endif
11: static PetscInt VecGetSubVectorSavedStateId = -1;
13: PETSC_EXTERN PetscErrorCode VecValidValues(Vec vec,PetscInt argnum,PetscBool begin)
14: {
15: #if defined(PETSC_USE_DEBUG)
16: PetscErrorCode ierr;
17: PetscInt n,i;
18: const PetscScalar *x;
21: #if defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_VIENNACL)
22: if ((vec->petscnative || vec->ops->getarray) && (vec->offloadmask & PETSC_OFFLOAD_CPU)) {
23: #else
24: if (vec->petscnative || vec->ops->getarray) {
25: #endif
26: VecGetLocalSize(vec,&n);
27: VecGetArrayRead(vec,&x);
28: for (i=0; i<n; i++) {
29: if (begin) {
30: if (PetscIsInfOrNanScalar(x[i])) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_FP,"Vec entry at local location %D is not-a-number or infinite at beginning of function: Parameter number %D",i,argnum);
31: } else {
32: if (PetscIsInfOrNanScalar(x[i])) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_FP,"Vec entry at local location %D is not-a-number or infinite at end of function: Parameter number %D",i,argnum);
33: }
34: }
35: VecRestoreArrayRead(vec,&x);
36: }
37: return(0);
38: #else
39: return 0;
40: #endif
41: }
43: /*@
44: VecMaxPointwiseDivide - Computes the maximum of the componentwise division max = max_i abs(x_i/y_i).
46: Logically Collective on Vec
48: Input Parameters:
49: . x, y - the vectors
51: Output Parameter:
52: . max - the result
54: Level: advanced
56: Notes:
57: x and y may be the same vector
58: if a particular y_i is zero, it is treated as 1 in the above formula
60: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs()
61: @*/
62: PetscErrorCode VecMaxPointwiseDivide(Vec x,Vec y,PetscReal *max)
63: {
73: VecCheckSameSize(x,1,y,2);
74: (*x->ops->maxpointwisedivide)(x,y,max);
75: return(0);
76: }
78: /*@
79: VecDot - Computes the vector dot product.
81: Collective on Vec
83: Input Parameters:
84: . x, y - the vectors
86: Output Parameter:
87: . val - the dot product
89: Performance Issues:
90: $ per-processor memory bandwidth
91: $ interprocessor latency
92: $ work load inbalance that causes certain processes to arrive much earlier than others
94: Notes for Users of Complex Numbers:
95: For complex vectors, VecDot() computes
96: $ val = (x,y) = y^H x,
97: where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
98: inner product where the SECOND argument gets the complex conjugate. Since the BLASdot() complex conjugates the first
99: first argument we call the BLASdot() with the arguments reversed.
101: Use VecTDot() for the indefinite form
102: $ val = (x,y) = y^T x,
103: where y^T denotes the transpose of y.
105: Level: intermediate
108: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDotRealPart()
109: @*/
110: PetscErrorCode VecDot(Vec x,Vec y,PetscScalar *val)
111: {
121: VecCheckSameSize(x,1,y,2);
123: PetscLogEventBegin(VEC_Dot,x,y,0,0);
124: (*x->ops->dot)(x,y,val);
125: PetscLogEventEnd(VEC_Dot,x,y,0,0);
126: return(0);
127: }
129: /*@
130: VecDotRealPart - Computes the real part of the vector dot product.
132: Collective on Vec
134: Input Parameters:
135: . x, y - the vectors
137: Output Parameter:
138: . val - the real part of the dot product;
140: Performance Issues:
141: $ per-processor memory bandwidth
142: $ interprocessor latency
143: $ work load inbalance that causes certain processes to arrive much earlier than others
145: Notes for Users of Complex Numbers:
146: See VecDot() for more details on the definition of the dot product for complex numbers
148: For real numbers this returns the same value as VecDot()
150: For complex numbers in C^n (that is a vector of n components with a complex number for each component) this is equal to the usual real dot product on the
151: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
153: Developer Note: This is not currently optimized to compute only the real part of the dot product.
155: Level: intermediate
158: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDot(), VecDotNorm2()
159: @*/
160: PetscErrorCode VecDotRealPart(Vec x,Vec y,PetscReal *val)
161: {
163: PetscScalar fdot;
166: VecDot(x,y,&fdot);
167: *val = PetscRealPart(fdot);
168: return(0);
169: }
171: /*@
172: VecNorm - Computes the vector norm.
174: Collective on Vec
176: Input Parameters:
177: + x - the vector
178: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
179: NORM_1_AND_2, which computes both norms and stores them
180: in a two element array.
182: Output Parameter:
183: . val - the norm
185: Notes:
186: $ NORM_1 denotes sum_i |x_i|
187: $ NORM_2 denotes sqrt(sum_i |x_i|^2)
188: $ NORM_INFINITY denotes max_i |x_i|
190: For complex numbers NORM_1 will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
191: norm of the absolutely values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
192: the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
193: people expect the former.
195: Level: intermediate
197: Performance Issues:
198: $ per-processor memory bandwidth
199: $ interprocessor latency
200: $ work load inbalance that causes certain processes to arrive much earlier than others
203: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNormAvailable(),
204: VecNormBegin(), VecNormEnd()
206: @*/
208: PetscErrorCode VecNorm(Vec x,NormType type,PetscReal *val)
209: {
210: PetscBool flg;
218: /*
219: * Cached data?
220: */
221: if (type!=NORM_1_AND_2) {
222: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,flg);
223: if (flg) return(0);
224: }
225: PetscLogEventBegin(VEC_Norm,x,0,0,0);
226: (*x->ops->norm)(x,type,val);
227: PetscLogEventEnd(VEC_Norm,x,0,0,0);
228: if (type!=NORM_1_AND_2) {
229: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[type],*val);
230: }
231: return(0);
232: }
234: /*@
235: VecNormAvailable - Returns the vector norm if it is already known.
237: Not Collective
239: Input Parameters:
240: + x - the vector
241: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
242: NORM_1_AND_2, which computes both norms and stores them
243: in a two element array.
245: Output Parameter:
246: + available - PETSC_TRUE if the val returned is valid
247: - val - the norm
249: Notes:
250: $ NORM_1 denotes sum_i |x_i|
251: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
252: $ NORM_INFINITY denotes max_i |x_i|
254: Level: intermediate
256: Performance Issues:
257: $ per-processor memory bandwidth
258: $ interprocessor latency
259: $ work load inbalance that causes certain processes to arrive much earlier than others
261: Compile Option:
262: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
263: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
264: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
267: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNorm()
268: VecNormBegin(), VecNormEnd()
270: @*/
271: PetscErrorCode VecNormAvailable(Vec x,NormType type,PetscBool *available,PetscReal *val)
272: {
280: *available = PETSC_FALSE;
281: if (type!=NORM_1_AND_2) {
282: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,*available);
283: }
284: return(0);
285: }
287: /*@
288: VecNormalize - Normalizes a vector by 2-norm.
290: Collective on Vec
292: Input Parameters:
293: + x - the vector
295: Output Parameter:
296: . x - the normalized vector
297: - val - the vector norm before normalization
299: Level: intermediate
302: @*/
303: PetscErrorCode VecNormalize(Vec x,PetscReal *val)
304: {
306: PetscReal norm;
311: PetscLogEventBegin(VEC_Normalize,x,0,0,0);
312: VecNorm(x,NORM_2,&norm);
313: if (norm == 0.0) {
314: PetscInfo(x,"Vector of zero norm can not be normalized; Returning only the zero norm\n");
315: } else if (norm != 1.0) {
316: PetscScalar tmp = 1.0/norm;
317: VecScale(x,tmp);
318: }
319: if (val) *val = norm;
320: PetscLogEventEnd(VEC_Normalize,x,0,0,0);
321: return(0);
322: }
324: /*@C
325: VecMax - Determines the vector component with maximum real part and its location.
327: Collective on Vec
329: Input Parameter:
330: . x - the vector
332: Output Parameters:
333: + p - the location of val (pass NULL if you don't want this)
334: - val - the maximum component
336: Notes:
337: Returns the value PETSC_MIN_REAL and p = -1 if the vector is of length 0.
339: Returns the smallest index with the maximum value
340: Level: intermediate
343: .seealso: VecNorm(), VecMin()
344: @*/
345: PetscErrorCode VecMax(Vec x,PetscInt *p,PetscReal *val)
346: {
353: PetscLogEventBegin(VEC_Max,x,0,0,0);
354: (*x->ops->max)(x,p,val);
355: PetscLogEventEnd(VEC_Max,x,0,0,0);
356: return(0);
357: }
359: /*@C
360: VecMin - Determines the vector component with minimum real part and its location.
362: Collective on Vec
364: Input Parameters:
365: . x - the vector
367: Output Parameter:
368: + p - the location of val (pass NULL if you don't want this location)
369: - val - the minimum component
371: Level: intermediate
373: Notes:
374: Returns the value PETSC_MAX_REAL and p = -1 if the vector is of length 0.
376: This returns the smallest index with the minumum value
379: .seealso: VecMax()
380: @*/
381: PetscErrorCode VecMin(Vec x,PetscInt *p,PetscReal *val)
382: {
389: PetscLogEventBegin(VEC_Min,x,0,0,0);
390: (*x->ops->min)(x,p,val);
391: PetscLogEventEnd(VEC_Min,x,0,0,0);
392: return(0);
393: }
395: /*@
396: VecTDot - Computes an indefinite vector dot product. That is, this
397: routine does NOT use the complex conjugate.
399: Collective on Vec
401: Input Parameters:
402: . x, y - the vectors
404: Output Parameter:
405: . val - the dot product
407: Notes for Users of Complex Numbers:
408: For complex vectors, VecTDot() computes the indefinite form
409: $ val = (x,y) = y^T x,
410: where y^T denotes the transpose of y.
412: Use VecDot() for the inner product
413: $ val = (x,y) = y^H x,
414: where y^H denotes the conjugate transpose of y.
416: Level: intermediate
418: .seealso: VecDot(), VecMTDot()
419: @*/
420: PetscErrorCode VecTDot(Vec x,Vec y,PetscScalar *val)
421: {
431: VecCheckSameSize(x,1,y,2);
433: PetscLogEventBegin(VEC_TDot,x,y,0,0);
434: (*x->ops->tdot)(x,y,val);
435: PetscLogEventEnd(VEC_TDot,x,y,0,0);
436: return(0);
437: }
439: /*@
440: VecScale - Scales a vector.
442: Not collective on Vec
444: Input Parameters:
445: + x - the vector
446: - alpha - the scalar
448: Output Parameter:
449: . x - the scaled vector
451: Note:
452: For a vector with n components, VecScale() computes
453: $ x[i] = alpha * x[i], for i=1,...,n.
455: Level: intermediate
458: @*/
459: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
460: {
461: PetscReal norms[4] = {0.0,0.0,0.0, 0.0};
462: PetscBool flgs[4];
464: PetscInt i;
469: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
470: PetscLogEventBegin(VEC_Scale,x,0,0,0);
471: if (alpha != (PetscScalar)1.0) {
472: VecSetErrorIfLocked(x,1);
473: /* get current stashed norms */
474: for (i=0; i<4; i++) {
475: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
476: }
477: (*x->ops->scale)(x,alpha);
478: PetscObjectStateIncrease((PetscObject)x);
479: /* put the scaled stashed norms back into the Vec */
480: for (i=0; i<4; i++) {
481: if (flgs[i]) {
482: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],PetscAbsScalar(alpha)*norms[i]);
483: }
484: }
485: }
486: PetscLogEventEnd(VEC_Scale,x,0,0,0);
487: return(0);
488: }
490: /*@
491: VecSet - Sets all components of a vector to a single scalar value.
493: Logically Collective on Vec
495: Input Parameters:
496: + x - the vector
497: - alpha - the scalar
499: Output Parameter:
500: . x - the vector
502: Note:
503: For a vector of dimension n, VecSet() computes
504: $ x[i] = alpha, for i=1,...,n,
505: so that all vector entries then equal the identical
506: scalar value, alpha. Use the more general routine
507: VecSetValues() to set different vector entries.
509: You CANNOT call this after you have called VecSetValues() but before you call
510: VecAssemblyBegin/End().
512: Level: beginner
514: .seealso VecSetValues(), VecSetValuesBlocked(), VecSetRandom()
516: @*/
517: PetscErrorCode VecSet(Vec x,PetscScalar alpha)
518: {
519: PetscReal val;
525: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"You cannot call this after you have called VecSetValues() but\n before you have called VecAssemblyBegin/End()");
527: VecSetErrorIfLocked(x,1);
529: PetscLogEventBegin(VEC_Set,x,0,0,0);
530: (*x->ops->set)(x,alpha);
531: PetscLogEventEnd(VEC_Set,x,0,0,0);
532: PetscObjectStateIncrease((PetscObject)x);
534: /* norms can be simply set (if |alpha|*N not too large) */
535: val = PetscAbsScalar(alpha);
536: if (x->map->N == 0) {
537: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],0.0l);
538: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],0.0);
539: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],0.0);
540: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],0.0);
541: } else if (val > PETSC_MAX_REAL/x->map->N) {
542: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
543: } else {
544: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],x->map->N * val);
545: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
546: val = PetscSqrtReal((PetscReal)x->map->N) * val;
547: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],val);
548: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],val);
549: }
550: return(0);
551: }
554: /*@
555: VecAXPY - Computes y = alpha x + y.
557: Logically Collective on Vec
559: Input Parameters:
560: + alpha - the scalar
561: - x, y - the vectors
563: Output Parameter:
564: . y - output vector
566: Level: intermediate
568: Notes:
569: x and y MUST be different vectors
570: This routine is optimized for alpha of 0.0, otherwise it calls the BLAS routine
572: $ VecAXPY(y,alpha,x) y = alpha x + y
573: $ VecAYPX(y,beta,x) y = x + beta y
574: $ VecAXPBY(y,alpha,beta,x) y = alpha x + beta y
575: $ VecWAXPY(w,alpha,x,y) w = alpha x + y
576: $ VecAXPBYPCZ(w,alpha,beta,gamma,x,y) z = alpha x + beta y + gamma z
577: $ VecMAXPY(y,nv,alpha[],x[]) y = sum alpha[i] x[i] + y
580: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPBYPCZ(), VecAXPBY()
581: @*/
582: PetscErrorCode VecAXPY(Vec y,PetscScalar alpha,Vec x)
583: {
592: VecCheckSameSize(x,1,y,3);
593: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
595: if (alpha == (PetscScalar)0.0) return(0);
596: VecSetErrorIfLocked(y,1);
598: VecLockReadPush(x);
599: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
600: (*y->ops->axpy)(y,alpha,x);
601: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
602: VecLockReadPop(x);
603: PetscObjectStateIncrease((PetscObject)y);
604: return(0);
605: }
607: /*@
608: VecAXPBY - Computes y = alpha x + beta y.
610: Logically Collective on Vec
612: Input Parameters:
613: + alpha,beta - the scalars
614: - x, y - the vectors
616: Output Parameter:
617: . y - output vector
619: Level: intermediate
621: Notes:
622: x and y MUST be different vectors
623: The implementation is optimized for alpha and/or beta values of 0.0 and 1.0
626: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ()
627: @*/
628: PetscErrorCode VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
629: {
638: VecCheckSameSize(y,1,x,4);
639: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
642: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)1.0) return(0);
643: VecSetErrorIfLocked(y,1);
644: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
645: (*y->ops->axpby)(y,alpha,beta,x);
646: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
647: PetscObjectStateIncrease((PetscObject)y);
648: return(0);
649: }
651: /*@
652: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
654: Logically Collective on Vec
656: Input Parameters:
657: + alpha,beta, gamma - the scalars
658: - x, y, z - the vectors
660: Output Parameter:
661: . z - output vector
663: Level: intermediate
665: Notes:
666: x, y and z must be different vectors
667: The implementation is optimized for alpha of 1.0 and gamma of 1.0 or 0.0
670: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBY()
671: @*/
672: PetscErrorCode VecAXPBYPCZ(Vec z,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec x,Vec y)
673: {
685: VecCheckSameSize(x,1,y,5);
686: VecCheckSameSize(x,1,z,6);
687: if (x == y || x == z) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
688: if (y == z) SETERRQ(PetscObjectComm((PetscObject)y),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
692: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)0.0 && gamma == (PetscScalar)1.0) return(0);
693: VecSetErrorIfLocked(z,1);
695: PetscLogEventBegin(VEC_AXPBYPCZ,x,y,z,0);
696: (*y->ops->axpbypcz)(z,alpha,beta,gamma,x,y);
697: PetscLogEventEnd(VEC_AXPBYPCZ,x,y,z,0);
698: PetscObjectStateIncrease((PetscObject)z);
699: return(0);
700: }
702: /*@
703: VecAYPX - Computes y = x + beta y.
705: Logically Collective on Vec
707: Input Parameters:
708: + beta - the scalar
709: - x, y - the vectors
711: Output Parameter:
712: . y - output vector
714: Level: intermediate
716: Notes:
717: x and y MUST be different vectors
718: The implementation is optimized for beta of -1.0, 0.0, and 1.0
721: .seealso: VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
722: @*/
723: PetscErrorCode VecAYPX(Vec y,PetscScalar beta,Vec x)
724: {
733: VecCheckSameSize(x,1,y,3);
734: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
736: VecSetErrorIfLocked(y,1);
738: PetscLogEventBegin(VEC_AYPX,x,y,0,0);
739: (*y->ops->aypx)(y,beta,x);
740: PetscLogEventEnd(VEC_AYPX,x,y,0,0);
741: PetscObjectStateIncrease((PetscObject)y);
742: return(0);
743: }
746: /*@
747: VecWAXPY - Computes w = alpha x + y.
749: Logically Collective on Vec
751: Input Parameters:
752: + alpha - the scalar
753: - x, y - the vectors
755: Output Parameter:
756: . w - the result
758: Level: intermediate
760: Notes:
761: w cannot be either x or y, but x and y can be the same
762: The implementation is optimzed for alpha of -1.0, 0.0, and 1.0
765: .seealso: VecAXPY(), VecAYPX(), VecAXPBY(), VecMAXPY(), VecAXPBYPCZ()
766: @*/
767: PetscErrorCode VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
768: {
780: VecCheckSameSize(x,3,y,4);
781: VecCheckSameSize(x,3,w,1);
782: if (w == y) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector y, suggest VecAXPY()");
783: if (w == x) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector x, suggest VecAYPX()");
785: VecSetErrorIfLocked(w,1);
787: PetscLogEventBegin(VEC_WAXPY,x,y,w,0);
788: (*w->ops->waxpy)(w,alpha,x,y);
789: PetscLogEventEnd(VEC_WAXPY,x,y,w,0);
790: PetscObjectStateIncrease((PetscObject)w);
791: return(0);
792: }
795: /*@C
796: VecSetValues - Inserts or adds values into certain locations of a vector.
798: Not Collective
800: Input Parameters:
801: + x - vector to insert in
802: . ni - number of elements to add
803: . ix - indices where to add
804: . y - array of values
805: - iora - either INSERT_VALUES or ADD_VALUES, where
806: ADD_VALUES adds values to any existing entries, and
807: INSERT_VALUES replaces existing entries with new values
809: Notes:
810: VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.
812: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
813: options cannot be mixed without intervening calls to the assembly
814: routines.
816: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
817: MUST be called after all calls to VecSetValues() have been completed.
819: VecSetValues() uses 0-based indices in Fortran as well as in C.
821: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
822: negative indices may be passed in ix. These rows are
823: simply ignored. This allows easily inserting element load matrices
824: with homogeneous Dirchlet boundary conditions that you don't want represented
825: in the vector.
827: Level: beginner
829: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
830: VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
831: @*/
832: PetscErrorCode VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
833: {
838: if (!ni) return(0);
843: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
844: (*x->ops->setvalues)(x,ni,ix,y,iora);
845: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
846: PetscObjectStateIncrease((PetscObject)x);
847: return(0);
848: }
850: /*@
851: VecGetValues - Gets values from certain locations of a vector. Currently
852: can only get values on the same processor
854: Not Collective
856: Input Parameters:
857: + x - vector to get values from
858: . ni - number of elements to get
859: - ix - indices where to get them from (in global 1d numbering)
861: Output Parameter:
862: . y - array of values
864: Notes:
865: The user provides the allocated array y; it is NOT allocated in this routine
867: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
869: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
871: VecGetValues() uses 0-based indices in Fortran as well as in C.
873: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
874: negative indices may be passed in ix. These rows are
875: simply ignored.
877: Level: beginner
879: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues()
880: @*/
881: PetscErrorCode VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
882: {
887: if (!ni) return(0);
891: (*x->ops->getvalues)(x,ni,ix,y);
892: return(0);
893: }
895: /*@C
896: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
898: Not Collective
900: Input Parameters:
901: + x - vector to insert in
902: . ni - number of blocks to add
903: . ix - indices where to add in block count, rather than element count
904: . y - array of values
905: - iora - either INSERT_VALUES or ADD_VALUES, where
906: ADD_VALUES adds values to any existing entries, and
907: INSERT_VALUES replaces existing entries with new values
909: Notes:
910: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
911: for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
913: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
914: options cannot be mixed without intervening calls to the assembly
915: routines.
917: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
918: MUST be called after all calls to VecSetValuesBlocked() have been completed.
920: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
922: Negative indices may be passed in ix, these rows are
923: simply ignored. This allows easily inserting element load matrices
924: with homogeneous Dirchlet boundary conditions that you don't want represented
925: in the vector.
927: Level: intermediate
929: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
930: VecSetValues()
931: @*/
932: PetscErrorCode VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
933: {
938: if (!ni) return(0);
943: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
944: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
945: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
946: PetscObjectStateIncrease((PetscObject)x);
947: return(0);
948: }
951: /*@C
952: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
953: using a local ordering of the nodes.
955: Not Collective
957: Input Parameters:
958: + x - vector to insert in
959: . ni - number of elements to add
960: . ix - indices where to add
961: . y - array of values
962: - iora - either INSERT_VALUES or ADD_VALUES, where
963: ADD_VALUES adds values to any existing entries, and
964: INSERT_VALUES replaces existing entries with new values
966: Level: intermediate
968: Notes:
969: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
971: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
972: options cannot be mixed without intervening calls to the assembly
973: routines.
975: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
976: MUST be called after all calls to VecSetValuesLocal() have been completed.
978: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
980: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
981: VecSetValuesBlockedLocal()
982: @*/
983: PetscErrorCode VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
984: {
986: PetscInt lixp[128],*lix = lixp;
990: if (!ni) return(0);
995: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
996: if (!x->ops->setvalueslocal) {
997: if (!x->map->mapping) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMapping()");
998: if (ni > 128) {
999: PetscMalloc1(ni,&lix);
1000: }
1001: ISLocalToGlobalMappingApply(x->map->mapping,ni,(PetscInt*)ix,lix);
1002: (*x->ops->setvalues)(x,ni,lix,y,iora);
1003: if (ni > 128) {
1004: PetscFree(lix);
1005: }
1006: } else {
1007: (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
1008: }
1009: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1010: PetscObjectStateIncrease((PetscObject)x);
1011: return(0);
1012: }
1014: /*@
1015: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
1016: using a local ordering of the nodes.
1018: Not Collective
1020: Input Parameters:
1021: + x - vector to insert in
1022: . ni - number of blocks to add
1023: . ix - indices where to add in block count, not element count
1024: . y - array of values
1025: - iora - either INSERT_VALUES or ADD_VALUES, where
1026: ADD_VALUES adds values to any existing entries, and
1027: INSERT_VALUES replaces existing entries with new values
1029: Level: intermediate
1031: Notes:
1032: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
1033: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
1035: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
1036: options cannot be mixed without intervening calls to the assembly
1037: routines.
1039: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1040: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
1042: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
1045: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(),
1046: VecSetLocalToGlobalMapping()
1047: @*/
1048: PetscErrorCode VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1049: {
1051: PetscInt lixp[128],*lix = lixp;
1055: if (!ni) return(0);
1059: if (ni > 128) {
1060: PetscMalloc1(ni,&lix);
1061: }
1063: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1064: ISLocalToGlobalMappingApplyBlock(x->map->mapping,ni,(PetscInt*)ix,lix);
1065: (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
1066: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1067: if (ni > 128) {
1068: PetscFree(lix);
1069: }
1070: PetscObjectStateIncrease((PetscObject)x);
1071: return(0);
1072: }
1074: /*@
1075: VecMTDot - Computes indefinite vector multiple dot products.
1076: That is, it does NOT use the complex conjugate.
1078: Collective on Vec
1080: Input Parameters:
1081: + x - one vector
1082: . nv - number of vectors
1083: - y - array of vectors. Note that vectors are pointers
1085: Output Parameter:
1086: . val - array of the dot products
1088: Notes for Users of Complex Numbers:
1089: For complex vectors, VecMTDot() computes the indefinite form
1090: $ val = (x,y) = y^T x,
1091: where y^T denotes the transpose of y.
1093: Use VecMDot() for the inner product
1094: $ val = (x,y) = y^H x,
1095: where y^H denotes the conjugate transpose of y.
1097: Level: intermediate
1100: .seealso: VecMDot(), VecTDot()
1101: @*/
1102: PetscErrorCode VecMTDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1103: {
1109: if (!nv) return(0);
1116: VecCheckSameSize(x,1,*y,3);
1118: PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
1119: (*x->ops->mtdot)(x,nv,y,val);
1120: PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
1121: return(0);
1122: }
1124: /*@
1125: VecMDot - Computes vector multiple dot products.
1127: Collective on Vec
1129: Input Parameters:
1130: + x - one vector
1131: . nv - number of vectors
1132: - y - array of vectors.
1134: Output Parameter:
1135: . val - array of the dot products (does not allocate the array)
1137: Notes for Users of Complex Numbers:
1138: For complex vectors, VecMDot() computes
1139: $ val = (x,y) = y^H x,
1140: where y^H denotes the conjugate transpose of y.
1142: Use VecMTDot() for the indefinite form
1143: $ val = (x,y) = y^T x,
1144: where y^T denotes the transpose of y.
1146: Level: intermediate
1149: .seealso: VecMTDot(), VecDot()
1150: @*/
1151: PetscErrorCode VecMDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1152: {
1158: if (!nv) return(0);
1159: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1166: VecCheckSameSize(x,1,*y,3);
1168: PetscLogEventBegin(VEC_MDot,x,*y,0,0);
1169: (*x->ops->mdot)(x,nv,y,val);
1170: PetscLogEventEnd(VEC_MDot,x,*y,0,0);
1171: return(0);
1172: }
1174: /*@
1175: VecMAXPY - Computes y = y + sum alpha[i] x[i]
1177: Logically Collective on Vec
1179: Input Parameters:
1180: + nv - number of scalars and x-vectors
1181: . alpha - array of scalars
1182: . y - one vector
1183: - x - array of vectors
1185: Level: intermediate
1187: Notes:
1188: y cannot be any of the x vectors
1190: .seealso: VecAYPX(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
1191: @*/
1192: PetscErrorCode VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec x[])
1193: {
1195: PetscInt i;
1196: PetscBool nonzero;
1201: if (!nv) return(0);
1202: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1209: VecCheckSameSize(y,1,*x,4);
1211: for (i=0, nonzero = PETSC_FALSE; i<nv && !nonzero; i++) nonzero = (PetscBool)(nonzero || alpha[i] != (PetscScalar)0.0);
1212: if (!nonzero) return(0);
1213: VecSetErrorIfLocked(y,1);
1214: PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
1215: (*y->ops->maxpy)(y,nv,alpha,x);
1216: PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
1217: PetscObjectStateIncrease((PetscObject)y);
1218: return(0);
1219: }
1221: /*@
1222: VecGetSubVector - Gets a vector representing part of another vector
1224: Collective on IS
1226: Input Arguments:
1227: + X - vector from which to extract a subvector
1228: - is - index set representing portion of X to extract
1230: Output Arguments:
1231: . Y - subvector corresponding to is
1233: Level: advanced
1235: Notes:
1236: The subvector Y should be returned with VecRestoreSubVector().
1238: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1239: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1241: .seealso: MatCreateSubMatrix()
1242: @*/
1243: PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
1244: {
1245: PetscErrorCode ierr;
1246: Vec Z;
1252: if (X->ops->getsubvector) {
1253: (*X->ops->getsubvector)(X,is,&Z);
1254: } else { /* Default implementation currently does no caching */
1255: PetscInt gstart,gend,start;
1256: PetscBool contiguous,gcontiguous;
1257: VecGetOwnershipRange(X,&gstart,&gend);
1258: ISContiguousLocal(is,gstart,gend,&start,&contiguous);
1259: MPIU_Allreduce(&contiguous,&gcontiguous,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1260: if (gcontiguous) { /* We can do a no-copy implementation */
1261: PetscInt n,N,bs;
1262: PetscInt state;
1264: ISGetSize(is,&N);
1265: ISGetLocalSize(is,&n);
1266: VecGetBlockSize(X,&bs);
1267: if (n%bs || bs == 1 || !n) bs = -1; /* Do not decide block size if we do not have to */
1268: VecLockGet(X,&state);
1269: if (state) {
1270: const PetscScalar *x;
1271: VecGetArrayRead(X,&x);
1272: VecCreate(PetscObjectComm((PetscObject)X),&Z);
1273: VecSetType(Z,((PetscObject)X)->type_name);
1274: VecSetSizes(Z,n,N);
1275: VecSetBlockSize(Z,bs);
1276: VecPlaceArray(Z,(PetscScalar*)x+start);
1277: VecLockReadPush(Z);
1278: VecRestoreArrayRead(X,&x);
1279: } else {
1280: PetscScalar *x;
1281: VecGetArray(X,&x);
1282: VecCreate(PetscObjectComm((PetscObject)X),&Z);
1283: VecSetType(Z,((PetscObject)X)->type_name);
1284: VecSetSizes(Z,n,N);
1285: VecSetBlockSize(Z,bs);
1286: VecPlaceArray(Z,(PetscScalar*)x+start);
1287: VecRestoreArray(X,&x);
1288: }
1289: } else { /* Have to create a scatter and do a copy */
1290: VecScatter scatter;
1291: PetscInt n,N;
1292: ISGetLocalSize(is,&n);
1293: ISGetSize(is,&N);
1294: VecCreate(PetscObjectComm((PetscObject)is),&Z);
1295: VecSetSizes(Z,n,N);
1296: VecSetType(Z,((PetscObject)X)->type_name);
1297: VecScatterCreate(X,is,Z,NULL,&scatter);
1298: VecScatterBegin(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1299: VecScatterEnd(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1300: PetscObjectCompose((PetscObject)Z,"VecGetSubVector_Scatter",(PetscObject)scatter);
1301: VecScatterDestroy(&scatter);
1302: }
1303: }
1304: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1305: if (VecGetSubVectorSavedStateId < 0) {PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);}
1306: PetscObjectComposedDataSetInt((PetscObject)Z,VecGetSubVectorSavedStateId,1);
1307: *Y = Z;
1308: return(0);
1309: }
1311: /*@
1312: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1314: Collective on IS
1316: Input Arguments:
1317: + X - vector from which subvector was obtained
1318: . is - index set representing the subset of X
1319: - Y - subvector being restored
1321: Level: advanced
1323: .seealso: VecGetSubVector()
1324: @*/
1325: PetscErrorCode VecRestoreSubVector(Vec X,IS is,Vec *Y)
1326: {
1334: if (X->ops->restoresubvector) {
1335: (*X->ops->restoresubvector)(X,is,Y);
1336: } else {
1337: PETSC_UNUSED PetscObjectState dummystate = 0;
1338: PetscBool valid;
1339: PetscObjectComposedDataGetInt((PetscObject)*Y,VecGetSubVectorSavedStateId,dummystate,valid);
1340: if (!valid) {
1341: VecScatter scatter;
1343: PetscObjectQuery((PetscObject)*Y,"VecGetSubVector_Scatter",(PetscObject*)&scatter);
1344: if (scatter) {
1345: VecScatterBegin(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1346: VecScatterEnd(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1347: }
1348: }
1349: VecDestroy(Y);
1350: }
1351: return(0);
1352: }
1354: /*@
1355: VecGetLocalVectorRead - Maps the local portion of a vector into a
1356: vector. You must call VecRestoreLocalVectorRead() when the local
1357: vector is no longer needed.
1359: Not collective.
1361: Input parameter:
1362: . v - The vector for which the local vector is desired.
1364: Output parameter:
1365: . w - Upon exit this contains the local vector.
1367: Level: beginner
1369: Notes:
1370: This function is similar to VecGetArrayRead() which maps the local
1371: portion into a raw pointer. VecGetLocalVectorRead() is usually
1372: almost as efficient as VecGetArrayRead() but in certain circumstances
1373: VecGetLocalVectorRead() can be much more efficient than
1374: VecGetArrayRead(). This is because the construction of a contiguous
1375: array representing the vector data required by VecGetArrayRead() can
1376: be an expensive operation for certain vector types. For example, for
1377: GPU vectors VecGetArrayRead() requires that the data between device
1378: and host is synchronized.
1380: Unlike VecGetLocalVector(), this routine is not collective and
1381: preserves cached information.
1383: .seealso: VecRestoreLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1384: @*/
1385: PetscErrorCode VecGetLocalVectorRead(Vec v,Vec w)
1386: {
1388: PetscScalar *a;
1393: VecCheckSameLocalSize(v,1,w,2);
1394: if (v->ops->getlocalvectorread) {
1395: (*v->ops->getlocalvectorread)(v,w);
1396: } else {
1397: VecGetArrayRead(v,(const PetscScalar**)&a);
1398: VecPlaceArray(w,a);
1399: }
1400: return(0);
1401: }
1403: /*@
1404: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1405: previously mapped into a vector using VecGetLocalVectorRead().
1407: Not collective.
1409: Input parameter:
1410: + v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1411: - w - The vector into which the local portion of v was mapped.
1413: Level: beginner
1415: .seealso: VecGetLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1416: @*/
1417: PetscErrorCode VecRestoreLocalVectorRead(Vec v,Vec w)
1418: {
1420: PetscScalar *a;
1425: if (v->ops->restorelocalvectorread) {
1426: (*v->ops->restorelocalvectorread)(v,w);
1427: } else {
1428: VecGetArrayRead(w,(const PetscScalar**)&a);
1429: VecRestoreArrayRead(v,(const PetscScalar**)&a);
1430: VecResetArray(w);
1431: }
1432: return(0);
1433: }
1435: /*@
1436: VecGetLocalVector - Maps the local portion of a vector into a
1437: vector.
1439: Collective on v, not collective on w.
1441: Input parameter:
1442: . v - The vector for which the local vector is desired.
1444: Output parameter:
1445: . w - Upon exit this contains the local vector.
1447: Level: beginner
1449: Notes:
1450: This function is similar to VecGetArray() which maps the local
1451: portion into a raw pointer. VecGetLocalVector() is usually about as
1452: efficient as VecGetArray() but in certain circumstances
1453: VecGetLocalVector() can be much more efficient than VecGetArray().
1454: This is because the construction of a contiguous array representing
1455: the vector data required by VecGetArray() can be an expensive
1456: operation for certain vector types. For example, for GPU vectors
1457: VecGetArray() requires that the data between device and host is
1458: synchronized.
1460: .seealso: VecRestoreLocalVector(), VecGetLocalVectorRead(), VecGetArrayRead(), VecGetArray()
1461: @*/
1462: PetscErrorCode VecGetLocalVector(Vec v,Vec w)
1463: {
1465: PetscScalar *a;
1470: VecCheckSameLocalSize(v,1,w,2);
1471: if (v->ops->getlocalvector) {
1472: (*v->ops->getlocalvector)(v,w);
1473: } else {
1474: VecGetArray(v,&a);
1475: VecPlaceArray(w,a);
1476: }
1477: return(0);
1478: }
1480: /*@
1481: VecRestoreLocalVector - Unmaps the local portion of a vector
1482: previously mapped into a vector using VecGetLocalVector().
1484: Logically collective.
1486: Input parameter:
1487: + v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1488: - w - The vector into which the local portion of v was mapped.
1490: Level: beginner
1492: .seealso: VecGetLocalVector(), VecGetLocalVectorRead(), VecRestoreLocalVectorRead(), LocalVectorRead(), VecGetArrayRead(), VecGetArray()
1493: @*/
1494: PetscErrorCode VecRestoreLocalVector(Vec v,Vec w)
1495: {
1497: PetscScalar *a;
1502: if (v->ops->restorelocalvector) {
1503: (*v->ops->restorelocalvector)(v,w);
1504: } else {
1505: VecGetArray(w,&a);
1506: VecRestoreArray(v,&a);
1507: VecResetArray(w);
1508: }
1509: return(0);
1510: }
1512: /*@C
1513: VecGetArray - Returns a pointer to a contiguous array that contains this
1514: processor's portion of the vector data. For the standard PETSc
1515: vectors, VecGetArray() returns a pointer to the local data array and
1516: does not use any copies. If the underlying vector data is not stored
1517: in a contiguous array this routine will copy the data to a contiguous
1518: array and return a pointer to that. You MUST call VecRestoreArray()
1519: when you no longer need access to the array.
1521: Logically Collective on Vec
1523: Input Parameter:
1524: . x - the vector
1526: Output Parameter:
1527: . a - location to put pointer to the array
1529: Fortran Note:
1530: This routine is used differently from Fortran 77
1531: $ Vec x
1532: $ PetscScalar x_array(1)
1533: $ PetscOffset i_x
1534: $ PetscErrorCode ierr
1535: $ call VecGetArray(x,x_array,i_x,ierr)
1536: $
1537: $ Access first local entry in vector with
1538: $ value = x_array(i_x + 1)
1539: $
1540: $ ...... other code
1541: $ call VecRestoreArray(x,x_array,i_x,ierr)
1542: For Fortran 90 see VecGetArrayF90()
1544: See the Fortran chapter of the users manual and
1545: petsc/src/snes/examples/tutorials/ex5f.F for details.
1547: Level: beginner
1549: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1550: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
1551: @*/
1552: PetscErrorCode VecGetArray(Vec x,PetscScalar **a)
1553: {
1555: #if defined(PETSC_HAVE_VIENNACL)
1556: PetscBool is_viennacltype = PETSC_FALSE;
1557: #endif
1561: VecSetErrorIfLocked(x,1);
1562: if (x->petscnative) {
1563: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1564: if (x->offloadmask == PETSC_OFFLOAD_GPU) {
1565: #if defined(PETSC_HAVE_VIENNACL)
1566: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1567: if (is_viennacltype) {
1568: VecViennaCLCopyFromGPU(x);
1569: } else
1570: #endif
1571: {
1572: #if defined(PETSC_HAVE_CUDA)
1573: VecCUDACopyFromGPU(x);
1574: #endif
1575: }
1576: } else if (x->offloadmask == PETSC_OFFLOAD_UNALLOCATED) {
1577: #if defined(PETSC_HAVE_VIENNACL)
1578: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1579: if (is_viennacltype) {
1580: VecViennaCLAllocateCheckHost(x);
1581: } else
1582: #endif
1583: {
1584: #if defined(PETSC_HAVE_CUDA)
1585: VecCUDAAllocateCheckHost(x);
1586: #endif
1587: }
1588: }
1589: #endif
1590: *a = *((PetscScalar**)x->data);
1591: } else {
1592: if (x->ops->getarray) {
1593: (*x->ops->getarray)(x,a);
1594: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array for vector type \"%s\"",((PetscObject)x)->type_name);
1595: }
1596: return(0);
1597: }
1599: /*@C
1600: VecGetArrayInPlace - Like VecGetArray(), but if this is a CUDA vector and it is currently offloaded to GPU,
1601: the returned pointer will be a GPU pointer to the GPU memory that contains this processor's portion of the
1602: vector data. Otherwise, it functions as VecGetArray().
1604: Logically Collective on Vec
1606: Input Parameter:
1607: . x - the vector
1609: Output Parameter:
1610: . a - location to put pointer to the array
1612: Level: beginner
1614: .seealso: VecRestoreArrayInPlace(), VecRestoreArrayInPlace(), VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(),
1615: VecPlaceArray(), VecGetArray2d(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
1616: @*/
1617: PetscErrorCode VecGetArrayInPlace(Vec x,PetscScalar **a)
1618: {
1623: VecSetErrorIfLocked(x,1);
1625: #if defined(PETSC_HAVE_CUDA)
1626: if (x->petscnative && (x->offloadmask & PETSC_OFFLOAD_GPU)) { /* Prefer working on GPU when offloadmask is PETSC_OFFLOAD_BOTH */
1627: PetscBool is_cudatype = PETSC_FALSE;
1628: PetscObjectTypeCompareAny((PetscObject)x,&is_cudatype,VECSEQCUDA,VECMPICUDA,VECCUDA,"");
1629: if (is_cudatype) {
1630: VecCUDAGetArray(x,a);
1631: x->offloadmask = PETSC_OFFLOAD_GPU; /* Change the mask once GPU gets write access, don't wait until restore array */
1632: return(0);
1633: }
1634: }
1635: #endif
1636: VecGetArray(x,a);
1637: return(0);
1638: }
1640: /*@C
1641: VecGetArrayWrite - Returns a pointer to a contiguous array that WILL contains this
1642: processor's portion of the vector data. The values in this array are NOT valid, the routine calling this
1643: routine is responsible for putting values into the array; any values it does not set will be invalid
1645: Logically Collective on Vec
1647: Input Parameter:
1648: . x - the vector
1650: Output Parameter:
1651: . a - location to put pointer to the array
1653: Level: intermediate
1655: This is for vectors associate with GPUs, the vector is not copied up before giving access. If you need correct
1656: values in the array use VecGetArray()
1658: Concepts: vector^accessing local values
1660: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1661: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArray(), VecRestoreArrayWrite()
1662: @*/
1663: PetscErrorCode VecGetArrayWrite(Vec x,PetscScalar **a)
1664: {
1669: VecSetErrorIfLocked(x,1);
1670: if (!x->ops->getarraywrite) {
1671: VecGetArray(x,a);
1672: } else {
1673: (*x->ops->getarraywrite)(x,a);
1674: }
1675: return(0);
1676: }
1678: /*@C
1679: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1681: Not Collective
1683: Input Parameters:
1684: . x - the vector
1686: Output Parameter:
1687: . a - the array
1689: Level: beginner
1691: Notes:
1692: The array must be returned using a matching call to VecRestoreArrayRead().
1694: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1696: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1697: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1698: only one copy is performed when this routine is called multiple times in sequence.
1700: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1701: @*/
1702: PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar **a)
1703: {
1705: #if defined(PETSC_HAVE_VIENNACL)
1706: PetscBool is_viennacltype = PETSC_FALSE;
1707: #endif
1711: if (x->petscnative) {
1712: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1713: if (x->offloadmask == PETSC_OFFLOAD_GPU) {
1714: #if defined(PETSC_HAVE_VIENNACL)
1715: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1716: if (is_viennacltype) {
1717: VecViennaCLCopyFromGPU(x);
1718: } else
1719: #endif
1720: {
1721: #if defined(PETSC_HAVE_CUDA)
1722: VecCUDACopyFromGPU(x);
1723: #endif
1724: }
1725: }
1726: #endif
1727: *a = *((PetscScalar **)x->data);
1728: } else if (x->ops->getarrayread) {
1729: (*x->ops->getarrayread)(x,a);
1730: } else {
1731: (*x->ops->getarray)(x,(PetscScalar**)a);
1732: }
1733: return(0);
1734: }
1736: /*@C
1737: VecGetArrayReadInPlace - Like VecGetArrayRead(), but if this is a CUDA vector and it is currently offloaded to GPU,
1738: the returned pointer will be a GPU pointer to the GPU memory that contains this processor's portion of the
1739: vector data. Otherwise, it functions as VecGetArrayRead().
1741: Not Collective
1743: Input Parameters:
1744: . x - the vector
1746: Output Parameter:
1747: . a - the array
1749: Level: beginner
1751: Notes:
1752: The array must be returned using a matching call to VecRestoreArrayReadInPlace().
1755: .seealso: VecRestoreArrayReadInPlace(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayInPlace()
1756: @*/
1757: PetscErrorCode VecGetArrayReadInPlace(Vec x,const PetscScalar **a)
1758: {
1763: #if defined(PETSC_HAVE_CUDA)
1764: if (x->petscnative && x->offloadmask & PETSC_OFFLOAD_GPU) {
1765: PetscBool is_cudatype = PETSC_FALSE;
1766: PetscObjectTypeCompareAny((PetscObject)x,&is_cudatype,VECSEQCUDA,VECMPICUDA,VECCUDA,"");
1767: if (is_cudatype) {
1768: VecCUDAGetArrayRead(x,a);
1769: return(0);
1770: }
1771: }
1772: #endif
1773: VecGetArrayRead(x,a);
1774: return(0);
1775: }
1777: /*@C
1778: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1779: that were created by a call to VecDuplicateVecs(). You MUST call
1780: VecRestoreArrays() when you no longer need access to the array.
1782: Logically Collective on Vec
1784: Input Parameter:
1785: + x - the vectors
1786: - n - the number of vectors
1788: Output Parameter:
1789: . a - location to put pointer to the array
1791: Fortran Note:
1792: This routine is not supported in Fortran.
1794: Level: intermediate
1796: .seealso: VecGetArray(), VecRestoreArrays()
1797: @*/
1798: PetscErrorCode VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1799: {
1801: PetscInt i;
1802: PetscScalar **q;
1808: if (n <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must get at least one array n = %D",n);
1809: PetscMalloc1(n,&q);
1810: for (i=0; i<n; ++i) {
1811: VecGetArray(x[i],&q[i]);
1812: }
1813: *a = q;
1814: return(0);
1815: }
1817: /*@C
1818: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
1819: has been called.
1821: Logically Collective on Vec
1823: Input Parameters:
1824: + x - the vector
1825: . n - the number of vectors
1826: - a - location of pointer to arrays obtained from VecGetArrays()
1828: Notes:
1829: For regular PETSc vectors this routine does not involve any copies. For
1830: any special vectors that do not store local vector data in a contiguous
1831: array, this routine will copy the data back into the underlying
1832: vector data structure from the arrays obtained with VecGetArrays().
1834: Fortran Note:
1835: This routine is not supported in Fortran.
1837: Level: intermediate
1839: .seealso: VecGetArrays(), VecRestoreArray()
1840: @*/
1841: PetscErrorCode VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1842: {
1844: PetscInt i;
1845: PetscScalar **q = *a;
1852: for (i=0; i<n; ++i) {
1853: VecRestoreArray(x[i],&q[i]);
1854: }
1855: PetscFree(q);
1856: return(0);
1857: }
1859: /*@C
1860: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1862: Logically Collective on Vec
1864: Input Parameters:
1865: + x - the vector
1866: - a - location of pointer to array obtained from VecGetArray()
1868: Level: beginner
1870: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1871: VecGetArrayPair(), VecRestoreArrayPair()
1872: @*/
1873: PetscErrorCode VecRestoreArray(Vec x,PetscScalar **a)
1874: {
1879: if (x->petscnative) {
1880: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1881: x->offloadmask = PETSC_OFFLOAD_CPU;
1882: #endif
1883: } else {
1884: (*x->ops->restorearray)(x,a);
1885: }
1886: if (a) *a = NULL;
1887: PetscObjectStateIncrease((PetscObject)x);
1888: return(0);
1889: }
1891: /*@C
1892: VecRestoreArrayInPlace - Restores a vector after VecGetArrayInPlace() has been called.
1894: Logically Collective on Vec
1896: Input Parameters:
1897: + x - the vector
1898: - a - location of pointer to array obtained from VecGetArrayInPlace()
1900: Level: beginner
1902: .seealso: VecGetArrayInPlace(), VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(),
1903: VecPlaceArray(), VecRestoreArray2d(), VecGetArrayPair(), VecRestoreArrayPair()
1904: @*/
1905: PetscErrorCode VecRestoreArrayInPlace(Vec x,PetscScalar **a)
1906: {
1911: #if defined(PETSC_HAVE_CUDA)
1912: if (x->petscnative && x->offloadmask == PETSC_OFFLOAD_GPU) {
1913: PetscBool is_cudatype = PETSC_FALSE;
1914: PetscObjectTypeCompareAny((PetscObject)x,&is_cudatype,VECSEQCUDA,VECMPICUDA,VECCUDA,"");
1915: if (is_cudatype) {
1916: VecCUDARestoreArray(x,a);
1917: return(0);
1918: }
1919: }
1920: #endif
1921: VecRestoreArray(x,a);
1922: return(0);
1923: }
1926: /*@C
1927: VecRestoreArrayWrite - Restores a vector after VecGetArrayWrite() has been called.
1929: Logically Collective on Vec
1931: Input Parameters:
1932: + x - the vector
1933: - a - location of pointer to array obtained from VecGetArray()
1935: Level: beginner
1937: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1938: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite()
1939: @*/
1940: PetscErrorCode VecRestoreArrayWrite(Vec x,PetscScalar **a)
1941: {
1946: if (x->petscnative) {
1947: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1948: x->offloadmask = PETSC_OFFLOAD_CPU;
1949: #endif
1950: } else {
1951: if (x->ops->restorearraywrite) {
1952: (*x->ops->restorearraywrite)(x,a);
1953: } else {
1954: (*x->ops->restorearray)(x,a);
1955: }
1956: }
1957: if (a) *a = NULL;
1958: PetscObjectStateIncrease((PetscObject)x);
1959: return(0);
1960: }
1962: /*@C
1963: VecRestoreArrayRead - Restore array obtained with VecGetArrayRead()
1965: Not Collective
1967: Input Parameters:
1968: + vec - the vector
1969: - array - the array
1971: Level: beginner
1973: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1974: @*/
1975: PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar **a)
1976: {
1981: if (x->petscnative) {
1982: /* nothing */
1983: } else if (x->ops->restorearrayread) {
1984: (*x->ops->restorearrayread)(x,a);
1985: } else {
1986: (*x->ops->restorearray)(x,(PetscScalar**)a);
1987: }
1988: if (a) *a = NULL;
1989: return(0);
1990: }
1992: /*@C
1993: VecRestoreArrayReadInPlace - Restore array obtained with VecGetArrayReadInPlace()
1995: Not Collective
1997: Input Parameters:
1998: + vec - the vector
1999: - array - the array
2001: Level: beginner
2003: .seealso: VecGetArrayReadInPlace(), VecRestoreArrayInPlace(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
2004: @*/
2005: PetscErrorCode VecRestoreArrayReadInPlace(Vec x,const PetscScalar **a)
2006: {
2010: VecRestoreArrayRead(x,a);
2011: return(0);
2012: }
2014: /*@
2015: VecPlaceArray - Allows one to replace the array in a vector with an
2016: array provided by the user. This is useful to avoid copying an array
2017: into a vector.
2019: Not Collective
2021: Input Parameters:
2022: + vec - the vector
2023: - array - the array
2025: Notes:
2026: You can return to the original array with a call to VecResetArray()
2028: Level: developer
2030: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()
2032: @*/
2033: PetscErrorCode VecPlaceArray(Vec vec,const PetscScalar array[])
2034: {
2041: if (vec->ops->placearray) {
2042: (*vec->ops->placearray)(vec,array);
2043: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot place array in this type of vector");
2044: PetscObjectStateIncrease((PetscObject)vec);
2045: return(0);
2046: }
2048: /*@C
2049: VecReplaceArray - Allows one to replace the array in a vector with an
2050: array provided by the user. This is useful to avoid copying an array
2051: into a vector.
2053: Not Collective
2055: Input Parameters:
2056: + vec - the vector
2057: - array - the array
2059: Notes:
2060: This permanently replaces the array and frees the memory associated
2061: with the old array.
2063: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
2064: freed by the user. It will be freed when the vector is destroy.
2066: Not supported from Fortran
2068: Level: developer
2070: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()
2072: @*/
2073: PetscErrorCode VecReplaceArray(Vec vec,const PetscScalar array[])
2074: {
2080: if (vec->ops->replacearray) {
2081: (*vec->ops->replacearray)(vec,array);
2082: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot replace array in this type of vector");
2083: PetscObjectStateIncrease((PetscObject)vec);
2084: return(0);
2085: }
2088: /*@C
2089: VecCUDAGetArray - Provides access to the CUDA buffer inside a vector.
2091: This function has semantics similar to VecGetArray(): the pointer
2092: returned by this function points to a consistent view of the vector
2093: data. This may involve a copy operation of data from the host to the
2094: device if the data on the device is out of date. If the device
2095: memory hasn't been allocated previously it will be allocated as part
2096: of this function call. VecCUDAGetArray() assumes that
2097: the user will modify the vector data. This is similar to
2098: intent(inout) in fortran.
2100: The CUDA device pointer has to be released by calling
2101: VecCUDARestoreArray(). Upon restoring the vector data
2102: the data on the host will be marked as out of date. A subsequent
2103: access of the host data will thus incur a data transfer from the
2104: device to the host.
2107: Input Parameter:
2108: . v - the vector
2110: Output Parameter:
2111: . a - the CUDA device pointer
2113: Fortran note:
2114: This function is not currently available from Fortran.
2116: Level: intermediate
2118: .seealso: VecCUDARestoreArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2119: @*/
2120: PETSC_EXTERN PetscErrorCode VecCUDAGetArray(Vec v, PetscScalar **a)
2121: {
2122: #if defined(PETSC_HAVE_CUDA)
2124: #endif
2128: #if defined(PETSC_HAVE_CUDA)
2129: *a = 0;
2130: VecCUDACopyToGPU(v);
2131: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2132: #endif
2133: return(0);
2134: }
2136: /*@C
2137: VecCUDARestoreArray - Restore a CUDA device pointer previously acquired with VecCUDAGetArray().
2139: This marks the host data as out of date. Subsequent access to the
2140: vector data on the host side with for instance VecGetArray() incurs a
2141: data transfer.
2143: Input Parameter:
2144: + v - the vector
2145: - a - the CUDA device pointer. This pointer is invalid after
2146: VecCUDARestoreArray() returns.
2148: Fortran note:
2149: This function is not currently available from Fortran.
2151: Level: intermediate
2153: .seealso: VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2154: @*/
2155: PETSC_EXTERN PetscErrorCode VecCUDARestoreArray(Vec v, PetscScalar **a)
2156: {
2161: #if defined(PETSC_HAVE_CUDA)
2162: v->offloadmask = PETSC_OFFLOAD_GPU;
2163: #endif
2165: PetscObjectStateIncrease((PetscObject)v);
2166: return(0);
2167: }
2169: /*@C
2170: VecCUDAGetArrayRead - Provides read access to the CUDA buffer inside a vector.
2172: This function is analogous to VecGetArrayRead(): The pointer
2173: returned by this function points to a consistent view of the vector
2174: data. This may involve a copy operation of data from the host to the
2175: device if the data on the device is out of date. If the device
2176: memory hasn't been allocated previously it will be allocated as part
2177: of this function call. VecCUDAGetArrayRead() assumes that the
2178: user will not modify the vector data. This is analgogous to
2179: intent(in) in Fortran.
2181: The CUDA device pointer has to be released by calling
2182: VecCUDARestoreArrayRead(). If the data on the host side was
2183: previously up to date it will remain so, i.e. data on both the device
2184: and the host is up to date. Accessing data on the host side does not
2185: incur a device to host data transfer.
2187: Input Parameter:
2188: . v - the vector
2190: Output Parameter:
2191: . a - the CUDA pointer.
2193: Fortran note:
2194: This function is not currently available from Fortran.
2196: Level: intermediate
2198: .seealso: VecCUDARestoreArrayRead(), VecCUDAGetArray(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2199: @*/
2200: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayRead(Vec v, const PetscScalar **a)
2201: {
2202: #if defined(PETSC_HAVE_CUDA)
2204: #endif
2208: #if defined(PETSC_HAVE_CUDA)
2209: *a = 0;
2210: VecCUDACopyToGPU(v);
2211: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2212: #endif
2213: return(0);
2214: }
2216: /*@C
2217: VecCUDARestoreArrayRead - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayRead().
2219: If the data on the host side was previously up to date it will remain
2220: so, i.e. data on both the device and the host is up to date.
2221: Accessing data on the host side e.g. with VecGetArray() does not
2222: incur a device to host data transfer.
2224: Input Parameter:
2225: + v - the vector
2226: - a - the CUDA device pointer. This pointer is invalid after
2227: VecCUDARestoreArrayRead() returns.
2229: Fortran note:
2230: This function is not currently available from Fortran.
2232: Level: intermediate
2234: .seealso: VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2235: @*/
2236: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayRead(Vec v, const PetscScalar **a)
2237: {
2240: *a = NULL;
2241: return(0);
2242: }
2244: /*@C
2245: VecCUDAGetArrayWrite - Provides write access to the CUDA buffer inside a vector.
2247: The data pointed to by the device pointer is uninitialized. The user
2248: may not read from this data. Furthermore, the entire array needs to
2249: be filled by the user to obtain well-defined behaviour. The device
2250: memory will be allocated by this function if it hasn't been allocated
2251: previously. This is analogous to intent(out) in Fortran.
2253: The device pointer needs to be released with
2254: VecCUDARestoreArrayWrite(). When the pointer is released the
2255: host data of the vector is marked as out of data. Subsequent access
2256: of the host data with e.g. VecGetArray() incurs a device to host data
2257: transfer.
2260: Input Parameter:
2261: . v - the vector
2263: Output Parameter:
2264: . a - the CUDA pointer
2266: Fortran note:
2267: This function is not currently available from Fortran.
2269: Level: advanced
2271: .seealso: VecCUDARestoreArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2272: @*/
2273: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayWrite(Vec v, PetscScalar **a)
2274: {
2275: #if defined(PETSC_HAVE_CUDA)
2277: #endif
2281: #if defined(PETSC_HAVE_CUDA)
2282: *a = 0;
2283: VecCUDAAllocateCheck(v);
2284: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2285: #endif
2286: return(0);
2287: }
2289: /*@C
2290: VecCUDARestoreArrayWrite - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayWrite().
2292: Data on the host will be marked as out of date. Subsequent access of
2293: the data on the host side e.g. with VecGetArray() will incur a device
2294: to host data transfer.
2296: Input Parameter:
2297: + v - the vector
2298: - a - the CUDA device pointer. This pointer is invalid after
2299: VecCUDARestoreArrayWrite() returns.
2301: Fortran note:
2302: This function is not currently available from Fortran.
2304: Level: intermediate
2306: .seealso: VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2307: @*/
2308: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayWrite(Vec v, PetscScalar **a)
2309: {
2314: #if defined(PETSC_HAVE_CUDA)
2315: v->offloadmask = PETSC_OFFLOAD_GPU;
2316: #endif
2318: PetscObjectStateIncrease((PetscObject)v);
2319: return(0);
2320: }
2322: /*@C
2323: VecCUDAPlaceArray - Allows one to replace the GPU array in a vector with a
2324: GPU array provided by the user. This is useful to avoid copying an
2325: array into a vector.
2327: Not Collective
2329: Input Parameters:
2330: + vec - the vector
2331: - array - the GPU array
2333: Notes:
2334: You can return to the original GPU array with a call to VecCUDAResetArray()
2335: It is not possible to use VecCUDAPlaceArray() and VecPlaceArray() at the
2336: same time on the same vector.
2338: Level: developer
2340: .seealso: VecPlaceArray(), VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAReplaceArray()
2342: @*/
2343: PetscErrorCode VecCUDAPlaceArray(Vec vin,PetscScalar *a)
2344: {
2349: #if defined(PETSC_HAVE_CUDA)
2350: VecCUDACopyToGPU(vin);
2351: if (((Vec_Seq*)vin->data)->unplacedarray) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"VecCUDAPlaceArray()/VecPlaceArray() was already called on this vector, without a call to VecCUDAResetArray()/VecResetArray()");
2352: ((Vec_Seq*)vin->data)->unplacedarray = (PetscScalar *) ((Vec_CUDA*)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2353: ((Vec_CUDA*)vin->spptr)->GPUarray = a;
2354: vin->offloadmask = PETSC_OFFLOAD_GPU;
2355: #endif
2356: PetscObjectStateIncrease((PetscObject)vin);
2357: return(0);
2358: }
2360: /*@C
2361: VecCUDAReplaceArray - Allows one to replace the GPU array in a vector
2362: with a GPU array provided by the user. This is useful to avoid copying
2363: a GPU array into a vector.
2365: Not Collective
2367: Input Parameters:
2368: + vec - the vector
2369: - array - the GPU array
2371: Notes:
2372: This permanently replaces the GPU array and frees the memory associated
2373: with the old GPU array.
2375: The memory passed in CANNOT be freed by the user. It will be freed
2376: when the vector is destroyed.
2378: Not supported from Fortran
2380: Level: developer
2382: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAPlaceArray(), VecReplaceArray()
2384: @*/
2385: PetscErrorCode VecCUDAReplaceArray(Vec vin,PetscScalar *a)
2386: {
2387: #if defined(PETSC_HAVE_CUDA)
2388: cudaError_t err;
2389: #endif
2394: #if defined(PETSC_HAVE_CUDA)
2395: err = cudaFree(((Vec_CUDA*)vin->spptr)->GPUarray);CHKERRCUDA(err);
2396: ((Vec_CUDA*)vin->spptr)->GPUarray = a;
2397: vin->offloadmask = PETSC_OFFLOAD_GPU;
2398: #endif
2399: PetscObjectStateIncrease((PetscObject)vin);
2400: return(0);
2401: }
2403: /*@C
2404: VecCUDAResetArray - Resets a vector to use its default memory. Call this
2405: after the use of VecCUDAPlaceArray().
2407: Not Collective
2409: Input Parameters:
2410: . vec - the vector
2412: Level: developer
2414: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray(), VecResetArray(), VecCUDAPlaceArray(), VecCUDAReplaceArray()
2416: @*/
2417: PetscErrorCode VecCUDAResetArray(Vec vin)
2418: {
2423: #if defined(PETSC_HAVE_CUDA)
2424: VecCUDACopyToGPU(vin);
2425: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar *) ((Vec_Seq*)vin->data)->unplacedarray;
2426: ((Vec_Seq*)vin->data)->unplacedarray = 0;
2427: vin->offloadmask = PETSC_OFFLOAD_GPU;
2428: #endif
2429: PetscObjectStateIncrease((PetscObject)vin);
2430: return(0);
2431: }
2433: /*MC
2434: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
2435: and makes them accessible via a Fortran90 pointer.
2437: Synopsis:
2438: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
2440: Collective on Vec
2442: Input Parameters:
2443: + x - a vector to mimic
2444: - n - the number of vectors to obtain
2446: Output Parameters:
2447: + y - Fortran90 pointer to the array of vectors
2448: - ierr - error code
2450: Example of Usage:
2451: .vb
2452: #include <petsc/finclude/petscvec.h>
2453: use petscvec
2455: Vec x
2456: Vec, pointer :: y(:)
2457: ....
2458: call VecDuplicateVecsF90(x,2,y,ierr)
2459: call VecSet(y(2),alpha,ierr)
2460: call VecSet(y(2),alpha,ierr)
2461: ....
2462: call VecDestroyVecsF90(2,y,ierr)
2463: .ve
2465: Notes:
2466: Not yet supported for all F90 compilers
2468: Use VecDestroyVecsF90() to free the space.
2470: Level: beginner
2472: .seealso: VecDestroyVecsF90(), VecDuplicateVecs()
2474: M*/
2476: /*MC
2477: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
2478: VecGetArrayF90().
2480: Synopsis:
2481: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2483: Logically Collective on Vec
2485: Input Parameters:
2486: + x - vector
2487: - xx_v - the Fortran90 pointer to the array
2489: Output Parameter:
2490: . ierr - error code
2492: Example of Usage:
2493: .vb
2494: #include <petsc/finclude/petscvec.h>
2495: use petscvec
2497: PetscScalar, pointer :: xx_v(:)
2498: ....
2499: call VecGetArrayF90(x,xx_v,ierr)
2500: xx_v(3) = a
2501: call VecRestoreArrayF90(x,xx_v,ierr)
2502: .ve
2504: Level: beginner
2506: .seealso: VecGetArrayF90(), VecGetArray(), VecRestoreArray(), UsingFortran, VecRestoreArrayReadF90()
2508: M*/
2510: /*MC
2511: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
2513: Synopsis:
2514: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
2516: Collective on Vec
2518: Input Parameters:
2519: + n - the number of vectors previously obtained
2520: - x - pointer to array of vector pointers
2522: Output Parameter:
2523: . ierr - error code
2525: Notes:
2526: Not yet supported for all F90 compilers
2528: Level: beginner
2530: .seealso: VecDestroyVecs(), VecDuplicateVecsF90()
2532: M*/
2534: /*MC
2535: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
2536: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2537: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
2538: when you no longer need access to the array.
2540: Synopsis:
2541: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2543: Logically Collective on Vec
2545: Input Parameter:
2546: . x - vector
2548: Output Parameters:
2549: + xx_v - the Fortran90 pointer to the array
2550: - ierr - error code
2552: Example of Usage:
2553: .vb
2554: #include <petsc/finclude/petscvec.h>
2555: use petscvec
2557: PetscScalar, pointer :: xx_v(:)
2558: ....
2559: call VecGetArrayF90(x,xx_v,ierr)
2560: xx_v(3) = a
2561: call VecRestoreArrayF90(x,xx_v,ierr)
2562: .ve
2564: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
2566: Level: beginner
2568: .seealso: VecRestoreArrayF90(), VecGetArray(), VecRestoreArray(), VecGetArrayReadF90(), UsingFortran
2570: M*/
2572: /*MC
2573: VecGetArrayReadF90 - Accesses a read only array from Fortran90. For default PETSc
2574: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2575: this routine is implementation dependent. You MUST call VecRestoreArrayReadF90()
2576: when you no longer need access to the array.
2578: Synopsis:
2579: VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2581: Logically Collective on Vec
2583: Input Parameter:
2584: . x - vector
2586: Output Parameters:
2587: + xx_v - the Fortran90 pointer to the array
2588: - ierr - error code
2590: Example of Usage:
2591: .vb
2592: #include <petsc/finclude/petscvec.h>
2593: use petscvec
2595: PetscScalar, pointer :: xx_v(:)
2596: ....
2597: call VecGetArrayReadF90(x,xx_v,ierr)
2598: a = xx_v(3)
2599: call VecRestoreArrayReadF90(x,xx_v,ierr)
2600: .ve
2602: If you intend to write entries into the array you must use VecGetArrayF90().
2604: Level: beginner
2606: .seealso: VecRestoreArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecGetArrayF90(), UsingFortran
2608: M*/
2610: /*MC
2611: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
2612: VecGetArrayReadF90().
2614: Synopsis:
2615: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2617: Logically Collective on Vec
2619: Input Parameters:
2620: + x - vector
2621: - xx_v - the Fortran90 pointer to the array
2623: Output Parameter:
2624: . ierr - error code
2626: Example of Usage:
2627: .vb
2628: #include <petsc/finclude/petscvec.h>
2629: use petscvec
2631: PetscScalar, pointer :: xx_v(:)
2632: ....
2633: call VecGetArrayReadF90(x,xx_v,ierr)
2634: a = xx_v(3)
2635: call VecRestoreArrayReadF90(x,xx_v,ierr)
2636: .ve
2638: Level: beginner
2640: .seealso: VecGetArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(),UsingFortran, VecRestoreArrayF90()
2642: M*/
2644: /*@C
2645: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
2646: processor's portion of the vector data. You MUST call VecRestoreArray2d()
2647: when you no longer need access to the array.
2649: Logically Collective
2651: Input Parameter:
2652: + x - the vector
2653: . m - first dimension of two dimensional array
2654: . n - second dimension of two dimensional array
2655: . mstart - first index you will use in first coordinate direction (often 0)
2656: - nstart - first index in the second coordinate direction (often 0)
2658: Output Parameter:
2659: . a - location to put pointer to the array
2661: Level: developer
2663: Notes:
2664: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2665: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2666: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2667: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2669: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2671: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2672: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2673: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2674: @*/
2675: PetscErrorCode VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2676: {
2678: PetscInt i,N;
2679: PetscScalar *aa;
2685: VecGetLocalSize(x,&N);
2686: if (m*n != N) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 2d array dimensions %D by %D",N,m,n);
2687: VecGetArray(x,&aa);
2689: PetscMalloc1(m,a);
2690: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
2691: *a -= mstart;
2692: return(0);
2693: }
2695: /*@C
2696: VecGetArray2dWrite - Returns a pointer to a 2d contiguous array that will contain this
2697: processor's portion of the vector data. You MUST call VecRestoreArray2dWrite()
2698: when you no longer need access to the array.
2700: Logically Collective
2702: Input Parameter:
2703: + x - the vector
2704: . m - first dimension of two dimensional array
2705: . n - second dimension of two dimensional array
2706: . mstart - first index you will use in first coordinate direction (often 0)
2707: - nstart - first index in the second coordinate direction (often 0)
2709: Output Parameter:
2710: . a - location to put pointer to the array
2712: Level: developer
2714: Notes:
2715: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2716: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2717: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2718: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2720: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2722: Concepts: vector^accessing local values as 2d array
2724: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2725: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2726: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2727: @*/
2728: PetscErrorCode VecGetArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2729: {
2731: PetscInt i,N;
2732: PetscScalar *aa;
2738: VecGetLocalSize(x,&N);
2739: if (m*n != N) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 2d array dimensions %D by %D",N,m,n);
2740: VecGetArrayWrite(x,&aa);
2742: PetscMalloc1(m,a);
2743: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
2744: *a -= mstart;
2745: return(0);
2746: }
2748: /*@C
2749: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
2751: Logically Collective
2753: Input Parameters:
2754: + x - the vector
2755: . m - first dimension of two dimensional array
2756: . n - second dimension of the two dimensional array
2757: . mstart - first index you will use in first coordinate direction (often 0)
2758: . nstart - first index in the second coordinate direction (often 0)
2759: - a - location of pointer to array obtained from VecGetArray2d()
2761: Level: developer
2763: Notes:
2764: For regular PETSc vectors this routine does not involve any copies. For
2765: any special vectors that do not store local vector data in a contiguous
2766: array, this routine will copy the data back into the underlying
2767: vector data structure from the array obtained with VecGetArray().
2769: This routine actually zeros out the a pointer.
2771: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2772: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2773: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2774: @*/
2775: PetscErrorCode VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2776: {
2778: void *dummy;
2784: dummy = (void*)(*a + mstart);
2785: PetscFree(dummy);
2786: VecRestoreArray(x,NULL);
2787: return(0);
2788: }
2790: /*@C
2791: VecRestoreArray2dWrite - Restores a vector after VecGetArray2dWrite() has been called.
2793: Logically Collective
2795: Input Parameters:
2796: + x - the vector
2797: . m - first dimension of two dimensional array
2798: . n - second dimension of the two dimensional array
2799: . mstart - first index you will use in first coordinate direction (often 0)
2800: . nstart - first index in the second coordinate direction (often 0)
2801: - a - location of pointer to array obtained from VecGetArray2d()
2803: Level: developer
2805: Notes:
2806: For regular PETSc vectors this routine does not involve any copies. For
2807: any special vectors that do not store local vector data in a contiguous
2808: array, this routine will copy the data back into the underlying
2809: vector data structure from the array obtained with VecGetArray().
2811: This routine actually zeros out the a pointer.
2813: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2814: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2815: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2816: @*/
2817: PetscErrorCode VecRestoreArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2818: {
2820: void *dummy;
2826: dummy = (void*)(*a + mstart);
2827: PetscFree(dummy);
2828: VecRestoreArrayWrite(x,NULL);
2829: return(0);
2830: }
2832: /*@C
2833: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
2834: processor's portion of the vector data. You MUST call VecRestoreArray1d()
2835: when you no longer need access to the array.
2837: Logically Collective
2839: Input Parameter:
2840: + x - the vector
2841: . m - first dimension of two dimensional array
2842: - mstart - first index you will use in first coordinate direction (often 0)
2844: Output Parameter:
2845: . a - location to put pointer to the array
2847: Level: developer
2849: Notes:
2850: For a vector obtained from DMCreateLocalVector() mstart are likely
2851: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2852: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2854: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2856: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2857: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2858: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2859: @*/
2860: PetscErrorCode VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2861: {
2863: PetscInt N;
2869: VecGetLocalSize(x,&N);
2870: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2871: VecGetArray(x,a);
2872: *a -= mstart;
2873: return(0);
2874: }
2876: /*@C
2877: VecGetArray1dWrite - Returns a pointer to a 1d contiguous array that will contain this
2878: processor's portion of the vector data. You MUST call VecRestoreArray1dWrite()
2879: when you no longer need access to the array.
2881: Logically Collective
2883: Input Parameter:
2884: + x - the vector
2885: . m - first dimension of two dimensional array
2886: - mstart - first index you will use in first coordinate direction (often 0)
2888: Output Parameter:
2889: . a - location to put pointer to the array
2891: Level: developer
2893: Notes:
2894: For a vector obtained from DMCreateLocalVector() mstart are likely
2895: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2896: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2898: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2900: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2901: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2902: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2903: @*/
2904: PetscErrorCode VecGetArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2905: {
2907: PetscInt N;
2913: VecGetLocalSize(x,&N);
2914: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2915: VecGetArrayWrite(x,a);
2916: *a -= mstart;
2917: return(0);
2918: }
2920: /*@C
2921: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
2923: Logically Collective
2925: Input Parameters:
2926: + x - the vector
2927: . m - first dimension of two dimensional array
2928: . mstart - first index you will use in first coordinate direction (often 0)
2929: - a - location of pointer to array obtained from VecGetArray21()
2931: Level: developer
2933: Notes:
2934: For regular PETSc vectors this routine does not involve any copies. For
2935: any special vectors that do not store local vector data in a contiguous
2936: array, this routine will copy the data back into the underlying
2937: vector data structure from the array obtained with VecGetArray1d().
2939: This routine actually zeros out the a pointer.
2941: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2942: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2943: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2944: @*/
2945: PetscErrorCode VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2946: {
2952: VecRestoreArray(x,NULL);
2953: return(0);
2954: }
2956: /*@C
2957: VecRestoreArray1dWrite - Restores a vector after VecGetArray1dWrite() has been called.
2959: Logically Collective
2961: Input Parameters:
2962: + x - the vector
2963: . m - first dimension of two dimensional array
2964: . mstart - first index you will use in first coordinate direction (often 0)
2965: - a - location of pointer to array obtained from VecGetArray21()
2967: Level: developer
2969: Notes:
2970: For regular PETSc vectors this routine does not involve any copies. For
2971: any special vectors that do not store local vector data in a contiguous
2972: array, this routine will copy the data back into the underlying
2973: vector data structure from the array obtained with VecGetArray1d().
2975: This routine actually zeros out the a pointer.
2977: Concepts: vector^accessing local values as 1d array
2979: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2980: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2981: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2982: @*/
2983: PetscErrorCode VecRestoreArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2984: {
2990: VecRestoreArrayWrite(x,NULL);
2991: return(0);
2992: }
2994: /*@C
2995: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
2996: processor's portion of the vector data. You MUST call VecRestoreArray3d()
2997: when you no longer need access to the array.
2999: Logically Collective
3001: Input Parameter:
3002: + x - the vector
3003: . m - first dimension of three dimensional array
3004: . n - second dimension of three dimensional array
3005: . p - third dimension of three dimensional array
3006: . mstart - first index you will use in first coordinate direction (often 0)
3007: . nstart - first index in the second coordinate direction (often 0)
3008: - pstart - first index in the third coordinate direction (often 0)
3010: Output Parameter:
3011: . a - location to put pointer to the array
3013: Level: developer
3015: Notes:
3016: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3017: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3018: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3019: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3021: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3023: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3024: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3025: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3026: @*/
3027: PetscErrorCode VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3028: {
3030: PetscInt i,N,j;
3031: PetscScalar *aa,**b;
3037: VecGetLocalSize(x,&N);
3038: if (m*n*p != N) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 3d array dimensions %D by %D by %D",N,m,n,p);
3039: VecGetArray(x,&aa);
3041: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3042: b = (PetscScalar**)((*a) + m);
3043: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3044: for (i=0; i<m; i++)
3045: for (j=0; j<n; j++)
3046: b[i*n+j] = aa + i*n*p + j*p - pstart;
3048: *a -= mstart;
3049: return(0);
3050: }
3052: /*@C
3053: VecGetArray3dWrite - Returns a pointer to a 3d contiguous array that will contain this
3054: processor's portion of the vector data. You MUST call VecRestoreArray3dWrite()
3055: when you no longer need access to the array.
3057: Logically Collective
3059: Input Parameter:
3060: + x - the vector
3061: . m - first dimension of three dimensional array
3062: . n - second dimension of three dimensional array
3063: . p - third dimension of three dimensional array
3064: . mstart - first index you will use in first coordinate direction (often 0)
3065: . nstart - first index in the second coordinate direction (often 0)
3066: - pstart - first index in the third coordinate direction (often 0)
3068: Output Parameter:
3069: . a - location to put pointer to the array
3071: Level: developer
3073: Notes:
3074: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3075: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3076: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3077: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3079: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3081: Concepts: vector^accessing local values as 3d array
3083: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3084: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3085: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3086: @*/
3087: PetscErrorCode VecGetArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3088: {
3090: PetscInt i,N,j;
3091: PetscScalar *aa,**b;
3097: VecGetLocalSize(x,&N);
3098: if (m*n*p != N) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 3d array dimensions %D by %D by %D",N,m,n,p);
3099: VecGetArrayWrite(x,&aa);
3101: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3102: b = (PetscScalar**)((*a) + m);
3103: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3104: for (i=0; i<m; i++)
3105: for (j=0; j<n; j++)
3106: b[i*n+j] = aa + i*n*p + j*p - pstart;
3108: *a -= mstart;
3109: return(0);
3110: }
3112: /*@C
3113: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
3115: Logically Collective
3117: Input Parameters:
3118: + x - the vector
3119: . m - first dimension of three dimensional array
3120: . n - second dimension of the three dimensional array
3121: . p - third dimension of the three dimensional array
3122: . mstart - first index you will use in first coordinate direction (often 0)
3123: . nstart - first index in the second coordinate direction (often 0)
3124: . pstart - first index in the third coordinate direction (often 0)
3125: - a - location of pointer to array obtained from VecGetArray3d()
3127: Level: developer
3129: Notes:
3130: For regular PETSc vectors this routine does not involve any copies. For
3131: any special vectors that do not store local vector data in a contiguous
3132: array, this routine will copy the data back into the underlying
3133: vector data structure from the array obtained with VecGetArray().
3135: This routine actually zeros out the a pointer.
3137: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3138: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3139: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3140: @*/
3141: PetscErrorCode VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3142: {
3144: void *dummy;
3150: dummy = (void*)(*a + mstart);
3151: PetscFree(dummy);
3152: VecRestoreArray(x,NULL);
3153: return(0);
3154: }
3156: /*@C
3157: VecRestoreArray3dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3159: Logically Collective
3161: Input Parameters:
3162: + x - the vector
3163: . m - first dimension of three dimensional array
3164: . n - second dimension of the three dimensional array
3165: . p - third dimension of the three dimensional array
3166: . mstart - first index you will use in first coordinate direction (often 0)
3167: . nstart - first index in the second coordinate direction (often 0)
3168: . pstart - first index in the third coordinate direction (often 0)
3169: - a - location of pointer to array obtained from VecGetArray3d()
3171: Level: developer
3173: Notes:
3174: For regular PETSc vectors this routine does not involve any copies. For
3175: any special vectors that do not store local vector data in a contiguous
3176: array, this routine will copy the data back into the underlying
3177: vector data structure from the array obtained with VecGetArray().
3179: This routine actually zeros out the a pointer.
3181: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3182: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3183: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3184: @*/
3185: PetscErrorCode VecRestoreArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3186: {
3188: void *dummy;
3194: dummy = (void*)(*a + mstart);
3195: PetscFree(dummy);
3196: VecRestoreArrayWrite(x,NULL);
3197: return(0);
3198: }
3200: /*@C
3201: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
3202: processor's portion of the vector data. You MUST call VecRestoreArray4d()
3203: when you no longer need access to the array.
3205: Logically Collective
3207: Input Parameter:
3208: + x - the vector
3209: . m - first dimension of four dimensional array
3210: . n - second dimension of four dimensional array
3211: . p - third dimension of four dimensional array
3212: . q - fourth dimension of four dimensional array
3213: . mstart - first index you will use in first coordinate direction (often 0)
3214: . nstart - first index in the second coordinate direction (often 0)
3215: . pstart - first index in the third coordinate direction (often 0)
3216: - qstart - first index in the fourth coordinate direction (often 0)
3218: Output Parameter:
3219: . a - location to put pointer to the array
3221: Level: beginner
3223: Notes:
3224: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3225: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3226: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3227: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3229: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3231: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3232: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3233: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3234: @*/
3235: PetscErrorCode VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3236: {
3238: PetscInt i,N,j,k;
3239: PetscScalar *aa,***b,**c;
3245: VecGetLocalSize(x,&N);
3246: if (m*n*p*q != N) SETERRQ5(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 4d array dimensions %D by %D by %D by %D",N,m,n,p,q);
3247: VecGetArray(x,&aa);
3249: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3250: b = (PetscScalar***)((*a) + m);
3251: c = (PetscScalar**)(b + m*n);
3252: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3253: for (i=0; i<m; i++)
3254: for (j=0; j<n; j++)
3255: b[i*n+j] = c + i*n*p + j*p - pstart;
3256: for (i=0; i<m; i++)
3257: for (j=0; j<n; j++)
3258: for (k=0; k<p; k++)
3259: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3260: *a -= mstart;
3261: return(0);
3262: }
3264: /*@C
3265: VecGetArray4dWrite - Returns a pointer to a 4d contiguous array that will contain this
3266: processor's portion of the vector data. You MUST call VecRestoreArray4dWrite()
3267: when you no longer need access to the array.
3269: Logically Collective
3271: Input Parameter:
3272: + x - the vector
3273: . m - first dimension of four dimensional array
3274: . n - second dimension of four dimensional array
3275: . p - third dimension of four dimensional array
3276: . q - fourth dimension of four dimensional array
3277: . mstart - first index you will use in first coordinate direction (often 0)
3278: . nstart - first index in the second coordinate direction (often 0)
3279: . pstart - first index in the third coordinate direction (often 0)
3280: - qstart - first index in the fourth coordinate direction (often 0)
3282: Output Parameter:
3283: . a - location to put pointer to the array
3285: Level: beginner
3287: Notes:
3288: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3289: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3290: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3291: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3293: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3295: Concepts: vector^accessing local values as 3d array
3297: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3298: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3299: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3300: @*/
3301: PetscErrorCode VecGetArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3302: {
3304: PetscInt i,N,j,k;
3305: PetscScalar *aa,***b,**c;
3311: VecGetLocalSize(x,&N);
3312: if (m*n*p*q != N) SETERRQ5(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 4d array dimensions %D by %D by %D by %D",N,m,n,p,q);
3313: VecGetArrayWrite(x,&aa);
3315: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3316: b = (PetscScalar***)((*a) + m);
3317: c = (PetscScalar**)(b + m*n);
3318: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3319: for (i=0; i<m; i++)
3320: for (j=0; j<n; j++)
3321: b[i*n+j] = c + i*n*p + j*p - pstart;
3322: for (i=0; i<m; i++)
3323: for (j=0; j<n; j++)
3324: for (k=0; k<p; k++)
3325: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3326: *a -= mstart;
3327: return(0);
3328: }
3330: /*@C
3331: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
3333: Logically Collective
3335: Input Parameters:
3336: + x - the vector
3337: . m - first dimension of four dimensional array
3338: . n - second dimension of the four dimensional array
3339: . p - third dimension of the four dimensional array
3340: . q - fourth dimension of the four dimensional array
3341: . mstart - first index you will use in first coordinate direction (often 0)
3342: . nstart - first index in the second coordinate direction (often 0)
3343: . pstart - first index in the third coordinate direction (often 0)
3344: . qstart - first index in the fourth coordinate direction (often 0)
3345: - a - location of pointer to array obtained from VecGetArray4d()
3347: Level: beginner
3349: Notes:
3350: For regular PETSc vectors this routine does not involve any copies. For
3351: any special vectors that do not store local vector data in a contiguous
3352: array, this routine will copy the data back into the underlying
3353: vector data structure from the array obtained with VecGetArray().
3355: This routine actually zeros out the a pointer.
3357: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3358: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3359: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3360: @*/
3361: PetscErrorCode VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3362: {
3364: void *dummy;
3370: dummy = (void*)(*a + mstart);
3371: PetscFree(dummy);
3372: VecRestoreArray(x,NULL);
3373: return(0);
3374: }
3376: /*@C
3377: VecRestoreArray4dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3379: Logically Collective
3381: Input Parameters:
3382: + x - the vector
3383: . m - first dimension of four dimensional array
3384: . n - second dimension of the four dimensional array
3385: . p - third dimension of the four dimensional array
3386: . q - fourth dimension of the four dimensional array
3387: . mstart - first index you will use in first coordinate direction (often 0)
3388: . nstart - first index in the second coordinate direction (often 0)
3389: . pstart - first index in the third coordinate direction (often 0)
3390: . qstart - first index in the fourth coordinate direction (often 0)
3391: - a - location of pointer to array obtained from VecGetArray4d()
3393: Level: beginner
3395: Notes:
3396: For regular PETSc vectors this routine does not involve any copies. For
3397: any special vectors that do not store local vector data in a contiguous
3398: array, this routine will copy the data back into the underlying
3399: vector data structure from the array obtained with VecGetArray().
3401: This routine actually zeros out the a pointer.
3403: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3404: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3405: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3406: @*/
3407: PetscErrorCode VecRestoreArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3408: {
3410: void *dummy;
3416: dummy = (void*)(*a + mstart);
3417: PetscFree(dummy);
3418: VecRestoreArrayWrite(x,NULL);
3419: return(0);
3420: }
3422: /*@C
3423: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
3424: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
3425: when you no longer need access to the array.
3427: Logically Collective
3429: Input Parameter:
3430: + x - the vector
3431: . m - first dimension of two dimensional array
3432: . n - second dimension of two dimensional array
3433: . mstart - first index you will use in first coordinate direction (often 0)
3434: - nstart - first index in the second coordinate direction (often 0)
3436: Output Parameter:
3437: . a - location to put pointer to the array
3439: Level: developer
3441: Notes:
3442: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3443: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3444: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3445: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3447: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3449: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3450: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3451: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3452: @*/
3453: PetscErrorCode VecGetArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3454: {
3455: PetscErrorCode ierr;
3456: PetscInt i,N;
3457: const PetscScalar *aa;
3463: VecGetLocalSize(x,&N);
3464: if (m*n != N) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 2d array dimensions %D by %D",N,m,n);
3465: VecGetArrayRead(x,&aa);
3467: PetscMalloc1(m,a);
3468: for (i=0; i<m; i++) (*a)[i] = (PetscScalar*) aa + i*n - nstart;
3469: *a -= mstart;
3470: return(0);
3471: }
3473: /*@C
3474: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
3476: Logically Collective
3478: Input Parameters:
3479: + x - the vector
3480: . m - first dimension of two dimensional array
3481: . n - second dimension of the two dimensional array
3482: . mstart - first index you will use in first coordinate direction (often 0)
3483: . nstart - first index in the second coordinate direction (often 0)
3484: - a - location of pointer to array obtained from VecGetArray2d()
3486: Level: developer
3488: Notes:
3489: For regular PETSc vectors this routine does not involve any copies. For
3490: any special vectors that do not store local vector data in a contiguous
3491: array, this routine will copy the data back into the underlying
3492: vector data structure from the array obtained with VecGetArray().
3494: This routine actually zeros out the a pointer.
3496: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3497: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3498: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3499: @*/
3500: PetscErrorCode VecRestoreArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3501: {
3503: void *dummy;
3509: dummy = (void*)(*a + mstart);
3510: PetscFree(dummy);
3511: VecRestoreArrayRead(x,NULL);
3512: return(0);
3513: }
3515: /*@C
3516: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
3517: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
3518: when you no longer need access to the array.
3520: Logically Collective
3522: Input Parameter:
3523: + x - the vector
3524: . m - first dimension of two dimensional array
3525: - mstart - first index you will use in first coordinate direction (often 0)
3527: Output Parameter:
3528: . a - location to put pointer to the array
3530: Level: developer
3532: Notes:
3533: For a vector obtained from DMCreateLocalVector() mstart are likely
3534: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3535: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3537: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3539: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3540: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3541: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3542: @*/
3543: PetscErrorCode VecGetArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3544: {
3546: PetscInt N;
3552: VecGetLocalSize(x,&N);
3553: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
3554: VecGetArrayRead(x,(const PetscScalar**)a);
3555: *a -= mstart;
3556: return(0);
3557: }
3559: /*@C
3560: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
3562: Logically Collective
3564: Input Parameters:
3565: + x - the vector
3566: . m - first dimension of two dimensional array
3567: . mstart - first index you will use in first coordinate direction (often 0)
3568: - a - location of pointer to array obtained from VecGetArray21()
3570: Level: developer
3572: Notes:
3573: For regular PETSc vectors this routine does not involve any copies. For
3574: any special vectors that do not store local vector data in a contiguous
3575: array, this routine will copy the data back into the underlying
3576: vector data structure from the array obtained with VecGetArray1dRead().
3578: This routine actually zeros out the a pointer.
3580: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3581: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3582: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3583: @*/
3584: PetscErrorCode VecRestoreArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3585: {
3591: VecRestoreArrayRead(x,NULL);
3592: return(0);
3593: }
3596: /*@C
3597: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
3598: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
3599: when you no longer need access to the array.
3601: Logically Collective
3603: Input Parameter:
3604: + x - the vector
3605: . m - first dimension of three dimensional array
3606: . n - second dimension of three dimensional array
3607: . p - third dimension of three dimensional array
3608: . mstart - first index you will use in first coordinate direction (often 0)
3609: . nstart - first index in the second coordinate direction (often 0)
3610: - pstart - first index in the third coordinate direction (often 0)
3612: Output Parameter:
3613: . a - location to put pointer to the array
3615: Level: developer
3617: Notes:
3618: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3619: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3620: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3621: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
3623: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3625: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3626: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3627: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3628: @*/
3629: PetscErrorCode VecGetArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3630: {
3631: PetscErrorCode ierr;
3632: PetscInt i,N,j;
3633: const PetscScalar *aa;
3634: PetscScalar **b;
3640: VecGetLocalSize(x,&N);
3641: if (m*n*p != N) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 3d array dimensions %D by %D by %D",N,m,n,p);
3642: VecGetArrayRead(x,&aa);
3644: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3645: b = (PetscScalar**)((*a) + m);
3646: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3647: for (i=0; i<m; i++)
3648: for (j=0; j<n; j++)
3649: b[i*n+j] = (PetscScalar *)aa + i*n*p + j*p - pstart;
3651: *a -= mstart;
3652: return(0);
3653: }
3655: /*@C
3656: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
3658: Logically Collective
3660: Input Parameters:
3661: + x - the vector
3662: . m - first dimension of three dimensional array
3663: . n - second dimension of the three dimensional array
3664: . p - third dimension of the three dimensional array
3665: . mstart - first index you will use in first coordinate direction (often 0)
3666: . nstart - first index in the second coordinate direction (often 0)
3667: . pstart - first index in the third coordinate direction (often 0)
3668: - a - location of pointer to array obtained from VecGetArray3dRead()
3670: Level: developer
3672: Notes:
3673: For regular PETSc vectors this routine does not involve any copies. For
3674: any special vectors that do not store local vector data in a contiguous
3675: array, this routine will copy the data back into the underlying
3676: vector data structure from the array obtained with VecGetArray().
3678: This routine actually zeros out the a pointer.
3680: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3681: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3682: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3683: @*/
3684: PetscErrorCode VecRestoreArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3685: {
3687: void *dummy;
3693: dummy = (void*)(*a + mstart);
3694: PetscFree(dummy);
3695: VecRestoreArrayRead(x,NULL);
3696: return(0);
3697: }
3699: /*@C
3700: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
3701: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
3702: when you no longer need access to the array.
3704: Logically Collective
3706: Input Parameter:
3707: + x - the vector
3708: . m - first dimension of four dimensional array
3709: . n - second dimension of four dimensional array
3710: . p - third dimension of four dimensional array
3711: . q - fourth dimension of four dimensional array
3712: . mstart - first index you will use in first coordinate direction (often 0)
3713: . nstart - first index in the second coordinate direction (often 0)
3714: . pstart - first index in the third coordinate direction (often 0)
3715: - qstart - first index in the fourth coordinate direction (often 0)
3717: Output Parameter:
3718: . a - location to put pointer to the array
3720: Level: beginner
3722: Notes:
3723: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3724: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3725: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3726: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3728: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3730: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3731: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3732: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3733: @*/
3734: PetscErrorCode VecGetArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3735: {
3736: PetscErrorCode ierr;
3737: PetscInt i,N,j,k;
3738: const PetscScalar *aa;
3739: PetscScalar ***b,**c;
3745: VecGetLocalSize(x,&N);
3746: if (m*n*p*q != N) SETERRQ5(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Local array size %D does not match 4d array dimensions %D by %D by %D by %D",N,m,n,p,q);
3747: VecGetArrayRead(x,&aa);
3749: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3750: b = (PetscScalar***)((*a) + m);
3751: c = (PetscScalar**)(b + m*n);
3752: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3753: for (i=0; i<m; i++)
3754: for (j=0; j<n; j++)
3755: b[i*n+j] = c + i*n*p + j*p - pstart;
3756: for (i=0; i<m; i++)
3757: for (j=0; j<n; j++)
3758: for (k=0; k<p; k++)
3759: c[i*n*p+j*p+k] = (PetscScalar*) aa + i*n*p*q + j*p*q + k*q - qstart;
3760: *a -= mstart;
3761: return(0);
3762: }
3764: /*@C
3765: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
3767: Logically Collective
3769: Input Parameters:
3770: + x - the vector
3771: . m - first dimension of four dimensional array
3772: . n - second dimension of the four dimensional array
3773: . p - third dimension of the four dimensional array
3774: . q - fourth dimension of the four dimensional array
3775: . mstart - first index you will use in first coordinate direction (often 0)
3776: . nstart - first index in the second coordinate direction (often 0)
3777: . pstart - first index in the third coordinate direction (often 0)
3778: . qstart - first index in the fourth coordinate direction (often 0)
3779: - a - location of pointer to array obtained from VecGetArray4dRead()
3781: Level: beginner
3783: Notes:
3784: For regular PETSc vectors this routine does not involve any copies. For
3785: any special vectors that do not store local vector data in a contiguous
3786: array, this routine will copy the data back into the underlying
3787: vector data structure from the array obtained with VecGetArray().
3789: This routine actually zeros out the a pointer.
3791: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3792: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3793: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3794: @*/
3795: PetscErrorCode VecRestoreArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3796: {
3798: void *dummy;
3804: dummy = (void*)(*a + mstart);
3805: PetscFree(dummy);
3806: VecRestoreArrayRead(x,NULL);
3807: return(0);
3808: }
3810: #if defined(PETSC_USE_DEBUG)
3812: /*@
3813: VecLockGet - Gets the current lock status of a vector
3815: Logically Collective on Vec
3817: Input Parameter:
3818: . x - the vector
3820: Output Parameter:
3821: . state - greater than zero indicates the vector is locked for read; less then zero indicates the vector is
3822: locked for write; equal to zero means the vector is unlocked, that is, it is free to read or write.
3824: Level: beginner
3826: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop()
3827: @*/
3828: PetscErrorCode VecLockGet(Vec x,PetscInt *state)
3829: {
3832: *state = x->lock;
3833: return(0);
3834: }
3836: /*@
3837: VecLockReadPush - Pushes a read-only lock on a vector to prevent it from writing
3839: Logically Collective on Vec
3841: Input Parameter:
3842: . x - the vector
3844: Notes:
3845: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
3847: The call can be nested, i.e., called multiple times on the same vector, but each VecLockReadPush(x) has to have one matching
3848: VecLockReadPop(x), which removes the latest read-only lock.
3850: Level: beginner
3852: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPop(), VecLockGet()
3853: @*/
3854: PetscErrorCode VecLockReadPush(Vec x)
3855: {
3858: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector is already locked for exclusive write access but you want to read it");
3859: x->lock++;
3860: return(0);
3861: }
3863: /*@
3864: VecLockReadPop - Pops a read-only lock from a vector
3866: Logically Collective on Vec
3868: Input Parameter:
3869: . x - the vector
3871: Level: beginner
3873: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockGet()
3874: @*/
3875: PetscErrorCode VecLockReadPop(Vec x)
3876: {
3879: x->lock--;
3880: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector has been unlocked from read-only access too many times");
3881: return(0);
3882: }
3884: /*@C
3885: VecLockWriteSet_Private - Lock or unlock a vector for exclusive read/write access
3887: Logically Collective on Vec
3889: Input Parameter:
3890: + x - the vector
3891: - flg - PETSC_TRUE to lock the vector for writing; PETSC_FALSE to unlock it.
3893: Notes:
3894: The function is usefull in split-phase computations, which usually have a begin phase and an end phase.
3895: One can call VecLockWriteSet_Private(x,PETSC_TRUE) in the begin phase to lock a vector for exclusive
3896: access, and call VecLockWriteSet_Private(x,PETSC_FALSE) in the end phase to unlock the vector from exclusive
3897: access. In this way, one is ensured no other operations can access the vector in between. The code may like
3900: VecGetArray(x,&xdata); // begin phase
3901: VecLockWriteSet_Private(v,PETSC_TRUE);
3903: Other operations, which can not acceess x anymore (they can access xdata, of course)
3905: VecRestoreArray(x,&vdata); // end phase
3906: VecLockWriteSet_Private(v,PETSC_FALSE);
3908: The call can not be nested on the same vector, in other words, one can not call VecLockWriteSet_Private(x,PETSC_TRUE)
3909: again before calling VecLockWriteSet_Private(v,PETSC_FALSE).
3911: Level: beginner
3913: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop(), VecLockGet()
3914: @*/
3915: PetscErrorCode VecLockWriteSet_Private(Vec x,PetscBool flg)
3916: {
3919: if (flg) {
3920: if (x->lock > 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector is already locked for read-only access but you want to write it");
3921: else if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector is already locked for exclusive write access but you want to write it");
3922: else x->lock = -1;
3923: } else {
3924: if (x->lock != -1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector is not locked for exclusive write access but you want to unlock it from that");
3925: x->lock = 0;
3926: }
3927: return(0);
3928: }
3930: /*@
3931: VecLockPush - Pushes a read-only lock on a vector to prevent it from writing
3933: Level: deprecated
3935: .seealso: VecLockReadPush()
3936: @*/
3937: PetscErrorCode VecLockPush(Vec x)
3938: {
3941: VecLockReadPush(x);
3942: return(0);
3943: }
3945: /*@
3946: VecLockPop - Pops a read-only lock from a vector
3948: Level: deprecated
3950: .seealso: VecLockReadPop()
3951: @*/
3952: PetscErrorCode VecLockPop(Vec x)
3953: {
3956: VecLockReadPop(x);
3957: return(0);
3958: }
3960: #endif