Actual source code: petsc.h

  1: /*
  2:    This is the main PETSc include file (for C and C++).  It is included by all
  3:    other PETSc include files, so it almost never has to be specifically included.
  4: */
  7: /* ========================================================================== */
  8: /* 
  9:    petscconf.h is contained in ${PETSC_ARCH}/include/petscconf.h it is 
 10:    found automatically by the compiler due to the -I${PETSC_DIR}/${PETSC_ARCH}/include
 11:    in the conf/variables definition of PETSC_INCLUDE
 12: */
 13: #include "petscconf.h"
 14: #include "petscfix.h"

 16: /* ========================================================================== */
 17: /* 
 18:    This facilitates using C version of PETSc from C++ and 
 19:    C++ version from C. Use --with-c-support --with-clanguage=c++ with config/configure.py for the latter)
 20: */
 22: #error "PETSc configured with --with-clanguage=c++ and NOT --with-c-support - it can be used only with a C++ compiler"
 23: #endif      

 28: #else
 31: #endif
 32: /* ========================================================================== */
 33: /* 
 34:    Current PETSc version number and release date. Also listed in
 35:     Web page
 36:     src/docs/tex/manual/intro.tex,
 37:     src/docs/tex/manual/manual.tex.
 38:     src/docs/website/index.html.
 39: */
 40:  #include petscversion.h
 41: #define PETSC_AUTHOR_INFO        "\
 42:        The PETSc Team\n\
 43:     petsc-maint@mcs.anl.gov\n\
 44:  http://www.mcs.anl.gov/petsc/\n"
 45: #if (PETSC_VERSION_RELEASE == 1)
 46: #define PetscGetVersion(version,len) (PetscSNPrintf(version,len,"Petsc Release Version %d.%d.%d, Patch %d, ", \
 47:                                          PETSC_VERSION_MAJOR,PETSC_VERSION_MINOR, PETSC_VERSION_SUBMINOR, \
 48:                                          PETSC_VERSION_PATCH),PetscStrcat(version,PETSC_VERSION_PATCH_DATE))
 49: #else
 50: #define PetscGetVersion(version,len) (PetscSNPrintf(version,len,"Petsc Development"), \
 51:                                          PetscStrcat(version," HG revision: "),PetscStrcat(version,PETSC_VERSION_HG), \
 52:                                          PetscStrcat(version," HG Date: "),PetscStrcat(version,PETSC_VERSION_DATE_HG))
 53: #endif

 55: /*MC
 56:     PetscGetVersion - Gets the PETSc version information in a string.

 58:     Input Parameter:
 59: .   len - length of the string

 61:     Output Parameter:
 62: .   version - version string

 64:     Level: developer

 66:     Usage:
 67:     char version[256];
 68:     PetscGetVersion(version,256);CHKERRQ(ierr)

 70:     Fortran Note:
 71:     This routine is not supported in Fortran.

 73: .seealso: PetscGetProgramName()

 75: M*/

 77: /* ========================================================================== */

 79: /*
 80:    Currently cannot check formatting for PETSc print statements because we have our
 81:    own format %D and %G
 82: */
 83: #undef  PETSC_PRINTF_FORMAT_CHECK
 84: #define PETSC_PRINTF_FORMAT_CHECK(a,b)
 85: #undef  PETSC_FPRINTF_FORMAT_CHECK
 86: #define PETSC_FPRINTF_FORMAT_CHECK(a,b)

 88: /*
 89:    Fixes for config/configure.py time choices which impact our interface. Currently only
 90:    calling conventions and extra compiler checking falls under this category.
 91: */
 92: #if !defined(PETSC_STDCALL)
 93: #define PETSC_STDCALL
 94: #endif
 95: #if !defined(PETSC_TEMPLATE)
 96: #define PETSC_TEMPLATE
 97: #endif
 98: #if !defined(PETSC_HAVE_DLL_EXPORT)
 99: #define PETSC_DLL_EXPORT
100: #define PETSC_DLL_IMPORT
101: #endif
102: #if !defined()
103: #define 
104: #endif
105: #if !defined()
106: #define 
107: #endif
108: #if !defined()
109: #define 
110: #endif
111: #if !defined()
112: #define 
113: #endif
114: #if !defined()
115: #define 
116: #endif
117: #if !defined()
118: #define 
119: #endif
120: #if !defined()
121: #define 
122: #endif
123: #if !defined()
124: #define 
125: #endif
126: /* ========================================================================== */

128: /*
129:     Defines the interface to MPI allowing the use of all MPI functions.

131:     PETSc does not use the C++ binding of MPI at ALL. The following flag
132:     makes sure the C++ bindings are not included. The C++ bindings REQUIRE
133:     putting mpi.h before ANY C++ include files, we cannot control this
134:     with all PETSc users. Users who want to use the MPI C++ bindings can include 
135:     mpicxx.h directly in their code
136: */
137: #define MPICH_SKIP_MPICXX 1
138: #define OMPI_SKIP_MPICXX 1
139: #include "mpi.h"
140: /*
141:     Yuck, we need to put stdio.h AFTER mpi.h for MPICH2 with C++ compiler 
142:     see the top of mpicxx.h in the MPICH2 distribution.

144:     The MPI STANDARD HAS TO BE CHANGED to prevent this nonsense.
145: */
146: #include <stdio.h>

148: /* MSMPI on 32bit windows requires this yukky hack - that breaks MPI standard compliance */
149: #if !defined(MPIAPI)
150: #define MPIAPI
151: #endif

153: /*MC
154:     PetscErrorCode - datatype used for return error code from all PETSc functions

156:     Level: beginner

158: .seealso: CHKERRQ, SETERRQ
159: M*/
160: typedef int PetscErrorCode;

162: /*MC

164:     PetscCookie - A unique id used to identify each PETSc object.
165:          (internal integer in the data structure used for error
166:          checking). These are all defined by an offset from the lowest
167:          one, PETSC_SMALLEST_COOKIE. 

169:     Level: advanced

171: .seealso: PetscCookieRegister(), PetscLogEventRegister(), PetscHeaderCreate()
172: M*/
173: typedef int PetscCookie;

175: /*MC
176:     PetscLogEvent - id used to identify PETSc or user events - primarily for logging

178:     Level: intermediate

180: .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscLogStage
181: M*/
182: typedef int PetscLogEvent;

184: /*MC
185:     PetscLogStage - id used to identify user stages of runs - for logging

187:     Level: intermediate

189: .seealso: PetscLogStageRegister(), PetscLogStageBegin(), PetscLogStageEnd(), PetscLogEvent
190: M*/
191: typedef int PetscLogStage;

193: /*MC
194:     PetscBLASInt - datatype used to represent 'int' parameters to BLAS/LAPACK functions.

196:     Level: intermediate

198:     Notes: usually this is the same as PetscInt, but if PETSc was built with --with-64-bit-indices but 
199:            standard C/Fortran integers are 32 bit then this is NOT the same as PetscInt

201: .seealso: PetscMPIInt, PetscInt

203: M*/
204: typedef int PetscBLASInt;

206: /*MC
207:     PetscMPIInt - datatype used to represent 'int' parameters to MPI functions.

209:     Level: intermediate

211:     Notes: usually this is the same as PetscInt, but if PETSc was built with --with-64-bit-indices but 
212:            standard C/Fortran integers are 32 bit then this is NOT the same as PetscInt

214:     PetscBLASIntCheck(a) checks if the given PetscInt a will fit in a PetscBLASInt, if not it generates a 
215:       PETSC_ERR_ARG_OUTOFRANGE.

217:     PetscBLASInt b = PetscBLASIntCast(a) checks if the given PetscInt a will fit in a PetscBLASInt, if not it 
218:       generates a PETSC_ERR_ARG_OUTOFRANGE

220: .seealso: PetscBLASInt, PetscInt

222: M*/
223: typedef int PetscMPIInt;

225: /*MC
226:     PetscEnum - datatype used to pass enum types within PETSc functions.

228:     Level: intermediate

230:     PetscMPIIntCheck(a) checks if the given PetscInt a will fit in a PetscMPIInt, if not it generates a 
231:       PETSC_ERR_ARG_OUTOFRANGE.

233:     PetscMPIInt b = PetscMPIIntCast(a) checks if the given PetscInt a will fit in a PetscMPIInt, if not it 
234:       generates a PETSC_ERR_ARG_OUTOFRANGE

236: .seealso: PetscOptionsGetEnum(), PetscOptionsEnum(), PetscBagRegisterEnum()
237: M*/
238: typedef enum { ENUM_DUMMY } PetscEnum;

240: /*MC
241:     PetscInt - PETSc type that represents integer - used primarily to
242:       represent size of objects. Its size can be configured with the option
243:       --with-64-bit-indices - to be either 32bit or 64bit [default 32 bit ints]

245:    Level: intermediate


248: .seealso: PetscScalar, PetscBLASInt, PetscMPIInt
249: M*/
250: #if defined(PETSC_USE_64BIT_INDICES)
251: typedef long long PetscInt;
252: #define MPIU_INT MPI_LONG_LONG_INT
253: #else
254: typedef int PetscInt;
255: #define MPIU_INT MPI_INT
256: #endif

258: /* add in MPIU type for size_t */
259: #if (PETSC_SIZEOF_SIZE_T) == (PETSC_SIZEOF_INT)
260: #define MPIU_SIZE_T MPI_INT
261: #elif  (PETSC_SIZEOF_SIZE_T) == (PETSC_SIZEOF_LONG)
262: #define MPIU_SIZE_T MPI_LONG
263: #elif  (PETSC_SIZEOF_SIZE_T) == (PETSC_SIZEOF_LONG_LONG)
264: #define MPIU_SIZE_T MPI_LONG_LONG_INT
265: #else
266: #error "Unknown size for size_t! Send us a bugreport at petsc-maint@mcs.anl.gov"
267: #endif


270: /*
271:       You can use PETSC_STDOUT as a replacement of stdout. You can also change
272:     the value of PETSC_STDOUT to redirect all standard output elsewhere
273: */


277: /*
278:       You can use PETSC_STDERR as a replacement of stderr. You can also change
279:     the value of PETSC_STDERR to redirect all standard error elsewhere
280: */

283: /*
284:       PETSC_ZOPEFD is used to send data to the PETSc webpage.  It can be used
285:     in conjunction with PETSC_STDOUT, or by itself.
286: */

290: /*MC
291:       PetscPolymorphicSubroutine - allows defining a C++ polymorphic version of 
292:             a PETSc function that remove certain optional arguments for a simplier user interface

294:      Not collective

296:    Synopsis:
297:    PetscPolymorphicSubroutine(Functionname,(arguments of C++ function),(arguments of C function))
298:  
299:    Level: developer

301:     Example:
302:       PetscPolymorphicSubroutine(VecNorm,(Vec x,PetscReal *r),(x,NORM_2,r)) generates the new routine
303:            PetscErrorCode VecNorm(Vec x,PetscReal *r) = VecNorm(x,NORM_2,r)

305: .seealso: PetscPolymorphicFunction()

307: M*/
308: #define PetscPolymorphicSubroutine(A,B,C) PETSC_STATIC_INLINE PetscErrorCode A B {return A C;}

