Actual source code: petscimpl.h

petsc-3.6.1 2015-07-22
Report Typos and Errors
  2: /*
  3:     Defines the basic header of all PETSc objects.
  4: */

  6: #if !defined(_PETSCHEAD_H)
  7: #define _PETSCHEAD_H
  8: #include <petscsys.h>

 10: /*
 11:    All major PETSc data structures have a common core; this is defined
 12:    below by PETSCHEADER.

 14:    PetscHeaderCreate() should be used whenever creating a PETSc structure.
 15: */

 17: /*
 18:    PetscOps: structure of core operations that all PETSc objects support.

 20:       getcomm()         - Gets the object's communicator.
 21:       view()            - Is the routine for viewing the entire PETSc object; for
 22:                           example, MatView() is the general matrix viewing routine.
 23:                           This is used by PetscObjectView((PetscObject)obj) to allow
 24:                           viewing any PETSc object.
 25:       destroy()         - Is the routine for destroying the entire PETSc object;
 26:                           for example,MatDestroy() is the general matrix
 27:                           destruction routine.
 28:                           This is used by PetscObjectDestroy((PetscObject*)&obj) to allow
 29:                           destroying any PETSc object.
 30:       compose()         - Associates a PETSc object with another PETSc object with a name
 31:       query()           - Returns a different PETSc object that has been associated
 32:                           with the first object using a name.
 33:       composefunction() - Attaches an a function to a PETSc object with a name.
 34:       queryfunction()   - Requests a registered function that has been attached to a PETSc object.
 35: */

 37: typedef struct {
 38:    PetscErrorCode (*getcomm)(PetscObject,MPI_Comm *);
 39:    PetscErrorCode (*view)(PetscObject,PetscViewer);
 40:    PetscErrorCode (*destroy)(PetscObject*);
 41:    PetscErrorCode (*compose)(PetscObject,const char[],PetscObject);
 42:    PetscErrorCode (*query)(PetscObject,const char[],PetscObject *);
 43:    PetscErrorCode (*composefunction)(PetscObject,const char[],void (*)(void));
 44:    PetscErrorCode (*queryfunction)(PetscObject,const char[],void (**)(void));
 45: } PetscOps;

 47: typedef enum {PETSC_FORTRAN_CALLBACK_CLASS,PETSC_FORTRAN_CALLBACK_SUBTYPE,PETSC_FORTRAN_CALLBACK_MAXTYPE} PetscFortranCallbackType;
 48: typedef int PetscFortranCallbackId;
 49: #define PETSC_SMALLEST_FORTRAN_CALLBACK ((PetscFortranCallbackId)1000)
 50: PETSC_EXTERN PetscErrorCode PetscFortranCallbackRegister(PetscClassId,const char*,PetscFortranCallbackId*);
 51: PETSC_EXTERN PetscErrorCode PetscFortranCallbackGetSizes(PetscClassId,PetscInt*,PetscInt*);

 53: typedef struct {
 54:   void (*func)(void);
 55:   void *ctx;
 56: } PetscFortranCallback;

 58: /*
 59:    All PETSc objects begin with the fields defined in PETSCHEADER.
 60:    The PetscObject is a way of examining these fields regardless of
 61:    the specific object. In C++ this could be a base abstract class
 62:    from which all objects are derived.
 63: */
 64: #define PETSC_MAX_OPTIONS_HANDLER 5
 65: typedef struct _p_PetscObject {
 66:   PetscClassId         classid;
 67:   PetscOps             bops[1];
 68:   MPI_Comm             comm;
 69:   PetscInt             type;
 70:   PetscLogDouble       flops,time,mem,memchildren;
 71:   PetscObjectId        id;
 72:   PetscInt             refct;
 73:   PetscMPIInt          tag;
 74:   PetscFunctionList    qlist;
 75:   PetscObjectList      olist;
 76:   char                 *class_name;    /*  for example, "Vec" */
 77:   char                 *description;
 78:   char                 *mansec;
 79:   char                 *type_name;     /*  this is the subclass, for example VECSEQ which equals "seq" */
 80:   PetscObject          parent;
 81:   PetscObjectId        parentid;
 82:   char*                name;
 83:   char                 *prefix;
 84:   PetscInt             tablevel;
 85:   void                 *cpp;
 86:   PetscObjectState     state;
 87:   PetscInt             int_idmax,        intstar_idmax;
 88:   PetscObjectState     *intcomposedstate,*intstarcomposedstate;
 89:   PetscInt             *intcomposeddata, **intstarcomposeddata;
 90:   PetscInt             real_idmax,        realstar_idmax;
 91:   PetscObjectState     *realcomposedstate,*realstarcomposedstate;
 92:   PetscReal            *realcomposeddata, **realstarcomposeddata;
 93:   PetscInt             scalar_idmax,        scalarstar_idmax;
 94:   PetscObjectState     *scalarcomposedstate,*scalarstarcomposedstate;
 95:   PetscScalar          *scalarcomposeddata, **scalarstarcomposeddata;
 96:   void                 (**fortran_func_pointers)(void);                  /* used by Fortran interface functions to stash user provided Fortran functions */
 97:   PetscInt             num_fortran_func_pointers;                        /* number of Fortran function pointers allocated */
 98:   PetscFortranCallback *fortrancallback[PETSC_FORTRAN_CALLBACK_MAXTYPE];
 99:   PetscInt             num_fortrancallback[PETSC_FORTRAN_CALLBACK_MAXTYPE];
100:   void                 *python_context;
101:   PetscErrorCode       (*python_destroy)(void*);

103:   PetscInt             noptionhandler;
104:   PetscErrorCode       (*optionhandler[PETSC_MAX_OPTIONS_HANDLER])(PetscObject,void*);
105:   PetscErrorCode       (*optiondestroy[PETSC_MAX_OPTIONS_HANDLER])(PetscObject,void*);
106:   void                 *optionctx[PETSC_MAX_OPTIONS_HANDLER];
107:   PetscPrecision       precision;
108:   PetscBool            optionsprinted;
109: #if defined(PETSC_HAVE_SAWS)
110:   PetscBool            amsmem;          /* if PETSC_TRUE then this object is registered with SAWs and visible to clients */
111:   PetscBool            amspublishblock; /* if PETSC_TRUE and publishing objects then will block at PetscObjectSAWsBlock() */
112: #endif
113: } _p_PetscObject;

