Actual source code: rvector.c
petsc-master 2019-10-13
2: /*
3: Provides the interface functions for vector operations that have PetscScalar/PetscReal in the signature
4: These are the vector functions the user calls.
5: */
6: #include <petsc/private/vecimpl.h>
7: #if defined(PETSC_HAVE_CUDA)
8: #include <../src/vec/vec/impls/dvecimpl.h>
9: #include <../src/vec/vec/impls/seq/seqcuda/cudavecimpl.h>
10: #endif
11: static PetscInt VecGetSubVectorSavedStateId = -1;
13: PETSC_EXTERN PetscErrorCode VecValidValues(Vec vec,PetscInt argnum,PetscBool begin)
14: {
15: #if defined(PETSC_USE_DEBUG)
16: PetscErrorCode ierr;
17: PetscInt n,i;
18: const PetscScalar *x;
21: #if defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_VIENNACL)
22: if ((vec->petscnative || vec->ops->getarray) && (vec->offloadmask & PETSC_OFFLOAD_CPU)) {
23: #else
24: if (vec->petscnative || vec->ops->getarray) {
25: #endif
26: VecGetLocalSize(vec,&n);
27: VecGetArrayRead(vec,&x);
28: for (i=0; i<n; i++) {
29: if (begin) {
30: if (PetscIsInfOrNanScalar(x[i])) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_FP,"Vec entry at local location %D is not-a-number or infinite at beginning of function: Parameter number %D",i,argnum);
31: } else {
32: if (PetscIsInfOrNanScalar(x[i])) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_FP,"Vec entry at local location %D is not-a-number or infinite at end of function: Parameter number %D",i,argnum);
33: }
34: }
35: VecRestoreArrayRead(vec,&x);
36: }
37: return(0);
38: #else
39: return 0;
40: #endif
41: }
43: /*@
44: VecMaxPointwiseDivide - Computes the maximum of the componentwise division max = max_i abs(x_i/y_i).
46: Logically Collective on Vec
48: Input Parameters:
49: . x, y - the vectors
51: Output Parameter:
52: . max - the result
54: Level: advanced
56: Notes:
57: x and y may be the same vector
58: if a particular y_i is zero, it is treated as 1 in the above formula
60: .seealso: VecPointwiseDivide(), VecPointwiseMult(), VecPointwiseMax(), VecPointwiseMin(), VecPointwiseMaxAbs()
61: @*/
62: PetscErrorCode VecMaxPointwiseDivide(Vec x,Vec y,PetscReal *max)
63: {
73: VecCheckSameSize(x,1,y,2);
74: (*x->ops->maxpointwisedivide)(x,y,max);
75: return(0);
76: }
78: /*@
79: VecDot - Computes the vector dot product.
81: Collective on Vec
83: Input Parameters:
84: . x, y - the vectors
86: Output Parameter:
87: . val - the dot product
89: Performance Issues:
90: $ per-processor memory bandwidth
91: $ interprocessor latency
92: $ work load inbalance that causes certain processes to arrive much earlier than others
94: Notes for Users of Complex Numbers:
95: For complex vectors, VecDot() computes
96: $ val = (x,y) = y^H x,
97: where y^H denotes the conjugate transpose of y. Note that this corresponds to the usual "mathematicians" complex
98: inner product where the SECOND argument gets the complex conjugate. Since the BLASdot() complex conjugates the first
99: first argument we call the BLASdot() with the arguments reversed.
101: Use VecTDot() for the indefinite form
102: $ val = (x,y) = y^T x,
103: where y^T denotes the transpose of y.
105: Level: intermediate
108: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDotRealPart()
109: @*/
110: PetscErrorCode VecDot(Vec x,Vec y,PetscScalar *val)
111: {
121: VecCheckSameSize(x,1,y,2);
123: PetscLogEventBegin(VEC_Dot,x,y,0,0);
124: (*x->ops->dot)(x,y,val);
125: PetscLogEventEnd(VEC_Dot,x,y,0,0);
126: return(0);
127: }
129: /*@
130: VecDotRealPart - Computes the real part of the vector dot product.
132: Collective on Vec
134: Input Parameters:
135: . x, y - the vectors
137: Output Parameter:
138: . val - the real part of the dot product;
140: Performance Issues:
141: $ per-processor memory bandwidth
142: $ interprocessor latency
143: $ work load inbalance that causes certain processes to arrive much earlier than others
145: Notes for Users of Complex Numbers:
146: See VecDot() for more details on the definition of the dot product for complex numbers
148: For real numbers this returns the same value as VecDot()
150: For complex numbers in C^n (that is a vector of n components with a complex number for each component) this is equal to the usual real dot product on the
151: the space R^{2n} (that is a vector of 2n components with the real or imaginary part of the complex numbers for components)
153: Developer Note: This is not currently optimized to compute only the real part of the dot product.
155: Level: intermediate
158: .seealso: VecMDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecDot(), VecDotNorm2()
159: @*/
160: PetscErrorCode VecDotRealPart(Vec x,Vec y,PetscReal *val)
161: {
163: PetscScalar fdot;
166: VecDot(x,y,&fdot);
167: *val = PetscRealPart(fdot);
168: return(0);
169: }
171: /*@
172: VecNorm - Computes the vector norm.
174: Collective on Vec
176: Input Parameters:
177: + x - the vector
178: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
179: NORM_1_AND_2, which computes both norms and stores them
180: in a two element array.
182: Output Parameter:
183: . val - the norm
185: Notes:
186: $ NORM_1 denotes sum_i |x_i|
187: $ NORM_2 denotes sqrt(sum_i |x_i|^2)
188: $ NORM_INFINITY denotes max_i |x_i|
190: For complex numbers NORM_1 will return the traditional 1 norm of the 2 norm of the complex numbers; that is the 1
191: norm of the absolutely values of the complex entries. In PETSc 3.6 and earlier releases it returned the 1 norm of
192: the 1 norm of the complex entries (what is returned by the BLAS routine asum()). Both are valid norms but most
193: people expect the former.
195: Level: intermediate
197: Performance Issues:
198: $ per-processor memory bandwidth
199: $ interprocessor latency
200: $ work load inbalance that causes certain processes to arrive much earlier than others
203: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNormAvailable(),
204: VecNormBegin(), VecNormEnd()
206: @*/
208: PetscErrorCode VecNorm(Vec x,NormType type,PetscReal *val)
209: {
210: PetscBool flg;
218: /*
219: * Cached data?
220: */
221: if (type!=NORM_1_AND_2) {
222: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,flg);
223: if (flg) return(0);
224: }
225: PetscLogEventBegin(VEC_Norm,x,0,0,0);
226: (*x->ops->norm)(x,type,val);
227: PetscLogEventEnd(VEC_Norm,x,0,0,0);
228: if (type!=NORM_1_AND_2) {
229: PetscObjectComposedDataSetReal((PetscObject)x,NormIds[type],*val);
230: }
231: return(0);
232: }
234: /*@
235: VecNormAvailable - Returns the vector norm if it is already known.
237: Not Collective
239: Input Parameters:
240: + x - the vector
241: - type - one of NORM_1, NORM_2, NORM_INFINITY. Also available
242: NORM_1_AND_2, which computes both norms and stores them
243: in a two element array.
245: Output Parameter:
246: + available - PETSC_TRUE if the val returned is valid
247: - val - the norm
249: Notes:
250: $ NORM_1 denotes sum_i |x_i|
251: $ NORM_2 denotes sqrt(sum_i (x_i)^2)
252: $ NORM_INFINITY denotes max_i |x_i|
254: Level: intermediate
256: Performance Issues:
257: $ per-processor memory bandwidth
258: $ interprocessor latency
259: $ work load inbalance that causes certain processes to arrive much earlier than others
261: Compile Option:
262: PETSC_HAVE_SLOW_BLAS_NORM2 will cause a C (loop unrolled) version of the norm to be used, rather
263: than the BLAS. This should probably only be used when one is using the FORTRAN BLAS routines
264: (as opposed to vendor provided) because the FORTRAN BLAS NRM2() routine is very slow.
267: .seealso: VecDot(), VecTDot(), VecNorm(), VecDotBegin(), VecDotEnd(), VecNorm()
268: VecNormBegin(), VecNormEnd()
270: @*/
271: PetscErrorCode VecNormAvailable(Vec x,NormType type,PetscBool *available,PetscReal *val)
272: {
280: *available = PETSC_FALSE;
281: if (type!=NORM_1_AND_2) {
282: PetscObjectComposedDataGetReal((PetscObject)x,NormIds[type],*val,*available);
283: }
284: return(0);
285: }
287: /*@
288: VecNormalize - Normalizes a vector by 2-norm.
290: Collective on Vec
292: Input Parameters:
293: + x - the vector
295: Output Parameter:
296: . x - the normalized vector
297: - val - the vector norm before normalization
299: Level: intermediate
302: @*/
303: PetscErrorCode VecNormalize(Vec x,PetscReal *val)
304: {
306: PetscReal norm;
311: PetscLogEventBegin(VEC_Normalize,x,0,0,0);
312: VecNorm(x,NORM_2,&norm);
313: if (norm == 0.0) {
314: PetscInfo(x,"Vector of zero norm can not be normalized; Returning only the zero norm\n");
315: } else if (norm != 1.0) {
316: PetscScalar tmp = 1.0/norm;
317: VecScale(x,tmp);
318: }
319: if (val) *val = norm;
320: PetscLogEventEnd(VEC_Normalize,x,0,0,0);
321: return(0);
322: }
324: /*@C
325: VecMax - Determines the vector component with maximum real part and its location.
327: Collective on Vec
329: Input Parameter:
330: . x - the vector
332: Output Parameters:
333: + p - the location of val (pass NULL if you don't want this)
334: - val - the maximum component
336: Notes:
337: Returns the value PETSC_MIN_REAL and p = -1 if the vector is of length 0.
339: Returns the smallest index with the maximum value
340: Level: intermediate
343: .seealso: VecNorm(), VecMin()
344: @*/
345: PetscErrorCode VecMax(Vec x,PetscInt *p,PetscReal *val)
346: {
353: PetscLogEventBegin(VEC_Max,x,0,0,0);
354: (*x->ops->max)(x,p,val);
355: PetscLogEventEnd(VEC_Max,x,0,0,0);
356: return(0);
357: }
359: /*@C
360: VecMin - Determines the vector component with minimum real part and its location.
362: Collective on Vec
364: Input Parameters:
365: . x - the vector
367: Output Parameter:
368: + p - the location of val (pass NULL if you don't want this location)
369: - val - the minimum component
371: Level: intermediate
373: Notes:
374: Returns the value PETSC_MAX_REAL and p = -1 if the vector is of length 0.
376: This returns the smallest index with the minumum value
379: .seealso: VecMax()
380: @*/
381: PetscErrorCode VecMin(Vec x,PetscInt *p,PetscReal *val)
382: {
389: PetscLogEventBegin(VEC_Min,x,0,0,0);
390: (*x->ops->min)(x,p,val);
391: PetscLogEventEnd(VEC_Min,x,0,0,0);
392: return(0);
393: }
395: /*@
396: VecTDot - Computes an indefinite vector dot product. That is, this
397: routine does NOT use the complex conjugate.
399: Collective on Vec
401: Input Parameters:
402: . x, y - the vectors
404: Output Parameter:
405: . val - the dot product
407: Notes for Users of Complex Numbers:
408: For complex vectors, VecTDot() computes the indefinite form
409: $ val = (x,y) = y^T x,
410: where y^T denotes the transpose of y.
412: Use VecDot() for the inner product
413: $ val = (x,y) = y^H x,
414: where y^H denotes the conjugate transpose of y.
416: Level: intermediate
418: .seealso: VecDot(), VecMTDot()
419: @*/
420: PetscErrorCode VecTDot(Vec x,Vec y,PetscScalar *val)
421: {
431: VecCheckSameSize(x,1,y,2);
433: PetscLogEventBegin(VEC_TDot,x,y,0,0);
434: (*x->ops->tdot)(x,y,val);
435: PetscLogEventEnd(VEC_TDot,x,y,0,0);
436: return(0);
437: }
439: /*@
440: VecScale - Scales a vector.
442: Not collective on Vec
444: Input Parameters:
445: + x - the vector
446: - alpha - the scalar
448: Output Parameter:
449: . x - the scaled vector
451: Note:
452: For a vector with n components, VecScale() computes
453: $ x[i] = alpha * x[i], for i=1,...,n.
455: Level: intermediate
458: @*/
459: PetscErrorCode VecScale(Vec x, PetscScalar alpha)
460: {
461: PetscReal norms[4] = {0.0,0.0,0.0, 0.0};
462: PetscBool flgs[4];
464: PetscInt i;
469: if (x->stash.insertmode != NOT_SET_VALUES) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled vector");
470: PetscLogEventBegin(VEC_Scale,x,0,0,0);
471: if (alpha != (PetscScalar)1.0) {
472: /* 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: VecSetErrorIfLocked(y,1);
596: VecLockReadPush(x);
597: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
598: (*y->ops->axpy)(y,alpha,x);
599: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
600: VecLockReadPop(x);
601: PetscObjectStateIncrease((PetscObject)y);
602: return(0);
603: }
605: /*@
606: VecAXPBY - Computes y = alpha x + beta y.
608: Logically Collective on Vec
610: Input Parameters:
611: + alpha,beta - the scalars
612: - x, y - the vectors
614: Output Parameter:
615: . y - output vector
617: Level: intermediate
619: Notes:
620: x and y MUST be different vectors
621: The implementation is optimized for alpha and/or beta values of 0.0 and 1.0
624: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ()
625: @*/
626: PetscErrorCode VecAXPBY(Vec y,PetscScalar alpha,PetscScalar beta,Vec x)
627: {
636: VecCheckSameSize(y,1,x,4);
637: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
640: PetscLogEventBegin(VEC_AXPY,x,y,0,0);
641: (*y->ops->axpby)(y,alpha,beta,x);
642: PetscLogEventEnd(VEC_AXPY,x,y,0,0);
643: PetscObjectStateIncrease((PetscObject)y);
644: return(0);
645: }
647: /*@
648: VecAXPBYPCZ - Computes z = alpha x + beta y + gamma z
650: Logically Collective on Vec
652: Input Parameters:
653: + alpha,beta, gamma - the scalars
654: - x, y, z - the vectors
656: Output Parameter:
657: . z - output vector
659: Level: intermediate
661: Notes:
662: x, y and z must be different vectors
663: The implementation is optimized for alpha of 1.0 and gamma of 1.0 or 0.0
666: .seealso: VecAYPX(), VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBY()
667: @*/
668: PetscErrorCode VecAXPBYPCZ(Vec z,PetscScalar alpha,PetscScalar beta,PetscScalar gamma,Vec x,Vec y)
669: {
681: VecCheckSameSize(x,1,y,5);
682: VecCheckSameSize(x,1,z,6);
683: if (x == y || x == z) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
684: if (y == z) SETERRQ(PetscObjectComm((PetscObject)y),PETSC_ERR_ARG_IDN,"x, y, and z must be different vectors");
689: PetscLogEventBegin(VEC_AXPBYPCZ,x,y,z,0);
690: (*y->ops->axpbypcz)(z,alpha,beta,gamma,x,y);
691: PetscLogEventEnd(VEC_AXPBYPCZ,x,y,z,0);
692: PetscObjectStateIncrease((PetscObject)z);
693: return(0);
694: }
696: /*@
697: VecAYPX - Computes y = x + beta y.
699: Logically Collective on Vec
701: Input Parameters:
702: + beta - the scalar
703: - x, y - the vectors
705: Output Parameter:
706: . y - output vector
708: Level: intermediate
710: Notes:
711: x and y MUST be different vectors
712: The implementation is optimized for beta of -1.0, 0.0, and 1.0
715: .seealso: VecMAXPY(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
716: @*/
717: PetscErrorCode VecAYPX(Vec y,PetscScalar beta,Vec x)
718: {
727: VecCheckSameSize(x,1,y,3);
728: if (x == y) SETERRQ(PetscObjectComm((PetscObject)x),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
731: PetscLogEventBegin(VEC_AYPX,x,y,0,0);
732: (*y->ops->aypx)(y,beta,x);
733: PetscLogEventEnd(VEC_AYPX,x,y,0,0);
734: PetscObjectStateIncrease((PetscObject)y);
735: return(0);
736: }
739: /*@
740: VecWAXPY - Computes w = alpha x + y.
742: Logically Collective on Vec
744: Input Parameters:
745: + alpha - the scalar
746: - x, y - the vectors
748: Output Parameter:
749: . w - the result
751: Level: intermediate
753: Notes:
754: w cannot be either x or y, but x and y can be the same
755: The implementation is optimzed for alpha of -1.0, 0.0, and 1.0
758: .seealso: VecAXPY(), VecAYPX(), VecAXPBY(), VecMAXPY(), VecAXPBYPCZ()
759: @*/
760: PetscErrorCode VecWAXPY(Vec w,PetscScalar alpha,Vec x,Vec y)
761: {
773: VecCheckSameSize(x,3,y,4);
774: VecCheckSameSize(x,3,w,1);
775: if (w == y) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector y, suggest VecAXPY()");
776: if (w == x) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Result vector w cannot be same as input vector x, suggest VecAYPX()");
779: PetscLogEventBegin(VEC_WAXPY,x,y,w,0);
780: (*w->ops->waxpy)(w,alpha,x,y);
781: PetscLogEventEnd(VEC_WAXPY,x,y,w,0);
782: PetscObjectStateIncrease((PetscObject)w);
783: return(0);
784: }
787: /*@C
788: VecSetValues - Inserts or adds values into certain locations of a vector.
790: Not Collective
792: Input Parameters:
793: + x - vector to insert in
794: . ni - number of elements to add
795: . ix - indices where to add
796: . y - array of values
797: - iora - either INSERT_VALUES or ADD_VALUES, where
798: ADD_VALUES adds values to any existing entries, and
799: INSERT_VALUES replaces existing entries with new values
801: Notes:
802: VecSetValues() sets x[ix[i]] = y[i], for i=0,...,ni-1.
804: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
805: options cannot be mixed without intervening calls to the assembly
806: routines.
808: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
809: MUST be called after all calls to VecSetValues() have been completed.
811: VecSetValues() uses 0-based indices in Fortran as well as in C.
813: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
814: negative indices may be passed in ix. These rows are
815: simply ignored. This allows easily inserting element load matrices
816: with homogeneous Dirchlet boundary conditions that you don't want represented
817: in the vector.
819: Level: beginner
821: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesLocal(),
822: VecSetValue(), VecSetValuesBlocked(), InsertMode, INSERT_VALUES, ADD_VALUES, VecGetValues()
823: @*/
824: PetscErrorCode VecSetValues(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
825: {
830: if (!ni) return(0);
834: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
835: (*x->ops->setvalues)(x,ni,ix,y,iora);
836: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
837: PetscObjectStateIncrease((PetscObject)x);
838: return(0);
839: }
841: /*@
842: VecGetValues - Gets values from certain locations of a vector. Currently
843: can only get values on the same processor
845: Not Collective
847: Input Parameters:
848: + x - vector to get values from
849: . ni - number of elements to get
850: - ix - indices where to get them from (in global 1d numbering)
852: Output Parameter:
853: . y - array of values
855: Notes:
856: The user provides the allocated array y; it is NOT allocated in this routine
858: VecGetValues() gets y[i] = x[ix[i]], for i=0,...,ni-1.
860: VecAssemblyBegin() and VecAssemblyEnd() MUST be called before calling this
862: VecGetValues() uses 0-based indices in Fortran as well as in C.
864: If you call VecSetOption(x, VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE),
865: negative indices may be passed in ix. These rows are
866: simply ignored.
868: Level: beginner
870: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues()
871: @*/
872: PetscErrorCode VecGetValues(Vec x,PetscInt ni,const PetscInt ix[],PetscScalar y[])
873: {
878: if (!ni) return(0);
882: (*x->ops->getvalues)(x,ni,ix,y);
883: return(0);
884: }
886: /*@C
887: VecSetValuesBlocked - Inserts or adds blocks of values into certain locations of a vector.
889: Not Collective
891: Input Parameters:
892: + x - vector to insert in
893: . ni - number of blocks to add
894: . ix - indices where to add in block count, rather than element count
895: . y - array of values
896: - iora - either INSERT_VALUES or ADD_VALUES, where
897: ADD_VALUES adds values to any existing entries, and
898: INSERT_VALUES replaces existing entries with new values
900: Notes:
901: VecSetValuesBlocked() sets x[bs*ix[i]+j] = y[bs*i+j],
902: for j=0,...,bs-1, for i=0,...,ni-1. where bs was set with VecSetBlockSize().
904: Calls to VecSetValuesBlocked() with the INSERT_VALUES and ADD_VALUES
905: options cannot be mixed without intervening calls to the assembly
906: routines.
908: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
909: MUST be called after all calls to VecSetValuesBlocked() have been completed.
911: VecSetValuesBlocked() uses 0-based indices in Fortran as well as in C.
913: Negative indices may be passed in ix, these rows are
914: simply ignored. This allows easily inserting element load matrices
915: with homogeneous Dirchlet boundary conditions that you don't want represented
916: in the vector.
918: Level: intermediate
920: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValuesBlockedLocal(),
921: VecSetValues()
922: @*/
923: PetscErrorCode VecSetValuesBlocked(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
924: {
929: if (!ni) return(0);
933: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
934: (*x->ops->setvaluesblocked)(x,ni,ix,y,iora);
935: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
936: PetscObjectStateIncrease((PetscObject)x);
937: return(0);
938: }
941: /*@C
942: VecSetValuesLocal - Inserts or adds values into certain locations of a vector,
943: using a local ordering of the nodes.
945: Not Collective
947: Input Parameters:
948: + x - vector to insert in
949: . ni - number of elements to add
950: . ix - indices where to add
951: . y - array of values
952: - iora - either INSERT_VALUES or ADD_VALUES, where
953: ADD_VALUES adds values to any existing entries, and
954: INSERT_VALUES replaces existing entries with new values
956: Level: intermediate
958: Notes:
959: VecSetValuesLocal() sets x[ix[i]] = y[i], for i=0,...,ni-1.
961: Calls to VecSetValues() with the INSERT_VALUES and ADD_VALUES
962: options cannot be mixed without intervening calls to the assembly
963: routines.
965: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
966: MUST be called after all calls to VecSetValuesLocal() have been completed.
968: VecSetValuesLocal() uses 0-based indices in Fortran as well as in C.
970: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetLocalToGlobalMapping(),
971: VecSetValuesBlockedLocal()
972: @*/
973: PetscErrorCode VecSetValuesLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
974: {
976: PetscInt lixp[128],*lix = lixp;
980: if (!ni) return(0);
985: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
986: if (!x->ops->setvalueslocal) {
987: if (!x->map->mapping) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Local to global never set with VecSetLocalToGlobalMapping()");
988: if (ni > 128) {
989: PetscMalloc1(ni,&lix);
990: }
991: ISLocalToGlobalMappingApply(x->map->mapping,ni,(PetscInt*)ix,lix);
992: (*x->ops->setvalues)(x,ni,lix,y,iora);
993: if (ni > 128) {
994: PetscFree(lix);
995: }
996: } else {
997: (*x->ops->setvalueslocal)(x,ni,ix,y,iora);
998: }
999: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1000: PetscObjectStateIncrease((PetscObject)x);
1001: return(0);
1002: }
1004: /*@
1005: VecSetValuesBlockedLocal - Inserts or adds values into certain locations of a vector,
1006: using a local ordering of the nodes.
1008: Not Collective
1010: Input Parameters:
1011: + x - vector to insert in
1012: . ni - number of blocks to add
1013: . ix - indices where to add in block count, not element count
1014: . y - array of values
1015: - iora - either INSERT_VALUES or ADD_VALUES, where
1016: ADD_VALUES adds values to any existing entries, and
1017: INSERT_VALUES replaces existing entries with new values
1019: Level: intermediate
1021: Notes:
1022: VecSetValuesBlockedLocal() sets x[bs*ix[i]+j] = y[bs*i+j],
1023: for j=0,..bs-1, for i=0,...,ni-1, where bs has been set with VecSetBlockSize().
1025: Calls to VecSetValuesBlockedLocal() with the INSERT_VALUES and ADD_VALUES
1026: options cannot be mixed without intervening calls to the assembly
1027: routines.
1029: These values may be cached, so VecAssemblyBegin() and VecAssemblyEnd()
1030: MUST be called after all calls to VecSetValuesBlockedLocal() have been completed.
1032: VecSetValuesBlockedLocal() uses 0-based indices in Fortran as well as in C.
1035: .seealso: VecAssemblyBegin(), VecAssemblyEnd(), VecSetValues(), VecSetValuesBlocked(),
1036: VecSetLocalToGlobalMapping()
1037: @*/
1038: PetscErrorCode VecSetValuesBlockedLocal(Vec x,PetscInt ni,const PetscInt ix[],const PetscScalar y[],InsertMode iora)
1039: {
1041: PetscInt lixp[128],*lix = lixp;
1045: if (!ni) return(0);
1049: if (ni > 128) {
1050: PetscMalloc1(ni,&lix);
1051: }
1053: PetscLogEventBegin(VEC_SetValues,x,0,0,0);
1054: ISLocalToGlobalMappingApplyBlock(x->map->mapping,ni,(PetscInt*)ix,lix);
1055: (*x->ops->setvaluesblocked)(x,ni,lix,y,iora);
1056: PetscLogEventEnd(VEC_SetValues,x,0,0,0);
1057: if (ni > 128) {
1058: PetscFree(lix);
1059: }
1060: PetscObjectStateIncrease((PetscObject)x);
1061: return(0);
1062: }
1064: /*@
1065: VecMTDot - Computes indefinite vector multiple dot products.
1066: That is, it does NOT use the complex conjugate.
1068: Collective on Vec
1070: Input Parameters:
1071: + x - one vector
1072: . nv - number of vectors
1073: - y - array of vectors. Note that vectors are pointers
1075: Output Parameter:
1076: . val - array of the dot products
1078: Notes for Users of Complex Numbers:
1079: For complex vectors, VecMTDot() computes the indefinite form
1080: $ val = (x,y) = y^T x,
1081: where y^T denotes the transpose of y.
1083: Use VecMDot() for the inner product
1084: $ val = (x,y) = y^H x,
1085: where y^H denotes the conjugate transpose of y.
1087: Level: intermediate
1090: .seealso: VecMDot(), VecTDot()
1091: @*/
1092: PetscErrorCode VecMTDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1093: {
1099: if (!nv) return(0);
1106: VecCheckSameSize(x,1,*y,3);
1108: PetscLogEventBegin(VEC_MTDot,x,*y,0,0);
1109: (*x->ops->mtdot)(x,nv,y,val);
1110: PetscLogEventEnd(VEC_MTDot,x,*y,0,0);
1111: return(0);
1112: }
1114: /*@
1115: VecMDot - Computes vector multiple dot products.
1117: Collective on Vec
1119: Input Parameters:
1120: + x - one vector
1121: . nv - number of vectors
1122: - y - array of vectors.
1124: Output Parameter:
1125: . val - array of the dot products (does not allocate the array)
1127: Notes for Users of Complex Numbers:
1128: For complex vectors, VecMDot() computes
1129: $ val = (x,y) = y^H x,
1130: where y^H denotes the conjugate transpose of y.
1132: Use VecMTDot() for the indefinite form
1133: $ val = (x,y) = y^T x,
1134: where y^T denotes the transpose of y.
1136: Level: intermediate
1139: .seealso: VecMTDot(), VecDot()
1140: @*/
1141: PetscErrorCode VecMDot(Vec x,PetscInt nv,const Vec y[],PetscScalar val[])
1142: {
1148: if (!nv) return(0);
1149: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1156: VecCheckSameSize(x,1,*y,3);
1158: PetscLogEventBegin(VEC_MDot,x,*y,0,0);
1159: (*x->ops->mdot)(x,nv,y,val);
1160: PetscLogEventEnd(VEC_MDot,x,*y,0,0);
1161: return(0);
1162: }
1164: /*@
1165: VecMAXPY - Computes y = y + sum alpha[i] x[i]
1167: Logically Collective on Vec
1169: Input Parameters:
1170: + nv - number of scalars and x-vectors
1171: . alpha - array of scalars
1172: . y - one vector
1173: - x - array of vectors
1175: Level: intermediate
1177: Notes:
1178: y cannot be any of the x vectors
1180: .seealso: VecAYPX(), VecWAXPY(), VecAXPY(), VecAXPBYPCZ(), VecAXPBY()
1181: @*/
1182: PetscErrorCode VecMAXPY(Vec y,PetscInt nv,const PetscScalar alpha[],Vec x[])
1183: {
1185: PetscInt i;
1190: if (!nv) return(0);
1191: if (nv < 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of vectors (given %D) cannot be negative",nv);
1198: VecCheckSameSize(y,1,*x,4);
1201: PetscLogEventBegin(VEC_MAXPY,*x,y,0,0);
1202: (*y->ops->maxpy)(y,nv,alpha,x);
1203: PetscLogEventEnd(VEC_MAXPY,*x,y,0,0);
1204: PetscObjectStateIncrease((PetscObject)y);
1205: return(0);
1206: }
1208: /*@
1209: VecGetSubVector - Gets a vector representing part of another vector
1211: Collective on IS
1213: Input Arguments:
1214: + X - vector from which to extract a subvector
1215: - is - index set representing portion of X to extract
1217: Output Arguments:
1218: . Y - subvector corresponding to is
1220: Level: advanced
1222: Notes:
1223: The subvector Y should be returned with VecRestoreSubVector().
1225: This function may return a subvector without making a copy, therefore it is not safe to use the original vector while
1226: modifying the subvector. Other non-overlapping subvectors can still be obtained from X using this function.
1228: .seealso: MatCreateSubMatrix()
1229: @*/
1230: PetscErrorCode VecGetSubVector(Vec X,IS is,Vec *Y)
1231: {
1232: PetscErrorCode ierr;
1233: Vec Z;
1239: if (X->ops->getsubvector) {
1240: (*X->ops->getsubvector)(X,is,&Z);
1241: } else { /* Default implementation currently does no caching */
1242: PetscInt gstart,gend,start;
1243: PetscBool contiguous,gcontiguous;
1244: VecGetOwnershipRange(X,&gstart,&gend);
1245: ISContiguousLocal(is,gstart,gend,&start,&contiguous);
1246: MPIU_Allreduce(&contiguous,&gcontiguous,1,MPIU_BOOL,MPI_LAND,PetscObjectComm((PetscObject)is));
1247: if (gcontiguous) { /* We can do a no-copy implementation */
1248: PetscInt n,N,bs;
1249: PetscInt state;
1251: ISGetSize(is,&N);
1252: ISGetLocalSize(is,&n);
1253: VecGetBlockSize(X,&bs);
1254: if (n%bs || bs == 1 || !n) bs = -1; /* Do not decide block size if we do not have to */
1255: VecLockGet(X,&state);
1256: if (state) {
1257: const PetscScalar *x;
1258: VecGetArrayRead(X,&x);
1259: VecCreate(PetscObjectComm((PetscObject)X),&Z);
1260: VecSetType(Z,((PetscObject)X)->type_name);
1261: VecSetSizes(Z,n,N);
1262: VecSetBlockSize(Z,bs);
1263: VecPlaceArray(Z,(PetscScalar*)x+start);
1264: VecLockReadPush(Z);
1265: VecRestoreArrayRead(X,&x);
1266: } else {
1267: PetscScalar *x;
1268: VecGetArray(X,&x);
1269: VecCreate(PetscObjectComm((PetscObject)X),&Z);
1270: VecSetType(Z,((PetscObject)X)->type_name);
1271: VecSetSizes(Z,n,N);
1272: VecSetBlockSize(Z,bs);
1273: VecPlaceArray(Z,(PetscScalar*)x+start);
1274: VecRestoreArray(X,&x);
1275: }
1276: } else { /* Have to create a scatter and do a copy */
1277: VecScatter scatter;
1278: PetscInt n,N;
1279: ISGetLocalSize(is,&n);
1280: ISGetSize(is,&N);
1281: VecCreate(PetscObjectComm((PetscObject)is),&Z);
1282: VecSetSizes(Z,n,N);
1283: VecSetType(Z,((PetscObject)X)->type_name);
1284: VecScatterCreate(X,is,Z,NULL,&scatter);
1285: VecScatterBegin(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1286: VecScatterEnd(scatter,X,Z,INSERT_VALUES,SCATTER_FORWARD);
1287: PetscObjectCompose((PetscObject)Z,"VecGetSubVector_Scatter",(PetscObject)scatter);
1288: VecScatterDestroy(&scatter);
1289: }
1290: }
1291: /* Record the state when the subvector was gotten so we know whether its values need to be put back */
1292: if (VecGetSubVectorSavedStateId < 0) {PetscObjectComposedDataRegister(&VecGetSubVectorSavedStateId);}
1293: PetscObjectComposedDataSetInt((PetscObject)Z,VecGetSubVectorSavedStateId,1);
1294: *Y = Z;
1295: return(0);
1296: }
1298: /*@
1299: VecRestoreSubVector - Restores a subvector extracted using VecGetSubVector()
1301: Collective on IS
1303: Input Arguments:
1304: + X - vector from which subvector was obtained
1305: . is - index set representing the subset of X
1306: - Y - subvector being restored
1308: Level: advanced
1310: .seealso: VecGetSubVector()
1311: @*/
1312: PetscErrorCode VecRestoreSubVector(Vec X,IS is,Vec *Y)
1313: {
1321: if (X->ops->restoresubvector) {
1322: (*X->ops->restoresubvector)(X,is,Y);
1323: } else {
1324: PETSC_UNUSED PetscObjectState dummystate = 0;
1325: PetscBool valid;
1326: PetscObjectComposedDataGetInt((PetscObject)*Y,VecGetSubVectorSavedStateId,dummystate,valid);
1327: if (!valid) {
1328: VecScatter scatter;
1330: PetscObjectQuery((PetscObject)*Y,"VecGetSubVector_Scatter",(PetscObject*)&scatter);
1331: if (scatter) {
1332: VecScatterBegin(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1333: VecScatterEnd(scatter,*Y,X,INSERT_VALUES,SCATTER_REVERSE);
1334: }
1335: }
1336: VecDestroy(Y);
1337: }
1338: return(0);
1339: }
1341: /*@
1342: VecGetLocalVectorRead - Maps the local portion of a vector into a
1343: vector. You must call VecRestoreLocalVectorRead() when the local
1344: vector is no longer needed.
1346: Not collective.
1348: Input parameter:
1349: . v - The vector for which the local vector is desired.
1351: Output parameter:
1352: . w - Upon exit this contains the local vector.
1354: Level: beginner
1356: Notes:
1357: This function is similar to VecGetArrayRead() which maps the local
1358: portion into a raw pointer. VecGetLocalVectorRead() is usually
1359: almost as efficient as VecGetArrayRead() but in certain circumstances
1360: VecGetLocalVectorRead() can be much more efficient than
1361: VecGetArrayRead(). This is because the construction of a contiguous
1362: array representing the vector data required by VecGetArrayRead() can
1363: be an expensive operation for certain vector types. For example, for
1364: GPU vectors VecGetArrayRead() requires that the data between device
1365: and host is synchronized.
1367: Unlike VecGetLocalVector(), this routine is not collective and
1368: preserves cached information.
1370: .seealso: VecRestoreLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1371: @*/
1372: PetscErrorCode VecGetLocalVectorRead(Vec v,Vec w)
1373: {
1375: PetscScalar *a;
1380: VecCheckSameLocalSize(v,1,w,2);
1381: if (v->ops->getlocalvectorread) {
1382: (*v->ops->getlocalvectorread)(v,w);
1383: } else {
1384: VecGetArrayRead(v,(const PetscScalar**)&a);
1385: VecPlaceArray(w,a);
1386: }
1387: return(0);
1388: }
1390: /*@
1391: VecRestoreLocalVectorRead - Unmaps the local portion of a vector
1392: previously mapped into a vector using VecGetLocalVectorRead().
1394: Not collective.
1396: Input parameter:
1397: + v - The local portion of this vector was previously mapped into w using VecGetLocalVectorRead().
1398: - w - The vector into which the local portion of v was mapped.
1400: Level: beginner
1402: .seealso: VecGetLocalVectorRead(), VecGetLocalVector(), VecGetArrayRead(), VecGetArray()
1403: @*/
1404: PetscErrorCode VecRestoreLocalVectorRead(Vec v,Vec w)
1405: {
1407: PetscScalar *a;
1412: if (v->ops->restorelocalvectorread) {
1413: (*v->ops->restorelocalvectorread)(v,w);
1414: } else {
1415: VecGetArrayRead(w,(const PetscScalar**)&a);
1416: VecRestoreArrayRead(v,(const PetscScalar**)&a);
1417: VecResetArray(w);
1418: }
1419: return(0);
1420: }
1422: /*@
1423: VecGetLocalVector - Maps the local portion of a vector into a
1424: vector.
1426: Collective on v, not collective on w.
1428: Input parameter:
1429: . v - The vector for which the local vector is desired.
1431: Output parameter:
1432: . w - Upon exit this contains the local vector.
1434: Level: beginner
1436: Notes:
1437: This function is similar to VecGetArray() which maps the local
1438: portion into a raw pointer. VecGetLocalVector() is usually about as
1439: efficient as VecGetArray() but in certain circumstances
1440: VecGetLocalVector() can be much more efficient than VecGetArray().
1441: This is because the construction of a contiguous array representing
1442: the vector data required by VecGetArray() can be an expensive
1443: operation for certain vector types. For example, for GPU vectors
1444: VecGetArray() requires that the data between device and host is
1445: synchronized.
1447: .seealso: VecRestoreLocalVector(), VecGetLocalVectorRead(), VecGetArrayRead(), VecGetArray()
1448: @*/
1449: PetscErrorCode VecGetLocalVector(Vec v,Vec w)
1450: {
1452: PetscScalar *a;
1457: VecCheckSameLocalSize(v,1,w,2);
1458: if (v->ops->getlocalvector) {
1459: (*v->ops->getlocalvector)(v,w);
1460: } else {
1461: VecGetArray(v,&a);
1462: VecPlaceArray(w,a);
1463: }
1464: return(0);
1465: }
1467: /*@
1468: VecRestoreLocalVector - Unmaps the local portion of a vector
1469: previously mapped into a vector using VecGetLocalVector().
1471: Logically collective.
1473: Input parameter:
1474: + v - The local portion of this vector was previously mapped into w using VecGetLocalVector().
1475: - w - The vector into which the local portion of v was mapped.
1477: Level: beginner
1479: .seealso: VecGetLocalVector(), VecGetLocalVectorRead(), VecRestoreLocalVectorRead(), LocalVectorRead(), VecGetArrayRead(), VecGetArray()
1480: @*/
1481: PetscErrorCode VecRestoreLocalVector(Vec v,Vec w)
1482: {
1484: PetscScalar *a;
1489: if (v->ops->restorelocalvector) {
1490: (*v->ops->restorelocalvector)(v,w);
1491: } else {
1492: VecGetArray(w,&a);
1493: VecRestoreArray(v,&a);
1494: VecResetArray(w);
1495: }
1496: return(0);
1497: }
1499: /*@C
1500: VecGetArray - Returns a pointer to a contiguous array that contains this
1501: processor's portion of the vector data. For the standard PETSc
1502: vectors, VecGetArray() returns a pointer to the local data array and
1503: does not use any copies. If the underlying vector data is not stored
1504: in a contiguous array this routine will copy the data to a contiguous
1505: array and return a pointer to that. You MUST call VecRestoreArray()
1506: when you no longer need access to the array.
1508: Logically Collective on Vec
1510: Input Parameter:
1511: . x - the vector
1513: Output Parameter:
1514: . a - location to put pointer to the array
1516: Fortran Note:
1517: This routine is used differently from Fortran 77
1518: $ Vec x
1519: $ PetscScalar x_array(1)
1520: $ PetscOffset i_x
1521: $ PetscErrorCode ierr
1522: $ call VecGetArray(x,x_array,i_x,ierr)
1523: $
1524: $ Access first local entry in vector with
1525: $ value = x_array(i_x + 1)
1526: $
1527: $ ...... other code
1528: $ call VecRestoreArray(x,x_array,i_x,ierr)
1529: For Fortran 90 see VecGetArrayF90()
1531: See the Fortran chapter of the users manual and
1532: petsc/src/snes/examples/tutorials/ex5f.F for details.
1534: Level: beginner
1536: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1537: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
1538: @*/
1539: PetscErrorCode VecGetArray(Vec x,PetscScalar **a)
1540: {
1542: #if defined(PETSC_HAVE_VIENNACL)
1543: PetscBool is_viennacltype = PETSC_FALSE;
1544: #endif
1548: VecSetErrorIfLocked(x,1);
1549: if (x->petscnative) {
1550: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1551: if (x->offloadmask == PETSC_OFFLOAD_GPU) {
1552: #if defined(PETSC_HAVE_VIENNACL)
1553: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1554: if (is_viennacltype) {
1555: VecViennaCLCopyFromGPU(x);
1556: } else
1557: #endif
1558: {
1559: #if defined(PETSC_HAVE_CUDA)
1560: VecCUDACopyFromGPU(x);
1561: #endif
1562: }
1563: } else if (x->offloadmask == PETSC_OFFLOAD_UNALLOCATED) {
1564: #if defined(PETSC_HAVE_VIENNACL)
1565: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1566: if (is_viennacltype) {
1567: VecViennaCLAllocateCheckHost(x);
1568: } else
1569: #endif
1570: {
1571: #if defined(PETSC_HAVE_CUDA)
1572: VecCUDAAllocateCheckHost(x);
1573: #endif
1574: }
1575: }
1576: #endif
1577: *a = *((PetscScalar**)x->data);
1578: } else {
1579: if (x->ops->getarray) {
1580: (*x->ops->getarray)(x,a);
1581: } else SETERRQ1(PetscObjectComm((PetscObject)x),PETSC_ERR_SUP,"Cannot get array for vector type \"%s\"",((PetscObject)x)->type_name);
1582: }
1583: return(0);
1584: }
1586: /*@C
1587: VecGetArrayInPlace - Like VecGetArray(), but if this is a CUDA vector and it is currently offloaded to GPU,
1588: the returned pointer will be a GPU pointer to the GPU memory that contains this processor's portion of the
1589: vector data. Otherwise, it functions as VecGetArray().
1591: Logically Collective on Vec
1593: Input Parameter:
1594: . x - the vector
1596: Output Parameter:
1597: . a - location to put pointer to the array
1599: Level: beginner
1601: .seealso: VecRestoreArrayInPlace(), VecRestoreArrayInPlace(), VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(),
1602: VecPlaceArray(), VecGetArray2d(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite(), VecRestoreArrayWrite()
1603: @*/
1604: PetscErrorCode VecGetArrayInPlace(Vec x,PetscScalar **a)
1605: {
1610: VecSetErrorIfLocked(x,1);
1612: #if defined(PETSC_HAVE_CUDA)
1613: if (x->petscnative && (x->offloadmask & PETSC_OFFLOAD_GPU)) { /* Prefer working on GPU when offloadmask is PETSC_OFFLOAD_BOTH */
1614: PetscBool is_cudatype = PETSC_FALSE;
1615: PetscObjectTypeCompareAny((PetscObject)x,&is_cudatype,VECSEQCUDA,VECMPICUDA,VECCUDA,"");
1616: if (is_cudatype) {
1617: VecCUDAGetArray(x,a);
1618: x->offloadmask = PETSC_OFFLOAD_GPU; /* Change the mask once GPU gets write access, don't wait until restore array */
1619: return(0);
1620: }
1621: }
1622: #endif
1623: VecGetArray(x,a);
1624: return(0);
1625: }
1627: /*@C
1628: VecGetArrayWrite - Returns a pointer to a contiguous array that WILL contains this
1629: processor's portion of the vector data. The values in this array are NOT valid, the routine calling this
1630: routine is responsible for putting values into the array; any values it does not set will be invalid
1632: Logically Collective on Vec
1634: Input Parameter:
1635: . x - the vector
1637: Output Parameter:
1638: . a - location to put pointer to the array
1640: Level: intermediate
1642: This is for vectors associate with GPUs, the vector is not copied up before giving access. If you need correct
1643: values in the array use VecGetArray()
1645: Concepts: vector^accessing local values
1647: .seealso: VecRestoreArray(), VecGetArrayRead(), VecGetArrays(), VecGetArrayF90(), VecGetArrayReadF90(), VecPlaceArray(), VecGetArray2d(),
1648: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArray(), VecRestoreArrayWrite()
1649: @*/
1650: PetscErrorCode VecGetArrayWrite(Vec x,PetscScalar **a)
1651: {
1656: VecSetErrorIfLocked(x,1);
1657: if (!x->ops->getarraywrite) {
1658: VecGetArray(x,a);
1659: } else {
1660: (*x->ops->getarraywrite)(x,a);
1661: }
1662: return(0);
1663: }
1665: /*@C
1666: VecGetArrayRead - Get read-only pointer to contiguous array containing this processor's portion of the vector data.
1668: Not Collective
1670: Input Parameters:
1671: . x - the vector
1673: Output Parameter:
1674: . a - the array
1676: Level: beginner
1678: Notes:
1679: The array must be returned using a matching call to VecRestoreArrayRead().
1681: Unlike VecGetArray(), this routine is not collective and preserves cached information like vector norms.
1683: Standard PETSc vectors use contiguous storage so that this routine does not perform a copy. Other vector
1684: implementations may require a copy, but must such implementations should cache the contiguous representation so that
1685: only one copy is performed when this routine is called multiple times in sequence.
1687: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1688: @*/
1689: PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar **a)
1690: {
1692: #if defined(PETSC_HAVE_VIENNACL)
1693: PetscBool is_viennacltype = PETSC_FALSE;
1694: #endif
1698: if (x->petscnative) {
1699: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1700: if (x->offloadmask == PETSC_OFFLOAD_GPU) {
1701: #if defined(PETSC_HAVE_VIENNACL)
1702: PetscObjectTypeCompareAny((PetscObject)x,&is_viennacltype,VECSEQVIENNACL,VECMPIVIENNACL,VECVIENNACL,"");
1703: if (is_viennacltype) {
1704: VecViennaCLCopyFromGPU(x);
1705: } else
1706: #endif
1707: {
1708: #if defined(PETSC_HAVE_CUDA)
1709: VecCUDACopyFromGPU(x);
1710: #endif
1711: }
1712: }
1713: #endif
1714: *a = *((PetscScalar **)x->data);
1715: } else if (x->ops->getarrayread) {
1716: (*x->ops->getarrayread)(x,a);
1717: } else {
1718: (*x->ops->getarray)(x,(PetscScalar**)a);
1719: }
1720: return(0);
1721: }
1723: /*@C
1724: VecGetArrayReadInPlace - Like VecGetArrayRead(), but if this is a CUDA vector and it is currently offloaded to GPU,
1725: the returned pointer will be a GPU pointer to the GPU memory that contains this processor's portion of the
1726: vector data. Otherwise, it functions as VecGetArrayRead().
1728: Not Collective
1730: Input Parameters:
1731: . x - the vector
1733: Output Parameter:
1734: . a - the array
1736: Level: beginner
1738: Notes:
1739: The array must be returned using a matching call to VecRestoreArrayReadInPlace().
1742: .seealso: VecRestoreArrayReadInPlace(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayInPlace()
1743: @*/
1744: PetscErrorCode VecGetArrayReadInPlace(Vec x,const PetscScalar **a)
1745: {
1750: #if defined(PETSC_HAVE_CUDA)
1751: if (x->petscnative && x->offloadmask & PETSC_OFFLOAD_GPU) {
1752: PetscBool is_cudatype = PETSC_FALSE;
1753: PetscObjectTypeCompareAny((PetscObject)x,&is_cudatype,VECSEQCUDA,VECMPICUDA,VECCUDA,"");
1754: if (is_cudatype) {
1755: VecCUDAGetArrayRead(x,a);
1756: return(0);
1757: }
1758: }
1759: #endif
1760: VecGetArrayRead(x,a);
1761: return(0);
1762: }
1764: /*@C
1765: VecGetArrays - Returns a pointer to the arrays in a set of vectors
1766: that were created by a call to VecDuplicateVecs(). You MUST call
1767: VecRestoreArrays() when you no longer need access to the array.
1769: Logically Collective on Vec
1771: Input Parameter:
1772: + x - the vectors
1773: - n - the number of vectors
1775: Output Parameter:
1776: . a - location to put pointer to the array
1778: Fortran Note:
1779: This routine is not supported in Fortran.
1781: Level: intermediate
1783: .seealso: VecGetArray(), VecRestoreArrays()
1784: @*/
1785: PetscErrorCode VecGetArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1786: {
1788: PetscInt i;
1789: PetscScalar **q;
1795: if (n <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must get at least one array n = %D",n);
1796: PetscMalloc1(n,&q);
1797: for (i=0; i<n; ++i) {
1798: VecGetArray(x[i],&q[i]);
1799: }
1800: *a = q;
1801: return(0);
1802: }
1804: /*@C
1805: VecRestoreArrays - Restores a group of vectors after VecGetArrays()
1806: has been called.
1808: Logically Collective on Vec
1810: Input Parameters:
1811: + x - the vector
1812: . n - the number of vectors
1813: - a - location of pointer to arrays obtained from VecGetArrays()
1815: Notes:
1816: For regular PETSc vectors this routine does not involve any copies. For
1817: any special vectors that do not store local vector data in a contiguous
1818: array, this routine will copy the data back into the underlying
1819: vector data structure from the arrays obtained with VecGetArrays().
1821: Fortran Note:
1822: This routine is not supported in Fortran.
1824: Level: intermediate
1826: .seealso: VecGetArrays(), VecRestoreArray()
1827: @*/
1828: PetscErrorCode VecRestoreArrays(const Vec x[],PetscInt n,PetscScalar **a[])
1829: {
1831: PetscInt i;
1832: PetscScalar **q = *a;
1839: for (i=0; i<n; ++i) {
1840: VecRestoreArray(x[i],&q[i]);
1841: }
1842: PetscFree(q);
1843: return(0);
1844: }
1846: /*@C
1847: VecRestoreArray - Restores a vector after VecGetArray() has been called.
1849: Logically Collective on Vec
1851: Input Parameters:
1852: + x - the vector
1853: - a - location of pointer to array obtained from VecGetArray()
1855: Level: beginner
1857: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1858: VecGetArrayPair(), VecRestoreArrayPair()
1859: @*/
1860: PetscErrorCode VecRestoreArray(Vec x,PetscScalar **a)
1861: {
1866: if (x->petscnative) {
1867: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1868: x->offloadmask = PETSC_OFFLOAD_CPU;
1869: #endif
1870: } else {
1871: (*x->ops->restorearray)(x,a);
1872: }
1873: if (a) *a = NULL;
1874: PetscObjectStateIncrease((PetscObject)x);
1875: return(0);
1876: }
1878: /*@C
1879: VecRestoreArrayInPlace - Restores a vector after VecGetArrayInPlace() has been called.
1881: Logically Collective on Vec
1883: Input Parameters:
1884: + x - the vector
1885: - a - location of pointer to array obtained from VecGetArrayInPlace()
1887: Level: beginner
1889: .seealso: VecGetArrayInPlace(), VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(),
1890: VecPlaceArray(), VecRestoreArray2d(), VecGetArrayPair(), VecRestoreArrayPair()
1891: @*/
1892: PetscErrorCode VecRestoreArrayInPlace(Vec x,PetscScalar **a)
1893: {
1898: #if defined(PETSC_HAVE_CUDA)
1899: if (x->petscnative && x->offloadmask == PETSC_OFFLOAD_GPU) {
1900: PetscBool is_cudatype = PETSC_FALSE;
1901: PetscObjectTypeCompareAny((PetscObject)x,&is_cudatype,VECSEQCUDA,VECMPICUDA,VECCUDA,"");
1902: if (is_cudatype) {
1903: VecCUDARestoreArray(x,a);
1904: return(0);
1905: }
1906: }
1907: #endif
1908: VecRestoreArray(x,a);
1909: return(0);
1910: }
1913: /*@C
1914: VecRestoreArrayWrite - Restores a vector after VecGetArrayWrite() has been called.
1916: Logically Collective on Vec
1918: Input Parameters:
1919: + x - the vector
1920: - a - location of pointer to array obtained from VecGetArray()
1922: Level: beginner
1924: .seealso: VecGetArray(), VecRestoreArrayRead(), VecRestoreArrays(), VecRestoreArrayF90(), VecRestoreArrayReadF90(), VecPlaceArray(), VecRestoreArray2d(),
1925: VecGetArrayPair(), VecRestoreArrayPair(), VecGetArrayWrite()
1926: @*/
1927: PetscErrorCode VecRestoreArrayWrite(Vec x,PetscScalar **a)
1928: {
1933: if (x->petscnative) {
1934: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA)
1935: x->offloadmask = PETSC_OFFLOAD_CPU;
1936: #endif
1937: } else {
1938: if (x->ops->restorearraywrite) {
1939: (*x->ops->restorearraywrite)(x,a);
1940: } else {
1941: (*x->ops->restorearray)(x,a);
1942: }
1943: }
1944: if (a) *a = NULL;
1945: PetscObjectStateIncrease((PetscObject)x);
1946: return(0);
1947: }
1949: /*@C
1950: VecRestoreArrayRead - Restore array obtained with VecGetArrayRead()
1952: Not Collective
1954: Input Parameters:
1955: + vec - the vector
1956: - array - the array
1958: Level: beginner
1960: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1961: @*/
1962: PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar **a)
1963: {
1968: if (x->petscnative) {
1969: /* nothing */
1970: } else if (x->ops->restorearrayread) {
1971: (*x->ops->restorearrayread)(x,a);
1972: } else {
1973: (*x->ops->restorearray)(x,(PetscScalar**)a);
1974: }
1975: if (a) *a = NULL;
1976: return(0);
1977: }
1979: /*@C
1980: VecRestoreArrayReadInPlace - Restore array obtained with VecGetArrayReadInPlace()
1982: Not Collective
1984: Input Parameters:
1985: + vec - the vector
1986: - array - the array
1988: Level: beginner
1990: .seealso: VecGetArrayReadInPlace(), VecRestoreArrayInPlace(), VecGetArray(), VecRestoreArray(), VecGetArrayPair(), VecRestoreArrayPair()
1991: @*/
1992: PetscErrorCode VecRestoreArrayReadInPlace(Vec x,const PetscScalar **a)
1993: {
1997: VecRestoreArrayRead(x,a);
1998: return(0);
1999: }
2001: /*@
2002: VecPlaceArray - Allows one to replace the array in a vector with an
2003: array provided by the user. This is useful to avoid copying an array
2004: into a vector.
2006: Not Collective
2008: Input Parameters:
2009: + vec - the vector
2010: - array - the array
2012: Notes:
2013: You can return to the original array with a call to VecResetArray()
2015: Level: developer
2017: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray()
2019: @*/
2020: PetscErrorCode VecPlaceArray(Vec vec,const PetscScalar array[])
2021: {
2028: if (vec->ops->placearray) {
2029: (*vec->ops->placearray)(vec,array);
2030: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot place array in this type of vector");
2031: PetscObjectStateIncrease((PetscObject)vec);
2032: return(0);
2033: }
2035: /*@C
2036: VecReplaceArray - Allows one to replace the array in a vector with an
2037: array provided by the user. This is useful to avoid copying an array
2038: into a vector.
2040: Not Collective
2042: Input Parameters:
2043: + vec - the vector
2044: - array - the array
2046: Notes:
2047: This permanently replaces the array and frees the memory associated
2048: with the old array.
2050: The memory passed in MUST be obtained with PetscMalloc() and CANNOT be
2051: freed by the user. It will be freed when the vector is destroy.
2053: Not supported from Fortran
2055: Level: developer
2057: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray()
2059: @*/
2060: PetscErrorCode VecReplaceArray(Vec vec,const PetscScalar array[])
2061: {
2067: if (vec->ops->replacearray) {
2068: (*vec->ops->replacearray)(vec,array);
2069: } else SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_SUP,"Cannot replace array in this type of vector");
2070: PetscObjectStateIncrease((PetscObject)vec);
2071: return(0);
2072: }
2075: /*@C
2076: VecCUDAGetArray - Provides access to the CUDA buffer inside a vector.
2078: This function has semantics similar to VecGetArray(): the pointer
2079: returned by this function points to a consistent view of the vector
2080: data. This may involve a copy operation of data from the host to the
2081: device if the data on the device is out of date. If the device
2082: memory hasn't been allocated previously it will be allocated as part
2083: of this function call. VecCUDAGetArray() assumes that
2084: the user will modify the vector data. This is similar to
2085: intent(inout) in fortran.
2087: The CUDA device pointer has to be released by calling
2088: VecCUDARestoreArray(). Upon restoring the vector data
2089: the data on the host will be marked as out of date. A subsequent
2090: access of the host data will thus incur a data transfer from the
2091: device to the host.
2094: Input Parameter:
2095: . v - the vector
2097: Output Parameter:
2098: . a - the CUDA device pointer
2100: Fortran note:
2101: This function is not currently available from Fortran.
2103: Level: intermediate
2105: .seealso: VecCUDARestoreArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2106: @*/
2107: PETSC_EXTERN PetscErrorCode VecCUDAGetArray(Vec v, PetscScalar **a)
2108: {
2109: #if defined(PETSC_HAVE_CUDA)
2111: #endif
2115: #if defined(PETSC_HAVE_CUDA)
2116: *a = 0;
2117: VecCUDACopyToGPU(v);
2118: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2119: #endif
2120: return(0);
2121: }
2123: /*@C
2124: VecCUDARestoreArray - Restore a CUDA device pointer previously acquired with VecCUDAGetArray().
2126: This marks the host data as out of date. Subsequent access to the
2127: vector data on the host side with for instance VecGetArray() incurs a
2128: data transfer.
2130: Input Parameter:
2131: + v - the vector
2132: - a - the CUDA device pointer. This pointer is invalid after
2133: VecCUDARestoreArray() returns.
2135: Fortran note:
2136: This function is not currently available from Fortran.
2138: Level: intermediate
2140: .seealso: VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2141: @*/
2142: PETSC_EXTERN PetscErrorCode VecCUDARestoreArray(Vec v, PetscScalar **a)
2143: {
2148: #if defined(PETSC_HAVE_CUDA)
2149: v->offloadmask = PETSC_OFFLOAD_GPU;
2150: #endif
2152: PetscObjectStateIncrease((PetscObject)v);
2153: return(0);
2154: }
2156: /*@C
2157: VecCUDAGetArrayRead - Provides read access to the CUDA buffer inside a vector.
2159: This function is analogous to VecGetArrayRead(): The pointer
2160: returned by this function points to a consistent view of the vector
2161: data. This may involve a copy operation of data from the host to the
2162: device if the data on the device is out of date. If the device
2163: memory hasn't been allocated previously it will be allocated as part
2164: of this function call. VecCUDAGetArrayRead() assumes that the
2165: user will not modify the vector data. This is analgogous to
2166: intent(in) in Fortran.
2168: The CUDA device pointer has to be released by calling
2169: VecCUDARestoreArrayRead(). If the data on the host side was
2170: previously up to date it will remain so, i.e. data on both the device
2171: and the host is up to date. Accessing data on the host side does not
2172: incur a device to host data transfer.
2174: Input Parameter:
2175: . v - the vector
2177: Output Parameter:
2178: . a - the CUDA pointer.
2180: Fortran note:
2181: This function is not currently available from Fortran.
2183: Level: intermediate
2185: .seealso: VecCUDARestoreArrayRead(), VecCUDAGetArray(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2186: @*/
2187: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayRead(Vec v, const PetscScalar **a)
2188: {
2189: #if defined(PETSC_HAVE_CUDA)
2191: #endif
2195: #if defined(PETSC_HAVE_CUDA)
2196: *a = 0;
2197: VecCUDACopyToGPU(v);
2198: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2199: #endif
2200: return(0);
2201: }
2203: /*@C
2204: VecCUDARestoreArrayRead - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayRead().
2206: If the data on the host side was previously up to date it will remain
2207: so, i.e. data on both the device and the host is up to date.
2208: Accessing data on the host side e.g. with VecGetArray() does not
2209: incur a device to host data transfer.
2211: Input Parameter:
2212: + v - the vector
2213: - a - the CUDA device pointer. This pointer is invalid after
2214: VecCUDARestoreArrayRead() returns.
2216: Fortran note:
2217: This function is not currently available from Fortran.
2219: Level: intermediate
2221: .seealso: VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2222: @*/
2223: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayRead(Vec v, const PetscScalar **a)
2224: {
2227: *a = NULL;
2228: return(0);
2229: }
2231: /*@C
2232: VecCUDAGetArrayWrite - Provides write access to the CUDA buffer inside a vector.
2234: The data pointed to by the device pointer is uninitialized. The user
2235: may not read from this data. Furthermore, the entire array needs to
2236: be filled by the user to obtain well-defined behaviour. The device
2237: memory will be allocated by this function if it hasn't been allocated
2238: previously. This is analogous to intent(out) in Fortran.
2240: The device pointer needs to be released with
2241: VecCUDARestoreArrayWrite(). When the pointer is released the
2242: host data of the vector is marked as out of data. Subsequent access
2243: of the host data with e.g. VecGetArray() incurs a device to host data
2244: transfer.
2247: Input Parameter:
2248: . v - the vector
2250: Output Parameter:
2251: . a - the CUDA pointer
2253: Fortran note:
2254: This function is not currently available from Fortran.
2256: Level: advanced
2258: .seealso: VecCUDARestoreArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecGetArrayRead()
2259: @*/
2260: PETSC_EXTERN PetscErrorCode VecCUDAGetArrayWrite(Vec v, PetscScalar **a)
2261: {
2262: #if defined(PETSC_HAVE_CUDA)
2264: #endif
2268: #if defined(PETSC_HAVE_CUDA)
2269: *a = 0;
2270: VecCUDAAllocateCheck(v);
2271: *a = ((Vec_CUDA*)v->spptr)->GPUarray;
2272: #endif
2273: return(0);
2274: }
2276: /*@C
2277: VecCUDARestoreArrayWrite - Restore a CUDA device pointer previously acquired with VecCUDAGetArrayWrite().
2279: Data on the host will be marked as out of date. Subsequent access of
2280: the data on the host side e.g. with VecGetArray() will incur a device
2281: to host data transfer.
2283: Input Parameter:
2284: + v - the vector
2285: - a - the CUDA device pointer. This pointer is invalid after
2286: VecCUDARestoreArrayWrite() returns.
2288: Fortran note:
2289: This function is not currently available from Fortran.
2291: Level: intermediate
2293: .seealso: VecCUDAGetArrayWrite(), VecCUDAGetArray(), VecCUDAGetArrayRead(), VecCUDAGetArrayWrite(), VecGetArray(), VecRestoreArray(), VecGetArrayRead()
2294: @*/
2295: PETSC_EXTERN PetscErrorCode VecCUDARestoreArrayWrite(Vec v, PetscScalar **a)
2296: {
2301: #if defined(PETSC_HAVE_CUDA)
2302: v->offloadmask = PETSC_OFFLOAD_GPU;
2303: #endif
2305: PetscObjectStateIncrease((PetscObject)v);
2306: return(0);
2307: }
2309: /*@C
2310: VecCUDAPlaceArray - Allows one to replace the GPU array in a vector with a
2311: GPU array provided by the user. This is useful to avoid copying an
2312: array into a vector.
2314: Not Collective
2316: Input Parameters:
2317: + vec - the vector
2318: - array - the GPU array
2320: Notes:
2321: You can return to the original GPU array with a call to VecCUDAResetArray()
2322: It is not possible to use VecCUDAPlaceArray() and VecPlaceArray() at the
2323: same time on the same vector.
2325: Level: developer
2327: .seealso: VecPlaceArray(), VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAReplaceArray()
2329: @*/
2330: PetscErrorCode VecCUDAPlaceArray(Vec vin,PetscScalar *a)
2331: {
2336: #if defined(PETSC_HAVE_CUDA)
2337: VecCUDACopyToGPU(vin);
2338: 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()");
2339: ((Vec_Seq*)vin->data)->unplacedarray = (PetscScalar *) ((Vec_CUDA*)vin->spptr)->GPUarray; /* save previous GPU array so reset can bring it back */
2340: ((Vec_CUDA*)vin->spptr)->GPUarray = a;
2341: vin->offloadmask = PETSC_OFFLOAD_GPU;
2342: #endif
2343: PetscObjectStateIncrease((PetscObject)vin);
2344: return(0);
2345: }
2347: /*@C
2348: VecCUDAReplaceArray - Allows one to replace the GPU array in a vector
2349: with a GPU array provided by the user. This is useful to avoid copying
2350: a GPU array into a vector.
2352: Not Collective
2354: Input Parameters:
2355: + vec - the vector
2356: - array - the GPU array
2358: Notes:
2359: This permanently replaces the GPU array and frees the memory associated
2360: with the old GPU array.
2362: The memory passed in CANNOT be freed by the user. It will be freed
2363: when the vector is destroyed.
2365: Not supported from Fortran
2367: Level: developer
2369: .seealso: VecGetArray(), VecRestoreArray(), VecPlaceArray(), VecResetArray(), VecCUDAResetArray(), VecCUDAPlaceArray(), VecReplaceArray()
2371: @*/
2372: PetscErrorCode VecCUDAReplaceArray(Vec vin,PetscScalar *a)
2373: {
2374: #if defined(PETSC_HAVE_CUDA)
2375: cudaError_t err;
2376: #endif
2381: #if defined(PETSC_HAVE_CUDA)
2382: err = cudaFree(((Vec_CUDA*)vin->spptr)->GPUarray);CHKERRCUDA(err);
2383: ((Vec_CUDA*)vin->spptr)->GPUarray = a;
2384: vin->offloadmask = PETSC_OFFLOAD_GPU;
2385: #endif
2386: PetscObjectStateIncrease((PetscObject)vin);
2387: return(0);
2388: }
2390: /*@C
2391: VecCUDAResetArray - Resets a vector to use its default memory. Call this
2392: after the use of VecCUDAPlaceArray().
2394: Not Collective
2396: Input Parameters:
2397: . vec - the vector
2399: Level: developer
2401: .seealso: VecGetArray(), VecRestoreArray(), VecReplaceArray(), VecPlaceArray(), VecResetArray(), VecCUDAPlaceArray(), VecCUDAReplaceArray()
2403: @*/
2404: PetscErrorCode VecCUDAResetArray(Vec vin)
2405: {
2410: #if defined(PETSC_HAVE_CUDA)
2411: VecCUDACopyToGPU(vin);
2412: ((Vec_CUDA*)vin->spptr)->GPUarray = (PetscScalar *) ((Vec_Seq*)vin->data)->unplacedarray;
2413: ((Vec_Seq*)vin->data)->unplacedarray = 0;
2414: vin->offloadmask = PETSC_OFFLOAD_GPU;
2415: #endif
2416: PetscObjectStateIncrease((PetscObject)vin);
2417: return(0);
2418: }
2423: /*MC
2424: VecDuplicateVecsF90 - Creates several vectors of the same type as an existing vector
2425: and makes them accessible via a Fortran90 pointer.
2427: Synopsis:
2428: VecDuplicateVecsF90(Vec x,PetscInt n,{Vec, pointer :: y(:)},integer ierr)
2430: Collective on Vec
2432: Input Parameters:
2433: + x - a vector to mimic
2434: - n - the number of vectors to obtain
2436: Output Parameters:
2437: + y - Fortran90 pointer to the array of vectors
2438: - ierr - error code
2440: Example of Usage:
2441: .vb
2442: #include <petsc/finclude/petscvec.h>
2443: use petscvec
2445: Vec x
2446: Vec, pointer :: y(:)
2447: ....
2448: call VecDuplicateVecsF90(x,2,y,ierr)
2449: call VecSet(y(2),alpha,ierr)
2450: call VecSet(y(2),alpha,ierr)
2451: ....
2452: call VecDestroyVecsF90(2,y,ierr)
2453: .ve
2455: Notes:
2456: Not yet supported for all F90 compilers
2458: Use VecDestroyVecsF90() to free the space.
2460: Level: beginner
2462: .seealso: VecDestroyVecsF90(), VecDuplicateVecs()
2464: M*/
2466: /*MC
2467: VecRestoreArrayF90 - Restores a vector to a usable state after a call to
2468: VecGetArrayF90().
2470: Synopsis:
2471: VecRestoreArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2473: Logically Collective on Vec
2475: Input Parameters:
2476: + x - vector
2477: - xx_v - the Fortran90 pointer to the array
2479: Output Parameter:
2480: . ierr - error code
2482: Example of Usage:
2483: .vb
2484: #include <petsc/finclude/petscvec.h>
2485: use petscvec
2487: PetscScalar, pointer :: xx_v(:)
2488: ....
2489: call VecGetArrayF90(x,xx_v,ierr)
2490: xx_v(3) = a
2491: call VecRestoreArrayF90(x,xx_v,ierr)
2492: .ve
2494: Level: beginner
2496: .seealso: VecGetArrayF90(), VecGetArray(), VecRestoreArray(), UsingFortran, VecRestoreArrayReadF90()
2498: M*/
2500: /*MC
2501: VecDestroyVecsF90 - Frees a block of vectors obtained with VecDuplicateVecsF90().
2503: Synopsis:
2504: VecDestroyVecsF90(PetscInt n,{Vec, pointer :: x(:)},PetscErrorCode ierr)
2506: Collective on Vec
2508: Input Parameters:
2509: + n - the number of vectors previously obtained
2510: - x - pointer to array of vector pointers
2512: Output Parameter:
2513: . ierr - error code
2515: Notes:
2516: Not yet supported for all F90 compilers
2518: Level: beginner
2520: .seealso: VecDestroyVecs(), VecDuplicateVecsF90()
2522: M*/
2524: /*MC
2525: VecGetArrayF90 - Accesses a vector array from Fortran90. For default PETSc
2526: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2527: this routine is implementation dependent. You MUST call VecRestoreArrayF90()
2528: when you no longer need access to the array.
2530: Synopsis:
2531: VecGetArrayF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2533: Logically Collective on Vec
2535: Input Parameter:
2536: . x - vector
2538: Output Parameters:
2539: + xx_v - the Fortran90 pointer to the array
2540: - ierr - error code
2542: Example of Usage:
2543: .vb
2544: #include <petsc/finclude/petscvec.h>
2545: use petscvec
2547: PetscScalar, pointer :: xx_v(:)
2548: ....
2549: call VecGetArrayF90(x,xx_v,ierr)
2550: xx_v(3) = a
2551: call VecRestoreArrayF90(x,xx_v,ierr)
2552: .ve
2554: If you ONLY intend to read entries from the array and not change any entries you should use VecGetArrayReadF90().
2556: Level: beginner
2558: .seealso: VecRestoreArrayF90(), VecGetArray(), VecRestoreArray(), VecGetArrayReadF90(), UsingFortran
2560: M*/
2562: /*MC
2563: VecGetArrayReadF90 - Accesses a read only array from Fortran90. For default PETSc
2564: vectors, VecGetArrayF90() returns a pointer to the local data array. Otherwise,
2565: this routine is implementation dependent. You MUST call VecRestoreArrayReadF90()
2566: when you no longer need access to the array.
2568: Synopsis:
2569: VecGetArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2571: Logically Collective on Vec
2573: Input Parameter:
2574: . x - vector
2576: Output Parameters:
2577: + xx_v - the Fortran90 pointer to the array
2578: - ierr - error code
2580: Example of Usage:
2581: .vb
2582: #include <petsc/finclude/petscvec.h>
2583: use petscvec
2585: PetscScalar, pointer :: xx_v(:)
2586: ....
2587: call VecGetArrayReadF90(x,xx_v,ierr)
2588: a = xx_v(3)
2589: call VecRestoreArrayReadF90(x,xx_v,ierr)
2590: .ve
2592: If you intend to write entries into the array you must use VecGetArrayF90().
2594: Level: beginner
2596: .seealso: VecRestoreArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(), VecGetArrayF90(), UsingFortran
2598: M*/
2600: /*MC
2601: VecRestoreArrayReadF90 - Restores a readonly vector to a usable state after a call to
2602: VecGetArrayReadF90().
2604: Synopsis:
2605: VecRestoreArrayReadF90(Vec x,{Scalar, pointer :: xx_v(:)},integer ierr)
2607: Logically Collective on Vec
2609: Input Parameters:
2610: + x - vector
2611: - xx_v - the Fortran90 pointer to the array
2613: Output Parameter:
2614: . ierr - error code
2616: Example of Usage:
2617: .vb
2618: #include <petsc/finclude/petscvec.h>
2619: use petscvec
2621: PetscScalar, pointer :: xx_v(:)
2622: ....
2623: call VecGetArrayReadF90(x,xx_v,ierr)
2624: a = xx_v(3)
2625: call VecRestoreArrayReadF90(x,xx_v,ierr)
2626: .ve
2628: Level: beginner
2630: .seealso: VecGetArrayReadF90(), VecGetArray(), VecRestoreArray(), VecGetArrayRead(), VecRestoreArrayRead(),UsingFortran, VecRestoreArrayF90()
2632: M*/
2634: /*@C
2635: VecGetArray2d - Returns a pointer to a 2d contiguous array that contains this
2636: processor's portion of the vector data. You MUST call VecRestoreArray2d()
2637: when you no longer need access to the array.
2639: Logically Collective
2641: Input Parameter:
2642: + x - the vector
2643: . m - first dimension of two dimensional array
2644: . n - second dimension of two dimensional array
2645: . mstart - first index you will use in first coordinate direction (often 0)
2646: - nstart - first index in the second coordinate direction (often 0)
2648: Output Parameter:
2649: . a - location to put pointer to the array
2651: Level: developer
2653: Notes:
2654: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2655: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2656: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2657: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2659: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2661: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2662: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2663: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2664: @*/
2665: PetscErrorCode VecGetArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2666: {
2668: PetscInt i,N;
2669: PetscScalar *aa;
2675: VecGetLocalSize(x,&N);
2676: 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);
2677: VecGetArray(x,&aa);
2679: PetscMalloc1(m,a);
2680: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
2681: *a -= mstart;
2682: return(0);
2683: }
2685: /*@C
2686: VecGetArray2dWrite - Returns a pointer to a 2d contiguous array that will contain this
2687: processor's portion of the vector data. You MUST call VecRestoreArray2dWrite()
2688: when you no longer need access to the array.
2690: Logically Collective
2692: Input Parameter:
2693: + x - the vector
2694: . m - first dimension of two dimensional array
2695: . n - second dimension of two dimensional array
2696: . mstart - first index you will use in first coordinate direction (often 0)
2697: - nstart - first index in the second coordinate direction (often 0)
2699: Output Parameter:
2700: . a - location to put pointer to the array
2702: Level: developer
2704: Notes:
2705: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
2706: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2707: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
2708: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
2710: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2712: Concepts: vector^accessing local values as 2d array
2714: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2715: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2716: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2717: @*/
2718: PetscErrorCode VecGetArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2719: {
2721: PetscInt i,N;
2722: PetscScalar *aa;
2728: VecGetLocalSize(x,&N);
2729: 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);
2730: VecGetArrayWrite(x,&aa);
2732: PetscMalloc1(m,a);
2733: for (i=0; i<m; i++) (*a)[i] = aa + i*n - nstart;
2734: *a -= mstart;
2735: return(0);
2736: }
2738: /*@C
2739: VecRestoreArray2d - Restores a vector after VecGetArray2d() has been called.
2741: Logically Collective
2743: Input Parameters:
2744: + x - the vector
2745: . m - first dimension of two dimensional array
2746: . n - second dimension of the two dimensional array
2747: . mstart - first index you will use in first coordinate direction (often 0)
2748: . nstart - first index in the second coordinate direction (often 0)
2749: - a - location of pointer to array obtained from VecGetArray2d()
2751: Level: developer
2753: Notes:
2754: For regular PETSc vectors this routine does not involve any copies. For
2755: any special vectors that do not store local vector data in a contiguous
2756: array, this routine will copy the data back into the underlying
2757: vector data structure from the array obtained with VecGetArray().
2759: This routine actually zeros out the a pointer.
2761: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2762: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2763: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2764: @*/
2765: PetscErrorCode VecRestoreArray2d(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2766: {
2768: void *dummy;
2774: dummy = (void*)(*a + mstart);
2775: PetscFree(dummy);
2776: VecRestoreArray(x,NULL);
2777: return(0);
2778: }
2780: /*@C
2781: VecRestoreArray2dWrite - Restores a vector after VecGetArray2dWrite() has been called.
2783: Logically Collective
2785: Input Parameters:
2786: + x - the vector
2787: . m - first dimension of two dimensional array
2788: . n - second dimension of the two dimensional array
2789: . mstart - first index you will use in first coordinate direction (often 0)
2790: . nstart - first index in the second coordinate direction (often 0)
2791: - a - location of pointer to array obtained from VecGetArray2d()
2793: Level: developer
2795: Notes:
2796: For regular PETSc vectors this routine does not involve any copies. For
2797: any special vectors that do not store local vector data in a contiguous
2798: array, this routine will copy the data back into the underlying
2799: vector data structure from the array obtained with VecGetArray().
2801: This routine actually zeros out the a pointer.
2803: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2804: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2805: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2806: @*/
2807: PetscErrorCode VecRestoreArray2dWrite(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
2808: {
2810: void *dummy;
2816: dummy = (void*)(*a + mstart);
2817: PetscFree(dummy);
2818: VecRestoreArrayWrite(x,NULL);
2819: return(0);
2820: }
2822: /*@C
2823: VecGetArray1d - Returns a pointer to a 1d contiguous array that contains this
2824: processor's portion of the vector data. You MUST call VecRestoreArray1d()
2825: when you no longer need access to the array.
2827: Logically Collective
2829: Input Parameter:
2830: + x - the vector
2831: . m - first dimension of two dimensional array
2832: - mstart - first index you will use in first coordinate direction (often 0)
2834: Output Parameter:
2835: . a - location to put pointer to the array
2837: Level: developer
2839: Notes:
2840: For a vector obtained from DMCreateLocalVector() mstart are likely
2841: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2842: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2844: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2846: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2847: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2848: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2849: @*/
2850: PetscErrorCode VecGetArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2851: {
2853: PetscInt N;
2859: VecGetLocalSize(x,&N);
2860: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2861: VecGetArray(x,a);
2862: *a -= mstart;
2863: return(0);
2864: }
2866: /*@C
2867: VecGetArray1dWrite - Returns a pointer to a 1d contiguous array that will contain this
2868: processor's portion of the vector data. You MUST call VecRestoreArray1dWrite()
2869: when you no longer need access to the array.
2871: Logically Collective
2873: Input Parameter:
2874: + x - the vector
2875: . m - first dimension of two dimensional array
2876: - mstart - first index you will use in first coordinate direction (often 0)
2878: Output Parameter:
2879: . a - location to put pointer to the array
2881: Level: developer
2883: Notes:
2884: For a vector obtained from DMCreateLocalVector() mstart are likely
2885: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
2886: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
2888: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
2890: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
2891: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
2892: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
2893: @*/
2894: PetscErrorCode VecGetArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2895: {
2897: PetscInt N;
2903: VecGetLocalSize(x,&N);
2904: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
2905: VecGetArrayWrite(x,a);
2906: *a -= mstart;
2907: return(0);
2908: }
2910: /*@C
2911: VecRestoreArray1d - Restores a vector after VecGetArray1d() has been called.
2913: Logically Collective
2915: Input Parameters:
2916: + x - the vector
2917: . m - first dimension of two dimensional array
2918: . mstart - first index you will use in first coordinate direction (often 0)
2919: - a - location of pointer to array obtained from VecGetArray21()
2921: Level: developer
2923: Notes:
2924: For regular PETSc vectors this routine does not involve any copies. For
2925: any special vectors that do not store local vector data in a contiguous
2926: array, this routine will copy the data back into the underlying
2927: vector data structure from the array obtained with VecGetArray1d().
2929: This routine actually zeros out the a pointer.
2931: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2932: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2933: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2934: @*/
2935: PetscErrorCode VecRestoreArray1d(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2936: {
2942: VecRestoreArray(x,NULL);
2943: return(0);
2944: }
2946: /*@C
2947: VecRestoreArray1dWrite - Restores a vector after VecGetArray1dWrite() has been called.
2949: Logically Collective
2951: Input Parameters:
2952: + x - the vector
2953: . m - first dimension of two dimensional array
2954: . mstart - first index you will use in first coordinate direction (often 0)
2955: - a - location of pointer to array obtained from VecGetArray21()
2957: Level: developer
2959: Notes:
2960: For regular PETSc vectors this routine does not involve any copies. For
2961: any special vectors that do not store local vector data in a contiguous
2962: array, this routine will copy the data back into the underlying
2963: vector data structure from the array obtained with VecGetArray1d().
2965: This routine actually zeros out the a pointer.
2967: Concepts: vector^accessing local values as 1d array
2969: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
2970: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
2971: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
2972: @*/
2973: PetscErrorCode VecRestoreArray1dWrite(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
2974: {
2980: VecRestoreArrayWrite(x,NULL);
2981: return(0);
2982: }
2984: /*@C
2985: VecGetArray3d - Returns a pointer to a 3d contiguous array that contains this
2986: processor's portion of the vector data. You MUST call VecRestoreArray3d()
2987: when you no longer need access to the array.
2989: Logically Collective
2991: Input Parameter:
2992: + x - the vector
2993: . m - first dimension of three dimensional array
2994: . n - second dimension of three dimensional array
2995: . p - third dimension of three dimensional array
2996: . mstart - first index you will use in first coordinate direction (often 0)
2997: . nstart - first index in the second coordinate direction (often 0)
2998: - pstart - first index in the third coordinate direction (often 0)
3000: Output Parameter:
3001: . a - location to put pointer to the array
3003: Level: developer
3005: Notes:
3006: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3007: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3008: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3009: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3011: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3013: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3014: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3015: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3016: @*/
3017: PetscErrorCode VecGetArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3018: {
3020: PetscInt i,N,j;
3021: PetscScalar *aa,**b;
3027: VecGetLocalSize(x,&N);
3028: 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);
3029: VecGetArray(x,&aa);
3031: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3032: b = (PetscScalar**)((*a) + m);
3033: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3034: for (i=0; i<m; i++)
3035: for (j=0; j<n; j++)
3036: b[i*n+j] = aa + i*n*p + j*p - pstart;
3038: *a -= mstart;
3039: return(0);
3040: }
3042: /*@C
3043: VecGetArray3dWrite - Returns a pointer to a 3d contiguous array that will contain this
3044: processor's portion of the vector data. You MUST call VecRestoreArray3dWrite()
3045: when you no longer need access to the array.
3047: Logically Collective
3049: Input Parameter:
3050: + x - the vector
3051: . m - first dimension of three dimensional array
3052: . n - second dimension of three dimensional array
3053: . p - third dimension of three dimensional array
3054: . mstart - first index you will use in first coordinate direction (often 0)
3055: . nstart - first index in the second coordinate direction (often 0)
3056: - pstart - first index in the third coordinate direction (often 0)
3058: Output Parameter:
3059: . a - location to put pointer to the array
3061: Level: developer
3063: Notes:
3064: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3065: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3066: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3067: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3069: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3071: Concepts: vector^accessing local values as 3d array
3073: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3074: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3075: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3076: @*/
3077: PetscErrorCode VecGetArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3078: {
3080: PetscInt i,N,j;
3081: PetscScalar *aa,**b;
3087: VecGetLocalSize(x,&N);
3088: 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);
3089: VecGetArrayWrite(x,&aa);
3091: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3092: b = (PetscScalar**)((*a) + m);
3093: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3094: for (i=0; i<m; i++)
3095: for (j=0; j<n; j++)
3096: b[i*n+j] = aa + i*n*p + j*p - pstart;
3098: *a -= mstart;
3099: return(0);
3100: }
3102: /*@C
3103: VecRestoreArray3d - Restores a vector after VecGetArray3d() has been called.
3105: Logically Collective
3107: Input Parameters:
3108: + x - the vector
3109: . m - first dimension of three dimensional array
3110: . n - second dimension of the three dimensional array
3111: . p - third dimension of the three dimensional array
3112: . mstart - first index you will use in first coordinate direction (often 0)
3113: . nstart - first index in the second coordinate direction (often 0)
3114: . pstart - first index in the third coordinate direction (often 0)
3115: - a - location of pointer to array obtained from VecGetArray3d()
3117: Level: developer
3119: Notes:
3120: For regular PETSc vectors this routine does not involve any copies. For
3121: any special vectors that do not store local vector data in a contiguous
3122: array, this routine will copy the data back into the underlying
3123: vector data structure from the array obtained with VecGetArray().
3125: This routine actually zeros out the a pointer.
3127: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3128: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3129: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3130: @*/
3131: PetscErrorCode VecRestoreArray3d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3132: {
3134: void *dummy;
3140: dummy = (void*)(*a + mstart);
3141: PetscFree(dummy);
3142: VecRestoreArray(x,NULL);
3143: return(0);
3144: }
3146: /*@C
3147: VecRestoreArray3dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3149: Logically Collective
3151: Input Parameters:
3152: + x - the vector
3153: . m - first dimension of three dimensional array
3154: . n - second dimension of the three dimensional array
3155: . p - third dimension of the three dimensional array
3156: . mstart - first index you will use in first coordinate direction (often 0)
3157: . nstart - first index in the second coordinate direction (often 0)
3158: . pstart - first index in the third coordinate direction (often 0)
3159: - a - location of pointer to array obtained from VecGetArray3d()
3161: Level: developer
3163: Notes:
3164: For regular PETSc vectors this routine does not involve any copies. For
3165: any special vectors that do not store local vector data in a contiguous
3166: array, this routine will copy the data back into the underlying
3167: vector data structure from the array obtained with VecGetArray().
3169: This routine actually zeros out the a pointer.
3171: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3172: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3173: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3174: @*/
3175: PetscErrorCode VecRestoreArray3dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3176: {
3178: void *dummy;
3184: dummy = (void*)(*a + mstart);
3185: PetscFree(dummy);
3186: VecRestoreArrayWrite(x,NULL);
3187: return(0);
3188: }
3190: /*@C
3191: VecGetArray4d - Returns a pointer to a 4d contiguous array that contains this
3192: processor's portion of the vector data. You MUST call VecRestoreArray4d()
3193: when you no longer need access to the array.
3195: Logically Collective
3197: Input Parameter:
3198: + x - the vector
3199: . m - first dimension of four dimensional array
3200: . n - second dimension of four dimensional array
3201: . p - third dimension of four dimensional array
3202: . q - fourth dimension of four dimensional array
3203: . mstart - first index you will use in first coordinate direction (often 0)
3204: . nstart - first index in the second coordinate direction (often 0)
3205: . pstart - first index in the third coordinate direction (often 0)
3206: - qstart - first index in the fourth coordinate direction (often 0)
3208: Output Parameter:
3209: . a - location to put pointer to the array
3211: Level: beginner
3213: Notes:
3214: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3215: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3216: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3217: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3219: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3221: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3222: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3223: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3224: @*/
3225: PetscErrorCode VecGetArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3226: {
3228: PetscInt i,N,j,k;
3229: PetscScalar *aa,***b,**c;
3235: VecGetLocalSize(x,&N);
3236: 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);
3237: VecGetArray(x,&aa);
3239: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3240: b = (PetscScalar***)((*a) + m);
3241: c = (PetscScalar**)(b + m*n);
3242: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3243: for (i=0; i<m; i++)
3244: for (j=0; j<n; j++)
3245: b[i*n+j] = c + i*n*p + j*p - pstart;
3246: for (i=0; i<m; i++)
3247: for (j=0; j<n; j++)
3248: for (k=0; k<p; k++)
3249: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3250: *a -= mstart;
3251: return(0);
3252: }
3254: /*@C
3255: VecGetArray4dWrite - Returns a pointer to a 4d contiguous array that will contain this
3256: processor's portion of the vector data. You MUST call VecRestoreArray4dWrite()
3257: when you no longer need access to the array.
3259: Logically Collective
3261: Input Parameter:
3262: + x - the vector
3263: . m - first dimension of four dimensional array
3264: . n - second dimension of four dimensional array
3265: . p - third dimension of four dimensional array
3266: . q - fourth dimension of four dimensional array
3267: . mstart - first index you will use in first coordinate direction (often 0)
3268: . nstart - first index in the second coordinate direction (often 0)
3269: . pstart - first index in the third coordinate direction (often 0)
3270: - qstart - first index in the fourth coordinate direction (often 0)
3272: Output Parameter:
3273: . a - location to put pointer to the array
3275: Level: beginner
3277: Notes:
3278: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3279: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3280: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3281: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3283: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3285: Concepts: vector^accessing local values as 3d array
3287: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3288: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3289: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3290: @*/
3291: PetscErrorCode VecGetArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3292: {
3294: PetscInt i,N,j,k;
3295: PetscScalar *aa,***b,**c;
3301: VecGetLocalSize(x,&N);
3302: 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);
3303: VecGetArrayWrite(x,&aa);
3305: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3306: b = (PetscScalar***)((*a) + m);
3307: c = (PetscScalar**)(b + m*n);
3308: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3309: for (i=0; i<m; i++)
3310: for (j=0; j<n; j++)
3311: b[i*n+j] = c + i*n*p + j*p - pstart;
3312: for (i=0; i<m; i++)
3313: for (j=0; j<n; j++)
3314: for (k=0; k<p; k++)
3315: c[i*n*p+j*p+k] = aa + i*n*p*q + j*p*q + k*q - qstart;
3316: *a -= mstart;
3317: return(0);
3318: }
3320: /*@C
3321: VecRestoreArray4d - Restores a vector after VecGetArray3d() has been called.
3323: Logically Collective
3325: Input Parameters:
3326: + x - the vector
3327: . m - first dimension of four dimensional array
3328: . n - second dimension of the four dimensional array
3329: . p - third dimension of the four dimensional array
3330: . q - fourth dimension of the four dimensional array
3331: . mstart - first index you will use in first coordinate direction (often 0)
3332: . nstart - first index in the second coordinate direction (often 0)
3333: . pstart - first index in the third coordinate direction (often 0)
3334: . qstart - first index in the fourth coordinate direction (often 0)
3335: - a - location of pointer to array obtained from VecGetArray4d()
3337: Level: beginner
3339: Notes:
3340: For regular PETSc vectors this routine does not involve any copies. For
3341: any special vectors that do not store local vector data in a contiguous
3342: array, this routine will copy the data back into the underlying
3343: vector data structure from the array obtained with VecGetArray().
3345: This routine actually zeros out the a pointer.
3347: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3348: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3349: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3350: @*/
3351: PetscErrorCode VecRestoreArray4d(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3352: {
3354: void *dummy;
3360: dummy = (void*)(*a + mstart);
3361: PetscFree(dummy);
3362: VecRestoreArray(x,NULL);
3363: return(0);
3364: }
3366: /*@C
3367: VecRestoreArray4dWrite - Restores a vector after VecGetArray3dWrite() has been called.
3369: Logically Collective
3371: Input Parameters:
3372: + x - the vector
3373: . m - first dimension of four dimensional array
3374: . n - second dimension of the four dimensional array
3375: . p - third dimension of the four dimensional array
3376: . q - fourth dimension of the four dimensional array
3377: . mstart - first index you will use in first coordinate direction (often 0)
3378: . nstart - first index in the second coordinate direction (often 0)
3379: . pstart - first index in the third coordinate direction (often 0)
3380: . qstart - first index in the fourth coordinate direction (often 0)
3381: - a - location of pointer to array obtained from VecGetArray4d()
3383: Level: beginner
3385: Notes:
3386: For regular PETSc vectors this routine does not involve any copies. For
3387: any special vectors that do not store local vector data in a contiguous
3388: array, this routine will copy the data back into the underlying
3389: vector data structure from the array obtained with VecGetArray().
3391: This routine actually zeros out the a pointer.
3393: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3394: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3395: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3396: @*/
3397: PetscErrorCode VecRestoreArray4dWrite(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3398: {
3400: void *dummy;
3406: dummy = (void*)(*a + mstart);
3407: PetscFree(dummy);
3408: VecRestoreArrayWrite(x,NULL);
3409: return(0);
3410: }
3412: /*@C
3413: VecGetArray2dRead - Returns a pointer to a 2d contiguous array that contains this
3414: processor's portion of the vector data. You MUST call VecRestoreArray2dRead()
3415: when you no longer need access to the array.
3417: Logically Collective
3419: Input Parameter:
3420: + x - the vector
3421: . m - first dimension of two dimensional array
3422: . n - second dimension of two dimensional array
3423: . mstart - first index you will use in first coordinate direction (often 0)
3424: - nstart - first index in the second coordinate direction (often 0)
3426: Output Parameter:
3427: . a - location to put pointer to the array
3429: Level: developer
3431: Notes:
3432: For a vector obtained from DMCreateLocalVector() mstart and nstart are likely
3433: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3434: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3435: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray2d().
3437: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3439: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3440: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3441: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3442: @*/
3443: PetscErrorCode VecGetArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3444: {
3445: PetscErrorCode ierr;
3446: PetscInt i,N;
3447: const PetscScalar *aa;
3453: VecGetLocalSize(x,&N);
3454: 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);
3455: VecGetArrayRead(x,&aa);
3457: PetscMalloc1(m,a);
3458: for (i=0; i<m; i++) (*a)[i] = (PetscScalar*) aa + i*n - nstart;
3459: *a -= mstart;
3460: return(0);
3461: }
3463: /*@C
3464: VecRestoreArray2dRead - Restores a vector after VecGetArray2dRead() has been called.
3466: Logically Collective
3468: Input Parameters:
3469: + x - the vector
3470: . m - first dimension of two dimensional array
3471: . n - second dimension of the two dimensional array
3472: . mstart - first index you will use in first coordinate direction (often 0)
3473: . nstart - first index in the second coordinate direction (often 0)
3474: - a - location of pointer to array obtained from VecGetArray2d()
3476: Level: developer
3478: Notes:
3479: For regular PETSc vectors this routine does not involve any copies. For
3480: any special vectors that do not store local vector data in a contiguous
3481: array, this routine will copy the data back into the underlying
3482: vector data structure from the array obtained with VecGetArray().
3484: This routine actually zeros out the a pointer.
3486: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3487: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3488: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3489: @*/
3490: PetscErrorCode VecRestoreArray2dRead(Vec x,PetscInt m,PetscInt n,PetscInt mstart,PetscInt nstart,PetscScalar **a[])
3491: {
3493: void *dummy;
3499: dummy = (void*)(*a + mstart);
3500: PetscFree(dummy);
3501: VecRestoreArrayRead(x,NULL);
3502: return(0);
3503: }
3505: /*@C
3506: VecGetArray1dRead - Returns a pointer to a 1d contiguous array that contains this
3507: processor's portion of the vector data. You MUST call VecRestoreArray1dRead()
3508: when you no longer need access to the array.
3510: Logically Collective
3512: Input Parameter:
3513: + x - the vector
3514: . m - first dimension of two dimensional array
3515: - mstart - first index you will use in first coordinate direction (often 0)
3517: Output Parameter:
3518: . a - location to put pointer to the array
3520: Level: developer
3522: Notes:
3523: For a vector obtained from DMCreateLocalVector() mstart are likely
3524: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3525: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners().
3527: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3529: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3530: VecRestoreArray2d(), DMDAVecGetArray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3531: VecGetArray2d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3532: @*/
3533: PetscErrorCode VecGetArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3534: {
3536: PetscInt N;
3542: VecGetLocalSize(x,&N);
3543: if (m != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Local array size %D does not match 1d array dimensions %D",N,m);
3544: VecGetArrayRead(x,(const PetscScalar**)a);
3545: *a -= mstart;
3546: return(0);
3547: }
3549: /*@C
3550: VecRestoreArray1dRead - Restores a vector after VecGetArray1dRead() has been called.
3552: Logically Collective
3554: Input Parameters:
3555: + x - the vector
3556: . m - first dimension of two dimensional array
3557: . mstart - first index you will use in first coordinate direction (often 0)
3558: - a - location of pointer to array obtained from VecGetArray21()
3560: Level: developer
3562: Notes:
3563: For regular PETSc vectors this routine does not involve any copies. For
3564: any special vectors that do not store local vector data in a contiguous
3565: array, this routine will copy the data back into the underlying
3566: vector data structure from the array obtained with VecGetArray1dRead().
3568: This routine actually zeros out the a pointer.
3570: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3571: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3572: VecGetArray1d(), VecRestoreArray2d(), VecGetArray4d(), VecRestoreArray4d()
3573: @*/
3574: PetscErrorCode VecRestoreArray1dRead(Vec x,PetscInt m,PetscInt mstart,PetscScalar *a[])
3575: {
3581: VecRestoreArrayRead(x,NULL);
3582: return(0);
3583: }
3586: /*@C
3587: VecGetArray3dRead - Returns a pointer to a 3d contiguous array that contains this
3588: processor's portion of the vector data. You MUST call VecRestoreArray3dRead()
3589: when you no longer need access to the array.
3591: Logically Collective
3593: Input Parameter:
3594: + x - the vector
3595: . m - first dimension of three dimensional array
3596: . n - second dimension of three dimensional array
3597: . p - third dimension of three dimensional array
3598: . mstart - first index you will use in first coordinate direction (often 0)
3599: . nstart - first index in the second coordinate direction (often 0)
3600: - pstart - first index in the third coordinate direction (often 0)
3602: Output Parameter:
3603: . a - location to put pointer to the array
3605: Level: developer
3607: Notes:
3608: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3609: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3610: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3611: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3dRead().
3613: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3615: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3616: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3617: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3618: @*/
3619: PetscErrorCode VecGetArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3620: {
3621: PetscErrorCode ierr;
3622: PetscInt i,N,j;
3623: const PetscScalar *aa;
3624: PetscScalar **b;
3630: VecGetLocalSize(x,&N);
3631: 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);
3632: VecGetArrayRead(x,&aa);
3634: PetscMalloc1(m*sizeof(PetscScalar**)+m*n,a);
3635: b = (PetscScalar**)((*a) + m);
3636: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3637: for (i=0; i<m; i++)
3638: for (j=0; j<n; j++)
3639: b[i*n+j] = (PetscScalar *)aa + i*n*p + j*p - pstart;
3641: *a -= mstart;
3642: return(0);
3643: }
3645: /*@C
3646: VecRestoreArray3dRead - Restores a vector after VecGetArray3dRead() has been called.
3648: Logically Collective
3650: Input Parameters:
3651: + x - the vector
3652: . m - first dimension of three dimensional array
3653: . n - second dimension of the three dimensional array
3654: . p - third dimension of the three dimensional array
3655: . mstart - first index you will use in first coordinate direction (often 0)
3656: . nstart - first index in the second coordinate direction (often 0)
3657: . pstart - first index in the third coordinate direction (often 0)
3658: - a - location of pointer to array obtained from VecGetArray3dRead()
3660: Level: developer
3662: Notes:
3663: For regular PETSc vectors this routine does not involve any copies. For
3664: any special vectors that do not store local vector data in a contiguous
3665: array, this routine will copy the data back into the underlying
3666: vector data structure from the array obtained with VecGetArray().
3668: This routine actually zeros out the a pointer.
3670: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3671: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3672: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3673: @*/
3674: PetscErrorCode VecRestoreArray3dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscScalar ***a[])
3675: {
3677: void *dummy;
3683: dummy = (void*)(*a + mstart);
3684: PetscFree(dummy);
3685: VecRestoreArrayRead(x,NULL);
3686: return(0);
3687: }
3689: /*@C
3690: VecGetArray4dRead - Returns a pointer to a 4d contiguous array that contains this
3691: processor's portion of the vector data. You MUST call VecRestoreArray4dRead()
3692: when you no longer need access to the array.
3694: Logically Collective
3696: Input Parameter:
3697: + x - the vector
3698: . m - first dimension of four dimensional array
3699: . n - second dimension of four dimensional array
3700: . p - third dimension of four dimensional array
3701: . q - fourth dimension of four dimensional array
3702: . mstart - first index you will use in first coordinate direction (often 0)
3703: . nstart - first index in the second coordinate direction (often 0)
3704: . pstart - first index in the third coordinate direction (often 0)
3705: - qstart - first index in the fourth coordinate direction (often 0)
3707: Output Parameter:
3708: . a - location to put pointer to the array
3710: Level: beginner
3712: Notes:
3713: For a vector obtained from DMCreateLocalVector() mstart, nstart, and pstart are likely
3714: obtained from the corner indices obtained from DMDAGetGhostCorners() while for
3715: DMCreateGlobalVector() they are the corner indices from DMDAGetCorners(). In both cases
3716: the arguments from DMDAGet[Ghost]Corners() are reversed in the call to VecGetArray3d().
3718: For standard PETSc vectors this is an inexpensive call; it does not copy the vector values.
3720: .seealso: VecGetArray(), VecRestoreArray(), VecGetArrays(), VecGetArrayF90(), VecPlaceArray(),
3721: VecRestoreArray2d(), DMDAVecGetarray(), DMDAVecRestoreArray(), VecGetArray3d(), VecRestoreArray3d(),
3722: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d()
3723: @*/
3724: PetscErrorCode VecGetArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3725: {
3726: PetscErrorCode ierr;
3727: PetscInt i,N,j,k;
3728: const PetscScalar *aa;
3729: PetscScalar ***b,**c;
3735: VecGetLocalSize(x,&N);
3736: 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);
3737: VecGetArrayRead(x,&aa);
3739: PetscMalloc1(m*sizeof(PetscScalar***)+m*n*sizeof(PetscScalar**)+m*n*p,a);
3740: b = (PetscScalar***)((*a) + m);
3741: c = (PetscScalar**)(b + m*n);
3742: for (i=0; i<m; i++) (*a)[i] = b + i*n - nstart;
3743: for (i=0; i<m; i++)
3744: for (j=0; j<n; j++)
3745: b[i*n+j] = c + i*n*p + j*p - pstart;
3746: for (i=0; i<m; i++)
3747: for (j=0; j<n; j++)
3748: for (k=0; k<p; k++)
3749: c[i*n*p+j*p+k] = (PetscScalar*) aa + i*n*p*q + j*p*q + k*q - qstart;
3750: *a -= mstart;
3751: return(0);
3752: }
3754: /*@C
3755: VecRestoreArray4dRead - Restores a vector after VecGetArray3d() has been called.
3757: Logically Collective
3759: Input Parameters:
3760: + x - the vector
3761: . m - first dimension of four dimensional array
3762: . n - second dimension of the four dimensional array
3763: . p - third dimension of the four dimensional array
3764: . q - fourth dimension of the four dimensional array
3765: . mstart - first index you will use in first coordinate direction (often 0)
3766: . nstart - first index in the second coordinate direction (often 0)
3767: . pstart - first index in the third coordinate direction (often 0)
3768: . qstart - first index in the fourth coordinate direction (often 0)
3769: - a - location of pointer to array obtained from VecGetArray4dRead()
3771: Level: beginner
3773: Notes:
3774: For regular PETSc vectors this routine does not involve any copies. For
3775: any special vectors that do not store local vector data in a contiguous
3776: array, this routine will copy the data back into the underlying
3777: vector data structure from the array obtained with VecGetArray().
3779: This routine actually zeros out the a pointer.
3781: .seealso: VecGetArray(), VecRestoreArray(), VecRestoreArrays(), VecRestoreArrayF90(), VecPlaceArray(),
3782: VecGetArray2d(), VecGetArray3d(), VecRestoreArray3d(), DMDAVecGetArray(), DMDAVecRestoreArray()
3783: VecGetArray1d(), VecRestoreArray1d(), VecGetArray4d(), VecRestoreArray4d(), VecGet
3784: @*/
3785: PetscErrorCode VecRestoreArray4dRead(Vec x,PetscInt m,PetscInt n,PetscInt p,PetscInt q,PetscInt mstart,PetscInt nstart,PetscInt pstart,PetscInt qstart,PetscScalar ****a[])
3786: {
3788: void *dummy;
3794: dummy = (void*)(*a + mstart);
3795: PetscFree(dummy);
3796: VecRestoreArrayRead(x,NULL);
3797: return(0);
3798: }
3800: #if defined(PETSC_USE_DEBUG)
3802: /*@
3803: VecLockGet - Gets the current lock status of a vector
3805: Logically Collective on Vec
3807: Input Parameter:
3808: . x - the vector
3810: Output Parameter:
3811: . state - greater than zero indicates the vector is locked for read; less then zero indicates the vector is
3812: locked for write; equal to zero means the vector is unlocked, that is, it is free to read or write.
3814: Level: beginner
3816: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop()
3817: @*/
3818: PetscErrorCode VecLockGet(Vec x,PetscInt *state)
3819: {
3822: *state = x->lock;
3823: return(0);
3824: }
3826: /*@
3827: VecLockReadPush - Pushes a read-only lock on a vector to prevent it from writing
3829: Logically Collective on Vec
3831: Input Parameter:
3832: . x - the vector
3834: Notes:
3835: If this is set then calls to VecGetArray() or VecSetValues() or any other routines that change the vectors values will fail.
3837: The call can be nested, i.e., called multiple times on the same vector, but each VecLockReadPush(x) has to have one matching
3838: VecLockReadPop(x), which removes the latest read-only lock.
3840: Level: beginner
3842: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPop(), VecLockGet()
3843: @*/
3844: PetscErrorCode VecLockReadPush(Vec x)
3845: {
3848: 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");
3849: x->lock++;
3850: return(0);
3851: }
3853: /*@
3854: VecLockReadPop - Pops a read-only lock from a vector
3856: Logically Collective on Vec
3858: Input Parameter:
3859: . x - the vector
3861: Level: beginner
3863: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockGet()
3864: @*/
3865: PetscErrorCode VecLockReadPop(Vec x)
3866: {
3869: x->lock--;
3870: if (x->lock < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Vector has been unlocked from read-only access too many times");
3871: return(0);
3872: }
3874: /*@C
3875: VecLockWriteSet_Private - Lock or unlock a vector for exclusive read/write access
3877: Logically Collective on Vec
3879: Input Parameter:
3880: + x - the vector
3881: - flg - PETSC_TRUE to lock the vector for writing; PETSC_FALSE to unlock it.
3883: Notes:
3884: The function is usefull in split-phase computations, which usually have a begin phase and an end phase.
3885: One can call VecLockWriteSet_Private(x,PETSC_TRUE) in the begin phase to lock a vector for exclusive
3886: access, and call VecLockWriteSet_Private(x,PETSC_FALSE) in the end phase to unlock the vector from exclusive
3887: access. In this way, one is ensured no other operations can access the vector in between. The code may like
3890: VecGetArray(x,&xdata); // begin phase
3891: VecLockWriteSet_Private(v,PETSC_TRUE);
3893: Other operations, which can not acceess x anymore (they can access xdata, of course)
3895: VecRestoreArray(x,&vdata); // end phase
3896: VecLockWriteSet_Private(v,PETSC_FALSE);
3898: The call can not be nested on the same vector, in other words, one can not call VecLockWriteSet_Private(x,PETSC_TRUE)
3899: again before calling VecLockWriteSet_Private(v,PETSC_FALSE).
3901: Level: beginner
3903: .seealso: VecRestoreArray(), VecGetArrayRead(), VecLockReadPush(), VecLockReadPop(), VecLockGet()
3904: @*/
3905: PetscErrorCode VecLockWriteSet_Private(Vec x,PetscBool flg)
3906: {
3909: if (flg) {
3910: 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");
3911: 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");
3912: else x->lock = -1;
3913: } else {
3914: 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");
3915: x->lock = 0;
3916: }
3917: return(0);
3918: }
3920: /*@
3921: VecLockPush - Pushes a read-only lock on a vector to prevent it from writing
3923: Level: deprecated
3925: .seealso: VecLockReadPush()
3926: @*/
3927: PetscErrorCode VecLockPush(Vec x)
3928: {
3931: VecLockReadPush(x);
3932: return(0);
3933: }
3935: /*@
3936: VecLockPop - Pops a read-only lock from a vector
3938: Level: deprecated
3940: .seealso: VecLockReadPop()
3941: @*/
3942: PetscErrorCode VecLockPop(Vec x)
3943: {
3946: VecLockReadPop(x);
3947: return(0);
3948: }
3950: #endif