slepc4py.SLEPc.PEP#

class slepc4py.SLEPc.PEP#

Bases: Object

PEP.

Enumerations

Basis

PEP basis type for the representation of the polynomial.

CISSExtraction

PEP CISS extraction technique.

Conv

PEP convergence test.

ConvergedReason

PEP convergence reasons.

ErrorType

PEP error type to assess accuracy of computed solutions.

Extract

PEP extraction strategy used.

JDProjection

PEP type of projection to be used in the Jacobi-Davidson solver.

ProblemType

PEP problem type.

Refine

PEP refinement strategy.

RefineScheme

PEP scheme for solving linear systems during iterative refinement.

Scale

PEP scaling strategy.

Stop

PEP stopping test.

Type

PEP type.

Which

PEP desired part of spectrum.

Methods Summary

appendOptionsPrefix([prefix])

Append to the prefix used for searching for all PEP options in the database.

cancelMonitor()

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.

getCISSExtraction()

Get the extraction technique used in the CISS solver.

getCISSKSPs()

Get the array of linear solver objects associated with the CISS solver.

getCISSRefinement()

Get the values of various refinement parameters in the CISS solver.

getCISSSizes()

Get the values of various size parameters in the CISS solver.

getCISSThreshold()

Get the values of various threshold parameters in the CISS solver.

getConverged()

Get the number of converged eigenpairs.

getConvergedReason()

Get the reason why the solve() iteration was stopped.

getConvergenceTest()

Get the method used to compute the error estimate used in the convergence test.

getDS()

Get the direct solver associated to the eigensolver.

getDimensions()

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().

getErrorEstimate(i)

Get the error estimate associated to the i-th computed eigenpair.

getExtract()

Get the extraction technique used by the PEP object.

getInterval()

Get the computational interval for spectrum slicing.

getIterationNumber()

Get the current iteration number.

getJDFix()

Get threshold for changing the target in the correction equation.

getJDMinimalityIndex()

Get the maximum allowed value of the minimality index.

getJDProjection()

Get the type of projection to be used in the Jacobi-Davidson solver.

getJDRestart()

Get the restart parameter used in the Jacobi-Davidson method.

getJDReusePreconditioner()

Get the flag for reusing the preconditioner.

getLinearEPS()

Get the eigensolver object associated to the polynomial eigenvalue solver.

getLinearExplicitMatrix()

Get if the matrices A and B for the linearization are built explicitly.

getLinearLinearization()

Get the coeffs.

getMonitor()

Get the list of monitor functions.

getOperators()

Get the matrices associated with the eigenvalue problem.

getOptionsPrefix()

Get the prefix used for searching for all PEP options in the database.

getProblemType()

Get the problem type from the PEP object.

getQArnoldiLocking()

Get the locking flag used in the Q-Arnoldi method.

getQArnoldiRestart()

Get the restart parameter used in the Q-Arnoldi method.

getRG()

Get the region object associated to the eigensolver.

getRefine()

Get the refinement strategy used by the PEP object.

getRefineKSP()

Get the KSP object used by the eigensolver in the refinement phase.

getST()

Get the ST object associated to the eigensolver object.

getSTOARCheckEigenvalueType()

Get the flag for the eigenvalue type check in spectrum slicing.

getSTOARDetectZeros()

Get the flag that enforces zero detection in spectrum slicing.

getSTOARDimensions()

Get the dimensions used for each subsolve step.

getSTOARInertias()

Get the values of the shifts and their corresponding inertias.

getSTOARLinearization()

Get the coefficients that define the linearization of a quadratic eigenproblem.

getSTOARLocking()

Get the locking flag used in the STOAR method.

getScale([Dl, Dr])

Get the strategy used for scaling the polynomial eigenproblem.

getStoppingTest()

Get the stopping function.

getTOARLocking()

Get the locking flag used in the TOAR method.

getTOARRestart()

Get the restart parameter used in the TOAR method.

getTarget()

Get the value of the target.

getTolerances()

Get the tolerance and maximum iteration count.

getTrackAll()

Get the flag indicating whether all residual norms must be computed.

getType()

Get the PEP type of this object.

getWhichEigenpairs()

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.

setFromOptions()

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.

setSTOARCheckEigenvalueType(flag)

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

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.

Methods Documentation

appendOptionsPrefix(prefix=None)#

Append to the prefix used for searching for all PEP options in the database.

Logically collective.

Parameters:

prefix (str | None) – The prefix string to prepend to all PEP option requests.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:353