115: #define PETSCHEADER(ObjectOps) \
116:   _p_PetscObject hdr;          \
117:   ObjectOps      ops[1]

119: #define  PETSCFREEDHEADER -1

121: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*PetscObjectDestroyFunction)(PetscObject*); /* force cast in next macro to NEVER use extern "C" style */
122: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*PetscObjectViewFunction)(PetscObject,PetscViewer);

124: /*@C
125:     PetscHeaderCreate - Creates a PETSc object of a particular class

127:     Input Parameters:
128: +   classid - the classid associated with this object (for example VEC_CLASSID)
129: .   class_name - string name of class; should be static (for example "Vec")
130: .   descr - string containing short description; should be static (for example "Vector")
131: .   mansec - string indicating section in manual pages; should be static (for example "Vec")
132: .   comm - the MPI Communicator
133: .   destroy - the destroy routine for this object (for example VecDestroy())
134: -   view - the view routine for this object (for example VecView())

136:     Output Parameter:
137: .   h - the newly created object

139:     Level: developer

141: .seealso: PetscHeaderDestroy(), PetscClassIdRegister()

143: @*/
144: #define PetscHeaderCreate(h,classid,class_name,descr,mansec,comm,destroy,view) \
145:   (PetscNew(&(h)) || \
146:    PetscHeaderCreate_Private((PetscObject)h,classid,class_name,descr,mansec,comm,(PetscObjectDestroyFunction)destroy,(PetscObjectViewFunction)view) || \
147:    PetscLogObjectCreate(h) || \
148:    PetscLogObjectMemory((PetscObject)h,sizeof(*(h))))

150: PETSC_EXTERN PetscErrorCode PetscComposedQuantitiesDestroy(PetscObject obj);
151: PETSC_EXTERN PetscErrorCode PetscHeaderCreate_Private(PetscObject,PetscClassId,const char[],const char[],const char[],MPI_Comm,PetscObjectDestroyFunction,PetscObjectViewFunction);

153: /*@C
154:     PetscHeaderDestroy - Final step in destroying a PetscObject

156:     Input Parameters:
157: .   h - the header created with PetscHeaderCreate()

159:     Level: developer

161: .seealso: PetscHeaderCreate()
162: @*/
163: #define PetscHeaderDestroy(h) (PetscHeaderDestroy_Private((PetscObject)(*h)) || PetscFree(*h))

165: PETSC_EXTERN PetscErrorCode PetscHeaderDestroy_Private(PetscObject);
166: PETSC_EXTERN PetscErrorCode PetscObjectCopyFortranFunctionPointers(PetscObject,PetscObject);
167: PETSC_EXTERN PetscErrorCode PetscObjectSetFortranCallback(PetscObject,PetscFortranCallbackType,PetscFortranCallbackId*,void(*)(void),void *ctx);
168: PETSC_EXTERN PetscErrorCode PetscObjectGetFortranCallback(PetscObject,PetscFortranCallbackType,PetscFortranCallbackId,void(**)(void),void **ctx);