310: /*MC
311:       PetscPolymorphicScalar - allows defining a C++ polymorphic version of 
312:             a PETSc function that replaces a PetscScalar * argument with a PetscScalar argument

314:      Not collective

316:    Synopsis:
317:    PetscPolymorphicScalar(Functionname,(arguments of C++ function),(arguments of C function))
318:  
319:    Level: developer

321:     Example:
322:       PetscPolymorphicScalar(VecAXPY,(PetscScalar _val,Vec x,Vec y),(&_Val,x,y)) generates the new routine
323:            PetscErrorCode VecAXPY(PetscScalar _val,Vec x,Vec y) = {PetscScalar _Val = _val; return VecAXPY(&_Val,x,y);}

325: .seealso: PetscPolymorphicFunction(),PetscPolymorphicSubroutine()

327: M*/
328: #define PetscPolymorphicScalar(A,B,C) PETSC_STATIC_INLINE PetscErrorCode A B {PetscScalar _Val = _val; return A C;}

330: /*MC
331:       PetscPolymorphicFunction - allows defining a C++ polymorphic version of 
332:             a PETSc function that remove certain optional arguments for a simplier user interface
333:             and returns the computed value (istead of an error code)

335:      Not collective

337:    Synopsis:
338:    PetscPolymorphicFunction(Functionname,(arguments of C++ function),(arguments of C function),return type,return variable name)
339:  
340:    Level: developer

342:     Example:
343:       PetscPolymorphicFunction(VecNorm,(Vec x,NormType t),(x,t,&r),PetscReal,r) generates the new routine
344:          PetscReal VecNorm(Vec x,NormType t) = {PetscReal r; VecNorm(x,t,&r); return r;}

346: .seealso: PetscPolymorphicSubroutine()

348: M*/
349: #define PetscPolymorphicFunction(A,B,C,D,E) PETSC_STATIC_INLINE D A B {D E; A C;return E;}

351: #else
352: #define PetscPolymorphicSubroutine(A,B,C)
353: #define PetscPolymorphicScalar(A,B,C)
354: #define PetscPolymorphicFunction(A,B,C,D,E)
355: #endif

357: /*
358:     Extern indicates a PETSc function defined elsewhere
359: */
360: #if !defined(EXTERN)
361: #define EXTERN extern
362: #endif

364: /*
365:     Defines some elementary mathematics functions and constants.
366: */
367:  #include petscmath.h

369: /*
371: */


375: /*
376:        Basic PETSc constants
377: */

379: /*E
380:     PetscTruth - Logical variable. Actually an integer

382:    Level: beginner

384: E*/
385: typedef enum { PETSC_FALSE,PETSC_TRUE } PetscTruth;

388: /*MC
389:     PETSC_FALSE - False value of PetscTruth

391:     Level: beginner

393:     Note: Zero integer

395: .seealso: PetscTruth, PETSC_TRUE
396: M*/

398: /*MC
399:     PETSC_TRUE - True value of PetscTruth

401:     Level: beginner

403:     Note: Nonzero integer

405: .seealso: PetscTruth, PETSC_FALSE
406: M*/

408: /*MC
409:     PETSC_YES - Alias for PETSC_TRUE

411:     Level: beginner

413:     Note: Zero integer

415: .seealso: PetscTruth, PETSC_TRUE, PETSC_FALSE, PETSC_NO
416: M*/
417: #define PETSC_YES            PETSC_TRUE

419: /*MC
420:     PETSC_NO - Alias for PETSC_FALSE

422:     Level: beginner

424:     Note: Nonzero integer

426: .seealso: PetscTruth, PETSC_TRUE, PETSC_FALSE, PETSC_YES
427: M*/
428: #define PETSC_NO             PETSC_FALSE

430: /*MC
431:     PETSC_NULL - standard way of passing in a null or array or pointer

433:    Level: beginner

435:    Notes: accepted by many PETSc functions to not set a parameter and instead use
436:           some default

438:           This macro does not exist in Fortran; you must use PETSC_NULL_INTEGER, 
439:           PETSC_NULL_DOUBLE_PRECISION etc

441: .seealso: PETSC_DECIDE, PETSC_DEFAULT, PETSC_IGNORE, PETSC_DETERMINE

443: M*/
444: #define PETSC_NULL           0

446: /*MC
447:     PETSC_DECIDE - standard way of passing in integer or floating point parameter
448:        where you wish PETSc to use the default.

450:    Level: beginner

452: .seealso: PETSC_NULL, PETSC_DEFAULT, PETSC_IGNORE, PETSC_DETERMINE

454: M*/
455: #define PETSC_DECIDE         -1

457: /*MC
458:     PETSC_DEFAULT - standard way of passing in integer or floating point parameter
459:        where you wish PETSc to use the default.

461:    Level: beginner

463:    Fortran Notes: You need to use PETSC_DEFAULT_INTEGER or PETSC_DEFAULT_DOUBLE_PRECISION.

465: .seealso: PETSC_DECIDE, PETSC_NULL, PETSC_IGNORE, PETSC_DETERMINE

467: M*/
468: #define PETSC_DEFAULT        -2


471: /*MC
472:     PETSC_IGNORE - same as PETSC_NULL, means PETSc will ignore this argument

474:    Level: beginner

476:    Note: accepted by many PETSc functions to not set a parameter and instead use
477:           some default

479:    Fortran Notes: This macro does not exist in Fortran; you must use PETSC_NULL_INTEGER, 
480:           PETSC_NULL_DOUBLE_PRECISION etc

482: .seealso: PETSC_DECIDE, PETSC_DEFAULT, PETSC_NULL, PETSC_DETERMINE

484: M*/
485: #define PETSC_IGNORE         PETSC_NULL

487: /*MC
488:     PETSC_DETERMINE - standard way of passing in integer or floating point parameter
489:        where you wish PETSc to compute the required value.

491:    Level: beginner

493: .seealso: PETSC_DECIDE, PETSC_DEFAULT, PETSC_IGNORE, PETSC_NULL, VecSetSizes()

495: M*/
496: #define PETSC_DETERMINE      PETSC_DECIDE

498: /*MC
499:     PETSC_COMM_WORLD - the equivalent of the MPI_COMM_WORLD communicator which represents
500:            all the processs that PETSc knows about. 

502:    Level: beginner

504:    Notes: By default PETSC_COMM_WORLD and MPI_COMM_WORLD are identical unless you wish to 
505:           run PETSc on ONLY a subset of MPI_COMM_WORLD. In that case create your new (smaller)
506:           communicator, call it, say comm, and set PETSC_COMM_WORLD = comm BEFORE calling
507:           PetscInitialize()

509: .seealso: PETSC_COMM_SELF

511: M*/

514: /*MC
515:     PETSC_COMM_SELF - a duplicate of the MPI_COMM_SELF communicator which represents
516:            the current process

518:    Level: beginner

520:    Notes: PETSC_COMM_SELF and MPI_COMM_SELF are equivalent.

522: .seealso: PETSC_COMM_WORLD

524: M*/
525: #define PETSC_COMM_SELF MPI_COMM_SELF


530: EXTERN PetscErrorCode  PetscSetHelpVersionFunctions(PetscErrorCode (*)(MPI_Comm),PetscErrorCode (*)(MPI_Comm));
531: EXTERN PetscErrorCode  PetscCommDuplicate(MPI_Comm,MPI_Comm*,int*);
532: EXTERN PetscErrorCode  PetscCommDestroy(MPI_Comm*);

534: /*MC
535:    PetscMalloc - Allocates memory

537:    Input Parameter:
538: .  m - number of bytes to allocate

540:    Output Parameter:
541: .  result - memory allocated

543:    Synopsis:
544:    PetscErrorCode PetscMalloc(size_t m,void **result)

546:    Level: beginner

548:    Notes: Memory is always allocated at least double aligned

550:           If you request memory of zero size it will allocate no space and assign the pointer to 0; PetscFree() will 
551:           properly handle not freeing the null pointer.

553: .seealso: PetscFree(), PetscNew()

555:   Concepts: memory allocation

557: M*/
558: #define PetscMalloc(a,b)  ((a != 0) ? (*PetscTrMalloc)((a),__LINE__,__FUNCT__,__FILE__,__SDIR__,(void**)(b)) : (*(b) = 0,0) )

560: /*MC
561:    PetscMalloc2 - Allocates 2 chunks of  memory

563:    Input Parameter:
564: +  m1 - number of elements to allocate in 1st chunk  (may be zero)
565: .  t1 - type of first memory elements 
566: .  m2 - number of elements to allocate in 2nd chunk  (may be zero)
567: -  t2 - type of second memory elements

569:    Output Parameter:
570: +  r1 - memory allocated in first chunk
571: -  r2 - memory allocated in second chunk

573:    Synopsis:
574:    PetscErrorCode PetscMalloc2(size_t m1,type, t1,void **r1,size_t m2,type t2,void **r2)

576:    Level: developer

578:    Notes: Memory of first chunk is always allocated at least double aligned

580: .seealso: PetscFree(), PetscNew(), PetscMalloc()

582:   Concepts: memory allocation

584: M*/
585: #if defined(PETSC_USE_DEBUG)
586: #define PetscMalloc2(m1,t1,r1,m2,t2,r2) (PetscMalloc((m1)*sizeof(t1),r1) || PetscMalloc((m2)*sizeof(t2),r2))
587: #else
588: #define PetscMalloc2(m1,t1,r1,m2,t2,r2) (PetscMalloc((m1)*sizeof(t1)+(m2)*sizeof(t2),r1) || (*(r2) = (t2*)(*(r1)+m1),0))
589: #endif

