当前位置: 首页>>代码示例>>C++>>正文

C++ PetscRandomSetFromOptions函数代码示例

本文整理汇总了C++中PetscRandomSetFromOptions函数的典型用法代码示例。如果您正苦于以下问题：C++ PetscRandomSetFromOptions函数的具体用法？C++ PetscRandomSetFromOptions怎么用？C++ PetscRandomSetFromOptions使用的例子？那么恭喜您, 这里精选的函数代码示例或许可以为您提供帮助。


在下文中一共展示了PetscRandomSetFromOptions函数的15个代码示例，这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞，您的评价将有助于系统推荐出更棒的C++代码示例。

示例1: SNESDiffParameterCreate_More

PetscErrorCode SNESDiffParameterCreate_More(SNES snes,Vec x,void **outneP)
{
  DIFFPAR_MORE   *neP;
  Vec            w;
  PetscRandom    rctx;  /* random number generator context */
  PetscErrorCode ierr;
  PetscBool      flg;
  char           noise_file[PETSC_MAX_PATH_LEN];

  PetscFunctionBegin;

  ierr = PetscNewLog(snes,DIFFPAR_MORE,&neP);CHKERRQ(ierr);
  neP->function_count = 0;
  neP->fnoise_min     = 1.0e-20;
  neP->hopt_min       = 1.0e-8;
  neP->h_first_try    = 1.0e-3;
  neP->fnoise_resets  = 0;
  neP->hopt_resets    = 0;

  /* Create work vectors */
  ierr = VecDuplicateVecs(x,3,&neP->workv);CHKERRQ(ierr);
  w = neP->workv[0];

  /* Set components of vector w to random numbers */
  ierr = PetscRandomCreate(((PetscObject)snes)->comm,&rctx);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(rctx);CHKERRQ(ierr);
  ierr = VecSetRandom(w,rctx);CHKERRQ(ierr);
  ierr = PetscRandomDestroy(&rctx);CHKERRQ(ierr);

  /* Open output file */
  ierr = PetscOptionsGetString(((PetscObject)snes)->prefix,"-snes_mf_noise_file",noise_file,PETSC_MAX_PATH_LEN,&flg);CHKERRQ(ierr);
  if (flg) neP->fp = fopen(noise_file,"w"); 
  else     neP->fp = fopen("noise.out","w"); 
  if (!neP->fp) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN,"Cannot open file");
  ierr = PetscInfo(snes,"Creating Jorge's differencing parameter context\n");CHKERRQ(ierr);

  *outneP = neP;
  PetscFunctionReturn(0);
}
开发者ID:Kun-Qu，项目名称:petsc，代码行数:39，代码来源:snesnoise.c

示例2: Ensure

Ensure(FFT, i_transform_is_inverse_of_transform)
{
    Vec v;
    PetscInt dim = 5;

    DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,dim,2,1,NULL,&da);

    DMCreateGlobalVector(da, &v);
    PetscRandom rdm;
    PetscRandomCreate(PETSC_COMM_WORLD, &rdm);
    PetscRandomSetFromOptions(rdm);
    VecSetRandom(v, rdm);
    PetscRandomDestroy(&rdm);
    Vec vIn;
    VecDuplicate(v, &vIn);
    VecCopy(v, vIn);

    scFftCreate(da, &fft);

    Vec x, y, z;
    scFftCreateVecsFFTW(fft, &x, &y, &z);

    int component = 0;
    scFftTransform(fft, v, component, y);
    scFftITransform(fft, v, component, y);
    VecAXPY(v, -dim, vIn);
    PetscReal error;
    VecStrideNorm(v, component, NORM_INFINITY, &error);
    assert_that(error < 1.0e-6, is_true);

    VecDestroy(&x);
    VecDestroy(&y);
    VecDestroy(&z);
    VecDestroy(&v);
    VecDestroy(&vIn);
    scFftDestroy(&fft);
    DMDestroy(&da);
}
开发者ID:ocramz，项目名称:SuperContinuum，代码行数:38，代码来源:test_FFT.c

示例3: main

PetscInt main(PetscInt argc,char **args)
{
  typedef enum {RANDOM, CONSTANT, TANH, NUM_FUNCS} FuncType;
  const char    *funcNames[NUM_FUNCS] = {"random", "constant", "tanh"};
  PetscMPIInt    size;
  PetscInt       n = 10,N,Ny,ndim=4,dim[4],DIM,i;
  Vec            x,y,z;
  PetscScalar    s;
  PetscRandom    rdm;
  PetscReal      enorm;
  PetscInt       func=RANDOM;
  FuncType       function = RANDOM;
  PetscBool      view = PETSC_FALSE;
  PetscErrorCode ierr;
  PetscScalar    *x_array,*y_array,*z_array;
  fftw_plan      fplan,bplan;
  const ptrdiff_t N0 = 20, N1 = 20;

  ptrdiff_t alloc_local, local_n0, local_0_start;
  ierr = PetscInitialize(&argc,&args,(char *)0,help);CHKERRQ(ierr);
#if defined(PETSC_USE_COMPLEX)
  SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP, "This example requires real numbers");
