int main(int argc, char **argv)
{
  PetscBool      transform = PETSC_FALSE;
  PetscErrorCode ierr;

  ierr = PetscInitialize(&argc, &argv, NULL, help);CHKERRQ(ierr);
  ierr = PetscOptionsGetBool(NULL, "-transform", &transform, NULL);CHKERRQ(ierr);
  if (transform) {ierr = PetscPrintf(PETSC_COMM_WORLD, "Using random transforms");CHKERRQ(ierr);}
  ierr = TestTriangle(PETSC_COMM_SELF, PETSC_FALSE, transform);CHKERRQ(ierr);
  ierr = TestTriangle(PETSC_COMM_SELF, PETSC_TRUE,  transform);CHKERRQ(ierr);
  ierr = TestQuadrilateral(PETSC_COMM_SELF, PETSC_FALSE, transform);CHKERRQ(ierr);
  ierr = TestQuadrilateral(PETSC_COMM_SELF, PETSC_TRUE,  transform);CHKERRQ(ierr);
  ierr = TestTetrahedron(PETSC_COMM_SELF, PETSC_FALSE, transform);CHKERRQ(ierr);
  ierr = TestTetrahedron(PETSC_COMM_SELF, PETSC_TRUE,  transform);CHKERRQ(ierr);
  ierr = TestHexahedron(PETSC_COMM_SELF, PETSC_FALSE, transform);CHKERRQ(ierr);
  ierr = TestHexahedron(PETSC_COMM_SELF, PETSC_TRUE, transform);CHKERRQ(ierr);
  ierr = PetscFinalize();
  return 0;
}

int main(int argc,char **args)
{
  PetscErrorCode ierr;
  PetscBool      usempiio = PETSC_FALSE;

  PetscInitialize(&argc,&args,(char*)0,help);

  ierr = PetscOptionsGetBool(NULL,"-usempiio",&usempiio,NULL);CHKERRQ(ierr);
  if (!usempiio) {
    ierr = TestDMDAVec(PETSC_FALSE);CHKERRQ(ierr);
  } else {
#if defined(PETSC_HAVE_MPIIO)
    ierr = TestDMDAVec(PETSC_TRUE);CHKERRQ(ierr);
#else
    ierr = PetscPrintf(PETSC_COMM_WORLD,"Warning: Executing TestDMDAVec(PETSC_TRUE) requires a working MPI-2 implementation\n");CHKERRQ(ierr);
#endif
  }
  ierr = PetscFinalize();
  return 0;
}

PetscErrorCode PetscDrawSetUpColormap_Shared(Display *display,int screen,Visual *visual,Colormap colormap)
{
  int            i;
  XColor         colordef,ecolordef;
  PetscBool      fast = PETSC_FALSE;
  PetscErrorCode ierr;

  PetscFunctionBegin;
  if (colormap) gColormap = colormap;
  else          gColormap = DefaultColormap(display,screen);

  /* set the basic colors into the color map */
  for (i=0; i<PETSC_DRAW_BASIC_COLORS; i++) {
    XAllocNamedColor(display,gColormap,colornames[i],&colordef,&ecolordef);
    gCmapping[i] = colordef.pixel;
  }

  /* set the contour colors into the colormap */
  ierr = PetscOptionsGetBool(NULL,NULL,"-draw_fast",&fast,NULL);CHKERRQ(ierr);
  if (!fast) {
    int ncolors = 256-PETSC_DRAW_BASIC_COLORS;
    unsigned char *red,*green,*blue;
    ierr = PetscMalloc3(ncolors,&red,ncolors,&green,ncolors,&blue);CHKERRQ(ierr);
    ierr = PetscDrawUtilitySetCmap(NULL,ncolors,red,green,blue);CHKERRQ(ierr);
    for (i=0; i<ncolors; i++) {
      colordef.red   = (unsigned short)(red[i]   << 8);
      colordef.green = (unsigned short)(green[i] << 8);
      colordef.blue  = (unsigned short)(blue[i]  << 8);
      colordef.flags = DoRed|DoGreen|DoBlue;
      XAllocColor(display,gColormap,&colordef);
      gCmapping[PETSC_DRAW_BASIC_COLORS + i] = colordef.pixel;
    }
    ierr = PetscFree3(red,green,blue);CHKERRQ(ierr);
  }

  ierr = PetscInfo(0,"Successfully allocated colors\n");CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

/*
   PetscSplitReductionCreate - Creates a data structure to contain the queued information.
*/
static PetscErrorCode  PetscSplitReductionCreate(MPI_Comm comm,PetscSplitReduction **sr)
{
  PetscErrorCode ierr;

  PetscFunctionBegin;
  ierr               = PetscNew(PetscSplitReduction,sr);CHKERRQ(ierr);
  (*sr)->numopsbegin = 0;
  (*sr)->numopsend   = 0;
  (*sr)->state       = STATE_BEGIN;
  (*sr)->maxops      = 32;
  ierr               = PetscMalloc(2*32*sizeof(PetscScalar),&(*sr)->lvalues);CHKERRQ(ierr);
  ierr               = PetscMalloc(2*32*sizeof(PetscScalar),&(*sr)->gvalues);CHKERRQ(ierr);
  ierr               = PetscMalloc(32*sizeof(void*),&(*sr)->invecs);CHKERRQ(ierr);
  (*sr)->comm        = comm;
  (*sr)->request     = MPI_REQUEST_NULL;
  ierr               = PetscMalloc(32*sizeof(PetscInt),&(*sr)->reducetype);CHKERRQ(ierr);
  (*sr)->async = PETSC_FALSE;
#if defined(PETSC_HAVE_MPIX_IALLREDUCE) || defined(PETSC_HAVE_PAMI) || defined(PETSC_HAVE_DCMF)
  (*sr)->async = PETSC_TRUE;    /* Enable by default */
  ierr = PetscOptionsGetBool(PETSC_NULL,"-splitreduction_async",&(*sr)->async,PETSC_NULL);CHKERRQ(ierr);
#endif
  PetscFunctionReturn(0);
}

PETSC_EXTERN PetscErrorCode PetscSFCreate_Window(PetscSF sf)
{
  PetscSF_Window *w = (PetscSF_Window*)sf->data;
  PetscErrorCode ierr;

  PetscFunctionBegin;
  sf->ops->SetUp           = PetscSFSetUp_Window;
  sf->ops->SetFromOptions  = PetscSFSetFromOptions_Window;
  sf->ops->Reset           = PetscSFReset_Window;
  sf->ops->Destroy         = PetscSFDestroy_Window;
  sf->ops->View            = PetscSFView_Window;
  sf->ops->Duplicate       = PetscSFDuplicate_Window;
  sf->ops->BcastBegin      = PetscSFBcastBegin_Window;
  sf->ops->BcastEnd        = PetscSFBcastEnd_Window;
  sf->ops->ReduceBegin     = PetscSFReduceBegin_Window;
  sf->ops->ReduceEnd       = PetscSFReduceEnd_Window;
  sf->ops->FetchAndOpBegin = PetscSFFetchAndOpBegin_Window;
  sf->ops->FetchAndOpEnd   = PetscSFFetchAndOpEnd_Window;

  ierr     = PetscNewLog(sf,&w);CHKERRQ(ierr);
  sf->data = (void*)w;
  w->sync  = PETSCSF_WINDOW_SYNC_FENCE;

  ierr = PetscObjectComposeFunction((PetscObject)sf,"PetscSFWindowSetSyncType_C",PetscSFWindowSetSyncType_Window);CHKERRQ(ierr);
  ierr = PetscObjectComposeFunction((PetscObject)sf,"PetscSFWindowGetSyncType_C",PetscSFWindowGetSyncType_Window);CHKERRQ(ierr);

#if defined(OMPI_MAJOR_VERSION) && (OMPI_MAJOR_VERSION < 1 || (OMPI_MAJOR_VERSION == 1 && OMPI_MINOR_VERSION <= 6))
  {
    PetscBool ackbug = PETSC_FALSE;
    ierr = PetscOptionsGetBool(NULL,NULL,"-acknowledge_ompi_onesided_bug",&ackbug,NULL);CHKERRQ(ierr);
    if (ackbug) {
      ierr = PetscInfo(sf,"Acknowledged Open MPI bug, proceeding anyway. Expect memory corruption.\n");CHKERRQ(ierr);
    } else SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_LIB,"Open MPI is known to be buggy (https://svn.open-mpi.org/trac/ompi/ticket/1905 and 2656), use -acknowledge_ompi_onesided_bug to proceed");
  }
#endif
  PetscFunctionReturn(0);
}