591: /*MC
592:    PetscMalloc3 - Allocates 3 chunks of  memory

594:    Input Parameter:
595: +  m1 - number of elements to allocate in 1st chunk  (may be zero)
596: .  t1 - type of first memory elements 
597: .  m2 - number of elements to allocate in 2nd chunk  (may be zero)
598: .  t2 - type of second memory elements
599: .  m3 - number of elements to allocate in 3rd chunk  (may be zero)
600: -  t3 - type of third memory elements

602:    Output Parameter:
603: +  r1 - memory allocated in first chunk
604: .  r2 - memory allocated in second chunk
605: -  r3 - memory allocated in third chunk

607:    Synopsis:
608:    PetscErrorCode PetscMalloc3(size_t m1,type, t1,void **r1,size_t m2,type t2,void **r2,size_t m3,type t3,void **r3)

610:    Level: developer

612:    Notes: Memory of first chunk is always allocated at least double aligned

614: .seealso: PetscFree(), PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree3()

616:   Concepts: memory allocation

618: M*/
619: #if defined(PETSC_USE_DEBUG)
620: #define PetscMalloc3(m1,t1,r1,m2,t2,r2,m3,t3,r3) (PetscMalloc((m1)*sizeof(t1),r1) || PetscMalloc((m2)*sizeof(t2),r2) || PetscMalloc((m3)*sizeof(t3),r3))
621: #else
622: #define PetscMalloc3(m1,t1,r1,m2,t2,r2,m3,t3,r3) (PetscMalloc((m1)*sizeof(t1)+(m2)*sizeof(t2)+(m3)*sizeof(t3),r1) || (*(r2) = (t2*)(*(r1)+m1),*(r3) = (t3*)(*(r2)+m2),0))
623: #endif

625: /*MC
626:    PetscMalloc4 - Allocates 4 chunks of  memory

628:    Input Parameter:
629: +  m1 - number of elements to allocate in 1st chunk  (may be zero)
630: .  t1 - type of first memory elements 
631: .  m2 - number of elements to allocate in 2nd chunk  (may be zero)
632: .  t2 - type of second memory elements
633: .  m3 - number of elements to allocate in 3rd chunk  (may be zero)
634: .  t3 - type of third memory elements
635: .  m4 - number of elements to allocate in 4th chunk  (may be zero)
636: -  t4 - type of fourth memory elements

638:    Output Parameter:
639: +  r1 - memory allocated in first chunk
640: .  r2 - memory allocated in second chunk
641: .  r3 - memory allocated in third chunk
642: -  r4 - memory allocated in fourth chunk

644:    Synopsis:
645:    PetscErrorCode PetscMalloc4(size_t m1,type, t1,void **r1,size_t m2,type t2,void **r2,size_t m3,type t3,void **r3,size_t m4,type t4,void **r4)

647:    Level: developer

649:    Notes: Memory of first chunk is always allocated at least double aligned

651: .seealso: PetscFree(), PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree3(), PetscFree4()

653:   Concepts: memory allocation

655: M*/
656: #if defined(PETSC_USE_DEBUG)
657: #define PetscMalloc4(m1,t1,r1,m2,t2,r2,m3,t3,r3,m4,t4,r4) (PetscMalloc((m1)*sizeof(t1),r1) || PetscMalloc((m2)*sizeof(t2),r2) || PetscMalloc((m3)*sizeof(t3),r3) || PetscMalloc((m4)*sizeof(t4),r4))
658: #else
659: #define PetscMalloc4(m1,t1,r1,m2,t2,r2,m3,t3,r3,m4,t4,r4) (PetscMalloc((m1)*sizeof(t1)+(m2)*sizeof(t2)+(m3)*sizeof(t3)+(m4)*sizeof(t4),r1) || (*(r2) = (t2*)(*(r1)+m1),*(r3) = (t3*)(*(r2)+m2),*(r4) = (t4*)(*(r3)+m3),0))
660: #endif

662: /*MC
663:    PetscMalloc5 - Allocates 5 chunks of  memory

665:    Input Parameter:
666: +  m1 - number of elements to allocate in 1st chunk  (may be zero)
667: .  t1 - type of first memory elements 
668: .  m2 - number of elements to allocate in 2nd chunk  (may be zero)
669: .  t2 - type of second memory elements
670: .  m3 - number of elements to allocate in 3rd chunk  (may be zero)
671: .  t3 - type of third memory elements
672: .  m4 - number of elements to allocate in 4th chunk  (may be zero)
673: .  t4 - type of fourth memory elements
674: .  m5 - number of elements to allocate in 5th chunk  (may be zero)
675: -  t5 - type of fifth memory elements

677:    Output Parameter:
678: +  r1 - memory allocated in first chunk
679: .  r2 - memory allocated in second chunk
680: .  r3 - memory allocated in third chunk
681: .  r4 - memory allocated in fourth chunk
682: -  r5 - memory allocated in fifth chunk

684:    Synopsis:
685:    PetscErrorCode PetscMalloc5(size_t m1,type, t1,void **r1,size_t m2,type t2,void **r2,size_t m3,type t3,void **r3,size_t m4,type t4,void **r4,size_t m5,type t5,void **r5)

687:    Level: developer

689:    Notes: Memory of first chunk is always allocated at least double aligned

691: .seealso: PetscFree(), PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree3(), PetscFree4(), PetscFree5()

693:   Concepts: memory allocation

695: M*/
696: #if defined(PETSC_USE_DEBUG)
697: #define PetscMalloc5(m1,t1,r1,m2,t2,r2,m3,t3,r3,m4,t4,r4,m5,t5,r5) (PetscMalloc((m1)*sizeof(t1),r1) || PetscMalloc((m2)*sizeof(t2),r2) || PetscMalloc((m3)*sizeof(t3),r3) || PetscMalloc((m4)*sizeof(t4),r4) || PetscMalloc((m5)*sizeof(t5),r5))
698: #else
699: #define PetscMalloc5(m1,t1,r1,m2,t2,r2,m3,t3,r3,m4,t4,r4,m5,t5,r5) (PetscMalloc((m1)*sizeof(t1)+(m2)*sizeof(t2)+(m3)*sizeof(t3)+(m4)*sizeof(t4)+(m5)*sizeof(t5),r1) || (*(r2) = (t2*)(*(r1)+m1),*(r3) = (t3*)(*(r2)+m2),*(r4) = (t4*)(*(r3)+m3),*(r5) = (t5*)(*(r4)+m4),0))
700: #endif


703: /*MC
704:    PetscMalloc6 - Allocates 6 chunks of  memory

706:    Input Parameter:
707: +  m1 - number of elements to allocate in 1st chunk  (may be zero)
708: .  t1 - type of first memory elements 
709: .  m2 - number of elements to allocate in 2nd chunk  (may be zero)
710: .  t2 - type of second memory elements
711: .  m3 - number of elements to allocate in 3rd chunk  (may be zero)
712: .  t3 - type of third memory elements
713: .  m4 - number of elements to allocate in 4th chunk  (may be zero)
714: .  t4 - type of fourth memory elements
715: .  m5 - number of elements to allocate in 5th chunk  (may be zero)
716: .  t5 - type of fifth memory elements
717: .  m6 - number of elements to allocate in 6th chunk  (may be zero)
718: -  t6 - type of sixth memory elements

720:    Output Parameter:
721: +  r1 - memory allocated in first chunk
722: .  r2 - memory allocated in second chunk
723: .  r3 - memory allocated in third chunk
724: .  r4 - memory allocated in fourth chunk
725: .  r5 - memory allocated in fifth chunk
726: -  r6 - memory allocated in sixth chunk

728:    Synopsis:
729:    PetscErrorCode PetscMalloc6(size_t m1,type, t1,void **r1,size_t m2,type t2,void **r2,size_t m3,type t3,void **r3,size_t m4,type t4,void **r4,size_t m5,type t5,void **r5,size_t m6,type t6,void **r6)

731:    Level: developer

733:    Notes: Memory of first chunk is always allocated at least double aligned

735: .seealso: PetscFree(), PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree3(), PetscFree4(), PetscFree5(), PetscFree6()

737:   Concepts: memory allocation

739: M*/
740: #if defined(PETSC_USE_DEBUG)
741: #define PetscMalloc6(m1,t1,r1,m2,t2,r2,m3,t3,r3,m4,t4,r4,m5,t5,r5,m6,t6,r6) (PetscMalloc((m1)*sizeof(t1),r1) || PetscMalloc((m2)*sizeof(t2),r2) || PetscMalloc((m3)*sizeof(t3),r3) || PetscMalloc((m4)*sizeof(t4),r4) || PetscMalloc((m5)*sizeof(t5),r5) || PetscMalloc((m6)*sizeof(t6),r6))
742: #else
743: #define PetscMalloc6(m1,t1,r1,m2,t2,r2,m3,t3,r3,m4,t4,r4,m5,t5,r5,m6,t6,r6) (PetscMalloc((m1)*sizeof(t1)+(m2)*sizeof(t2)+(m3)*sizeof(t3)+(m4)*sizeof(t4)+(m5)*sizeof(t5)+(m6)*sizeof(t6),r1) || (*(r2) = (t2*)(*(r1)+m1),*(r3) = (t3*)(*(r2)+m2),*(r4) = (t4*)(*(r3)+m3),*(r5) = (t5*)(*(r4)+m4),*(r6) = (t6*)(*(r5)+m5),0))
744: #endif

746: /*MC
747:    PetscMalloc7 - Allocates 7 chunks of  memory

749:    Input Parameter:
750: +  m1 - number of elements to allocate in 1st chunk  (may be zero)
751: .  t1 - type of first memory elements 
752: .  m2 - number of elements to allocate in 2nd chunk  (may be zero)
753: .  t2 - type of second memory elements
754: .  m3 - number of elements to allocate in 3rd chunk  (may be zero)
755: .  t3 - type of third memory elements
756: .  m4 - number of elements to allocate in 4th chunk  (may be zero)
757: .  t4 - type of fourth memory elements
758: .  m5 - number of elements to allocate in 5th chunk  (may be zero)
759: .  t5 - type of fifth memory elements
760: .  m6 - number of elements to allocate in 6th chunk  (may be zero)
761: .  t6 - type of sixth memory elements
762: .  m7 - number of elements to allocate in 7th chunk  (may be zero)
763: -  t7 - type of sixth memory elements

765:    Output Parameter:
766: +  r1 - memory allocated in first chunk
767: .  r2 - memory allocated in second chunk
768: .  r3 - memory allocated in third chunk
769: .  r4 - memory allocated in fourth chunk
770: .  r5 - memory allocated in fifth chunk
771: .  r6 - memory allocated in sixth chunk
772: -  r7 - memory allocated in sixth chunk

774:    Synopsis:
775:    PetscErrorCode PetscMalloc7(size_t m1,type, t1,void **r1,size_t m2,type t2,void **r2,size_t m3,type t3,void **r3,size_t m4,type t4,void **r4,size_t m5,type t5,void **r5,size_t m6,type t6,void **r6,size_t m7,type t7,void **r7)

777:    Level: developer

779:    Notes: Memory of first chunk is always allocated at least double aligned

781: .seealso: PetscFree(), PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree3(), PetscFree4(), PetscFree5(), PetscFree6(), PetscFree7()

783:   Concepts: memory allocation

785: M*/
786: #if defined(PETSC_USE_DEBUG)
787: #define PetscMalloc7(m1,t1,r1,m2,t2,r2,m3,t3,r3,m4,t4,r4,m5,t5,r5,m6,t6,r6,m7,t7,r7) (PetscMalloc((m1)*sizeof(t1),r1) || PetscMalloc((m2)*sizeof(t2),r2) || PetscMalloc((m3)*sizeof(t3),r3) || PetscMalloc((m4)*sizeof(t4),r4) || PetscMalloc((m5)*sizeof(t5),r5) || PetscMalloc((m6)*sizeof(t6),r6) || PetscMalloc((m7)*sizeof(t7),r7))
788: #else
789: #define PetscMalloc7(m1,t1,r1,m2,t2,r2,m3,t3,r3,m4,t4,r4,m5,t5,r5,m6,t6,r6,m7,t7,r7) (PetscMalloc((m1)*sizeof(t1)+(m2)*sizeof(t2)+(m3)*sizeof(t3)+(m4)*sizeof(t4)+(m5)*sizeof(t5)+(m6)*sizeof(t6)+(m7)*sizeof(t7),r1) || (*(r2) = (t2*)(*(r1)+m1),*(r3) = (t3*)(*(r2)+m2),*(r4) = (t4*)(*(r3)+m3),*(r5) = (t5*)(*(r4)+m4),*(r6) = (t6*)(*(r5)+m5),*(r7) = (t7*)(*(r6)+m6),0))
790: #endif

