-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathreference.html
16 lines (16 loc) · 67.6 KB
/
reference.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Reference · StandardPacked.jl</title><meta name="title" content="Reference · StandardPacked.jl"/><meta property="og:title" content="Reference · StandardPacked.jl"/><meta property="twitter:title" content="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/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 class="is-active"><a class="tocitem" href="reference.html">Reference</a><ul class="internal"><li><a class="tocitem" href="#The-SPMatrix-type"><span>The SPMatrix type</span></a></li><li><a class="tocitem" href="#Creating-an-SPMatrix"><span>Creating an SPMatrix</span></a></li><li><a class="tocitem" href="#Formats-of-an-SPMatrix"><span>Formats of an SPMatrix</span></a></li><li><a class="tocitem" href="#Accessing-the-data"><span>Accessing the data</span></a></li><li><a class="tocitem" href="#Diagonal-and-off-diagonal-elements"><span>Diagonal and off-diagonal elements</span></a></li><li><a class="tocitem" href="#Extensions-in-LinearAlgebra"><span>Extensions in LinearAlgebra</span></a></li></ul></li><li><a class="tocitem" href="lapack.html">BLAS and LAPACK reference</a></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="reference.html">Reference</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href="reference.html">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/reference.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="Reference"><a class="docs-heading-anchor" href="#Reference">Reference</a><a id="Reference-1"></a><a class="docs-heading-anchor-permalink" href="#Reference" title="Permalink"></a></h1><h2 id="The-SPMatrix-type"><a class="docs-heading-anchor" href="#The-SPMatrix-type">The SPMatrix type</a><a id="The-SPMatrix-type-1"></a><a class="docs-heading-anchor-permalink" href="#The-SPMatrix-type" 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.SPMatrix" href="#StandardPacked.SPMatrix"><code>StandardPacked.SPMatrix</code></a> — <span class="docstring-category">Type</span></header><section><div><pre><code class="language-julia hljs">SPMatrix{R,V,Fmt}</code></pre><p>Two-dimensional dense Hermitian square matrix in packed storage. Depending on Fmt, we only store the upper or lower triangle in col-major vectorized or scaled vectorized form:</p><ul><li><code>:U</code>: <code>PM[i, j] = PM[i + (j -1) * j ÷ 2]</code> for <code>1 ≤ i ≤ j ≤ n</code></li><li><code>:L</code>: <code>PM[i, j] = PM[i + (j -1) * (2n - j) ÷ 2]</code> for <code>1 ≤ j ≤ i ≤ n</code></li><li><code>:US</code>: <code>PM[i, j] = s(i, j) PM[i + (j -1) * j ÷ 2]</code> for <code>1 ≤ i ≤ j ≤ n</code></li><li><code>:LS</code>: <code>PM[i, j] = s(i, j) PM[i + (j -1) * (2n - j) ÷ 2]</code> for <code>1 ≤ j ≤ i ≤ n</code></li></ul><p>where <code>n</code> is the side dimension of the packed matrix <code>PM</code>, <code>s(i, i) = 1</code>, and <code>s(i, j) = inv(sqrt(2))</code> for <code>i ≠ j</code>.</p><div class="admonition is-warning"><header class="admonition-header">Warning</header><div class="admonition-body"><p>The SPMatrix can be effciently broadcast to other matrices, which works on the vector representation. It can also receive values from generic matrices by <code>copyto!</code> or broadcasting. Using it in a mixed chain of broadcastings of different type is not implemented and will potentially lead to logical errors (in the sense that the types will not match) or even segfaults (as the correct index mapping is not implemented).</p></div></div></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L4-L20">source</a></section></article><h2 id="Creating-an-SPMatrix"><a class="docs-heading-anchor" href="#Creating-an-SPMatrix">Creating an SPMatrix</a><a id="Creating-an-SPMatrix-1"></a><a class="docs-heading-anchor-permalink" href="#Creating-an-SPMatrix" 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.SPMatrix-Union{Tuple{R}, Tuple{Integer, AbstractVector{R}, Symbol}} where R" href="#StandardPacked.SPMatrix-Union{Tuple{R}, Tuple{Integer, AbstractVector{R}, Symbol}} where R"><code>StandardPacked.SPMatrix</code></a> — <span class="docstring-category">Method</span></header><section><div><pre><code class="language-julia hljs">SPMatrix(dim::Integer, data::AbstractVector{R}, type::Symbol=:U)</code></pre><p>Creates a packed matrix from already existing vectorized data. <code>type</code> must be one of <code>:U</code>, <code>:L</code>, <code>:US</code>, or <code>:LS</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L25-L29">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.SPMatrix-Union{Tuple{R}, Tuple{UndefInitializer, Integer, Symbol}} where R" href="#StandardPacked.SPMatrix-Union{Tuple{R}, Tuple{UndefInitializer, Integer, Symbol}} where R"><code>StandardPacked.SPMatrix</code></a> — <span class="docstring-category">Method</span></header><section><div><pre><code class="language-julia hljs">SPMatrix{R}(undef, dim, type::Symbol=:U)</code></pre><p>Creates a packed matrix with undefined data. <code>type</code> must be one of <code>:U</code>, <code>:L</code>, <code>:US</code>, or <code>:LS</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L36-L40">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.SPMatrix-Union{Tuple{LinearAlgebra.Symmetric{R, <:AbstractMatrix{R}}}, Tuple{R}} where R" href="#StandardPacked.SPMatrix-Union{Tuple{LinearAlgebra.Symmetric{R, <:AbstractMatrix{R}}}, Tuple{R}} where R"><code>StandardPacked.SPMatrix</code></a> — <span class="docstring-category">Method</span></header><section><div><pre><code class="language-julia hljs">SPMatrix(A::Union{<:Hermitian,<:Symmetric})</code></pre><p>Construct a packed matrix from a Hermitian or symmetric wrapper of any other matrix. Note that in case a symmetric wrapper is used, the element type must be invariant under conjugation (but this is not checked)!</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L353-L358">source</a></section></article><h2 id="Formats-of-an-SPMatrix"><a class="docs-heading-anchor" href="#Formats-of-an-SPMatrix">Formats of an SPMatrix</a><a id="Formats-of-an-SPMatrix-1"></a><a class="docs-heading-anchor-permalink" href="#Formats-of-an-SPMatrix" 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.SPMatrixUpper" href="#StandardPacked.SPMatrixUpper"><code>StandardPacked.SPMatrixUpper</code></a> — <span class="docstring-category">Type</span></header><section><div><pre><code class="language-julia hljs">SPMatrixUpper{R,V}</code></pre><p>This is a union type for scaled and unscaled packed matrices with upper triangular storage.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L66-L70">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.SPMatrixLower" href="#StandardPacked.SPMatrixLower"><code>StandardPacked.SPMatrixLower</code></a> — <span class="docstring-category">Type</span></header><section><div><pre><code class="language-julia hljs">SPMatrixLower{R,V}</code></pre><p>This is a union type for scaled and unscaled packed matrices with lower triangular storage.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L72-L76">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.SPMatrixUnscaled" href="#StandardPacked.SPMatrixUnscaled"><code>StandardPacked.SPMatrixUnscaled</code></a> — <span class="docstring-category">Type</span></header><section><div><pre><code class="language-julia hljs">SPMatrixUnscaled{R,V}</code></pre><p>This is a union type for unscaled packed matrices with upper or lower triangular storage.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L78-L82">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.SPMatrixScaled" href="#StandardPacked.SPMatrixScaled"><code>StandardPacked.SPMatrixScaled</code></a> — <span class="docstring-category">Type</span></header><section><div><pre><code class="language-julia hljs">SPMatrixScaled{R,V}</code></pre><p>This is a union type for scaled packed matrices with upper or lower triangular storage.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L84-L88">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.packed_isupper" href="#StandardPacked.packed_isupper"><code>StandardPacked.packed_isupper</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">packed_isupper(::SPMatrix)
packed_isupper(::Type{<:SPMatrix})</code></pre><p>Returns <code>true</code> iff the given packed matrix has upper triangular storage.</p><p>See also <a href="reference.html#StandardPacked.packed_islower"><code>packed_islower</code></a>, <a href="reference.html#StandardPacked.packed_isscaled"><code>packed_isscaled</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L92-L99">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.packed_islower" href="#StandardPacked.packed_islower"><code>StandardPacked.packed_islower</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">packed_islower(::SPMatrix)</code></pre><p>Returns <code>true</code> iff the given packed matrix has lower triangular storage.</p><p>See also <a href="reference.html#StandardPacked.packed_isupper"><code>packed_isupper</code></a>, <a href="reference.html#StandardPacked.packed_isscaled"><code>packed_isscaled</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L103-L109">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.packed_isscaled" href="#StandardPacked.packed_isscaled"><code>StandardPacked.packed_isscaled</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">packed_isscaled(::SPMatrix)</code></pre><p>Returns <code>true</code> iff the given packed matrix has scaled storage, i.e., if the off-diagonal elements are internally stored with a scaling of <span>$\sqrt2$</span>.</p><p>See also <a href="reference.html#StandardPacked.packed_isupper"><code>packed_isupper</code></a>, <a href="reference.html#StandardPacked.packed_islower"><code>packed_islower</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L113-L120">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.packed_scale!" href="#StandardPacked.packed_scale!"><code>StandardPacked.packed_scale!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">packed_scale!(P::SPMatrix)</code></pre><p>Ensures that <code>P</code> is a scaled packed matrix by multiplying off-diagonals by <span>$\sqrt2$</span> if necessary. Returns a scaled packed matrix. <code>P</code> itself should not be referenced any more, only the result of this function.</p><p>See also <a href="reference.html#StandardPacked.packed_unscale!"><code>packed_unscale!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L126-L133">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.packed_unscale!" href="#StandardPacked.packed_unscale!"><code>StandardPacked.packed_unscale!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">packed_unscale!(P::SPMatrix)</code></pre><p>Ensures that <code>P</code> is an unscaled packed matrix by dividing off-diagonals by <span>$\sqrt2$</span> if necessary. Returns an unscaled packed matrix. <code>P</code> itself should not be referenced any more, only the result of this function.</p><p>See also <a href="reference.html#StandardPacked.packed_scale!"><code>packed_scale!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L137-L144">source</a></section></article><h2 id="Accessing-the-data"><a class="docs-heading-anchor" href="#Accessing-the-data">Accessing the data</a><a id="Accessing-the-data-1"></a><a class="docs-heading-anchor-permalink" href="#Accessing-the-data" 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="Base.getindex" href="#Base.getindex"><code>Base.getindex</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">P[idx]</code></pre><p>Returns the value that is stored at the position <code>idx</code> in the vectorized (possibly scaled) representation of <code>P</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L149-L153">source</a></section><section><div><pre><code class="language-julia hljs">P[row, col]</code></pre><p>Returns the value that is stored in the given row and column of <code>P</code>. This corresponds to the value in the matrix, so even if <code>P</code> is of a scaled type, this does not affect the result.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L155-L160">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="Base.setindex!" href="#Base.setindex!"><code>Base.setindex!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">P[idx] = X</code></pre><p>Sets the value that is stored at the position <code>idx</code> in the vectorized (possibly scaled) representation of <code>P</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L183-L187">source</a></section><section><div><pre><code class="language-julia hljs">P[row, col] = X</code></pre><p>Sets the value that is stored in the given row and column of <code>P</code>. This corresponds to the value in the matrix, so if <code>P</code> is of a scaled type, <code>X</code> will internally be multiplied by <span>$\sqrt2$</span>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L189-L194">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="Base.vec-Tuple{SPMatrix}" href="#Base.vec-Tuple{SPMatrix}"><code>Base.vec</code></a> — <span class="docstring-category">Method</span></header><section><div><pre><code class="language-julia hljs">vec(P::SPMatrix)</code></pre><p>Returns the vectorized data associated with <code>P</code>. Note that this returns the actual vector, not a copy.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L317-L321">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="Base.Matrix-Union{Tuple{SPMatrixUnscaled{R}}, Tuple{R}} where R" href="#Base.Matrix-Union{Tuple{SPMatrixUnscaled{R}}, Tuple{R}} where R"><code>Base.Matrix</code></a> — <span class="docstring-category">Method</span></header><section><div><pre><code class="language-julia hljs">Matrix{R}(P::SPMatrix{R}, viewtype=R isa Complex ? Hermitian : Symmetric) where {R}</code></pre><p>Construct a dense matrix from a packed matrix. The parameter <code>symmetric</code> determines which kind of view is returned. Use <code>convert(Matrix{R}, SPMatrix{R})</code> to return the matrix itself.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/spmatrix.jl#L345-L350">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.packedsize" href="#StandardPacked.packedsize"><code>StandardPacked.packedsize</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">packedsize(dim::Integer)
packedsize(A::AbstractMatrix)</code></pre><p>Calculates the number of unique entries in a symmetric matrix <code>A</code> with side dimension <code>dim</code>, <span>$\frac{\mathit{dim}(\mathit{dim} +1)}{2}$</span>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/utils.jl#L3-L9">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.packedside" href="#StandardPacked.packedside"><code>StandardPacked.packedside</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">packedside(fullsize::Integer)
packedside(AP::AbstractVector)
packedside(P::SPMatrix)</code></pre><p>Calculates the side dimension of a vector <code>AP</code> of size <code>fullside</code> representing a symmetric packed matrix, <span>$\frac{\sqrt{8\mathit{fullsize} -1}}{2}$</span>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/utils.jl#L17-L24">source</a></section></article><h2 id="Diagonal-and-off-diagonal-elements"><a class="docs-heading-anchor" href="#Diagonal-and-off-diagonal-elements">Diagonal and off-diagonal elements</a><a id="Diagonal-and-off-diagonal-elements-1"></a><a class="docs-heading-anchor-permalink" href="#Diagonal-and-off-diagonal-elements" 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.SPDiagonalIterator" href="#StandardPacked.SPDiagonalIterator"><code>StandardPacked.SPDiagonalIterator</code></a> — <span class="docstring-category">Type</span></header><section><div><pre><code class="language-julia hljs">SPDiagonalIterator(P::SPMatrix, k=0)
SPDiagonalIterator(fmt::Symbol, dim, k=0)</code></pre><p>Creates an iterator that returns the linear indices for iterating through the <code>k</code>th diagonal of a packed matrix.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/utils.jl#L43-L48">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.rmul_diags!" href="#StandardPacked.rmul_diags!"><code>StandardPacked.rmul_diags!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">rmul_diags!(P::SPMatrix, factor)</code></pre><p>Right-multiplies all diagonal entries in <code>P</code> by <code>factor</code>. Returns <code>P</code>.</p><p>See also <a href="reference.html#StandardPacked.lmul_diags!"><code>lmul_diags!</code></a>, <a href="reference.html#StandardPacked.rmul_offdiags!"><code>rmul_offdiags!</code></a>, <a href="reference.html#StandardPacked.lmul_offdiags!"><code>lmul_offdiags!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/utils.jl#L100-L106">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.rmul_offdiags!" href="#StandardPacked.rmul_offdiags!"><code>StandardPacked.rmul_offdiags!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">rmul_offdiags!(P::SPMatrix, factor)</code></pre><p>Right-multiplies all off-diagonal entries in <code>P</code> by <code>factor</code>. Returns <code>P</code>.</p><p>See also <a href="reference.html#StandardPacked.rmul_diags!"><code>rmul_diags!</code></a>, <a href="reference.html#StandardPacked.lmul_diags!"><code>lmul_diags!</code></a>, <a href="reference.html#StandardPacked.lmul_offdiags!"><code>lmul_offdiags!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/utils.jl#L114-L120">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.lmul_diags!" href="#StandardPacked.lmul_diags!"><code>StandardPacked.lmul_diags!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">lmul_diags!(factor, P::SPMatrix)</code></pre><p>Left-multiplies all diagonal entries in <code>P</code> by <code>factor</code>. Returns <code>P</code>.</p><p>See also <a href="reference.html#StandardPacked.rmul_diags!"><code>rmul_diags!</code></a>, <a href="reference.html#StandardPacked.rmul_offdiags!"><code>rmul_offdiags!</code></a>, <a href="reference.html#StandardPacked.lmul_offdiags!"><code>lmul_offdiags!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/utils.jl#L129-L135">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.lmul_offdiags!" href="#StandardPacked.lmul_offdiags!"><code>StandardPacked.lmul_offdiags!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">lmul_offdiags!(factor, P::SPMatrix)</code></pre><p>Left-multiplies all diagonal entries in <code>P</code> by <code>factor</code>. Returns <code>P</code>.</p><p>See also <a href="reference.html#StandardPacked.rmul_diags!"><code>rmul_diags!</code></a>, <a href="reference.html#StandardPacked.rmul_offdiags!"><code>rmul_offdiags!</code></a>, <a href="reference.html#StandardPacked.lmul_diags!"><code>lmul_diags!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/utils.jl#L143-L149">source</a></section></article><h2 id="Extensions-in-LinearAlgebra"><a class="docs-heading-anchor" href="#Extensions-in-LinearAlgebra">Extensions in LinearAlgebra</a><a id="Extensions-in-LinearAlgebra-1"></a><a class="docs-heading-anchor-permalink" href="#Extensions-in-LinearAlgebra" title="Permalink"></a></h2><p>This section documents functions for which the interface provided by <code>LinearAlgebra</code> has been extended to cover the <code>SPMatrix</code> type. However, since <code>SPMatrix</code> is an <code>AbstractMatrix</code>, many more helper and convenience functions are actually available that will (hopefully) fall back to the efficient implementations documented here.</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.dot" href="#LinearAlgebra.dot"><code>LinearAlgebra.dot</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">dot(::SPMatrix, ::SPMatrix)</code></pre><p>Gives the scalar product between two packed matrices, which corresponds to the Frobenius/Hilbert-Schmidt scalar product for matrices. This is fastest for scaled matrices.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L5-L10">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.axpy!" href="#LinearAlgebra.axpy!"><code>LinearAlgebra.axpy!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">axpy!(α, x::Union{<:SPMatrix,<:AbstractVector},
y::Union{<:SPMatrix,<:AbstractVector})</code></pre><p>Overwrite <code>y</code> with <code>x * α + y</code> and return <code>y</code>, where both <code>x</code> and <code>y</code> are interpreted as their corresponding vectorizations.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L110-L115">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.axpby!" href="#LinearAlgebra.axpby!"><code>LinearAlgebra.axpby!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">axpby!(α, x::Union{<:SPMatrix,<:AbstractVector}, β,
y::Union{<:SPMatrix,<:AbstractVector})</code></pre><p>Overwrite <code>y</code> with <code>x * α + y * β</code> and return <code>y</code>, where both <code>x</code> and <code>y</code> are interpreted as their corresponding vectorizations.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L120-L126">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.factorize" href="#LinearAlgebra.factorize"><code>LinearAlgebra.factorize</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">factorize(A::SPMatrix)</code></pre><p>Compute a Bunch-Kaufman factorization of <code>A</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L131-L135">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.cholesky" href="#LinearAlgebra.cholesky"><code>LinearAlgebra.cholesky</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">cholesky(P::SPMatrix, NoPivot(); shift=0, check=true) -> Cholesky</code></pre><p>Compute the Cholesky factorization of a packed positive definite matrix <code>P</code> and return a <code>Cholesky</code> factorization.</p><p>The triangular Cholesky factor can be obtained from the factorization <code>F</code> via <code>F.L</code> and <code>F.U</code>, where <code>P ≈ F.U' * F.U ≈ F.L * F.L'</code>.</p><p>The following functions are available for Cholesky objects from packed matrices: <code>size</code>, <code>\</code>, <code>inv</code>, <code>det</code>, <code>logdet</code>, <code>isposdef</code>.</p><p>When <code>check = true</code>, an error is thrown if the decomposition fails. When <code>check = false</code>, responsibility for checking the decomposition's validity (via issuccess) lies with the user.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L143-L156">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.cholesky!" href="#LinearAlgebra.cholesky!"><code>LinearAlgebra.cholesky!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">cholesky!(P::SPMatrix, NoPivot(); shift=0, check=true) -> Cholesky</code></pre><p>The same as <a href="reference.html#LinearAlgebra.cholesky"><code>cholesky</code></a>, but saves space by overwriting the input <code>P</code>, instead of creating a copy.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L159-L163">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.bunchkaufman" href="#LinearAlgebra.bunchkaufman"><code>LinearAlgebra.bunchkaufman</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">bunchkaufman(P::SPMatrix, ipiv=missing; check = true) -> S::BunchKaufman</code></pre><p>Compute the Bunch-Kaufman factorization of a packed matrix <code>P</code> <code>P'*U*D*U'*P</code> or <code>P'*L*D*L'*P</code>, depending on which triangle is stored in <code>P</code>, and return a <code>BunchKaufman</code> object.</p><p>Iterating the decomposition produces the components <code>S.D</code>, <code>S.U</code> or <code>S.L</code> as appropriate given <code>S.uplo</code>, and <code>S.p</code>.</p><p>When <code>check = true</code>, an error is thrown if the decomposition fails. When <code>check = false</code>, responsibility for checking the decomposition's validity (via issuccess) lies with the user.</p><p>The following functions are available for BunchKaufman objects: <code>size</code>, <code>\</code>, <code>inv</code>, <code>getindex</code>.</p><p>The pivoting vector <code>ipiv</code> may be passed beforehand to save allocations.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L180-L194">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.bunchkaufman!" href="#LinearAlgebra.bunchkaufman!"><code>LinearAlgebra.bunchkaufman!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">bunchkaufman!(P::SPMatrix, ipiv=missing; check = true) -> S::BunchKaufman</code></pre><p>The same as <a href="reference.html#LinearAlgebra.bunchkaufman"><code>bunchkaufman</code></a>, but saves space by overwriting the input <code>P</code>, instead of creating a copy.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L198-L202">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.eigvals" href="#LinearAlgebra.eigvals"><code>LinearAlgebra.eigvals</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigvals(P::SPMatrix, args...)</code></pre><p>Return the eigenvalues of <code>P</code>.</p><p>This function calls <a href="lapack.html#StandardPacked.spevd!"><code>spevd!</code></a> and will forward all additional arguments. However, consider using <a href="reference.html#LinearAlgebra.eigvals!"><code>eigvals!</code></a> for a non-allocationg version.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L342-L349">source</a></section><section><div><pre><code class="language-julia hljs">eigvals(P::SPMatrix, vl::Real, vu::Real, args...)</code></pre><p>Return the eigenvalues of <code>P</code>. It is possible to calculate only a subset of the eigenvalues by specifying a pair <code>vl</code> and <code>vu</code> for the lower and upper boundaries of the eigenvalues.</p><p>This function calls <a href="lapack.html#StandardPacked.spevx!"><code>spevx!</code></a> and will forward all additional arguments. However, consider using <a href="reference.html#LinearAlgebra.eigvals!"><code>eigvals!</code></a> for a non-allocationg version.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L351-L359">source</a></section><section><div><pre><code class="language-julia hljs">eigvals(P::SPMatrix, range::UnitRange, args...)</code></pre><p>Return the eigenvalues of <code>P</code>, overwriting <code>P</code> in the progress. It is possible to calculate only a subset of the eigenvalues by specifying a <code>UnitRange</code> <code>range</code> covering indices of the sorted eigenvalues, e.g. the 2nd to 8th eigenvalues.</p><p>This function calls <a href="lapack.html#StandardPacked.spevx!"><code>spevx!</code></a> and will forward all additional arguments. However, consider using <a href="reference.html#LinearAlgebra.eigvals!"><code>eigvals!</code></a> for a non-allocationg version.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L361-L369">source</a></section><section><div><pre><code class="language-julia hljs">eigvals(AP::SPMatrix, BP::SPMatrix, args...)</code></pre><p>Compute the generalized eigenvalues of <code>AP</code> and <code>BP</code>.</p><p>See also <a href="reference.html#LinearAlgebra.eigen"><code>eigen</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L371-L377">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.eigvals!" href="#LinearAlgebra.eigvals!"><code>LinearAlgebra.eigvals!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigvals!(P::SPMatrix, args...)
eigvals!(P::SPMatrix, vl::Real, vu::Real, args...)
eigvals!(P::SPMatrix, range::UnitRange, args...)</code></pre><p>The same as <a href="reference.html#LinearAlgebra.eigvals"><code>eigvals</code></a>, but saves space by overwriting the input <code>P</code>, instead of creating a copy.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L388-L394">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.eigmax" href="#LinearAlgebra.eigmax"><code>LinearAlgebra.eigmax</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigmax(P::SPMatrix, args...)</code></pre><p>Return the largest eigenvalue of <code>P</code>.</p><p>See also <a href="reference.html#StandardPacked.eigmax!"><code>eigmax!</code></a>, <a href="reference.html#LinearAlgebra.eigvals"><code>eigvals</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L397-L403">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.eigmax!" href="#StandardPacked.eigmax!"><code>StandardPacked.eigmax!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigmax!(P::SPMatrix, args...)</code></pre><p>Return the largest eigenvalue of <code>P</code>, overwriting <code>P</code> in the progress.</p><p>See also <a href="reference.html#LinearAlgebra.eigvals!"><code>eigvals!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L405-L411">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.eigmin" href="#LinearAlgebra.eigmin"><code>LinearAlgebra.eigmin</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigmin(P::SPMatrix, args...)</code></pre><p>Return the smallest eigenvalue of <code>P</code></p><p>See also <a href="reference.html#StandardPacked.eigmin!"><code>eigmin!</code></a>, <a href="reference.html#LinearAlgebra.eigvals"><code>eigvals</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L413-L419">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.eigmin!" href="#StandardPacked.eigmin!"><code>StandardPacked.eigmin!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigmin!(P::SPMatrix, args...)</code></pre><p>Return the smallest eigenvalue of <code>P</code>, overwriting <code>P</code> in the progress.</p><p>See also <a href="reference.html#LinearAlgebra.eigvals!"><code>eigvals!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L421-L427">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.eigvecs" href="#LinearAlgebra.eigvecs"><code>LinearAlgebra.eigvecs</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigvecs(P::SPMatrix, args...)</code></pre><p>Return a matrix <code>M</code> whose columns are the eigenvectors of <code>P</code>. (The <code>k</code>th eigenvector can be obtained from the slice <code>M[:, k]</code>.)</p><p>See also <a href="reference.html#LinearAlgebra.eigen"><code>eigen</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L430-L437">source</a></section><section><div><pre><code class="language-julia hljs">eigvecs(A::SPMatrix, B::SPMatrix, args...)</code></pre><p>Return a matrix <code>M</code> whose columns are the generalized eigenvectors of <code>A</code> and <code>B</code>. (The <code>k</code>th eigenvector can be obtained from the slice <code>M[:, k]</code>.)</p><p>See also <a href="reference.html#LinearAlgebra.eigen"><code>eigen</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L439-L446">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.eigvecs!" href="#StandardPacked.eigvecs!"><code>StandardPacked.eigvecs!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigvecs!(P::SPMatrix, args...)
eigvecs(A::SPMatrix, B::SPMatrix, args...)</code></pre><p>The same as <a href="reference.html#LinearAlgebra.eigvecs"><code>eigvecs</code></a>, but saves space by overwriting the input <code>P</code> (or <code>A</code> and <code>B</code>), instead of creating a copy.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L448-L453">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.eigen" href="#LinearAlgebra.eigen"><code>LinearAlgebra.eigen</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigen(P::SPMatrix, args...) -> Eigen</code></pre><p>Compute the eigenvalue decomposition of <code>P</code>, returning an <code>Eigen</code> factorization object <code>F</code> which contains the eigenvalues in <code>F.values</code> and the eigenvectors in the columns of the matrix <code>F.vectors</code>. (The <code>k</code>th eigenvector can be obtained from the slice <code>F.vectors[:, k]</code>.)</p><p>This function calls <a href="lapack.html#StandardPacked.spevd!"><code>spevd!</code></a> and will forward all arguments. However, consider using <a href="reference.html#LinearAlgebra.eigen!"><code>eigen!</code></a> for a non-allocating version.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L456-L465">source</a></section><section><div><pre><code class="language-julia hljs">eigen(A::SPMatrix, B::SPMatrix, args...) -> GeneralizedEigen</code></pre><p>Compute the generalized eigenvalue decomposition of <code>A</code> and <code>B</code>, returning a <code>GeneralizedEigen</code> factorization object <code>F</code> which contains the generalized eigenvalues in <code>F.values</code> and the generalized eigenvectors in the columns of the matrix <code>F.vectors</code>. (The <code>k</code>th generalized eigenvector can be obtained from the slice <code>F.vectors[:, k]</code>.)</p><p>This function calls <a href="lapack.html#StandardPacked.spgvd!"><code>spgvd!</code></a> and will forward all arguments. However, consider using <a href="reference.html#LinearAlgebra.eigen!"><code>eigen!</code></a> for a non-allocating version.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L467-L476">source</a></section><section><div><pre><code class="language-julia hljs">eigen(P::SPMatrix, irange::UnitRange, args...)</code></pre><p>Compute the eigenvalue decomposition of <code>P</code>, returning an <code>Eigen</code> factorization object <code>F</code> which contains the eigenvalues in <code>F.values</code> and the eigenvectors in the columns of the matrix <code>F.vectors</code>. (The <code>k</code>th eigenvector can be obtained from the slice <code>F.vectors[:, k]</code>.)</p><p>The <code>UnitRange</code> <code>irange</code> specifies indices of the sorted eigenvalues to search for.</p><p>This function calls <a href="lapack.html#StandardPacked.spevx!"><code>spevx!</code></a> and will forward all arguments. However, consider using <a href="reference.html#LinearAlgebra.eigen!"><code>eigen!</code></a> for a non-allocating version.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L478-L489">source</a></section><section><div><pre><code class="language-julia hljs">eigen(P::SPMatrix, vl::Real, vu::Real, args...)</code></pre><p>Compute the eigenvalue decomposition of <code>P</code>, returning an <code>Eigen</code> factorization object <code>F</code> which contains the eigenvalues in <code>F.values</code> and the eigenvectors in the columns of the matrix <code>F.vectors</code>. (The <code>k</code>th eigenvector can be obtained from the slice <code>F.vectors[:, k]</code>.)</p><p><code>vl</code> is the lower bound of the window of eigenvalues to search for, and <code>vu</code> is the upper bound.</p><p>This function calls <a href="lapack.html#StandardPacked.spevx!"><code>spevx!</code></a> and will forward all arguments. However, consider using <a href="reference.html#LinearAlgebra.eigen!"><code>eigen!</code></a> for a non-allocating version.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L491-L502">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.eigen!" href="#LinearAlgebra.eigen!"><code>LinearAlgebra.eigen!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">eigen!(P::SPMatrix, args...)
eigen!(A::SPMatrix, B::SPMatrix, args...)</code></pre><p>Same as <a href="reference.html#LinearAlgebra.eigen"><code>eigen</code></a>, but saves space by overwriting the input <code>A</code> (and <code>B</code>) instead of creating a copy.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L505-L510">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.diagind" href="#LinearAlgebra.diagind"><code>LinearAlgebra.diagind</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">diagind(P::SPMatrix, k::Integer=0)</code></pre><p>A vector containing the indices of the <code>k</code>the diagonal of the packed matrix <code>P</code>.</p><p>See also: <a href="reference.html#LinearAlgebra.diag"><code>diag</code></a>, <a href="reference.html#StandardPacked.SPDiagonalIterator"><code>SPDiagonalIterator</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L526-L532">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.diag" href="#LinearAlgebra.diag"><code>LinearAlgebra.diag</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">diag(P::SPMatrix, k::Integer=0)</code></pre><p>The <code>k</code>th diagonal of a packed matrix, as a vector.</p><p>See also: <a href="reference.html#LinearAlgebra.diagind"><code>diagind</code></a>, <a href="reference.html#StandardPacked.SPDiagonalIterator"><code>SPDiagonalIterator</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L535-L541">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.norm" href="#LinearAlgebra.norm"><code>LinearAlgebra.norm</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">norm(P::SPMatrix, p::Real=2)</code></pre><p>For any packed matrix <code>P</code>, compute the <code>p</code>-norm as if <code>P</code> was interpreted as a full vector (containing the off-diagonals twice without scaling).</p><p>The <code>p</code>-norm is defined as</p><p class="math-container">\[ \lVert P\rVert_p = \Biggl( \sum_{i, j = 1}^n \lvert P_{i, j}\rvert^p \Biggr)^{1/p}\]</p><p>with <span>$P_{i, j}$</span> the entries of <span>$P$</span>, <span>$\lvert P_{i, j}\rvert$</span> the <code>norm</code> of <span>$P_{i, j}$</span>, and <span>$n$</span> its side dimension.</p><p><code>p</code> can assume any numeric value (even though not all values produce a mathematically valid vector norm). In particular, <code>norm(P, Inf)</code> returns the largest value in <code>abs.(P)</code>, whereas <code>norm(P, -Inf)</code> returns the smallest. If <code>p=2</code>, then this is equivalent to the Frobenius norm.</p><p>The second argument <code>p</code> is not necessarily a part of the interface for <code>norm</code>, i.e. a custom type may only implement <code>norm(P)</code> without second argument.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L609-L627">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.tr" href="#LinearAlgebra.tr"><code>LinearAlgebra.tr</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">tr(P::SPMatrix)</code></pre><p>Matrix trace. Sums the diagonal elements of <code>P</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L630-L634">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.isposdef" href="#LinearAlgebra.isposdef"><code>LinearAlgebra.isposdef</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">isposdef(P::SPMatrix, tol=0)</code></pre><p>Test whether a matrix is positive definite by trying to perform a Cholesky factorization of <code>P</code>.</p><p>See also <a href="reference.html#LinearAlgebra.isposdef!"><code>isposdef!</code></a>, <a href="reference.html#LinearAlgebra.cholesky"><code>cholesky</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L646-L652">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.isposdef!" href="#LinearAlgebra.isposdef!"><code>LinearAlgebra.isposdef!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">isposdef!(P::SPMatrix, tol=0)</code></pre><p>Test whether a matrix is positive definite by trying to perform a Cholesky factorization of <code>P</code>, overwriting <code>P</code> in the process. See also <a href="reference.html#LinearAlgebra.isposdef"><code>isposdef</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L654-L659">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="Base.transpose" href="#Base.transpose"><code>Base.transpose</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">transpose(P::SPMatrix{<:Real})</code></pre><p>Identity operation for real-valued packed matrices.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L666-L670">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="Base.adjoint" href="#Base.adjoint"><code>Base.adjoint</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">P'
adjoint(P::SPMatrix)</code></pre><p>Identity operation for packed matrices</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L672-L677">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.checksquare" href="#LinearAlgebra.checksquare"><code>LinearAlgebra.checksquare</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">checksquare(P::SPMatrix)</code></pre><p>Returns the side dimension of <code>P</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L680-L684">source</a></section></article><h3 id="Low-level-matrix-operations"><a class="docs-heading-anchor" href="#Low-level-matrix-operations">Low-level matrix operations</a><a id="Low-level-matrix-operations-1"></a><a class="docs-heading-anchor-permalink" href="#Low-level-matrix-operations" 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="LinearAlgebra.mul!" href="#LinearAlgebra.mul!"><code>LinearAlgebra.mul!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">mul!(::SPMatrix, B::AbstractMatrix, ::AbstractVector, α, β)
mul!(::SPMatrix, B::AbstractMatrix, ::SPMatrix, α, β)
mul!(::AbstractVector, B::AbstractMatrix, ::SPMatrix, α, β)</code></pre><p>Combined inplace matrix-vector multiply-add <span>$A B α + C β$</span>. The result is stored in <code>C</code> by overwriting it. Note that <code>C</code> must not be aliased with either <code>A</code> or <code>B</code>. This method will automatically interpret all <code>SPMatrix</code> arguments as their (scaled) vectorizations (so this is simply an alias to calling the standard <code>mul!</code> method on the <a href="reference.html#Base.vec-Tuple{SPMatrix}"><code>vec</code></a> of the packed matrices).</p><p>Note that <code>B</code> must not be a packed matrix.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L689-L700">source</a></section><section><div><pre><code class="language-julia hljs">mul!(C::AbstractVector, A::SPMatrixUnscaled, B::AbstractVector, α, β)</code></pre><p>Combined inplace matrix-vector multiply-add <code>A B α + C β</code>. The result is stored in <code>C</code> by overwriting it. Note that <code>C</code> must not be aliased with either <code>A</code> or <code>B</code>. This method requires an unscaled packed matrix and only works on native floating point types supported by BLAS.</p><p>See also <a href="lapack.html#LinearAlgebra.BLAS.spmv!"><code>spmv!</code></a>, <a href="lapack.html#LinearAlgebra.BLAS.hpmv!"><code>hpmv!</code></a>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L712-L720">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.lmul!" href="#LinearAlgebra.lmul!"><code>LinearAlgebra.lmul!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">lmul!(a::Number, P::SPMatrix)</code></pre><p>Scale a packed matrix <code>P</code> by a scalar <code>a</code> overwriting <code>P</code> in-place. Use <a href="reference.html#LinearAlgebra.rmul!"><code>rmul!</code></a> to multiply scalar from right. The scaling operation respects the semantics of the multiplication <code>*</code> between <code>a</code> and an element of <code>P</code>. In particular, this also applies to multiplication involving non-finite numbers such as <code>NaN</code> and <code>±Inf</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L726-L732">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.rmul!" href="#LinearAlgebra.rmul!"><code>LinearAlgebra.rmul!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">rmul!(P::SPMatrix, b::Number)</code></pre><p>Scale a packed matrix <code>P</code> by a scalar <code>b</code> overwriting <code>P</code> in-place. Use <a href="reference.html#LinearAlgebra.lmul!"><code>lmul!</code></a> to multiply scalar from left. The scaling operation respects the semantics of the multiplication <code>*</code> between an element of <code>P</code> and <code>b</code>. In particular, this also applies to multiplication involving non-finite numbers such as <code>NaN</code> and <code>±Inf</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L735-L741">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.ldiv!" href="#LinearAlgebra.ldiv!"><code>LinearAlgebra.ldiv!</code></a> — <span class="docstring-category">Function</span></header><section><div><pre><code class="language-julia hljs">ldiv!(P, B)</code></pre><p>Compute <code>P \ B</code> in-place and overwriting <code>B</code> to store the result.</p><p>The argument <code>P</code> should not be a matrix. Rather, instead of matrices it should be a factorization object (e.g. produced by <a href="reference.html#LinearAlgebra.cholesky"><code>cholesky</code></a> or <a href="reference.html#LinearAlgebra.bunchkaufman"><code>bunchkaufman</code></a> from a <code>SPMatrix</code>). The reason for this is that factorization itself is both expensive and typically allocates memory (although it can also be done in-place via, e.g., <a href="reference.html#LinearAlgebra.cholesky!"><code>cholesky!</code></a>), and performance-critical situations requiring <code>ldiv!</code> usually also require fine-grained control over the factorization of <code>P</code>.</p></div><a class="docs-sourcelink" target="_blank" href="https://github.com/projekter/StandardPacked.jl/blob/7df8a3cc7e2908b8ac2d4597c67a706ecdf2bd6c/src/linalg.jl#L744-L753">source</a></section></article></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="index.html">« Introduction</a><a class="docs-footer-nextpage" href="lapack.html">BLAS and LAPACK 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="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option><option value="auto">Automatic (OS)</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> version 1.2.1 on <span class="colophon-date" title="Saturday 2 December 2023 17:21">Saturday 2 December 2023</span>. Using Julia version 1.9.4.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>