#endif
  ierr = MPI_Comm_size(PETSC_COMM_WORLD, &size);CHKERRQ(ierr);
  alloc_local=fftw_mpi_local_size_2d(N0, N1, PETSC_COMM_WORLD,
                                              &local_n0, &local_0_start);

  if (size != 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP, "This is a uniprocessor example only!");
  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, PETSC_NULL, "FFTW Options", "ex142");CHKERRQ(ierr);
    ierr = PetscOptionsEList("-function", "Function type", "ex142", funcNames, NUM_FUNCS, funcNames[function], &func, PETSC_NULL);CHKERRQ(ierr);
    ierr = PetscOptionsBool("-vec_view draw", "View the functions", "ex112", view, &view, PETSC_NULL);CHKERRQ(ierr);
    function = (FuncType) func;
  ierr = PetscOptionsEnd();CHKERRQ(ierr);

  for (DIM = 0; DIM < ndim; DIM++){
    dim[DIM]  = n; /* size of real space vector in DIM-dimension */
  }
  ierr = PetscRandomCreate(PETSC_COMM_SELF, &rdm);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(rdm);CHKERRQ(ierr);

  for (DIM = 1; DIM < 5; DIM++){
    /* create vectors of length N=dim[0]*dim[1]* ...*dim[DIM-1] */
    /*----------------------------------------------------------*/
    N = Ny = 1;
    for (i = 0; i < DIM-1; i++) {
      N *= dim[i];
    }
    Ny = N; Ny *= 2*(dim[DIM-1]/2 + 1); /* add padding elements to output vector y */
    N *= dim[DIM-1];


    ierr = PetscPrintf(PETSC_COMM_SELF, "\n %d-D: FFTW on vector of size %d \n",DIM,N);CHKERRQ(ierr);
    ierr = VecCreateSeq(PETSC_COMM_SELF,N,&x);CHKERRQ(ierr);
    ierr = PetscObjectSetName((PetscObject) x, "Real space vector");CHKERRQ(ierr);

    ierr = VecCreateSeq(PETSC_COMM_SELF,Ny,&y);CHKERRQ(ierr);
    ierr = PetscObjectSetName((PetscObject) y, "Frequency space vector");CHKERRQ(ierr);

    ierr = VecDuplicate(x,&z);CHKERRQ(ierr);
    ierr = PetscObjectSetName((PetscObject) z, "Reconstructed vector");CHKERRQ(ierr);


    /* Set fftw plan                    */
    /*----------------------------------*/
    ierr = VecGetArray(x,&x_array);CHKERRQ(ierr);
    ierr = VecGetArray(y,&y_array);CHKERRQ(ierr);
    ierr = VecGetArray(z,&z_array);CHKERRQ(ierr);

    unsigned int flags = FFTW_ESTIMATE; //or FFTW_MEASURE
    /* The data in the in/out arrays is overwritten during FFTW_MEASURE planning, so such planning
     should be done before the input is initialized by the user. */
    printf("DIM: %d, N %d, Ny %d\n",DIM,N,Ny);

    switch (DIM){
    case 1:
      fplan = fftw_plan_dft_r2c_1d(dim[0], (double *)x_array, (fftw_complex*)y_array, flags);
      bplan = fftw_plan_dft_c2r_1d(dim[0], (fftw_complex*)y_array, (double *)z_array, flags);
      break;
    case 2:
      fplan = fftw_plan_dft_r2c_2d(dim[0],dim[1],(double *)x_array, (fftw_complex*)y_array,flags);
      bplan = fftw_plan_dft_c2r_2d(dim[0],dim[1],(fftw_complex*)y_array,(double *)z_array,flags);
      break;
    case 3:
      fplan = fftw_plan_dft_r2c_3d(dim[0],dim[1],dim[2],(double *)x_array, (fftw_complex*)y_array,flags);
      bplan = fftw_plan_dft_c2r_3d(dim[0],dim[1],dim[2],(fftw_complex*)y_array,(double *)z_array,flags);
      break;
    default:
      fplan = fftw_plan_dft_r2c(DIM,dim,(double *)x_array, (fftw_complex*)y_array,flags);
      bplan = fftw_plan_dft_c2r(DIM,dim,(fftw_complex*)y_array,(double *)z_array,flags);
      break;
    }

    ierr = VecRestoreArray(x,&x_array);CHKERRQ(ierr);
    ierr = VecRestoreArray(y,&y_array);CHKERRQ(ierr);
    ierr = VecRestoreArray(z,&z_array);CHKERRQ(ierr);

    /* Initialize Real space vector x:
       The data in the in/out arrays is overwritten during FFTW_MEASURE planning, so planning
       should be done before the input is initialized by the user.
    --------------------------------------------------------*/
//.........这里部分代码省略.........
开发者ID:erdc-cm，项目名称:petsc-dev，代码行数:101，代码来源:ex142.c

示例4: TestTriangle

PetscErrorCode TestTriangle(MPI_Comm comm, PetscBool interpolate, PetscBool transform)
{
  DM             dm;
  PetscRandom    r, ang, ang2;
  PetscInt       dim, t;
  PetscErrorCode ierr;

  PetscFunctionBegin;
  /* Create reference triangle */
  dim  = 2;
  ierr = DMCreate(comm, &dm);CHKERRQ(ierr);
  ierr = PetscObjectSetName((PetscObject) dm, "triangle");CHKERRQ(ierr);
  ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr);
  ierr = DMPlexSetDimension(dm, dim);CHKERRQ(ierr);
  {
    PetscInt    numPoints[2]        = {3, 1};
    PetscInt    coneSize[4]         = {3, 0, 0, 0};
    PetscInt    cones[3]            = {1, 2, 3};
    PetscInt    coneOrientations[3] = {0, 0, 0};
    PetscScalar vertexCoords[6]     = {-1.0, -1.0, 1.0, -1.0, -1.0, 1.0};

    ierr = DMPlexCreateFromDAG(dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr);
    if (interpolate) {
      DM idm;

      ierr = DMPlexInterpolate(dm, &idm);CHKERRQ(ierr);
      ierr = PetscObjectSetName((PetscObject) idm, "triangle");CHKERRQ(ierr);
      ierr = DMPlexCopyCoordinates(dm, idm);CHKERRQ(ierr);
      ierr = DMDestroy(&dm);CHKERRQ(ierr);
      dm   = idm;
    }
    ierr = DMSetFromOptions(dm);CHKERRQ(ierr);
  }
  /* Check reference geometry: determinant is scaled by reference volume (2.0) */
  {
    PetscReal v0Ex[2]       = {-1.0, -1.0};
    PetscReal JEx[4]        = {1.0, 0.0, 0.0, 1.0};
    PetscReal invJEx[4]     = {1.0, 0.0, 0.0, 1.0};
    PetscReal detJEx        = 1.0;
    PetscReal centroidEx[2] = {-0.333333333333, -0.333333333333};
    PetscReal normalEx[2]   = {0.0, 0.0};
    PetscReal volEx         = 2.0;

    ierr = CheckFEMGeometry(dm, 0, dim, v0Ex, JEx, invJEx, detJEx);CHKERRQ(ierr);
    if (interpolate) {ierr = CheckFVMGeometry(dm, 0, dim, centroidEx, normalEx, volEx);CHKERRQ(ierr);}
  }
  /* Check random triangles: rotate, scale, then translate */
  if (transform) {
    ierr = PetscRandomCreate(PETSC_COMM_SELF, &r);CHKERRQ(ierr);
    ierr = PetscRandomSetFromOptions(r);CHKERRQ(ierr);
    ierr = PetscRandomSetInterval(r, 0.0, 10.0);CHKERRQ(ierr);
    ierr = PetscRandomCreate(PETSC_COMM_SELF, &ang);CHKERRQ(ierr);
    ierr = PetscRandomSetFromOptions(ang);CHKERRQ(ierr);
    ierr = PetscRandomSetInterval(ang, 0.0, 2*PETSC_PI);CHKERRQ(ierr);
    for (t = 0; t < 100; ++t) {
      PetscScalar vertexCoords[6] = {-1.0, -1.0, 1.0, -1.0, -1.0, 1.0}, trans[2];
      PetscReal   v0Ex[2]         = {-1.0, -1.0};
      PetscReal   JEx[4]          = {1.0, 0.0, 0.0, 1.0}, R[4], rot[2], rotM[4];
      PetscReal   invJEx[4]       = {1.0, 0.0, 0.0, 1.0};
      PetscReal   detJEx          = 1.0, scale, phi;
      PetscReal   centroidEx[2]   = {-0.333333333333, -0.333333333333};
      PetscReal   normalEx[2]     = {0.0, 0.0};
      PetscReal   volEx           = 2.0;
      PetscInt    d, e, f, p;

      ierr = PetscRandomGetValueReal(r, &scale);CHKERRQ(ierr);
      ierr = PetscRandomGetValueReal(ang, &phi);CHKERRQ(ierr);
      R[0] = cos(phi); R[1] = -sin(phi);
      R[2] = sin(phi); R[3] =  cos(phi);
      for (p = 0; p < 3; ++p) {
        for (d = 0; d < dim; ++d) {
          for (e = 0, rot[d] = 0.0; e < dim; ++e) {
            rot[d] += R[d*dim+e] * vertexCoords[p*dim+e];
          }
        }
        for (d = 0; d < dim; ++d) vertexCoords[p*dim+d] = rot[d];
      }
      for (d = 0; d < dim; ++d) {
        for (e = 0, rot[d] = 0.0; e < dim; ++e) {
          rot[d] += R[d*dim+e] * centroidEx[e];
        }
      }
      for (d = 0; d < dim; ++d) centroidEx[d] = rot[d];
      for (d = 0; d < dim; ++d) {
        for (e = 0; e < dim; ++e) {
          for (f = 0, rotM[d*dim+e] = 0.0; f < dim; ++f) {
            rotM[d*dim+e] += R[d*dim+f] * JEx[f*dim+e];
          }
        }
      }
      for (d = 0; d < dim; ++d) {
        for (e = 0; e < dim; ++e) {
          JEx[d*dim+e] = rotM[d*dim+e];
        }
      }
      for (d = 0; d < dim; ++d) {
        for (e = 0; e < dim; ++e) {
          for (f = 0, rotM[d*dim+e] = 0.0; f < dim; ++f) {
            rotM[d*dim+e] += invJEx[d*dim+f] * R[e*dim+f];
          }
//.........这里部分代码省略.........
开发者ID:feelpp，项目名称:debian-petsc，代码行数:101，代码来源:ex8.c

示例5: main

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;
}
开发者ID:Kun-Qu，项目名称:petsc，代码行数:80，代码来源:ex17.c

示例6: main

PetscInt main(PetscInt argc,char **args)
{
  typedef enum {RANDOM, CONSTANT, TANH, NUM_FUNCS} FuncType;
  const char    *funcNames[NUM_FUNCS] = {"random", "constant", "tanh"};
  Mat            A;
  PetscMPIInt    size;
  PetscInt       n = 10,N,ndim=4,dim[4],DIM,i;
  Vec            x,y,z;
  PetscScalar    s;
  PetscRandom    rdm;
  PetscReal      enorm;
  PetscInt       func;
  FuncType       function = RANDOM;
  PetscBool      view = PETSC_FALSE;
  PetscErrorCode ierr;

  ierr = PetscInitialize(&argc,&args,(char *)0,help);CHKERRQ(ierr);
#if !defined(PETSC_USE_COMPLEX)
  SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP, "This example requires complex numbers");
#endif
  ierr = MPI_Comm_size(PETSC_COMM_WORLD, &size);CHKERRQ(ierr);
  if (size != 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP, "This is a uniprocessor example only!");
  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, PETSC_NULL, "FFTW Options", "ex112");CHKERRQ(ierr);
    ierr = PetscOptionsEList("-function", "Function type", "ex112", funcNames, NUM_FUNCS, funcNames[function], &func, PETSC_NULL);CHKERRQ(ierr);
    ierr = PetscOptionsBool("-vec_view draw", "View the functions", "ex112", view, &view, PETSC_NULL);CHKERRQ(ierr);
    function = (FuncType) func;
  ierr = PetscOptionsEnd();CHKERRQ(ierr);

  for (DIM = 0; DIM < ndim; DIM++){
    dim[DIM] = n;  /* size of transformation in DIM-dimension */
  }
  ierr = PetscRandomCreate(PETSC_COMM_SELF, &rdm);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(rdm);CHKERRQ(ierr);

  for (DIM = 1; DIM < 5; DIM++){
    for (i = 0, N = 1; i < DIM; i++) N *= dim[i];
    ierr = PetscPrintf(PETSC_COMM_SELF, "\n %d-D: FFTW on vector of size %d \n",DIM,N);CHKERRQ(ierr);

    /* create FFTW object */
    ierr = MatCreateFFT(PETSC_COMM_SELF,DIM,dim,MATFFTW,&A);CHKERRQ(ierr);

    /* create vectors of length N=n^DIM */
    ierr = MatGetVecs(A,&x,&y);CHKERRQ(ierr);
    ierr = MatGetVecs(A,&z,PETSC_NULL);CHKERRQ(ierr);
    ierr = PetscObjectSetName((PetscObject) x, "Real space vector");CHKERRQ(ierr);
    ierr = PetscObjectSetName((PetscObject) y, "Frequency space vector");CHKERRQ(ierr);
    ierr = PetscObjectSetName((PetscObject) z, "Reconstructed vector");CHKERRQ(ierr);

    /* set values of space vector x */
    if (function == RANDOM) {
      ierr = VecSetRandom(x, rdm);CHKERRQ(ierr);
    } else if (function == CONSTANT) {
      ierr = VecSet(x, 1.0);CHKERRQ(ierr);
    } else if (function == TANH) {
      PetscScalar *a;
      ierr = VecGetArray(x, &a);CHKERRQ(ierr);
      for (i = 0; i < N; ++i) {
        a[i] = tanh((i - N/2.0)*(10.0/N));
      }
      ierr = VecRestoreArray(x, &a);CHKERRQ(ierr);
    }
    if (view) {ierr = VecView(x, PETSC_VIEWER_DRAW_WORLD);CHKERRQ(ierr);}

    /* apply FFTW_FORWARD and FFTW_BACKWARD several times on same x, y, and z */
    for (i=0; i<3; i++){
      ierr = MatMult(A,x,y);CHKERRQ(ierr);
      if (view && i == 0) {ierr = VecView(y, PETSC_VIEWER_DRAW_WORLD);CHKERRQ(ierr);}
      ierr = MatMultTranspose(A,y,z);CHKERRQ(ierr);

      /* compare x and z. FFTW computes an unnormalized DFT, thus z = N*x */
      s = 1.0/(PetscReal)N;
      ierr = VecScale(z,s);CHKERRQ(ierr);
      if (view && i == 0) {ierr = VecView(z, PETSC_VIEWER_DRAW_WORLD);CHKERRQ(ierr);}
      ierr = VecAXPY(z,-1.0,x);CHKERRQ(ierr);
      ierr = VecNorm(z,NORM_1,&enorm);CHKERRQ(ierr);
      if (enorm > 1.e-11){
        ierr = PetscPrintf(PETSC_COMM_SELF,"  Error norm of |x - z| %G\n",enorm);CHKERRQ(ierr);
      }
    }

    /* apply FFTW_FORWARD and FFTW_BACKWARD several times on different x */
    for (i=0; i<3; i++){
      ierr = VecDestroy(&x);CHKERRQ(ierr);
      ierr = VecCreateSeq(PETSC_COMM_SELF,N,&x);CHKERRQ(ierr);
      ierr = VecSetRandom(x, rdm);CHKERRQ(ierr);

      ierr = MatMult(A,x,y);CHKERRQ(ierr);
      ierr = MatMultTranspose(A,y,z);CHKERRQ(ierr);

      /* compare x and z. FFTW computes an unnormalized DFT, thus z = N*x */
      s = 1.0/(PetscReal)N;
      ierr = VecScale(z,s);CHKERRQ(ierr);
      if (view && i == 0) {ierr = VecView(z, PETSC_VIEWER_DRAW_WORLD);CHKERRQ(ierr);}
      ierr = VecAXPY(z,-1.0,x);CHKERRQ(ierr);
      ierr = VecNorm(z,NORM_1,&enorm);CHKERRQ(ierr);
      if (enorm > 1.e-11){
        ierr = PetscPrintf(PETSC_COMM_SELF,"  Error norm of new |x - z| %G\n",enorm);CHKERRQ(ierr);
      }
    }

//.........这里部分代码省略.........
开发者ID:erdc-cm，项目名称:petsc-dev，代码行数:101，代码来源:ex112.c

示例7: main

PetscInt main(PetscInt argc,char **args)
{
  PetscErrorCode  ierr;
  PetscMPIInt     rank,size;
  PetscInt        N0=3,N1=3,N2=3,N=N0*N1*N2;
  PetscRandom     rdm;
  PetscScalar     a;
  PetscReal       enorm;
  Vec             x,y,z,input,output;
  PetscBool       view=PETSC_FALSE,use_interface=PETSC_TRUE;
  Mat             A;
  PetscInt        DIM, dim[3],vsize;
  PetscReal       fac;

  ierr = PetscInitialize(&argc,&args,(char *)0,help);CHKERRQ(ierr);
#if defined(PETSC_USE_COMPLEX)
 SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP, "This example requires real numbers");
#endif

  ierr = MPI_Comm_size(PETSC_COMM_WORLD, &size);CHKERRQ(ierr);
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank);CHKERRQ(ierr);


  ierr = PetscRandomCreate(PETSC_COMM_WORLD, &rdm);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(rdm);CHKERRQ(ierr);

  ierr = VecCreate(PETSC_COMM_WORLD,&input);CHKERRQ(ierr);
  ierr = VecSetSizes(input,PETSC_DECIDE,N);CHKERRQ(ierr);
  ierr = VecSetFromOptions(input);CHKERRQ(ierr);
  ierr = VecSetRandom(input,rdm);CHKERRQ(ierr);
  ierr = VecDuplicate(input,&output);
//  ierr = VecGetSize(input,&vsize);CHKERRQ(ierr);
//  printf("Size of the input Vector is %d\n",vsize);

  DIM = 3;
  dim[0] = N0; dim[1] = N1; dim[2] = N2;

  ierr = MatCreateFFT(PETSC_COMM_WORLD,DIM,dim,MATFFTW,&A);CHKERRQ(ierr);
  ierr = MatGetVecs(A,&x,&y);CHKERRQ(ierr);
  ierr = MatGetVecs(A,&z,PETSC_NULL);CHKERRQ(ierr);
  ierr = VecGetSize(y,&vsize);CHKERRQ(ierr);
  printf("The vector size from the main routine is %d\n",vsize);

  ierr = InputTransformFFT(A,input,x);CHKERRQ(ierr);
  ierr = MatMult(A,x,y);CHKERRQ(ierr);
  ierr = MatMultTranspose(A,y,z);CHKERRQ(ierr);
  ierr = OutputTransformFFT(A,z,output);CHKERRQ(ierr);

  fac = 1.0/(PetscReal)N;
  ierr = VecScale(output,fac);CHKERRQ(ierr);

  ierr = VecAssemblyBegin(input);CHKERRQ(ierr);
  ierr = VecAssemblyEnd(input);CHKERRQ(ierr);
  ierr = VecAssemblyBegin(output);CHKERRQ(ierr);
  ierr = VecAssemblyEnd(output);CHKERRQ(ierr);

  ierr = VecView(input,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  ierr = VecView(output,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);

  ierr = VecAXPY(output,-1.0,input);CHKERRQ(ierr);
  ierr = VecNorm(output,NORM_1,&enorm);CHKERRQ(ierr);
//  if (enorm > 1.e-14){
    if (!rank)
      ierr = PetscPrintf(PETSC_COMM_SELF,"  Error norm of |x - z| %e\n",enorm);CHKERRQ(ierr);
//      }




// ierr = MatGetVecs(A,&z,PETSC_NULL);CHKERRQ(ierr);
//  printf("Vector size from ex148 %d\n",vsize);
//  ierr = PetscObjectSetName((PetscObject) x, "Real space vector");CHKERRQ(ierr);
//      ierr = PetscObjectSetName((PetscObject) y, "Frequency space vector");CHKERRQ(ierr);
//      ierr = PetscObjectSetName((PetscObject) z, "Reconstructed vector");CHKERRQ(ierr);

  ierr = PetscFinalize();
  return 0;

}
开发者ID:erdc-cm，项目名称:petsc-dev，代码行数:79，代码来源:ex149.c

示例8: main


//.........这里部分代码省略.........
  rnorm = PetscAbsReal(norm1-norm2)/norm2;
  if (rnorm > tol) {
    ierr = PetscPrintf(PETSC_COMM_SELF,"Error: MatNorm_INFINITY(), NormA=%16.14e NormsB=%16.14e\n",norm1,norm2);CHKERRQ(ierr);
  }

  /* Test MatGetInfo(), MatGetSize(), MatGetBlockSize() */
  ierr = MatGetInfo(A,MAT_LOCAL,&minfo1);CHKERRQ(ierr);
  ierr = MatGetInfo(sB,MAT_LOCAL,&minfo2);CHKERRQ(ierr);
  /*
  printf("matrix nonzeros (BAIJ format) = %d, allocated nonzeros= %d\n", (int)minfo1.nz_used,(int)minfo1.nz_allocated);
  printf("matrix nonzeros(SBAIJ format) = %d, allocated nonzeros= %d\n", (int)minfo2.nz_used,(int)minfo2.nz_allocated);
  */
  i  = (int) (minfo1.nz_used - minfo2.nz_used);
  j  = (int) (minfo1.nz_allocated - minfo2.nz_allocated);
  k1 = (int) (minfo1.nz_allocated - minfo1.nz_used);
  k2 = (int) (minfo2.nz_allocated - minfo2.nz_used);
  if (i < 0 || j < 0 || k1 < 0 || k2 < 0) {
    ierr = PetscPrintf(PETSC_COMM_SELF,"Error(compare A and sB): MatGetInfo()\n");CHKERRQ(ierr);
  }

  ierr = MatGetSize(A,&Ii,&J);CHKERRQ(ierr);
  ierr = MatGetSize(sB,&i,&j);CHKERRQ(ierr);
  if (i-Ii || j-J) {
    PetscPrintf(PETSC_COMM_SELF,"Error: MatGetSize()\n");CHKERRQ(ierr);
  }

  ierr = MatGetBlockSize(A, &Ii);CHKERRQ(ierr);
  ierr = MatGetBlockSize(sB, &i);CHKERRQ(ierr);
  if (i-Ii) {
    ierr = PetscPrintf(PETSC_COMM_SELF,"Error: MatGetBlockSize()\n");CHKERRQ(ierr);
  }

  ierr = PetscRandomCreate(PETSC_COMM_SELF,&rdm);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(rdm);CHKERRQ(ierr);
  ierr = VecCreateSeq(PETSC_COMM_SELF,n,&x);CHKERRQ(ierr);
  ierr = VecDuplicate(x,&s1);CHKERRQ(ierr);
  ierr = VecDuplicate(x,&s2);CHKERRQ(ierr);
  ierr = VecDuplicate(x,&y);CHKERRQ(ierr);
  ierr = VecDuplicate(x,&b);CHKERRQ(ierr);
  ierr = VecSetRandom(x,rdm);CHKERRQ(ierr);

  /* Test MatDiagonalScale(), MatGetDiagonal(), MatScale() */
#if !defined(PETSC_USE_COMPLEX)
  /* Scaling matrix with complex numbers results non-spd matrix,
     causing crash of MatForwardSolve() and MatBackwardSolve() */
  ierr = MatDiagonalScale(A,x,x);CHKERRQ(ierr);
  ierr = MatDiagonalScale(sB,x,x);CHKERRQ(ierr);
  ierr = MatMultEqual(A,sB,10,&equal);CHKERRQ(ierr);
  if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in MatDiagonalScale");

  ierr = MatGetDiagonal(A,s1);CHKERRQ(ierr);
  ierr = MatGetDiagonal(sB,s2);CHKERRQ(ierr);
  ierr = VecAXPY(s2,neg_one,s1);CHKERRQ(ierr);
  ierr = VecNorm(s2,NORM_1,&norm1);CHKERRQ(ierr);
  if (norm1>tol) {
    ierr = PetscPrintf(PETSC_COMM_SELF,"Error:MatGetDiagonal(), ||s1-s2||=%G\n",norm1);CHKERRQ(ierr);
  }

  {
    PetscScalar alpha=0.1;
    ierr = MatScale(A,alpha);CHKERRQ(ierr);
    ierr = MatScale(sB,alpha);CHKERRQ(ierr);
  }
#endif

  /* Test MatGetRowMaxAbs() */
开发者ID:ZJLi2013，项目名称:petsc，代码行数:67，代码来源:ex74.c

示例9: main

int main(int Argc,char **Args)
{
  PetscBool      flg;
  PetscInt       n   = -6;
  PetscScalar    rho = 1.0;
  PetscReal      h;
  PetscReal      beta = 1.0;
  DM             da;
  PetscRandom    rctx;
  PetscMPIInt    comm_size;
  Mat            H,HtH;
  PetscInt       x, y, xs, ys, xm, ym;
  PetscReal      r1, r2;
  PetscScalar    uxy1, uxy2;
  MatStencil     sxy, sxy_m;
  PetscScalar    val, valconj;
  Vec            b, Htb,xvec;
  KSP            kspmg;
  PC             pcmg;
  PetscErrorCode ierr;
  PetscInt       ix[1]   = {0};
  PetscScalar    vals[1] = {1.0};

  PetscInitialize(&Argc,&Args,(char*)0,help);
  ierr = PetscOptionsGetInt(NULL,"-size",&n,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetReal(NULL,"-beta",&beta,&flg);CHKERRQ(ierr);
  ierr = PetscOptionsGetScalar(NULL,"-rho",&rho,&flg);CHKERRQ(ierr);

  /* Set the fudge parameters, we scale the whole thing by 1/(2*h) later */
  h    = 1.;
  rho *= 1./(2.*h);

  /* Geometry info */
  ierr = DMDACreate2d(PETSC_COMM_WORLD, DMDA_BOUNDARY_PERIODIC,DMDA_BOUNDARY_PERIODIC, DMDA_STENCIL_STAR, n, n,
                      PETSC_DECIDE, PETSC_DECIDE, 2 /* this is the # of dof's */,
                      1, NULL, NULL, &da);CHKERRQ(ierr);

  /* Random numbers */
  ierr = PetscRandomCreate(PETSC_COMM_WORLD,&rctx);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(rctx);CHKERRQ(ierr);

  /* Single or multi processor ? */
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&comm_size);CHKERRQ(ierr);

  /* construct matrix */
  ierr = DMSetMatType(da,MATAIJ);CHKERRQ(ierr);
  ierr = DMCreateMatrix(da, &H);CHKERRQ(ierr);

  /* get local corners for this processor */
  ierr = DMDAGetCorners(da,&xs,&ys,0,&xm,&ym,0);CHKERRQ(ierr);

  /* Assemble the matrix */
  for (x=xs; x<xs+xm; x++) {
    for (y=ys; y<ys+ym; y++) {
      /* each lattice point sets only the *forward* pointing parameters (right, down),
         i.e. Nabla_1^+ and Nabla_2^+.
         In this way we can use only local random number creation. That means
         we also have to set the corresponding backward pointing entries. */
      /* Compute some normally distributed random numbers via Box-Muller */
      ierr = PetscRandomGetValueReal(rctx, &r1);CHKERRQ(ierr);
      r1   = 1.-r1; /* to change from [0,1) to (0,1], which we need for the log */
      ierr = PetscRandomGetValueReal(rctx, &r2);CHKERRQ(ierr);
      PetscReal R = PetscSqrtReal(-2.*PetscLogReal(r1));
      PetscReal c = PetscCosReal(2.*PETSC_PI*r2);
      PetscReal s = PetscSinReal(2.*PETSC_PI*r2);

      /* use those to set the field */
      uxy1 = PetscExpScalar(((PetscScalar) (R*c/beta))*PETSC_i);
      uxy2 = PetscExpScalar(((PetscScalar) (R*s/beta))*PETSC_i);

      sxy.i = x; sxy.j = y; /* the point where we are */

      /* center action */
      sxy.c = 0; /* spin 0, 0 */
      ierr  = MatSetValuesStencil(H, 1, &sxy, 1, &sxy, &rho, ADD_VALUES);CHKERRQ(ierr);
      sxy.c = 1; /* spin 1, 1 */
      val   = -rho;
      ierr  = MatSetValuesStencil(H, 1, &sxy, 1, &sxy, &val, ADD_VALUES);CHKERRQ(ierr);

      sxy_m.i = x+1; sxy_m.j = y; /* right action */
      sxy.c   = 0; sxy_m.c = 0; /* spin 0, 0 */
      val     = -uxy1; valconj = PetscConj(val);
      ierr    = MatSetValuesStencil(H, 1, &sxy_m, 1, &sxy, &val, ADD_VALUES);CHKERRQ(ierr);
      ierr    = MatSetValuesStencil(H, 1, &sxy, 1, &sxy_m, &valconj, ADD_VALUES);CHKERRQ(ierr);
      sxy.c   = 0; sxy_m.c = 1; /* spin 0, 1 */
      val     = -uxy1; valconj = PetscConj(val);
      ierr    = MatSetValuesStencil(H, 1, &sxy_m, 1, &sxy, &val, ADD_VALUES);CHKERRQ(ierr);
      ierr    = MatSetValuesStencil(H, 1, &sxy, 1, &sxy_m, &valconj, ADD_VALUES);CHKERRQ(ierr);
      sxy.c   = 1; sxy_m.c = 0; /* spin 1, 0 */
      val     = uxy1; valconj = PetscConj(val);
      ierr    = MatSetValuesStencil(H, 1, &sxy_m, 1, &sxy, &val, ADD_VALUES);CHKERRQ(ierr);
      ierr    = MatSetValuesStencil(H, 1, &sxy, 1, &sxy_m, &valconj, ADD_VALUES);CHKERRQ(ierr);
      sxy.c   = 1; sxy_m.c = 1; /* spin 1, 1 */
      val     = uxy1; valconj = PetscConj(val);
      ierr    = MatSetValuesStencil(H, 1, &sxy_m, 1, &sxy, &val, ADD_VALUES);CHKERRQ(ierr);
      ierr    = MatSetValuesStencil(H, 1, &sxy, 1, &sxy_m, &valconj, ADD_VALUES);CHKERRQ(ierr);

      sxy_m.i = x; sxy_m.j = y+1; /* down action */
      sxy.c   = 0; sxy_m.c = 0; /* spin 0, 0 */
      val     = -uxy2; valconj = PetscConj(val);
//.........这里部分代码省略.........
开发者ID:hsahasra，项目名称:petsc-magma-dense-mat，代码行数:101，代码来源:ex39.c

示例10: main

PetscInt main(PetscInt argc,char **args)
{
  PetscErrorCode  ierr;
  PetscMPIInt     rank,size;
  PetscInt        N0=10,N1=10,N2=10,N3=10,N4=10,N=N0*N1*N2*N3*N4;
  PetscRandom     rdm;
  PetscReal       enorm;
  Vec             x,y,z,input,output;
  Mat             A;
  PetscInt        DIM, dim[5],vsize;
  PetscReal       fac;

  ierr = PetscInitialize(&argc,&args,(char *)0,help);CHKERRQ(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_WORLD,PETSC_ERR_SUP, "This is a uni-processor example only");
  ierr = PetscRandomCreate(PETSC_COMM_SELF, &rdm);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(rdm);CHKERRQ(ierr);
  ierr = VecCreate(PETSC_COMM_SELF,&input);CHKERRQ(ierr);
  ierr = VecSetSizes(input,N,N);CHKERRQ(ierr);
  ierr = VecSetFromOptions(input);CHKERRQ(ierr);
  ierr = VecSetRandom(input,rdm);CHKERRQ(ierr);
  ierr = VecDuplicate(input,&output);

  DIM = 5; dim[0] = N0; dim[1] = N1; dim[2] = N2; dim[3] = N3; dim[4] = N4;
  ierr = MatCreateFFT(PETSC_COMM_SELF,DIM,dim,MATFFTW,&A);CHKERRQ(ierr);
  ierr = MatGetVecs(A,&x,&y);CHKERRQ(ierr);
  ierr = MatGetVecs(A,&z,PETSC_NULL);CHKERRQ(ierr);

  ierr = VecGetSize(x,&vsize);CHKERRQ(ierr);
  printf("The vector size  of input from the main routine is %d\n",vsize);

  ierr = VecGetSize(z,&vsize);CHKERRQ(ierr);
  printf("The vector size of output from the main routine is %d\n",vsize);

  ierr = InputTransformFFT(A,input,x);CHKERRQ(ierr);

  ierr = MatMult(A,x,y);CHKERRQ(ierr);
  ierr = MatMultTranspose(A,y,z);CHKERRQ(ierr);

  ierr = OutputTransformFFT(A,z,output);CHKERRQ(ierr);
  fac = 1.0/(PetscReal)N;
  ierr = VecScale(output,fac);CHKERRQ(ierr);
/*
  ierr = VecAssemblyBegin(input);CHKERRQ(ierr);
  ierr = VecAssemblyEnd(input);CHKERRQ(ierr);
  ierr = VecAssemblyBegin(output);CHKERRQ(ierr);
  ierr = VecAssemblyEnd(output);CHKERRQ(ierr);

  ierr = VecView(input,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  ierr = VecView(output,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
*/
  ierr = VecAXPY(output,-1.0,input);CHKERRQ(ierr);
  ierr = VecNorm(output,NORM_1,&enorm);CHKERRQ(ierr);
//  if (enorm > 1.e-14){
      ierr = PetscPrintf(PETSC_COMM_SELF,"  Error norm of |x - z| %e\n",enorm);CHKERRQ(ierr);
//      }

  ierr = VecDestroy(&output);CHKERRQ(ierr);
  ierr = VecDestroy(&input);CHKERRQ(ierr);
  ierr = VecDestroy(&x);CHKERRQ(ierr);
  ierr = VecDestroy(&y);CHKERRQ(ierr);
  ierr = VecDestroy(&z);CHKERRQ(ierr);
  ierr = MatDestroy(&A);CHKERRQ(ierr);
  ierr = PetscRandomDestroy(&rdm);CHKERRQ(ierr);
  PetscFinalize();
  return 0;

}
开发者ID:erdc-cm，项目名称:petsc-dev，代码行数:71，代码来源:ex153.c

示例11: main

int main(int argc,char **argv)
{
  Mat            A,B,C,D;
  PetscInt       i,M=10,N=5,j,nrows,ncols,am,an,rstart,rend;
  PetscErrorCode ierr;
  PetscRandom    r;
  PetscBool      equal,iselemental;
  PetscReal      fill = 1.0;
  IS             isrows,iscols;
  const PetscInt *rows,*cols;
  PetscScalar    *v,rval;
#if defined(PETSC_HAVE_ELEMENTAL)
  PetscBool      Test_MatMatMult=PETSC_TRUE;
#else
  PetscBool      Test_MatMatMult=PETSC_FALSE;
#endif
  PetscMPIInt    size;

  ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);

  ierr = PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL);CHKERRQ(ierr);
  ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
  ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);CHKERRQ(ierr);
  ierr = MatSetType(A,MATDENSE);CHKERRQ(ierr);
  ierr = MatSetFromOptions(A);CHKERRQ(ierr);
  ierr = MatSetUp(A);CHKERRQ(ierr);
  ierr = PetscRandomCreate(PETSC_COMM_WORLD,&r);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(r);CHKERRQ(ierr);

  /* Set local matrix entries */
  ierr = MatGetOwnershipIS(A,&isrows,&iscols);CHKERRQ(ierr);
  ierr = ISGetLocalSize(isrows,&nrows);CHKERRQ(ierr);
  ierr = ISGetIndices(isrows,&rows);CHKERRQ(ierr);
  ierr = ISGetLocalSize(iscols,&ncols);CHKERRQ(ierr);
  ierr = ISGetIndices(iscols,&cols);CHKERRQ(ierr);
  ierr = PetscMalloc1(nrows*ncols,&v);CHKERRQ(ierr);
  for (i=0; i<nrows; i++) {
    for (j=0; j<ncols; j++) {
      ierr         = PetscRandomGetValue(r,&rval);CHKERRQ(ierr);
      v[i*ncols+j] = rval;
    }
  }
  ierr = MatSetValues(A,nrows,rows,ncols,cols,v,INSERT_VALUES);CHKERRQ(ierr);
  ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = ISRestoreIndices(isrows,&rows);CHKERRQ(ierr);
  ierr = ISRestoreIndices(iscols,&cols);CHKERRQ(ierr);
  ierr = ISDestroy(&isrows);CHKERRQ(ierr);
  ierr = ISDestroy(&iscols);CHKERRQ(ierr);
  ierr = PetscRandomDestroy(&r);CHKERRQ(ierr);

  /* Test MatTranspose() */
  ierr = MatCreateTranspose(A,&C);CHKERRQ(ierr);
  ierr = MatTranspose(A,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr); /* B = A^T */
  ierr = MatMultEqual(C,B,10,&equal);CHKERRQ(ierr);
  if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"A^T*x != (x^T*A)^T");
  ierr = MatTranspose(A,MAT_REUSE_MATRIX,&B);CHKERRQ(ierr); /* B = A^T */
  ierr = MatMultEqual(C,B,10,&equal);CHKERRQ(ierr);
  if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"A^T*x != (x^T*A)^T");
  ierr = MatDestroy(&B);CHKERRQ(ierr);
  ierr = MatDuplicate(A,MAT_COPY_VALUES,&B);CHKERRQ(ierr);
  ierr = MatTranspose(B,MAT_INPLACE_MATRIX,&B);CHKERRQ(ierr);
  ierr = MatMultEqual(C,B,10,&equal);CHKERRQ(ierr);
  if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"A^T*x != (x^T*A)^T");
  ierr = MatDestroy(&B);CHKERRQ(ierr);
  ierr = MatDestroy(&C);CHKERRQ(ierr);

  /* Test MatMatMult() */
  if (Test_MatMatMult) {
#if !defined(PETSC_HAVE_ELEMENTAL)
    if (size > 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"This test requires ELEMENTAL");
#endif
    ierr = MatTranspose(A,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr); /* B = A^T */
    ierr = MatMatMult(B,A,MAT_INITIAL_MATRIX,fill,&C);CHKERRQ(ierr); /* C = B*A = A^T*A */
    ierr = MatMatMult(B,A,MAT_REUSE_MATRIX,fill,&C);CHKERRQ(ierr);

    /* Test MatDuplicate for matrix product */
    ierr = MatDuplicate(C,MAT_COPY_VALUES,&D);CHKERRQ(ierr);
    ierr = MatDestroy(&D);CHKERRQ(ierr);

    /* Test B*A*x = C*x for n random vector x */
    ierr = MatMatMultEqual(B,A,C,10,&equal);CHKERRQ(ierr);
    if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"B*A*x != C*x");
    ierr = MatDestroy(&C);CHKERRQ(ierr);

    ierr = MatMatMultSymbolic(B,A,fill,&C);CHKERRQ(ierr);
    for (i=0; i<2; i++) {
      /* Repeat the numeric product to test reuse of the previous symbolic product */
      ierr = MatMatMultNumeric(B,A,C);CHKERRQ(ierr);

      ierr = MatMatMultEqual(B,A,C,10,&equal);CHKERRQ(ierr);
      if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"B*A*x != C*x");
    }
    ierr = MatDestroy(&C);CHKERRQ(ierr);
    ierr = MatDestroy(&B);CHKERRQ(ierr);
  }

  /* Test MatTransposeMatMult() */