170: PETSC_INTERN PetscErrorCode PetscCitationsInitialize(void);
171: PETSC_INTERN PetscErrorCode PetscOptionsFindPair_Private(const char[],const char[],char**,PetscBool*);


175: /*
176:     Macros to test if a PETSc object is valid and if pointers are valid
177: */
178: #if !defined(PETSC_USE_DEBUG)


189: #else

192:   do {                                                                  \
193:     if (!h) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"Null Object: Parameter # %d",arg); \
195:     if (((PetscObject)(h))->classid != ck) {                            \
196:       if (((PetscObject)(h))->classid == PETSCFREEDHEADER) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"Object already free: Parameter # %d",arg); \
197:       else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Wrong type of object: Parameter # %d",arg); \
198:     }                                                                   \
199:   } while (0)

202:   do {                                                                  \
203:     if (!h) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"Null Object: Parameter # %d",arg); \
205:     if (((PetscObject)(h))->classid == PETSCFREEDHEADER) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"Object already free: Parameter # %d",arg); \
206:     else if (((PetscObject)(h))->classid < PETSC_SMALLEST_CLASSID || ((PetscObject)(h))->classid > PETSC_LARGEST_CLASSID) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"Invalid type of object: Parameter # %d",arg); \
207:   } while (0)

210:   do {                                                                  \
211:     if (!h) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"Null Pointer: Parameter # %d",arg); \
213:   } while (0)

216:   do {                                                                  \
217:     if (!h) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"Null Pointer: Parameter # %d",arg);\
219:   } while (0)

222:   do {                                                                  \
223:     if (!h) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_BADPTR,"Null Pointer: Parameter # %d",arg); \
225:   } while (0)

228:   do {                                                                  \
229:     if (!h) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"Null Pointer: Parameter # %d",arg); \
231:   } while (0)

234:   do {                                                                  \
235:     if (!h) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"Null Pointer: Parameter # %d",arg); \
237:   } while (0)

240:   do {                                                                  \
241:     if (!f) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"Null Function Pointer: Parameter # %d",arg); \
242:   } while (0)

244: #endif

246: #if !defined(PETSC_USE_DEBUG)


258: #else

260: /*
261:     For example, in the dot product between two vectors,
262:   both vectors must be either Seq or MPI, not one of each
263: */
265:   if (((PetscObject)a)->type != ((PetscObject)b)->type) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Objects not of same type: Argument # %d and %d",arga,argb);
266: /*
267:    Use this macro to check if the type is set
268: */
270:   if (!((PetscObject)a)->type_name) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"%s object's type is not set: Argument # %d",((PetscObject)a)->class_name,arg);
271: /*
272:    Sometimes object must live on same communicator to inter-operate
273: */
275:   do {                                                                  \
276:     PetscErrorCode _6_ierr,__flag;                                      \
277:     _6_MPI_Comm_compare(PetscObjectComm((PetscObject)a),PetscObjectComm((PetscObject)b),&__flag);CHKERRQ(_6_ierr);                                                   \
278:     if (__flag != MPI_CONGRUENT && __flag != MPI_IDENT) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMECOMM,"Different communicators in the two objects: Argument # %d and %d flag %d",arga,argb,__flag); \
279:   } while (0)

282:   do {                                                  \
285:   } while (0)

288:   do {                                                                  \
289:     PetscErrorCode _7_ierr;                                             \
290:     PetscReal b1[2],b2[2];                                              \
291:     b1[0] = -PetscRealPart(b); b1[1] = PetscRealPart(b);                \
292:     _7_MPI_Allreduce(b1,b2,2,MPIU_REAL,MPIU_MAX,PetscObjectComm((PetscObject)a));CHKERRQ(_7_ierr); \
293:     if (-b2[0] != b2[1]) SETERRQ1(PetscObjectComm((PetscObject)a),PETSC_ERR_ARG_WRONG,"Scalar value must be same on all processes, argument # %d",c); \
294:   } while (0)

297:   do {                                                                  \
298:     PetscErrorCode _7_ierr;                                             \
299:     PetscReal b1[2],b2[2];                                              \
300:     b1[0] = -b; b1[1] = b;                                              \
301:     _7_MPI_Allreduce(b1,b2,2,MPIU_REAL,MPIU_MAX,PetscObjectComm((PetscObject)a));CHKERRQ(_7_ierr); \
302:     if (-b2[0] != b2[1]) SETERRQ1(PetscObjectComm((PetscObject)a),PETSC_ERR_ARG_WRONG,"Real value must be same on all processes, argument # %d",c); \
303:   } while (0)

