Actual source code: rvector.c
petsc-master 2020-07-05
1: /*
2: Provides the interface functions for vector operations that have PetscScalar/PetscReal in the signature
3: These are the vector functions the user calls.
4: */
5: #include <petsc/private/vecimpl.h>
6: #if defined(PETSC_HAVE_CUDA)
7: #include <../src/vec/vec/impls/dvecimpl.h>
8: #include <../src/vec/vec/impls/seq/seqcuda/cudavecimpl.h>
9: #endif
10: static PetscInt VecGetSubVectorSavedStateId = -1;
12: PETSC_EXTERN PetscErrorCode VecValidValues(Vec vec,PetscInt argnum,PetscBool begin)
13: {
14: #if defined(PETSC_USE_DEBUG)
15: PetscErrorCode ierr;
16: PetscInt n,i;
17: const PetscScalar *x;
20: #if defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_VIENNACL)
21: if ((vec->petscnative || vec->ops->getarray) && (vec->offloadmask & PETSC_OFFLOAD_CPU)) {
22: #else
23: if (vec->petscnative || vec->ops->getarray) {
24: #endif
25: VecGetLocalSize(vec,&n);
26: VecGetArrayRead(vec,&x);
27: for (i=0; i<n; i++) {
28: if (begin) {
29: 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);
30: } else {
31: 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);
32: }
33: }
34: VecRestoreArrayRead(vec,&x);
35: }
36: return(0);
37: #else
38: return 0;
39: #endif
40: }
42: /*@
43: VecMaxPointwiseDivide - Computes the maximum of the componentwise division max = max_i abs(x_i/y_i).
45: Logically Collective on Vec
47: Input Parameters:
48: . x, y - the vectors
50: Output Parameter:
51: . max - the result
53: Level: advanced
55: Notes:
56: x and y may be the same vector
57: if a particular y_i is zero, it is treated as 1 in the above formula
59: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs()
60: @*/
61: PetscErrorCode VecMaxPointwiseDivide(Vec x,Vec y,PetscReal *max)
62: {
72: VecCheckSameSize(x,1,y,2);
73: (*x->ops->maxpointwisedivide)(x,y,max);
74: return(0);
75: }
77: /*@
78: VecDot - Computes the vector dot product.
80: Collective on Vec
82: Input Parameters:
83: . x, y - the vectors
85: Output Parameter:
86: . val - the dot product
88: Performance Issues:
89: $ per-processor memory bandwidth
90: $ interprocessor latency
91: $ work load inbalance that causes certain processes to arrive much earlier than others
93: Notes for Users of Complex Numbers:
94: For complex vectors, VecDot() computes
95: $ val = (x,y) = y^H x,
96: where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
97: inner product where the SECOND argument gets the complex conjugate. Since the BLASdot() complex conjugates the first
98: first argument we call the BLASdot() with the arguments reversed.
100: Use VecTDot() for the indefinite form
101: $ val = (x,y) = y^T x,
102: where y^T denotes the transpose of y.
104: Level: intermediate
107: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDotRealPart()
108: @*/
109: PetscErrorCode VecDot(Vec x,Vec y,PetscScalar *val)
110: {
120: VecCheckSameSize(x,1,y,2);
122: PetscLogEventBegin(VEC_Dot,x,y,0,0);
123: (*x->ops->dot)(x,y,val);
124: PetscLogEventEnd(VEC_Dot,x,y,0,0);
125: return(0);
126: }
128: /*@
129: VecDotRealPart - Computes the real part of the vector dot product.
131: Collective on Vec
133: Input Parameters:
134: . x, y - the vectors
136: Output Parameter:
137: . val - the real part of the dot product;
139: Performance Issues:
140: $ per-processor memory bandwidth
141: $ interprocessor latency
142: $ work load inbalance that causes certain processes to arrive much earlier than others
144: Notes for Users of Complex Numbers:
145: See VecDot() for more details on the definition of the dot product for complex numbers
147: For real numbers this returns the same value as VecDot()
149: 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
150: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
152: Developer Note: This is not currently optimized to compute only the real part of the dot product.
154: Level: intermediate
157: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDot(), VecDotNorm2()
158: @*/
159: PetscErrorCode VecDotRealPart(Vec x,Vec y,PetscReal *val)
160: {
162: PetscScalar fdot;
165: VecDot(x,y,&fdot);
166: *val = PetscRealPart(fdot);
167: return(0);
168: }
170: /*@
171: VecNorm - Computes the vector norm.
173: Collective on Vec
175: Input Parameters:
176: + x - the vector
177: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
178: NORM_1_AND_2, which computes both norms and stores them
179: in a two element array.
181: Output Parameter:
182: . val - the norm
184: Notes:
185: $ NORM_1 denotes sum_i |x_i|
186: $ NORM_2 denotes sqrt(sum_i |x_i|^2)
187: $ NORM_INFINITY denotes max_i |x_i|
189: For complex numbers NORM_1 will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
190: norm of the absolutely values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
191: the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
192: people expect the former.
194: Level: intermediate
196: Performance Issues:
197: $ per-processor memory bandwidth
198: $ interprocessor latency
199: $ work load inbalance that causes certain processes to arrive much earlier than others
202: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNormAvailable(),
203: VecNormBegin(), VecNormEnd()
205: @*/
207: PetscErrorCode VecNorm(Vec x,NormType type,PetscReal *val)
208: {
209: PetscBool flg;
217: /*
218: * Cached data?
219: */
220: if (type!=NORM_1_AND_2) {
221: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,flg);
222: if (flg) return(0);
223: }
224: PetscLogEventBegin(VEC_Norm,x,0,0,0);
225: (*x->ops->norm)(x,type,val);
226: PetscLogEventEnd(VEC_Norm,x,0,0,0);
227: if (type!=NORM_1_AND_2) {
228: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[type],*val);
229: }
230: return(0);
231: }
233: /*@
234: VecNormAvailable - Returns the vector norm if it is already known.
236: Not Collective
238: Input Parameters:
239: + x - the vector
240: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
241: NORM_1_AND_2, which computes both norms and stores them
242: in a two element array.
244: Output Parameter:
245: + available - PETSC_TRUE if the val returned is valid
246: - val - the norm
248: Notes:
249: $ NORM_1 denotes sum_i |x_i|
250: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
251: $ NORM_INFINITY denotes max_i |x_i|
253: Level: intermediate
255: Performance Issues:
256: $ per-processor memory bandwidth
257: $ interprocessor latency
258: $ work load inbalance that causes certain processes to arrive much earlier than others
260: Compile Option:
261: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
262: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
263: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
266: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNorm()
267: VecNormBegin(), VecNormEnd()
269: @*/
270: PetscErrorCode VecNormAvailable(Vec x,NormType type,PetscBool *available,PetscReal *val)
271: {
279: *available = PETSC_FALSE;
280: if (type!=NORM_1_AND_2) {
281: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,*available);
282: }
283: return(0);
284: }
286: /*@
287: VecNormalize - Normalizes a vector by 2-norm.
289: Collective on Vec
291: Input Parameters:
292: + x - the vector
294: Output Parameter:
295: . x - the normalized vector
296: - val - the vector norm before normalization
298: Level: intermediate
301: @*/
302: PetscErrorCode VecNormalize(Vec x,PetscReal *val)
303: {
305: PetscReal norm;
310: PetscLogEventBegin(VEC_Normalize,x,0,0,0);
311: VecNorm(x,NORM_2,&norm);
312: if (norm == 0.0) {
313: PetscInfo(x,"Vector of zero norm can not be normalized; Returning only the zero norm\n");
314: } else if (norm != 1.0) {
315: PetscScalar tmp = 1.0/norm;
316: VecScale(x,tmp);
317: }
318: if (val) *val = norm;
319: PetscLogEventEnd(VEC_Normalize,x,0,0,0);
320: return(0);
321: }
323: /*@C
324: VecMax - Determines the vector component with maximum real part and its location.
326: Collective on Vec
328: Input Parameter:
329: . x - the vector
331: Output Parameters:
332: + p - the location of val (pass NULL if you don't want this)
333: - val - the maximum component
335: Notes:
336: Returns the value PETSC_MIN_REAL and p = -1 if the vector is of length 0.
338: Returns the smallest index with the maximum value
339: Level: intermediate
342: .seealso: VecNorm(), VecMin()
343: @*/
344: PetscErrorCode VecMax(Vec x,PetscInt *p,PetscReal *val)
345: {
352: PetscLogEventBegin(VEC_Max,x,0,0,0);
353: (*x->ops->max)(x,p,val);
354: PetscLogEventEnd(VEC_Max,x,0,0,0);
355: return(0);
356: }
358: /*@C
359: VecMin - Determines the vector component with minimum real part and its location.
361: Collective on Vec
363: Input Parameters:
364: . x - the vector
366: Output Parameter:
367: + p - the location of val (pass NULL if you don't want this location)
368: - val - the minimum component
370: Level: intermediate
372: Notes:
373: Returns the value PETSC_MAX_REAL and p = -1 if the vector is of length 0.
375: This returns the smallest index with the minumum value
378: .seealso: VecMax()
379: @*/
380: PetscErrorCode VecMin(Vec x,PetscInt *p,PetscReal *val)
381: {
388: PetscLogEventBegin(VEC_Min,x,0,0,0);
389: (*x->ops->min)(x,p,val);
390: PetscLogEventEnd(VEC_Min,x,0,0,0);
391: return(0);
392: }
394: /*@
395: VecTDot - Computes an indefinite vector dot product. That is, this
396: routine does NOT use the complex conjugate.
398: Collective on Vec
400: Input Parameters:
401: . x, y - the vectors
403: Output Parameter:
404: . val - the dot product
406: Notes for Users of Complex Numbers:
407: For complex vectors, VecTDot() computes the indefinite form
408: $ val = (x,y) = y^T x,
409: where y^T denotes the transpose of y.
411: Use VecDot() for the inner product
412: $ val = (x,y) = y^H x,
413: where y^H denotes the conjugate transpose of y.
415: Level: intermediate
417: .seealso: VecDot(), VecMTDot()
418: @*/
419: PetscErrorCode VecTDot(Vec x,Vec y,PetscScalar *val)
420: {
430: VecCheckSameSize(x,1,y,2);
432: PetscLogEventBegin(VEC_TDot,x,y,0,0);
433: (*x->ops->tdot)(x,y,val);
434: PetscLogEventEnd(VEC_TDot,x,y,0,0);
435: return(0);
436: }
438: /*@
439: VecScale - Scales a vector.
441: Not collective on Vec
443: Input Parameters:
444: + x - the vector
445: - alpha - the scalar
447: Output Parameter:
448: . x - the scaled vector
450: Note:
451: For a vector with n components, VecScale() computes
452: $ x[i] = alpha * x[i], for i=1,...,n.
454: Level: intermediate
457: @*/
458: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
459: {
460: PetscReal norms[4] = {0.0,0.0,0.0, 0.0};
461: PetscBool flgs[4];
463: PetscInt i;
468: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
469: PetscLogEventBegin(VEC_Scale,x,0,0,0);
470: if (alpha != (PetscScalar)1.0) {
471: VecSetErrorIfLocked(x,1);
472: /* get current stashed norms */
473: for (i=0; i<4; i++) {
474: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[i],norms[i],flgs[i]);
475: }
476: (*x->ops->scale)(x,alpha);
477: PetscObjectStateIncrease((PetscObject)x);
478: /* put the scaled stashed norms back into the Vec */
479: for (i=0; i<4; i++) {
480: if (flgs[i]) {
481: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[i],PetscAbsScalar(alpha)*norms[i]);
482: }
483: }
484: }
485: PetscLogEventEnd(VEC_Scale,x,0,0,0);
486: return(0);
487: }
489: /*@
490: VecSet - Sets all components of a vector to a single scalar value.
492: Logically Collective on Vec
494: Input Parameters:
495: + x - the vector
496: - alpha - the scalar
498: Output Parameter:
499: . x - the vector
501: Note:
502: For a vector of dimension n, VecSet() computes
503: $ x[i] = alpha, for i=1,...,n,
504: so that all vector entries then equal the identical
505: scalar value, alpha. Use the more general routine
506: VecSetValues() to set different vector entries.
508: You CANNOT call this after you have called VecSetValues() but before you call
509: VecAssemblyBegin/End().
511: Level: beginner
513: .seealso VecSetValues(), VecSetValuesBlocked(), VecSetRandom()
515: @*/
516: PetscErrorCode VecSet(Vec x,PetscScalar alpha)
517: {
518: PetscReal val;
524: 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()");
526: VecSetErrorIfLocked(x,1);
528: PetscLogEventBegin(VEC_Set,x,0,0,0);
529: (*x->ops->set)(x,alpha);
530: PetscLogEventEnd(VEC_Set,x,0,0,0);
531: PetscObjectStateIncrease((PetscObject)x);
533: /* norms can be simply set (if |alpha|*N not too large) */
534: val = PetscAbsScalar(alpha);
535: if (x->map->N == 0) {
536: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],0.0l);
537: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],0.0);
538: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],0.0);
539: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],0.0);
540: } else if (val > PETSC_MAX_REAL/x->map->N) {
541: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
542: } else {
543: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_1],x->map->N * val);
544: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_INFINITY],val);
545: val = PetscSqrtReal((PetscReal)x->map->N) * val;
546: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_2],val);
547: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[NORM_FROBENIUS],val);
548: }
549: return(0);
550: }
553: /*@
554: VecAXPY - Computes y = alpha x + y.
556: Logically Collective on Vec
558: Input Parameters:
559: + alpha - the scalar
560: - x, y - the vectors
562: Output Parameter:
563: . y - output vector
565: Level: intermediate
567: Notes:
568: x and y MUST be different vectors
569: This routine is optimized for alpha of 0.0, otherwise it calls the BLAS routine
571: $ VecAXPY(y,alpha,x) y = alpha x + y
572: $ VecAYPX(y,beta,x) y = x + beta y
573: $ VecAXPBY(y,alpha,beta,x) y = alpha x + beta y
574: $ VecWAXPY(w,alpha,x,y) w = alpha x + y
575: $ VecAXPBYPCZ(w,alpha,beta,gamma,x,y) z = alpha x + beta y + gamma z
576: $ VecMAXPY(y,nv,alpha[],x[]) y = sum alpha[i] x[i] + y
579: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPBYPCZ(), VecAXPBY()
580: @*/
581: PetscErrorCode VecAXPY(Vec y,PetscScalar alpha,Vec x)
582: {
591: VecCheckSameSize(x,1,y,3);
592: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
594: if (alpha == (PetscScalar)0.0) return(0);
595: VecSetErrorIfLocked(y,1);
597: VecLockReadPush(x);
598: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
599: (*y->ops->axpy)(y,alpha,x);
600: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
601: VecLockReadPop(x);
602: PetscObjectStateIncrease((PetscObject)y);
603: return(0);
604: }
606: /*@
607: VecAXPBY - Computes y = alpha x + beta y.
609: Logically Collective on Vec
611: Input Parameters:
612: + alpha,beta - the scalars
613: - x, y - the vectors
615: Output Parameter:
616: . y - output vector
618: Level: intermediate
620: Notes:
621: x and y MUST be different vectors
622: The implementation is optimized for alpha and/or beta values of 0.0 and 1.0
625: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ()
626: @*/
627: PetscErrorCode VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
628: {
637: VecCheckSameSize(y,1,x,4);
638: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
641: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)1.0) return(0);
642: VecSetErrorIfLocked(y,1);
643: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
644: (*y->ops->axpby)(y,alpha,beta,x);
645: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
646: PetscObjectStateIncrease((PetscObject)y);
647: return(0);
648: }
650: /*@
651: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
653: Logically Collective on Vec
655: Input Parameters:
656: + alpha,beta, gamma - the scalars
657: - x, y, z - the vectors
659: Output Parameter:
660: . z - output vector
662: Level: intermediate
664: Notes:
665: x, y and z must be different vectors
666: The implementation is optimized for alpha of 1.0 and gamma of 1.0 or 0.0
669: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBY()
670: @*/
671: PetscErrorCode VecAXPBYPCZ(Vec z,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec x,Vec y)
672: {
684: VecCheckSameSize(x,1,y,5);
685: VecCheckSameSize(x,1,z,6);
686: if (x == y || x == z) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
687: if (y == z) SETERRQ(PetscObjectComm((PetscObject)y),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
691: if (alpha == (PetscScalar)0.0 && beta == (PetscScalar)0.0 && gamma == (PetscScalar)1.0) return(0);
692: VecSetErrorIfLocked(z,1);
694: PetscLogEventBegin(VEC_AXPBYPCZ,x,y,z,0);
695: (*y->ops->axpbypcz)(z,alpha,beta,gamma,x,y);
696: PetscLogEventEnd(VEC_AXPBYPCZ,x,y,z,0);
697: PetscObjectStateIncrease((PetscObject)z);
698: return(0);
699: }
701: /*@
702: VecAYPX - Computes y = x + beta y.
704: Logically Collective on Vec
706: Input Parameters:
707: + beta - the scalar
708: - x, y - the vectors
710: Output Parameter:
711: . y - output vector
713: Level: intermediate
715: Notes:
716: x and y MUST be different vectors
717: The implementation is optimized for beta of -1.0, 0.0, and 1.0
720: .seealso: VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
721: @*/
722: PetscErrorCode VecAYPX(Vec y,PetscScalar beta,Vec x)
723: {
732: VecCheckSameSize(x,1,y,3);
733: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
735: VecSetErrorIfLocked(y,1);
737: PetscLogEventBegin(VEC_AYPX,x,y,0,0);
738: (*y->ops->aypx)(y,beta,x);
739: PetscLogEventEnd(VEC_AYPX,x,y,0,0);
740: PetscObjectStateIncrease((PetscObject)y);
741: return(0);
742: }
745: /*@
746: VecWAXPY - Computes w = alpha x + y.
748: Logically Collective on Vec
750: Input Parameters:
751: + alpha - the scalar
752: - x, y - the vectors
754: Output Parameter:
755: . w - the result
757: Level: intermediate
759: Notes:
760: w cannot be either x or y, but x and y can be the same
761: The implementation is optimzed for alpha of -1.0, 0.0, and 1.0
764: .seealso: VecAXPY(), VecAYPX(), VecAXPBY(), VecMAXPY(), VecAXPBYPCZ()
765: @*/
766: PetscErrorCode VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
767: {
779: VecCheckSameSize(x,3,y,4);
780: VecCheckSameSize(x,3,w,1);
781: if (w == y) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector y, suggest VecAXPY()");
782: if (w == x) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector x, suggest VecAYPX()");
784: VecSetErrorIfLocked(w,1);
786: PetscLogEventBegin(VEC_WAXPY,x,y,w,0);
787: (*w->ops->waxpy)(w,alpha,x,y);
788: PetscLogEventEnd(VEC_WAXPY,x,y,w,0);
789: PetscObjectStateIncrease((PetscObject)w);
790: return(0);
791: }
794: /*@C
795: VecSetValues - Inserts or adds values into certain locations of a vector.
797: Not Collective
799: Input Parameters:
800: + x - vector to insert in
801: . ni - number of elements to add
802: . ix - indices where to add
803: . y - array of values
804: - iora - either INSERT_VALUES or ADD_VALUES, where
805: ADD_VALUES adds values to any existing entries, and
806: INSERT_VALUES replaces existing entries with new values
808: Notes:
809: VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.
811: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
812: options cannot be mixed without intervening calls to the assembly
813: routines.
815: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
816: MUST be called after all calls to VecSetValues() have been completed.
818: VecSetValues() uses 0-based indices in Fortran as well as in C.
820: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
821: negative indices may be passed in ix. These rows are
822: simply ignored. This allows easily inserting element load matrices
823: with homogeneous Dirchlet boundary conditions that you don't want represented
824: in the vector.
826: Level: beginner
828: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
829: VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
830: @*/
831: PetscErrorCode VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
832: {
837: if (!ni) return(0);
842: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
843: (*x->ops->setvalues)(x,ni,ix,y,iora);
844: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
845: PetscObjectStateIncrease((PetscObject)x);
846: return(0);
847: }
849: /*@C
850: VecGetValues - Gets values from certain locations of a vector. Currently
851: can only get values on the same processor
853: Not Collective
855: Input Parameters:
856: + x - vector to get values from
857: . ni - number of elements to get
858: - ix - indices where to get them from (in global 1d numbering)
860: Output Parameter:
861: . y - array of values
863: Notes:
864: The user provides the allocated array y; it is NOT allocated in this routine
866: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
868: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
870: VecGetValues() uses 0-based indices in Fortran as well as in C.
872: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
873: negative indices may be passed in ix. These rows are
874: simply ignored.
876: Level: beginner
878: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues()
879: @*/
880: PetscErrorCode VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
881: {
886: if (!ni) return(0);
890: (*x->ops->getvalues)(x,ni,ix,y);
891: return(0);
892: }
894: /*@C
895: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
897: Not Collective
899: Input Parameters:
900: + x - vector to insert in
901: . ni - number of blocks to add
902: . ix - indices where to add in block count, rather than element count
903: . y - array of values
904: - iora - either INSERT_VALUES or ADD_VALUES, where
905: ADD_VALUES adds values to any existing entries, and
906: INSERT_VALUES replaces existing entries with new values
908: Notes:
909: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
910: for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
912: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
913: options cannot be mixed without intervening calls to the assembly
914: routines.
916: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
917: MUST be called after all calls to VecSetValuesBlocked() have been completed.
919: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
921: Negative indices may be passed in ix, these rows are
922: simply ignored. This allows easily inserting element load matrices
923: with homogeneous Dirchlet boundary conditions that you don't want represented
924: in the vector.
926: Level: intermediate
928: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
929: VecSetValues()
930: @*/
931: PetscErrorCode VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
932: {
937: if (!ni) return(0);
942: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
943: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
944: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
945: PetscObjectStateIncrease((PetscObject)x);
946: return(0);
947: }
950: /*@C
951: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
952: using a local ordering of the nodes.
954: Not Collective
956: Input Parameters:
957: + x - vector to insert in
958: . ni - number of elements to add
959: . ix - indices where to add
960: . y - array of values
961: - iora - either INSERT_VALUES or ADD_VALUES, where
962: ADD_VALUES adds values to any existing entries, and
963: INSERT_VALUES replaces existing entries with new values
965: Level: intermediate
967: Notes:
968: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
970: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
971: options cannot be mixed without intervening calls to the assembly
972: routines.
974: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
975: MUST be called after all calls to VecSetValuesLocal() have been completed.
977: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
979: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
980: VecSetValuesBlockedLocal()
981: @*/
982: PetscErrorCode VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
983: {
985: PetscInt lixp[128],*lix = lixp;
989: if (!ni) return(0);
994: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
995: if (!x->ops->setvalueslocal) {
996: if (!x->map->mapping) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMapping()");
997: if (ni > 128) {
998: PetscMalloc1(ni,&lix);
999: }
1000: ISLocalToGlobalMappingApply(x->map->mapping,ni,(PetscInt*)ix,lix);
1001: (*x->ops->setvalues)(x,ni,lix,y,iora);
1002: if (ni > 128) {
1003: PetscFree(lix);
1004: }
1005: } else {
1006: (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
1007: }
1008: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1009: PetscObjectStateIncrease((PetscObject)x);
1010: return(0);
1011: }
1013: /*@
1014: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
1015: using a local ordering of the nodes.
1017: Not Collective
1019: Input Parameters:
1020: + x - vector to insert in
1021: . ni - number of blocks to add
1022: . ix - indices where to add in block count, not element count
1023: . y - array of values
1024: - iora - either INSERT_VALUES or ADD_VALUES, where
1025: ADD_VALUES adds values to any existing entries, and
1026: INSERT_VALUES replaces existing entries with new values
1028: Level: intermediate
1030: Notes:
1031: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
1032: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
1034: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
1035: options cannot be mixed without intervening calls to the assembly
1036: routines.
1038: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1039: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
1041: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
1044: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(),
1045: VecSetLocalToGlobalMapping()
1046: @*/
1047: PetscErrorCode VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1048: {
1050: PetscInt lixp[128],*lix = lixp;
1054: if (!ni) return(0);
1058: if (ni > 128) {
1059: PetscMalloc1(ni,&lix);
1060: }
1062: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1063: ISLocalToGlobalMappingApplyBlock(x->map->mapping,ni,(PetscInt*)ix,lix);
1064: (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
1065: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1066: if (ni > 128) {
1067: PetscFree(lix);
1068: }
1069: PetscObjectStateIncrease((PetscObject)x);
1070: return(0);
1071: }
1073: /*@
1074: VecMTDot - Computes indefinite vector multiple dot products.
1075: That is, it does NOT use the complex conjugate.
1077: Collective on Vec
1079: Input Parameters:
1080: + x - one vector
1081: . nv - number of vectors
1082: - y - array of vectors. Note that vectors are pointers
1084: Output Parameter:
1085: . val - array of the dot products
1087: Notes for Users of Complex Numbers:
1088: For complex vectors, VecMTDot() computes the indefinite form
1089: $ val = (x,y) = y^T x,
1090: where y^T denotes the transpose of y.
1092: Use VecMDot() for the inner product
1093: $ val = (x,y) = y^H x,
1094: where y^H denotes the conjugate transpose of y.
1096: Level: intermediate
1099: .seealso: VecMDot(), VecTDot()
1100: @*/
1101: PetscErrorCode VecMTDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1102: {
1108: if (!nv) return(0);
1115: VecCheckSameSize(x,1,*y,3);
1117: PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
1118: (*x->ops->mtdot)(x,nv,y,val);
1119: PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
1120: return(0);
1121: }
1123: /*@
1124: VecMDot - Computes vector multiple dot products.
1126: Collective on Vec
1128: Input Parameters:
1129: + x - one vector
1130: . nv - number of vectors
1131: - y - array of vectors.
1133: Output Parameter:
1134: . val - array of the dot products (does not allocate the array)
1136: Notes for Users of Complex Numbers:
1137: For complex vectors, VecMDot() computes
1138: $ val = (x,y) = y^H x,
1139: where y^H denotes the conjugate transpose of y.
1141: Use VecMTDot() for the indefinite form
1142: $ val = (x,y) = y^T x,
1143: where y^T denotes the transpose of y.
1145: Level: intermediate
1148: .seealso: VecMTDot(), VecDot()
1149: @*/
1150: PetscErrorCode VecMDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1151: {
1157: if (!nv) return(0);
1158: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1165: VecCheckSameSize(x,1,*y,3);
1167: PetscLogEventBegin(VEC_MDot,x,*y,0,0);
1168: (*x->ops->mdot)(x,nv,y,val);
1169: PetscLogEventEnd(VEC_MDot,x,*y,0,0);
1170: return(0);
1171: }
1173: /*@
1174: VecMAXPY - Computes y = y + sum alpha[i] x[i]
1176: Logically Collective on Vec
1178: Input Parameters:
1179: + nv - number of scalars and x-vectors
1180: . alpha - array of scalars
1181: . y - one vector
1182: - x - array of vectors
1184: Level: intermediate
1186: Notes:
1187: y cannot be any of the x vectors
1189: .seealso: VecAYPX(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
1190: @*/
1191: PetscErrorCode VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec x[])
1192: {
1194: PetscInt i;
1195: PetscBool nonzero;
1200: if (!nv) return(0);
1201: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1208: VecCheckSameSize(y,1,*x,4);
1210: for (i=0, nonzero = PETSC_FALSE; i<nv && !nonzero; i++) nonzero = (PetscBool)(nonzero || alpha[i] != (PetscScalar)0.0);
1211: if (!nonzero) return(0);
1212: VecSetErrorIfLocked(y,1);
1213: PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
1214: (*y->ops->maxpy)(y,nv,alpha,x);
1215: PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
1216: PetscObjectStateIncrease((PetscObject)y);
1217: return(0);
1218: }
1220: /*@
1221: VecGetSubVector - Gets a vector representing part of another vector
1223: Collective on IS
1225: Input Arguments:
1226: + X - vector from which to extract a subvector
1227: - is - index set representing portion of X to extract
1229: Output Arguments:
1230: . Y - subvector corresponding to is
1232: Level: advanced
1234: Notes:
1235: The subvector Y should be returned with VecRestoreSubVector().
1237: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1238: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1240: .seealso: MatCreateSubMatrix()
1241: @*/
1242: PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
1243: {
1244: PetscErrorCode ierr;
1245: Vec Z;
1251: if (X->ops->getsubvector) {
1252: (*X->ops->getsubvector)(X,is,&Z);
1253: } else { /* Default implementation currently does no caching */
1254: PetscInt gstart,gend,start;
1255: PetscBool contiguous,gcontiguous;
1256: VecGetOwnershipRange(X,&gstart,&gend);
1257: ISContiguousLocal(is,gstart,gend,&start,&contiguous);
1258: MPIU_Allreduce(&contiguous,&gcontiguous,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1259: if (gcontiguous) { /* We can do a no-copy implementation */
1260: PetscInt n,N,bs;
1261: PetscInt state;
1263: ISGetSize(is,&N);
1264: ISGetLocalSize(is,&n);
1265: VecGetBlockSize(X,&bs);
1266: if (n%bs || bs == 1 || !n) bs = -1; /* Do not decide block size if we do not have to */
1267: VecLockGet(X,&state);
1268: if (state) {
1269: const PetscScalar *x;
1270: VecGetArrayRead(X,&x);
1271: VecCreate(PetscObjectComm((PetscObject)X),&Z);
1272: VecSetType(Z,((PetscObject)X)->type_name);
1273: VecSetSizes(Z,n,N);
1274: VecSetBlockSize(Z,bs);
1275: VecPlaceArray(Z,(PetscScalar*)x+start);
1276: VecLockReadPush(Z);
1277: VecRestoreArrayRead(X,&x);
1278: } else {
1279: PetscScalar *x;
1280: VecGetArray(X,&x);
1281: VecCreate(PetscObjectComm((PetscObject)X),&Z);
1282: VecSetType(Z,((PetscObject)X)->type_name);
1283: VecSetSizes(Z,n,N);
1284: VecSetBlockSize(Z,bs);
1285: VecPlaceArray(Z,(PetscScalar*)x+start);
1286: VecRestoreArray(X,&x);
1287: }
1288: } else { /* Have to create a scatter and do a copy */
1289: VecScatter scatter;
1290: PetscInt n,N;
1291: ISGetLocalSize(is,&n);
1292: ISGetSize(is,&N);
1293: VecCreate(PetscObjectComm((PetscObject)is),&Z);
1294: VecSetSizes(Z,n,N);
1295: VecSetType(Z,((PetscObject)X)->type_name);
1296: VecScatterCreate(X,is,Z,NULL,&scatter);
1297: VecScatterBegin(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1298: VecScatterEnd(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1299: PetscObjectCompose((PetscObject)Z,"VecGetSubVector_Scatter",(PetscObject)scatter);
1300: VecScatterDestroy(&scatter);
1301: }
1302: }
1303: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1304: if (VecGetSubVectorSavedStateId < 0) {PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);}
1305: PetscObjectComposedDataSetInt((PetscObject)Z,VecGetSubVectorSavedStateId,1);
1306: *Y = Z;
1307: return(0);
1308: }
1310: /*@
1311: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1313: Collective on IS
1315: Input Arguments:
1316: + X - vector from which subvector was obtained
1317: . is - index set representing the subset of X
1318: - Y - subvector being restored
1320: Level: advanced
1322: .seealso: VecGetSubVector()
1323: @*/
1324: PetscErrorCode VecRestoreSubVector(Vec X,IS is,Vec *Y)
1325: {
1333: if (X->ops->restoresubvector) {
1334: (*X->ops->restoresubvector)(X,is,Y);
1335: } else {
1336: PETSC_UNUSED PetscObjectState dummystate = 0;
1337: PetscBool valid;
1338: PetscObjectComposedDataGetInt((PetscObject)*Y,VecGetSubVectorSavedStateId,dummystate,valid);
1339: if (!valid) {
1340: VecScatter scatter;
1342: PetscObjectQuery((PetscObject)*Y,"VecGetSubVector_Scatter",(PetscObject*)&scatter);
1343: if (scatter) {
1344: VecScatterBegin(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1345: VecScatterEnd(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1346: }
1347: }
1348: VecDestroy(Y);
1349: }
1350: return(0);
1351: }
1353: /*@
1354: VecGetLocalVectorRead - Maps the local portion of a vector into a
1355: vector. You must call VecRestoreLocalVectorRead() when the local
1356: vector is no longer needed.
1358: Not collective.
1360: Input parameter:
1361: . v - The vector for which the local vector is desired.
1363: Output parameter:
1364: . w - Upon exit this contains the local vector.
1366: Level: beginner
1368: Notes:
1369: This function is similar to VecGetArrayRead() which maps the local
1370: portion into a raw pointer. VecGetLocalVectorRead() is usually
1371: almost as efficient as VecGetArrayRead() but in certain circumstances
1372: VecGetLocalVectorRead() can be much more efficient than
1373: VecGetArrayRead(). This is because the construction of a contiguous
1374: array representing the vector data required by VecGetArrayRead() can
1375: be an expensive operation for certain vector types. For example, for
1376: GPU vectors VecGetArrayRead() requires that the data between device
1377: and host is synchronized.
1379: Unlike VecGetLocalVector(), this routine is not collective and
1380: preserves cached information.
1382: .seealso: VecRestoreLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1383: @*/
1384: PetscErrorCode VecGetLocalVectorRead(Vec v,Vec w)
1385: {
1387: PetscScalar *a;
1392: VecCheckSameLocalSize(v,1,w,2);
1393: if (v->ops->getlocalvectorread) {
1394: (*v->ops->getlocalvectorread)(v,w);
1395: } else {
1396: VecGetArrayRead(v,(const PetscScalar**)&a);
1397: VecPlaceArray(w,a);
1398: }
1399: return(0);
1400: }
1402: /*@
1403: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1404: previously mapped into a vector using VecGetLocalVectorRead().
1406: Not collective.
1408: Input parameter:
1409: + v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1410: - w - The vector into which the local portion of v was mapped.
1412: Level: beginner
1414: .seealso: VecGetLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1415: @*/
1416: PetscErrorCode VecRestoreLocalVectorRead(Vec v,Vec w)
1417: {
1419: PetscScalar *a;
1424: if (v->ops->restorelocalvectorread) {
1425: (*v->ops->restorelocalvectorread)(v,w);
1426: } else {
1427: VecGetArrayRead(w,(const PetscScalar**)&a);
1428: VecRestoreArrayRead(v,(const PetscScalar**)&a);
1429: VecResetArray(w);
1430: }
1431: return(0);
1432: }
1434: /*@
1435: VecGetLocalVector - Maps the local portion of a vector into a
1436: vector.
1438: Collective on v, not collective on w.
1440: Input parameter:
1441: . v - The vector for which the local vector is desired.
1443: Output parameter:
1444: . w - Upon exit this contains the local vector.
1446: Level: beginner
1448: Notes:
1449: This function is similar to VecGetArray() which maps the local
1450: portion into a raw pointer. VecGetLocalVector() is usually about as
1451: efficient as VecGetArray() but in certain circumstances
1452: VecGetLocalVector() can be much more efficient than VecGetArray().
1453: This is because the construction of a contiguous array representing
1454: the vector data required by VecGetArray() can be an expensive
1455: operation for certain vector types. For example, for GPU vectors
1456: VecGetArray() requires that the data between device and host is
1457: synchronized.
1459: .seealso: VecRestoreLocalVector(), VecGetLocalVectorRead(), VecGetArrayRead(), VecGetArray()
1460: @*/
1461: PetscErrorCode VecGetLocalVector(Vec v,Vec w)
1462: {
1464: PetscScalar *a;
1469: VecCheckSameLocalSize(v,1,w,2);
1470: if (v->ops->getlocalvector) {
1471: (*v->ops->getlocalvector)(v,w);
1472: } else {
1473: VecGetArray(v,&a);
1474: VecPlaceArray(w,a);
1475: }
1476: return(0);
1477: }
1479: /*@
1480: VecRestoreLocalVector - Unmaps the local portion of a vector
1481: previously mapped into a vector using VecGetLocalVector().
1483: Logically collective.
1485: Input parameter:
1486: + v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1487: - w - The vector into which the local portion of v was mapped.
1489: Level: beginner
1491: .seealso: VecGetLocalVector(), VecGetLocalVectorRead(), VecRestoreLocalVectorRead(), LocalVectorRead(), VecGetArrayRead(), VecGetArray()
1492: @*/
1493: PetscErrorCode VecRestoreLocalVector(Vec v,Vec w)
1494: {
1496: PetscScalar *a;
1501: if (v->ops->restorelocalvector) {
1502: (*v->ops->restorelocalvector)(v,w);
1503: } else {
1504: VecGetArray(w,&a);
1505: VecRestoreArray(v,&a);
1506: VecResetArray(w);
1507: }
1508: return(0);
1509: }
1511: /*@C
1512: VecGetArray - Returns a pointer to a contiguous array that contains this
1513: processor's portion of the vector data. For the standard PETSc
1514: vectors, VecGetArray() returns a pointer to the local data array and
1515: does not use any copies. If the underlying vector data is not stored
1516: in a contiguous array this routine will copy the data to a contiguous
1517: array and return a pointer to that. You MUST call VecRestoreArray()
1518: when you no longer need access to the array.
1520: Logically Collective on Vec
1522: Input Parameter:
1523: . x - the vector
1525: Output Parameter:
1526: . a - location to put pointer to the array
1528: Fortran Note:
1529: This routine is used differently from Fortran 77
1530: $ Vec x
1531: $ PetscScalar x_array(1)
1532: $ PetscOffset i_x
1533: $ PetscErrorCode ierr
1534: $ call VecGetArray(x,x_array,i_x,ierr)
1535: $
1536: $ Access first local entry in vector with
1537: $ value = x_array(i_x + 1)
1538: $
1539: $ ...... other code
1540: $ call VecRestoreArray(x,x_array,i_x,ierr)
1541: For Fortran 90 see VecGetArrayF90()
1543: See the Fortran chapter of the users manual and
1544: petsc/src/snes/tutorials/ex5f.F for details.
1546: Level: beginner
1548: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1549: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
1550: @*/
1551: PetscErrorCode VecGetArray(Vec x,PetscScalar **a)
1552: {
1554: #if defined(PETSC_HAVE_VIENNACL)
1555: PetscBool is_viennacltype = PETSC_FALSE;
1556: #endif
1560: VecSetErrorIfLocked(x,1);
1561: if (x->petscnative) {
1562: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1563: if (x->offloadmask == PETSC_OFFLOAD_GPU) {
1564: #if defined(PETSC_HAVE_VIENNACL)
1565: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1566: if (is_viennacltype) {
1567: VecViennaCLCopyFromGPU(x);
1568: } else
1569: #endif
1570: {
1571: #if defined(PETSC_HAVE_CUDA)
1572: VecCUDACopyFromGPU(x);
1573: #endif
1574: }
1575: } else if (x->offloadmask == PETSC_OFFLOAD_UNALLOCATED) {
1576: #if defined(PETSC_HAVE_VIENNACL)
1577: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1578: if (is_viennacltype) {
1579: VecViennaCLAllocateCheckHost(x);
1580: } else
1581: #endif
1582: {
1583: #if defined(PETSC_HAVE_CUDA)
1584: VecCUDAAllocateCheckHost(x);
1585: #endif
1586: }
1587: }
1588: #endif
1589: *a = *((PetscScalar**)x->data);
1590: } else {
1591: if (x->ops->getarray) {
1592: (*x->ops->getarray)(x,a);
1593: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array for vector type \"%s\"",((PetscObject)x)->type_name);
1594: }
1595: return(0);
1596: }
1598: /*@C
1599: VecGetArrayInPlace - Like VecGetArray(), but if this is a CUDA vector and it is currently offloaded to GPU,
1600: the returned pointer will be a GPU pointer to the GPU memory that contains this processor's portion of the
1601: vector data. Otherwise, it functions as VecGetArray().
1603: Logically Collective on Vec
1605: Input Parameter:
1606: . x - the vector
1608: Output Parameter:
1609: . a - location to put pointer to the array
1611: Level: beginner
1613: .seealso: VecRestoreArrayInPlace(), VecRestoreArrayInPlace(), VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(),
1614: VecPlaceArray(), VecGetArray2d(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
1615: @*/
1616: PetscErrorCode VecGetArrayInPlace(Vec x,PetscScalar **a)
1617: {
1622: VecSetErrorIfLocked(x,1);
1624: #if defined(PETSC_HAVE_CUDA)
1625: if (x->petscnative && (x->offloadmask & PETSC_OFFLOAD_GPU)) { /* Prefer working on GPU when offloadmask is PETSC_OFFLOAD_BOTH */
1626: PetscBool is_cudatype = PETSC_FALSE;
1627: PetscObjectTypeCompareAny((PetscObject)x,&is_cudatype,VECSEQCUDA,VECMPICUDA,VECCUDA,"");
1628: if (is_cudatype) {
1629: VecCUDAGetArray(x,a);
1630: x->offloadmask = PETSC_OFFLOAD_GPU; /* Change the mask once GPU gets write access, don't wait until restore array */
1631: return(0);
1632: }
1633: }
1634: #endif
1635: VecGetArray(x,a);
1636: return(0);
1637: }
1639: /*@C
1640: VecGetArrayWrite - Returns a pointer to a contiguous array that WILL contains this
1641: processor's portion of the vector data. The values in this array are NOT valid, the routine calling this
1642: routine is responsible for putting values into the array; any values it does not set will be invalid
1644: Logically Collective on Vec
1646: Input Parameter:
1647: . x - the vector
1649: Output Parameter:
1650: . a - location to put pointer to the array
1652: Level: intermediate
1654: This is for vectors associate with GPUs, the vector is not copied up before giving access. If you need correct
1655: values in the array use VecGetArray()
1657: Concepts: vector^accessing local values
1659: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1660: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArray(), VecRestoreArrayWrite()
1661: @*/
1662: PetscErrorCode VecGetArrayWrite(Vec x,PetscScalar **a)
1663: {
1668: VecSetErrorIfLocked(x,1);
1669: if (!x->ops->getarraywrite) {
1670: VecGetArray(x,a);
1671: } else {
1672: (*x->ops->getarraywrite)(x,a);
1673: }
1674: return(0);
1675: }
1677: /*@C
1678: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1680: Not Collective
1682: Input Parameters:
1683: . x - the vector
1685: Output Parameter:
1686: . a - the array
1688: Level: beginner
1690: Notes:
1691: The array must be returned using a matching call to VecRestoreArrayRead().
1693: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1695: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1696: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1697: only one copy is performed when this routine is called multiple times in sequence.
1699: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1700: @*/
1701: PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar **a)
1702: {
1704: #if defined(PETSC_HAVE_VIENNACL)
1705: PetscBool is_viennacltype = PETSC_FALSE;
1706: #endif
1710: if (x->petscnative) {
1711: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1712: if (x->offloadmask == PETSC_OFFLOAD_GPU) {
1713: #if defined(PETSC_HAVE_VIENNACL)
1714: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1715: if (is_viennacltype) {
1716: VecViennaCLCopyFromGPU(x);
1717: } else
1718: #endif
1719: {
1720: #if defined(PETSC_HAVE_CUDA)
1721: VecCUDACopyFromGPU(x);
1722: #endif
1723: }
1724: }
1725: #endif
1726: *a = *((PetscScalar **)x->data);
1727: } else if (x->ops->getarrayread) {
1728: (*x->ops->getarrayread)(x,a);
1729: } else {
1730: (*x->ops->getarray)(x,(PetscScalar**)a);
1731: }
1732: return(0);
1733: }
1735: /*@C
1736: VecGetArrayReadInPlace - Like VecGetArrayRead(), but if this is a CUDA vector and it is currently offloaded to GPU,
1737: the returned pointer will be a GPU pointer to the GPU memory that contains this processor's portion of the
1738: vector data. Otherwise, it functions as VecGetArrayRead().
1740: Not Collective
1742: Input Parameters:
1743: . x - the vector
1745: Output Parameter:
1746: . a - the array
1748: Level: beginner
1750: Notes:
1751: The array must be returned using a matching call to VecRestoreArrayReadInPlace().
1754: .seealso: VecRestoreArrayReadInPlace(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayInPlace()
1755: @*/
1756: PetscErrorCode VecGetArrayReadInPlace(Vec x,const PetscScalar **a)
1757: {
1762: #if defined(PETSC_HAVE_CUDA)
1763: if (x->petscnative && x->offloadmask & PETSC_OFFLOAD_GPU) {
1764: PetscBool is_cudatype = PETSC_FALSE;
1765: PetscObjectTypeCompareAny((PetscObject)x,&is_cudatype,VECSEQCUDA,VECMPICUDA,VECCUDA,"");
1766: if (is_cudatype) {
1767: VecCUDAGetArrayRead(x,a);
1768: return(0);
1769: }
1770: }
1771: #endif
1772: VecGetArrayRead(x,a);
1773: return(0);
1774: }
1776: /*@C
1777: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1778: that were created by a call to VecDuplicateVecs(). You MUST call
1779: VecRestoreArrays() when you no longer need access to the array.
1781: Logically Collective on Vec
1783: Input Parameter:
1784: + x - the vectors
1785: - n - the number of vectors
1787: Output Parameter:
1788: . a - location to put pointer to the array
1790: Fortran Note:
1791: This routine is not supported in Fortran.
1793: Level: intermediate
1795: .seealso: VecGetArray(), VecRestoreArrays()
1796: @*/
1797: PetscErrorCode VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1798: {
1800: PetscInt i;
1801: PetscScalar **q;
1807: if (n <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must get at least one array n = %D",n);
1808: PetscMalloc1(n,&q);
1809: for (i=0; i<n; ++i) {
1810: VecGetArray(x[i],&q[i]);
1811: }
1812: *a = q;
1813: return(0);
1814: }
1816: /*@C
1817: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
1818: has been called.
1820: Logically Collective on Vec
1822: Input Parameters:
1823: + x - the vector
1824: . n - the number of vectors
1825: - a - location of pointer to arrays obtained from VecGetArrays()
1827: Notes:
1828: For regular PETSc vectors this routine does not involve any copies. For
1829: any special vectors that do not store local vector data in a contiguous
1830: array, this routine will copy the data back into the underlying
1831: vector data structure from the arrays obtained with VecGetArrays().
1833: Fortran Note:
1834: This routine is not supported in Fortran.
1836: Level: intermediate
1838: .seealso: VecGetArrays(), VecRestoreArray()
1839: @*/
1840: PetscErrorCode VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1841: {
1843: PetscInt i;
1844: PetscScalar **q = *a;
1851: for (i=0; i<n; ++i) {
1852: VecRestoreArray(x[i],&q[i]);
1853: }
1854: PetscFree(q);
1855: return(0);
1856: }
1858: /*@C
1859: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1861: Logically Collective on Vec
1863: Input Parameters:
1864: + x - the vector
1865: - a - location of pointer to array obtained from VecGetArray()
1867: Level: beginner
1869: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1870: VecGetArrayPair(), VecRestoreArrayPair()
1871: @*/
1872: PetscErrorCode VecRestoreArray(Vec x,PetscScalar **a)
1873: {
1878: if (x->petscnative) {
1879: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1880: x->offloadmask = PETSC_OFFLOAD_CPU;
1881: #endif
1882: } else {
1883: (*x->ops->restorearray)(x,a);
1884: }
1885: if (a) *a = NULL;
1886: PetscObjectStateIncrease((PetscObject)x);
1887: return(0);
1888: }
1890: /*@C
1891: VecRestoreArrayInPlace - Restores a vector after VecGetArrayInPlace() has been called.
1893: Logically Collective on Vec
1895: Input Parameters:
1896: + x - the vector
1897: - a - location of pointer to array obtained from VecGetArrayInPlace()
1899: Level: beginner
1901: .seealso: VecGetArrayInPlace(), VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(),
1902: VecPlaceArray(), VecRestoreArray2d(), VecGetArrayPair(), VecRestoreArrayPair()
1903: @*/
1904: PetscErrorCode VecRestoreArrayInPlace(Vec x,PetscScalar **a)
1905: {
1910: #if defined(PETSC_HAVE_CUDA)
1911: if (x->petscnative && x->offloadmask == PETSC_OFFLOAD_GPU) {
1912: PetscBool is_cudatype = PETSC_FALSE;
1913: PetscObjectTypeCompareAny((PetscObject)x,&is_cudatype,VECSEQCUDA,VECMPICUDA,VECCUDA,"");
1914: if (is_cudatype) {
1915: VecCUDARestoreArray(x,a);
1916: return(0);
1917: }
1918: }
1919: #endif
1920: VecRestoreArray(x,a);
1921: return(0);
1922: }
1925: /*@C
1926: VecRestoreArrayWrite - Restores a vector after VecGetArrayWrite() has been called.
1928: Logically Collective on Vec
1930: Input Parameters:
1931: + x - the vector
1932: - a - location of pointer to array obtained from VecGetArray()
1934: Level: beginner
1936: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1937: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite()
1938: @*/
1939: PetscErrorCode VecRestoreArrayWrite(Vec x,PetscScalar **a)
1940: {
1945: if (x->petscnative) {
1946: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1947: x->offloadmask = PETSC_OFFLOAD_CPU;
1948: #endif
1949: } else {
1950: if (x->ops->restorearraywrite) {
1951: (*x->ops->restorearraywrite)(x,a);
1952: } else {
1953: (*x->ops->restorearray)(x,a);
1954: }
1955: }
1956: if (a) *a = NULL;
1957: PetscObjectStateIncrease((PetscObject)x);
1958: return(0);
1959: }
1961: /*@C
1962: VecRestoreArrayRead - Restore array obtained with VecGetArrayRead()
1964: Not Collective
1966: Input Parameters:
1967: + vec - the vector
1968: - array - the array
1970: Level: beginner
1972: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1973: @*/
1974: PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar **a)
1975: {
1980: if (x->petscnative) {
1981: /* nothing */
1982: } else if (x->ops->restorearrayread) {
1983: (*x->ops->restorearrayread)(x,a);
1984: } else {
1985: (*x->ops->restorearray)(x,(PetscScalar**)a);
1986: }
1987: if (a) *a = NULL;
1988: return(0);
1989: }
1991: /*@C
1992: VecRestoreArrayReadInPlace - Restore array obtained with VecGetArrayReadInPlace()
1994: Not Collective
1996: Input Parameters:
1997: + vec - the vector
1998: - array - the array
2000: Level: beginner
2002: .seealso: VecGetArrayReadInPlace(), VecRestoreArrayInPlace(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
2003: @*/
2004: PetscErrorCode VecRestoreArrayReadInPlace(Vec x,const PetscScalar **a)
2005: {
2009: VecRestoreArrayRead(x,a);
2010: return(0);
2011: }
2013: /*@
2014: VecPlaceArray - Allows one to replace the array in a vector with an
2015: array provided by the user. This is useful to avoid copying an array
2016: into a vector.
2018: Not Collective
2020: Input Parameters:
2021: + vec - the vector
2022: - array - the array
2024: Notes:
2025: You can return to the original array with a call to VecResetArray()
2027: Level: developer
2029: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()
2031: @*/
2032: PetscErrorCode VecPlaceArray(Vec vec,const PetscScalar array[])
2033: {
2040: if (vec->ops->placearray) {
2041: (*vec->ops->placearray)(vec,array);
2042: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot place array in this type of vector");
2043: PetscObjectStateIncrease((PetscObject)vec);
2044: return(0);
2045: }
2047: /*@C
2048: VecReplaceArray - Allows one to replace the array in a vector with an
2049: array provided by the user. This is useful to avoid copying an array
2050: into a vector.
2052: Not Collective
2054: Input Parameters:
2055: + vec - the vector
2056: - array - the array
2058: Notes:
2059: This permanently replaces the array and frees the memory associated
2060: with the old array.
2062: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
2063: freed by the user. It will be freed when the vector is destroyed.
2065: Not supported from Fortran
2067: Level: developer
2069: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()
2071: @*/
2072: PetscErrorCode VecReplaceArray(Vec vec,const PetscScalar array[])
2073: {
2079: if (vec->ops->replacearray) {
2080: (*vec->ops->replacearray)(vec,array);
2081: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot replace array in this type of vector");
2082: PetscObjectStateIncrease((PetscObject)vec);
2083: return(0);
2084: }
2087: /*@C
2088: VecCUDAGetArray - Provides access to the CUDA buffer inside a vector.
2090: This function has semantics similar to VecGetArray(): the pointer
2091: returned by this function points to a consistent view of the vector
2092: data. This may involve a copy operation of data from the host to the
2093: device if the data on the device is out of date. If the device
2094: memory hasn't been allocated previously it will be allocated as part
2095: of this function call. VecCUDAGetArray() assumes that
2096: the user will modify the vector data. This is similar to
2097: intent(inout) in fortran.
2099: The CUDA device pointer has to be released by calling
2100: VecCUDARestoreArray(). Upon restoring the vector data
2101: the data on the host will be marked as out of date. A subsequent
2102: access of the host data will thus incur a data transfer from the
2103: device to the host.
2106: Input Parameter:
2107: . v - the vector
2109: Output Parameter:
2110: . a - the CUDA device pointer
2112: Fortran note:
2113: This function is not currently available from Fortran.
2115: Level: intermediate
2117: .seealso: VecCUDARestoreArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2118: @*/
2119: PETSC_EXTERN PetscErrorCode VecCUDAGetArray(Vec v, PetscScalar **a)
2120: {
2121: #if defined(PETSC_HAVE_CUDA)
2123: #endif
2127: #if defined(PETSC_HAVE_CUDA)
2128: *a = 0;
2129: VecCUDACopyToGPU(v);
2130: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2131: #endif
2132: return(0);
2133: }
2135: /*@C
2136: VecCUDARestoreArray - Restore a CUDA device pointer previously acquired with VecCUDAGetArray().
2138: This marks the host data as out of date. Subsequent access to the
2139: vector data on the host side with for instance VecGetArray() incurs a
2140: data transfer.
2142: Input Parameter:
2143: + v - the vector
2144: - a - the CUDA device pointer. This pointer is invalid after
2145: VecCUDARestoreArray() returns.
2147: Fortran note:
2148: This function is not currently available from Fortran.
2150: Level: intermediate
2152: .seealso: VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2153: @*/
2154: PETSC_EXTERN PetscErrorCode VecCUDARestoreArray(Vec v, PetscScalar **a)
2155: {
2160: #if defined(PETSC_HAVE_CUDA)
2161: v->offloadmask = PETSC_OFFLOAD_GPU;
2162: #endif
2164: PetscObjectStateIncrease((PetscObject)v);
2165: return(0);
2166: }
2168: /*@C
2169: VecCUDAGetArrayRead - Provides read access to the CUDA buffer inside a vector.
2171: This function is analogous to VecGetArrayRead(): The pointer
2172: returned by this function points to a consistent view of the vector
2173: data. This may involve a copy operation of data from the host to the
2174: device if the data on the device is out of date. If the device
2175: memory hasn't been allocated previously it will be allocated as part
2176: of this function call. VecCUDAGetArrayRead() assumes that the
2177: user will not modify the vector data. This is analgogous to
2178: intent(in) in Fortran.
2180: The CUDA device pointer has to be released by calling
2181: VecCUDARestoreArrayRead(). If the data on the host side was
2182: previously up to date it will remain so, i.e. data on both the device
2183: and the host is up to date. Accessing data on the host side does not
2184: incur a device to host data transfer.
2186: Input Parameter:
2187: . v - the vector
2189: Output Parameter:
2190: . a - the CUDA pointer.
2192: Fortran note:
2193: This function is not currently available from Fortran.
2195: Level: intermediate
2197: .seealso: VecCUDARestoreArrayRead(), VecCUDAGetArray(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2198: @*/
2199: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayRead(Vec v, const PetscScalar **a)
2200: {
2201: #if defined(PETSC_HAVE_CUDA)
2203: #endif
2207: #if defined(PETSC_HAVE_CUDA)
2208: *a = 0;
2209: VecCUDACopyToGPU(v);
2210: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2211: #endif
2212: return(0);
2213: }
2215: /*@C
2216: VecCUDARestoreArrayRead - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayRead().
2218: If the data on the host side was previously up to date it will remain
2219: so, i.e. data on both the device and the host is up to date.
2220: Accessing data on the host side e.g. with VecGetArray() does not
2221: incur a device to host data transfer.
2223: Input Parameter:
2224: + v - the vector
2225: - a - the CUDA device pointer. This pointer is invalid after
2226: VecCUDARestoreArrayRead() returns.
2228: Fortran note:
2229: This function is not currently available from Fortran.
2231: Level: intermediate
2233: .seealso: VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2234: @*/
2235: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayRead(Vec v, const PetscScalar **a)
2236: {
2239: *a = NULL;
2240: return(0);
2241: }
2243: /*@C
2244: VecCUDAGetArrayWrite - Provides write access to the CUDA buffer inside a vector.
2246: The data pointed to by the device pointer is uninitialized. The user
2247: may not read from this data. Furthermore, the entire array needs to
2248: be filled by the user to obtain well-defined behaviour. The device
2249: memory will be allocated by this function if it hasn't been allocated
2250: previously. This is analogous to intent(out) in Fortran.
2252: The device pointer needs to be released with
2253: VecCUDARestoreArrayWrite(). When the pointer is released the
2254: host data of the vector is marked as out of data. Subsequent access
2255: of the host data with e.g. VecGetArray() incurs a device to host data
2256: transfer.
2259: Input Parameter:
2260: . v - the vector
2262: Output Parameter:
2263: . a - the CUDA pointer
2265: Fortran note:
2266: This function is not currently available from Fortran.
2268: Level: advanced
2270: .seealso: VecCUDARestoreArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2271: @*/
2272: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayWrite(Vec v, PetscScalar **a)
2273: {
2274: #if defined(PETSC_HAVE_CUDA)
2276: #endif
2280: #if defined(PETSC_HAVE_CUDA)
2281: *a = 0;
2282: VecCUDAAllocateCheck(v);
2283: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2284: #endif
2285: return(0);
2286: }
2288: /*@C
2289: VecCUDARestoreArrayWrite - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayWrite().
2291: Data on the host will be marked as out of date. Subsequent access of
2292: the data on the host side e.g. with VecGetArray() will incur a device
2293: to host data transfer.
2295: Input Parameter:
2296: + v - the vector
2297: - a - the CUDA device pointer. This pointer is invalid after
2298: VecCUDARestoreArrayWrite() returns.
2300: Fortran note:
2301: This function is not currently available from Fortran.
2303: Level: intermediate
2305: .seealso: VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2306: @*/
2307: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayWrite(Vec v, PetscScalar **a)
2308: {
2313: #if defined(PETSC_HAVE_CUDA)
2314: v->offloadmask = PETSC_OFFLOAD_GPU;
2315: #endif
2317: PetscObjectStateIncrease((PetscObject)v);
2318: return(0);
2319: }
2321: /*@C
2322: VecCUDAPlaceArray - Allows one to replace the GPU array in a vector with a
2323: GPU array provided by the user. This is useful to avoid copying an
2324: array into a vector.
2326: Not Collective
2328: Input Parameters:
2329: + vec - the vector
2330: - array - the GPU array
2332: Notes:
2333: You can return to the original GPU array with a call to VecCUDAResetArray()
2334: It is not possible to use VecCUDAPlaceArray() and VecPlaceArray() at the
2335: same time on the same vector.
2337: Level: developer
2339: .seealso: VecPlaceArray(), VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAReplaceArray()
2341: @*/
2342: PetscErrorCode VecCUDAPlaceArray(Vec vin,const PetscScalar a[])
2343: {
2348: #if defined(PETSC_HAVE_CUDA)
2349: VecCUDACopyToGPU(vin);
2350: 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()");
2351: ((Vec_Seq*)vin->data)->unplacedarray = (PetscScalar *) ((Vec_CUDA*)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2352: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar*)a;
2353: vin->offloadmask = PETSC_OFFLOAD_GPU;
2354: #endif
2355: PetscObjectStateIncrease((PetscObject)vin);
2356: return(0);
2357: }
2359: /*@C
2360: VecCUDAReplaceArray - Allows one to replace the GPU array in a vector
2361: with a GPU array provided by the user. This is useful to avoid copying
2362: a GPU array into a vector.
2364: Not Collective
2366: Input Parameters:
2367: + vec - the vector
2368: - array - the GPU array
2370: Notes:
2371: This permanently replaces the GPU array and frees the memory associated
2372: with the old GPU array.
2374: The memory passed in CANNOT be freed by the user. It will be freed
2375: when the vector is destroyed.
2377: Not supported from Fortran
2379: Level: developer
2381: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAPlaceArray(), VecReplaceArray()
2383: @*/
2384: PetscErrorCode VecCUDAReplaceArray(Vec vin,const PetscScalar a[])
2385: {
2386: #if defined(PETSC_HAVE_CUDA)
2387: cudaError_t err;
2388: #endif
2393: #if defined(PETSC_HAVE_CUDA)
2394: err = cudaFree(((Vec_CUDA*)vin->spptr)->GPUarray);CHKERRCUDA(err);
2395: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar*)a;
2396: vin->offloadmask = PETSC_OFFLOAD_GPU;
2397: #endif
2398: PetscObjectStateIncrease((PetscObject)vin);
2399: return(0);
2400: }
2402: /*@C
2403: VecCUDAResetArray - Resets a vector to use its default memory. Call this
2404: after the use of VecCUDAPlaceArray().
2406: Not Collective
2408: Input Parameters:
2409: . vec - the vector
2411: Level: developer
2413: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray(), VecResetArray(), VecCUDAPlaceArray(), VecCUDAReplaceArray()
2415: @*/
2416: PetscErrorCode VecCUDAResetArray(Vec vin)
2417: {
2422: #if defined(PETSC_HAVE_CUDA)
2423: VecCUDACopyToGPU(vin);
2424: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar *) ((Vec_Seq*)vin->data)->unplacedarray;
2425: ((Vec_Seq*)vin->data)->unplacedarray = 0;
2426: vin->offloadmask = PETSC_OFFLOAD_GPU;
2427: #endif
2428: PetscObjectStateIncrease((PetscObject)vin);
2429: return(0);
2430: }
2432: /*MC
2433: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
2434: and makes them accessible via a Fortran90 pointer.
2436: Synopsis:
2437: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
2439: Collective on Vec
2441: Input Parameters:
2442: + x - a vector to mimic
2443: - n - the number of vectors to obtain
2445: Output Parameters:
2446: + y - Fortran90 pointer to the array of vectors
2447: - ierr - error code
2449: Example of Usage:
2450: .vb
2451: #include <petsc/finclude/petscvec.h>
2452: use petscvec
2454: Vec x
2455: Vec, pointer :: y(:)
2456: ....
2457: call VecDuplicateVecsF90(x,2,y,ierr)
2458: call VecSet(y(2),alpha,ierr)
2459: call VecSet(y(2),alpha,ierr)
2460: ....
2461: call VecDestroyVecsF90(2,y,ierr)
2462: .ve
2464: Notes:
2465: Not yet supported for all F90 compilers
2467: Use VecDestroyVecsF90() to free the space.
2469: Level: beginner
2471: .seealso: VecDestroyVecsF90(), VecDuplicateVecs()
2473: M*/
2475: /*MC
2476: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
2477: VecGetArrayF90().
2479: Synopsis:
2480: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2482: Logically Collective on Vec
2484: Input Parameters:
2485: + x - vector
2486: - xx_v - the Fortran90 pointer to the array
2488: Output Parameter:
2489: . ierr - error code
2491: Example of Usage:
2492: .vb
2493: #include <petsc/finclude/petscvec.h>
2494: use petscvec
2496: PetscScalar, pointer :: xx_v(:)
2497: ....
2498: call VecGetArrayF90(x,xx_v,ierr)
2499: xx_v(3) = a
2500: call VecRestoreArrayF90(x,xx_v,ierr)
2501: .ve
2503: Level: beginner
2505: .seealso: VecGetArrayF90(), VecGetArray(), VecRestoreArray(), UsingFortran, VecRestoreArrayReadF90()
2507: M*/
2509: /*MC
2510: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
2512: Synopsis:
2513: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
2515: Collective on Vec
2517: Input Parameters:
2518: + n - the number of vectors previously obtained
2519: - x - pointer to array of vector pointers
2521: Output Parameter:
2522: . ierr - error code
2524: Notes:
2525: Not yet supported for all F90 compilers
2527: Level: beginner
2529: .seealso: VecDestroyVecs(), VecDuplicateVecsF90()
2531: M*/
2533: /*MC
2534: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
2535: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2536: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
2537: when you no longer need access to the array.
2539: Synopsis:
2540: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2542: Logically Collective on Vec
2544: Input Parameter:
2545: . x - vector
2547: Output Parameters:
2548: + xx_v - the Fortran90 pointer to the array
2549: - ierr - error code
2551: Example of Usage:
2552: .vb
2553: #include <petsc/finclude/petscvec.h>
2554: use petscvec
2556: PetscScalar, pointer :: xx_v(:)
2557: ....
2558: call VecGetArrayF90(x,xx_v,ierr)
2559: xx_v(3) = a
2560: call VecRestoreArrayF90(x,xx_v,ierr)
2561: .ve
2563: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
2565: Level: beginner
2567: .seealso: VecRestoreArrayF90(), VecGetArray(), VecRestoreArray(), VecGetArrayReadF90(), UsingFortran
2569: M*/
2571: /*MC
2572: VecGetArrayReadF90 - Accesses a read only array from Fortran90. For default PETSc
2573: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2574: this routine is implementation dependent. You MUST call VecRestoreArrayReadF90()
2575: when you no longer need access to the array.
2577: Synopsis:
2578: VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2580: Logically Collective on Vec
2582: Input Parameter:
2583: . x - vector
2585: Output Parameters:
2586: + xx_v - the Fortran90 pointer to the array
2587: - ierr - error code
2589: Example of Usage:
2590: .vb
2591: #include <petsc/finclude/petscvec.h>
2592: use petscvec
2594: PetscScalar, pointer :: xx_v(:)
2595: ....
2596: call VecGetArrayReadF90(x,xx_v,ierr)
2597: a = xx_v(3)
2598: call VecRestoreArrayReadF90(x,xx_v,ierr)
2599: .ve
2601: If you intend to write entries into the array you must use VecGetArrayF90().
2603: Level: beginner
2605: .seealso: VecRestoreArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecGetArrayF90(), UsingFortran
2607: M*/
2609: /*MC
2610: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
2611: VecGetArrayReadF90().
2613: Synopsis:
2614: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2616: Logically Collective on Vec
2618: Input Parameters:
2619: + x - vector
2620: - xx_v - the Fortran90 pointer to the array
2622: Output Parameter:
2623: . ierr - error code
2625: Example of Usage:
2626: .vb
2627: #include <petsc/finclude/petscvec.h>
2628: use petscvec
2630: PetscScalar, pointer :: xx_v(:)
2631: ....
2632: call VecGetArrayReadF90(x,xx_v,ierr)
2633: a = xx_v(3)
2634: call VecRestoreArrayReadF90(x,xx_v,ierr)
2635: .ve
2637: Level: beginner
2639: .seealso: VecGetArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(),UsingFortran, VecRestoreArrayF90()
2641: M*/
2643: /*@C
2644: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
2645: processor's portion of the vector data. You MUST call VecRestoreArray2d()
2646: when you no longer need access to the array.
2648: Logically Collective
2650: Input Parameter:
2651: + x - the vector
2652: . m - first dimension of two dimensional array
2653: . n - second dimension of two dimensional array
2654: . mstart - first index you will use in first coordinate direction (often 0)
2655: - nstart - first index in the second coordinate direction (often 0)
2657: Output Parameter:
2658: . a - location to put pointer to the array
2660: Level: developer
2662: Notes:
2663: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2664: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2665: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2666: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2668: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2670: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2671: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2672: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2673: @*/
2674: PetscErrorCode VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2675: {
2677: PetscInt i,N;
2678: PetscScalar *aa;
2684: VecGetLocalSize(x,&N);
2685: 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);
2686: VecGetArray(x,&aa);
2688: PetscMalloc1(m,a);
2689: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
2690: *a -= mstart;
2691: return(0);
2692: }
2694: /*@C
2695: VecGetArray2dWrite - Returns a pointer to a 2d contiguous array that will contain this
2696: processor's portion of the vector data. You MUST call VecRestoreArray2dWrite()
2697: when you no longer need access to the array.
2699: Logically Collective
2701: Input Parameter:
2702: + x - the vector
2703: . m - first dimension of two dimensional array
2704: . n - second dimension of two dimensional array
2705: . mstart - first index you will use in first coordinate direction (often 0)
2706: - nstart - first index in the second coordinate direction (often 0)
2708: Output Parameter:
2709: . a - location to put pointer to the array
2711: Level: developer
2713: Notes:
2714: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2715: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2716: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2717: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2719: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2721: Concepts: vector^accessing local values as 2d array
2723: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2724: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2725: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2726: @*/
2727: PetscErrorCode VecGetArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2728: {
2730: PetscInt i,N;
2731: PetscScalar *aa;
2737: VecGetLocalSize(x,&N);
2738: 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);
2739: VecGetArrayWrite(x,&aa);
2741: PetscMalloc1(m,a);
2742: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
2743: *a -= mstart;
2744: return(0);
2745: }
2747: /*@C
2748: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
2750: Logically Collective
2752: Input Parameters:
2753: + x - the vector
2754: . m - first dimension of two dimensional array
2755: . n - second dimension of the two dimensional array
2756: . mstart - first index you will use in first coordinate direction (often 0)
2757: . nstart - first index in the second coordinate direction (often 0)
2758: - a - location of pointer to array obtained from VecGetArray2d()
2760: Level: developer
2762: Notes:
2763: For regular PETSc vectors this routine does not involve any copies. For
2764: any special vectors that do not store local vector data in a contiguous
2765: array, this routine will copy the data back into the underlying
2766: vector data structure from the array obtained with VecGetArray().
2768: This routine actually zeros out the a pointer.
2770: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2771: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2772: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2773: @*/
2774: PetscErrorCode VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2775: {
2777: void *dummy;
2783: dummy = (void*)(*a + mstart);
2784: PetscFree(dummy);
2785: VecRestoreArray(x,NULL);
2786: return(0);
2787: }
2789: /*@C
2790: VecRestoreArray2dWrite - Restores a vector after VecGetArray2dWrite() has been called.
2792: Logically Collective
2794: Input Parameters:
2795: + x - the vector
2796: . m - first dimension of two dimensional array
2797: . n - second dimension of the two dimensional array
2798: . mstart - first index you will use in first coordinate direction (often 0)
2799: . nstart - first index in the second coordinate direction (often 0)
2800: - a - location of pointer to array obtained from VecGetArray2d()
2802: Level: developer
2804: Notes:
2805: For regular PETSc vectors this routine does not involve any copies. For
2806: any special vectors that do not store local vector data in a contiguous
2807: array, this routine will copy the data back into the underlying
2808: vector data structure from the array obtained with VecGetArray().
2810: This routine actually zeros out the a pointer.
2812: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2813: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2814: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2815: @*/
2816: PetscErrorCode VecRestoreArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2817: {
2819: void *dummy;
2825: dummy = (void*)(*a + mstart);
2826: PetscFree(dummy);
2827: VecRestoreArrayWrite(x,NULL);
2828: return(0);
2829: }
2831: /*@C
2832: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
2833: processor's portion of the vector data. You MUST call VecRestoreArray1d()
2834: when you no longer need access to the array.
2836: Logically Collective
2838: Input Parameter:
2839: + x - the vector
2840: . m - first dimension of two dimensional array
2841: - mstart - first index you will use in first coordinate direction (often 0)
2843: Output Parameter:
2844: . a - location to put pointer to the array
2846: Level: developer
2848: Notes:
2849: For a vector obtained from DMCreateLocalVector() mstart are likely
2850: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2851: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2853: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2855: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2856: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2857: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2858: @*/
2859: PetscErrorCode VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2860: {
2862: PetscInt N;
2868: VecGetLocalSize(x,&N);
2869: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2870: VecGetArray(x,a);
2871: *a -= mstart;
2872: return(0);
2873: }
2875: /*@C
2876: VecGetArray1dWrite - Returns a pointer to a 1d contiguous array that will contain this
2877: processor's portion of the vector data. You MUST call VecRestoreArray1dWrite()
2878: when you no longer need access to the array.
2880: Logically Collective
2882: Input Parameter:
2883: + x - the vector
2884: . m - first dimension of two dimensional array
2885: - mstart - first index you will use in first coordinate direction (often 0)
2887: Output Parameter:
2888: . a - location to put pointer to the array
2890: Level: developer
2892: Notes:
2893: For a vector obtained from DMCreateLocalVector() mstart are likely
2894: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2895: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2897: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2899: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2900: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2901: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2902: @*/
2903: PetscErrorCode VecGetArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2904: {
2906: PetscInt N;
2912: VecGetLocalSize(x,&N);
2913: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2914: VecGetArrayWrite(x,a);
2915: *a -= mstart;
2916: return(0);
2917: }
2919: /*@C
2920: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
2922: Logically Collective
2924: Input Parameters:
2925: + x - the vector
2926: . m - first dimension of two dimensional array
2927: . mstart - first index you will use in first coordinate direction (often 0)
2928: - a - location of pointer to array obtained from VecGetArray21()
2930: Level: developer
2932: Notes:
2933: For regular PETSc vectors this routine does not involve any copies. For
2934: any special vectors that do not store local vector data in a contiguous
2935: array, this routine will copy the data back into the underlying
2936: vector data structure from the array obtained with VecGetArray1d().
2938: This routine actually zeros out the a pointer.
2940: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2941: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2942: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2943: @*/
2944: PetscErrorCode VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2945: {
2951: VecRestoreArray(x,NULL);
2952: return(0);
2953: }
2955: /*@C
2956: VecRestoreArray1dWrite - Restores a vector after VecGetArray1dWrite() has been called.
2958: Logically Collective
2960: Input Parameters:
2961: + x - the vector
2962: . m - first dimension of two dimensional array
2963: . mstart - first index you will use in first coordinate direction (often 0)
2964: - a - location of pointer to array obtained from VecGetArray21()
2966: Level: developer
2968: Notes:
2969: For regular PETSc vectors this routine does not involve any copies. For
2970: any special vectors that do not store local vector data in a contiguous
2971: array, this routine will copy the data back into the underlying
2972: vector data structure from the array obtained with VecGetArray1d().
2974: This routine actually zeros out the a pointer.
2976: Concepts: vector^accessing local values as 1d array
2978: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2979: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2980: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2981: @*/
2982: PetscErrorCode VecRestoreArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2983: {
2989: VecRestoreArrayWrite(x,NULL);
2990: return(0);
2991: }
2993: /*@C
2994: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
2995: processor's portion of the vector data. You MUST call VecRestoreArray3d()
2996: when you no longer need access to the array.
2998: Logically Collective
3000: Input Parameter:
3001: + x - the vector
3002: . m - first dimension of three dimensional array
3003: . n - second dimension of three dimensional array
3004: . p - third dimension of three dimensional array
3005: . mstart - first index you will use in first coordinate direction (often 0)
3006: . nstart - first index in the second coordinate direction (often 0)
3007: - pstart - first index in the third coordinate direction (often 0)
3009: Output Parameter:
3010: . a - location to put pointer to the array
3012: Level: developer
3014: Notes:
3015: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3016: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3017: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3018: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3020: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3022: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3023: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3024: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3025: @*/
3026: PetscErrorCode VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3027: {
3029: PetscInt i,N,j;
3030: PetscScalar *aa,**b;
3036: VecGetLocalSize(x,&N);
3037: 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);
3038: VecGetArray(x,&aa);
3040: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3041: b = (PetscScalar**)((*a) + m);
3042: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3043: for (i=0; i<m; i++)
3044: for (j=0; j<n; j++)
3045: b[i*n+j] = aa + i*n*p + j*p - pstart;
3047: *a -= mstart;
3048: return(0);
3049: }
3051: /*@C
3052: VecGetArray3dWrite - Returns a pointer to a 3d contiguous array that will contain this
3053: processor's portion of the vector data. You MUST call VecRestoreArray3dWrite()
3054: when you no longer need access to the array.
3056: Logically Collective
3058: Input Parameter:
3059: + x - the vector
3060: . m - first dimension of three dimensional array
3061: . n - second dimension of three dimensional array
3062: . p - third dimension of three dimensional array
3063: . mstart - first index you will use in first coordinate direction (often 0)
3064: . nstart - first index in the second coordinate direction (often 0)
3065: - pstart - first index in the third coordinate direction (often 0)
3067: Output Parameter:
3068: . a - location to put pointer to the array
3070: Level: developer
3072: Notes:
3073: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3074: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3075: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3076: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3078: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3080: Concepts: vector^accessing local values as 3d array
3082: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3083: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3084: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3085: @*/
3086: PetscErrorCode VecGetArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3087: {
3089: PetscInt i,N,j;
3090: PetscScalar *aa,**b;
3096: VecGetLocalSize(x,&N);
3097: 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);
3098: VecGetArrayWrite(x,&aa);
3100: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3101: b = (PetscScalar**)((*a) + m);
3102: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3103: for (i=0; i<m; i++)
3104: for (j=0; j<n; j++)
3105: b[i*n+j] = aa + i*n*p + j*p - pstart;
3107: *a -= mstart;
3108: return(0);
3109: }
3111: /*@C
3112: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
3114: Logically Collective
3116: Input Parameters:
3117: + x - the vector
3118: . m - first dimension of three dimensional array
3119: . n - second dimension of the three dimensional array
3120: . p - third dimension of the three dimensional array
3121: . mstart - first index you will use in first coordinate direction (often 0)
3122: . nstart - first index in the second coordinate direction (often 0)
3123: . pstart - first index in the third coordinate direction (often 0)
3124: - a - location of pointer to array obtained from VecGetArray3d()
3126: Level: developer
3128: Notes:
3129: For regular PETSc vectors this routine does not involve any copies. For
3130: any special vectors that do not store local vector data in a contiguous
3131: array, this routine will copy the data back into the underlying
3132: vector data structure from the array obtained with VecGetArray().
3134: This routine actually zeros out the a pointer.
3136: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3137: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3138: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3139: @*/
3140: PetscErrorCode VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3141: {
3143: void *dummy;
3149: dummy = (void*)(*a + mstart);
3150: PetscFree(dummy);
3151: VecRestoreArray(x,NULL);
3152: return(0);
3153: }
3155: /*@C
3156: VecRestoreArray3dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3158: Logically Collective
3160: Input Parameters:
3161: + x - the vector
3162: . m - first dimension of three dimensional array
3163: . n - second dimension of the three dimensional array
3164: . p - third dimension of the three dimensional array
3165: . mstart - first index you will use in first coordinate direction (often 0)
3166: . nstart - first index in the second coordinate direction (often 0)
3167: . pstart - first index in the third coordinate direction (often 0)
3168: - a - location of pointer to array obtained from VecGetArray3d()
3170: Level: developer
3172: Notes:
3173: For regular PETSc vectors this routine does not involve any copies. For
3174: any special vectors that do not store local vector data in a contiguous
3175: array, this routine will copy the data back into the underlying
3176: vector data structure from the array obtained with VecGetArray().
3178: This routine actually zeros out the a pointer.
3180: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3181: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3182: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3183: @*/
3184: PetscErrorCode VecRestoreArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3185: {
3187: void *dummy;
3193: dummy = (void*)(*a + mstart);
3194: PetscFree(dummy);
3195: VecRestoreArrayWrite(x,NULL);
3196: return(0);
3197: }
3199: /*@C
3200: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
3201: processor's portion of the vector data. You MUST call VecRestoreArray4d()
3202: when you no longer need access to the array.
3204: Logically Collective
3206: Input Parameter:
3207: + x - the vector
3208: . m - first dimension of four dimensional array
3209: . n - second dimension of four dimensional array
3210: . p - third dimension of four dimensional array
3211: . q - fourth dimension of four dimensional array
3212: . mstart - first index you will use in first coordinate direction (often 0)
3213: . nstart - first index in the second coordinate direction (often 0)
3214: . pstart - first index in the third coordinate direction (often 0)
3215: - qstart - first index in the fourth coordinate direction (often 0)
3217: Output Parameter:
3218: . a - location to put pointer to the array
3220: Level: beginner
3222: Notes:
3223: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3224: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3225: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3226: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3228: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3230: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3231: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3232: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3233: @*/
3234: PetscErrorCode VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3235: {
3237: PetscInt i,N,j,k;
3238: PetscScalar *aa,***b,**c;
3244: VecGetLocalSize(x,&N);
3245: 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);
3246: VecGetArray(x,&aa);
3248: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3249: b = (PetscScalar***)((*a) + m);
3250: c = (PetscScalar**)(b + m*n);
3251: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3252: for (i=0; i<m; i++)
3253: for (j=0; j<n; j++)
3254: b[i*n+j] = c + i*n*p + j*p - pstart;
3255: for (i=0; i<m; i++)
3256: for (j=0; j<n; j++)
3257: for (k=0; k<p; k++)
3258: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3259: *a -= mstart;
3260: return(0);
3261: }
3263: /*@C
3264: VecGetArray4dWrite - Returns a pointer to a 4d contiguous array that will contain this
3265: processor's portion of the vector data. You MUST call VecRestoreArray4dWrite()
3266: when you no longer need access to the array.
3268: Logically Collective
3270: Input Parameter:
3271: + x - the vector
3272: . m - first dimension of four dimensional array
3273: . n - second dimension of four dimensional array
3274: . p - third dimension of four dimensional array
3275: . q - fourth dimension of four dimensional array
3276: . mstart - first index you will use in first coordinate direction (often 0)
3277: . nstart - first index in the second coordinate direction (often 0)
3278: . pstart - first index in the third coordinate direction (often 0)
3279: - qstart - first index in the fourth coordinate direction (often 0)
3281: Output Parameter:
3282: . a - location to put pointer to the array
3284: Level: beginner
3286: Notes:
3287: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3288: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3289: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3290: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3292: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3294: Concepts: vector^accessing local values as 3d array
3296: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3297: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3298: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3299: @*/
3300: PetscErrorCode VecGetArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3301: {
3303: PetscInt i,N,j,k;
3304: PetscScalar *aa,***b,**c;
3310: VecGetLocalSize(x,&N);
3311: 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);
3312: VecGetArrayWrite(x,&aa);
3314: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3315: b = (PetscScalar***)((*a) + m);
3316: c = (PetscScalar**)(b + m*n);
3317: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3318: for (i=0; i<m; i++)
3319: for (j=0; j<n; j++)
3320: b[i*n+j] = c + i*n*p + j*p - pstart;
3321: for (i=0; i<m; i++)
3322: for (j=0; j<n; j++)
3323: for (k=0; k<p; k++)
3324: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3325: *a -= mstart;
3326: return(0);
3327: }
3329: /*@C
3330: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
3332: Logically Collective
3334: Input Parameters:
3335: + x - the vector
3336: . m - first dimension of four dimensional array
3337: . n - second dimension of the four dimensional array
3338: . p - third dimension of the four dimensional array
3339: . q - fourth dimension of the four dimensional array
3340: . mstart - first index you will use in first coordinate direction (often 0)
3341: . nstart - first index in the second coordinate direction (often 0)
3342: . pstart - first index in the third coordinate direction (often 0)
3343: . qstart - first index in the fourth coordinate direction (often 0)
3344: - a - location of pointer to array obtained from VecGetArray4d()
3346: Level: beginner
3348: Notes:
3349: For regular PETSc vectors this routine does not involve any copies. For
3350: any special vectors that do not store local vector data in a contiguous
3351: array, this routine will copy the data back into the underlying
3352: vector data structure from the array obtained with VecGetArray().
3354: This routine actually zeros out the a pointer.
3356: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3357: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3358: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3359: @*/
3360: PetscErrorCode VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3361: {
3363: void *dummy;
3369: dummy = (void*)(*a + mstart);
3370: PetscFree(dummy);
3371: VecRestoreArray(x,NULL);
3372: return(0);
3373: }
3375: /*@C
3376: VecRestoreArray4dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3378: Logically Collective
3380: Input Parameters:
3381: + x - the vector
3382: . m - first dimension of four dimensional array
3383: . n - second dimension of the four dimensional array
3384: . p - third dimension of the four dimensional array
3385: . q - fourth dimension of the four dimensional array
3386: . mstart - first index you will use in first coordinate direction (often 0)
3387: . nstart - first index in the second coordinate direction (often 0)
3388: . pstart - first index in the third coordinate direction (often 0)
3389: . qstart - first index in the fourth coordinate direction (often 0)
3390: - a - location of pointer to array obtained from VecGetArray4d()
3392: Level: beginner
3394: Notes:
3395: For regular PETSc vectors this routine does not involve any copies. For
3396: any special vectors that do not store local vector data in a contiguous
3397: array, this routine will copy the data back into the underlying
3398: vector data structure from the array obtained with VecGetArray().
3400: This routine actually zeros out the a pointer.
3402: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3403: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3404: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3405: @*/
3406: PetscErrorCode VecRestoreArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3407: {
3409: void *dummy;
3415: dummy = (void*)(*a + mstart);
3416: PetscFree(dummy);
3417: VecRestoreArrayWrite(x,NULL);
3418: return(0);
3419: }
3421: /*@C
3422: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
3423: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
3424: when you no longer need access to the array.
3426: Logically Collective
3428: Input Parameter:
3429: + x - the vector
3430: . m - first dimension of two dimensional array
3431: . n - second dimension of two dimensional array
3432: . mstart - first index you will use in first coordinate direction (often 0)
3433: - nstart - first index in the second coordinate direction (often 0)
3435: Output Parameter:
3436: . a - location to put pointer to the array
3438: Level: developer
3440: Notes:
3441: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3442: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3443: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3444: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3446: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3448: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3449: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3450: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3451: @*/
3452: PetscErrorCode VecGetArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3453: {
3454: PetscErrorCode ierr;
3455: PetscInt i,N;
3456: const PetscScalar *aa;
3462: VecGetLocalSize(x,&N);
3463: 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);
3464: VecGetArrayRead(x,&aa);
3466: PetscMalloc1(m,a);
3467: for (i=0; i<m; i++) (*a)[i] = (PetscScalar*) aa + i*n - nstart;
3468: *a -= mstart;
3469: return(0);
3470: }
3472: /*@C
3473: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
3475: Logically Collective
3477: Input Parameters:
3478: + x - the vector
3479: . m - first dimension of two dimensional array
3480: . n - second dimension of the two dimensional array
3481: . mstart - first index you will use in first coordinate direction (often 0)
3482: . nstart - first index in the second coordinate direction (often 0)
3483: - a - location of pointer to array obtained from VecGetArray2d()
3485: Level: developer
3487: Notes:
3488: For regular PETSc vectors this routine does not involve any copies. For
3489: any special vectors that do not store local vector data in a contiguous
3490: array, this routine will copy the data back into the underlying
3491: vector data structure from the array obtained with VecGetArray().
3493: This routine actually zeros out the a pointer.
3495: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3496: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3497: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3498: @*/
3499: PetscErrorCode VecRestoreArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3500: {
3502: void *dummy;
3508: dummy = (void*)(*a + mstart);
3509: PetscFree(dummy);
3510: VecRestoreArrayRead(x,NULL);
3511: return(0);
3512: }
3514: /*@C
3515: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
3516: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
3517: when you no longer need access to the array.
3519: Logically Collective
3521: Input Parameter:
3522: + x - the vector
3523: . m - first dimension of two dimensional array
3524: - mstart - first index you will use in first coordinate direction (often 0)
3526: Output Parameter:
3527: . a - location to put pointer to the array
3529: Level: developer
3531: Notes:
3532: For a vector obtained from DMCreateLocalVector() mstart are likely
3533: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3534: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3536: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3538: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3539: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3540: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3541: @*/
3542: PetscErrorCode VecGetArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3543: {
3545: PetscInt N;
3551: VecGetLocalSize(x,&N);
3552: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
3553: VecGetArrayRead(x,(const PetscScalar**)a);
3554: *a -= mstart;
3555: return(0);
3556: }
3558: /*@C
3559: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
3561: Logically Collective
3563: Input Parameters:
3564: + x - the vector
3565: . m - first dimension of two dimensional array
3566: . mstart - first index you will use in first coordinate direction (often 0)
3567: - a - location of pointer to array obtained from VecGetArray21()
3569: Level: developer
3571: Notes:
3572: For regular PETSc vectors this routine does not involve any copies. For
3573: any special vectors that do not store local vector data in a contiguous
3574: array, this routine will copy the data back into the underlying
3575: vector data structure from the array obtained with VecGetArray1dRead().
3577: This routine actually zeros out the a pointer.
3579: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3580: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3581: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3582: @*/
3583: PetscErrorCode VecRestoreArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3584: {
3590: VecRestoreArrayRead(x,NULL);
3591: return(0);
3592: }
3595: /*@C
3596: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
3597: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
3598: when you no longer need access to the array.
3600: Logically Collective
3602: Input Parameter:
3603: + x - the vector
3604: . m - first dimension of three dimensional array
3605: . n - second dimension of three dimensional array
3606: . p - third dimension of three dimensional array
3607: . mstart - first index you will use in first coordinate direction (often 0)
3608: . nstart - first index in the second coordinate direction (often 0)
3609: - pstart - first index in the third coordinate direction (often 0)
3611: Output Parameter:
3612: . a - location to put pointer to the array
3614: Level: developer
3616: Notes:
3617: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3618: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3619: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3620: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
3622: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3624: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3625: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3626: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3627: @*/
3628: PetscErrorCode VecGetArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3629: {
3630: PetscErrorCode ierr;
3631: PetscInt i,N,j;
3632: const PetscScalar *aa;
3633: PetscScalar **b;
3639: VecGetLocalSize(x,&N);
3640: 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);
3641: VecGetArrayRead(x,&aa);
3643: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3644: b = (PetscScalar**)((*a) + m);
3645: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3646: for (i=0; i<m; i++)
3647: for (j=0; j<n; j++)
3648: b[i*n+j] = (PetscScalar *)aa + i*n*p + j*p - pstart;
3650: *a -= mstart;
3651: return(0);
3652: }
3654: /*@C
3655: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
3657: Logically Collective
3659: Input Parameters:
3660: + x - the vector
3661: . m - first dimension of three dimensional array
3662: . n - second dimension of the three dimensional array
3663: . p - third dimension of the three dimensional array
3664: . mstart - first index you will use in first coordinate direction (often 0)
3665: . nstart - first index in the second coordinate direction (often 0)
3666: . pstart - first index in the third coordinate direction (often 0)
3667: - a - location of pointer to array obtained from VecGetArray3dRead()
3669: Level: developer
3671: Notes:
3672: For regular PETSc vectors this routine does not involve any copies. For
3673: any special vectors that do not store local vector data in a contiguous
3674: array, this routine will copy the data back into the underlying
3675: vector data structure from the array obtained with VecGetArray().
3677: This routine actually zeros out the a pointer.
3679: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3680: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3681: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3682: @*/
3683: PetscErrorCode VecRestoreArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3684: {
3686: void *dummy;
3692: dummy = (void*)(*a + mstart);
3693: PetscFree(dummy);
3694: VecRestoreArrayRead(x,NULL);
3695: return(0);
3696: }
3698: /*@C
3699: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
3700: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
3701: when you no longer need access to the array.
3703: Logically Collective
3705: Input Parameter:
3706: + x - the vector
3707: . m - first dimension of four dimensional array
3708: . n - second dimension of four dimensional array
3709: . p - third dimension of four dimensional array
3710: . q - fourth dimension of four dimensional array
3711: . mstart - first index you will use in first coordinate direction (often 0)
3712: . nstart - first index in the second coordinate direction (often 0)
3713: . pstart - first index in the third coordinate direction (often 0)
3714: - qstart - first index in the fourth coordinate direction (often 0)
3716: Output Parameter:
3717: . a - location to put pointer to the array
3719: Level: beginner
3721: Notes:
3722: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3723: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3724: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3725: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3727: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3729: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3730: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3731: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3732: @*/
3733: PetscErrorCode VecGetArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3734: {
3735: PetscErrorCode ierr;
3736: PetscInt i,N,j,k;
3737: const PetscScalar *aa;
3738: PetscScalar ***b,**c;
3744: VecGetLocalSize(x,&N);
3745: 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);
3746: VecGetArrayRead(x,&aa);
3748: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3749: b = (PetscScalar***)((*a) + m);
3750: c = (PetscScalar**)(b + m*n);
3751: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3752: for (i=0; i<m; i++)
3753: for (j=0; j<n; j++)
3754: b[i*n+j] = c + i*n*p + j*p - pstart;
3755: for (i=0; i<m; i++)
3756: for (j=0; j<n; j++)
3757: for (k=0; k<p; k++)
3758: c[i*n*p+j*p+k] = (PetscScalar*) aa + i*n*p*q + j*p*q + k*q - qstart;
3759: *a -= mstart;
3760: return(0);
3761: }
3763: /*@C
3764: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
3766: Logically Collective
3768: Input Parameters:
3769: + x - the vector
3770: . m - first dimension of four dimensional array
3771: . n - second dimension of the four dimensional array
3772: . p - third dimension of the four dimensional array
3773: . q - fourth dimension of the four dimensional array
3774: . mstart - first index you will use in first coordinate direction (often 0)
3775: . nstart - first index in the second coordinate direction (often 0)
3776: . pstart - first index in the third coordinate direction (often 0)
3777: . qstart - first index in the fourth coordinate direction (often 0)
3778: - a - location of pointer to array obtained from VecGetArray4dRead()
3780: Level: beginner
3782: Notes:
3783: For regular PETSc vectors this routine does not involve any copies. For
3784: any special vectors that do not store local vector data in a contiguous
3785: array, this routine will copy the data back into the underlying
3786: vector data structure from the array obtained with VecGetArray().
3788: This routine actually zeros out the a pointer.
3790: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3791: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3792: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3793: @*/
3794: PetscErrorCode VecRestoreArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3795: {
3797: void *dummy;
3803: dummy = (void*)(*a + mstart);
3804: PetscFree(dummy);
3805: VecRestoreArrayRead(x,NULL);
3806: return(0);
3807: }
3809: #if defined(PETSC_USE_DEBUG)
3811: /*@
3812: VecLockGet - Gets the current lock status of a vector
3814: Logically Collective on Vec
3816: Input Parameter:
3817: . x - the vector
3819: Output Parameter:
3820: . state - greater than zero indicates the vector is locked for read; less then zero indicates the vector is
3821: locked for write; equal to zero means the vector is unlocked, that is, it is free to read or write.
3823: Level: beginner
3825: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop()
3826: @*/
3827: PetscErrorCode VecLockGet(Vec x,PetscInt *state)
3828: {
3831: *state = x->lock;
3832: return(0);
3833: }
3835: /*@
3836: VecLockReadPush - Pushes a read-only lock on a vector to prevent it from writing
3838: Logically Collective on Vec
3840: Input Parameter:
3841: . x - the vector
3843: Notes:
3844: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
3846: The call can be nested, i.e., called multiple times on the same vector, but each VecLockReadPush(x) has to have one matching
3847: VecLockReadPop(x), which removes the latest read-only lock.
3849: Level: beginner
3851: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPop(), VecLockGet()
3852: @*/
3853: PetscErrorCode VecLockReadPush(Vec x)
3854: {
3857: 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");
3858: x->lock++;
3859: return(0);
3860: }
3862: /*@
3863: VecLockReadPop - Pops a read-only lock from a vector
3865: Logically Collective on Vec
3867: Input Parameter:
3868: . x - the vector
3870: Level: beginner
3872: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockGet()
3873: @*/
3874: PetscErrorCode VecLockReadPop(Vec x)
3875: {
3878: x->lock--;
3879: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector has been unlocked from read-only access too many times");
3880: return(0);
3881: }
3883: /*@C
3884: VecLockWriteSet_Private - Lock or unlock a vector for exclusive read/write access
3886: Logically Collective on Vec
3888: Input Parameter:
3889: + x - the vector
3890: - flg - PETSC_TRUE to lock the vector for writing; PETSC_FALSE to unlock it.
3892: Notes:
3893: The function is usefull in split-phase computations, which usually have a begin phase and an end phase.
3894: One can call VecLockWriteSet_Private(x,PETSC_TRUE) in the begin phase to lock a vector for exclusive
3895: access, and call VecLockWriteSet_Private(x,PETSC_FALSE) in the end phase to unlock the vector from exclusive
3896: access. In this way, one is ensured no other operations can access the vector in between. The code may like
3899: VecGetArray(x,&xdata); // begin phase
3900: VecLockWriteSet_Private(v,PETSC_TRUE);
3902: Other operations, which can not acceess x anymore (they can access xdata, of course)
3904: VecRestoreArray(x,&vdata); // end phase
3905: VecLockWriteSet_Private(v,PETSC_FALSE);
3907: The call can not be nested on the same vector, in other words, one can not call VecLockWriteSet_Private(x,PETSC_TRUE)
3908: again before calling VecLockWriteSet_Private(v,PETSC_FALSE).
3910: Level: beginner
3912: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop(), VecLockGet()
3913: @*/
3914: PetscErrorCode VecLockWriteSet_Private(Vec x,PetscBool flg)
3915: {
3918: if (flg) {
3919: 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");
3920: 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");
3921: else x->lock = -1;
3922: } else {
3923: 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");
3924: x->lock = 0;
3925: }
3926: return(0);
3927: }
3929: /*@
3930: VecLockPush - Pushes a read-only lock on a vector to prevent it from writing
3932: Level: deprecated
3934: .seealso: VecLockReadPush()
3935: @*/
3936: PetscErrorCode VecLockPush(Vec x)
3937: {
3940: VecLockReadPush(x);
3941: return(0);
3942: }
3944: /*@
3945: VecLockPop - Pops a read-only lock from a vector
3947: Level: deprecated
3949: .seealso: VecLockReadPop()
3950: @*/
3951: PetscErrorCode VecLockPop(Vec x)
3952: {
3955: VecLockReadPop(x);
3956: return(0);
3957: }
3959: #endif