Actual source code: slepcinit.c
2: #include slepc.h
3: #include slepceps.h
4: #include slepcst.h
8: /*
9: SlepcPrintVersion - Prints SLEPc version info.
11: Collective on MPI_Comm
12: */
13: PetscErrorCode SlepcPrintVersion(MPI_Comm comm)
14: {
15: int info = 0;
16:
19: info = (*PetscHelpPrintf)(comm,"--------------------------------------------\
20: ------------------------------\n"); CHKERRQ(info);
21: info = (*PetscHelpPrintf)(comm,"\t %s\n",SLEPC_VERSION_NUMBER); CHKERRQ(info);
22: info = (*PetscHelpPrintf)(comm,"%s",SLEPC_AUTHOR_INFO); CHKERRQ(info);
23: info = (*PetscHelpPrintf)(comm,"See docs/index.html for help. \n"); CHKERRQ(info);
24: #if !defined(PARCH_win32)
25: info = (*PetscHelpPrintf)(comm,"SLEPc libraries linked from %s\n",SLEPC_LIB_DIR); CHKERRQ(info);
26: #endif
27: info = (*PetscHelpPrintf)(comm,"--------------------------------------------\
28: ------------------------------\n"); CHKERRQ(info);
30: PetscFunctionReturn(info);
31: }
35: /*
36: SlepcPrintHelpIntro - Prints introductory SLEPc help info.
38: Collective on MPI_Comm
39: */
40: PetscErrorCode SlepcPrintHelpIntro(MPI_Comm comm)
41: {
42: int info = 0;
43:
46: info = (*PetscHelpPrintf)(comm,"--------------------------------------------\
47: ------------------------------\n"); CHKERRQ(info);
48: info = (*PetscHelpPrintf)(comm,"SLEPc help information includes that for the PETSc libraries, which provide\n"); CHKERRQ(info);
49: info = (*PetscHelpPrintf)(comm,"low-level system infrastructure and linear algebra tools.\n"); CHKERRQ(info);
50: info = (*PetscHelpPrintf)(comm,"--------------------------------------------\
51: ------------------------------\n"); CHKERRQ(info);
53: PetscFunctionReturn(info);
54: }
58: PetscEvent EPS_SetUp, EPS_Solve, ST_SetUp, ST_Apply, ST_ApplyB, ST_ApplyNoB, EPS_Orthogonalization, ST_InnerProduct;
60: /*
61: SlepcRegisterEvents - Registers SLEPc events for use in performance logging.
62: */
63: PetscErrorCode SlepcRegisterEvents()
64: {
65: int info = 0;
66:
69: info = PetscLogEventRegister(&EPS_SetUp,"EPSSetUp",PETSC_NULL); CHKERRQ(info);
70: info = PetscLogEventRegister(&EPS_Solve,"EPSSolve",PETSC_NULL); CHKERRQ(info);
71: info = PetscLogEventRegister(&ST_SetUp,"STSetUp",PETSC_NULL); CHKERRQ(info);
72: info = PetscLogEventRegister(&ST_Apply,"STApply",PETSC_NULL); CHKERRQ(info);
73: info = PetscLogEventRegister(&ST_ApplyB,"STApplyB",PETSC_NULL); CHKERRQ(info);
74: info = PetscLogEventRegister(&ST_ApplyNoB,"STApplyNoB",PETSC_NULL); CHKERRQ(info);
75: info = PetscLogEventRegister(&EPS_Orthogonalization,"EPSOrthogonalization",PETSC_NULL); CHKERRQ(info);
76: info = PetscLogEventRegister(&ST_InnerProduct,"STInnerProduct",PETSC_NULL); CHKERRQ(info);
78: PetscFunctionReturn(info);
79: }
81: /* ------------------------Nasty global variables -------------------------------*/
82: /*
83: Indicates whether SLEPc started PETSc, or whether it was
84: already started before SLEPc was initialized.
85: */
86: PetscTruth SlepcBeganPetsc = PETSC_FALSE;
87: PetscTruth SlepcInitializeCalled = PETSC_FALSE;
88: PetscCookie EPS_COOKIE = 0;
89: PetscCookie ST_COOKIE = 0;
91: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
92: extern PetscDLLibraryList DLLibrariesLoaded;
93: #endif
97: /*@C
98: SlepcInitialize - Initializes the SLEPc library. SlepcInitialize() calls
99: PetscInitialize() if that has not been called yet, so this routine should
100: always be called near the beginning of your program.
102: Collective on MPI_COMM_WORLD or PETSC_COMM_WORLD if it has been set
104: Input Parameters:
105: + argc - count of number of command line arguments
106: . args - the command line arguments
107: . file - [optional] PETSc database file, defaults to ~username/.petscrc
108: (use PETSC_NULL for default)
109: - help - [optional] Help message to print, use PETSC_NULL for no message
111: Fortran Note:
112: Fortran syntax is very similar to that of PetscInitialize()
113:
114: Level: beginner
116: .seealso: SlepcInitializeFortran(), SlepcFinalize(), PetscInitialize()
117: @*/
118: PetscErrorCode SlepcInitialize(int *argc,char ***args,char file[],const char help[])
119: {
120: PetscErrorCode ierr,info=0;
121: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
122: char libs[PETSC_MAX_PATH_LEN],dlib[PETSC_MAX_PATH_LEN];
123: PetscTruth found;
124: #endif
128: if (SlepcInitializeCalled==PETSC_TRUE) {
129: return(0);
130: }
132: #if !defined(PARCH_t3d)
133: info = PetscSetHelpVersionFunctions(SlepcPrintHelpIntro,SlepcPrintVersion);CHKERRQ(info);
134: #endif
136: if (!PetscInitializeCalled) {
137: info = PetscInitialize(argc,args,file,help);CHKERRQ(info);
138: SlepcBeganPetsc = PETSC_TRUE;
139: }
141: EPS_COOKIE = 0;
142: PetscLogClassRegister(&EPS_COOKIE,"Eigenproblem Solver");
143: ST_COOKIE = 0;
144: PetscLogClassRegister(&ST_COOKIE,"Spectral Transform");
146: /*
147: Load the dynamic libraries
148: */
150: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
151: PetscStrcpy(libs,SLEPC_LIB_DIR);
152: PetscStrcat(libs,"/libslepc");
153: PetscDLLibraryRetrieve(PETSC_COMM_WORLD,libs,dlib,1024,&found);
154: if (found) {
155: PetscDLLibraryAppend(PETSC_COMM_WORLD,&DLLibrariesLoaded,libs);
156: } else {
157: SETERRQ1(1,"Unable to locate SLEPc dynamic library %s \n You cannot move the dynamic libraries!\n or remove USE_DYNAMIC_LIBRARIES from ${PETSC_DIR}/bmake/$PETSC_ARCH/petscconf.h\n and rebuild libraries before moving",libs);
158: }
159: #else
161: EPSRegisterAll(PETSC_NULL);
162: STRegisterAll(PETSC_NULL);
164: #endif
166: /*
167: Register SLEPc events
168: */
169: info = SlepcRegisterEvents();CHKERRQ(info);
170: SlepcInitializeCalled = PETSC_TRUE;
171: PetscLogInfo(0,"SlepcInitialize: SLEPc successfully started\n");
172: PetscFunctionReturn(info);
173: }
177: /*@
178: SlepcFinalize - Checks for options to be called at the conclusion
179: of the SLEPc program and calls PetscFinalize().
181: Collective on PETSC_COMM_WORLD
183: Level: beginner
185: .seealso: SlepcInitialize(), PetscFinalize()
186: @*/
187: PetscErrorCode SlepcFinalize(void)
188: {
189: PetscErrorCode info=0;
190:
192: PetscLogInfo(0,"SlepcFinalize: SLEPc successfully ended!\n");
194: if (SlepcBeganPetsc) {
195: info = PetscFinalize();CHKERRQ(info);
196: }
198: SlepcInitializeCalled = PETSC_FALSE;
200: PetscFunctionReturn(info);
201: }