//.........这里部分代码省略.........
开发者ID:firedrakeproject，项目名称:petsc，代码行数:101，代码来源:ex104.c

示例12: main

int main(int argc,char **argv)
{
  Mat            A,B,C,D;
  PetscInt       i,M=10,N=5,j,nrows,ncols;
  PetscErrorCode ierr;
  PetscRandom    r;
  PetscBool      equal,iselemental;
  PetscReal      fill = 1.0;
  IS             isrows,iscols;
  const PetscInt *rows,*cols;
  PetscScalar    *v,rval;

  PetscInitialize(&argc,&argv,(char*)0,help);
  ierr = PetscOptionsGetInt(NULL,"-M",&M,NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(NULL,"-N",&N,NULL);CHKERRQ(ierr);
  ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
  ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,M,N);CHKERRQ(ierr);
  ierr = MatSetType(A,MATDENSE);CHKERRQ(ierr);
  ierr = MatSetFromOptions(A);CHKERRQ(ierr);
  ierr = MatSetUp(A);CHKERRQ(ierr);
  ierr = PetscRandomCreate(PETSC_COMM_WORLD,&r);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(r);CHKERRQ(ierr);

  /* Set local matrix entries */
  ierr = MatGetOwnershipIS(A,&isrows,&iscols);CHKERRQ(ierr);
  ierr = ISGetLocalSize(isrows,&nrows);CHKERRQ(ierr);
  ierr = ISGetIndices(isrows,&rows);CHKERRQ(ierr);
  ierr = ISGetLocalSize(iscols,&ncols);CHKERRQ(ierr);
  ierr = ISGetIndices(iscols,&cols);CHKERRQ(ierr);
  ierr = PetscMalloc(nrows*ncols*sizeof(*v),&v);CHKERRQ(ierr);
  for (i=0; i<nrows; i++) {
    for (j=0; j<ncols; j++) {
      ierr         = PetscRandomGetValue(r,&rval);CHKERRQ(ierr);
      v[i*ncols+j] = rval;
    }
  }
  ierr = MatSetValues(A,nrows,rows,ncols,cols,v,INSERT_VALUES);CHKERRQ(ierr);
  ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = ISRestoreIndices(isrows,&rows);CHKERRQ(ierr);
  ierr = ISRestoreIndices(iscols,&cols);CHKERRQ(ierr);
  ierr = ISDestroy(&isrows);CHKERRQ(ierr);
  ierr = ISDestroy(&iscols);CHKERRQ(ierr);
  ierr = PetscFree(v);CHKERRQ(ierr);
  ierr = PetscRandomDestroy(&r);CHKERRQ(ierr);

  /* Test MatMatMult() */
  ierr = MatTranspose(A,MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr); /* B = A^T */
  ierr = MatMatMult(B,A,MAT_INITIAL_MATRIX,fill,&C);CHKERRQ(ierr); /* C = B*A = A^T*A */

  ierr = MatMatMultSymbolic(B,A,fill,&D);CHKERRQ(ierr); /* D = B*A = A^T*A */
  for (i=0; i<2; i++) {
    /* Repeat the numeric product to test reuse of the previous symbolic product */
    ierr = MatMatMultNumeric(B,A,D);CHKERRQ(ierr);
  }
  ierr = MatMultEqual(C,D,10,&equal);CHKERRQ(ierr);
  if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"C != D");
  ierr = MatDestroy(&D);CHKERRQ(ierr);

  /* Test MatTransposeMatMult() */
  ierr = PetscObjectTypeCompare((PetscObject)A,MATELEMENTAL,&iselemental);CHKERRQ(ierr);
  if (!iselemental) {
    ierr = MatTransposeMatMult(A,A,MAT_INITIAL_MATRIX,fill,&D);CHKERRQ(ierr); /* D = A^T*A */
    ierr = MatEqual(C,D,&equal);CHKERRQ(ierr);
    if (!equal) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"C != D");
    ierr = MatDestroy(&D);CHKERRQ(ierr);
  }

  ierr = MatDestroy(&C);CHKERRQ(ierr);
  ierr = MatDestroy(&B);
  ierr = MatDestroy(&A);
  PetscFinalize();
  return(0);
}
开发者ID:hansec，项目名称:petsc，代码行数:74，代码来源:ex104.c

