-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathlapack.html
64 lines (64 loc) · 59.4 KB
/
lapack.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>BLAS and LAPACK reference · StandardPacked.jl</title><meta name="title" content="BLAS and LAPACK reference · StandardPacked.jl"/><meta property="og:title" content="BLAS and LAPACK reference · StandardPacked.jl"/><meta property="twitter:title" content="BLAS and LAPACK reference · StandardPacked.jl"/><meta name="description" content="Documentation for StandardPacked.jl."/><meta property="og:description" content="Documentation for StandardPacked.jl."/><meta property="twitter:description" content="Documentation for StandardPacked.jl."/><script data-outdated-warner src="assets/warner.js"></script><link href="https://cdnjs.cloudflare.com/ajax/libs/lato-font/3.0.0/css/lato-font.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/juliamono/0.050/juliamono.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/6.4.2/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/6.4.2/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/6.4.2/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.16.8/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="assets/documenter.js"></script><script src="search_index.js"></script><script src="siteinfo.js"></script><script src="../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/catppuccin-mocha.css" data-theme-name="catppuccin-mocha"/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/catppuccin-macchiato.css" data-theme-name="catppuccin-macchiato"/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/catppuccin-frappe.css" data-theme-name="catppuccin-frappe"/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/catppuccin-latte.css" data-theme-name="catppuccin-latte"/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="assets/themeswap.js"></script></head><body><div id="documenter"><nav class="docs-sidebar"><div class="docs-package-name"><span class="docs-autofit"><a href="index.html">StandardPacked.jl</a></span></div><button class="docs-search-query input is-rounded is-small is-clickable my-2 mx-auto py-1 px-2" id="documenter-search-query">Search docs (Ctrl + /)</button><ul class="docs-menu"><li><a class="tocitem" href="index.html">Introduction</a></li><li><a class="tocitem" href="reference.html">Reference</a></li><li class="is-active"><a class="tocitem" href="lapack.html">BLAS and LAPACK reference</a><ul class="internal"><li><a class="tocitem" href="#BLAS-Level-2"><span>BLAS Level 2</span></a></li><li><a class="tocitem" href="#BLAS-Level-3"><span>BLAS Level 3</span></a></li><li><a class="tocitem" href="#LAPACK"><span>LAPACK</span></a></li></ul></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><a class="docs-sidebar-button docs-navbar-link fa-solid fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a><nav class="breadcrumb"><ul class="is-hidden-mobile"><li class="is-active"><a href="lapack.html">BLAS and LAPACK reference</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="lapack.html">BLAS and LAPACK reference</a></li></ul></nav><div class="docs-right"><a class="docs-navbar-link" href="https://github.com/projekter/StandardPacked.jl" title="View the repository on GitHub"><span class="docs-icon fa-brands"></span><span class="docs-label is-hidden-touch">GitHub</span></a><a class="docs-navbar-link" href="https://github.com/projekter/StandardPacked.jl/blob/main/docs/src/lapack.md" title="Edit source on GitHub"><span class="docs-icon fa-solid"></span></a><a class="docs-settings-button docs-navbar-link fa-solid fa-gear" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-article-toggle-button fa-solid fa-chevron-up" id="documenter-article-toggle-button" href="javascript:;" title="Collapse all docstrings"></a></div></header><article class="content" id="documenter-page"><h1 id="BLAS-and-LAPACK-reference"><a class="docs-heading-anchor" href="#BLAS-and-LAPACK-reference">BLAS and LAPACK reference</a><a id="BLAS-and-LAPACK-reference-1"></a><a class="docs-heading-anchor-permalink" href="#BLAS-and-LAPACK-reference" title="Permalink"></a></h1><p><code>StandardPacked.jl</code> makes several BLAS and LAPACK functions related to packed matrices available that are not provided by Base.</p><div class="admonition is-info"><header class="admonition-header">Coverage</header><div class="admonition-body"><p>The general rule of thumb is: Whenever a function is provided in Base for general matrices, this package makes its packed counterpart available. Note that there are a lot of additional BLAS and LAPACK functions with respect to packed matrices for which no interface is provided because the corresponding general function is not made available in Base.</p><p>If you want to have one of these functions included, just open an issue or pull request.</p></div></div><div class="admonition is-warning"><header class="admonition-header">Naming</header><div class="admonition-body"><p>It is sometimes rather mysterious why a function has the <code>s</code> (for symmetric) or <code>h</code> (for Hermitian) prefix: for example,</p><ul><li><a href="lapack.html#LinearAlgebra.BLAS.spr!"><code>spr!</code></a> is the symmetric rank-1 update, and since it is also implemented (though not made available in Base) for symmetric complex-valued matrices, there is a separate instance <a href="lapack.html#StandardPacked.hpr!"><code>hpr!</code></a> that does the Hermitian equivalent.</li><li><a href="lapack.html#StandardPacked.spev!"><code>spev!</code></a> is the symmetric <em>or</em> Hermitian eigensolver (depending on the type), and it does not change its name to <code>hpev!</code> even though LAPACK does. There is no symmetric complex eigensolver.</li><li><a href="lapack.html#StandardPacked.hptrd!"><code>hptrd!</code></a> is the Hermitian <em>or</em> symmetric tridiagonalized (depending on the type), and it does not change its name to <code>sptrd!</code> even though LAPACK does.</li></ul><p>This mystery comes from following Base's convention of names.</p></div></div><div class="admonition is-info"><header class="admonition-header">Variants</header><div class="admonition-body"><p>Note that for every function, this package provides an undocumented "raw" variant that does almost no checking apart from whether the <code>uplo</code> parameter is correct plus a check of the return value. This raw version is usually not made available directly by the Julia Base.</p><p>Additionally, <code>StandardPacked.jl</code> provide a convenience wrapper that works with the <code>AbstractArray</code> interface, automatically infers dimensions and checks for the compatibility of all the sizes - this is what is commonly provided by Julia. And finally, the vector representing the packed matrix can also be replaced by a <a href="reference.html#StandardPacked.SPMatrixUnscaled"><code>SPMatrixUnscaled</code></a> - then, the <code>uplo</code> argument is not present, as it is determined from the type.</p></div></div><div class="admonition is-success"><header class="admonition-header">Performance</header><div class="admonition-body"><p>Even the convenience wrappers give more control than you might be used to. Whenever some allocations might occur due to the need for temporaries or additional output arrays, these wrappers will allow to pass preallocated arrays if they are sized appropriately (where output arguments must match exactly in size, but temporaries may also be larger than necessary).</p><p>Whenever a parameter is documented with a <code>missing</code> default, this is such a parameter that might be preallocated if necessary. The parameters are all passed by position, not by keyword, which allows Julia to remove unnecessary checks by compiling the correct version.</p></div></div><h2 id="BLAS-Level-2"><a class="docs-heading-anchor" href="#BLAS-Level-2">BLAS Level 2</a><a id="BLAS-Level-2-1"></a><a class="docs-heading-anchor-permalink" href="#BLAS-Level-2" title="Permalink"></a></h2><p>Some of these functions are already present in <code>Base.BLAS</code> and are re-exported in <code>StandardPacked</code>. They are listed here for completeness only.</p><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="LinearAlgebra.BLAS.spmv!" href="#LinearAlgebra.BLAS.spmv!"><code>LinearAlgebra.BLAS.spmv!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">spmv!(α, AP::SPMatrixUnscaled, x, β)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L791-L793">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="LinearAlgebra.BLAS.hpmv!" href="#LinearAlgebra.BLAS.hpmv!"><code>LinearAlgebra.BLAS.hpmv!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">hpmv!(α, AP::SPMatrixUnscaled, x, β, y)</code></pre></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L802-L804">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="LinearAlgebra.BLAS.spr!" href="#LinearAlgebra.BLAS.spr!"><code>LinearAlgebra.BLAS.spr!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">spr!(α, x, AP::SPMatrix)</code></pre><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L810-L812">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.hpr!" href="#StandardPacked.hpr!"><code>StandardPacked.hpr!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">hpr!(uplo, α, x, AP::AbstractVector)
hpr!(α, x, AP::SPMatrix)</code></pre><p>Update matrix <span>$A$</span> as <span>$A + \alpha x x'$</span>, where <span>$A$</span> is a Hermitian matrix provided in packed format <code>AP</code> and <code>x</code> is a vector.</p><p>With <code>uplo = 'U'</code>, the vector <code>AP</code> must contain the upper triangular part of the Hermitian matrix packed sequentially, column by column, so that <code>AP[1]</code> contains <code>A[1, 1]</code>, <code>AP[2]</code> and <code>AP[3]</code> contain <code>A[1, 2]</code> and <code>A[2, 2]</code> respectively, and so on.</p><p>With <code>uplo = 'L'</code>, the vector <code>AP</code> must contain the lower triangular part of the symmetric matrix packed sequentially, column by column, so that <code>AP[1]</code> contains <code>A[1, 1]</code>, <code>AP[2]</code> and <code>AP[3]</code> contain <code>A[2, 1]</code> and <code>A[3, 1]</code> respectively, and so on.</p><p>The scalar input <code>α</code> must be real.</p><p>The array inputs <code>x</code> and <code>AP</code> must all be of <code>ComplexF32</code> or <code>ComplexF64</code> type. Return the updated <code>AP</code>.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L834-L851">source</a></section></article><h2 id="BLAS-Level-3"><a class="docs-heading-anchor" href="#BLAS-Level-3">BLAS Level 3</a><a id="BLAS-Level-3-1"></a><a class="docs-heading-anchor-permalink" href="#BLAS-Level-3" title="Permalink"></a></h2><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.gemmt!" href="#StandardPacked.gemmt!"><code>StandardPacked.gemmt!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">gemmt!(uplo, transA, transB, alpha, A, B, beta, C)</code></pre><p><code>gemmt!</code> computes a matrix-matrix product with general matrices but updates only the upper or lower triangular part of the result matrix.</p><div class="admonition is-info"><header class="admonition-header">Info</header><div class="admonition-body"><p>This function is a recent BLAS extension; for OpenBLAS, it requires at least version 0.3.22 (which is not yet shipped with Julia). If the currently available BLAS does not offer <code>gemmt</code>, the function falls back to <code>gemm</code>.</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L854-L863">source</a></section></article><h2 id="LAPACK"><a class="docs-heading-anchor" href="#LAPACK">LAPACK</a><a id="LAPACK-1"></a><a class="docs-heading-anchor-permalink" href="#LAPACK" title="Permalink"></a></h2><h3 id="Conversion"><a class="docs-heading-anchor" href="#Conversion">Conversion</a><a id="Conversion-1"></a><a class="docs-heading-anchor-permalink" href="#Conversion" title="Permalink"></a></h3><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.trttp!" href="#StandardPacked.trttp!"><code>StandardPacked.trttp!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">trttp!(uplo, A, AP::AbstractVector)
trttp!(A, AP::SPMatrixUnscaled)</code></pre><p><code>trttp!</code> copies a triangular matrix from the standard full format (TR) to the standard packed format (TP).</p><div class="admonition is-info"><header class="admonition-header">Info</header><div class="admonition-body"><p>This function is also implemented in plain Julia and therefore works with arbitrary element types.</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L894-L902">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.tpttr!" href="#StandardPacked.tpttr!"><code>StandardPacked.tpttr!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">tpttr!(uplo, AP::AbstractVector, A)
tpttr!(AP::SPMatrixUnscaled, A)</code></pre><p><code>tpttr!</code> copies a triangular matrix from the standard packed format (TP) to the standard full format (TR).</p><div class="admonition is-info"><header class="admonition-header">Info</header><div class="admonition-body"><p>This function is also implemented in plain Julia and therefore works with arbitrary element types.</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L933-L941">source</a></section></article><h3 id="Linear-systems-Cholesky"><a class="docs-heading-anchor" href="#Linear-systems-Cholesky">Linear systems - Cholesky</a><a id="Linear-systems-Cholesky-1"></a><a class="docs-heading-anchor-permalink" href="#Linear-systems-Cholesky" title="Permalink"></a></h3><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.pptrf!" href="#StandardPacked.pptrf!"><code>StandardPacked.pptrf!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">pptrf!(uplo, AP::AbstractVector) -> (AP, info)
pptrf!(AP::SPMatrix) -> (AP, info)</code></pre><p><code>pptrf!</code> computes the Cholesky factorization of a real symmetric or complex Hermitian positive definite matrix <span>$A$</span> stored in packed format <code>AP</code>.</p><p>The factorization has the form</p><ul><li><span>$A = U' U$</span>, if <code>uplo == 'U'</code>, or</li><li><span>$A = L L'$</span>, if <code>uplo == 'L'</code>,</li></ul><p>where <span>$U$</span> is an upper triangular matrix and <span>$L$</span> is lower triangular.</p><p>If <code>info > 0</code>, the leading minor of order <code>info</code> is not positive definite, and the factorization could not be completed.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L952-L966">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.pptrs!" href="#StandardPacked.pptrs!"><code>StandardPacked.pptrs!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">pptrs!(uplo, AP::AbstractVector, B)
pptrs!(AP::SPMatrixUnscaled, B)</code></pre><p><code>pptrs!</code> solves a system of linear equations <span>$A X = B$</span> with a symmetric or Hermitian positive definite matrix <span>$A$</span> in packed storage <code>AP</code> using the Cholesky factorization <span>$A = U' U$</span> or <span>$A = L L'$</span> computed by <a href="lapack.html#StandardPacked.pptrf!"><code>pptrf!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L979-L985">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.pptri!" href="#StandardPacked.pptri!"><code>StandardPacked.pptri!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">pptri!(uplo, AP::AbstractVector)
pptri!(AP::SPMatrixUnscaled)</code></pre><p><code>pptri!</code> computes the inverse of a real symmetric or complex Hermitian positive definite matrix <code>A</code> using the Cholesky factorization <span>$A = U' U$</span> or <span>$A = L L'$</span> computed by <a href="lapack.html#StandardPacked.pptrf!"><code>pptrf!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L995-L1001">source</a></section></article><h3 id="Linear-systems-symmetric-indefinite"><a class="docs-heading-anchor" href="#Linear-systems-symmetric-indefinite">Linear systems - symmetric indefinite</a><a id="Linear-systems-symmetric-indefinite-1"></a><a class="docs-heading-anchor-permalink" href="#Linear-systems-symmetric-indefinite" title="Permalink"></a></h3><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.spsv!" href="#StandardPacked.spsv!"><code>StandardPacked.spsv!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">spsv!(uplo, AP::AbstractVector, ipiv=missing, B) -> (B, AP, ipiv)
spsv!(AP::SPMatrix, ipiv=missing, B) -> (B, AP, ipiv)</code></pre><p><code>spsv!</code> computes the solution to a real or complex system of linear equations <span>$A X = B$</span>, where <span>$A$</span> is an <span>$N$</span>-by-<span>$N$</span> symmetric matrix stored in packed format <code>AP</code> and <span>$X$</span> and <code>B</code> are <span>$N$</span>-by-<span>$N_{\mathrm{rhs}}$</span> matrices.</p><p>The diagonal pivoting method is used to factor <span>$A$</span> as</p><ul><li><span>$A = U D U^\top$</span>, if <code>UPLO = 'U'</code>, or</li><li><span>$A = L D L^\top$</span>, if <code>UPLO = 'L'</code>,</li></ul><p>where <code>U</code> (or <code>L</code>) is a product of permutation and unit upper (lower) triangular matrices, <code>D</code> is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The factored form of <span>$A$</span> is then used to solve the system of equations <span>$A X = B$</span>.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1026-L1039">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.hpsv!" href="#StandardPacked.hpsv!"><code>StandardPacked.hpsv!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">hpsv!(uplo, AP::AbstractVector, ipiv=missing, B) -> (B, AP, ipiv)
hpsv!(AP::SPMatrix, ipiv=missing, B) -> (B, AP, ipiv)</code></pre><p><code>hpsv!</code> computes the solution to a complex system of linear equations <span>$A X = B$</span>, where <span>$A$</span> is an <span>$N$</span>-by-<span>$N$</span> Hermitian matrix stored in packed format <code>AP</code> and <span>$X$</span> and <code>B</code> are <span>$N$</span>-by-<span>$N_{\mathrm{rhs}}$</span> matrices.</p><p>The diagonal pivoting method is used to factor <span>$A$</span> as</p><ul><li><span>$A = U D U'$</span>, if <code>UPLO = 'U'</code>, or</li><li><span>$A = L D L'$</span>, if <code>UPLO = 'L'</code>,</li></ul><p>where <code>U</code> (or <code>L</code>) is a product of permutation and unit upper (lower) triangular matrices, <code>D</code> is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks. The factored form of <span>$A$</span> is then used to solve the system of equations <span>$A X = B$</span>.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1026-L1039">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.sptrf!" href="#StandardPacked.sptrf!"><code>StandardPacked.sptrf!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">sptrf!(uplo, AP::AbstractVector, ipiv=missing) -> (AP, ipiv, info)
sptrf!(AP::SPMatrix, ipiv=missing) -> (AP, ipiv, info)</code></pre><p><code>sptrf!</code> computes the factorization of a real or complex symmetric matrix <span>$A$</span> stored in packed format <code>AP</code> using the Bunch-Kaufman diagonal pivoting method:</p><p><span>$A = U D U^\top$</span> or <span>$A = L D L^\top$</span></p><p>where <span>$U$</span> (or <span>$L$</span>) is a product of permutation and unit upper (lower) triangular matrices, and <span>$D$</span> is symmetric and block diagonal with 1-by-1 and 2-by-2 diagonal blocks.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1060-L1072">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.hptrf!" href="#StandardPacked.hptrf!"><code>StandardPacked.hptrf!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">hptrf!(uplo, AP::AbstractVector, ipiv=missing) -> (AP, ipiv, info)
hptrf!(AP::SPMatrix, ipiv=missing) -> (AP, ipiv, info)</code></pre><p><code>hptrf!</code> computes the factorization of a complex Hermitian matrix <span>$A$</span> stored in packed format <code>AP</code> using the Bunch-Kaufman diagonal pivoting method:</p><p><span>$A = U D U'$</span> or <span>$A = L D L'$</span></p><p>where <span>$U$</span> (or <span>$L$</span>) is a product of permutation and unit upper (lower) triangular matrices, and <span>$D$</span> is Hermitian and block diagonal with 1-by-1 and 2-by-2 diagonal blocks.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1060-L1072">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.sptrs!" href="#StandardPacked.sptrs!"><code>StandardPacked.sptrs!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">sptrs!(uplo, AP::AbstractVector, ipiv, B)
sptrs!(AP::SPMatrixUnscaled, ipiv, B)</code></pre><p><code>sptrs!</code> solves a system of linear equations <span>$A X = B$</span> with a real or complex symmetric matrix <span>$A$</span> stored in packed format using the factorization <span>$A = U D U^\top$</span> or <span>$A = L D L^\top$</span> computed by <a href="lapack.html#StandardPacked.sptrf!"><code>sptrf!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1086-L1092">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.hptrs!" href="#StandardPacked.hptrs!"><code>StandardPacked.hptrs!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">hptrs!(uplo, AP::AbstractVector, ipiv, B)
hptrs!(AP::SPMatrixUnscaled, ipiv, B)</code></pre><p><code>hptrs!</code> solves a system of linear equations <span>$A X = B$</span> with a complex Hermitian matrix <span>$A$</span> stored in packed format using the factorization <span>$A = U D U'$</span> or <span>$A = L D L'$</span> computed by <a href="lapack.html#StandardPacked.hptrf!"><code>hptrf!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1086-L1092">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.sptri!" href="#StandardPacked.sptri!"><code>StandardPacked.sptri!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">sptri!(uplo, AP::AbstractVector, ipiv, work=missing)
sptri!(AP::SPMatrixUnscaled, ipiv, work=missing)</code></pre><p><code>sptri!</code> computes the inverse of a real or complex symmetric indefinite matrix <span>$A$</span> in packed storage <code>AP</code> using the factorization <span>$A = U D U^\top$</span> or <span>$A = L D L^\top$</span> computed by <a href="lapack.html#StandardPacked.sptrf!"><code>sptrf!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1114-L1120">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.hptri!" href="#StandardPacked.hptri!"><code>StandardPacked.hptri!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">hptri!(uplo, AP::AbstractVector, ipiv, work=missing)
hptri!(AP::SPMatrixUnscaled, ipiv, work=missing)</code></pre><p><code>hptri!</code> computes the inverse of a complex Hermitian indefinite matrix <span>$A$</span> in packed storage <code>AP</code> using the factorization <span>$A = U D U'$</span> or <span>$A = L D L'$</span> computed by <a href="lapack.html#StandardPacked.hptrf!"><code>hptrf!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1114-L1120">source</a></section></article><h3 id="Eigenvalue"><a class="docs-heading-anchor" href="#Eigenvalue">Eigenvalue</a><a id="Eigenvalue-1"></a><a class="docs-heading-anchor-permalink" href="#Eigenvalue" title="Permalink"></a></h3><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.spev!" href="#StandardPacked.spev!"><code>StandardPacked.spev!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">spev!('N', uplo, AP::AbstractVector, W=missing, work=missing[, rwork=missing]) -> W
spev!('N', AP::SPMatrix, W=missing, work=missing[, rwork=missing]) -> W
spev!('V', uplo, AP::AbstractVector, W=missing, Z=missing, work=missing
[, rwork=missing]) -> (W, Z)
spev!('V', AP::SPMatrix, W=missing, Z=missing, work=missing[, rwork=missing])
-> (W, Z)</code></pre><p>Finds the eigenvalues (first parameter <code>'N'</code>) or eigenvalues and eigenvectors (first parameter <code>'V'</code>) of a Hermitian matrix <span>$A$</span> in packed storage <code>AP</code>. If <code>uplo = 'U'</code>, <code>AP</code> is the upper triangle of <span>$A$</span>. If <code>uplo = 'L'</code>, <code>AP</code> is the lower triangle of <span>$A$</span>.</p><p>The output vector <code>W</code> and matrix <code>Z</code>, as well as the temporaries <code>work</code> and <code>rwork</code> (in the complex case) may be specified to save allocations.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1256-L1270">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.spevx!" href="#StandardPacked.spevx!"><code>StandardPacked.spevx!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">spevx!('N', range, uplo, AP::AbstractVector, vl, vu, il, iu, abstol, W=missing,
work=missing[, rwork=missing], iwork=missing) -> view(W)
spevx!('N', range, AP::SPMatrix, vl, vu, il, iu, abstol, W=missing, work=missing
[, rwork=missing], iwork=missing) -> view(W)
spevx!('V', range, uplo, AP::AbstractVector, vl, vu, il, iu, abstol, W=missing,
Z=missing, work=missing[, rwork=missing], iwork=missing, ifail=missing)
-> (view(W), view(Z), info, view(ifail))
spevx!('V', range, AP::SPMatrix, vl, vu, il, iu, abstol, W=missing, Z=missing,
work=missing[, rwork=missing], iwork=missing, ifail=missing)
-> (view(W), view(Z), info, view(ifail))</code></pre><p>Finds the eigenvalues (first parameter <code>'N'</code>) or eigenvalues and eigenvectors (first parameter <code>'V'</code>) of a Hermitian matrix <span>$A$</span> in packed storage <code>AP</code>. If <code>uplo = 'U'</code>, <code>AP</code> is the upper triangle of <span>$A$</span>. If <code>uplo = 'L'</code>, <code>AP</code> is the lower triangle of <span>$A$</span>. If <code>range = 'A'</code>, all the eigenvalues are found. If <code>range = 'V'</code>, the eigenvalues in the half-open interval <code>(vl, vu]</code> are found. If <code>range = 'I'</code>, the eigenvalues with indices between <code>il</code> and <code>iu</code> are found. <code>abstol</code> can be set as a tolerance for convergence.</p><p>The eigenvalues are returned in <code>W</code> and the eigenvectors in <code>Z</code>. <code>info</code> is zero upon success or the number of eigenvalues that failed to converge; their indices are stored in <code>ifail</code>. The resulting views to the original data are sized according to the number of eigenvalues that fulfilled the bounds given by the range.</p><p>The output vector <code>W</code>, <code>ifail</code> and matrix <code>Z</code>, as well as the temporaries <code>work</code>, <code>rwork</code> (in the complex case), and <code>iwork</code> may be specified to save allocations.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1515-L1539">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.spevd!" href="#StandardPacked.spevd!"><code>StandardPacked.spevd!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">spevd!('N', uplo, AP::AbstractVector, W=missing, work=missing[, rwork=missing],
iwork=missing) -> W
spevd!('N', AP::SPMatrix, W=missing, work=missing[, rwork=missing],
iwork=missing) -> W
spevd!('V', uplo, AP::AbstractVector, W=missing, Z=missing, work=missing
[, rwork=missing], iwork=missing) -> (W, Z)
spevd!('V', AP::SPMatrix, W=missing, Z=missing, work=missing[, rwork=missing],
iwork=missing) -> (W, Z)</code></pre><p>Finds the eigenvalues (first parameter <code>'N'</code>) or eigenvalues and eigenvectors (first parameter <code>'V'</code>) of a Hermitian matrix <span>$A$</span> in packed storage <code>AP</code>. If <code>uplo = 'U'</code>, <code>AP</code> is the upper triangle of <span>$A$</span>. If <code>uplo = 'L'</code>, <code>AP</code> is the lower triangle of <span>$A$</span>.</p><p>The output vector <code>W</code> and matrix <code>Z</code>, as well as the temporaries <code>work</code>, <code>rwork</code> (in the complex case), and <code>iwork</code> may be specified to save allocations.</p><p>If eigenvectors are desired, it uses a divide and conquer algorithm.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1706-L1724">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.spgv!" href="#StandardPacked.spgv!"><code>StandardPacked.spgv!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">spgv!(itype, 'N', uplo, AP::AbstractVector, BP::AbstractVector, W=missing,
work=missing[, rwork=missing]) -> (W, BP)
spgv!(itype, 'N', AP::SPMatrix, BP::SPMatrix, W=missing, work=missing
[, rwork=missing]) -> (W, BP)
spgv!(itype, 'V', uplo, AP::AbstractVector, BP::AbstractVector, W=missing, Z=missing,
work=missing[, rwork=missing]) -> (W, Z, BP)
spgv!(itype, 'V', AP::SPMatrix, BP::SPMatrix, W=missing, Z=missing,
work=missing[, rwork=missing]) -> (W, Z, BP)</code></pre><p>Finds the generalized eigenvalues (second parameter <code>'N'</code>) or eigenvalues and eigenvectors (second parameter <code>'V'</code>) of a Hermitian matrix <span>$A$</span> in packed storage <code>AP</code> and Hermitian positive-definite matrix <span>$B$</span> in packed storage <code>BP</code>. If <code>uplo = 'U'</code>, <code>AP</code> and <code>BP</code> are the upper triangles of <span>$A$</span> and <span>$B$</span>, respectively. If <code>uplo = 'L'</code>, they are the lower triangles.</p><ul><li>If <code>itype = 1</code>, the problem to solve is <span>$A x = \lambda B x$</span>.</li><li>If <code>itype = 2</code>, the problem to solve is <span>$A B x = \lambda x$</span>.</li><li>If <code>itype = 3</code>, the problem to solve is <span>$B A x = \lambda x$</span>.</li></ul><p>On exit, <span>$B$</span> is overwritten with the triangular factor <span>$U$</span> or <span>$L$</span> from the Cholesky factorization <span>$B = U' U$</span> or <span>$B = L L'$</span>.</p><p>The output vector <code>W</code> and matrix <code>Z</code>, as well as the temporaries <code>work</code> and <code>rwork</code> (in the complex case) may be specified to save allocations.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L1877-L1902">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.spgvx!" href="#StandardPacked.spgvx!"><code>StandardPacked.spgvx!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">spgvx!(itype, 'N', range, uplo, AP::AbstractVector, BP::AbstractVector, vl, vu,
il, iu, abstol, W=missing, work=missing[, rwork=missing], iwork=missing)
-> (view(W), BP)
spgvx!(itype, 'N', range, AP::SPMatrix, BP::SPMatrix, vl, vu, il, iu, abstol,
W=missing, work=missing[, rwork=missing], iwork=missing) -> (view(W), BP)
spgvx!(itype, 'V', range, uplo, AP::AbstractVector, BP::AbstractVector, vl, vu,
il, iu, abstol, W=missing, Z=missing, work=missing[, rwork=missing],
iwork=missing, ifail=missing) -> (view(W), view(Z), info, view(ifail), BP)
spgvx!(itype, 'V', range, AP::SPMatrix, BP::SPMatrix, vl, vu, il, iu, abstol,
W=missing, Z=missing, work=missing[, rwork=missing], iwork=missing, ifail=missing)
-> (view(W), view(Z), info, view(ifail), BP)</code></pre><p>Finds the generalized eigenvalues (second parameter <code>'N'</code>) or eigenvalues and eigenvectors (second parameter <code>'V'</code>) of a Hermitian matrix <span>$A$</span> in packed storage <code>AP</code> and Hermitian positive-definite matrix <span>$B$</span> in packed storage <code>BP</code>. If <code>uplo = 'U'</code>, <code>AP</code> and <code>BP</code> are the upper triangles of <span>$A$</span> and <span>$B$</span>, respectively. If <code>uplo = 'L'</code>, they are the lower triangles. If <code>range = 'A'</code>, all the eigenvalues are found. If <code>range = 'V'</code>, the eigenvalues in the half-open interval <code>(vl, vu]</code> are found. If <code>range = 'I'</code>, the eigenvalues with indices between <code>il</code> and <code>iu</code> are found. <code>abstol</code> can be set as a tolerance for convergence.</p><ul><li>If <code>itype = 1</code>, the problem to solve is <span>$A x = \lambda B x$</span>.</li><li>If <code>itype = 2</code>, the problem to solve is <span>$A B x = \lambda x$</span>.</li><li>If <code>itype = 3</code>, the problem to solve is <span>$B A x = \lambda x$</span>.</li></ul><p>On exit, <span>$B$</span> is overwritten with the triangular factor <span>$U$</span> or <span>$L$</span> from the Cholesky factorization <span>$B = U' U$</span> or <span>$B = L L'$</span>.</p><p>The eigenvalues are returned in <code>W</code> and the eigenvectors in <code>Z</code>. <code>info</code> is zero upon success or the number of eigenvalues that failed to converge; their indices are stored in <code>ifail</code>. The resulting views to the original data are sized according to the number of eigenvalues that fulfilled the bounds given by the range.</p><p>The output vector <code>W</code>, <code>ifail</code> and matrix <code>Z</code>, as well as the temporaries <code>work</code>, <code>rwork</code> (in the complex case), and <code>iwork</code> may be specified to save allocations.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L2146-L2180">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.spgvd!" href="#StandardPacked.spgvd!"><code>StandardPacked.spgvd!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">spgvd!(itype, 'N', uplo, AP::AbstractVector, BP::AbstractVector, W=missing,
work=missing[, rwork=missing], iwork=missing) -> (W, BP)
spgvd!(itype, 'N', AP::SPMatrix, BP::SPMatrix, W=missing, work=missing
[, rwork=missing], iwork=missing) -> (W, BP)
spgvd!(itype, 'V', uplo, AP::AbstractVector, BP::AbstractVector, W=missing, Z=missing,
work=missing[, rwork=missing], iwork=missing) -> (W, Z, BP)
spgvd!(itype, 'V', AP::SPMatrix, BP::SPMatrix, W=missing, Z=missing,
work=missing[, rwork=missing], iwork=missing) -> (W, Z, BP)</code></pre><p>Finds the generalized eigenvalues (second parameter <code>'N'</code>) or eigenvalues and eigenvectors (second parameter <code>'V'</code>) of a Hermitian matrix <span>$A$</span> in packed storage <code>AP</code> and Hermitian positive-definite matrix <span>$B$</span> in packed storage <code>BP</code>. If <code>uplo = 'U'</code>, <code>AP</code> and <code>BP</code> are the upper triangles of <span>$A$</span>, and <span>$B$</span>, respectively. If <code>uplo = 'L'</code>, they are the lower triangles.</p><ul><li>If <code>itype = 1</code>, the problem to solve is <span>$A x = \lambda B x$</span>.</li><li>If <code>itype = 2</code>, the problem to solve is <span>$A B x = \lambda x$</span>.</li><li>If <code>itype = 3</code>, the problem to solve is <span>$B A x = \lambda x$</span>.</li></ul><p>On exit, <span>$B$</span> is overwritten with the triangular factor <span>$U$</span> or <span>$L$</span> from the Cholesky factorization <span>$B = U' U$</span> or <span>$B = L L'$</span>.</p><p>The output vector <code>W</code> and matrix <code>Z</code>, as well as the temporaries <code>work</code>, <code>rwork</code> (in the complex case), and <code>iwork</code> may be specified to save allocations.</p><p>If eigenvectors are desired, it uses a divide and conquer algorithm.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L2366-L2393">source</a></section></article><h3 id="Tridiagonal"><a class="docs-heading-anchor" href="#Tridiagonal">Tridiagonal</a><a id="Tridiagonal-1"></a><a class="docs-heading-anchor-permalink" href="#Tridiagonal" title="Permalink"></a></h3><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.hptrd!" href="#StandardPacked.hptrd!"><code>StandardPacked.hptrd!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">hptrd!(uplo, AP::AbstractVector, D=missing, E=missing, τ=missing) -> (A, τ, D, E)
hptrd!(AP::SPMatrix, D=missing, E=missing, τ=missing) -> (A, τ, D, E)</code></pre><p><code>hptrd!</code> reduces a Hermitian matrix <span>$A$</span> stored in packed form <code>AP</code> to real-symmetric tridiagonal form <span>$T$</span> by an orthogonal similarity transformation: <span>$Q' A Q = T$</span>. If <code>uplo = 'U'</code>, <code>AP</code> is the upper triangle of <span>$A$</span>, else it is the lower triangle.</p><p>On exit, if <code>uplo = 'U'</code>, the diagonal and first superdiagonal of <span>$A$</span> are overwritten by the corresponding elements of the tridiagonal matrix <span>$T$</span>, and the elements above the first superdiagonal, with the array <code>τ</code>, represent the orthogonal matrix <span>$Q$</span> as a product of elementary reflectors. If <code>uplo = 'L'</code> the diagonal and first subdiagonal of <span>$A$</span> are over-written by the corresponding elements of the tridiagonal matrix <span>$T$</span>, and the elements below the first subdiagonal, with the array <code>τ</code>, represent the orthogonal matrix <span>$Q$</span> as a product of elementary reflectors.</p><p><code>tau</code> contains the scalar factors of the elementary reflectors, <code>D</code> contains the diagonal elements of <span>$T$</span> and <code>E</code> contains the off-diagonal elements of <span>$T$</span>.</p><p>The output vectors <code>D</code>, <code>E</code>, and <code>τ</code> may be specified to save allocations.</p><div class="admonition is-warning"><header class="admonition-header">Scaled matrices</header><div class="admonition-body"><p>The variant of this function that takes a <a href="reference.html#StandardPacked.SPMatrix"><code>SPMatrix</code></a> also allows for scaled packed matrices. It will automatically call <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a> on the matrix and return the unscaled result. Do not use the reference to the scaled matrix any more, only the result of this function!</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L2429-L2448">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.opgtr!" href="#StandardPacked.opgtr!"><code>StandardPacked.opgtr!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">opgtr!(uplo, AP::AbstractVector, τ, Q=missing, work=missing)
opgtr!(AP::SPMatrixUnscaled, τ, Q=missing, work=missing)</code></pre><p>Explicitly finds <code>Q</code>, the orthogonal/unitary matrix from <a href="lapack.html#StandardPacked.hptrd!"><code>hptrd!</code></a>. <code>uplo</code>, <code>AP</code>, and <code>τ</code> must correspond to the input/output to <code>hptrd!</code>.</p><p>The output matrix <code>Q</code> and the temporary <code>work</code> may be specified to save allocations.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L2474-L2482">source</a></section></article><article class="docstring"><header><a class="docstring-article-toggle-button fa-solid fa-chevron-down" href="javascript:;" title="Collapse docstring"></a><a class="docstring-binding" id="StandardPacked.opmtr!" href="#StandardPacked.opmtr!"><code>StandardPacked.opmtr!</code></a> — <span class="docstring-category">Function</span><span class="is-flex-grow-1 docstring-article-toggle-button" title="Collapse docstring"></span></header><section><div><pre><code class="language-julia hljs">opmtr!(side, uplo, trans, AP::AbstractVector, τ, C, work=missing)
opmtr!(side, trans, AP::SPMatrixUnscaled, τ, C, work=missing)</code></pre><p><code>opmtr!</code> overwrites the general <span>$m \times n$</span> matrix <code>C</code> with</p><table><tr><th style="text-align: right"></th><th style="text-align: center"><code>SIDE = 'L'</code></th><th style="text-align: center"><code>SIDE = 'R'</code></th></tr><tr><td style="text-align: right"><code>TRANS = 'N'</code></td><td style="text-align: center"><span>$Q C$</span></td><td style="text-align: center"><span>$C Q$</span></td></tr><tr><td style="text-align: right"><code>TRANS = 'T'</code> or <code>'C'</code> (complex case)</td><td style="text-align: center"><span>$Q' C$</span></td><td style="text-align: center"><span>$C Q'$</span></td></tr></table><p>where <span>$Q$</span> is a real orthogonal or complex unitary matrix of order <span>$n_q$</span>, with <span>$n_q = m$</span> if <code>SIDE = 'L'</code> and <span>$n_q = n$</span> if <code>SIDE = 'R'</code>. <span>$Q$</span> is defined as the product of <span>$n_q -1$</span> elementary reflectors, as returned by <a href="lapack.html#StandardPacked.hptrd!"><code>hptrd!</code></a> using packed storage:</p><ul><li>if <code>UPLO = 'U'</code>, <span>$Q = H_{n_q -1} \dotsm H_2 H_1$</span>;</li><li>if <code>UPLO = 'L'</code>, <span>$Q = H_1 H_2 \dotsm H_{n_q -1}$</span>.</li></ul><p>The temporary <code>work</code> may be specified to save allocations.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/f6495d88df16fab91e1c1715576a52bc34e00b07/src/lapack.jl#L2514-L2532">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="reference.html">« Reference</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="auto">Automatic (OS)</option><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option><option value="catppuccin-latte">catppuccin-latte</option><option value="catppuccin-frappe">catppuccin-frappe</option><option value="catppuccin-macchiato">catppuccin-macchiato</option><option value="catppuccin-mocha">catppuccin-mocha</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> version 1.8.0 on <span class="colophon-date" title="Friday 27 December 2024 16:03">Friday 27 December 2024</span>. Using Julia version 1.11.2.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>