cancelMonitor()#

Clear all monitors for a PEP object.

Logically collective.

Source code at slepc4py/SLEPc/PEP.pyx:1162

Return type:

None

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:
  • i (int) – Index of the solution to be considered.

  • etype (ErrorType | None) – The error type to compute.

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:

float

Notes

The index i should be a value between 0 and nconv-1 (see getConverged()).

Source code at slepc4py/SLEPc/PEP.pyx:1290

create(comm=None)#

Create the PEP object.

Collective.

Parameters:

comm (Comm | None) – MPI communicator. If not provided, it defaults to all processes.

Return type:

Self

Source code at slepc4py/SLEPc/PEP.pyx:276

destroy()#

Destroy the PEP object.

Collective.

Source code at slepc4py/SLEPc/PEP.pyx:258

Return type:

Self

errorView(etype=None, viewer=None)#

Display the errors associated with the computed solution.

Collective.

Display the errors and the eigenvalues.

Parameters:
Return type:

None

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.

Source code at slepc4py/SLEPc/PEP.pyx:1324

getBV()#

Get the basis vectors object associated to the eigensolver.

Not collective.

Returns:

The basis vectors context.

Return type:

BV

Source code at slepc4py/SLEPc/PEP.pyx:954

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:

Basis

Source code at slepc4py/SLEPc/PEP.pyx:379

getCISSExtraction()#

Get the extraction technique used in the CISS solver.

Not collective.

Returns:

The extraction technique.

Return type:

CISSExtraction

Source code at slepc4py/SLEPc/PEP.pyx:2049

getCISSKSPs()#

Get the array of linear solver objects associated with the CISS solver.

Collective.

Returns:

The linear solver objects.

Return type:

list of petsc4py.PETSc.KSP

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.

Source code at slepc4py/SLEPc/PEP.pyx:2218

getCISSRefinement()#

Get the values of various refinement parameters in the CISS solver.

Not collective.

Returns:

  • inner (int) – Number of iterative refinement iterations (inner loop).

  • blsize (int) – Number of iterative refinement iterations (blocksize loop).

Return type:

tuple[int, int]

Source code at slepc4py/SLEPc/PEP.pyx:2200

getCISSSizes()#

Get the values of various size parameters in the CISS solver.

Not collective.

Returns:

  • ip (int) – Number of integration points.

  • bs (int) – Block size.

  • ms (int) – Moment size.

  • npart (int) – Number of partitions when splitting the communicator.

  • bsmax (int) – Maximum block size.

  • realmats (bool) – True if A and B are real.

Return type:

tuple[int, int, int, int, int, bool]

Source code at slepc4py/SLEPc/PEP.pyx:2114

getCISSThreshold()#

Get the values of various threshold parameters in the CISS solver.

Not collective.