306:   do {                                                                  \
307:     PetscErrorCode _7_ierr;                                             \
308:     PetscInt b1[2],b2[2];                                               \
309:     b1[0] = -b; b1[1] = b;                                              \
310:     _7_MPI_Allreduce(b1,b2,2,MPIU_INT,MPI_MAX,PetscObjectComm((PetscObject)a));CHKERRQ(_7_ierr); \
311:     if (-b2[0] != b2[1]) SETERRQ1(PetscObjectComm((PetscObject)a),PETSC_ERR_ARG_WRONG,"Int value must be same on all processes, argument # %d",c); \
312:   } while (0)

315:   do {                                                                  \
316:     PetscErrorCode _7_ierr;                                             \
317:     PetscMPIInt b1[2],b2[2];                                            \
318:     b1[0] = -(PetscMPIInt)b; b1[1] = (PetscMPIInt)b;                    \
319:     _7_MPI_Allreduce(b1,b2,2,MPI_INT,MPI_MAX,PetscObjectComm((PetscObject)a));CHKERRQ(_7_ierr); \
320:     if (-b2[0] != b2[1]) SETERRQ1(PetscObjectComm((PetscObject)a),PETSC_ERR_ARG_WRONG,"Bool value must be same on all processes, argument # %d",c); \
321:   } while (0)

324:   do {                                                                  \
325:     PetscErrorCode _7_ierr;                                             \
326:     PetscMPIInt b1[2],b2[2];                                            \
327:     b1[0] = -(PetscMPIInt)b; b1[1] = (PetscMPIInt)b;                    \
328:     _7_MPI_Allreduce(b1,b2,2,MPI_INT,MPI_MAX,PetscObjectComm((PetscObject)a));CHKERRQ(_7_ierr); \
329:     if (-b2[0] != b2[1]) SETERRQ1(PetscObjectComm((PetscObject)a),PETSC_ERR_ARG_WRONG,"Enum value must be same on all processes, argument # %d",c); \
330:   } while (0)

332: #endif

334: /*
335:    PetscTryMethod - Queries an object for a method, if it exists then calls it.
336:               These are intended to be used only inside PETSc functions.

338:    Level: developer

340: .seealso: PetscUseMethod()
341: */
342: #define  PetscTryMethod(obj,A,B,C) \
343:   0;{ PetscErrorCode (*f)B, __ierr; \
344:     __PetscObjectQueryFunction((PetscObject)obj,A,&f);CHKERRQ(__ierr); \
345:     if (f) {__(*f)C;CHKERRQ(__ierr);}\
346:   }

348: /*
349:    PetscUseMethod - Queries an object for a method, if it exists then calls it, otherwise generates an error.
350:               These are intended to be used only inside PETSc functions.

352:    Level: developer

354: .seealso: PetscTryMethod()
355: */
356: #define  PetscUseMethod(obj,A,B,C) \
357:   0;{ PetscErrorCode (*f)B, __ierr; \
358:     __PetscObjectQueryFunction((PetscObject)obj,A,&f);CHKERRQ(__ierr); \
359:     if (f) {__(*f)C;CHKERRQ(__ierr);}\
360:     else SETERRQ1(PetscObjectComm((PetscObject)obj),PETSC_ERR_SUP,"Cannot locate function %s in object",A); \
361:   }

363: /*MC
364:    PetscObjectStateIncrease - Increases the state of any PetscObject

366:    Synopsis:
367:    #include "petsc/private/petscimpl.h"
368:    PetscErrorCode PetscObjectStateIncrease(PetscObject obj)

370:    Logically Collective

372:    Input Parameter:
373: .  obj - any PETSc object, for example a Vec, Mat or KSP. This must be
374:          cast with a (PetscObject), for example,
375:          PetscObjectStateIncrease((PetscObject)mat);

377:    Notes: object state is an integer which gets increased every time
378:    the object is changed internally. By saving and later querying the object state
379:    one can determine whether information about the object is still current.
380:    Currently, state is maintained for Vec and Mat objects.

382:    This routine is mostly for internal use by PETSc; a developer need only
383:    call it after explicit access to an object's internals. Routines such
384:    as VecSet() or MatScale() already call this routine. It is also called, as a
385:    precaution, in VecRestoreArray(), MatRestoreRow(), MatDenseRestoreArray().

387:    This routine is logically collective because state equality comparison needs to be possible without communication.

389:    Level: developer

391:    seealso: PetscObjectStateGet()

393:    Concepts: state

395: M*/
396: #define PetscObjectStateIncrease(obj) ((obj)->state++,0)