792: /*MC
793:    PetscNew - Allocates memory of a particular type, zeros the memory!

795:    Input Parameter:
796: .  type - structure name of space to be allocated. Memory of size sizeof(type) is allocated

798:    Output Parameter:
799: .  result - memory allocated

801:    Synopsis:
802:    PetscErrorCode PetscNew(struct type,((type *))result)

804:    Level: beginner

806: .seealso: PetscFree(), PetscMalloc()

808:   Concepts: memory allocation

810: M*/
811: #define PetscNew(A,b)      (PetscMalloc(sizeof(A),(b)) || PetscMemzero(*(b),sizeof(A)))
812: #define PetscNewLog(o,A,b) (PetscNew(A,b) || ((o) ? PetscLogObjectMemory(o,sizeof(A)) : 0))

814: /*MC
815:    PetscFree - Frees memory

817:    Input Parameter:
818: .   memory - memory to free (the pointer is ALWAYS set to 0 upon sucess)

820:    Synopsis:
821:    PetscErrorCode PetscFree(void *memory)

823:    Level: beginner

825:    Notes: Memory must have been obtained with PetscNew() or PetscMalloc()

827: .seealso: PetscNew(), PetscMalloc(), PetscFreeVoid()

829:   Concepts: memory allocation

831: M*/
832: #define PetscFree(a)   ((a) ? ((*PetscTrFree)((void*)(a),__LINE__,__FUNCT__,__FILE__,__SDIR__) || ((a = 0),0)) : 0)

834: /*MC
835:    PetscFreeVoid - Frees memory

837:    Input Parameter:
838: .   memory - memory to free

840:    Synopsis:
841:    void PetscFreeVoid(void *memory)

843:    Level: beginner

845:    Notes: This is different from PetscFree() in that no error code is returned

847: .seealso: PetscFree(), PetscNew(), PetscMalloc()

849:   Concepts: memory allocation

851: M*/
852: #define PetscFreeVoid(a) ((*PetscTrFree)((a),__LINE__,__FUNCT__,__FILE__,__SDIR__),(a) = 0)


855: /*MC
856:    PetscFree2 - Frees 2 chunks of memory obtained with PetscMalloc2()

858:    Input Parameter:
859: +   memory1 - memory to free
860: -   memory2 - 2nd memory to free


863:    Synopsis:
864:    PetscErrorCode PetscFree2(void *memory1,void *memory2)

866:    Level: developer

868:    Notes: Memory must have been obtained with PetscMalloc2()

870: .seealso: PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree()

872:   Concepts: memory allocation

874: M*/
875: #if defined(PETSC_USE_DEBUG)
876: #define PetscFree2(m1,m2)   (PetscFree(m2) || PetscFree(m1))
877: #else
878: #define PetscFree2(m1,m2)   (PetscFree(m1))
879: #endif

881: /*MC
882:    PetscFree3 - Frees 3 chunks of memory obtained with PetscMalloc3()

884:    Input Parameter:
885: +   memory1 - memory to free
886: .   memory2 - 2nd memory to free
887: -   memory3 - 3rd memory to free


890:    Synopsis:
891:    PetscErrorCode PetscFree3(void *memory1,void *memory2,void *memory3)

893:    Level: developer

895:    Notes: Memory must have been obtained with PetscMalloc3()

897: .seealso: PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree(), PetscMalloc3()

899:   Concepts: memory allocation

901: M*/
902: #if defined(PETSC_USE_DEBUG)
903: #define PetscFree3(m1,m2,m3)   (PetscFree(m3) || PetscFree(m2) || PetscFree(m1))
904: #else
905: #define PetscFree3(m1,m2,m3)   (PetscFree(m1))
906: #endif

908: /*MC
909:    PetscFree4 - Frees 4 chunks of memory obtained with PetscMalloc4()

911:    Input Parameter:
912: +   m1 - memory to free
913: .   m2 - 2nd memory to free
914: .   m3 - 3rd memory to free
915: -   m4 - 4th memory to free


918:    Synopsis:
919:    PetscErrorCode PetscFree4(void *m1,void *m2,void *m3,void *m4)

921:    Level: developer

923:    Notes: Memory must have been obtained with PetscMalloc4()

925: .seealso: PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree(), PetscMalloc3(), PetscMalloc4()

927:   Concepts: memory allocation

929: M*/
930: #if defined(PETSC_USE_DEBUG)
931: #define PetscFree4(m1,m2,m3,m4)   (PetscFree(m4) || PetscFree(m3) || PetscFree(m2) || PetscFree(m1))
932: #else
933: #define PetscFree4(m1,m2,m3,m4)   (PetscFree(m1))
934: #endif

936: /*MC
937:    PetscFree5 - Frees 5 chunks of memory obtained with PetscMalloc5()

939:    Input Parameter:
940: +   m1 - memory to free
941: .   m2 - 2nd memory to free
942: .   m3 - 3rd memory to free
943: .   m4 - 4th memory to free
944: -   m5 - 5th memory to free


947:    Synopsis:
948:    PetscErrorCode PetscFree5(void *m1,void *m2,void *m3,void *m4,void *m5)

950:    Level: developer

952:    Notes: Memory must have been obtained with PetscMalloc5()

954: .seealso: PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree(), PetscMalloc3(), PetscMalloc4(), PetscMalloc5()

956:   Concepts: memory allocation

958: M*/
959: #if defined(PETSC_USE_DEBUG)
960: #define PetscFree5(m1,m2,m3,m4,m5)   (PetscFree(m5) || PetscFree(m4) || PetscFree(m3) || PetscFree(m2) || PetscFree(m1))
961: #else
962: #define PetscFree5(m1,m2,m3,m4,m5)   (PetscFree(m1))
963: #endif


966: /*MC
967:    PetscFree6 - Frees 6 chunks of memory obtained with PetscMalloc6()

969:    Input Parameter:
970: +   m1 - memory to free
971: .   m2 - 2nd memory to free
972: .   m3 - 3rd memory to free
973: .   m4 - 4th memory to free
974: .   m5 - 5th memory to free
975: -   m6 - 6th memory to free


978:    Synopsis:
979:    PetscErrorCode PetscFree6(void *m1,void *m2,void *m3,void *m4,void *m5,void *m6)

981:    Level: developer

983:    Notes: Memory must have been obtained with PetscMalloc6()

985: .seealso: PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree(), PetscMalloc3(), PetscMalloc4(), PetscMalloc5(), PetscMalloc6()

987:   Concepts: memory allocation

989: M*/
990: #if defined(PETSC_USE_DEBUG)
991: #define PetscFree6(m1,m2,m3,m4,m5,m6)   (PetscFree(m6) || PetscFree(m5) || PetscFree(m4) || PetscFree(m3) || PetscFree(m2) || PetscFree(m1))
992: #else
993: #define PetscFree6(m1,m2,m3,m4,m5,m6)   (PetscFree(m1))
994: #endif

996: /*MC
997:    PetscFree7 - Frees 7 chunks of memory obtained with PetscMalloc7()

999:    Input Parameter:
1000: +   m1 - memory to free
1001: .   m2 - 2nd memory to free
1002: .   m3 - 3rd memory to free
1003: .   m4 - 4th memory to free
1004: .   m5 - 5th memory to free
1005: .   m6 - 6th memory to free
1006: -   m7 - 7th memory to free


1009:    Synopsis:
1010:    PetscErrorCode PetscFree7(void *m1,void *m2,void *m3,void *m4,void *m5,void *m6,void *m7)

1012:    Level: developer

1014:    Notes: Memory must have been obtained with PetscMalloc7()

1016: .seealso: PetscNew(), PetscMalloc(), PetscMalloc2(), PetscFree(), PetscMalloc3(), PetscMalloc4(), PetscMalloc5(), PetscMalloc6(),
1017:           PetscMalloc7()

1019:   Concepts: memory allocation

1021: M*/
1022: #if defined(PETSC_USE_DEBUG)
1023: #define PetscFree7(m1,m2,m3,m4,m5,m6,m7)   (PetscFree(m7) || PetscFree(m6) || PetscFree(m5) || PetscFree(m4) || PetscFree(m3) || PetscFree(m2) || PetscFree(m1))
1024: #else
1025: #define PetscFree7(m1,m2,m3,m4,m5,m6,m7)   (PetscFree(m1))
1026: #endif

1028: EXTERN  PetscErrorCode (*PetscTrMalloc)(size_t,int,const char[],const char[],const char[],void**);
1029: EXTERN  PetscErrorCode (*PetscTrFree)(void*,int,const char[],const char[],const char[]);
1030: EXTERN PetscErrorCode   PetscMallocSet(PetscErrorCode (*)(size_t,int,const char[],const char[],const char[],void**),PetscErrorCode (*)(void*,int,const char[],const char[],const char[]));
1031: EXTERN PetscErrorCode   PetscMallocClear(void);

1033: /*
1034:    Routines for tracing memory corruption/bleeding with default PETSc 
1035:    memory allocation
1036: */
1037: EXTERN PetscErrorCode    PetscMallocDump(FILE *);
1038: EXTERN PetscErrorCode    PetscMallocDumpLog(FILE *);
1039: EXTERN PetscErrorCode    PetscMallocGetCurrentUsage(PetscLogDouble *);
1040: EXTERN PetscErrorCode    PetscMallocGetMaximumUsage(PetscLogDouble *);
1041: EXTERN PetscErrorCode    PetscMallocDebug(PetscTruth);
1042: EXTERN PetscErrorCode    PetscMallocValidate(int,const char[],const char[],const char[]);
1043: EXTERN PetscErrorCode    PetscMallocSetDumpLog(void);


