Actual source code: ex17.c

  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2013, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.

  8:    SLEPc is free software: you can redistribute it and/or modify it under  the
  9:    terms of version 3 of the GNU Lesser General Public License as published by
 10:    the Free Software Foundation.

 12:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY
 13:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS
 14:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for
 15:    more details.

 17:    You  should have received a copy of the GNU Lesser General  Public  License
 18:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 19:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 20: */

 22: static char help[] = "Solves a quadratic eigenproblem (l^2*M + l*C + K)*x = 0 with matrices loaded from a file.\n\n"
 23:   "The command line options are:\n"
 24:   "  -M <filename>, where <filename> = matrix (M) file in PETSc binary form.\n"
 25:   "  -C <filename>, where <filename> = matrix (C) file in PETSc binary form.\n"
 26:   "  -K <filename>, where <filename> = matrix (K) file in PETSc binary form.\n\n";

 28: #include <slepcqep.h>

 32: int main(int argc,char **argv)
 33: {
 34:   Mat            M,C,K;           /* problem matrices */
 35:   QEP            qep;             /* quadratic eigenproblem solver context */
 36:   QEPType        type;
 37:   PetscReal      tol;
 38:   PetscInt       nev,maxit,its;
 39:   char           filename[PETSC_MAX_PATH_LEN];
 40:   PetscViewer    viewer;
 41:   PetscBool      flg;

 44:   SlepcInitialize(&argc,&argv,(char*)0,help);

 46:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 47:         Load the matrices that define the quadratic eigenproblem
 48:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

 50:   PetscPrintf(PETSC_COMM_WORLD,"\nQuadratic eigenproblem stored in file.\n\n");
 51: #if defined(PETSC_USE_COMPLEX)
 52:   PetscPrintf(PETSC_COMM_WORLD," Reading COMPLEX matrices from binary files...\n");
 53: #else
 54:   PetscPrintf(PETSC_COMM_WORLD," Reading REAL matrices from binary files...\n");
 55: #endif

 57:   PetscOptionsGetString(NULL,"-M",filename,PETSC_MAX_PATH_LEN,&flg);
 58:   if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate a file name for matrix M with the -M option");
 59:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
 60:   MatCreate(PETSC_COMM_WORLD,&M);
 61:   MatSetFromOptions(M);
 62:   MatLoad(M,viewer);
 63:   PetscViewerDestroy(&viewer);

 65:   PetscOptionsGetString(NULL,"-C",filename,PETSC_MAX_PATH_LEN,&flg);
 66:   if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate a file name for matrix C with the -C option");
 67:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
 68:   MatCreate(PETSC_COMM_WORLD,&C);
 69:   MatSetFromOptions(C);
 70:   MatLoad(C,viewer);
 71:   PetscViewerDestroy(&viewer);

 73:   PetscOptionsGetString(NULL,"-K",filename,PETSC_MAX_PATH_LEN,&flg);
 74:   if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate a file name for matrix K with the -K option");
 75:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);
 76:   MatCreate(PETSC_COMM_WORLD,&K);
 77:   MatSetFromOptions(K);
 78:   MatLoad(K,viewer);
 79:   PetscViewerDestroy(&viewer);

 81:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 82:                 Create the eigensolver and set various options
 83:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

 85:   /*
 86:      Create eigensolver context
 87:   */
 88:   QEPCreate(PETSC_COMM_WORLD,&qep);

 90:   /*
 91:      Set matrices
 92:   */
 93:   QEPSetOperators(qep,M,C,K);

 95:   /*
 96:      Set solver parameters at runtime
 97:   */
 98:   QEPSetFromOptions(qep);

100:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
101:                       Solve the eigensystem
102:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

104:   QEPSolve(qep);
105:   QEPGetIterationNumber(qep,&its);
106:   PetscPrintf(PETSC_COMM_WORLD," Number of iterations of the method: %D\n",its);

108:   /*
109:      Optional: Get some information from the solver and display it
110:   */
111:   QEPGetType(qep,&type);
112:   PetscPrintf(PETSC_COMM_WORLD," Solution method: %s\n\n",type);
113:   QEPGetDimensions(qep,&nev,NULL,NULL);
114:   PetscPrintf(PETSC_COMM_WORLD," Number of requested eigenvalues: %D\n",nev);
115:   QEPGetTolerances(qep,&tol,&maxit);
116:   PetscPrintf(PETSC_COMM_WORLD," Stopping condition: tol=%.4G, maxit=%D\n",tol,maxit);

118:   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
119:                     Display solution and clean up
120:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

122:   QEPPrintSolution(qep,NULL);
123:   QEPDestroy(&qep);
124:   MatDestroy(&M);
125:   MatDestroy(&C);
126:   MatDestroy(&K);
127:   SlepcFinalize();
128:   return 0;
129: }