710 likes | 1.27k Views
Use PETSC and SLEPC package to solve large scale linear system ---quick start. Sheng Yi Wang , Department of Mathematics, National Taiwan University 2011/07/29. Outline. PETSC ---Vector ---Matrix ---linear system solver SLEPC --- eigenvalue problem solver PBS .
E N D
Use PETSC and SLEPC package to solve large scale linear system---quick start Sheng Yi Wang, Department of Mathematics, National Taiwan University2011/07/29
Outline • PETSC ---Vector ---Matrix ---linear system solver • SLEPC ---eigenvalue problem solver • PBS PETSC & SLEPC
What is PETSC ?? • PETSC = Portable Extensible Toolkit for Scientific Computation • PETSC is intended for use in large-scale application projects. • PETSC support MPI, that can parallel execute. • PETSC can interfaces many external software ex: Matlab, Mathematica, MUMPS…etc PETSC & SLEPC
What can PETSC do?? • Vector operation • Matrix operation • Linear system ( sparse or dense) Ax=b • Nonlinear solver • ODE & PDE solve ( steady state or time dependent) PETSC & SLEPC
User’s background • Some knowledge of parallel (MPI command) • You are familiar with C/Fortran language. • You are familiar with Linux environment. PS: PETSC can be installed in Windows (but you have to install many other packages) PETSC & SLEPC
How install?? • Step1 : Set environment variables PETSC_DIR ex: PETSC_DIR=/home/sywang/petsc-3.1-p8; export PETSC_DIR • Step 2. Configure (to check your environment) ---BLAS, C compiler, (MPI) must be installed ex1: ./configure PETSC_ARCH=intel-64-complex --with-cc=icc --with-fc=ifort --with-debugging=0 --download-c-blas-lapack=1 --download-mpich=1 --with-scalar-type=complex PETSC & SLEPC
Ex2: ./configure PETSC_ARCH=intel-64-O3-real --with-cc=gcc --with-fc=gfortran --with-debugging=1 --with-blas-lapack-dir=/opt/intel/mkl/10.2.5.035 --with-mpi-dir=/opt/mpich2 --with-scalar-type=real--download-mumps=1 PETSC & SLEPC
Step 2: make all • Step 3: make test • Set PETSC_DIR and PETSC_ARCH to ~/.bashrc • If you install PETSC successfully, install SLEPC is very easy PETSC & SLEPC
Don’t be afraid, install the PETSC is the most difficult process, after that , all you need to do is the following: 1. Call function 2. Set parameter 3. Show the output and explain the results PETSC & SLEPC
Start using PETSC PETSC & SLEPC
Declare Variables Vec sol, rhs; MatMtx_A; KSPksp; PetscInt ii= 10; PetscScalarvalue[3]; PetscScalarval_rhs; PetscErrorCodeierr; PETSC & SLEPC
Rough coding PetscInitialize(); ObjCreate(MPI_comm,&obj); ObjSetType(obj, ); ObjSetFromOptions(obj, ); ObjSolve(obj, ); ObjGetxxx(obj, ); ObjDestroy(obj); PetscFinalize() PETSC & SLEPC
Vector PETSC & SLEPC
I. Vector • Vec a • VecCreate(MPI_Commcomm,Vec* x) • VecCreate(PETSC_COMM_WORLD, &a) • VecSetSizes(Vec v, PetscInt n, PetscInt N) • VecSetSizes(a,PETSC_DECIDE,20); PETSC & SLEPC
I. Vector • VecSet(Vecx,PetscScalar alpha) • VecSet(a,1.0) • VecSetValues(Vecx,PetscIntni,constPetscInt ix[],const PetscScalar y[],InsertModeiora) • VecSetValues(a,1,0,-3, INSERT_VALUES) InsertMode : INSERT_VALUES, ADD_VALUES PETSC & SLEPC
I. Vector • VecSetRandom(Vecx,PetscRandomrctx) • VecSetRandom(b,r) • VecSetFromOptions(Vecvec) • VecSetFromOptions(a) • VecDuplicate(Vecv,Vec *newv) • VecDuplicate(a,&b) PETSC & SLEPC
I. Vector • VecView(a,PETSC_VIEWER_STDOUT_WORLD); • VecAssemblyBegin(a); • VecAssemblyEnd(a); • VecDestroy(a); PETSC & SLEPC
Example presentation • Example 1: Create Vec and use some basic vector operation PETSC & SLEPC
Matrix PETSC & SLEPC
2. Matrix • Mat mtx_a • MatCreate(MPI_Commcomm,Mat* A) • MatCreate(PETSC_SOMM_WORLD,&mtx_a) • MatSetSizes(Mat A,intm,intn,intM,int N) • MatSetSizes(mtx_a,PETSC_DECIDE, PETSC_DECIDE,10,10) • MatSetFromOptions(Mat A) • MatSetFromOptions(mtx_a) PETSC & SLEPC
2. Matrix • MatSetValues(Mat A,PetscIntm,const PetscIntidxm[],PetscIntn,constPetscInt idxn[],const PetscScalar v[],InsertMode addv) • MatSetValues(mtx_a,1,1,1,1,2.0,INSERT_VALUE) PETSC & SLEPC
Example • I have a small matrix s_a [1 2 ;3 4] but I want to insert s_a into global matrix mtx_a and the index is (1,2) MatSetValues(mtx_a,2,1,2,2,s_a,INSERT_VALUE) 0 0 0 0 0 0 0 1 2 0 0 0 3 4 0 0 0 0 0 0 PETSC & SLEPC
2. Matrix • MatAssemblyBegin(Mat, MAT_FINAL_ASSEMBLY) • MatAssemblyEnd(Mat ,MAT_FINAL_ASSEMBLY) • MatDestroy(Mat ) PETSC & SLEPC
**2. Matrix: Matrix-Free • If you don’t want to create all matrix ( the matrix is too large),PETSC allow users write their own matrix operator • extern intmult(Mat,Vec,Vec); MatCreateShell(comm,m,n,M,N,ctx,&mat); MatShellSetOperation(mat,MATOP_MULT,(void(*)(void))mult); MatDestroy(mat); PETSC & SLEPC
2. Read matrix from file • In many case, someone has created the matrix, so we don’t rewrite the matrix, we just want to read the matrix and solve them quickly, how to do that?? • If you have matrix which is ASCII format(just store the nonzero element , i , j, aij) we can call PETSC function to transfer ASCII to binary file ( PETSC can read them very quickly) PETSC & SLEPC
Binary file • PetscViewerBinaryOpen(MPI_Commcomm,const char name[],PetscViewerFileTypetype,PetscViewer *binv) • PetscViewerBinaryOpen(PETSC_COMM_WORLD,fileout,FILE_MODE_WRITE,&view) PETSC & SLEPC
Example presentation • Example 2: Read matrix in ASCII format and transfer them to PETSC binary file • Example 3: Read matrix from PETSC binary file and create a vector and use a basic matrix operation PETSC & SLEPC
Linear System (KSP) PETSC & SLEPC
3.Linear system solver : KSP • Linear system problem: give matrix A and vector b solve Ax=b • The dimension is too large to find inverse, and the matrix is sparse, so we need to use iterative method to solve them. • The basic method to solve linear system is Krylov subspace methods. PETSC & SLEPC
Linear Solvers in PETSC PETSC & SLEPC
3. Linear system solver : KSP First of all, set the matrix A and rhs vector b Declare Variables KSP ksp • KSPCreate(MPI_Commcomm,KSP *ksp) • KSPCreate(PETSC_COMM_WORLD,&ksp) • KSPSetOperators(KSP ksp,MatAmat,MatPmat,MatStructure flag) • KSPSetOperators(ksp,A,A, SAME_NONZERO_PATTERN) PETSC & SLEPC
3. Linear system solver : KSP • KSPSetType(KSP ksp, const KSPType type) • KSPSetType(ksp,KSPGMRES) -ksp_type • KSPSetFromOptions(KSP ksp) • KSPSetFromOptions(ksp) • KSPSolve(ksp,Vecb,Vec x) • KSPSolve(ksp,b,x) • KSPGetIterationNumber(KSP ksp,PetscInt *its) • KSPGetIterationNumber(ksp,&it) PETSC & SLEPC
Preconditioner • In many case, the matrix A has large condition number, so we need use a proconditoner to reduce condition number • PETSC provide many preconditioner, and some of them can use MPI to parallel. PETSC & SLEPC
Precondioner in PETSC PETSC & SLEPC
3. Linear system solver : KSP • KSPSetFromOptions(KSP ksp) • KSPSetFromOptions(ksp) • KSPGetPC(KSP ksp,PC *pc) • KSPGetPC(ksp,&pc) • PCSetType(PC pc, PCType type) • PCSetType(pc,PCBJACOBI) -pc_type PETSC & SLEPC
How to Parallel ?? • MPI- Message Passing Interface PETSC and SLEPC support MPI and user don’t have to call MPI function • mpiexec -np 4 yourprogram • mpirun -np 8 -machinefile mf yourprogram • PETSC_COMM_SELF PETSC_COMM_WORLD • PETSC and SLEPC initial will include mpi.h, so if you want to use MPI function,you can use them, too PETSC & SLEPC
Assign parameters at run time • Serial : • ./ex4 -ksp_typebcgs -pc_typelu -ksp_rtol 1e-4 • ./ex4 -ksp_typebcgs -pc_typelu -ksp_rtol 1e-4 • ./ex4 -ksp_typegmres -pc_typeasm -ksp_rtol 1e-8 -ksp_max_it 20 • Parallel: In one node : • mpiexec –np 4 ./ex4 -ksp_type cg -pc_typesor - pc_sor_omega 1.8 -ksp_rtol 1e-4 -ksp_view In multi nodes: • mpiexec –np 16 –machinefile mf ./ex4 -ksp_type cg -pc_typeasm -ksp_rtol 1e-4 -ksp_view -ksp_monitor_true_residual PETSC & SLEPC
Example presentation • Example 4 : Read matrix from a binary file, and call PETSC function to solve linear system. And show how to assign parameters at run time PETSC & SLEPC
makefile include ${PETSC_DIR}/conf/base linearsym : linearsym.ochkopts -${CLINKER} –o linearsymlinearsym.o{PETSC_LIB} ${RM} linearsym.o PETSC & SLEPC
Bash script #!/bin/bash for((i=10;i<=30;i=i+5)) do ./linearsym -n $i -ksp_typegmres -pc_typejacobi -ksp_max_it 100 -ksp_view > result_gmres_$i done for((i=10;i<=30;i=i+5)) do ./linearsym -n $i -ksp_type cg -pc_typejacobi -ksp_max_it 100 -ksp_view > result_cg_$i done PETSC & SLEPC
Perl (bash script like) • #!/usr/bin/perl for ($i=1; $i<=199; $i++) { $sor_omega=0.01*$i; system("~/program/openmpi2/bin/mpiexec -np $i -machinefile mf10 ex4 -fin m22103 –ksp_type cg –pc_typesor –pc_sor_omega$sor_omega-log_summary > m22103_cg_sor_omega$sor_omega"); system("echo $i"); sleep(3);} PETSC & SLEPC
Eigenvalue (EPS) PETSC & SLEPC
4.Eigenvalue Problem Solver: SLEPC • SLEPC the Scalable Library for Eigenvalue Problem Computations • Standard Eigenvalue problem : give matrix A, want to find unknown vector x and value k Ax=kx • General Eigenvalue problem : give matrix A and matrix B, want to find unknown vector x and value k Ax=kBx • Still , the matrix is very large and sparse. PETSC & SLEPC
How to install SLEPC • Export SLEPC_DIR =/home/sywang/petsc-3.1-p6; • ./configure (they will follow the PETSC environmental setting) • Make all • Make test PETSC & SLEPC
Declare Variables EPS eps; Mat A,B ; PetscIntii,nn = 10,col[3]; PetscScalar value[3]; PetscScalarkr,ki; PetscErrorCodeierr; PETSC & SLEPC
4.Eigenvalue Problem Solver: SLEPC • EPSCreate(MPI_Commcomm,EPS *eps) • EPSCreate(PETSC_COMM_WORLD,eps) • EPSSetOperators(EPS eps,MatA,Mat B) • EPSSetOperators(eps,A,B) Ax=kBx • EPSSetOperators(eps,A,PETSC_NULL) Ax=kx • EPSSetProblemType(EPS eps,EPSProblemType type) • EPSSetProblemType(eps,EPS_HEP) -eps_hermitian • EPSSetProblemType(eps,EPS_NHEP) -eps_non_hermitian PETSC & SLEPC