1046: /*
1047:     Variable type where we stash PETSc object pointers in Fortran.
1048:     On most machines size(pointer) == sizeof(long) - except windows
1049:     where its sizeof(long long)
1050: */

1052: #if (PETSC_SIZEOF_VOID_P) == (PETSC_SIZEOF_LONG)
1053: #define PetscFortranAddr   long
1054: #elif  (PETSC_SIZEOF_VOID_P) == (PETSC_SIZEOF_LONG_LONG)
1055: #define PetscFortranAddr   long long
1056: #else
1057: #error "Unknown size for PetscFortranAddr! Send us a bugreport at petsc-maint@mcs.anl.gov"
1058: #endif

1060: /*E
1061:     PetscDataType - Used for handling different basic data types.

1063:    Level: beginner

1065:    Developer comment: It would be nice if we could always just use MPI Datatypes, why can we not?

1067: .seealso: PetscBinaryRead(), PetscBinaryWrite(), PetscDataTypeToMPIDataType(),
1068:           PetscDataTypeGetSize()

1070: E*/
1071: typedef enum {PETSC_INT = 0,PETSC_DOUBLE = 1,PETSC_COMPLEX = 2, PETSC_LONG = 3 ,PETSC_SHORT = 4,PETSC_FLOAT = 5,
1072:               PETSC_CHAR = 6,PETSC_LOGICAL = 7,PETSC_ENUM = 8,PETSC_TRUTH=9, PETSC_LONG_DOUBLE = 10} PetscDataType;

1075: #if defined(PETSC_USE_COMPLEX)
1076: #define PETSC_SCALAR PETSC_COMPLEX
1077: #else
1078: #if defined(PETSC_USE_SINGLE)
1079: #define PETSC_SCALAR PETSC_FLOAT
1080: #elif defined(PETSC_USE_LONG_DOUBLE)
1081: #define PETSC_SCALAR PETSC_LONG_DOUBLE
1082: #elif defined(PETSC_USE_INT)
1083: #define PETSC_SCALAR PETSC_INT
1084: #else
1085: #define PETSC_SCALAR PETSC_DOUBLE
1086: #endif
1087: #endif
1088: #if defined(PETSC_USE_SINGLE)
1089: #define PETSC_REAL PETSC_FLOAT
1090: #elif defined(PETSC_USE_LONG_DOUBLE)
1091: #define PETSC_REAL PETSC_LONG_DOUBLE
1092: #elif defined(PETSC_USE_INT)
1093: #define PETSC_REAL PETSC_INT
1094: #else
1095: #define PETSC_REAL PETSC_DOUBLE
1096: #endif
1097: #define PETSC_FORTRANADDR PETSC_LONG

1099: EXTERN PetscErrorCode  PetscDataTypeToMPIDataType(PetscDataType,MPI_Datatype*);
1100: EXTERN PetscErrorCode  PetscMPIDataTypeToPetscDataType(MPI_Datatype,PetscDataType*);
1101: EXTERN PetscErrorCode  PetscDataTypeGetSize(PetscDataType,size_t*);

1103: /*
1104:     Basic memory and string operations. These are usually simple wrappers
1105:    around the basic Unix system calls, but a few of them have additional
1106:    functionality and/or error checking.
1107: */
1108: EXTERN PetscErrorCode    PetscMemcpy(void*,const void *,size_t);
1109: EXTERN PetscErrorCode    PetscBitMemcpy(void*,PetscInt,const void*,PetscInt,PetscInt,PetscDataType);
1110: EXTERN PetscErrorCode    PetscMemmove(void*,void *,size_t);
1111: EXTERN PetscErrorCode    PetscMemzero(void*,size_t);
1112: EXTERN PetscErrorCode    PetscMemcmp(const void*,const void*,size_t,PetscTruth *);
1113: EXTERN PetscErrorCode    PetscStrlen(const char[],size_t*);
1114: EXTERN PetscErrorCode    PetscStrcmp(const char[],const char[],PetscTruth *);
1115: EXTERN PetscErrorCode    PetscStrgrt(const char[],const char[],PetscTruth *);
1116: EXTERN PetscErrorCode    PetscStrcasecmp(const char[],const char[],PetscTruth*);
1117: EXTERN PetscErrorCode    PetscStrncmp(const char[],const char[],size_t,PetscTruth*);
1118: EXTERN PetscErrorCode    PetscStrcpy(char[],const char[]);
1119: EXTERN PetscErrorCode    PetscStrcat(char[],const char[]);
1120: EXTERN PetscErrorCode    PetscStrncat(char[],const char[],size_t);
1121: EXTERN PetscErrorCode    PetscStrncpy(char[],const char[],size_t);
1122: EXTERN PetscErrorCode    PetscStrchr(const char[],char,char *[]);
1123: EXTERN PetscErrorCode    PetscStrtolower(char[]);
1124: EXTERN PetscErrorCode    PetscStrrchr(const char[],char,char *[]);
1125: EXTERN PetscErrorCode    PetscStrstr(const char[],const char[],char *[]);
1126: EXTERN PetscErrorCode    PetscStrrstr(const char[],const char[],char *[]);
1127: EXTERN PetscErrorCode    PetscStrallocpy(const char[],char *[]);
1128: EXTERN PetscErrorCode    PetscStrreplace(MPI_Comm,const char[],char[],size_t);
1129: #define      PetscStrfree(a) ((a) ? PetscFree(a) : 0) 

1131: /*S
1132:     PetscToken - 'Token' used for managing tokenizing strings

1134:   Level: intermediate

1136: .seealso: PetscTokenCreate(), PetscTokenFind(), PetscTokenDestroy()
1137: S*/
1138: typedef struct _p_PetscToken* PetscToken;

1140: EXTERN PetscErrorCode    PetscTokenCreate(const char[],const char,PetscToken*);
1141: EXTERN PetscErrorCode    PetscTokenFind(PetscToken,char *[]);
1142: EXTERN PetscErrorCode    PetscTokenDestroy(PetscToken);

1144: /*
1145:    These are  MPI operations for MPI_Allreduce() etc
1146: */
1147: EXTERN  MPI_Op PetscMaxSum_Op;
1148: #if defined(PETSC_USE_COMPLEX)
1149: EXTERN  MPI_Op PetscSum_Op;
1150: #else
1151: #define PetscSum_Op MPI_SUM
1152: #endif
1153: EXTERN PetscErrorCode  PetscMaxSum(MPI_Comm,const PetscInt[],PetscInt*,PetscInt*);

1155: /*S
1156:      PetscObject - any PETSc object, PetscViewer, Mat, Vec, KSP etc

1158:    Level: beginner

1160:    Note: This is the base class from which all objects appear.

1162: .seealso:  PetscObjectDestroy(), PetscObjectView(), PetscObjectGetName(), PetscObjectSetName()
1163: S*/
1164: typedef struct _p_PetscObject* PetscObject;

1166: /*S
1167:      PetscFList - Linked list of functions, possibly stored in dynamic libraries, accessed
1168:       by string name

1170:    Level: advanced

1172: .seealso:  PetscFListAdd(), PetscFListDestroy()
1173: S*/
1174: typedef struct _n_PetscFList *PetscFList;

1176: /*E
1177:   PetscFileMode - Access mode for a file.

1179:   Level: beginner

1181:   FILE_MODE_READ - open a file at its beginning for reading

1183:   FILE_MODE_WRITE - open a file at its beginning for writing (will create if the file does not exist)

1185:   FILE_MODE_APPEND - open a file at end for writing

1187:   FILE_MODE_UPDATE - open a file for updating, meaning for reading and writing

1189:   FILE_MODE_APPEND_UPDATE - open a file for updating, meaning for reading and writing, at the end

1191: .seealso: PetscViewerFileSetMode()
1192: E*/
1193: typedef enum {FILE_MODE_READ, FILE_MODE_WRITE, FILE_MODE_APPEND, FILE_MODE_UPDATE, FILE_MODE_APPEND_UPDATE} PetscFileMode;

1195:  #include petscviewer.h
1196:  #include petscoptions.h

1198: #define PETSC_SMALLEST_COOKIE 1211211
1201: EXTERN PetscErrorCode  PetscCookieRegister(const char[],PetscCookie *);

1203: /*
1204:    Routines that get memory usage information from the OS
1205: */
1206: EXTERN PetscErrorCode  PetscMemoryGetCurrentUsage(PetscLogDouble *);
1207: EXTERN PetscErrorCode  PetscMemoryGetMaximumUsage(PetscLogDouble *);
1208: EXTERN PetscErrorCode  PetscMemorySetGetMaximumUsage(void);
1209: EXTERN PetscErrorCode  PetscMemoryShowUsage(PetscViewer,const char[]);

1211: EXTERN PetscErrorCode  PetscInfoAllow(PetscTruth,const char []);
1212: EXTERN PetscErrorCode  PetscGetTime(PetscLogDouble*);
1213: EXTERN PetscErrorCode  PetscGetCPUTime(PetscLogDouble*);
1214: EXTERN PetscErrorCode  PetscSleep(int);

1216: /*
1217:     Initialization of PETSc
1218: */
1219: EXTERN PetscErrorCode  PetscInitialize(int*,char***,const char[],const char[]);
1220: PetscPolymorphicSubroutine(PetscInitialize,(int *argc,char ***args),(argc,args,PETSC_NULL,PETSC_NULL))
1221: EXTERN PetscErrorCode  PetscInitializeNoArguments(void);
1222: EXTERN PetscErrorCode  PetscInitialized(PetscTruth *);
1223: EXTERN PetscErrorCode  PetscFinalized(PetscTruth *);
1224: EXTERN PetscErrorCode  PetscFinalize(void);
1225: EXTERN PetscErrorCode PetscInitializeFortran(void);
1226: EXTERN PetscErrorCode  PetscGetArgs(int*,char ***);
1227: EXTERN PetscErrorCode  PetscGetArguments(char ***);
1228: EXTERN PetscErrorCode  PetscFreeArguments(char **);

1230: EXTERN PetscErrorCode  PetscEnd(void);
1231: EXTERN PetscErrorCode  PetscInitializePackage(const char[]);

1233: EXTERN PetscErrorCode  PetscOpenMPMerge(PetscMPIInt);
1234: EXTERN PetscErrorCode  PetscOpenMPSpawn(PetscMPIInt);
1235: EXTERN PetscErrorCode  PetscOpenMPFinalize(void);
1236: EXTERN PetscErrorCode  PetscOpenMPRun(MPI_Comm,PetscErrorCode (*)(MPI_Comm,void *),void*);
1237: EXTERN PetscErrorCode  PetscOpenMPFree(MPI_Comm,void*);
1238: EXTERN PetscErrorCode  PetscOpenMPNew(MPI_Comm,PetscInt,void**);