398: PETSC_EXTERN PetscErrorCode PetscObjectStateGet(PetscObject,PetscObjectState*);
399: PETSC_EXTERN PetscErrorCode PetscObjectStateSet(PetscObject,PetscObjectState);
400: PETSC_EXTERN PetscErrorCode PetscObjectComposedDataRegister(PetscInt*);
401: PETSC_EXTERN PetscErrorCode PetscObjectComposedDataIncreaseInt(PetscObject);
402: PETSC_EXTERN PetscErrorCode PetscObjectComposedDataIncreaseIntstar(PetscObject);
403: PETSC_EXTERN PetscErrorCode PetscObjectComposedDataIncreaseReal(PetscObject);
404: PETSC_EXTERN PetscErrorCode PetscObjectComposedDataIncreaseRealstar(PetscObject);
405: PETSC_EXTERN PetscErrorCode PetscObjectComposedDataIncreaseScalar(PetscObject);
406: PETSC_EXTERN PetscErrorCode PetscObjectComposedDataIncreaseScalarstar(PetscObject);
407: PETSC_EXTERN PetscInt         PetscObjectComposedDataMax;
408: /*MC
409:    PetscObjectComposedDataSetInt - attach integer data to a PetscObject

411:    Synopsis:
412:    #include "petsc/private/petscimpl.h"
413:    PetscErrorCode PetscObjectComposedDataSetInt(PetscObject obj,int id,int data)

415:    Not collective

417:    Input parameters:
418: +  obj - the object to which data is to be attached
419: .  id - the identifier for the data
420: -  data - the data to  be attached

422:    Notes
423:    The data identifier can best be created through a call to  PetscObjectComposedDataRegister()

425:    Level: developer
426: M*/
427: #define PetscObjectComposedDataSetInt(obj,id,data)                                      \
428:   ((((obj)->int_idmax < PetscObjectComposedDataMax) && PetscObjectComposedDataIncreaseInt(obj)) ||  \
429:    ((obj)->intcomposeddata[id] = data,(obj)->intcomposedstate[id] = (obj)->state, 0))

431: /*MC
432:    PetscObjectComposedDataGetInt - retrieve integer data attached to an object

434:    Synopsis:
435:    #include "petsc/private/petscimpl.h"
436:    PetscErrorCode PetscObjectComposedDataGetInt(PetscObject obj,int id,int data,PetscBool  flag)

438:    Not collective

440:    Input parameters:
441: +  obj - the object from which data is to be retrieved
442: -  id - the identifier for the data

444:    Output parameters
445: +  data - the data to be retrieved
446: -  flag - PETSC_TRUE if the data item exists and is valid, PETSC_FALSE otherwise

448:    The 'data' and 'flag' variables are inlined, so they are not pointers.

450:    Level: developer
451: M*/
452: #define PetscObjectComposedDataGetInt(obj,id,data,flag)                            \
453:   ((((obj)->intcomposedstate && ((obj)->intcomposedstate[id] == (obj)->state)) ?   \
454:    (data = (obj)->intcomposeddata[id],flag = PETSC_TRUE) : (flag = PETSC_FALSE)),0)

456: /*MC
457:    PetscObjectComposedDataSetIntstar - attach integer array data to a PetscObject

459:    Synopsis:
460:    #include "petsc/private/petscimpl.h"
461:    PetscErrorCode PetscObjectComposedDataSetIntstar(PetscObject obj,int id,int *data)

463:    Not collective

465:    Input parameters:
466: +  obj - the object to which data is to be attached
467: .  id - the identifier for the data
468: -  data - the data to  be attached

470:    Notes
471:    The data identifier can best be determined through a call to
472:    PetscObjectComposedDataRegister()

474:    Level: developer
475: M*/
476: #define PetscObjectComposedDataSetIntstar(obj,id,data)                                          \
477:   ((((obj)->intstar_idmax < PetscObjectComposedDataMax) && PetscObjectComposedDataIncreaseIntstar(obj)) ||  \
478:    ((obj)->intstarcomposeddata[id] = data,(obj)->intstarcomposedstate[id] = (obj)->state, 0))

