Actual source code: slepcimpl.h
slepc-main 2023-10-18
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: #pragma once
13: #include <slepcsys.h>
14: #include <petsc/private/petscimpl.h>
16: /* SUBMANSEC = sys */
18: SLEPC_INTERN PetscBool SlepcBeganPetsc;
20: /*@C
21: SlepcHeaderCreate - Creates a SLEPc object
23: Input Parameters:
24: + classid - the classid associated with this object
25: . class_name - string name of class; should be static
26: . descr - string containing short description; should be static
27: . mansec - string indicating section in manual pages; should be static
28: . comm - the MPI Communicator
29: . destroy - the destroy routine for this object
30: - view - the view routine for this object
32: Output Parameter:
33: . h - the newly created object
35: Note:
36: This is equivalent to PetscHeaderCreate but makes sure that SlepcInitialize
37: has been called.
39: Level: developer
40: @*/
41: #define SlepcHeaderCreate(h,classid,class_name,descr,mansec,comm,destroy,view) \
42: ((PetscErrorCode)((!SlepcInitializeCalled && \
43: PetscError(comm,__LINE__,PETSC_FUNCTION_NAME,__FILE__,PETSC_ERR_ORDER,PETSC_ERROR_INITIAL, \
44: "Must call SlepcInitialize instead of PetscInitialize to use SLEPc classes")) || \
45: PetscHeaderCreate(h,classid,class_name,descr,mansec,comm,destroy,view)))
47: /* context for monitors of type XXXMonitorConverged */
48: struct _n_SlepcConvMon {
49: void *ctx;
50: PetscInt oldnconv; /* previous value of nconv */
51: };
53: /*
54: SlepcPrintEigenvalueASCII - Print an eigenvalue on an ASCII viewer.
55: */
56: static inline PetscErrorCode SlepcPrintEigenvalueASCII(PetscViewer viewer,PetscScalar eigr,PetscScalar eigi)
57: {
58: PetscReal re,im;
60: PetscFunctionBegin;
61: #if defined(PETSC_USE_COMPLEX)
62: re = PetscRealPart(eigr);
63: im = PetscImaginaryPart(eigr);
64: (void)eigi;
65: #else
66: re = eigr;
67: im = eigi;
68: #endif
69: /* print zero instead of tiny value */
70: if (PetscAbs(im) && PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
71: if (PetscAbs(re) && PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
72: /* print as real if imaginary part is zero */
73: if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im));
74: else PetscCall(PetscViewerASCIIPrintf(viewer,"%.5f",(double)re));
75: PetscFunctionReturn(PETSC_SUCCESS);
76: }
78: /*
79: SlepcViewEigenvector - Outputs an eigenvector xr,xi to a viewer.
80: In complex scalars only xr is written.
81: The name of xr,xi is set before writing, based on the label, the index, and the name of obj.
82: */
83: static inline PetscErrorCode SlepcViewEigenvector(PetscViewer viewer,Vec xr,Vec xi,const char *label,PetscInt index,PetscObject obj)
84: {
85: size_t count;
86: char vname[30];
87: const char *pname;
89: PetscFunctionBegin;
90: PetscCall(PetscObjectGetName(obj,&pname));
91: PetscCall(PetscSNPrintfCount(vname,sizeof(vname),"%s%s",&count,label,PetscDefined(USE_COMPLEX)?"":"r"));
92: count--;
93: PetscCall(PetscSNPrintf(vname+count,sizeof(vname)-count,"%" PetscInt_FMT "_%s",index,pname));
94: PetscCall(PetscObjectSetName((PetscObject)xr,vname));
95: PetscCall(VecView(xr,viewer));
96: #if !defined(PETSC_USE_COMPLEX)
97: vname[count-1] = 'i';
98: PetscCall(PetscObjectSetName((PetscObject)xi,vname));
99: PetscCall(VecView(xi,viewer));
100: #else
101: (void)xi;
102: #endif
103: PetscFunctionReturn(PETSC_SUCCESS);
104: }
106: /* Macros for strings with different value in real and complex */
107: #if defined(PETSC_USE_COMPLEX)
108: #define SLEPC_STRING_HERMITIAN "hermitian"
109: #else
110: #define SLEPC_STRING_HERMITIAN "symmetric"
111: #endif
113: /* Private functions that are shared by several classes */
114: SLEPC_SINGLE_LIBRARY_INTERN PetscErrorCode SlepcBasisReference_Private(PetscInt,Vec*,PetscInt*,Vec**);
115: SLEPC_SINGLE_LIBRARY_INTERN PetscErrorCode SlepcBasisDestroy_Private(PetscInt*,Vec**);
116: SLEPC_SINGLE_LIBRARY_INTERN PetscErrorCode SlepcMonitorMakeKey_Internal(const char[],PetscViewerType,PetscViewerFormat,char[]);
117: SLEPC_SINGLE_LIBRARY_INTERN PetscErrorCode PetscViewerAndFormatCreate_Internal(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**);
119: SLEPC_INTERN PetscErrorCode SlepcCitationsInitialize(void);
120: SLEPC_INTERN PetscErrorCode SlepcInitialize_DynamicLibraries(void);
121: SLEPC_INTERN PetscErrorCode SlepcInitialize_Packages(void);
123: /* Definitions needed to work with CUDA kernels */
124: #if defined(PETSC_HAVE_CUDA)
125: #include <petscdevice_cuda.h>
127: #define X_AXIS 0
128: #define Y_AXIS 1
130: #define SLEPC_TILE_SIZE_X 32
131: #define SLEPC_BLOCK_SIZE_X 128
132: #define SLEPC_TILE_SIZE_Y 32
133: #define SLEPC_BLOCK_SIZE_Y 128
135: static inline PetscErrorCode SlepcKernelSetGrid1D(PetscInt rows,dim3 *dimGrid,dim3 *dimBlock,PetscInt *dimGrid_xcount)
136: {
137: int card;
138: struct cudaDeviceProp devprop;
140: PetscFunctionBegin;
141: PetscCallCUDA(cudaGetDevice(&card));
142: PetscCallCUDA(cudaGetDeviceProperties(&devprop,card));
143: *dimGrid_xcount = 1;
145: /* X axis */
146: dimGrid->x = 1;
147: dimBlock->x = SLEPC_BLOCK_SIZE_X;
148: if (rows>SLEPC_BLOCK_SIZE_X) dimGrid->x = (rows+SLEPC_BLOCK_SIZE_X-1)/SLEPC_BLOCK_SIZE_X;
149: else dimBlock->x = rows;
150: if (dimGrid->x>(unsigned)devprop.maxGridSize[X_AXIS]) {
151: *dimGrid_xcount = (dimGrid->x+(devprop.maxGridSize[X_AXIS]-1))/devprop.maxGridSize[X_AXIS];
152: dimGrid->x = devprop.maxGridSize[X_AXIS];
153: }
154: PetscFunctionReturn(PETSC_SUCCESS);
155: }
157: static inline PetscErrorCode SlepcKernelSetGrid2DTiles(PetscInt rows,PetscInt cols,dim3 *dimGrid,dim3 *dimBlock,PetscInt *dimGrid_xcount,PetscInt *dimGrid_ycount)
158: {
159: int card;
160: struct cudaDeviceProp devprop;
162: PetscFunctionBegin;
163: PetscCallCUDA(cudaGetDevice(&card));
164: PetscCallCUDA(cudaGetDeviceProperties(&devprop,card));
165: *dimGrid_xcount = *dimGrid_ycount = 1;
167: /* X axis */
168: dimGrid->x = 1;
169: dimBlock->x = SLEPC_BLOCK_SIZE_X;
170: if (rows>SLEPC_BLOCK_SIZE_X*SLEPC_TILE_SIZE_X) dimGrid->x = (rows+SLEPC_BLOCK_SIZE_X*SLEPC_TILE_SIZE_X-1)/(SLEPC_BLOCK_SIZE_X*SLEPC_TILE_SIZE_X);
171: else dimBlock->x = (rows+SLEPC_TILE_SIZE_X-1)/SLEPC_TILE_SIZE_X;
172: if (dimGrid->x>(unsigned)devprop.maxGridSize[X_AXIS]) {
173: *dimGrid_xcount = (dimGrid->x+(devprop.maxGridSize[X_AXIS]-1))/devprop.maxGridSize[X_AXIS];
174: dimGrid->x = devprop.maxGridSize[X_AXIS];
175: }
177: /* Y axis */
178: dimGrid->y = 1;
179: dimBlock->y = SLEPC_BLOCK_SIZE_Y;
180: if (cols>SLEPC_BLOCK_SIZE_Y*SLEPC_TILE_SIZE_Y) dimGrid->y = (cols+SLEPC_BLOCK_SIZE_Y*SLEPC_TILE_SIZE_Y-1)/(SLEPC_BLOCK_SIZE_Y*SLEPC_TILE_SIZE_Y);
181: else dimBlock->y = (cols+SLEPC_TILE_SIZE_Y-1)/SLEPC_TILE_SIZE_Y;
182: if (dimGrid->y>(unsigned)devprop.maxGridSize[Y_AXIS]) {
183: *dimGrid_ycount = (dimGrid->y+(devprop.maxGridSize[Y_AXIS]-1))/devprop.maxGridSize[Y_AXIS];
184: dimGrid->y = devprop.maxGridSize[Y_AXIS];
185: }
186: PetscFunctionReturn(PETSC_SUCCESS);
187: }
188: #undef X_AXIS
189: #undef Y_AXIS
190: #endif