petsc toy

      petsc 结构化数据使用 是今年初的认识。今天重新看 solving coupled nonlinear partial differential equations using petsc, 完善下。

      1 Distributed Arrays(DA)

      -- when solving problems on logically regular rectangular grids. DA object will partition the grid among the processors and wiki keep track of which grid points ( including ghost points) each processor will need access to. 

      -- For a processor to obtain the information about the global grid and its local grid calling  DAGetLocalInfo(DA da, DALocalInfo *info) 

  

     2 Vectors (Vec)

      -- Create global/local vectors based on DA

      DACreateGlobalVector(DA da, Vec *global_vec)

      DACreateLocalVector(DA da, Vec *local_vec)

     -- Assigning/accessing Vectors

     To assign global DA vector, usually (for 2d/3d) turn a DA vector into a multicomponent structure of multidimensional array (e.g. Field[][]) , then elements of vector can be assigned or accessed through this array. When it is done, the vector should restore calling  DAVecRestoreArray(DA da, Vec vex, void* array) 

     To assign local DA vector(containing ghost point information), calling  DAGlobalToLocalBegin(DA da, Vec global, InsertMode mode, Vec local), followed by DAGlobalToLocalEnd(DA da, Vec global, InsertMode mode, Vec local)

     3 Matrices(Mat)

     -- Create Mat based on DA

     DAGetMatrix(DA da, MatType mtype, Mat *J)

    -- Assembling matrices

    PETSC handle global indexing and the arrangement of the matrix, by providing PETSC with the grid indices and the component indices of the rows and columns to be assigned, it can determined where in the matrix to place the corresponding values. 

   Multiple rows and columns can be inserted at once calling

  MatSetValuesStencil(Mat mat, PetscInt m, const MatStencil row[], PetscInt n, const MatStencil col[], const PetscScalar val[], InsertMode mode)

  ## row, col are MatStencil structures that contain grid indices and component indices of the rows and columns being inserted.

  A common way to build a matrix is row by row. the nonzero elements in a row are inserted at once, the "col" MatStructure will be an array with as many entries as there are non zeros in the row being inserted.

 Once the matrix values have been assigned, the matrix needs to be constructed. calling 

  MatAssemblyBegin(Mat mat, MatAssemblyType type)

 MatAssemblyEnd(Mat mat, MatAssemblyType type)

 4  Preconditioning and solving

  -- Scalable Nonlinear Equation Solvers(SNES)

  Create the SNES object  SNESCreate(MPI_Comm comm, SNES *snes)

   -- associating SNES with nonlinear equations

  SNESSetFunction(SNES snes, Vec soln, PetscErrorCode(*func) (SNES, Vec, Vec, void*), void* ctx)

    ## ctx, a user defined data structure

  ##  func a C function that evaluates the nonlinear equations

        # SNESSetFromOptions(SNES snes)

  --solving

  SNESSolve(SNES snes, Vec b, Vec soln)

   --KSP

   Solving   nonlinear equations via Newton's method requires the solution of linear systems, and the linear system  are often preconditioned. (a hierarchy)

   thus, each SNES  object is linked to a KSP object, and each KSP object is linked to a PC object.

    -- linking a ksp

    SNESGetKSP(SNES snes, KSP *ksp)

    KSPSetType(KSP ksp, KSPType type)  # set the iterative method used (e.g. cg, gmres, bcd)

         --additional options can be set on certain ksp type e.g. KSPGMRESSetRestart(KSP kip, PetscInt restart) 

       KSPSetFromOptions(KSP ksp)

        --linking a PC

        KSPGetPC(KSP ksp, PC* pc)

        PCSetType(PC pc, PCType type)  #set pc types(e.g.  Jacobi, block Jacobi, additive Schwarz)

        PCSetFromOptions(PC pc)