480: /*MC
481:    PetscObjectComposedDataGetIntstar - retrieve integer array data
482:    attached to an object

484:    Synopsis:
485:    #include "petsc/private/petscimpl.h"
486:    PetscErrorCode PetscObjectComposedDataGetIntstar(PetscObject obj,int id,int *data,PetscBool  flag)

488:    Not collective

490:    Input parameters:
491: +  obj - the object from which data is to be retrieved
492: -  id - the identifier for the data

494:    Output parameters
495: +  data - the data to be retrieved
496: -  flag - PETSC_TRUE if the data item exists and is valid, PETSC_FALSE otherwise

498:    The 'data' and 'flag' variables are inlined, so they are not pointers.

500:    Level: developer
501: M*/
502: #define PetscObjectComposedDataGetIntstar(obj,id,data,flag)                               \
503:   ((((obj)->intstarcomposedstate && ((obj)->intstarcomposedstate[id] == (obj)->state)) ?  \
504:    (data = (obj)->intstarcomposeddata[id],flag = PETSC_TRUE) : (flag = PETSC_FALSE)),0)

506: /*MC
507:    PetscObjectComposedDataSetReal - attach real data to a PetscObject

509:    Synopsis:
510:    #include "petsc/private/petscimpl.h"
511:    PetscErrorCode PetscObjectComposedDataSetReal(PetscObject obj,int id,PetscReal data)

513:    Not collective

515:    Input parameters:
516: +  obj - the object to which data is to be attached
517: .  id - the identifier for the data
518: -  data - the data to  be attached

520:    Notes
521:    The data identifier can best be determined through a call to
522:    PetscObjectComposedDataRegister()

524:    Level: developer
525: M*/
526: #define PetscObjectComposedDataSetReal(obj,id,data)                                       \
527:   ((((obj)->real_idmax < PetscObjectComposedDataMax) && PetscObjectComposedDataIncreaseReal(obj)) ||  \
528:    ((obj)->realcomposeddata[id] = data,(obj)->realcomposedstate[id] = (obj)->state, 0))

530: /*MC
531:    PetscObjectComposedDataGetReal - retrieve real data attached to an object

533:    Synopsis:
534:    #include "petsc/private/petscimpl.h"
535:    PetscErrorCode PetscObjectComposedDataGetReal(PetscObject obj,int id,PetscReal data,PetscBool  flag)

537:    Not collective

539:    Input parameters:
540: +  obj - the object from which data is to be retrieved
541: -  id - the identifier for the data

543:    Output parameters
544: +  data - the data to be retrieved
545: -  flag - PETSC_TRUE if the data item exists and is valid, PETSC_FALSE otherwise

547:    The 'data' and 'flag' variables are inlined, so they are not pointers.

549:    Level: developer
550: M*/
551: #define PetscObjectComposedDataGetReal(obj,id,data,flag)                            \
552:   ((((obj)->realcomposedstate && ((obj)->realcomposedstate[id] == (obj)->state)) ?  \
553:    (data = (obj)->realcomposeddata[id],flag = PETSC_TRUE) : (flag = PETSC_FALSE)),0)

555: /*MC
556:    PetscObjectComposedDataSetRealstar - attach real array data to a PetscObject

558:    Synopsis:
559:    #include "petsc/private/petscimpl.h"
560:    PetscErrorCode PetscObjectComposedDataSetRealstar(PetscObject obj,int id,PetscReal *data)

562:    Not collective

564:    Input parameters:
565: +  obj - the object to which data is to be attached
566: .  id - the identifier for the data
567: -  data - the data to  be attached

569:    Notes
570:    The data identifier can best be determined through a call to
571:    PetscObjectComposedDataRegister()

573:    Level: developer
574: M*/
575: #define PetscObjectComposedDataSetRealstar(obj,id,data)                                           \
576:   ((((obj)->realstar_idmax < PetscObjectComposedDataMax) && PetscObjectComposedDataIncreaseRealstar(obj)) ||  \
577:    ((obj)->realstarcomposeddata[id] = data, (obj)->realstarcomposedstate[id] = (obj)->state, 0))

579: /*MC
580:    PetscObjectComposedDataGetRealstar - retrieve real array data
581:    attached to an object

583:    Synopsis:
584:    #include "petsc/private/petscimpl.h"
585:    PetscErrorCode PetscObjectComposedDataGetRealstar(PetscObject obj,int id,PetscReal *data,PetscBool  flag)

587:    Not collective

589:    Input parameters:
590: +  obj - the object from which data is to be retrieved
591: -  id - the identifier for the data

593:    Output parameters
594: +  data - the data to be retrieved
595: -  flag - PETSC_TRUE if the data item exists and is valid, PETSC_FALSE otherwise

597:    The 'data' and 'flag' variables are inlined, so they are not pointers.

599:    Level: developer
600: M*/
601: #define PetscObjectComposedDataGetRealstar(obj,id,data,flag)                                \
602:   ((((obj)->realstarcomposedstate && ((obj)->realstarcomposedstate[id] == (obj)->state)) ?  \
603:    (data = (obj)->realstarcomposeddata[id],flag = PETSC_TRUE) : (flag = PETSC_FALSE)),0)