示例13: main

int main(int argc,char **argv)
{
  PetscMPIInt      rank,size;
  PetscInt         N            = 6,M=8,P=5,dof=1;
  PetscInt         stencil_width=1,pt=0,st=0;
  PetscErrorCode   ierr;
  PetscBool        flg2,flg3,isbinary,mpiio;
  DMBoundaryType   bx           = DM_BOUNDARY_NONE,by = DM_BOUNDARY_NONE,bz = DM_BOUNDARY_NONE;
  DMDAStencilType  stencil_type = DMDA_STENCIL_STAR;
  DM               da,da2;
  Vec              global1,global2;
  PetscScalar      mone = -1.0;
  PetscReal        norm;
  PetscViewer      viewer;
  PetscRandom      rdm;
#if defined(PETSC_HAVE_HDF5)
  PetscBool ishdf5;
#endif

  ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);

  ierr = PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(NULL,NULL,"-P",&P,NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(NULL,NULL,"-dof",&dof,NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(NULL,NULL,"-stencil_width",&stencil_width,NULL);CHKERRQ(ierr);
  ierr = PetscOptionsGetInt(NULL,NULL,"-periodic",&pt,NULL);CHKERRQ(ierr);
  if (pt == 1) bx = DM_BOUNDARY_PERIODIC;
  if (pt == 2) by = DM_BOUNDARY_PERIODIC;
  if (pt == 4) {bx = DM_BOUNDARY_PERIODIC; by = DM_BOUNDARY_PERIODIC;}

  ierr         = PetscOptionsGetInt(NULL,NULL,"-stencil_type",&st,NULL);CHKERRQ(ierr);
  stencil_type = (DMDAStencilType) st;

  ierr = PetscOptionsHasName(NULL,NULL,"-oned",&flg2);CHKERRQ(ierr);
  ierr = PetscOptionsHasName(NULL,NULL,"-twod",&flg2);CHKERRQ(ierr);
  ierr = PetscOptionsHasName(NULL,NULL,"-threed",&flg3);CHKERRQ(ierr);

  ierr = PetscOptionsHasName(NULL,NULL,"-binary",&isbinary);CHKERRQ(ierr);
#if defined(PETSC_HAVE_HDF5)
  ierr = PetscOptionsHasName(NULL,NULL,"-hdf5",&ishdf5);CHKERRQ(ierr);
#endif
  ierr = PetscOptionsHasName(NULL,NULL,"-mpiio",&mpiio);CHKERRQ(ierr);
  if (flg2) {
    ierr = DMDACreate2d(PETSC_COMM_WORLD,bx,by,stencil_type,M,N,PETSC_DECIDE,PETSC_DECIDE,dof,stencil_width,0,0,&da);CHKERRQ(ierr);
  } else if (flg3) {
    ierr = DMDACreate3d(PETSC_COMM_WORLD,bx,by,bz,stencil_type,M,N,P,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,dof,stencil_width,0,0,0,&da);CHKERRQ(ierr);
  } else {
    ierr = DMDACreate1d(PETSC_COMM_WORLD,bx,M,dof,stencil_width,NULL,&da);CHKERRQ(ierr);
  }

  ierr = DMCreateGlobalVector(da,&global1);CHKERRQ(ierr);
  ierr = PetscObjectSetName((PetscObject)global1,"Test_Vec");CHKERRQ(ierr);
  ierr = PetscRandomCreate(PETSC_COMM_WORLD,&rdm);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(rdm);CHKERRQ(ierr);
  ierr = VecSetRandom(global1,rdm);CHKERRQ(ierr);
  if (isbinary) {
    if (mpiio) {
      ierr = PetscOptionsSetValue(NULL,"-viewer_binary_mpiio","");CHKERRQ(ierr);
    }
    ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_WRITE,&viewer);CHKERRQ(ierr);
#if defined(PETSC_HAVE_HDF5)
  } else if (ishdf5) {
    ierr = PetscViewerHDF5Open(PETSC_COMM_WORLD,"temp",FILE_MODE_WRITE,&viewer);CHKERRQ(ierr);
#endif
  } else SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Invalid Viewer : Run with -binary or -hdf5 option\n");
  ierr = VecView(global1,viewer);CHKERRQ(ierr);
  ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);

  ierr = PetscPrintf(PETSC_COMM_WORLD,"Global vector written to temp file is \n");CHKERRQ(ierr);
  ierr = VecView(global1,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);

  if (flg2) {
    ierr = DMDACreate2d(PETSC_COMM_WORLD,bx,by,stencil_type,M,N,PETSC_DECIDE,PETSC_DECIDE,dof,stencil_width,0,0,&da2);CHKERRQ(ierr);
  } else if (flg3) {
    ierr = DMDACreate3d(PETSC_COMM_WORLD,bx,by,bz,stencil_type,M,N,P,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,dof,stencil_width,0,0,0,&da2);CHKERRQ(ierr);
  } else {
    ierr = DMDACreate1d(PETSC_COMM_WORLD,bx,M,dof,stencil_width,NULL,&da2);CHKERRQ(ierr);
  }

  if (isbinary) {
    ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_READ,&viewer);CHKERRQ(ierr);
