Actual source code: ex1.c
petsc-3.4.2 2013-07-02

      2: static char help[] = "Basic vector routines.\n\n";

      4: /*T
      5:    Concepts: vectors^basic routines;
      6:    Processors: n
      7: T*/

      9: /*
     10:   Include "petscvec.h" so that we can use vectors.  Note that this file
     11:   automatically includes:
     12:      petscsys.h       - base PETSc routines   petscis.h     - index sets
     13:      petscviewer.h - viewers
     14: */

     16: #include <petscvec.h>

     20: int main(int argc,char **argv)
     21: {
     22:   Vec            x,y,w;               /* vectors */
     23:   Vec            *z;                    /* array of vectors */
     24:   PetscReal      norm,v,v1,v2,maxval;
     25:   PetscInt       n = 20,maxind;
     27:   PetscScalar    one = 1.0,two = 2.0,three = 3.0,dots[3],dot;

     29:   PetscInitialize(&argc,&argv,(char*)0,help);

     32:   PetscOptionsGetInt(NULL,"-n",&n,NULL);

     34:   /*
     35:      Create a vector, specifying only its global dimension.
     36:      When using VecCreate(), VecSetSizes() and VecSetFromOptions(), the vector format
     37:      (currently parallel, shared, or sequential) is determined at runtime.  Also, the
     38:      parallel partitioning of the vector is determined by PETSc at runtime.

     40:      Routines for creating particular vector types directly are:
     41:         VecCreateSeq() - uniprocessor vector
     42:         VecCreateMPI() - distributed vector, where the user can
     43:                          determine the parallel partitioning
     44:         VecCreateShared() - parallel vector that uses shared memory
     45:                             (available only on the SGI); otherwise,
     46:                             is the same as VecCreateMPI()

     48:      With VecCreate(), VecSetSizes() and VecSetFromOptions() the option -vec_type mpi or
     49:      -vec_type shared causes the particular type of vector to be formed.
     50: y

     52:   */

     54:   VecCreate(PETSC_COMM_WORLD,&x);
     55:   VecSetSizes(x,PETSC_DECIDE,n);
     56:   VecSetFromOptions(x);
     57:   /*
     58:      Duplicate some work vectors (of the same format and
     59:      partitioning as the initial vector).
     60:   */
     61:   VecDuplicate(x,&y);
     62:   VecDuplicate(x,&w);

     64:   /*
     65:      Duplicate more work vectors (of the same format and
     66:      partitioning as the initial vector).  Here we duplicate
     67:      an array of vectors, which is often more convenient than
     68:      duplicating individual ones.
     69:   */
     70:   VecDuplicateVecs(x,3,&z);
     71:   /*
     72:      Set the vectors to entries to a constant value.
     73:   */
     74:   VecSet(x,one);
     75:   VecSet(y,two);
     76:   VecSet(z[0],one);
     77:   VecSet(z[1],two);
     78:   VecSet(z[2],three);
     79:   /*
     80:      Demonstrate various basic vector routines.
     81:   */
     82:   MPI_Barrier(PETSC_COMM_WORLD);
     83:   VecDot(x,y,&dot);
     84:   VecMDot(x,3,z,dots);

     86:   /*
     87:      Note: If using a complex numbers version of PETSc, then
     88:      PETSC_USE_COMPLEX is defined in the makefiles; otherwise,
     89:      (when using real numbers) it is undefined.
     90:   */

     92:   PetscPrintf(PETSC_COMM_WORLD,"Vector length %D\n",n);
     93:   VecMax(x,&maxind,&maxval);
     94:   PetscPrintf(PETSC_COMM_WORLD,"VecMax %g, VecInd %D\n",(double)maxval,maxind);

     96:   VecMin(x,&maxind,&maxval);
     97:   PetscPrintf(PETSC_COMM_WORLD,"VecMin %g, VecInd %D\n",(double)maxval,maxind);
     98:   PetscPrintf(PETSC_COMM_WORLD,"All other values should be near zero\n");

    101:   VecScale(x,two);
    102:   VecNorm(x,NORM_2,&norm);
    103:   v    = norm-2.0*sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    104:   PetscPrintf(PETSC_COMM_WORLD,"VecScale %g\n",(double)v);

    107:   VecCopy(x,w);
    108:   VecNorm(w,NORM_2,&norm);
    109:   v    = norm-2.0*sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    110:   PetscPrintf(PETSC_COMM_WORLD,"VecCopy  %g\n",(double)v);

    112:   VecAXPY(y,three,x);
    113:   VecNorm(y,NORM_2,&norm);
    114:   v    = norm-8.0*sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    115:   PetscPrintf(PETSC_COMM_WORLD,"VecAXPY %g\n",(double)v);

    117:   VecAYPX(y,two,x);
    118:   VecNorm(y,NORM_2,&norm);
    119:   v    = norm-18.0*sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    120:   PetscPrintf(PETSC_COMM_WORLD,"VecAYPX %g\n",(double)v);

    122:   VecSwap(x,y);
    123:   VecNorm(y,NORM_2,&norm);
    124:   v    = norm-2.0*sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    125:   PetscPrintf(PETSC_COMM_WORLD,"VecSwap  %g\n",(double)v);
    126:   VecNorm(x,NORM_2,&norm);
    127:   v = norm-18.0*sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    128:   PetscPrintf(PETSC_COMM_WORLD,"VecSwap  %g\n",(double)v);

    130:   VecWAXPY(w,two,x,y);
    131:   VecNorm(w,NORM_2,&norm);
    132:   v    = norm-38.0*sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    133:   PetscPrintf(PETSC_COMM_WORLD,"VecWAXPY %g\n",(double)v);

    135:   VecPointwiseMult(w,y,x);
    136:   VecNorm(w,NORM_2,&norm);
    137:   v    = norm-36.0*sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    138:   PetscPrintf(PETSC_COMM_WORLD,"VecPointwiseMult %g\n",(double)v);

    140:   VecPointwiseDivide(w,x,y);
    141:   VecNorm(w,NORM_2,&norm);
    142:   v    = norm-9.0*sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    143:   PetscPrintf(PETSC_COMM_WORLD,"VecPointwiseDivide %g\n",(double)v);

    145:   dots[0] = one;
    146:   dots[1] = three;
    147:   dots[2] = two;

    149:   VecSet(x,one);
    150:   VecMAXPY(x,3,dots,z);
    151:   VecNorm(z[0],NORM_2,&norm);
    152:   v    = norm-sqrt((double)n); if (v > -PETSC_SMALL && v < PETSC_SMALL) v = 0.0;
    153:   VecNorm(z[1],NORM_2,&norm);
    154:   v1   = norm-2.0*sqrt((double)n); if (v1 > -PETSC_SMALL && v1 < PETSC_SMALL) v1 = 0.0;
    155:   VecNorm(z[2],NORM_2,&norm);
    156:   v2   = norm-3.0*sqrt((double)n); if (v2 > -PETSC_SMALL && v2 < PETSC_SMALL) v2 = 0.0;
    157:   PetscPrintf(PETSC_COMM_WORLD,"VecMAXPY %g %g %g \n",(double)v,(double)v1,(double)v2);

    159:   /*
    160:      Free work space.  All PETSc objects should be destroyed when they
    161:      are no longer needed.
    162:   */
    163:   VecDestroy(&x);
    164:   VecDestroy(&y);
    165:   VecDestroy(&w);
    166:   VecDestroyVecs(3,&z);
    167:   PetscFinalize();
    168:   return 0;
    169: }