Actual source code: vecimpl.h

  2: /* $Id: vecimpl.h,v 1.89 2001/09/19 16:07:46 bsmith Exp $ */

  4: /* 
  5:    This private file should not be included in users' code.
  6:    Defines the fields shared by all vector implementations.
  7: */

  9: #ifndef __VECIMPL_H

 12:  #include petscvec.h

 14: struct _PetscMapOps {
 15:   int (*setfromoptions)(PetscMap),
 16:       (*destroy)(PetscMap);
 17: };

 19: struct _p_PetscMap {
 20:   PETSCHEADER(struct _PetscMapOps)
 21:   int  n,N;         /* local, global vector size */
 22:   int  rstart,rend; /* local start, local end + 1 */
 23:   int *range;       /* the offset of each processor */
 24: };

 26: /* ----------------------------------------------------------------------------*/

 28: typedef struct _VecOps *VecOps;
 29: struct _VecOps {
 30:   int  (*duplicate)(Vec,Vec*),              /* get single vector */
 31:        (*duplicatevecs)(Vec,int,Vec**),     /* get array of vectors */
 32:        (*destroyvecs)(const Vec[],int),     /* free array of vectors */
 33:        (*dot)(Vec,Vec,PetscScalar*),             /* z = x^H * y */
 34:        (*mdot)(int,Vec,const Vec[],PetscScalar*), /* z[j] = x dot y[j] */
 35:        (*norm)(Vec,NormType,PetscReal*),        /* z = sqrt(x^H * x) */
 36:        (*tdot)(Vec,Vec,PetscScalar*),             /* x'*y */
 37:        (*mtdot)(int,Vec,const Vec[],PetscScalar*),/* z[j] = x dot y[j] */
 38:        (*scale)(const PetscScalar*,Vec),          /* x = alpha * x   */
 39:        (*copy)(Vec,Vec),                     /* y = x */
 40:        (*set)(const PetscScalar*,Vec),            /* y = alpha  */
 41:        (*swap)(Vec,Vec),                     /* exchange x and y */
 42:        (*axpy)(const PetscScalar*,Vec,Vec),       /* y = y + alpha * x */
 43:        (*axpby)(const PetscScalar*,const PetscScalar*,Vec,Vec), /* y = y + alpha * x + beta * y*/
 44:        (*maxpy)(int,const PetscScalar*,Vec,Vec*), /* y = y + alpha[j] x[j] */
 45:        (*aypx)(const PetscScalar*,Vec,Vec),       /* y = x + alpha * y */
 46:        (*waxpy)(const PetscScalar*,Vec,Vec,Vec),  /* w = y + alpha * x */
 47:        (*pointwisemult)(Vec,Vec,Vec),        /* w = x .* y */
 48:        (*pointwisedivide)(Vec,Vec,Vec),      /* w = x ./ y */
 49:        (*setvalues)(Vec,int,const int[],const PetscScalar[],InsertMode),
 50:        (*assemblybegin)(Vec),                /* start global assembly */
 51:        (*assemblyend)(Vec),                  /* end global assembly */
 52:        (*getarray)(Vec,PetscScalar**),            /* get data array */
 53:        (*getsize)(Vec,int*),
 54:        (*getlocalsize)(Vec,int*),
 55:        (*restorearray)(Vec,PetscScalar**),        /* restore data array */
 56:        (*max)(Vec,int*,PetscReal*),      /* z = max(x); idx=index of max(x) */
 57:        (*min)(Vec,int*,PetscReal*),      /* z = min(x); idx=index of min(x) */
 58:        (*setrandom)(PetscRandom,Vec),        /* set y[j] = random numbers */
 59:        (*setoption)(Vec,VecOption),
 60:        (*setvaluesblocked)(Vec,int,const int[],const PetscScalar[],InsertMode),
 61:        (*destroy)(Vec),
 62:        (*view)(Vec,PetscViewer),
 63:        (*placearray)(Vec,const PetscScalar*),     /* place data array */
 64:        (*replacearray)(Vec,const PetscScalar*),     /* replace data array */
 65:        (*dot_local)(Vec,Vec,PetscScalar*),
 66:        (*tdot_local)(Vec,Vec,PetscScalar*),
 67:        (*norm_local)(Vec,NormType,PetscReal*),
 68:        (*loadintovector)(PetscViewer,Vec),
 69:        (*reciprocal)(Vec),
 70:        (*viewnative)(Vec,PetscViewer),
 71:        (*conjugate)(Vec),
 72:        (*setlocaltoglobalmapping)(Vec,ISLocalToGlobalMapping),
 73:        (*setvalueslocal)(Vec,int,const int *,const PetscScalar *,InsertMode),
 74:        (*resetarray)(Vec),      /* vector points to its original array, i.e. undoes any VecPlaceArray() */
 75:        (*setfromoptions)(Vec),
 76:        (*maxpointwisedivide)(Vec,Vec,PetscReal*),      /* m = max abs(x ./ y) */
 77:        (*load)(PetscViewer,const VecType,Vec*);
 78: };

 80: /* 
 81:     The stash is used to temporarily store inserted vec values that 
 82:   belong to another processor. During the assembly phase the stashed 
 83:   values are moved to the correct processor and 
 84: */

 86: typedef struct {
 87:   int           nmax;                   /* maximum stash size */
 88:   int           umax;                   /* max stash size user wants */
 89:   int           oldnmax;                /* the nmax value used previously */
 90:   int           n;                      /* stash size */
 91:   int           bs;                     /* block size of the stash */
 92:   int           reallocs;               /* preserve the no of mallocs invoked */
 93:   int           *idx;                   /* global row numbers in stash */
 94:   PetscScalar   *array;                 /* array to hold stashed values */
 95:   /* The following variables are used for communication */
 96:   MPI_Comm      comm;
 97:   int           size,rank;
 98:   int           tag1,tag2;
 99:   MPI_Request   *send_waits;            /* array of send requests */
100:   MPI_Request   *recv_waits;            /* array of receive requests */
101:   MPI_Status    *send_status;           /* array of send status */
102:   int           nsends,nrecvs;          /* numbers of sends and receives */
103:   PetscScalar   *svalues,*rvalues;      /* sending and receiving data */
104:   int           rmax;                   /* maximum message length */
105:   int           *nprocs;                /* tmp data used both duiring scatterbegin and end */
106:   int           nprocessed;             /* number of messages already processed */
107:   PetscTruth    donotstash;
108:   InsertMode    insertmode;
109:   int           *bowners;
110: } VecStash;