605: /*MC
606:    PetscObjectComposedDataSetScalar - attach scalar data to a PetscObject

608:    Synopsis:
609:    #include "petsc/private/petscimpl.h"
610:    PetscErrorCode PetscObjectComposedDataSetScalar(PetscObject obj,int id,PetscScalar data)

612:    Not collective

614:    Input parameters:
615: +  obj - the object to which data is to be attached
616: .  id - the identifier for the data
617: -  data - the data to  be attached

619:    Notes
620:    The data identifier can best be determined through a call to
621:    PetscObjectComposedDataRegister()

623:    Level: developer
624: M*/
625: #if defined(PETSC_USE_COMPLEX)
626: #define PetscObjectComposedDataSetScalar(obj,id,data)                                        \
627:   ((((obj)->scalar_idmax < PetscObjectComposedDataMax) && PetscObjectComposedDataIncreaseScalar(obj)) || \
628:    ((obj)->scalarcomposeddata[id] = data,(obj)->scalarcomposedstate[id] = (obj)->state, 0))
629: #else
630: #define PetscObjectComposedDataSetScalar(obj,id,data) \
631:         PetscObjectComposedDataSetReal(obj,id,data)
632: #endif
633: /*MC
634:    PetscObjectComposedDataGetScalar - retrieve scalar data attached to an object

636:    Synopsis:
637:    #include "petsc/private/petscimpl.h"
638:    PetscErrorCode PetscObjectComposedDataGetScalar(PetscObject obj,int id,PetscScalar data,PetscBool  flag)

640:    Not collective

642:    Input parameters:
643: +  obj - the object from which data is to be retrieved
644: -  id - the identifier for the data

646:    Output parameters
647: +  data - the data to be retrieved
648: -  flag - PETSC_TRUE if the data item exists and is valid, PETSC_FALSE otherwise

650:    The 'data' and 'flag' variables are inlined, so they are not pointers.

652:    Level: developer
653: M*/
654: #if defined(PETSC_USE_COMPLEX)
655: #define PetscObjectComposedDataGetScalar(obj,id,data,flag)                              \
656:   ((((obj)->scalarcomposedstate && ((obj)->scalarcomposedstate[id] == (obj)->state) ) ? \
657:    (data = (obj)->scalarcomposeddata[id],flag = PETSC_TRUE) : (flag = PETSC_FALSE)),0)
658: #else
659: #define PetscObjectComposedDataGetScalar(obj,id,data,flag)                             \
660:         PetscObjectComposedDataGetReal(obj,id,data,flag)
661: #endif

663: /*MC
664:    PetscObjectComposedDataSetScalarstar - attach scalar array data to a PetscObject

666:    Synopsis:
667:    #include "petsc/private/petscimpl.h"
668:    PetscErrorCode PetscObjectComposedDataSetScalarstar(PetscObject obj,int id,PetscScalar *data)

670:    Not collective

672:    Input parameters:
673: +  obj - the object to which data is to be attached
674: .  id - the identifier for the data
675: -  data - the data to  be attached

677:    Notes
678:    The data identifier can best be determined through a call to
679:    PetscObjectComposedDataRegister()

681:    Level: developer
682: M*/
683: #if defined(PETSC_USE_COMPLEX)
684: #define PetscObjectComposedDataSetScalarstar(obj,id,data)                                             \
685:   ((((obj)->scalarstar_idmax < PetscObjectComposedDataMax) && PetscObjectComposedDataIncreaseScalarstar(obj)) ||  \
686:    ((obj)->scalarstarcomposeddata[id] = data,(obj)->scalarstarcomposedstate[id] = (obj)->state, 0))
687: #else
688: #define PetscObjectComposedDataSetScalarstar(obj,id,data) \
689:         PetscObjectComposedDataSetRealstar(obj,id,data)
690: #endif
691: /*MC
692:    PetscObjectComposedDataGetScalarstar - retrieve scalar array data
693:    attached to an object

695:    Synopsis:
696:    #include "petsc/private/petscimpl.h"
697:    PetscErrorCode PetscObjectComposedDataGetScalarstar(PetscObject obj,int id,PetscScalar *data,PetscBool  flag)

699:    Not collective

701:    Input parameters:
702: +  obj - the object from which data is to be retrieved
703: -  id - the identifier for the data

705:    Output parameters
706: +  data - the data to be retrieved
707: -  flag - PETSC_TRUE if the data item exists and is valid, PETSC_FALSE otherwise

709:    The 'data' and 'flag' variables are inlined, so they are not pointers.

711:    Level: developer
712: M*/
713: #if defined(PETSC_USE_COMPLEX)
714: #define PetscObjectComposedDataGetScalarstar(obj,id,data,flag)                                 \
715:   ((((obj)->scalarstarcomposedstate && ((obj)->scalarstarcomposedstate[id] == (obj)->state)) ? \
716:        (data = (obj)->scalarstarcomposeddata[id],flag = PETSC_TRUE) : (flag = PETSC_FALSE)),0)
717: #else
718: #define PetscObjectComposedDataGetScalarstar(obj,id,data,flag)         \
719:         PetscObjectComposedDataGetRealstar(obj,id,data,flag)
720: #endif