Returns:

  • delta (float) – Threshold for numerical rank.

  • spur (float) – Spurious threshold (to discard spurious eigenpairs.

Return type:

tuple[float, float]

Source code at slepc4py/SLEPc/PEP.pyx:2163

getConverged()#

Get the number of converged eigenpairs.

Not collective.

Returns:

Number of converged eigenpairs.

Return type:

int

Source code at slepc4py/SLEPc/PEP.pyx:1226

getConvergedReason()#

Get the reason why the solve() iteration was stopped.

Not collective.

Returns:

Negative value indicates diverged, positive value converged.

Return type:

ConvergedReason

Source code at slepc4py/SLEPc/PEP.pyx:1210

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:

Conv

Source code at slepc4py/SLEPc/PEP.pyx:601

getDS()#

Get the direct solver associated to the eigensolver.

Not collective.

Returns:

The direct solver context.

Return type:

DS

Source code at slepc4py/SLEPc/PEP.pyx:1012

getDimensions()#

Get the number of eigenvalues to compute and the dimension of the subspace.

Not collective.

Returns:

  • nev (int) – Number of eigenvalues to compute.

  • ncv (int) – Maximum dimension of the subspace to be used by the solver.

  • mpd (int) – Maximum dimension allowed for the projected problem.

Return type:

tuple[int, int, int]

Source code at slepc4py/SLEPc/PEP.pyx:779

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:
  • i (int) – Index of the solution to be obtained.

  • Vr (Vec | None) – Placeholder for the returned eigenvector (real part).

  • Vi (Vec | None) – Placeholder for the returned eigenvector (imaginary part).

Returns:

The computed eigenvalue.

Return type:

complex

Source code at slepc4py/SLEPc/PEP.pyx:1241

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:

float

Source code at slepc4py/SLEPc/PEP.pyx:1270

getExtract()#

Get the extraction technique used by the PEP object.

Not collective.

Returns:

The extraction strategy.

Return type:

Extract

Source code at slepc4py/SLEPc/PEP.pyx:732

getInterval()#

Get the computational interval for spectrum slicing.

Not collective.

Returns:

  • inta (float) – The left end of the interval.

  • intb (float) – The right end of the interval.

Return type:

tuple[float, float]

Notes

If the interval was not set by the user, then zeros are returned.

Source code at slepc4py/SLEPc/PEP.pyx:554

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:

int

Source code at slepc4py/SLEPc/PEP.pyx:1192

getJDFix()#

Get threshold for changing the target in the correction equation.

Not collective.

Returns:

The threshold for changing the target.

Return type:

float

Source code at slepc4py/SLEPc/PEP.pyx:1931

getJDMinimalityIndex()#

Get the maximum allowed value of the minimality index.

Not collective.

Returns:

The maximum minimality index.

Return type:

int

Source code at slepc4py/SLEPc/PEP.pyx:1989

getJDProjection()#

Get the type of projection to be used in the Jacobi-Davidson solver.

Not collective.

Returns:

The type of projection.

Return type:

JDProjection

Source code at slepc4py/SLEPc/PEP.pyx:2018

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:

float

Source code at slepc4py/SLEPc/PEP.pyx:1896

getJDReusePreconditioner()#

Get the flag for reusing the preconditioner.

Not collective.

Returns:

The reuse flag.

Return type:

bool

Source code at slepc4py/SLEPc/PEP.pyx:1960

getLinearEPS()#

Get the eigensolver object associated to the polynomial eigenvalue solver.

Collective.

Returns:

The linear eigensolver.

Return type:

EPS

Source code at slepc4py/SLEPc/PEP.pyx:1398

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:

bool

Source code at slepc4py/SLEPc/PEP.pyx:1472

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:

  • alpha (float) – First parameter of the linearization.

  • beta (float) – Second parameter of the linearization.

Return type:

tuple[float, float]

Source code at slepc4py/SLEPc/PEP.pyx:1434

getMonitor()#

Get the list of monitor functions.

Not collective.

Source code at slepc4py/SLEPc/PEP.pyx:1154

Return type:

PEPMonitorFunction

getOperators()#

Get the matrices associated with the eigenvalue problem.

Collective.

Returns:

The matrices associated with the eigensystem.

Return type:

list of petsc4py.PETSc.Mat

Source code at slepc4py/SLEPc/PEP.pyx:1041

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:

str

Source code at slepc4py/SLEPc/PEP.pyx:323

getProblemType()#

Get the problem type from the PEP object.

Not collective.

Returns:

The problem type that was previously set.

Return type:

ProblemType

Source code at slepc4py/SLEPc/PEP.pyx:414

getQArnoldiLocking()#

Get the locking flag used in the Q-Arnoldi method.

Not collective.

Returns:

The locking flag.

Return type:

bool

Source code at slepc4py/SLEPc/PEP.pyx:1550

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:

float

Source code at slepc4py/SLEPc/PEP.pyx:1514

getRG()#

Get the region object associated to the eigensolver.

Not collective.

Returns:

The region context.

Return type:

RG

Source code at slepc4py/SLEPc/PEP.pyx:983

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:

tuple[Refine, int, float, int, RefineScheme]

Source code at slepc4py/SLEPc/PEP.pyx:632

getRefineKSP()#

Get the KSP object used by the eigensolver in the refinement phase.

Collective.

Returns:

The linear solver object.

Return type:

petsc4py.PETSc.KSP

Source code at slepc4py/SLEPc/PEP.pyx:702

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:

ST

Source code at slepc4py/SLEPc/PEP.pyx:828

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:

bool

Source code at slepc4py/SLEPc/PEP.pyx:1857

getSTOARDetectZeros()#

Get the flag that enforces zero detection in spectrum slicing.

Not collective.

Returns:

The zero detection flag.

Return type:

bool

Source code at slepc4py/SLEPc/PEP.pyx:1739

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:

  • nev (int) – Number of eigenvalues to compute.

  • ncv (int) – Maximum dimension of the subspace to be used by the solver.

  • mpd (int) – Maximum dimension allowed for the projected problem.

Return type:

tuple[int, int, int]

Source code at slepc4py/SLEPc/PEP.pyx:1786

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:

  • shifts (ArrayReal) – The values of the shifts used internally in the solver.

  • inertias (ArrayInt) – The values of the inertia in each shift.

Return type:

tuple[ArrayReal, ArrayInt]

Source code at slepc4py/SLEPc/PEP.pyx:1810

getSTOARLinearization()#

Get the coefficients that define the linearization of a quadratic eigenproblem.

Not collective.

Returns:

  • alpha (float) – First parameter of the linearization.

  • beta (float) – Second parameter of the linearization.

Return type:

tuple[float, float]

Source code at slepc4py/SLEPc/PEP.pyx:1659

getSTOARLocking()#

Get the locking flag used in the STOAR method.

Not collective.

Returns:

The locking flag.

Return type:

bool

Source code at slepc4py/SLEPc/PEP.pyx:1698

getScale(Dl=None, Dr=None)#

Get the strategy used for scaling the polynomial eigenproblem.

Not collective.

Parameters:
Returns:

  • scale (Scale) – The scaling strategy.

  • alpha (float) – The scaling factor.

  • its (int) – The number of iteration of diagonal scaling.

  • lbda (float) – Approximation of the wanted eigenvalues (modulus).

Return type:

tuple[Scale, float, int, float]

Source code at slepc4py/SLEPc/PEP.pyx:860

getStoppingTest()#

Get the stopping function.

Not collective.

Source code at slepc4py/SLEPc/PEP.pyx:1123

Return type:

PEPStoppingFunction

getTOARLocking()#

Get the locking flag used in the TOAR method.

Not collective.

Returns:

The locking flag.

Return type:

bool

Source code at slepc4py/SLEPc/PEP.pyx:1625

getTOARRestart()#

Get the restart parameter used in the TOAR method.

Not collective.

Returns:

The number of vectors to be kept at restart.

Return type:

float

Source code at slepc4py/SLEPc/PEP.pyx:1589

getTarget()#

Get the value of the target.

Not collective.

Returns:

The value of the target.

Return type:

Scalar

Notes

If the target was not set by the user, then zero is returned.

Source code at slepc4py/SLEPc/PEP.pyx:472

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:

  • tol (float) – The convergence tolerance.

  • max_it (int) – The maximum number of iterations

Return type:

tuple[float, int]

Source code at slepc4py/SLEPc/PEP.pyx:511

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:

bool

Source code at slepc4py/SLEPc/PEP.pyx:747

getType()#

Get the PEP type of this object.

Not collective.

Returns:

The solver currently being used.

Return type:

str

Source code at slepc4py/SLEPc/PEP.pyx:308

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:

Which

Source code at slepc4py/SLEPc/PEP.pyx:443

reset()#

Reset the PEP object.

Collective.

Source code at slepc4py/SLEPc/PEP.pyx:268

Return type:

None

setBV(bv)#

Set a basis vectors object associated to the eigensolver.

Collective.

Parameters:

bv (BV) – The basis vectors context.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:970

setBasis(basis)#

Set the type of polynomial basis used.

Logically collective.

Set the type of polynomial basis used to describe the polynomial eigenvalue problem.

Parameters:

basis (Basis) – The basis to be set.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:397

setCISSExtraction(extraction)#

Set the extraction technique used in the CISS solver.

Logically collective.

Parameters:

extraction (CISSExtraction) – The extraction technique.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:2035

setCISSRefinement(inner=None, blsize=None)#

Set the values of various refinement parameters in the CISS solver.

Logically collective.

Parameters:
  • inner (int | None) – Number of iterative refinement iterations (inner loop).

  • blsize (int | None) – Number of iterative refinement iterations (blocksize loop).

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:2181

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:
  • ip (int | None) – Number of integration points.

  • bs (int | None) – Block size.

  • ms (int | None) – Moment size.

  • npart (int | None) – Number of partitions when splitting the communicator.

  • bsmax (int | None) – Maximum block size.

  • realmats (bool) – True if A and B are real.

Return type:

None

Notes

The default number of partitions is 1. This means the internal petsc4py.PETSc.KSP object is shared among all processes of the PEP communicator. Otherwise, the communicator is split into npart communicators, so that npart petsc4py.PETSc.KSP solves proceed simultaneously.

Source code at slepc4py/SLEPc/PEP.pyx:2064

setCISSThreshold(delta=None, spur=None)#

Set the values of various threshold parameters in the CISS solver.

Logically collective.

Parameters:
  • delta (float | None) – Threshold for numerical rank.

  • spur (float | None) – Spurious threshold (to discard spurious eigenpairs).

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:2144

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:

None

Source code at slepc4py/SLEPc/PEP.pyx:617

setDS(ds)#

Set a direct solver object associated to the eigensolver.

Collective.

Parameters:

ds (DS) – The direct solver context.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1028

setDimensions(nev=None, ncv=None, mpd=None)#

Set the number of eigenvalues to compute and the dimension of the subspace.

Logically collective.

Parameters:
  • nev (int | None) – Number of eigenvalues to compute.

  • ncv (int | None) – Maximum dimension of the subspace to be used by the solver.

  • mpd (int | None) – Maximum dimension allowed for the projected problem.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:800

setExtract(extract)#

Set the extraction strategy to be used.

Logically collective.

Parameters:

extract (Extract) – The extraction strategy.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:718

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:

None

setInitialSpace(space)#

Set the initial space from which the eigensolver starts to iterate.

Collective.

Parameters:

space (Vec | list[Vec]) – The initial space

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1083

setInterval(inta, intb)#

Set the computational interval for spectrum slicing.

Logically collective.

Parameters:
  • inta (float) – The left end of the interval.

  • intb (float) – The right end of the interval.

Return type:

None

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, see setWhichEigenpairs().

Source code at slepc4py/SLEPc/PEP.pyx:576

setJDFix(fix)#

Set the threshold for changing the target in the correction equation.

Logically collective.

Parameters:

fix (float) – Threshold for changing the target.

Return type:

None

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.

Source code at slepc4py/SLEPc/PEP.pyx:1911

setJDMinimalityIndex(flag)#

Set the maximum allowed value for the minimality index.

Logically collective.

Parameters:

flag (int) – The maximum minimality index.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1975

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:

None

Source code at slepc4py/SLEPc/PEP.pyx:2004

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.

Parameters:

keep (float) – The number of vectors to be kept at restart.

Return type:

None

Notes

Allowed values are in the range [0.1,0.9]. The default is 0.5.

Source code at slepc4py/SLEPc/PEP.pyx:1874

setJDReusePreconditioner(flag)#

Set a flag indicating whether the preconditioner must be reused or not.

Logically collective.

Parameters:

flag (bool) – The reuse flag.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1946

setLinearEPS(eps)#

Set an eigensolver object associated to the polynomial eigenvalue solver.

Collective.

Parameters:

eps (EPS) – The linear eigensolver.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1385

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:

None

Source code at slepc4py/SLEPc/PEP.pyx:1455

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:
  • alpha (float) – First parameter of the linearization.

  • beta (float) – Second parameter of the linearization.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1414

setMonitor(monitor, args=None, kargs=None)#

Append a monitor function to the list of monitors.

Logically collective.

Source code at slepc4py/SLEPc/PEP.pyx:1133

Parameters:
Return type:

None

setOperators(operators)#

Set the matrices associated with the eigenvalue problem.

Collective.

Parameters:

operators (list[Mat]) – The matrices associated with the eigensystem.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1063

setOptionsPrefix(prefix=None)#

Set the prefix used for searching for all PEP options in the database.

Logically collective.

Parameters:

prefix (str | None) – The prefix string to prepend to all PEP option requests.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:338

setProblemType(problem_type)#

Set the type of the eigenvalue problem.

Logically collective.

Parameters:

problem_type (ProblemType) – The problem type to be set.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:429

setQArnoldiLocking(lock)#

Toggle between locking and non-locking variants of the Q-Arnoldi method.

Logically collective.

Parameters:

lock (bool) – True if the locking variant must be selected.

Return type:

None

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).