/*@C
   PFDestroy - Destroys PF context that was created with PFCreate().

   Collective on PF

   Input Parameter:
.  pf - the function context

   Level: beginner

.keywords: PF, destroy

.seealso: PFCreate(), PFSet(), PFSetType()
@*/
PetscErrorCode  PFDestroy(PF *pf)
{
  PetscErrorCode ierr;
  PetscBool      flg = PETSC_FALSE;

  PetscFunctionBegin;
  if (!*pf) PetscFunctionReturn(0);
  PetscValidHeaderSpecific((*pf),PF_CLASSID,1);
  if (--((PetscObject)(*pf))->refct > 0) PetscFunctionReturn(0);

  ierr = PetscOptionsGetBool(((PetscObject)(*pf))->prefix,"-pf_view",&flg,PETSC_NULL);CHKERRQ(ierr);
  if (flg) {
    PetscViewer viewer;
    ierr = PetscViewerASCIIGetStdout(((PetscObject)(*pf))->comm,&viewer);CHKERRQ(ierr);
    ierr = PFView((*pf),viewer);CHKERRQ(ierr);
  }

  /* if memory was published with AMS then destroy it */
  ierr = PetscObjectDepublish((*pf));CHKERRQ(ierr);

  if ((*pf)->ops->destroy) {ierr =  (*(*pf)->ops->destroy)((*pf)->data);CHKERRQ(ierr);}
  ierr = PetscHeaderDestroy(pf);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

int main(int argc,char **argv)
{
  PetscInt       M = -10,N = -8;
  PetscErrorCode ierr;
  PetscBool      flg = PETSC_FALSE;
  DM             da;
  Vec            global1,global2,global3;
  DMDABoundaryType bx = DMDA_BOUNDARY_NONE,by = DMDA_BOUNDARY_NONE;
  DMDAStencilType  stype = DMDA_STENCIL_BOX;

  ierr = PetscInitialize(&argc,&argv,(char*)0,help);CHKERRQ(ierr); 
  ierr = PetscOptionsGetBool(PETSC_NULL,"-star_stencil",&flg,PETSC_NULL);CHKERRQ(ierr);
  if (flg) stype = DMDA_STENCIL_STAR;
      
  /* Create distributed array and get vectors */
  ierr = DMDACreate2d(PETSC_COMM_WORLD,bx,by,stype,M,N,PETSC_DECIDE,PETSC_DECIDE,1,1,PETSC_NULL,PETSC_NULL,&da);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da,&global1);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da,&global2);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da,&global1);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da,&global2);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da,&global1);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da,&global3);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da,&global2);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da,&global1);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da,&global3);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da,&global2);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da,&global1);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da,&global3);CHKERRQ(ierr);
  ierr = DMGetGlobalVector(da,&global2);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da,&global1);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da,&global3);CHKERRQ(ierr);
  ierr = DMRestoreGlobalVector(da,&global2);CHKERRQ(ierr);
  ierr = DMDestroy(&da);CHKERRQ(ierr);
  ierr = PetscFinalize();
  return 0;
}