1240: EXTERN PetscErrorCode  PetscPythonInitialize(const char[],const char[]);
1241: EXTERN PetscErrorCode  PetscPythonFinalize(void);

1243: /*
1245:    function pointers. Since the regular C++ code expects its function pointers to be 
1246:    C++.
1247: */
1248: typedef void (**PetscVoidStarFunction)(void);
1249: typedef void (*PetscVoidFunction)(void);
1250: typedef PetscErrorCode (*PetscErrorCodeFunction)(void);

1252: /*
1253:    PetscTryMethod - Queries an object for a method, if it exists then calls it.
1254:               These are intended to be used only inside PETSc functions.
1255: */
1256: #define  PetscTryMethod(obj,A,B,C) \
1257:   0;{ PetscErrorCode (*f)B, __ierr; \
1258:     __PetscObjectQueryFunction((PetscObject)obj,A,(PetscVoidStarFunction)&f);CHKERRQ(__ierr); \
1259:     if (f) {__(*f)C;CHKERRQ(__ierr);}\
1260:   }
1261: #define  PetscUseMethod(obj,A,B,C) \
1262:   0;{ PetscErrorCode (*f)B, __ierr; \
1263:     __PetscObjectQueryFunction((PetscObject)obj,A,(PetscVoidStarFunction)&f);CHKERRQ(__ierr); \
1264:     if (f) {__(*f)C;CHKERRQ(__ierr);}\
1265:     else {SETERRQ1(PETSC_ERR_SUP,"Cannot locate function %s in object",A);} \
1266:   }
1267: /*
1268:     Functions that can act on any PETSc object.
1269: */
1270: EXTERN PetscErrorCode  PetscObjectCreate(MPI_Comm,PetscObject*);
1271: EXTERN PetscErrorCode  PetscObjectCreateGeneric(MPI_Comm, PetscCookie, const char [], PetscObject *);
1272: EXTERN PetscErrorCode  PetscObjectDestroy(PetscObject);
1273: EXTERN PetscErrorCode  PetscObjectExists(PetscObject,PetscTruth*);
1274: EXTERN PetscErrorCode  PetscObjectGetComm(PetscObject,MPI_Comm *);
1275: EXTERN PetscErrorCode  PetscObjectGetCookie(PetscObject,PetscCookie *);
1276: EXTERN PetscErrorCode  PetscObjectSetType(PetscObject,const char []);
1277: EXTERN PetscErrorCode  PetscObjectGetType(PetscObject,const char *[]);
1278: EXTERN PetscErrorCode  PetscObjectSetName(PetscObject,const char[]);
1279: EXTERN PetscErrorCode  PetscObjectGetName(PetscObject,const char*[]);
1280: EXTERN PetscErrorCode  PetscObjectSetTabLevel(PetscObject,PetscInt);
1281: EXTERN PetscErrorCode  PetscObjectGetTabLevel(PetscObject,PetscInt*);
1282: EXTERN PetscErrorCode  PetscObjectIncrementTabLevel(PetscObject,PetscObject,PetscInt);
1283: EXTERN PetscErrorCode  PetscObjectReference(PetscObject);
1284: EXTERN PetscErrorCode  PetscObjectGetReference(PetscObject,PetscInt*);
1285: EXTERN PetscErrorCode  PetscObjectDereference(PetscObject);
1286: EXTERN PetscErrorCode  PetscObjectGetNewTag(PetscObject,PetscMPIInt *);
1287: EXTERN PetscErrorCode  PetscObjectView(PetscObject,PetscViewer);
1288: EXTERN PetscErrorCode  PetscObjectCompose(PetscObject,const char[],PetscObject);
1289: EXTERN PetscErrorCode  PetscObjectQuery(PetscObject,const char[],PetscObject *);
1290: EXTERN PetscErrorCode  PetscObjectComposeFunction(PetscObject,const char[],const char[],void (*)(void));
1291: EXTERN PetscErrorCode  PetscObjectSetFromOptions(PetscObject);
1292: EXTERN PetscErrorCode  PetscObjectSetUp(PetscObject);
1293: EXTERN PetscErrorCode  PetscCommGetNewTag(MPI_Comm,PetscMPIInt *);

1295: /*MC
1296:    PetscObjectComposeFunctionDynamic - Associates a function with a given PETSc object. 
1297:                        
1298:    Collective on PetscObject

1300:    Input Parameters:
1301: +  obj - the PETSc object; this must be cast with a (PetscObject), for example, 
1302:          PetscObjectCompose((PetscObject)mat,...);
1303: .  name - name associated with the child function
1304: .  fname - name of the function
1305: -  ptr - function pointer (or PETSC_NULL if using dynamic libraries)

1307:    Level: advanced

1309:     Synopsis:
1310:     PetscErrorCode PetscObjectComposeFunctionDynamic(PetscObject obj,const char name[],const char fname[],void *ptr)

1312:    Notes:
1313:    To remove a registered routine, pass in a PETSC_NULL rname and fnc().

1315:    PetscObjectComposeFunctionDynamic() can be used with any PETSc object (such as
1316:    Mat, Vec, KSP, SNES, etc.) or any user-provided object. 

1319:    work in C++/complex with dynamic link libraries (config/configure.py options --with-shared --with-dynamic)
1320:    enabled.

1322:    Concepts: objects^composing functions
1323:    Concepts: composing functions
1324:    Concepts: functions^querying
1325:    Concepts: objects^querying
1326:    Concepts: querying objects

1328: .seealso: PetscObjectQueryFunction()
1329: M*/
1330: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
1331: #define PetscObjectComposeFunctionDynamic(a,b,c,d) PetscObjectComposeFunction(a,b,c,0)
1332: #else
1333: #define PetscObjectComposeFunctionDynamic(a,b,c,d) PetscObjectComposeFunction(a,b,c,(PetscVoidFunction)(d))
1334: #endif

1336: EXTERN PetscErrorCode  PetscObjectQueryFunction(PetscObject,const char[],void (**)(void));
1337: EXTERN PetscErrorCode  PetscObjectSetOptionsPrefix(PetscObject,const char[]);
1338: EXTERN PetscErrorCode  PetscObjectAppendOptionsPrefix(PetscObject,const char[]);
1339: EXTERN PetscErrorCode  PetscObjectPrependOptionsPrefix(PetscObject,const char[]);
1340: EXTERN PetscErrorCode  PetscObjectGetOptionsPrefix(PetscObject,const char*[]);
1341: EXTERN PetscErrorCode  PetscObjectPublish(PetscObject);
1342: EXTERN PetscErrorCode  PetscObjectChangeTypeName(PetscObject,const char[]);
1343: EXTERN PetscErrorCode  PetscObjectRegisterDestroy(PetscObject);
1344: EXTERN PetscErrorCode  PetscObjectRegisterDestroyAll(void);
1345: EXTERN PetscErrorCode  PetscObjectName(PetscObject);
1346: EXTERN PetscErrorCode  PetscTypeCompare(PetscObject,const char[],PetscTruth*);
1347: EXTERN PetscErrorCode  PetscRegisterFinalize(PetscErrorCode (*)(void));
1348: EXTERN PetscErrorCode  PetscRegisterFinalizeAll(void);

1350: /*
1351:     Defines PETSc error handling.
1352: */
1353:  #include petscerror.h

1355: /*S
1356:      PetscOList - Linked list of PETSc objects, accessable by string name

1358:    Level: advanced

1360: .seealso:  PetscOListAdd(), PetscOListDestroy(), PetscOListFind()
1361: S*/
1362: typedef struct _n_PetscOList *PetscOList;

1364: EXTERN PetscErrorCode  PetscOListDestroy(PetscOList);
1365: EXTERN PetscErrorCode  PetscOListFind(PetscOList,const char[],PetscObject*);
1366: EXTERN PetscErrorCode  PetscOListReverseFind(PetscOList,PetscObject,char**);
1367: EXTERN PetscErrorCode  PetscOListAdd(PetscOList *,const char[],PetscObject);
1368: EXTERN PetscErrorCode  PetscOListDuplicate(PetscOList,PetscOList *);

1370: /*
1371:     Dynamic library lists. Lists of names of routines in dynamic 
1372:   link libraries that will be loaded as needed.
1373: */
1374: EXTERN PetscErrorCode  PetscFListAdd(PetscFList*,const char[],const char[],void (*)(void));
1375: EXTERN PetscErrorCode  PetscFListDestroy(PetscFList*);
1376: EXTERN PetscErrorCode  PetscFListFind(PetscFList,MPI_Comm,const char[],void (**)(void));
1377: EXTERN PetscErrorCode  PetscFListPrintTypes(PetscFList,MPI_Comm,FILE*,const char[],const char[],const char[],const char[]);
1378: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
1379: #define    PetscFListAddDynamic(a,b,p,c) PetscFListAdd(a,b,p,0)
1380: #else
1381: #define    PetscFListAddDynamic(a,b,p,c) PetscFListAdd(a,b,p,(void (*)(void))c)
1382: #endif
1383: EXTERN PetscErrorCode  PetscFListDuplicate(PetscFList,PetscFList *);
1384: EXTERN PetscErrorCode  PetscFListView(PetscFList,PetscViewer);
1385: EXTERN PetscErrorCode  PetscFListConcat(const char [],const char [],char []);
1386: EXTERN PetscErrorCode  PetscFListGet(PetscFList,char ***,int*);

1388: /*S
1389:      PetscDLLibrary - Linked list of dynamics libraries to search for functions

1391:    Level: advanced

1393:    --with-shared --with-dynamic must be used with config/configure.py to use dynamic libraries

1395: .seealso:  PetscDLLibraryOpen()
1396: S*/
1397: typedef struct _n_PetscDLLibrary *PetscDLLibrary;
1399: EXTERN PetscErrorCode  PetscDLLibraryAppend(MPI_Comm,PetscDLLibrary *,const char[]);
1400: EXTERN PetscErrorCode  PetscDLLibraryPrepend(MPI_Comm,PetscDLLibrary *,const char[]);
1401: EXTERN PetscErrorCode  PetscDLLibrarySym(MPI_Comm,PetscDLLibrary *,const char[],const char[],void **);
1402: EXTERN PetscErrorCode  PetscDLLibraryPrintPath(PetscDLLibrary);
1403: EXTERN PetscErrorCode  PetscDLLibraryRetrieve(MPI_Comm,const char[],char *,size_t,PetscTruth *);
1404: EXTERN PetscErrorCode  PetscDLLibraryOpen(MPI_Comm,const char[],PetscDLLibrary *);
1405: EXTERN PetscErrorCode  PetscDLLibraryClose(PetscDLLibrary);
1406: EXTERN PetscErrorCode  PetscDLLibraryCCAAppend(MPI_Comm,PetscDLLibrary *,const char[]);