Source code at slepc4py/SLEPc/PEP.pyx:1529

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.

Parameters:

keep (float) – The number of vectors to be kept at restart.

Return type:

None

Notes

Allowed values are in the range [0.1,0.9]. The default is 0.5.

Source code at slepc4py/SLEPc/PEP.pyx:1492

setRG(rg)#

Set a region object associated to the eigensolver.

Collective.

Parameters:

rg (RG) – The region context.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:999

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:
  • ref (Refine) – The refinement type.

  • npart (int | None) – The number of partitions of the communicator.

  • tol (float | None) – The convergence tolerance.

  • its (int | None) – The maximum number of refinement iterations.

  • scheme (RefineScheme | None) – Scheme for linear system solves

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:662

setST(st)#

Set a spectral transformation object associated to the eigensolver.

Collective.

Parameters:

st (ST) – The spectral transformation.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:847

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.

Parameters:

flag (bool) – Whether the eigenvalue type is checked or not.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1840

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.

Parameters:

detect (bool) – True if zeros must checked for.

Return type:

None

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.

Source code at slepc4py/SLEPc/PEP.pyx:1713

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:
  • nev (int | None) – Number of eigenvalues to compute.

  • ncv (int | None) – Maximum dimension of the subspace to be used by the solver.

  • mpd (int | None) – Maximum dimension allowed for the projected problem.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1754