722: PETSC_EXTERN PetscErrorCode PetscObjectGetId(PetscObject,PetscObjectId*);

724: PETSC_EXTERN PetscMPIInt Petsc_Counter_keyval;
725: PETSC_EXTERN PetscMPIInt Petsc_InnerComm_keyval;
726: PETSC_EXTERN PetscMPIInt Petsc_OuterComm_keyval;

728: /*
729:   PETSc communicators have this attribute, see
730:   PetscCommDuplicate(), PetscCommDestroy(), PetscCommGetNewTag(), PetscObjectGetName()
731: */
732: typedef struct {
733:   PetscMPIInt tag;              /* next free tag value */
734:   PetscInt    refcount;         /* number of references, communicator can be freed when this reaches 0 */
735:   PetscInt    namecount;        /* used to generate the next name, as in Vec_0, Mat_1, ... */
736: } PetscCommCounter;

738: #if defined(PETSC_HAVE_CUSP)
739: /*E
740:     PetscCUSPFlag - indicates which memory (CPU, GPU, or none contains valid vector

742:    PETSC_CUSP_UNALLOCATED  - no memory contains valid matrix entries; NEVER used for vectors
743:    PETSC_CUSP_GPU - GPU has valid vector/matrix entries
744:    PETSC_CUSP_CPU - CPU has valid vector/matrix entries
745:    PETSC_CUSP_BOTH - Both GPU and CPU have valid vector/matrix entries and they match

747:    Level: developer
748: E*/
749: typedef enum {PETSC_CUSP_UNALLOCATED,PETSC_CUSP_GPU,PETSC_CUSP_CPU,PETSC_CUSP_BOTH} PetscCUSPFlag;
750: #endif

752: #if defined(PETSC_HAVE_VIENNACL)
753: /*E
754:     PetscViennaCLFlag - indicates which memory (CPU, GPU, or none contains valid vector

756:    PETSC_VIENNACL_UNALLOCATED  - no memory contains valid matrix entries; NEVER used for vectors
757:    PETSC_VIENNACL_GPU - GPU has valid vector/matrix entries
758:    PETSC_VIENNACL_CPU - CPU has valid vector/matrix entries
759:    PETSC_VIENNACL_BOTH - Both GPU and CPU have valid vector/matrix entries and they match

761:    Level: developer
762: E*/
763: typedef enum {PETSC_VIENNACL_UNALLOCATED,PETSC_VIENNACL_GPU,PETSC_VIENNACL_CPU,PETSC_VIENNACL_BOTH} PetscViennaCLFlag;
764: #endif

766: typedef enum {STATE_BEGIN, STATE_PENDING, STATE_END} SRState;

768: #define REDUCE_SUM  0
769: #define REDUCE_MAX  1
770: #define REDUCE_MIN  2

772: typedef struct {
773:   MPI_Comm    comm;
774:   MPI_Request request;
775:   PetscBool   async;
776:   PetscScalar *lvalues;     /* this are the reduced values before call to MPI_Allreduce() */
777:   PetscScalar *gvalues;     /* values after call to MPI_Allreduce() */
778:   void        **invecs;     /* for debugging only, vector/memory used with each op */
779:   PetscInt    *reducetype;  /* is particular value to be summed or maxed? */
780:   SRState     state;        /* are we calling xxxBegin() or xxxEnd()? */
781:   PetscInt    maxops;       /* total amount of space we have for requests */
782:   PetscInt    numopsbegin;  /* number of requests that have been queued in */
783:   PetscInt    numopsend;    /* number of requests that have been gotten by user */
784: } PetscSplitReduction;

786: PETSC_EXTERN PetscErrorCode PetscSplitReductionGet(MPI_Comm,PetscSplitReduction**);
787: PETSC_EXTERN PetscErrorCode PetscSplitReductionEnd(PetscSplitReduction*);
788: PETSC_EXTERN PetscErrorCode PetscSplitReductionExtend(PetscSplitReduction*);

790: #endif /* _PETSCHEAD_H */