1408: /*
1409:      Useful utility routines
1410: */
1411: EXTERN PetscErrorCode  PetscSplitOwnership(MPI_Comm,PetscInt*,PetscInt*);
1412: EXTERN PetscErrorCode  PetscSplitOwnershipBlock(MPI_Comm,PetscInt,PetscInt*,PetscInt*);
1413: EXTERN PetscErrorCode  PetscSequentialPhaseBegin(MPI_Comm,PetscMPIInt);
1414: PetscPolymorphicSubroutine(PetscSequentialPhaseBegin,(MPI_Comm comm),(comm,1))
1415: PetscPolymorphicSubroutine(PetscSequentialPhaseBegin,(void),(PETSC_COMM_WORLD,1))
1416: EXTERN PetscErrorCode  PetscSequentialPhaseEnd(MPI_Comm,PetscMPIInt);
1417: PetscPolymorphicSubroutine(PetscSequentialPhaseEnd,(MPI_Comm comm),(comm,1))
1418: PetscPolymorphicSubroutine(PetscSequentialPhaseEnd,(void),(PETSC_COMM_WORLD,1))
1419: EXTERN PetscErrorCode  PetscBarrier(PetscObject);
1420: EXTERN PetscErrorCode  PetscMPIDump(FILE*);

1422: #define PetscNot(a) ((a) ? PETSC_FALSE : PETSC_TRUE)
1423: /*
1424:     Defines basic graphics available from PETSc.
1425: */
1426:  #include petscdraw.h

1428: /*
1429:     Defines the base data structures for all PETSc objects
1430: */
1431:  #include private/petscimpl.h
1432: /*
1433:      Defines PETSc profiling.
1434: */
1435:  #include petsclog.h

1437: /*
1438:           For locking, unlocking and destroying AMS memories associated with 
1439:     PETSc objects. Not currently used.
1440: */
1441: #define PetscPublishAll(v)           0
1442: #define PetscObjectTakeAccess(obj)   0
1443: #define PetscObjectGrantAccess(obj)  0
1444: #define PetscObjectDepublish(obj)    0

1446: /*
1447:       Simple PETSc parallel IO for ASCII printing
1448: */
1449: EXTERN PetscErrorCode   PetscFixFilename(const char[],char[]);
1450: EXTERN PetscErrorCode   PetscFOpen(MPI_Comm,const char[],const char[],FILE**);
1451: EXTERN PetscErrorCode   PetscFClose(MPI_Comm,FILE*);
1452: EXTERN PetscErrorCode   PetscFPrintf(MPI_Comm,FILE*,const char[],...) PETSC_PRINTF_FORMAT_CHECK(3,4);
1453: EXTERN PetscErrorCode   PetscPrintf(MPI_Comm,const char[],...)  PETSC_PRINTF_FORMAT_CHECK(2,3);
1454: EXTERN PetscErrorCode   PetscSNPrintf(char*,size_t,const char [],...);

1456: /* These are used internally by PETSc ASCII IO routines*/
1457: #include <stdarg.h>
1458: EXTERN PetscErrorCode   PetscVSNPrintf(char*,size_t,const char[],int*,va_list);
1459: EXTERN PetscErrorCode   (*PetscVFPrintf)(FILE*,const char[],va_list);
1460: EXTERN PetscErrorCode   PetscVFPrintfDefault(FILE*,const char[],va_list);

1462: /*MC
1463:     PetscErrorPrintf - Prints error messages.

1465:     Not Collective

1467:    Synopsis:
1468:      PetscErrorCode (*PetscErrorPrintf)(const char format[],...);

1470:     Input Parameters:
1471: .   format - the usual printf() format string 

1473:    Options Database Keys:
1474: +    -error_output_stdout - cause error messages to be printed to stdout instead of the
1475:          (default) stderr
1476: -    -error_output_none to turn off all printing of error messages (does not change the way the 
1477:           error is handled.)

1479:    Notes: Use
1480: $     PetscErrorPrintf = PetscErrorPrintfNone; to turn off all printing of error messages (does not change the way the 
1481: $                        error is handled.) and
1482: $     PetscErrorPrintf = PetscErrorPrintfDefault; to turn it back on

1484:           Use
1485:      PETSC_STDERR = FILE* obtained from a file open etc. to have stderr printed to the file. 
1486:      PETSC_STDOUT = FILE* obtained from a file open etc. to have stdout printed to the file. 

1488:           Use
1489:       PetscPushErrorHandler() to provide your own error handler that determines what kind of messages to print

1491:    Level: developer

1493:     Fortran Note:
1494:     This routine is not supported in Fortran.

1496:     Concepts: error messages^printing
1497:     Concepts: printing^error messages

1499: .seealso: PetscFPrintf(), PetscSynchronizedPrintf(), PetscHelpPrintf(), PetscPrintf(), PetscErrorHandlerPush()
1500: M*/
1501: EXTERN  PetscErrorCode (*PetscErrorPrintf)(const char[],...);

1503: /*MC
1504:     PetscHelpPrintf - Prints help messages.

1506:     Not Collective

1508:    Synopsis:
1509:      PetscErrorCode (*PetscHelpPrintf)(const char format[],...);

1511:     Input Parameters:
1512: .   format - the usual printf() format string 

1514:    Level: developer

1516:     Fortran Note:
1517:     This routine is not supported in Fortran.

1519:     Concepts: help messages^printing
1520:     Concepts: printing^help messages

1522: .seealso: PetscFPrintf(), PetscSynchronizedPrintf(), PetscErrorPrintf()
1523: M*/
1524: EXTERN  PetscErrorCode  (*PetscHelpPrintf)(MPI_Comm,const char[],...);

1526: EXTERN PetscErrorCode  PetscErrorPrintfDefault(const char [],...);
1527: EXTERN PetscErrorCode  PetscErrorPrintfNone(const char [],...);
1528: EXTERN PetscErrorCode  PetscHelpPrintfDefault(MPI_Comm,const char [],...);

1530: #if defined(PETSC_HAVE_POPEN)
1531: EXTERN PetscErrorCode   PetscPOpen(MPI_Comm,const char[],const char[],const char[],FILE **);
1532: EXTERN PetscErrorCode   PetscPClose(MPI_Comm,FILE*);
1533: #endif

1535: EXTERN PetscErrorCode   PetscSynchronizedPrintf(MPI_Comm,const char[],...) PETSC_PRINTF_FORMAT_CHECK(2,3);
1536: EXTERN PetscErrorCode   PetscSynchronizedFPrintf(MPI_Comm,FILE*,const char[],...) PETSC_PRINTF_FORMAT_CHECK(3,4);
1537: EXTERN PetscErrorCode   PetscSynchronizedFlush(MPI_Comm);
1538: EXTERN PetscErrorCode   PetscSynchronizedFGets(MPI_Comm,FILE*,size_t,char[]);
1539: EXTERN PetscErrorCode   PetscStartMatlab(MPI_Comm,const char[],const char[],FILE**);
1540: EXTERN PetscErrorCode   PetscStartJava(MPI_Comm,const char[],const char[],FILE**);
1541: EXTERN PetscErrorCode   PetscGetPetscDir(const char*[]);

1543: EXTERN PetscErrorCode   PetscPopUpSelect(MPI_Comm,const char*,const char*,int,const char**,int*);

1545: /*S
1546:      PetscContainer - Simple PETSc object that contains a pointer to any required data

1548:    Level: advanced

1550: .seealso:  PetscObject, PetscContainerCreate()
1551: S*/
1553: typedef struct _p_PetscContainer*  PetscContainer;
1554: EXTERN PetscErrorCode  PetscContainerGetPointer(PetscContainer,void **);
1555: EXTERN PetscErrorCode  PetscContainerSetPointer(PetscContainer,void *);
1556: EXTERN PetscErrorCode  PetscContainerDestroy(PetscContainer);
1557: EXTERN PetscErrorCode  PetscContainerCreate(MPI_Comm,PetscContainer *);
1558: EXTERN PetscErrorCode  PetscContainerSetUserDestroy(PetscContainer, PetscErrorCode (*)(void*));

1560: /*
1561:    For use in debuggers 
1562: */

1566: EXTERN PetscErrorCode  PetscIntView(PetscInt,PetscInt[],PetscViewer);
1567: EXTERN PetscErrorCode  PetscRealView(PetscInt,PetscReal[],PetscViewer);
1568: EXTERN PetscErrorCode  PetscScalarView(PetscInt,PetscScalar[],PetscViewer);

1570: /*
1571:     Allows accessing Matlab Engine
1572: */
1573:  #include petscmatlab.h

1575: /*
1576:       Determine if some of the kernel computation routines use
1577:    Fortran (rather than C) for the numerical calculations. On some machines
1578:    and compilers (like complex numbers) the Fortran version of the routines
1579:    is faster than the C/C++ versions. The flag --with-fortran-kernels
1580:    should be used with config/configure.py to turn these on.
1581: */
1582: #if defined(PETSC_USE_FORTRAN_KERNELS)

1584: #if !defined(PETSC_USE_FORTRAN_KERNEL_MULTCRL)
1585: #define PETSC_USE_FORTRAN_KERNEL_MULTCRL
1586: #endif

1588: #if !defined(PETSC_USE_FORTRAN_KERNEL_MULTCSRPERM)
1589: #define PETSC_USE_FORTRAN_KERNEL_MULTCSRPERM
1590: #endif

1592: #if !defined(PETSC_USE_FORTRAN_KERNEL_MULTAIJ)
1593: #define PETSC_USE_FORTRAN_KERNEL_MULTAIJ
1594: #endif

1596: #if !defined(PETSC_USE_FORTRAN_KERNEL_MULTTRANSPOSEAIJ)
1597: #define PETSC_USE_FORTRAN_KERNEL_MULTTRANSPOSEAIJ
1598: #endif

1600: #if !defined(PETSC_USE_FORTRAN_KERNEL_NORM)
1601: #define PETSC_USE_FORTRAN_KERNEL_NORM
1602: #endif

1604: #if !defined(PETSC_USE_FORTRAN_KERNEL_MAXPY)
1605: #define PETSC_USE_FORTRAN_KERNEL_MAXPY
1606: #endif

1608: #if !defined(PETSC_USE_FORTRAN_KERNEL_SOLVEAIJ)
1609: #define PETSC_USE_FORTRAN_KERNEL_SOLVEAIJ
1610: #endif

1612: #if !defined(PETSC_USE_FORTRAN_KERNEL_RELAXAIJ)
1613: #define PETSC_USE_FORTRAN_KERNEL_RELAXAIJ
1614: #endif

1616: #if !defined(PETSC_USE_FORTRAN_KERNEL_SOLVEBAIJ)
1617: #define PETSC_USE_FORTRAN_KERNEL_SOLVEBAIJ
1618: #endif