#if defined(PETSC_HAVE_HDF5)
  } else if (ishdf5) {
    ierr = PetscViewerHDF5Open(PETSC_COMM_WORLD,"temp",FILE_MODE_READ,&viewer);CHKERRQ(ierr);
#endif
  } else SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Invalid Viewer : Run with -binary or -hdf5 option\n");

  ierr = DMCreateGlobalVector(da2,&global2);CHKERRQ(ierr);
  ierr = PetscObjectSetName((PetscObject)global2,"Test_Vec");CHKERRQ(ierr);
  ierr = VecLoad(global2,viewer);CHKERRQ(ierr);
  ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr);

  ierr = PetscPrintf(PETSC_COMM_WORLD,"Global vector read from temp file is \n");CHKERRQ(ierr);
  ierr = VecView(global2,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  ierr = VecAXPY(global2,mone,global1);CHKERRQ(ierr);
  ierr = VecNorm(global2,NORM_MAX,&norm);CHKERRQ(ierr);
  if (norm != 0.0) {
//.........这里部分代码省略.........
开发者ID:tom-klotz，项目名称:petsc，代码行数:101，代码来源:ex33.c

示例14: main

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
//.........这里部分代码省略.........
开发者ID:00liujj，项目名称:petsc，代码行数:101，代码来源:ex11.c

示例15: main

int main(int argc,char **argv)
{
  TS             ts;            /* ODE integrator */
  Vec            U;             /* solution will be stored here */
  Mat            A;             /* Jacobian matrix */
  PetscErrorCode ierr;
  PetscMPIInt    size;
  PetscInt       n = 2;
  AppCtx         ctx;
  PetscScalar    *u;

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     Initialize program
     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
  if (size > 1) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Only for sequential runs");

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    Create necessary matrix and vectors
    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
  ierr = MatSetSizes(A,n,n,PETSC_DETERMINE,PETSC_DETERMINE);CHKERRQ(ierr);
  ierr = MatSetFromOptions(A);CHKERRQ(ierr);
  ierr = MatSetUp(A);CHKERRQ(ierr);

  ierr = MatCreateVecs(A,&U,NULL);CHKERRQ(ierr);

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    Set runtime options
    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  ierr = PetscOptionsBegin(PETSC_COMM_WORLD,NULL,"Reaction options","");CHKERRQ(ierr);
  {
    ctx.omega_s = 1.0;
    ierr        = PetscOptionsScalar("-omega_s","","",ctx.omega_s,&ctx.omega_s,NULL);CHKERRQ(ierr);
    ctx.H       = 1.0;
    ierr        = PetscOptionsScalar("-H","","",ctx.H,&ctx.H,NULL);CHKERRQ(ierr);
    ctx.E       = 1.0;
    ierr        = PetscOptionsScalar("-E","","",ctx.E,&ctx.E,NULL);CHKERRQ(ierr);
    ctx.V       = 1.0;
    ierr        = PetscOptionsScalar("-V","","",ctx.V,&ctx.V,NULL);CHKERRQ(ierr);
    ctx.X       = 1.0;
    ierr        = PetscOptionsScalar("-X","","",ctx.X,&ctx.X,NULL);CHKERRQ(ierr);

    ierr = VecGetArray(U,&u);CHKERRQ(ierr);
    u[0] = 1;
    u[1] = .7;
    ierr = VecRestoreArray(U,&u);CHKERRQ(ierr);
    ierr = PetscOptionsGetVec(NULL,NULL,"-initial",U,NULL);CHKERRQ(ierr);
  }
  ierr = PetscOptionsEnd();CHKERRQ(ierr);

  ierr = PetscRandomCreate(PETSC_COMM_WORLD,&ctx.rand);CHKERRQ(ierr);
  ierr = PetscRandomSetFromOptions(ctx.rand);CHKERRQ(ierr);

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     Create timestepping solver context
     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  ierr = TSCreate(PETSC_COMM_WORLD,&ts);CHKERRQ(ierr);
  ierr = TSSetProblemType(ts,TS_NONLINEAR);CHKERRQ(ierr);
  ierr = TSSetType(ts,TSROSW);CHKERRQ(ierr);
  ierr = TSSetIFunction(ts,NULL,(TSIFunction) IFunction,&ctx);CHKERRQ(ierr);
  ierr = TSSetIJacobian(ts,A,A,(TSIJacobian)IJacobian,&ctx);CHKERRQ(ierr);

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     Set initial conditions
   - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  ierr = TSSetSolution(ts,U);CHKERRQ(ierr);

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     Set solver options
   - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  ierr = TSSetDuration(ts,100000,2000.0);CHKERRQ(ierr);
  ierr = TSSetExactFinalTime(ts,TS_EXACTFINALTIME_STEPOVER);CHKERRQ(ierr);
  ierr = TSSetInitialTimeStep(ts,0.0,.001);CHKERRQ(ierr);
  ierr = TSSetFromOptions(ts);CHKERRQ(ierr);

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     Solve nonlinear system
     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  ierr = TSSolve(ts,U);CHKERRQ(ierr);

  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
     Free work space.  All PETSc objects should be destroyed when they are no longer needed.
   - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
  ierr = MatDestroy(&A);CHKERRQ(ierr);
  ierr = VecDestroy(&U);CHKERRQ(ierr);
  ierr = TSDestroy(&ts);CHKERRQ(ierr);
  ierr = PetscRandomDestroy(&ctx.rand);CHKERRQ(ierr);
  ierr = PetscFinalize();
  return ierr;
}
开发者ID:tom-klotz，项目名称:petsc，代码行数:92，代码来源:ex1.c


注：本文中的PetscRandomSetFromOptions函数示例由纯净天空整理自Github/MSDocs等开源代码及文档管理平台，相关代码片段筛选自各路编程大神贡献的开源项目，源码版权归原作者所有，传播和使用请参考对应项目的License；未经允许，请勿转载。