Actual source code: petsc-vecimpl.h
petsc-3.3-p7 2013-05-11
2: /*
3: This private file should not be included in users' code.
4: Defines the fields shared by all vector implementations.
6: */
8: #ifndef __VECIMPL_H
11: #include <petscvec.h>
13: /*S PetscThreadsLayout - defines layout of vectors and matrices across threads (which rows are assigned to which threads)
15: Level: developer
17: S*/
18: /* PetscThreadsLayout will be removed eventually with the removal of the old pthread interface. */
19: typedef struct _n_PetscThreadsLayout* PetscThreadsLayout;
20: struct _n_PetscThreadsLayout{
21: PetscInt nthreads; /* Number of threads used for vector/matrix operations */
22: PetscInt N; /* global size = sum(local sizes of all threads) */
23: PetscInt *trstarts; /* local start for each thread */
24: PetscInt *affinity; /* Core affinity of each thread */
25: };
27: /*S
28: PetscLayout - defines layout of vectors and matrices across processes (which rows are owned by which processes)
30: Level: developer
33: .seealso: PetscLayoutCreate(), PetscLayoutDestroy()
34: S*/
35: typedef struct _n_PetscLayout* PetscLayout;
36: struct _n_PetscLayout{
37: MPI_Comm comm;
38: PetscInt n,N; /* local, global vector size */
39: PetscInt rstart,rend; /* local start, local end + 1 */
40: PetscInt *range; /* the offset of each processor */
41: PetscInt bs; /* number of elements in each block (generally for multi-component problems) Do NOT multiply above numbers by bs */
42: PetscInt refcnt; /* MPI Vecs obtained with VecDuplicate() and from MatGetVecs() reuse map of input object */
43: ISLocalToGlobalMapping mapping; /* mapping used in Vec/MatSetValuesLocal() */
44: ISLocalToGlobalMapping bmapping; /* mapping used in Vec/MatSetValuesBlockedLocal() */
45: PetscThreadsLayout tmap; /* threads specific layout info */
46: PetscInt *trstarts; /* local start for each thread */
47: };
49: PETSC_EXTERN PetscErrorCode PetscLayoutCreate(MPI_Comm,PetscLayout*);
50: PETSC_EXTERN PetscErrorCode PetscLayoutSetUp(PetscLayout);
51: PETSC_EXTERN PetscErrorCode PetscLayoutDestroy(PetscLayout*);
52: PETSC_EXTERN PetscErrorCode PetscLayoutDuplicate(PetscLayout,PetscLayout*);
53: PETSC_EXTERN PetscErrorCode PetscLayoutReference(PetscLayout,PetscLayout*);
54: PETSC_EXTERN PetscErrorCode PetscLayoutSetLocalSize(PetscLayout,PetscInt);
55: PETSC_EXTERN PetscErrorCode PetscLayoutGetLocalSize(PetscLayout,PetscInt *);
56: PETSC_EXTERN PetscErrorCode PetscLayoutSetSize(PetscLayout,PetscInt);
57: PETSC_EXTERN PetscErrorCode PetscLayoutGetSize(PetscLayout,PetscInt *);
58: PETSC_EXTERN PetscErrorCode PetscLayoutSetBlockSize(PetscLayout,PetscInt);
59: PETSC_EXTERN PetscErrorCode PetscLayoutGetBlockSize(PetscLayout,PetscInt*);
60: PETSC_EXTERN PetscErrorCode PetscLayoutGetRange(PetscLayout,PetscInt *,PetscInt *);
61: PETSC_EXTERN PetscErrorCode PetscLayoutGetRanges(PetscLayout,const PetscInt *[]);
62: PETSC_EXTERN PetscErrorCode PetscLayoutSetISLocalToGlobalMapping(PetscLayout,ISLocalToGlobalMapping);
63: PETSC_EXTERN PetscErrorCode PetscLayoutSetISLocalToGlobalMappingBlock(PetscLayout,ISLocalToGlobalMapping);
65: PETSC_EXTERN PetscErrorCode PetscThreadsLayoutCreate(PetscThreadsLayout*);
66: PETSC_EXTERN PetscErrorCode PetscThreadsLayoutDestroy(PetscThreadsLayout*);
67: PETSC_EXTERN PetscErrorCode PetscThreadsLayoutSetNThreads(PetscThreadsLayout,PetscInt);
68: PETSC_EXTERN PetscErrorCode PetscThreadsLayoutSetThreadAffinities(PetscThreadsLayout, PetscInt[]);
69: PETSC_EXTERN PetscErrorCode PetscThreadsLayoutSetUp(PetscThreadsLayout);
73: /*@C
74: PetscLayoutFindOwner - Find the owning rank for a global index
76: Not Collective
78: Input Parameters:
79: + map - the layout
80: - idx - global index to find the owner of
82: Output Parameter:
83: . owner - the owning rank
85: Level: developer
87: Fortran Notes:
88: Not available from Fortran
90: @*/
91: PETSC_STATIC_INLINE PetscErrorCode PetscLayoutFindOwner(PetscLayout map,PetscInt idx,PetscInt *owner)
92: {
94: PetscMPIInt lo = 0,hi,t;
97: if (!((map->n >= 0) && (map->N >= 0) && (map->range))) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"PetscLayoutSetUp() must be called first");
98: if (idx < 0 || idx > map->N) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Index %D is out of range",idx);
99: MPI_Comm_size(map->comm,&hi);
100: while (hi - lo > 1) {
101: t = lo + (hi - lo) / 2;
102: if (idx < map->range[t]) hi = t;
103: else lo = t;
104: }
105: *owner = lo;
106: return(0);
107: }
109: /* ----------------------------------------------------------------------------*/
110: typedef struct _n_PetscUniformSection *PetscUniformSection;
111: struct _n_PetscUniformSection {
112: MPI_Comm comm;
113: PetscInt pStart, pEnd; /* The chart: all points are contained in [pStart, pEnd) */
114: PetscInt numDof; /* Describes layout of storage, point --> (constant # of values, (p - pStart)*constant # of values) */
115: };
117: #if 0
118: /* Should I protect these for C++? */
119: PETSC_EXTERN PetscErrorCode PetscSectionGetDof(PetscUniformSection, PetscInt, PetscInt*);
120: PETSC_EXTERN PetscErrorCode PetscSectionGetOffset(PetscUniformSection, PetscInt, PetscInt*);
121: #endif
123: /*S
124: PetscSection - This is a mapping from DMMESH points to sets of values, which is
125: our presentation of a fibre bundle.
127: Level: developer
129: .seealso: PetscSectionCreate(), PetscSectionDestroy()
130: S*/
131: typedef struct _n_PetscSection *PetscSection;
132: struct _n_PetscSection {
133: struct _n_PetscUniformSection atlasLayout; /* Layout for the atlas */
134: PetscInt *atlasDof; /* Describes layout of storage, point --> # of values */
135: PetscInt *atlasOff; /* Describes layout of storage, point --> offset into storage */
136: PetscSection bc; /* Describes constraints, point --> # local dofs which are constrained */
137: PetscInt *bcIndices; /* Local indices for constrained dofs */
138: PetscInt refcnt; /* Vecs obtained with VecDuplicate() and from MatGetVecs() reuse map of input object */
139: PetscBool setup;
141: PetscInt numFields; /* The number of fields making up the degrees of freedom */
142: const char **fieldNames; /* The field names */
143: PetscInt *numFieldComponents; /* The number of components in each field */
144: PetscSection *field; /* A section describing the layout and constraints for each field */
145: };
147: PETSC_EXTERN PetscErrorCode PetscSectionCreate(MPI_Comm,PetscSection*);
148: PETSC_EXTERN PetscErrorCode PetscSectionGetNumFields(PetscSection, PetscInt *);
149: PETSC_EXTERN PetscErrorCode PetscSectionSetNumFields(PetscSection, PetscInt);
150: PETSC_EXTERN PetscErrorCode PetscSectionGetFieldName(PetscSection, PetscInt, const char *[]);
151: PETSC_EXTERN PetscErrorCode PetscSectionSetFieldName(PetscSection, PetscInt, const char []);
152: PETSC_EXTERN PetscErrorCode PetscSectionGetFieldComponents(PetscSection, PetscInt, PetscInt *);
153: PETSC_EXTERN PetscErrorCode PetscSectionSetFieldComponents(PetscSection, PetscInt, PetscInt);
154: PETSC_EXTERN PetscErrorCode PetscSectionGetChart(PetscSection, PetscInt *, PetscInt *);
155: PETSC_EXTERN PetscErrorCode PetscSectionSetChart(PetscSection, PetscInt, PetscInt);
156: PETSC_EXTERN PetscErrorCode PetscSectionGetDof(PetscSection, PetscInt, PetscInt*);
157: PETSC_EXTERN PetscErrorCode PetscSectionSetDof(PetscSection, PetscInt, PetscInt);
158: PETSC_EXTERN PetscErrorCode PetscSectionAddDof(PetscSection, PetscInt, PetscInt);
159: PETSC_EXTERN PetscErrorCode PetscSectionGetFieldDof(PetscSection, PetscInt, PetscInt, PetscInt*);
160: PETSC_EXTERN PetscErrorCode PetscSectionSetFieldDof(PetscSection, PetscInt, PetscInt, PetscInt);
161: PETSC_EXTERN PetscErrorCode PetscSectionGetConstraintDof(PetscSection, PetscInt, PetscInt*);
162: PETSC_EXTERN PetscErrorCode PetscSectionSetConstraintDof(PetscSection, PetscInt, PetscInt);
163: PETSC_EXTERN PetscErrorCode PetscSectionAddConstraintDof(PetscSection, PetscInt, PetscInt);
164: PETSC_EXTERN PetscErrorCode PetscSectionGetFieldConstraintDof(PetscSection, PetscInt, PetscInt, PetscInt*);
165: PETSC_EXTERN PetscErrorCode PetscSectionSetFieldConstraintDof(PetscSection, PetscInt, PetscInt, PetscInt);
166: PETSC_EXTERN PetscErrorCode PetscSectionGetConstraintIndices(PetscSection, PetscInt, PetscInt**);
167: PETSC_EXTERN PetscErrorCode PetscSectionSetConstraintIndices(PetscSection, PetscInt, PetscInt*);
168: PETSC_EXTERN PetscErrorCode PetscSectionGetFieldConstraintIndices(PetscSection, PetscInt, PetscInt, PetscInt**);
169: PETSC_EXTERN PetscErrorCode PetscSectionSetFieldConstraintIndices(PetscSection, PetscInt, PetscInt, PetscInt*);
170: PETSC_EXTERN PetscErrorCode PetscSectionSetUp(PetscSection);
171: PETSC_EXTERN PetscErrorCode PetscSectionGetStorageSize(PetscSection, PetscInt*);
172: PETSC_EXTERN PetscErrorCode PetscSectionGetConstrainedStorageSize(PetscSection, PetscInt*);
173: PETSC_EXTERN PetscErrorCode PetscSectionGetOffset(PetscSection, PetscInt, PetscInt*);
174: PETSC_EXTERN PetscErrorCode PetscSectionGetFieldOffset(PetscSection, PetscInt, PetscInt, PetscInt*);
175: PETSC_EXTERN PetscErrorCode PetscSectionGetOffsetRange(PetscSection, PetscInt *, PetscInt *);
176: PETSC_EXTERN PetscErrorCode PetscSectionView(PetscSection, PetscViewer);
177: PETSC_EXTERN PetscErrorCode PetscSectionVecView(PetscSection, Vec, PetscViewer);
178: PETSC_EXTERN PetscErrorCode PetscSectionDestroy(PetscSection*);
179: PETSC_EXTERN PetscErrorCode PetscSectionCreateGlobalSection(PetscSection, PetscSF, PetscSection *);
180: PETSC_EXTERN PetscErrorCode PetscSectionGetPointLayout(MPI_Comm, PetscSection, PetscLayout *);
181: PETSC_EXTERN PetscErrorCode PetscSectionGetValueLayout(MPI_Comm, PetscSection, PetscLayout *);
183: /* Sieve support */
184: PETSC_EXTERN PetscErrorCode PetscSFConvertPartition(PetscSF, PetscSection, IS, ISLocalToGlobalMapping *, PetscSF *);
185: PETSC_EXTERN PetscErrorCode PetscSFDistributeSection(PetscSF, PetscSection, PetscInt **, PetscSection);
186: PETSC_EXTERN PetscErrorCode PetscSFCreateSectionSF(PetscSF, PetscSection, PetscInt [], PetscSection, PetscSF *);
188: PETSC_EXTERN PetscErrorCode VecGetValuesSection(Vec, PetscSection, PetscInt, PetscScalar **);
189: PETSC_EXTERN PetscErrorCode VecSetValuesSection(Vec, PetscSection, PetscInt, PetscScalar [], InsertMode);
191: /* ----------------------------------------------------------------------------*/
193: typedef struct _VecOps *VecOps;
194: struct _VecOps {
195: PetscErrorCode (*duplicate)(Vec,Vec*); /* get single vector */
196: PetscErrorCode (*duplicatevecs)(Vec,PetscInt,Vec**); /* get array of vectors */
197: PetscErrorCode (*destroyvecs)(PetscInt,Vec[]); /* free array of vectors */
198: PetscErrorCode (*dot)(Vec,Vec,PetscScalar*); /* z = x^H * y */
199: PetscErrorCode (*mdot)(Vec,PetscInt,const Vec[],PetscScalar*); /* z[j] = x dot y[j] */
200: PetscErrorCode (*norm)(Vec,NormType,PetscReal*); /* z = sqrt(x^H * x) */
201: PetscErrorCode (*tdot)(Vec,Vec,PetscScalar*); /* x'*y */
202: PetscErrorCode (*mtdot)(Vec,PetscInt,const Vec[],PetscScalar*);/* z[j] = x dot y[j] */
203: PetscErrorCode (*scale)(Vec,PetscScalar); /* x = alpha * x */
204: PetscErrorCode (*copy)(Vec,Vec); /* y = x */
205: PetscErrorCode (*set)(Vec,PetscScalar); /* y = alpha */
206: PetscErrorCode (*swap)(Vec,Vec); /* exchange x and y */
207: PetscErrorCode (*axpy)(Vec,PetscScalar,Vec); /* y = y + alpha * x */
208: PetscErrorCode (*axpby)(Vec,PetscScalar,PetscScalar,Vec); /* y = alpha * x + beta * y*/
209: PetscErrorCode (*maxpy)(Vec,PetscInt,const PetscScalar*,Vec*); /* y = y + alpha[j] x[j] */
210: PetscErrorCode (*aypx)(Vec,PetscScalar,Vec); /* y = x + alpha * y */
211: PetscErrorCode (*waxpy)(Vec,PetscScalar,Vec,Vec); /* w = y + alpha * x */
212: PetscErrorCode (*axpbypcz)(Vec,PetscScalar,PetscScalar,PetscScalar,Vec,Vec); /* z = alpha * x + beta *y + gamma *z*/
213: PetscErrorCode (*pointwisemult)(Vec,Vec,Vec); /* w = x .* y */
214: PetscErrorCode (*pointwisedivide)(Vec,Vec,Vec); /* w = x ./ y */
215: PetscErrorCode (*setvalues)(Vec,PetscInt,const PetscInt[],const PetscScalar[],InsertMode);
216: PetscErrorCode (*assemblybegin)(Vec); /* start global assembly */
217: PetscErrorCode (*assemblyend)(Vec); /* end global assembly */
218: PetscErrorCode (*getarray)(Vec,PetscScalar**); /* get data array */
219: PetscErrorCode (*getsize)(Vec,PetscInt*);
220: PetscErrorCode (*getlocalsize)(Vec,PetscInt*);
221: PetscErrorCode (*restorearray)(Vec,PetscScalar**); /* restore data array */
222: PetscErrorCode (*max)(Vec,PetscInt*,PetscReal*); /* z = max(x); idx=index of max(x) */
223: PetscErrorCode (*min)(Vec,PetscInt*,PetscReal*); /* z = min(x); idx=index of min(x) */
224: PetscErrorCode (*setrandom)(Vec,PetscRandom); /* set y[j] = random numbers */
225: PetscErrorCode (*setoption)(Vec,VecOption,PetscBool );
226: PetscErrorCode (*setvaluesblocked)(Vec,PetscInt,const PetscInt[],const PetscScalar[],InsertMode);
227: PetscErrorCode (*destroy)(Vec);
228: PetscErrorCode (*view)(Vec,PetscViewer);
229: PetscErrorCode (*placearray)(Vec,const PetscScalar*); /* place data array */
230: PetscErrorCode (*replacearray)(Vec,const PetscScalar*); /* replace data array */
231: PetscErrorCode (*dot_local)(Vec,Vec,PetscScalar*);
232: PetscErrorCode (*tdot_local)(Vec,Vec,PetscScalar*);
233: PetscErrorCode (*norm_local)(Vec,NormType,PetscReal*);
234: PetscErrorCode (*mdot_local)(Vec,PetscInt,const Vec[],PetscScalar*);
235: PetscErrorCode (*mtdot_local)(Vec,PetscInt,const Vec[],PetscScalar*);
236: PetscErrorCode (*load)(Vec,PetscViewer);
237: PetscErrorCode (*reciprocal)(Vec);
238: PetscErrorCode (*conjugate)(Vec);
239: PetscErrorCode (*setlocaltoglobalmapping)(Vec,ISLocalToGlobalMapping);
240: PetscErrorCode (*setvalueslocal)(Vec,PetscInt,const PetscInt *,const PetscScalar *,InsertMode);
241: PetscErrorCode (*resetarray)(Vec); /* vector points to its original array, i.e. undoes any VecPlaceArray() */
242: PetscErrorCode (*setfromoptions)(Vec);
243: PetscErrorCode (*maxpointwisedivide)(Vec,Vec,PetscReal*); /* m = max abs(x ./ y) */
244: PetscErrorCode (*pointwisemax)(Vec,Vec,Vec);
245: PetscErrorCode (*pointwisemaxabs)(Vec,Vec,Vec);
246: PetscErrorCode (*pointwisemin)(Vec,Vec,Vec);
247: PetscErrorCode (*getvalues)(Vec,PetscInt,const PetscInt[],PetscScalar[]);
248: PetscErrorCode (*sqrt)(Vec);
249: PetscErrorCode (*abs)(Vec);
250: PetscErrorCode (*exp)(Vec);
251: PetscErrorCode (*log)(Vec);
252: PetscErrorCode (*shift)(Vec);
253: PetscErrorCode (*create)(Vec);
254: PetscErrorCode (*stridegather)(Vec,PetscInt,Vec,InsertMode);
255: PetscErrorCode (*stridescatter)(Vec,PetscInt,Vec,InsertMode);
256: PetscErrorCode (*dotnorm2)(Vec,Vec,PetscScalar*,PetscScalar*);
257: PetscErrorCode (*getsubvector)(Vec,IS,Vec*);
258: PetscErrorCode (*restoresubvector)(Vec,IS,Vec*);
259: };
261: /*
262: The stash is used to temporarily store inserted vec values that
263: belong to another processor. During the assembly phase the stashed
264: values are moved to the correct processor and
265: */
267: typedef struct {
268: PetscInt nmax; /* maximum stash size */
269: PetscInt umax; /* max stash size user wants */
270: PetscInt oldnmax; /* the nmax value used previously */
271: PetscInt n; /* stash size */
272: PetscInt bs; /* block size of the stash */
273: PetscInt reallocs; /* preserve the no of mallocs invoked */
274: PetscInt *idx; /* global row numbers in stash */
275: PetscScalar *array; /* array to hold stashed values */
276: /* The following variables are used for communication */
277: MPI_Comm comm;
278: PetscMPIInt size,rank;
279: PetscMPIInt tag1,tag2;
280: MPI_Request *send_waits; /* array of send requests */
281: MPI_Request *recv_waits; /* array of receive requests */
282: MPI_Status *send_status; /* array of send status */
283: PetscInt nsends,nrecvs; /* numbers of sends and receives */
284: PetscScalar *svalues,*rvalues; /* sending and receiving data */
285: PetscInt *sindices,*rindices;
286: PetscInt rmax; /* maximum message length */
287: PetscInt *nprocs; /* tmp data used both during scatterbegin and end */
288: PetscInt nprocessed; /* number of messages already processed */
289: PetscBool donotstash;
290: PetscBool ignorenegidx; /* ignore negative indices passed into VecSetValues/VetGetValues */
291: InsertMode insertmode;
292: PetscInt *bowners;
293: } VecStash;
295: #if defined(PETSC_HAVE_CUSP)
296: /*E
297: PetscCUSPFlag - indicates which memory (CPU, GPU, or none contains valid vector
299: PETSC_CUSP_UNALLOCATED - no memory contains valid matrix entries; NEVER used for vectors
300: PETSC_CUSP_GPU - GPU has valid vector/matrix entries
301: PETSC_CUSP_CPU - CPU has valid vector/matrix entries
302: PETSC_CUSP_BOTH - Both GPU and CPU have valid vector/matrix entries and they match
304: Level: developer
305: E*/
306: typedef enum {PETSC_CUSP_UNALLOCATED,PETSC_CUSP_GPU,PETSC_CUSP_CPU,PETSC_CUSP_BOTH} PetscCUSPFlag;
307: #endif
309: struct _p_Vec {
310: PETSCHEADER(struct _VecOps);
311: PetscLayout map;
312: void *data; /* implementation-specific data */
313: PetscBool array_gotten;
314: VecStash stash,bstash; /* used for storing off-proc values during assembly */
315: PetscBool petscnative; /* means the ->data starts with VECHEADER and can use VecGetArrayFast()*/
316: #if defined(PETSC_HAVE_CUSP)
317: PetscCUSPFlag valid_GPU_array; /* indicates where the most recently modified vector data is (GPU or CPU) */
318: void *spptr; /* if we're using CUSP, then this is the special pointer to the array on the GPU */
319: #endif
320: };
322: PETSC_EXTERN PetscLogEvent VEC_View, VEC_Max, VEC_Min, VEC_DotBarrier, VEC_Dot, VEC_MDotBarrier, VEC_MDot, VEC_TDot, VEC_MTDot;
323: PETSC_EXTERN PetscLogEvent VEC_Norm, VEC_Normalize, VEC_Scale, VEC_Copy, VEC_Set, VEC_AXPY, VEC_AYPX, VEC_WAXPY, VEC_MAXPY;
324: PETSC_EXTERN PetscLogEvent VEC_AssemblyEnd, VEC_PointwiseMult, VEC_SetValues, VEC_Load, VEC_ScatterBarrier, VEC_ScatterBegin, VEC_ScatterEnd;
325: PETSC_EXTERN PetscLogEvent VEC_SetRandom, VEC_ReduceArithmetic, VEC_ReduceBarrier, VEC_ReduceCommunication;
326: PETSC_EXTERN PetscLogEvent VEC_ReduceBegin,VEC_ReduceEnd;
327: PETSC_EXTERN PetscLogEvent VEC_Swap, VEC_AssemblyBegin, VEC_NormBarrier, VEC_DotNormBarrier, VEC_DotNorm, VEC_AXPBYPCZ, VEC_Ops;
328: PETSC_EXTERN PetscLogEvent VEC_CUSPCopyToGPU, VEC_CUSPCopyFromGPU;
329: PETSC_EXTERN PetscLogEvent VEC_CUSPCopyToGPUSome, VEC_CUSPCopyFromGPUSome;
331: #if defined(PETSC_HAVE_CUSP)
332: PETSC_EXTERN PetscErrorCode VecCUSPCopyFromGPU(Vec v);
333: #endif
337: PETSC_STATIC_INLINE PetscErrorCode VecGetArrayRead(Vec x,const PetscScalar *a[])
338: {
342: if (x->petscnative){
343: #if defined(PETSC_HAVE_CUSP)
344: if (x->valid_GPU_array == PETSC_CUSP_GPU || !*((PetscScalar**)x->data)){
345: VecCUSPCopyFromGPU(x);
346: }
347: #endif
348: *a = *((PetscScalar **)x->data);
349: } else {
350: (*x->ops->getarray)(x,(PetscScalar**)a);
351: }
352: return(0);
353: }
357: PETSC_STATIC_INLINE PetscErrorCode VecRestoreArrayRead(Vec x,const PetscScalar *a[])
358: {
362: if (!x->petscnative){
363: (*x->ops->restorearray)(x,(PetscScalar**)a);
364: }
365: if (a) *a = PETSC_NULL;
366: return(0);
367: }
371: PETSC_STATIC_INLINE PetscErrorCode VecGetArray(Vec x,PetscScalar *a[])
372: {
377: if (x->petscnative){
378: #if defined(PETSC_HAVE_CUSP)
379: if (x->valid_GPU_array == PETSC_CUSP_GPU || !*((PetscScalar**)x->data)){
380: VecCUSPCopyFromGPU(x);
381: }
382: #endif
383: *a = *((PetscScalar **)x->data);
384: } else {
385: (*x->ops->getarray)(x,a);
386: }
387: return(0);
388: }
392: PETSC_STATIC_INLINE PetscErrorCode VecRestoreArray(Vec x,PetscScalar *a[])
393: {
397: if (x->petscnative){
398: #if defined(PETSC_HAVE_CUSP)
399: x->valid_GPU_array = PETSC_CUSP_CPU;
400: #endif
401: } else {
402: (*x->ops->restorearray)(x,a);
403: }
404: PetscObjectStateIncrease((PetscObject)x);
405: if (a) *a = PETSC_NULL;
406: return(0);
407: }
410: /*
411: Common header shared by array based vectors,
412: currently Vec_Seq and Vec_MPI
413: */
414: #define VECHEADER \
415: PetscScalar *array; \
416: PetscScalar *array_allocated; /* if the array was allocated by PETSc this is its pointer */ \
417: PetscScalar *unplacedarray; /* if one called VecPlaceArray(), this is where it stashed the original */
419: /* Default obtain and release vectors; can be used by any implementation */
420: PETSC_EXTERN PetscErrorCode VecDuplicateVecs_Default(Vec,PetscInt,Vec *[]);
421: PETSC_EXTERN PetscErrorCode VecDestroyVecs_Default(PetscInt,Vec []);
422: PETSC_EXTERN PetscErrorCode VecLoad_Binary(Vec, PetscViewer);
423: PETSC_EXTERN PetscErrorCode VecLoad_Default(Vec, PetscViewer);
425: PETSC_EXTERN PetscInt NormIds[7]; /* map from NormType to IDs used to cache/retreive values of norms */
427: /* --------------------------------------------------------------------*/
428: /* */
429: /* Defines the data structures used in the Vec Scatter operations */
431: typedef enum { VEC_SCATTER_SEQ_GENERAL,VEC_SCATTER_SEQ_STRIDE,
432: VEC_SCATTER_MPI_GENERAL,VEC_SCATTER_MPI_TOALL,
433: VEC_SCATTER_MPI_TOONE} VecScatterType;
435: /*
436: These scatters are for the purely local case.
437: */
438: typedef struct {
439: VecScatterType type;
440: PetscInt n; /* number of components to scatter */
441: PetscInt *vslots; /* locations of components */
442: /*
443: The next three fields are used in parallel scatters, they contain
444: optimization in the special case that the "to" vector and the "from"
445: vector are the same, so one only needs copy components that truly
446: copies instead of just y[idx[i]] = y[jdx[i]] where idx[i] == jdx[i].
447: */
448: PetscBool nonmatching_computed;
449: PetscInt n_nonmatching; /* number of "from"s != "to"s */
450: PetscInt *slots_nonmatching; /* locations of "from"s != "to"s */
451: PetscBool is_copy;
452: PetscInt copy_start; /* local scatter is a copy starting at copy_start */
453: PetscInt copy_length;
454: } VecScatter_Seq_General;
456: typedef struct {
457: VecScatterType type;
458: PetscInt n;
459: PetscInt first;
460: PetscInt step;
461: } VecScatter_Seq_Stride;
463: /*
464: This scatter is for a global vector copied (completely) to each processor (or all to one)
465: */
466: typedef struct {
467: VecScatterType type;
468: PetscMPIInt *count; /* elements of vector on each processor */
469: PetscMPIInt *displx;
470: PetscScalar *work1;
471: PetscScalar *work2;
472: } VecScatter_MPI_ToAll;
474: /*
475: This is the general parallel scatter
476: */
477: typedef struct {
478: VecScatterType type;
479: PetscInt n; /* number of processors to send/receive */
480: PetscInt *starts; /* starting point in indices and values for each proc*/
481: PetscInt *indices; /* list of all components sent or received */
482: PetscMPIInt *procs; /* processors we are communicating with in scatter */
483: MPI_Request *requests,*rev_requests;
484: PetscScalar *values; /* buffer for all sends or receives */
485: VecScatter_Seq_General local; /* any part that happens to be local */
486: MPI_Status *sstatus,*rstatus;
487: PetscBool use_readyreceiver;
488: PetscInt bs;
489: PetscBool sendfirst;
490: PetscBool contiq;
491: /* for MPI_Alltoallv() approach */
492: PetscBool use_alltoallv;
493: PetscMPIInt *counts,*displs;
494: /* for MPI_Alltoallw() approach */
495: PetscBool use_alltoallw;
496: #if defined(PETSC_HAVE_MPI_ALLTOALLW)
497: PetscMPIInt *wcounts,*wdispls;
498: MPI_Datatype *types;
499: #endif
500: PetscBool use_window;
501: #if defined(PETSC_HAVE_MPI_WIN_CREATE)
502: MPI_Win window;
503: PetscInt *winstarts; /* displacements in the processes I am putting to */
504: #endif
505: } VecScatter_MPI_General;
507: struct _p_VecScatter {
508: PETSCHEADER(int);
509: PetscInt to_n,from_n;
510: PetscBool inuse; /* prevents corruption from mixing two scatters */
511: PetscBool beginandendtogether; /* indicates that the scatter begin and end function are called together, VecScatterEnd()
512: is then treated as a nop */
513: PetscBool packtogether; /* packs all the messages before sending, same with receive */
514: PetscBool reproduce; /* always receive the ghost points in the same order of processes */
515: PetscErrorCode (*begin)(VecScatter,Vec,Vec,InsertMode,ScatterMode);
516: PetscErrorCode (*end)(VecScatter,Vec,Vec,InsertMode,ScatterMode);
517: PetscErrorCode (*copy)(VecScatter,VecScatter);
518: PetscErrorCode (*destroy)(VecScatter);
519: PetscErrorCode (*view)(VecScatter,PetscViewer);
520: void *fromdata,*todata;
521: void *spptr;
522: };
524: PETSC_EXTERN PetscErrorCode VecStashCreate_Private(MPI_Comm,PetscInt,VecStash*);
525: PETSC_EXTERN PetscErrorCode VecStashDestroy_Private(VecStash*);
526: PETSC_EXTERN PetscErrorCode VecStashExpand_Private(VecStash*,PetscInt);
527: PETSC_EXTERN PetscErrorCode VecStashScatterEnd_Private(VecStash*);
528: PETSC_EXTERN PetscErrorCode VecStashSetInitialSize_Private(VecStash*,PetscInt);
529: PETSC_EXTERN PetscErrorCode VecStashGetInfo_Private(VecStash*,PetscInt*,PetscInt*);
530: PETSC_EXTERN PetscErrorCode VecStashScatterBegin_Private(VecStash*,PetscInt*);
531: PETSC_EXTERN PetscErrorCode VecStashScatterGetMesg_Private(VecStash*,PetscMPIInt*,PetscInt**,PetscScalar**,PetscInt*);
533: /*
534: VecStashValue_Private - inserts a single value into the stash.
536: Input Parameters:
537: stash - the stash
538: idx - the global of the inserted value
539: values - the value inserted
540: */
541: PETSC_STATIC_INLINE PetscErrorCode VecStashValue_Private(VecStash *stash,PetscInt row,PetscScalar value)
542: {
544: /* Check and see if we have sufficient memory */
545: if (((stash)->n + 1) > (stash)->nmax) {
546: VecStashExpand_Private(stash,1);
547: }
548: (stash)->idx[(stash)->n] = row;
549: (stash)->array[(stash)->n] = value;
550: (stash)->n++;
551: return 0;
552: }
554: /*
555: VecStashValuesBlocked_Private - inserts 1 block of values into the stash.
557: Input Parameters:
558: stash - the stash
559: idx - the global block index
560: values - the values inserted
561: */
562: PETSC_STATIC_INLINE PetscErrorCode VecStashValuesBlocked_Private(VecStash *stash,PetscInt row,PetscScalar *values)
563: {
564: PetscInt jj,stash_bs=(stash)->bs;
565: PetscScalar *array;
567: if (((stash)->n+1) > (stash)->nmax) {
568: VecStashExpand_Private(stash,1);
569: }
570: array = (stash)->array + stash_bs*(stash)->n;
571: (stash)->idx[(stash)->n] = row;
572: for (jj=0; jj<stash_bs; jj++) { array[jj] = values[jj];}
573: (stash)->n++;
574: return 0;
575: }
577: PETSC_EXTERN PetscErrorCode VecStrideGather_Default(Vec,PetscInt,Vec,InsertMode);
578: PETSC_EXTERN PetscErrorCode VecStrideScatter_Default(Vec,PetscInt,Vec,InsertMode);
579: PETSC_EXTERN PetscErrorCode VecReciprocal_Default(Vec);
581: #if defined(PETSC_HAVE_MATLAB_ENGINE)
582: EXTERN_C_BEGIN
583: PETSC_EXTERN PetscErrorCode VecMatlabEnginePut_Default(PetscObject,void*);
584: PETSC_EXTERN PetscErrorCode VecMatlabEngineGet_Default(PetscObject,void*);
585: EXTERN_C_END
586: #endif
589: /* Reset __FUNCT__ in case the user does not define it themselves */
593: #endif