1620: #if !defined(PETSC_USE_FORTRAN_KERNEL_MULTADDAIJ)
1621: #define PETSC_USE_FORTRAN_KERNEL_MULTADDAIJ
1622: #endif

1624: #if !defined(PETSC_USE_FORTRAN_KERNEL_MDOT)
1625: #define PETSC_USE_FORTRAN_KERNEL_MDOT
1626: #endif

1628: #if !defined(PETSC_USE_FORTRAN_KERNEL_XTIMESY)
1629: #define PETSC_USE_FORTRAN_KERNEL_XTIMESY
1630: #endif

1632: #if !defined(PETSC_USE_FORTRAN_KERNEL_AYPX)
1633: #define PETSC_USE_FORTRAN_KERNEL_AYPX
1634: #endif

1636: #if !defined(PETSC_USE_FORTRAN_KERNEL_WAXPY)
1637: #define PETSC_USE_FORTRAN_KERNEL_WAXPY
1638: #endif

1640: #endif

1642: /*
1643:     Macros for indicating code that should be compiled with a C interface,
1644:    rather than a C++ interface. Any routines that are dynamically loaded
1645:    (such as the PCCreate_XXX() routines) must be wrapped so that the name
1646:    mangler does not change the functions symbol name. This just hides the 
1648: */
1649: #if defined(__cplusplus)
1652: #else
1655: #endif

1657: /* --------------------------------------------------------------------*/

1659: /*MC
1660:     size - integer variable used to contain the number of processors in
1661:            the relevent MPI_Comm

1663:    Level: beginner

1665: .seealso: rank, comm
1666: M*/

1668: /*MC
1669:     rank - integer variable used to contain the number of this processor relative
1670:            to all in the relevent MPI_Comm

1672:    Level: beginner

1674: .seealso: size, comm
1675: M*/

1677: /*MC
1678:     comm - MPI_Comm used in the current routine or object

1680:    Level: beginner

1682: .seealso: size, rank
1683: M*/

1685: /*MC
1686:     MPI_Comm - the basic object used by MPI to determine which processes are involved in a 
1687:         communication

1689:    Level: beginner

1691:    Note: This manual page is a place-holder because MPICH does not have a manual page for MPI_Comm

1693: .seealso: size, rank, comm, PETSC_COMM_WORLD, PETSC_COMM_SELF
1694: M*/

1696: /*MC
1697:     PetscScalar - PETSc type that represents either a double precision real number, a double precision
1698:        complex number, a single precision real number, a long double or an int - if the code is configured 
1699:        with --with-scalar-type=real,complex --with-precision=single,double,longdouble,int,matsingle


1702:    Level: beginner

1704: .seealso: PetscReal, PassiveReal, PassiveScalar, MPIU_SCALAR, PetscInt
1705: M*/

1707: /*MC
1708:     PetscReal - PETSc type that represents a real number version of PetscScalar

1710:    Level: beginner

1712: .seealso: PetscScalar, PassiveReal, PassiveScalar
1713: M*/

1715: /*MC
1716:     PassiveScalar - PETSc type that represents a PetscScalar
1717:    Level: beginner

1719:     This is the same as a PetscScalar except in code that is automatically differentiated it is
1720:    treated as a constant (not an indendent or dependent variable)

1722: .seealso: PetscReal, PassiveReal, PetscScalar
1723: M*/

1725: /*MC
1726:     PassiveReal - PETSc type that represents a PetscReal

1728:    Level: beginner

1730:     This is the same as a PetscReal except in code that is automatically differentiated it is
1731:    treated as a constant (not an indendent or dependent variable)

1733: .seealso: PetscScalar, PetscReal, PassiveScalar
1734: M*/

1736: /*MC
1737:     MPIU_SCALAR - MPI datatype corresponding to PetscScalar

1739:    Level: beginner

1741:     Note: In MPI calls that require an MPI datatype that matches a PetscScalar or array of PetscScalars
1742:           pass this value

1744: .seealso: PetscReal, PassiveReal, PassiveScalar, PetscScalar, MPIU_INT
1745: M*/

1747: #if defined(PETSC_HAVE_MPIIO)
1748: #if !defined(PETSC_WORDS_BIGENDIAN)
1751: #else
1752: #define MPIU_File_write_all(a,b,c,d,e) MPI_File_write_all(a,b,c,d,e) 
1753: #define MPIU_File_read_all(a,b,c,d,e) MPI_File_read_all(a,b,c,d,e) 
1754: #endif
1755: #endif

1757: /* the following petsc_static_inline require petscerror.h */

1759: /* Limit MPI to 32-bits */
1760: #define PETSC_MPI_INT_MAX  2147483647
1761: #define PETSC_MPI_INT_MIN -2147483647
1762: /* Limit BLAS to 32-bits */
1763: #define PETSC_BLAS_INT_MAX  2147483647
1764: #define PETSC_BLAS_INT_MIN -2147483647

1766: #if defined(PETSC_USE_64BIT_INDICES)
1767: #define PetscMPIIntCheck(a)  if ((a) > PETSC_MPI_INT_MAX) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Message too long for MPI")
1768: #define PetscBLASIntCheck(a)  if ((a) > PETSC_BLAS_INT_MAX) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Array too long for BLAS/LAPACK")
1769: #define PetscMPIIntCast(a) (a);PetscMPIIntCheck(a)
1770: #define PetscBLASIntCast(a) (a);PetscBLASIntCheck(a)
1771: #else
1772: #define PetscMPIIntCheck(a) 
1773: #define PetscBLASIntCheck(a) 
1774: #define PetscMPIIntCast(a) a
1775: #define PetscBLASIntCast(a) a
1776: #endif  


1779: /*
1780:      The IBM include files define hz, here we hide it so that it may be used
1781:    as a regular user variable.
1782: */
1783: #if defined(hz)
1784: #undef hz
1785: #endif

1787: /*  For arrays that contain filenames or paths */


1790: #if defined(PETSC_HAVE_LIMITS_H)
1791: #include <limits.h>
1792: #endif
1793: #if defined(PETSC_HAVE_SYS_PARAM_H)
1794: #include <sys/param.h>
1795: #endif
1796: #if defined(PETSC_HAVE_SYS_TYPES_H)
1797: #include <sys/types.h>
1798: #endif
1799: #if defined(MAXPATHLEN)
1800: #  define PETSC_MAX_PATH_LEN     MAXPATHLEN
1801: #elif defined(MAX_PATH)
1802: #  define PETSC_MAX_PATH_LEN     MAX_PATH
1803: #elif defined(_MAX_PATH)
1804: #  define PETSC_MAX_PATH_LEN     _MAX_PATH
1805: #else
1806: #  define PETSC_MAX_PATH_LEN     4096
1807: #endif

1809: /* Special support for C++ */
1810:  #include petsc.hh

1812: /*MC

1814:     UsingFortran - Fortran can be used with PETSc in four distinct approaches

1816: $    1) classic Fortran 77 style
1817: $#include "petscXXX.h" to work with material from the XXX component of PETSc
1818: $       XXX variablename
1819: $      You cannot use this approach if you wish to use the Fortran 90 specific PETSc routines
1820: $      which end in F90; such as VecGetArrayF90()
1821: $
1822: $    2) classic Fortran 90 style
1823: $#include "petscXXX.h" 
1824: $#include "petscXXX.h90" to work with material from the XXX component of PETSc
1825: $       XXX variablename
1826: $
1827: $    3) Using Fortran modules
1828: $#include "petscXXXdef.h" 
1829: $         use petscXXXX
1830: $       XXX variablename
1831: $
1832: $    4) Use Fortran modules and Fortran data types for PETSc types
1833: $#include "petscXXXdef.h" 
1834: $         use petscXXXX
1835: $       type(XXX) variablename
1836: $      To use this approach you must config/configure.py PETSc with the additional
1837: $      option --with-fortran-datatypes You cannot use the type(XXX) declaration approach without using Fortran modules

1839:     Finally if you absolutely do not want to use any #include you can use either 

1841: $    3a) skip the #include BUT you cannot use any PETSc data type names like Vec, Mat, PetscInt, PetscErrorCode etc
1842: $        and you must declare the variables as integer, for example 
1843: $        integer variablename
1844: $
1845: $    4a) skip the #include, you use the object types like type(Vec) type(Mat) but cannot use the data type
1846: $        names like PetscErrorCode, PetscInt etc. again for those you must use integer

1848:    We recommend either 2 or 3. Approaches 2 and 3 provide type checking for most PETSc function calls; 4 has type checking 
1849: for only a few PETSc functions.

1851:    Fortran type checking with interfaces is strick, this means you cannot pass a scalar value when an array value
1852: is expected (even though it is legal Fortran). For example when setting a single value in a matrix with MatSetValues()
1853: you cannot have something like
1854: $      PetscInt row,col
1855: $      PetscScalar val
1856: $        ...
1857: $      call MatSetValues(mat,1,row,1,col,val,INSERT_VALUES,ierr)
1858: You must instead have 
1859: $      PetscInt row(1),col(1)
1860: $      PetscScalar val(1)
1861: $        ...
1862: $      call MatSetValues(mat,1,row,1,col,val,INSERT_VALUES,ierr)


1865:     See the example src/vec/vec/examples/tutorials/ex20f90.F90 for an example that can use all four approaches

1867:     Developer Notes: The finclude/petscXXXdef.h contain all the #defines (would be typedefs in C code) these
1868:      automatically include their predecessors; for example finclude/petscvecdef.h includes finclude/petscisdef.h

1870:      The finclude/petscXXXX.h contain all the parameter statements for that package. These automatically include
1871:      their finclude/petscXXXdef.h file but DO NOT automatically include their predecessors;  for example 
1872:      finclude/petscvec.h does NOT automatically include finclude/petscis.h

1874:      The finclude/ftn-custom/petscXXXdef.h90 are not intended to be used directly in code, they define the
1875:      Fortran data type type(XXX) (for example type(Vec)) when PETSc is config/configure.py with the --with-fortran-datatypes option.

1877:      The finclude/ftn-custom/petscXXX.h90 (not included directly by code) contain interface definitions for
1878:      the PETSc Fortran stubs that have different bindings then their C version (for example VecGetArrayF90).

1880:      The finclude/ftn-auto/petscXXX.h90 (not included directly by code) contain interface definitions generated
1881:      automatically by "make allfortranstubs".

1883:      The finclude/petscXXX.h90 includes the custom finclude/ftn-custom/petscXXX.h90 and if config/configure.py 
1884:      was run with --with-fortran-interfaces it also includes the finclude/ftn-auto/petscXXX.h90 These DO NOT automatically
1885:      include their predecessors

1887:     Level: beginner

1889: M*/
1891: #endif