112: struct _p_Vec {
113:   PETSCHEADER(struct _VecOps)
114:   PetscMap               map;
115:   void                   *data;     /* implementation-specific data */
116:   int                    N,n;      /* global, local vector size */
117:   int                    bs;
118:   ISLocalToGlobalMapping mapping;   /* mapping used in VecSetValuesLocal() */
119:   ISLocalToGlobalMapping bmapping;  /* mapping used in VecSetValuesBlockedLocal() */
120:   PetscTruth             array_gotten;
121:   VecStash               stash,bstash; /* used for storing off-proc values during assembly */
122:   PetscTruth             petscnative;  /* means the ->data starts with VECHEADER and can use VecGetArrayFast()*/
123:   void                   *esivec;      /* ESI wrapper of vector */
124: };

126: #define VecGetArray(x,a)     ((x)->petscnative ? (*(a) = *((PetscScalar **)(x)->data),0) : VecGetArray_Private((x),(a)))
127: #define VecRestoreArray(x,a) ((x)->petscnative ? PetscObjectIncreaseState((PetscObject)x) : VecRestoreArray_Private((x),(a)))

129: /*
130:      Common header shared by array based vectors, 
131:    currently Vec_Seq and Vec_MPI
132: */
133: #define VECHEADER                         \
134:   PetscScalar *array;                          \
135:   PetscScalar *array_allocated;            

137: /* Default obtain and release vectors; can be used by any implementation */
138: EXTERN int VecDuplicateVecs_Default(Vec,int,Vec *[]);
139: EXTERN int VecDestroyVecs_Default(const Vec [],int);

141: EXTERN int VecLoadIntoVector_Default(PetscViewer,Vec);

143: /* --------------------------------------------------------------------*/
144: /*                                                                     */
145: /* Defines the data structures used in the Vec Scatter operations      */

147: typedef enum { VEC_SCATTER_SEQ_GENERAL,VEC_SCATTER_SEQ_STRIDE,
148:                VEC_SCATTER_MPI_GENERAL,VEC_SCATTER_MPI_TOALL,
149:                VEC_SCATTER_MPI_TOONE} VecScatterType;

151: /* 
152:    These scatters are for the purely local case.
153: */
154: typedef struct {
155:   VecScatterType type;
156:   int            n;                    /* number of components to scatter */
157:   int            *slots;               /* locations of components */
158:   /*
159:        The next three fields are used in parallel scatters, they contain 
160:        optimization in the special case that the "to" vector and the "from" 
161:        vector are the same, so one only needs copy components that truly 
162:        copies instead of just y[idx[i]] = y[jdx[i]] where idx[i] == jdx[i].
163:   */
164:   PetscTruth     nonmatching_computed;
165:   int            n_nonmatching;        /* number of "from"s  != "to"s */
166:   int            *slots_nonmatching;   /* locations of "from"s  != "to"s */
167:   PetscTruth     is_copy;
168:   int            copy_start;   /* local scatter is a copy starting at copy_start */
169:   int            copy_length;
170: } VecScatter_Seq_General;

172: typedef struct {
173:   VecScatterType type;
174:   int            n;
175:   int            first;
176:   int            step;
177: } VecScatter_Seq_Stride;