setSTOARLinearization(alpha=1.0, beta=0.0)#

Set the coefficients that define the linearization of a quadratic eigenproblem.

Logically collective.

Parameters:
  • alpha (float) – First parameter of the linearization.

  • beta (float) – Second parameter of the linearization.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1642

setSTOARLocking(lock)#

Toggle between locking and non-locking variants of the STOAR method.

Logically collective.

Parameters:

lock (bool) – True if the locking variant must be selected.

Return type:

None

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).

Source code at slepc4py/SLEPc/PEP.pyx:1677

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:
Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:909

setStoppingTest(stopping, args=None, kargs=None)#

Set a function to decide when to stop the outer iteration of the eigensolver.

Logically collective.

Source code at slepc4py/SLEPc/PEP.pyx:1103

Parameters:
Return type:

None

setTOARLocking(lock)#

Toggle between locking and non-locking variants of the TOAR method.

Logically collective.

Parameters:

lock (bool) – True if the locking variant must be selected.

Return type:

None

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).

Source code at slepc4py/SLEPc/PEP.pyx:1604

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.

Parameters:

keep (float) – The number of vectors to be kept at restart.

Return type:

None

Notes

Allowed values are in the range [0.1,0.9]. The default is 0.5.

Source code at slepc4py/SLEPc/PEP.pyx:1567

