slepc4py.SLEPc.PEP#
- class slepc4py.SLEPc.PEP#
Bases:
Object
PEP.
Enumerations
PEP basis type for the representation of the polynomial.
PEP CISS extraction technique.
PEP convergence test.
PEP convergence reasons.
PEP error type to assess accuracy of computed solutions.
PEP extraction strategy used.
PEP type of projection to be used in the Jacobi-Davidson solver.
PEP problem type.
PEP refinement strategy.
PEP scheme for solving linear systems during iterative refinement.
PEP scaling strategy.
PEP stopping test.
PEP type.
PEP desired part of spectrum.
Methods Summary
appendOptionsPrefix
([prefix])Append to the prefix used for searching for all PEP options in the database.
Clear all monitors for a
PEP
object.computeError
(i[, etype])Compute the error associated with the i-th computed eigenpair.
create
([comm])Create the PEP object.
destroy
()Destroy the PEP object.
errorView
([etype, viewer])Display the errors associated with the computed solution.
getBV
()Get the basis vectors object associated to the eigensolver.
getBasis
()Get the type of polynomial basis used.
Get the extraction technique used in the CISS solver.
Get the array of linear solver objects associated with the CISS solver.
Get the values of various refinement parameters in the CISS solver.
Get the values of various size parameters in the CISS solver.
Get the values of various threshold parameters in the CISS solver.
Get the number of converged eigenpairs.
Get the reason why the
solve()
iteration was stopped.Get the method used to compute the error estimate used in the convergence test.
getDS
()Get the direct solver associated to the eigensolver.
Get the number of eigenvalues to compute and the dimension of the subspace.
getEigenpair
(i[, Vr, Vi])Get the i-th solution of the eigenproblem as computed by
solve()
.Get the error estimate associated to the i-th computed eigenpair.
Get the extraction technique used by the
PEP
object.Get the computational interval for spectrum slicing.
Get the current iteration number.
getJDFix
()Get threshold for changing the target in the correction equation.
Get the maximum allowed value of the minimality index.
Get the type of projection to be used in the Jacobi-Davidson solver.
Get the restart parameter used in the Jacobi-Davidson method.
Get the flag for reusing the preconditioner.
Get the eigensolver object associated to the polynomial eigenvalue solver.
Get if the matrices A and B for the linearization are built explicitly.
Get the coeffs.
Get the list of monitor functions.
Get the matrices associated with the eigenvalue problem.
Get the prefix used for searching for all PEP options in the database.
Get the problem type from the PEP object.
Get the locking flag used in the Q-Arnoldi method.
Get the restart parameter used in the Q-Arnoldi method.
getRG
()Get the region object associated to the eigensolver.
Get the refinement strategy used by the PEP object.
Get the
KSP
object used by the eigensolver in the refinement phase.getST
()Get the
ST
object associated to the eigensolver object.Get the flag for the eigenvalue type check in spectrum slicing.
Get the flag that enforces zero detection in spectrum slicing.
Get the dimensions used for each subsolve step.
Get the values of the shifts and their corresponding inertias.
Get the coefficients that define the linearization of a quadratic eigenproblem.
Get the locking flag used in the STOAR method.
getScale
([Dl, Dr])Get the strategy used for scaling the polynomial eigenproblem.
Get the stopping function.
Get the locking flag used in the TOAR method.
Get the restart parameter used in the TOAR method.
Get the value of the target.
Get the tolerance and maximum iteration count.
Get the flag indicating whether all residual norms must be computed.
getType
()Get the PEP type of this object.
Get which portion of the spectrum is to be sought.
reset
()Reset the PEP object.
setBV
(bv)Set a basis vectors object associated to the eigensolver.
setBasis
(basis)Set the type of polynomial basis used.
setCISSExtraction
(extraction)Set the extraction technique used in the CISS solver.
setCISSRefinement
([inner, blsize])Set the values of various refinement parameters in the CISS solver.
setCISSSizes
([ip, bs, ms, npart, bsmax, ...])Set the values of various size parameters in the CISS solver.
setCISSThreshold
([delta, spur])Set the values of various threshold parameters in the CISS solver.
setConvergenceTest
(conv)Set how to compute the error estimate used in the convergence test.
setDS
(ds)Set a direct solver object associated to the eigensolver.
setDimensions
([nev, ncv, mpd])Set the number of eigenvalues to compute and the dimension of the subspace.
setExtract
(extract)Set the extraction strategy to be used.
Set PEP options from the options database.
setInitialSpace
(space)Set the initial space from which the eigensolver starts to iterate.
setInterval
(inta, intb)Set the computational interval for spectrum slicing.
setJDFix
(fix)Set the threshold for changing the target in the correction equation.
setJDMinimalityIndex
(flag)Set the maximum allowed value for the minimality index.
setJDProjection
(proj)Set the type of projection to be used in the Jacobi-Davidson solver.
setJDRestart
(keep)Set the restart parameter for the Jacobi-Davidson method.
setJDReusePreconditioner
(flag)Set a flag indicating whether the preconditioner must be reused or not.
setLinearEPS
(eps)Set an eigensolver object associated to the polynomial eigenvalue solver.
setLinearExplicitMatrix
(flag)Set flag to explicitly build the matrices A and B.
setLinearLinearization
([alpha, beta])Set the coefficients that define the linearization of a quadratic eigenproblem.
setMonitor
(monitor[, args, kargs])Append a monitor function to the list of monitors.
setOperators
(operators)Set the matrices associated with the eigenvalue problem.
setOptionsPrefix
([prefix])Set the prefix used for searching for all PEP options in the database.
setProblemType
(problem_type)Set the type of the eigenvalue problem.
setQArnoldiLocking
(lock)Toggle between locking and non-locking variants of the Q-Arnoldi method.
setQArnoldiRestart
(keep)Set the restart parameter for the Q-Arnoldi method.
setRG
(rg)Set a region object associated to the eigensolver.
setRefine
(ref[, npart, tol, its, scheme])Set the refinement strategy used by the PEP object.
setST
(st)Set a spectral transformation object associated to the eigensolver.
Set flag to check if all eigenvalues have the same definite type.
setSTOARDetectZeros
(detect)Set flag to enforce detection of zeros during the factorizations.
setSTOARDimensions
([nev, ncv, mpd])Set the dimensions used for each subsolve step.
setSTOARLinearization
([alpha, beta])Set the coefficients that define the linearization of a quadratic eigenproblem.
setSTOARLocking
(lock)Toggle between locking and non-locking variants of the STOAR method.
setScale
(scale[, alpha, Dl, Dr, its, lbda])Set the scaling strategy to be used.
setStoppingTest
(stopping[, args, kargs])Set a function to decide when to stop the outer iteration of the eigensolver.
setTOARLocking
(lock)Toggle between locking and non-locking variants of the TOAR method.
setTOARRestart
(keep)Set the restart parameter for the TOAR method.
setTarget
(target)Set the value of the target.
setTolerances
([tol, max_it])Set the tolerance and maximum iteration count.
setTrackAll
(trackall)Set flag to compute the residual of all approximate eigenpairs.
setType
(pep_type)Set the particular solver to be used in the PEP object.
setUp
()Set up all the necessary internal data structures.
setWhichEigenpairs
(which)Set which portion of the spectrum is to be sought.
solve
()Solve the eigensystem.
valuesView
([viewer])Display the computed eigenvalues in a viewer.
vectorsView
([viewer])Output computed eigenvectors to a viewer.
view
([viewer])Print the PEP data structure.
Attributes Summary
The basis vectors (BV) object associated.
The direct solver (DS) object associated.
The type of extraction technique to be employed.
The maximum iteration count.
The type of the eigenvalue problem.
The region (RG) object associated.
The spectral transformation (ST) object associated.
The value of the target.
The tolerance.
Compute the residual norm of all approximate eigenpairs.
The portion of the spectrum to be sought.
Methods Documentation
- appendOptionsPrefix(prefix=None)#
Append to the prefix used for searching for all PEP options in the database.
Logically collective.
- cancelMonitor()#
Clear all monitors for a
PEP
object.Logically collective.
Source code at slepc4py/SLEPc/PEP.pyx:1162
- Return type:
- computeError(i, etype=None)#
Compute the error associated with the i-th computed eigenpair.
Collective.
Compute the error (based on the residual norm) associated with the i-th computed eigenpair.
- Parameters:
- Returns:
The error bound, computed in various ways from the residual norm \(||P(l)x||_2\) where \(l\) is the eigenvalue and \(x\) is the eigenvector.
- Return type:
Notes
The index
i
should be a value between0
andnconv-1
(seegetConverged()
).
- create(comm=None)#
Create the PEP object.
Collective.
- destroy()#
Destroy the PEP object.
Collective.
Source code at slepc4py/SLEPc/PEP.pyx:258
- Return type:
- errorView(etype=None, viewer=None)#
Display the errors associated with the computed solution.
Collective.
Display the errors and the eigenvalues.
- Parameters:
viewer (petsc4py.PETSc.Viewer | None) – Visualization context; if not provided, the standard output is used.
- Return type:
Notes
By default, this function checks the error of all eigenpairs and prints the eigenvalues if all of them are below the requested tolerance. If the viewer has format
ASCII_INFO_DETAIL
then a table with eigenvalues and corresponding errors is printed.
- getBV()#
Get the basis vectors object associated to the eigensolver.
Not collective.
- Returns:
The basis vectors context.
- Return type:
- getBasis()#
Get the type of polynomial basis used.
Not collective.
Get the type of polynomial basis used to describe the polynomial eigenvalue problem.
- Returns:
The basis that was previously set.
- Return type:
- getCISSExtraction()#
Get the extraction technique used in the CISS solver.
Not collective.
- Returns:
The extraction technique.
- Return type:
- getCISSKSPs()#
Get the array of linear solver objects associated with the CISS solver.
Collective.
- Returns:
The linear solver objects.
- Return type:
Notes
The number of
petsc4py.PETSc.KSP
solvers is equal to the number of integration points divided by the number of partitions. This value is halved in the case of real matrices with a region centered at the real axis.
- getCISSRefinement()#
Get the values of various refinement parameters in the CISS solver.
Not collective.
- Returns:
- Return type:
- getCISSSizes()#
Get the values of various size parameters in the CISS solver.
Not collective.
- Returns:
- Return type:
- getCISSThreshold()#
Get the values of various threshold parameters in the CISS solver.
Not collective.
- Returns:
- Return type:
- getConverged()#
Get the number of converged eigenpairs.
Not collective.
- Returns:
Number of converged eigenpairs.
- Return type:
- getConvergedReason()#
Get the reason why the
solve()
iteration was stopped.Not collective.
- Returns:
Negative value indicates diverged, positive value converged.
- Return type:
- getConvergenceTest()#
Get the method used to compute the error estimate used in the convergence test.
Not collective.
- Returns:
The method used to compute the error estimate used in the convergence test.
- Return type:
- getDS()#
Get the direct solver associated to the eigensolver.
Not collective.
- Returns:
The direct solver context.
- Return type:
- getDimensions()#
Get the number of eigenvalues to compute and the dimension of the subspace.
Not collective.
- Returns:
- Return type:
- getEigenpair(i, Vr=None, Vi=None)#
Get the i-th solution of the eigenproblem as computed by
solve()
.Collective.
The solution consists of both the eigenvalue and the eigenvector.
- Parameters:
- Returns:
The computed eigenvalue.
- Return type:
- getErrorEstimate(i)#
Get the error estimate associated to the i-th computed eigenpair.
Not collective.
- Parameters:
i (int) – Index of the solution to be considered.
- Returns:
Error estimate.
- Return type:
- getExtract()#
Get the extraction technique used by the
PEP
object.Not collective.
- Returns:
The extraction strategy.
- Return type:
- getInterval()#
Get the computational interval for spectrum slicing.
Not collective.
- Returns:
- Return type:
Notes
If the interval was not set by the user, then zeros are returned.
- getIterationNumber()#
Get the current iteration number.
Not collective.
If the call to
solve()
is complete, then it returns the number of iterations carried out by the solution method.- Returns:
Iteration number.
- Return type:
- getJDFix()#
Get threshold for changing the target in the correction equation.
Not collective.
- Returns:
The threshold for changing the target.
- Return type:
- getJDMinimalityIndex()#
Get the maximum allowed value of the minimality index.
Not collective.
- Returns:
The maximum minimality index.
- Return type:
- getJDProjection()#
Get the type of projection to be used in the Jacobi-Davidson solver.
Not collective.
- Returns:
The type of projection.
- Return type:
- getJDRestart()#
Get the restart parameter used in the Jacobi-Davidson method.
Not collective.
- Returns:
The number of vectors to be kept at restart.
- Return type:
- getJDReusePreconditioner()#
Get the flag for reusing the preconditioner.
Not collective.
- Returns:
The reuse flag.
- Return type:
- getLinearEPS()#
Get the eigensolver object associated to the polynomial eigenvalue solver.
Collective.
- Returns:
The linear eigensolver.
- Return type:
- getLinearExplicitMatrix()#
Get if the matrices A and B for the linearization are built explicitly.
Not collective.
Get the flag indicating if the matrices A and B for the linearization are built explicitly.
- Returns:
Boolean flag indicating if the matrices are built explicitly.
- Return type:
- getLinearLinearization()#
Get the coeffs. defining the linearization of a quadratic eigenproblem.
Not collective.
Return the coefficients that define the linearization of a quadratic eigenproblem.
- Returns:
- Return type:
- getMonitor()#
Get the list of monitor functions.
Not collective.
Source code at slepc4py/SLEPc/PEP.pyx:1154
- Return type:
- getOperators()#
Get the matrices associated with the eigenvalue problem.
Collective.
- Returns:
The matrices associated with the eigensystem.
- Return type:
- getOptionsPrefix()#
Get the prefix used for searching for all PEP options in the database.
Not collective.
- Returns:
The prefix string set for this PEP object.
- Return type:
- getProblemType()#
Get the problem type from the PEP object.
Not collective.
- Returns:
The problem type that was previously set.
- Return type:
- getQArnoldiLocking()#
Get the locking flag used in the Q-Arnoldi method.
Not collective.
- Returns:
The locking flag.
- Return type:
- getQArnoldiRestart()#
Get the restart parameter used in the Q-Arnoldi method.
Not collective.
- Returns:
The number of vectors to be kept at restart.
- Return type:
- getRG()#
Get the region object associated to the eigensolver.
Not collective.
- Returns:
The region context.
- Return type:
- getRefine()#
Get the refinement strategy used by the PEP object.
Not collective.
Get the refinement strategy used by the PEP object, and the associated parameters.
- Returns:
ref (
Refine
) – The refinement type.npart (
int
) – The number of partitions of the communicator.tol (
float
) – The convergence tolerance.its (
int
) – The maximum number of refinement iterations.scheme (
RefineScheme
) – Scheme for solving linear systems
- Return type:
- getRefineKSP()#
Get the
KSP
object used by the eigensolver in the refinement phase.Collective.
- Returns:
The linear solver object.
- Return type:
- getST()#
Get the
ST
object associated to the eigensolver object.Not collective.
Get the spectral transformation object associated to the eigensolver object.
- Returns:
The spectral transformation.
- Return type:
- getSTOARCheckEigenvalueType()#
Get the flag for the eigenvalue type check in spectrum slicing.
Not collective.
- Returns:
Whether the eigenvalue type is checked or not.
- Return type:
- getSTOARDetectZeros()#
Get the flag that enforces zero detection in spectrum slicing.
Not collective.
- Returns:
The zero detection flag.
- Return type:
- getSTOARDimensions()#
Get the dimensions used for each subsolve step.
Not collective.
Get the dimensions used for each subsolve step in case of doing spectrum slicing for a computational interval.
- Returns:
- Return type:
- getSTOARInertias()#
Get the values of the shifts and their corresponding inertias.
Not collective.
Get the values of the shifts and their corresponding inertias in case of doing spectrum slicing for a computational interval.
- Returns:
- Return type:
- getSTOARLinearization()#
Get the coefficients that define the linearization of a quadratic eigenproblem.
Not collective.
- Returns:
- Return type:
- getSTOARLocking()#
Get the locking flag used in the STOAR method.
Not collective.
- Returns:
The locking flag.
- Return type:
- getScale(Dl=None, Dr=None)#
Get the strategy used for scaling the polynomial eigenproblem.
Not collective.
- Parameters:
Dl (petsc4py.PETSc.Vec | None) – Placeholder for the returned left diagonal matrix.
Dr (petsc4py.PETSc.Vec | None) – Placeholder for the returned right diagonal matrix.
- Returns:
- Return type:
- getStoppingTest()#
Get the stopping function.
Not collective.
Source code at slepc4py/SLEPc/PEP.pyx:1123
- Return type:
- getTOARLocking()#
Get the locking flag used in the TOAR method.
Not collective.
- Returns:
The locking flag.
- Return type:
- getTOARRestart()#
Get the restart parameter used in the TOAR method.
Not collective.
- Returns:
The number of vectors to be kept at restart.
- Return type:
- getTarget()#
Get the value of the target.
Not collective.
- Returns:
The value of the target.
- Return type:
Notes
If the target was not set by the user, then zero is returned.
- getTolerances()#
Get the tolerance and maximum iteration count.
Not collective.
Get the tolerance and maximum iteration count used by the default PEP convergence tests.
- Returns:
- Return type:
- getTrackAll()#
Get the flag indicating whether all residual norms must be computed.
Not collective.
- Returns:
Whether the solver compute all residuals or not.
- Return type:
- getType()#
Get the PEP type of this object.
Not collective.
- Returns:
The solver currently being used.
- Return type:
- getWhichEigenpairs()#
Get which portion of the spectrum is to be sought.
Not collective.
- Returns:
The portion of the spectrum to be sought by the solver.
- Return type:
- reset()#
Reset the PEP object.
Collective.
Source code at slepc4py/SLEPc/PEP.pyx:268
- Return type:
- setBV(bv)#
Set a basis vectors object associated to the eigensolver.
Collective.
- setBasis(basis)#
Set the type of polynomial basis used.
Logically collective.
Set the type of polynomial basis used to describe the polynomial eigenvalue problem.
- setCISSExtraction(extraction)#
Set the extraction technique used in the CISS solver.
Logically collective.
- Parameters:
extraction (CISSExtraction) – The extraction technique.
- Return type:
- setCISSRefinement(inner=None, blsize=None)#
Set the values of various refinement parameters in the CISS solver.
Logically collective.
- Parameters:
- Return type:
- setCISSSizes(ip=None, bs=None, ms=None, npart=None, bsmax=None, realmats=False)#
Set the values of various size parameters in the CISS solver.
Logically collective.
- Parameters:
- Return type:
Notes
The default number of partitions is 1. This means the internal
petsc4py.PETSc.KSP
object is shared among all processes of thePEP
communicator. Otherwise, the communicator is split into npart communicators, so thatnpart
petsc4py.PETSc.KSP
solves proceed simultaneously.
- setCISSThreshold(delta=None, spur=None)#
Set the values of various threshold parameters in the CISS solver.
Logically collective.
- Parameters:
- Return type:
- setConvergenceTest(conv)#
Set how to compute the error estimate used in the convergence test.
Logically collective.
- Parameters:
conv (Conv) – The method used to compute the error estimate used in the convergence test.
- Return type:
- setDS(ds)#
Set a direct solver object associated to the eigensolver.
Collective.
- setDimensions(nev=None, ncv=None, mpd=None)#
Set the number of eigenvalues to compute and the dimension of the subspace.
Logically collective.
- Parameters:
- Return type:
- setExtract(extract)#
Set the extraction strategy to be used.
Logically collective.
- setFromOptions()#
Set PEP options from the options database.
Collective.
This routine must be called before
setUp()
if the user is to be allowed to set the solver type.Source code at slepc4py/SLEPc/PEP.pyx:368
- Return type:
- setInitialSpace(space)#
Set the initial space from which the eigensolver starts to iterate.
Collective.
- setInterval(inta, intb)#
Set the computational interval for spectrum slicing.
Logically collective.
- Parameters:
- Return type:
Notes
Spectrum slicing is a technique employed for computing all eigenvalues of symmetric quadratic eigenproblems in a given interval. This function provides the interval to be considered. It must be used in combination with
PEP.Which.ALL
, seesetWhichEigenpairs()
.
- setJDFix(fix)#
Set the threshold for changing the target in the correction equation.
Logically collective.
Notes
The target in the correction equation is fixed at the first iterations. When the norm of the residual vector is lower than the fix value, the target is set to the corresponding eigenvalue.
- setJDMinimalityIndex(flag)#
Set the maximum allowed value for the minimality index.
Logically collective.
- setJDProjection(proj)#
Set the type of projection to be used in the Jacobi-Davidson solver.
Logically collective.
- Parameters:
proj (JDProjection) – The type of projection.
- Return type:
- setJDRestart(keep)#
Set the restart parameter for the Jacobi-Davidson method.
Logically collective.
Set the restart parameter for the Jacobi-Davidson method, in particular the proportion of basis vectors that must be kept after restart.
Notes
Allowed values are in the range [0.1,0.9]. The default is 0.5.
- setJDReusePreconditioner(flag)#
Set a flag indicating whether the preconditioner must be reused or not.
Logically collective.
- setLinearEPS(eps)#
Set an eigensolver object associated to the polynomial eigenvalue solver.
Collective.
- setLinearExplicitMatrix(flag)#
Set flag to explicitly build the matrices A and B.
Logically collective.
Toggle if the matrices A and B for the linearization of the problem must be built explicitly.
- Parameters:
flag (bool) – Boolean flag indicating if the matrices are built explicitly.
- Return type:
- setLinearLinearization(alpha=1.0, beta=0.0)#
Set the coefficients that define the linearization of a quadratic eigenproblem.
Logically collective.
Set the coefficients that define the linearization of a quadratic eigenproblem.
- Parameters:
- Return type:
- setMonitor(monitor, args=None, kargs=None)#
Append a monitor function to the list of monitors.
Logically collective.
- setOperators(operators)#
Set the matrices associated with the eigenvalue problem.
Collective.
- setOptionsPrefix(prefix=None)#
Set the prefix used for searching for all PEP options in the database.
Logically collective.
- setProblemType(problem_type)#
Set the type of the eigenvalue problem.
Logically collective.
- Parameters:
problem_type (ProblemType) – The problem type to be set.
- Return type:
- setQArnoldiLocking(lock)#
Toggle between locking and non-locking variants of the Q-Arnoldi method.
Logically collective.
Notes
The default is to lock converged eigenpairs when the method restarts. This behaviour can be changed so that all directions are kept in the working subspace even if already converged to working accuracy (the non-locking variant).
- setQArnoldiRestart(keep)#
Set the restart parameter for the Q-Arnoldi method.
Logically collective.
Set the restart parameter for the Q-Arnoldi method, in particular the proportion of basis vectors that must be kept after restart.
Notes
Allowed values are in the range [0.1,0.9]. The default is 0.5.
- setRG(rg)#
Set a region object associated to the eigensolver.
Collective.
- setRefine(ref, npart=None, tol=None, its=None, scheme=None)#
Set the refinement strategy used by the PEP object.
Logically collective.
Set the refinement strategy used by the PEP object, and the associated parameters.
- Parameters:
- Return type:
- setST(st)#
Set a spectral transformation object associated to the eigensolver.
Collective.
- setSTOARCheckEigenvalueType(flag)#
Set flag to check if all eigenvalues have the same definite type.
Logically collective.
Set a flag to check that all the eigenvalues obtained throughout the spectrum slicing computation have the same definite type.
- setSTOARDetectZeros(detect)#
Set flag to enforce detection of zeros during the factorizations.
Logically collective.
Set a flag to enforce detection of zeros during the factorizations throughout the spectrum slicing computation.
Notes
A zero in the factorization indicates that a shift coincides with an eigenvalue.
This flag is turned off by default, and may be necessary in some cases. This feature currently requires an external package for factorizations with support for zero detection, e.g. MUMPS.
- setSTOARDimensions(nev=None, ncv=None, mpd=None)#
Set the dimensions used for each subsolve step.
Logically collective.
Set the dimensions used for each subsolve step in case of doing spectrum slicing for a computational interval. The meaning of the parameters is the same as in
setDimensions()
.- Parameters:
- Return type:
- setSTOARLinearization(alpha=1.0, beta=0.0)#
Set the coefficients that define the linearization of a quadratic eigenproblem.
Logically collective.
- Parameters:
- Return type:
- setSTOARLocking(lock)#
Toggle between locking and non-locking variants of the STOAR method.
Logically collective.
Notes
The default is to lock converged eigenpairs when the method restarts. This behaviour can be changed so that all directions are kept in the working subspace even if already converged to working accuracy (the non-locking variant).
- setScale(scale, alpha=None, Dl=None, Dr=None, its=None, lbda=None)#
Set the scaling strategy to be used.
Collective.
Set the scaling strategy to be used for scaling the polynomial problem before attempting to solve.
- Parameters:
scale (Scale) – The scaling strategy.
Dl (petsc4py.PETSc.Vec | None) – The left diagonal matrix.
Dr (petsc4py.PETSc.Vec | None) – The right diagonal matrix.
its (int | None) – The number of iteration of diagonal scaling.
lbda (float | None) – Approximation of the wanted eigenvalues (modulus).
- Return type:
- setStoppingTest(stopping, args=None, kargs=None)#
Set a function to decide when to stop the outer iteration of the eigensolver.
Logically collective.
- setTOARLocking(lock)#
Toggle between locking and non-locking variants of the TOAR method.
Logically collective.
Notes
The default is to lock converged eigenpairs when the method restarts. This behaviour can be changed so that all directions are kept in the working subspace even if already converged to working accuracy (the non-locking variant).
- setTOARRestart(keep)#
Set the restart parameter for the TOAR method.
Logically collective.
Set the restart parameter for the TOAR method, in particular the proportion of basis vectors that must be kept after restart.
Notes
Allowed values are in the range [0.1,0.9]. The default is 0.5.
- setTarget(target)#
Set the value of the target.
Logically collective.
Notes
The target is a scalar value used to determine the portion of the spectrum of interest. It is used in combination with
setWhichEigenpairs()
.
- setTolerances(tol=None, max_it=None)#
Set the tolerance and maximum iteration count.
Logically collective.
Set the tolerance and maximum iteration count used by the default PEP convergence tests.
- Parameters:
- Return type:
- setTrackAll(trackall)#
Set flag to compute the residual of all approximate eigenpairs.
Logically collective.
Set if the solver must compute the residual of all approximate eigenpairs or not.
- setType(pep_type)#
Set the particular solver to be used in the PEP object.
Logically collective.
- setUp()#
Set up all the necessary internal data structures.
Collective.
Set up all the internal data structures necessary for the execution of the eigensolver.
Source code at slepc4py/SLEPc/PEP.pyx:1173
- Return type:
- setWhichEigenpairs(which)#
Set which portion of the spectrum is to be sought.
Logically collective.
- Parameters:
which (Which) – The portion of the spectrum to be sought by the solver.
- Return type:
- solve()#
Solve the eigensystem.
Collective.
Source code at slepc4py/SLEPc/PEP.pyx:1184
- Return type:
- valuesView(viewer=None)#
Display the computed eigenvalues in a viewer.
Collective.
- vectorsView(viewer=None)#
Output computed eigenvectors to a viewer.
Collective.
- view(viewer=None)#
Print the PEP data structure.
Collective.
Attributes Documentation
- bv#
The basis vectors (BV) object associated.
- ds#
The direct solver (DS) object associated.
- extract#
The type of extraction technique to be employed.
- max_it#
The maximum iteration count.
- problem_type#
The type of the eigenvalue problem.
- rg#
The region (RG) object associated.
- st#
The spectral transformation (ST) object associated.
- target#
The value of the target.
- tol#
The tolerance.
- track_all#
Compute the residual norm of all approximate eigenpairs.
- which#
The portion of the spectrum to be sought.