179: /*
180:    This scatter is for a global vector copied (completely) to each processor (or all to one)
181: */
182: typedef struct {
183:   VecScatterType type;
184:   int            *count;        /* elements of vector on each processor */
185:   PetscScalar    *work1;
186:   PetscScalar    *work2;
187: } VecScatter_MPI_ToAll;

189: /*
190:    This is the general parallel scatter
191: */
192: typedef struct {
193:   VecScatterType         type;
194:   int                    n;        /* number of processors to send/receive */
195:   int                    *starts;  /* starting point in indices and values for each proc*/
196:   int                    *indices; /* list of all components sent or received */
197:   int                    *procs;   /* processors we are communicating with in scatter */
198:   MPI_Request            *requests,*rev_requests;
199:   PetscScalar            *values;  /* buffer for all sends or receives */
200:   VecScatter_Seq_General local;    /* any part that happens to be local */
201:   MPI_Status             *sstatus,*rstatus;
202:   PetscTruth             use_readyreceiver;
203:   int                    bs;
204:   PetscTruth             sendfirst;
205: } VecScatter_MPI_General;

207: struct _p_VecScatter {
208:   PETSCHEADER(int)
209:   int        to_n,from_n;
210:   PetscTruth inuse;   /* prevents corruption from mixing two scatters */
211:   PetscTruth beginandendtogether;         /* indicates that the scatter begin and end
212:                                           function are called together, VecScatterEnd()
213:                                           is then treated as a nop */
214:   PetscTruth packtogether; /* packs all the messages before sending, same with receive */
215:   int        (*postrecvs)(Vec,Vec,InsertMode,ScatterMode,VecScatter);
216:   int        (*begin)(Vec,Vec,InsertMode,ScatterMode,VecScatter);
217:   int        (*end)(Vec,Vec,InsertMode,ScatterMode,VecScatter);
218:   int        (*copy)(VecScatter,VecScatter);
219:   int        (*destroy)(VecScatter);
220:   int        (*view)(VecScatter,PetscViewer);
221:   void       *fromdata,*todata;
222: };

224: EXTERN int VecStashCreate_Private(MPI_Comm,int,VecStash*);
225: EXTERN int VecStashDestroy_Private(VecStash*);
226: EXTERN int VecStashExpand_Private(VecStash*,int);
227: EXTERN int VecStashScatterEnd_Private(VecStash*);
228: EXTERN int VecStashSetInitialSize_Private(VecStash*,int);
229: EXTERN int VecStashGetInfo_Private(VecStash*,int*,int*);
230: EXTERN int VecStashScatterBegin_Private(VecStash*,int*);
231: EXTERN int VecStashScatterGetMesg_Private(VecStash*,int*,int**,PetscScalar**,int*);

233: /* 
234:    The following are implemented as macros to avoid the function
235:    call overhead.

237:    extern int VecStashValue_Private(VecStash*,int,PetscScalar);
238:    extern int VecStashValuesBlocked_Private(VecStash*,int,PetscScalar*);
239: */

241: /*
242:   VecStashValue_Private - inserts a single value into the stash.

244:   Input Parameters:
245:   stash  - the stash
246:   idx    - the global of the inserted value
247:   values - the value inserted
248: */
249: #define VecStashValue_Private(stash,row,value) \
250: {  \
251:   /* Check and see if we have sufficient memory */ \
252:   if (((stash)->n + 1) > (stash)->nmax) { \
253:     VecStashExpand_Private(stash,1); \
254:   } \
255:   (stash)->idx[(stash)->n]   = row; \
256:   (stash)->array[(stash)->n] = value; \
257:   (stash)->n++; \
258: }

260: /*
261:   VecStashValuesBlocked_Private - inserts 1 block of values into the stash. 

263:   Input Parameters:
264:   stash  - the stash
265:   idx    - the global block index
266:   values - the values inserted
267: */
268: #define VecStashValuesBlocked_Private(stash,row,values) \
269: { \
270:   int    jj,stash_bs=(stash)->bs; \
271:   PetscScalar *array; \
272:   if (((stash)->n+1) > (stash)->nmax) { \
273:     VecStashExpand_Private(stash,1); \
274:   } \
275:   array = (stash)->array + stash_bs*(stash)->n; \
276:   (stash)->idx[(stash)->n]   = row; \
277:   for (jj=0; jj<stash_bs; jj++) { array[jj] = values[jj];} \
278:   (stash)->n++; \
279: }

281: EXTERN int VecReciprocal_Default(Vec);

283: #if defined(PETSC_HAVE_MATLAB) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
284: EXTERN_C_BEGIN
285: EXTERN int VecMatlabEnginePut_Default(PetscObject,void*);
286: EXTERN int VecMatlabEngineGet_Default(PetscObject,void*);
287: EXTERN_C_END
288: #endif


291: #endif