Actual source code: test42.c
slepc-3.20.2 2024-03-15
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: static char help[] = "Block-diagonal orthogonal eigenproblem.\n"
12: "The command line options are:\n"
13: " -n <n>, where <n> = matrix dimension.\n"
14: " -seed <s>, where <s> = seed for random number generation.\n\n";
16: #include <slepceps.h>
18: int main(int argc,char **argv)
19: {
20: EPS eps;
21: Mat A;
22: PetscRandom rand;
23: PetscScalar val,c,s;
24: PetscInt n=30,i,seed=0x12345678;
25: PetscMPIInt rank;
27: PetscFunctionBeginUser;
28: PetscCall(SlepcInitialize(&argc,&argv,(char*)0,help));
29: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank));
30: PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL));
31: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\nOrthogonal eigenproblem, n=%" PetscInt_FMT "\n\n",n));
33: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
34: Generate the matrix
35: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
37: PetscCall(PetscRandomCreate(PETSC_COMM_WORLD,&rand));
38: PetscCall(PetscRandomSetFromOptions(rand));
39: PetscCall(PetscOptionsGetInt(NULL,NULL,"-seed",&seed,NULL));
40: PetscCall(PetscRandomSetSeed(rand,seed));
41: PetscCall(PetscRandomSeed(rand));
42: PetscCall(PetscRandomSetInterval(rand,0,2*PETSC_PI));
44: PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
45: PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n));
46: PetscCall(MatSetFromOptions(A));
47: PetscCall(MatSetUp(A));
49: if (!rank) {
50: for (i=0;i<n/2;i++) {
51: PetscCall(PetscRandomGetValue(rand,&val));
52: c = PetscCosReal(PetscRealPart(val));
53: s = PetscSinReal(PetscRealPart(val));
54: PetscCall(MatSetValue(A,2*i,2*i,c,INSERT_VALUES));
55: PetscCall(MatSetValue(A,2*i+1,2*i+1,c,INSERT_VALUES));
56: PetscCall(MatSetValue(A,2*i,2*i+1,s,INSERT_VALUES));
57: PetscCall(MatSetValue(A,2*i+1,2*i,-s,INSERT_VALUES));
58: }
59: if (n%2) PetscCall(MatSetValue(A,n-1,n-1,-1.0,INSERT_VALUES));
60: }
61: PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
62: PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
64: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
65: Create the eigensolver and solve the problem
66: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
68: PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
69: PetscCall(EPSSetOperators(eps,A,NULL));
70: PetscCall(EPSSetProblemType(eps,EPS_NHEP));
71: PetscCall(EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL));
72: PetscCall(EPSSetFromOptions(eps));
73: PetscCall(EPSSolve(eps));
75: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
76: Display solution and clean up
77: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
79: PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL));
80: PetscCall(EPSDestroy(&eps));
81: PetscCall(MatDestroy(&A));
82: PetscCall(PetscRandomDestroy(&rand));
83: PetscCall(SlepcFinalize());
84: return 0;
85: }
87: /*TEST
89: testset:
90: requires: complex double
91: args: -eps_type ciss -eps_all -rg_type ring -rg_ring_center 0 -rg_ring_radius 1 -rg_ring_width 0.05 -rg_ring_startangle .93 -rg_ring_endangle .07
92: filter: sed -e "s/[+-]\([0-9]\.[0-9]*i\)/+-\\1/g"
93: test:
94: suffix: 1_ring
95: args: -eps_ciss_extraction {{ritz hankel}}
97: TEST*/