setTarget(target)#

Set the value of the target.

Logically collective.

Parameters:

target (Scalar) – The value of the target.

Return type:

None

Notes

The target is a scalar value used to determine the portion of the spectrum of interest. It is used in combination with setWhichEigenpairs().

Source code at slepc4py/SLEPc/PEP.pyx:491

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:
  • tol (float | None) – The convergence tolerance.

  • max_it (int | None) – The maximum number of iterations

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:532

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.

Parameters:

trackall (bool) – Whether compute all residuals or not.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:762

setType(pep_type)#

Set the particular solver to be used in the PEP object.

Logically collective.

Parameters:

pep_type (Type | str) – The solver to be used.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:293

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:

None

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:

None

Source code at slepc4py/SLEPc/PEP.pyx:458

solve()#

Solve the eigensystem.

Collective.

Source code at slepc4py/SLEPc/PEP.pyx:1184

Return type:

None

valuesView(viewer=None)#

Display the computed eigenvalues in a viewer.

Collective.

Parameters:

viewer (Viewer | None) – Visualization context; if not provided, the standard output is used.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1353

vectorsView(viewer=None)#

Output computed eigenvectors to a viewer.

Collective.

Parameters:

viewer (Viewer | None) – Visualization context; if not provided, the standard output is used.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:1368

view(viewer=None)#

Print the PEP data structure.

Collective.

Parameters:

viewer (Viewer | None) – Visualization context; if not provided, the standard output is used.

Return type:

None

Source code at slepc4py/SLEPc/PEP.pyx:243

Attributes Documentation

bv#

The basis vectors (BV) object associated.

Source code at slepc4py/SLEPc/PEP.pyx:2299

ds#

The direct solver (DS) object associated.

Source code at slepc4py/SLEPc/PEP.pyx:2313

extract#

The type of extraction technique to be employed.

Source code at slepc4py/SLEPc/PEP.pyx:2264

max_it#

The maximum iteration count.

Source code at slepc4py/SLEPc/PEP.pyx:2278

problem_type#

The type of the eigenvalue problem.

Source code at slepc4py/SLEPc/PEP.pyx:2243

rg#

The region (RG) object associated.

Source code at slepc4py/SLEPc/PEP.pyx:2306

st#

The spectral transformation (ST) object associated.

Source code at slepc4py/SLEPc/PEP.pyx:2292

target#

The value of the target.

Source code at slepc4py/SLEPc/PEP.pyx:2257

tol#

The tolerance.

Source code at slepc4py/SLEPc/PEP.pyx:2271

track_all#

Compute the residual norm of all approximate eigenpairs.

Source code at slepc4py/SLEPc/PEP.pyx:2285

which#

The portion of the spectrum to be sought.

Source code at slepc4py/SLEPc/PEP.pyx:2250