/*@
   MatPartitioningApply - Gets a partitioning for a matrix.

   Collective on Mat

   Input Parameters:
.  matp - the matrix partitioning object

   Output Parameters:
.   partitioning - the partitioning. For each local node this tells the processor
                   number that that node is assigned to.

   Options Database Keys:
   To specify the partitioning through the options database, use one of
   the following
$    -mat_partitioning_type parmetis, -mat_partitioning current
   To see the partitioning result
$    -mat_partitioning_view

   Level: beginner

   The user can define additional partitionings; see MatPartitioningRegister().

.keywords: matrix, get, partitioning

.seealso:  MatPartitioningRegister(), MatPartitioningCreate(),
           MatPartitioningDestroy(), MatPartitioningSetAdjacency(), ISPartitioningToNumbering(),
           ISPartitioningCount()
@*/
PetscErrorCode  MatPartitioningApply(MatPartitioning matp,IS *partitioning)
{
  PetscErrorCode ierr;
  PetscBool      flag = PETSC_FALSE;

  PetscFunctionBegin;
  PetscValidHeaderSpecific(matp,MAT_PARTITIONING_CLASSID,1);
  PetscValidPointer(partitioning,2);
  if (!matp->adj->assembled) SETERRQ(PetscObjectComm((PetscObject)matp),PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled matrix");
  if (matp->adj->factortype) SETERRQ(PetscObjectComm((PetscObject)matp),PETSC_ERR_ARG_WRONGSTATE,"Not for factored matrix");
  if (!matp->ops->apply) SETERRQ(PetscObjectComm((PetscObject)matp),PETSC_ERR_ARG_WRONGSTATE,"Must set type with MatPartitioningSetFromOptions() or MatPartitioningSetType()");
  ierr = PetscLogEventBegin(MAT_Partitioning,matp,0,0,0);CHKERRQ(ierr);
  ierr = (*matp->ops->apply)(matp,partitioning);CHKERRQ(ierr);
  ierr = PetscLogEventEnd(MAT_Partitioning,matp,0,0,0);CHKERRQ(ierr);

  ierr = PetscOptionsGetBool(NULL,"-mat_partitioning_view",&flag,NULL);CHKERRQ(ierr);
  if (flag) {
    PetscViewer viewer;
    ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)matp),&viewer);CHKERRQ(ierr);
    ierr = MatPartitioningView(matp,viewer);CHKERRQ(ierr);
    ierr = ISView(*partitioning,viewer);CHKERRQ(ierr);
  }
  PetscFunctionReturn(0);
}

PetscErrorCode  PetscSSEIsEnabled(MPI_Comm comm,PetscBool  *lflag,PetscBool  *gflag)
{
  PetscErrorCode ierr;
  PetscBool      disabled_option;

  PetscFunctionBegin;
  if (petsc_sse_local_is_untested && petsc_sse_global_is_untested) {
    disabled_option = PETSC_FALSE;

    ierr = PetscOptionsGetBool(NULL,"-disable_sse",&disabled_option,NULL);CHKERRQ(ierr);
    if (disabled_option) {
      petsc_sse_local_is_untested  = PETSC_FALSE;
      petsc_sse_enabled_local      = PETSC_FALSE;
      petsc_sse_global_is_untested = PETSC_FALSE;
      petsc_sse_enabled_global     = PETSC_FALSE;
    }

    if (petsc_sse_local_is_untested) {
      ierr = PetscSSEHardwareTest(&petsc_sse_enabled_local);CHKERRQ(ierr);
      if (petsc_sse_enabled_local) {
        ierr = PetscSSEOSEnabledTest(&petsc_sse_enabled_local);CHKERRQ(ierr);
      }
      petsc_sse_local_is_untested = PETSC_FALSE;
    }

    if (gflag && petsc_sse_global_is_untested) {
      ierr = MPI_Allreduce(&petsc_sse_enabled_local,&petsc_sse_enabled_global,1,MPIU_BOOL,MPI_LAND,comm);CHKERRQ(ierr);

      petsc_sse_global_is_untested = PETSC_FALSE;
    }
  }

  if (lflag) *lflag = petsc_sse_enabled_local;
  if (gflag) *gflag = petsc_sse_enabled_global;
  PetscFunctionReturn(0);
}

