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: };

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

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

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

124: #define VecGetArrayFast(x,a)     ((x)->petscnative ? (*(a) = *((PetscScalar **)(x)->data),0) : VecGetArray((x),(a)))
125: #define VecRestoreArrayFast(x,a) ((x)->petscnative ? 0 : VecRestoreArray((x),(a)))

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

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

139: EXTERN int VecLoadIntoVector_Default(PetscViewer,Vec);

141: /* --------------------------------------------------------------------*/
142: /*                                                                     */
143: /* Defines the data structures used in the Vec Scatter operations      */

145: typedef enum { VEC_SCATTER_SEQ_GENERAL,VEC_SCATTER_SEQ_STRIDE,
146:                VEC_SCATTER_MPI_GENERAL,VEC_SCATTER_MPI_TOALL,
147:                VEC_SCATTER_MPI_TOONE} VecScatterType;

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

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

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

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

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

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

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

235:    extern int VecStashValue_Private(VecStash*,int,PetscScalar);
236:    extern int VecStashValuesBlocked_Private(VecStash*,int,PetscScalar*);
237: */

239: /*
240:   VecStashValue_Private - inserts a single values into the stash.

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

258: /*
259:   VecStashValuesBlocked_Private - inserts 1 block of values into the stash. 

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

279: EXTERN int VecReciprocal_Default(Vec);

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


289: #endif