/*@C
   MatGetColoring - Gets a coloring for a matrix, from its sparsity structure,
      to reduce the number of function evaluations needed to compute a sparse Jacobian via differencing.

   Collective on Mat

   Input Parameters:
.  mat - the matrix
.  type - type of coloring, one of the following:
$      MATCOLORINGNATURAL - natural (one color for each column, very slow)
$      MATCOLORINGSL - smallest-last
$      MATCOLORINGLF - largest-first
$      MATCOLORINGID - incidence-degree

   Output Parameters:
.   iscoloring - the coloring

   Options Database Keys:
   To specify the coloring through the options database, use one of
   the following
$    -mat_coloring_type natural, -mat_coloring_type sl, -mat_coloring_type lf,
$    -mat_coloring_type id
   To see the coloring use
$    -mat_coloring_view

   Level: intermediate

   Notes:
$    A graph coloring C(A)  is a division of vertices  so that two vertices of the same color do not share any common edges.
$    A suitable coloring for a  smoother  is simply C(A).
$    A suitable coloring for efficient Jacobian computation is a division of the columns so that two columns of the same color do not share any common rows.
$         This corresponds to C(A^{T} A).  This is what MatGetColoring() computes.

   The user can define additional colorings; see MatColoringRegisterDynamic().

   For parallel matrices currently converts to sequential matrix and uses the sequential coloring
   on that.

   The colorings SL, LF, and ID are obtained via the Minpack software that was
   converted to C using f2c.

   For BAIJ matrices this colors the blocks. The true number of colors would be block size times the number of colors
   returned here.

   References:
$     Thomas F. Coleman and Jorge J. More, Estimation of Sparse {J}acobian Matrices and Graph Coloring Problems,
$         SIAM Journal on Numerical Analysis, 1983, pages 187-209, volume 20
$     Jorge J. Mor\'{e} and Danny C. Sorenson and  Burton S. Garbow and Kenneth E. Hillstrom, The {MINPACK} Project,
$         Sources and Development of Mathematical Software, Wayne R. Cowell editor, 1984, pages 88-111

.keywords: matrix, get, coloring

.seealso:  MatGetColoringTypeFromOptions(), MatColoringRegisterDynamic(), MatFDColoringCreate(),
           SNESDefaultComputeJacobianColor()
@*/
PetscErrorCode  MatGetColoring(Mat mat,MatColoringType type,ISColoring *iscoloring)
{
  PetscBool      flag;
  PetscErrorCode ierr,(*r)(Mat,MatColoringType,ISColoring *);
  char           tname[PETSC_MAX_PATH_LEN];
  MPI_Comm       comm;

  PetscFunctionBegin;
  PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
  PetscValidPointer(iscoloring,3);
  if (!mat->assembled) SETERRQ(((PetscObject)mat)->comm,PETSC_ERR_ARG_WRONGSTATE,"Not for unassembled matrix");
  if (mat->factortype) SETERRQ(((PetscObject)mat)->comm,PETSC_ERR_ARG_WRONGSTATE,"Not for factored matrix");

  /* look for type on command line */
  if (!MatColoringRegisterAllCalled) {ierr = MatColoringRegisterAll(PETSC_NULL);CHKERRQ(ierr);}
  ierr = PetscOptionsGetString(((PetscObject)mat)->prefix,"-mat_coloring_type",tname,256,&flag);CHKERRQ(ierr);
  if (flag) { type = tname; }

  ierr = PetscObjectGetComm((PetscObject)mat,&comm);CHKERRQ(ierr);
  ierr = PetscFunctionListFind(comm,MatColoringList,type,PETSC_TRUE,(void (**)(void)) &r);CHKERRQ(ierr);
  if (!r) SETERRQ1(((PetscObject)mat)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Unknown or unregistered type: %s",type);

  ierr = PetscLogEventBegin(MAT_GetColoring,mat,0,0,0);CHKERRQ(ierr);
  ierr = (*r)(mat,type,iscoloring);CHKERRQ(ierr);
  ierr = PetscLogEventEnd(MAT_GetColoring,mat,0,0,0);CHKERRQ(ierr);

  ierr = PetscInfo1(mat,"Number of colors %d\n",(*iscoloring)->n);CHKERRQ(ierr);
  flag = PETSC_FALSE;
  ierr = PetscOptionsGetBool(PETSC_NULL,"-mat_coloring_view",&flag,PETSC_NULL);CHKERRQ(ierr);
  if (flag) {
    PetscViewer viewer;
    ierr = PetscViewerASCIIGetStdout((*iscoloring)->comm,&viewer);CHKERRQ(ierr);
    ierr = ISColoringView(*iscoloring,viewer);CHKERRQ(ierr);
  }
  PetscFunctionReturn(0);
}

/*
   PetscSplitReductionCreate - Creates a data structure to contain the queued information.
*/
static PetscErrorCode  PetscSplitReductionCreate(MPI_Comm comm,PetscSplitReduction **sr)
{
  PetscErrorCode ierr;

  PetscFunctionBegin;
  ierr               = PetscNew(sr);CHKERRQ(ierr);
  (*sr)->numopsbegin = 0;
  (*sr)->numopsend   = 0;
  (*sr)->state       = STATE_BEGIN;
  (*sr)->maxops      = 32;
  ierr               = PetscMalloc1(2*32,&(*sr)->lvalues);CHKERRQ(ierr);
  ierr               = PetscMalloc1(2*32,&(*sr)->gvalues);CHKERRQ(ierr);
  ierr               = PetscMalloc1(32,&(*sr)->invecs);CHKERRQ(ierr);
  (*sr)->comm        = comm;
  (*sr)->request     = MPI_REQUEST_NULL;
  ierr               = PetscMalloc1(32,&(*sr)->reducetype);CHKERRQ(ierr);
  (*sr)->async       = PETSC_FALSE;
#if defined(PETSC_HAVE_MPI_IALLREDUCE) || defined(PETSC_HAVE_MPIX_IALLREDUCE)
  (*sr)->async = PETSC_TRUE;    /* Enable by default */
#endif
  /* always check for option; so that tests that run on systems without support don't warn about unhandled options */
  ierr = PetscOptionsGetBool(NULL,NULL,"-splitreduction_async",&(*sr)->async,NULL);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

PetscErrorCode GetParams(AppCtx *user)
{
  PetscErrorCode ierr;
  PetscBool      flg;

  PetscFunctionBeginUser;
  /* Set default parameters */
  user->tsmonitor = PETSC_FALSE;
  user->xmin      = 0.0;    user->xmax    = 1.0;
  user->ymin      = 0.0;    user->ymax    = 1.0;
  user->T         = 0.0002; user->dt      = 0.0001;
  user->gamma     = 3.2E-4; user->theta_c = 0;

  ierr = PetscOptionsGetBool(NULL,"-ts_monitor",&user->tsmonitor,NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetReal(NULL,"-xmin",&user->xmin,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetReal(NULL,"-xmax",&user->xmax,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetReal(NULL,"-ymin",&user->ymin,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetReal(NULL,"-ymax",&user->ymax,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetReal(NULL,"-T",&user->T,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetReal(NULL,"-dt",&user->dt,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetScalar(NULL,"-gamma",&user->gamma,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetScalar(NULL,"-theta_c",&user->theta_c,&flg);CHKERRQ(ierr);
  PetscFunctionReturn(0);
}

/*@C
   AOCreate - Creates an application ordering.

   Collective on MPI_Comm

   Input Parameters:
.  comm - MPI communicator that is to share AO

   Output Parameter:
.  ao - the new application ordering

   Options Database Key:
+   -ao_type <aotype> - create ao with particular format
-   -ao_view - call AOView() at the conclusion of AOCreate()

   Level: beginner

.keywords: AO, create

.seealso: AOSetIS(), AODestroy(), AOPetscToApplication(), AOApplicationToPetsc()
@*/
PetscErrorCode  AOCreate(MPI_Comm comm,AO *ao)
{
  PetscErrorCode ierr;
  AO             aonew;
  PetscBool      opt;

  PetscFunctionBegin;
  PetscValidPointer(ao,2);
  *ao = PETSC_NULL;
#ifndef PETSC_USE_DYNAMIC_LIBRARIES
  ierr = AOInitializePackage(PETSC_NULL);CHKERRQ(ierr);
#endif

  ierr = PetscHeaderCreate(aonew,_p_AO,struct _AOOps,AO_CLASSID,-1,"AO","Application Ordering","AO",comm,AODestroy,AOView);CHKERRQ(ierr);
  ierr = PetscMemzero(aonew->ops, sizeof(struct _AOOps));CHKERRQ(ierr);
  *ao = aonew;

  opt = PETSC_FALSE;
  ierr = PetscOptionsGetBool(PETSC_NULL, "-ao_view", &opt,PETSC_NULL);CHKERRQ(ierr);
  if (opt) {
    ierr = AOView(aonew, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  }
  PetscFunctionReturn(0);
}

int main(int argc, char* argv[])
{ 
  int status = 0;
#ifdef CH_USE_PETSC
  PetscInt nlev;
  PetscErrorCode ierr;
  char string[64]; 
  PetscBool set;
  PetscInitialize(&argc,&argv,(char*)0,help);

  PetscOptionsGetInt(PETSC_NULL,"-n",&s_nCells0,PETSC_NULL);
  PetscOptionsGetBool(PETSC_NULL,"-bc_debug",&s_debug,PETSC_NULL);
  PetscOptionsGetBool(PETSC_NULL,"-plot_fas",&s_amrfas,PETSC_NULL);
  PetscOptionsGetBool(PETSC_NULL,"-plot_mg",&s_amrmg,PETSC_NULL);
  PetscOptionsGetBool(PETSC_NULL,"-corner_stencil",&s_corner_stencil,PETSC_NULL);
  PetscOptionsGetBool(PETSC_NULL,"-matlab",&s_matlab,PETSC_NULL);
  PetscOptionsGetBool(PETSC_NULL,"-plot",&s_plot,PETSC_NULL);
  PetscOptionsGetReal(PETSC_NULL,"-error_thresh",&s_error_thresh,PETSC_NULL);
  PetscOptionsGetInt(PETSC_NULL,"-blocking_factor",&s_blockingfactor,PETSC_NULL);
  PetscOptionsGetInt(PETSC_NULL,"-order",&s_order,PETSC_NULL);
  PetscOptionsGetInt(PETSC_NULL,"-nesting",&s_nesting,PETSC_NULL);
  PetscOptionsGetInt(PETSC_NULL,"-max_box_size",&s_maxboxsz,PETSC_NULL);
  PetscOptionsGetInt(PETSC_NULL,"-refinement_ratio",&s_refRatio,PETSC_NULL);
  PetscOptionsGetString(PETSC_NULL,"-amr_type",string,64,&set);
  if (set && strcmp(string,"error")==0) s_amr_type_iserror = true;
  else s_amr_type_iserror = false;

  nlev = 4;
  PetscOptionsGetInt(PETSC_NULL,"-nlevels",&nlev,PETSC_NULL);
  ierr = go(nlev,status); CHKERRQ(ierr);

  CH_TIMER_REPORT();

  ierr = PetscFinalize(); CHKERRQ(ierr);
#endif
  return status;
}

int main(int argc,char **argv)
{
  PetscErrorCode     ierr;                  /* used to check for functions returning nonzeros */
  PetscReal          zero=0.0;
  Vec                x;                     /* solution vector */
  Mat                H;
  Tao                tao;                   /* Tao solver context */
  PetscBool          flg, test_lmvm = PETSC_FALSE;
  PetscMPIInt        size,rank;                  /* number of processes running */
  AppCtx             user;                  /* user-defined application context */
  TaoConvergedReason reason;
  PetscInt           its, recycled_its=0, oneshot_its=0;

  /* Initialize TAO and PETSc */
  ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
  if (size >1) SETERRQ(PETSC_COMM_SELF,1,"Incorrect number of processors");

  /* Initialize problem parameters */
  user.n = 2; user.alpha = 99.0; user.chained = PETSC_FALSE;
  /* Check for command line arguments to override defaults */
  ierr = PetscOptionsGetInt(NULL,NULL,"-n",&user.n,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetReal(NULL,NULL,"-alpha",&user.alpha,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetBool(NULL,NULL,"-chained",&user.chained,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetBool(NULL,NULL,"-test_lmvm",&test_lmvm,&flg);CHKERRQ(ierr);

  /* Allocate vectors for the solution and gradient */
  ierr = VecCreateSeq(PETSC_COMM_SELF,user.n,&x);CHKERRQ(ierr);
  ierr = MatCreateSeqBAIJ(PETSC_COMM_SELF,2,user.n,user.n,1,NULL,&H);CHKERRQ(ierr);

  /* The TAO code begins here */

  /* Create TAO solver with desired solution method */
  ierr = TaoCreate(PETSC_COMM_SELF,&tao);CHKERRQ(ierr);
  ierr = TaoSetType(tao,TAOLMVM);CHKERRQ(ierr);

  /* Set solution vec and an initial guess */
  ierr = VecSet(x, zero);CHKERRQ(ierr);
  ierr = TaoSetInitialVector(tao,x);CHKERRQ(ierr);

  /* Set routines for function, gradient, hessian evaluation */
  ierr = TaoSetObjectiveAndGradientRoutine(tao,FormFunctionGradient,&user);CHKERRQ(ierr);
  ierr = TaoSetHessianRoutine(tao,H,H,FormHessian,&user);CHKERRQ(ierr);

  /* Check for TAO command line options */
  ierr = TaoSetFromOptions(tao);CHKERRQ(ierr);

  /* Solve the problem */
  ierr = TaoSetTolerances(tao, 1.e-5, 0.0, 0.0);CHKERRQ(ierr);
  ierr = TaoSetMaximumIterations(tao, 5);CHKERRQ(ierr);
  ierr = TaoLMVMRecycle(tao, PETSC_TRUE);CHKERRQ(ierr);
  reason = TAO_CONTINUE_ITERATING;
  while (reason != TAO_CONVERGED_GATOL) {
    ierr = TaoSolve(tao);CHKERRQ(ierr);
    ierr = TaoGetConvergedReason(tao, &reason);CHKERRQ(ierr);
    ierr = TaoGetIterationNumber(tao, &its);CHKERRQ(ierr);
    recycled_its += its;
    ierr = PetscPrintf(PETSC_COMM_SELF, "-----------------------\n");CHKERRQ(ierr);
  }
  
  /* Disable recycling and solve again! */
  ierr = TaoSetMaximumIterations(tao, 100);CHKERRQ(ierr);
  ierr = TaoLMVMRecycle(tao, PETSC_FALSE);CHKERRQ(ierr);
  ierr = VecSet(x, zero);CHKERRQ(ierr);
  ierr = TaoSolve(tao);CHKERRQ(ierr);
  ierr = TaoGetConvergedReason(tao, &reason);CHKERRQ(ierr);
  if (reason != TAO_CONVERGED_GATOL) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_NOT_CONVERGED, "Solution failed to converge!");
  ierr = TaoGetIterationNumber(tao, &oneshot_its);CHKERRQ(ierr);
  ierr = PetscPrintf(PETSC_COMM_SELF, "-----------------------\n");CHKERRQ(ierr);
  ierr = PetscPrintf(PETSC_COMM_SELF, "recycled its: %D | oneshot its: %D\n", recycled_its, oneshot_its);CHKERRQ(ierr);
  if (recycled_its != oneshot_its) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_NOT_CONVERGED, "LMVM recycling does not work!");
  
  ierr = TaoDestroy(&tao);CHKERRQ(ierr);
  ierr = VecDestroy(&x);CHKERRQ(ierr);
  ierr = MatDestroy(&H);CHKERRQ(ierr);

  ierr = PetscFinalize();
  return ierr;
}

int main(int argc,char **args)
{
  Vec            x,b,u;      /* approx solution, RHS, exact solution */
  Mat            A;            /* linear system matrix */
  KSP            ksp;         /* KSP context */
  PetscErrorCode ierr;
  PetscInt       n = 10,its, dim,p = 1,use_random;
  PetscScalar    none = -1.0,pfive = 0.5;
  PetscReal      norm;
  PetscRandom    rctx;
  TestType       type;
  PetscBool      flg;

  PetscInitialize(&argc,&args,(char *)0,help);
  ierr = PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(PETSC_NULL,"-p",&p,PETSC_NULL);CHKERRQ(ierr);
  switch (p) {
    case 1:  type = TEST_1;      dim = n;   break;
    case 2:  type = TEST_2;      dim = n;   break;
    case 3:  type = TEST_3;      dim = n;   break;
    case 4:  type = HELMHOLTZ_1; dim = n*n; break;
    case 5:  type = HELMHOLTZ_2; dim = n*n; break;
    default: type = TEST_1;      dim = n;
  }

  /* Create vectors */
  ierr = VecCreate(PETSC_COMM_WORLD,&x);CHKERRQ(ierr);
  ierr = VecSetSizes(x,PETSC_DECIDE,dim);CHKERRQ(ierr);
  ierr = VecSetFromOptions(x);CHKERRQ(ierr);
  ierr = VecDuplicate(x,&b);CHKERRQ(ierr);
  ierr = VecDuplicate(x,&u);CHKERRQ(ierr);

  use_random = 1;
  flg        = PETSC_FALSE;
  ierr = PetscOptionsGetBool(PETSC_NULL,"-norandom",&flg,PETSC_NULL);CHKERRQ(ierr);
  if (flg) {
    use_random = 0;
    ierr = VecSet(u,pfive);CHKERRQ(ierr);
  } else {
    ierr = PetscRandomCreate(PETSC_COMM_WORLD,&rctx);CHKERRQ(ierr);
    ierr = PetscRandomSetFromOptions(rctx);CHKERRQ(ierr);
    ierr = VecSetRandom(u,rctx);CHKERRQ(ierr);
  }

  /* Create and assemble matrix */
  ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
  ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,dim,dim);CHKERRQ(ierr);
  ierr = MatSetFromOptions(A);CHKERRQ(ierr);
  ierr = MatSetUp(A);CHKERRQ(ierr);
  ierr = FormTestMatrix(A,n,type);CHKERRQ(ierr);
  ierr = MatMult(A,u,b);CHKERRQ(ierr);
  flg  = PETSC_FALSE;
  ierr = PetscOptionsGetBool(PETSC_NULL,"-printout",&flg,PETSC_NULL);CHKERRQ(ierr);
  if (flg) {
    ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
    ierr = VecView(u,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
    ierr = VecView(b,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  }

  /* Create KSP context; set operators and options; solve linear system */
  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp);CHKERRQ(ierr);
  ierr = KSPSetOperators(ksp,A,A,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr);
  ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
  ierr = KSPSolve(ksp,b,x);CHKERRQ(ierr);
  ierr = KSPView(ksp,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);

  /* Check error */
  ierr = VecAXPY(x,none,u);CHKERRQ(ierr);
  ierr  = VecNorm(x,NORM_2,&norm);CHKERRQ(ierr);
  ierr = KSPGetIterationNumber(ksp,&its);CHKERRQ(ierr);
  ierr = PetscPrintf(PETSC_COMM_WORLD,"Norm of error %G,Iterations %D\n",norm,its);CHKERRQ(ierr);

  /* Free work space */
  ierr = VecDestroy(&x);CHKERRQ(ierr); ierr = VecDestroy(&u);CHKERRQ(ierr);
  ierr = VecDestroy(&b);CHKERRQ(ierr); ierr = MatDestroy(&A);CHKERRQ(ierr);
  if (use_random) {ierr = PetscRandomDestroy(&rctx);CHKERRQ(ierr);}
  ierr = KSPDestroy(&ksp);CHKERRQ(ierr);
  ierr = PetscFinalize();
  return 0;
}

int main(int argc,char **args)
{
  Vec            x1,b1,x2,b2; /* solution and RHS vectors for systems #1 and #2 */
  Vec            u;              /* exact solution vector */
  Mat            C1,C2;         /* matrices for systems #1 and #2 */
  KSP            ksp1,ksp2;   /* KSP contexts for systems #1 and #2 */
  PetscInt       ntimes = 3;     /* number of times to solve the linear systems */
  PetscLogEvent  CHECK_ERROR;    /* event number for error checking */
  PetscInt       ldim,low,high,iglobal,Istart,Iend,Istart2,Iend2;
  PetscInt       Ii,J,i,j,m = 3,n = 2,its,t;
  PetscErrorCode ierr;
  PetscBool      flg = PETSC_FALSE;
  PetscScalar    v;
  PetscMPIInt    rank,size;
#if defined (PETSC_USE_LOG)
  PetscLogStage  stages[3];
#endif

  PetscInitialize(&argc,&args,(char *)0,help);
  ierr = PetscOptionsGetInt(PETSC_NULL,"-m",&m,PETSC_NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(PETSC_NULL,"-t",&ntimes,PETSC_NULL);CHKERRQ(ierr);
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
  n = 2*size;

  /*
     Register various stages for profiling
  */
  ierr = PetscLogStageRegister("Prelim setup",&stages[0]);CHKERRQ(ierr);
  ierr = PetscLogStageRegister("Linear System 1",&stages[1]);CHKERRQ(ierr);
  ierr = PetscLogStageRegister("Linear System 2",&stages[2]);CHKERRQ(ierr);

  /*
     Register a user-defined event for profiling (error checking).
  */
  CHECK_ERROR = 0;
  ierr = PetscLogEventRegister("Check Error",KSP_CLASSID,&CHECK_ERROR);CHKERRQ(ierr);

  /* - - - - - - - - - - - - Stage 0: - - - - - - - - - - - - - -
                        Preliminary Setup
     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

  ierr = PetscLogStagePush(stages[0]);CHKERRQ(ierr);

  /*
     Create data structures for first linear system.
      - Create parallel matrix, specifying only its global dimensions.
        When using MatCreate(), the matrix format can be specified at
        runtime. Also, the parallel partitioning of the matrix is
        determined by PETSc at runtime.
      - Create parallel vectors.
        - When using VecSetSizes(), we specify only the vector's global
          dimension; the parallel partitioning is determined at runtime.
        - Note: We form 1 vector from scratch and then duplicate as needed.
  */
  ierr = MatCreate(PETSC_COMM_WORLD,&C1);CHKERRQ(ierr);
  ierr = MatSetSizes(C1,PETSC_DECIDE,PETSC_DECIDE,m*n,m*n);CHKERRQ(ierr);
  ierr = MatSetFromOptions(C1);CHKERRQ(ierr);
  ierr = MatSetUp(C1);CHKERRQ(ierr);
  ierr = MatGetOwnershipRange(C1,&Istart,&Iend);CHKERRQ(ierr);
  ierr = VecCreate(PETSC_COMM_WORLD,&u);CHKERRQ(ierr);
  ierr = VecSetSizes(u,PETSC_DECIDE,m*n);CHKERRQ(ierr);
  ierr = VecSetFromOptions(u);CHKERRQ(ierr);
  ierr = VecDuplicate(u,&b1);CHKERRQ(ierr);
  ierr = VecDuplicate(u,&x1);CHKERRQ(ierr);

  /*
     Create first linear solver context.
     Set runtime options (e.g., -pc_type <type>).
     Note that the first linear system uses the default option
     names, while the second linear systme uses a different
     options prefix.
  */
  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp1);CHKERRQ(ierr);
  ierr = KSPSetFromOptions(ksp1);CHKERRQ(ierr);

  /*
     Set user-defined monitoring routine for first linear system.
  */
  ierr = PetscOptionsGetBool(PETSC_NULL,"-my_ksp_monitor",&flg,PETSC_NULL);CHKERRQ(ierr);
  if (flg) {ierr = KSPMonitorSet(ksp1,MyKSPMonitor,PETSC_NULL,0);CHKERRQ(ierr);}

  /*
     Create data structures for second linear system.
  */
  ierr = MatCreate(PETSC_COMM_WORLD,&C2);CHKERRQ(ierr);
  ierr = MatSetSizes(C2,PETSC_DECIDE,PETSC_DECIDE,m*n,m*n);CHKERRQ(ierr);
  ierr = MatSetFromOptions(C2);CHKERRQ(ierr);
  ierr = MatSetUp(C2);CHKERRQ(ierr);
  ierr = MatGetOwnershipRange(C2,&Istart2,&Iend2);CHKERRQ(ierr);
  ierr = VecDuplicate(u,&b2);CHKERRQ(ierr);
  ierr = VecDuplicate(u,&x2);CHKERRQ(ierr);

  /*
     Create second linear solver context
  */
  ierr = KSPCreate(PETSC_COMM_WORLD,&ksp2);CHKERRQ(ierr);

  /*
     Set different options prefix for second linear system.
//.........这里部分代码省略.........

int main(int argc,char **args)
{
  Vec            x,b,u;      /* approx solution, RHS, exact solution */
  Mat            A;            /* linear system matrix */
  KSP            ksp;         /* linear solver context */
  PetscReal      norm;         /* norm of solution error */
  PetscInt       dim,i,j,Ii,J,Istart,Iend,n = 6,its,use_random;
  PetscErrorCode ierr;
  PetscScalar    v,none = -1.0,sigma2,pfive = 0.5,*xa;
  PetscRandom    rctx;
  PetscReal      h2,sigma1 = 100.0;
  PetscBool      flg = PETSC_FALSE;
  PetscScalar    a   = 1.0+PETSC_i;

  PetscInitialize(&argc,&args,(char*)0,help);
#if !defined(PETSC_USE_COMPLEX)
  SETERRQ(PETSC_COMM_WORLD,1,"This example requires complex numbers");
#endif

  a=1.0+PETSC_i;
  printf("%g+%gi\n",(double)PetscRealPart(a),(double)PetscImaginaryPart(a));

  ierr = PetscOptionsGetReal(NULL,"-sigma1",&sigma1,NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(NULL,"-n",&n,NULL);CHKERRQ(ierr);
  dim  = n*n;

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
         Compute the matrix and right-hand-side vector that define
         the linear system, Ax = b.
     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  /*
     Create parallel matrix, specifying only its global dimensions.
     When using MatCreate(), the matrix format can be specified at
     runtime. Also, the parallel partitioning of the matrix is
     determined by PETSc at runtime.
  */
  ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
  ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,dim,dim);CHKERRQ(ierr);
  ierr = MatSetFromOptions(A);CHKERRQ(ierr);
  ierr = MatSetUp(A);CHKERRQ(ierr);

  /*
     Currently, all PETSc parallel matrix formats are partitioned by
     contiguous chunks of rows across the processors.  Determine which
     rows of the matrix are locally owned.
  */
  ierr = MatGetOwnershipRange(A,&Istart,&Iend);CHKERRQ(ierr);

  /*
     Set matrix elements in parallel.
      - Each processor needs to insert only elements that it owns
        locally (but any non-local elements will be sent to the
        appropriate processor during matrix assembly).
      - Always specify global rows and columns of matrix entries.
  */

  ierr = PetscOptionsGetBool(NULL,"-norandom",&flg,NULL);CHKERRQ(ierr);
  if (flg) use_random = 0;
  else use_random = 1;
  if (use_random) {
    ierr = PetscRandomCreate(PETSC_COMM_WORLD,&rctx);CHKERRQ(ierr);
    ierr = PetscRandomSetFromOptions(rctx);CHKERRQ(ierr);
    ierr = PetscRandomSetInterval(rctx,0.0,PETSC_i);CHKERRQ(ierr);
  } else {
    sigma2 = 10.0*PETSC_i;
  }
  h2 = 1.0/((n+1)*(n+1));
  for (Ii=Istart; Ii<Iend; Ii++) {
    v = -1.0; i = Ii/n; j = Ii - i*n;
    if (i>0) {
      J = Ii-n; ierr = MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);
    }
    if (i<n-1) {
      J = Ii+n; ierr = MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);
    }
    if (j>0) {
      J = Ii-1; ierr = MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);
    }
    if (j<n-1) {
      J = Ii+1; ierr = MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);
    }
    if (use_random) {ierr = PetscRandomGetValue(rctx,&sigma2);CHKERRQ(ierr);}
    v    = 4.0 - sigma1*h2 + sigma2*h2;
    ierr = MatSetValues(A,1,&Ii,1,&Ii,&v,ADD_VALUES);CHKERRQ(ierr);
  }
  if (use_random) {ierr = PetscRandomDestroy(&rctx);CHKERRQ(ierr);}

  /*
     Assemble matrix, using the 2-step process:
       MatAssemblyBegin(), MatAssemblyEnd()
     Computations can be done while messages are in transition
     by placing code between these two statements.
  */
  ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);

  /*
     Create parallel vectors.
      - When using VecCreate(), VecSetSizes() and VecSetFromOptions(),
      we specify only the vector's global
//.........这里部分代码省略.........

int main(int argc,char **argv)
{
  PetscMPIInt      rank;
  PetscInt         M=8,dof=1,stencil_width=1,i,start,end,P=5,N = 6,m=PETSC_DECIDE,n=PETSC_DECIDE,p=PETSC_DECIDE,pt = 0,st = 0;
  PetscErrorCode   ierr;