singe_cusolver/
eigen.rs

1#[allow(unused_imports)]
2use crate::error::Status;
3
4use std::{mem::size_of, ptr};
5
6use singe_cuda::{
7    data_type::{DataType, DataTypeLike},
8    memory::DeviceMemory,
9    types::{Complex32, Complex64},
10};
11
12use crate::{
13    context::Context,
14    error::{Error, Result},
15    layout::{ByteWorkspaceMut, MatrixMut, MatrixRef, SelectionWorkspaceSizes, WorkspaceSizes},
16    params::Params,
17    sys, try_ffi,
18    types::{EigenMode, EigenRange, EigenType, FillMode},
19    utility::{to_i32, to_i64, to_usize},
20};
21
22#[derive(Debug, Clone, Copy, PartialEq)]
23pub enum EigenSelection<T> {
24    All,
25    ByValue { lower: T, upper: T },
26    ByIndex { start: usize, end: usize },
27}
28
29#[derive(Debug, Clone, Copy, PartialEq, Eq)]
30pub struct Syevd {
31    pub mode: EigenMode,
32    pub fill_mode: FillMode,
33    pub n: usize,
34    pub leading_dimension: usize,
35}
36
37impl Syevd {
38    pub fn new(mode: EigenMode, fill_mode: FillMode, n: usize, leading_dimension: usize) -> Self {
39        Self {
40            mode,
41            fill_mode,
42            n,
43            leading_dimension,
44        }
45    }
46
47    pub fn workspace_size<TA: DataTypeLike, TW: DataTypeLike>(
48        self,
49        ctx: &Context,
50        params: &Params,
51        input: SyevdInput<'_, TA, TW>,
52    ) -> Result<WorkspaceSizes> {
53        xsyevd_buffer_size(
54            ctx,
55            params,
56            self.mode,
57            self.fill_mode,
58            self.n,
59            input.a,
60            self.leading_dimension,
61            input.eigenvalues,
62        )
63    }
64
65    pub fn execute<TA: DataTypeLike, TW: DataTypeLike>(
66        self,
67        ctx: &Context,
68        params: &Params,
69        bindings: SyevdBindings<'_, TA, TW>,
70    ) -> Result<()> {
71        xsyevd(
72            ctx,
73            params,
74            self.mode,
75            self.fill_mode,
76            self.n,
77            bindings.a,
78            self.leading_dimension,
79            bindings.eigenvalues,
80            bindings.workspace,
81            bindings.dev_info,
82        )
83    }
84}
85
86#[derive(Debug, Clone, Copy)]
87pub struct SyevdInput<'a, TA, TW> {
88    pub a: &'a DeviceMemory<TA>,
89    pub eigenvalues: &'a DeviceMemory<TW>,
90}
91
92#[derive(Debug)]
93pub struct SyevdBindings<'a, TA, TW> {
94    pub a: &'a mut DeviceMemory<TA>,
95    pub eigenvalues: &'a mut DeviceMemory<TW>,
96    pub workspace: ByteWorkspaceMut<'a>,
97    pub dev_info: &'a mut DeviceMemory<i32>,
98}
99
100#[derive(Debug, Clone, Copy, PartialEq)]
101pub struct Syevdx<T> {
102    pub mode: EigenMode,
103    pub fill_mode: FillMode,
104    pub selection: EigenSelection<T>,
105    pub n: usize,
106    pub leading_dimension: usize,
107}
108
109impl<T> Syevdx<T> {
110    pub fn new(
111        mode: EigenMode,
112        fill_mode: FillMode,
113        selection: EigenSelection<T>,
114        n: usize,
115        leading_dimension: usize,
116    ) -> Self {
117        Self {
118            mode,
119            fill_mode,
120            selection,
121            n,
122            leading_dimension,
123        }
124    }
125}
126
127impl<T: Copy + Default> Syevdx<T> {
128    pub fn workspace_size<TA: DataTypeLike, TW: DataTypeLike>(
129        self,
130        ctx: &Context,
131        params: &Params,
132        input: SyevdInput<'_, TA, TW>,
133    ) -> Result<SelectionWorkspaceSizes> {
134        xsyevdx_buffer_size(
135            ctx,
136            params,
137            self.mode,
138            self.fill_mode,
139            self.selection,
140            self.n,
141            input.a,
142            self.leading_dimension,
143            input.eigenvalues,
144        )
145    }
146
147    pub fn execute<TA: DataTypeLike, TW: DataTypeLike>(
148        self,
149        ctx: &Context,
150        params: &Params,
151        bindings: SyevdBindings<'_, TA, TW>,
152    ) -> Result<usize> {
153        xsyevdx(
154            ctx,
155            params,
156            self.mode,
157            self.fill_mode,
158            self.selection,
159            self.n,
160            bindings.a,
161            self.leading_dimension,
162            bindings.eigenvalues,
163            bindings.workspace,
164            bindings.dev_info,
165        )
166    }
167}
168
169#[derive(Debug)]
170pub struct SyevjInfo {
171    handle: sys::syevjInfo_t,
172}
173
174// syevj info handles store solver options and expose mutation only through
175// &mut self, so immutable sharing is allowed.
176unsafe impl Send for SyevjInfo {}
177unsafe impl Sync for SyevjInfo {}
178
179impl SyevjInfo {
180    /// Creates `syevj`, `syevj_batched`, and `sygvj` parameter storage with default values.
181    ///
182    /// The returned [`SyevjInfo`] owns the cuSOLVER parameter handle and destroys
183    /// it when dropped.
184    ///
185    /// # Errors
186    ///
187    /// Returns an error if cuSOLVER cannot allocate the parameter storage or if
188    /// it does not return a valid handle.
189    pub fn create() -> Result<Self> {
190        let mut handle = ptr::null_mut();
191        unsafe {
192            try_ffi!(sys::cusolverDnCreateSyevjInfo(&raw mut handle))?;
193        }
194
195        if handle.is_null() {
196            return Err(Error::NullHandle);
197        }
198
199        Ok(Self { handle })
200    }
201
202    /// Configures the `syevj` tolerance.
203    ///
204    /// # Errors
205    ///
206    /// Returns an error if cuSOLVER rejects the tolerance value.
207    pub fn set_tolerance(&mut self, tolerance: f64) -> Result<()> {
208        unsafe {
209            try_ffi!(sys::cusolverDnXsyevjSetTolerance(self.as_raw(), tolerance,))?;
210        }
211        Ok(())
212    }
213
214    /// Configures the maximum number of `syevj` sweeps.
215    /// The default value is 100.
216    ///
217    /// # Errors
218    ///
219    /// Returns an error if cuSOLVER rejects the sweep count.
220    pub fn set_max_sweeps(&mut self, max_sweeps: i32) -> Result<()> {
221        unsafe {
222            try_ffi!(sys::cusolverDnXsyevjSetMaxSweeps(self.as_raw(), max_sweeps,))?;
223        }
224        Ok(())
225    }
226
227    /// If `sort_eigenvalues` is false, the eigenvalues are not sorted.
228    /// This setting only applies to `syevj_batched`.
229    /// `syevj` and `sygvj` always sort eigenvalues in ascending order.
230    /// By default, eigenvalues are always sorted in ascending order.
231    ///
232    /// # Errors
233    ///
234    /// Returns an error if cuSOLVER rejects the sort setting.
235    pub fn set_sort_eigenvalues(&mut self, sort_eigenvalues: bool) -> Result<()> {
236        unsafe {
237            try_ffi!(sys::cusolverDnXsyevjSetSortEig(
238                self.as_raw(),
239                i32::from(sort_eigenvalues),
240            ))?;
241        }
242        Ok(())
243    }
244
245    /// Returns the residual reported by `syevj` or `sygvj`.
246    /// This accessor does not support `syevj_batched`.
247    /// Calling this after `syevj_batched` returns [`Status::NotSupported`].
248    ///
249    /// # Errors
250    ///
251    /// Returns an error if the info handle was used with `syevj_batched`,
252    /// which does not report a residual.
253    pub fn residual(&self, ctx: &Context) -> Result<f64> {
254        ctx.bind()?;
255
256        let mut residual = 0.0;
257        unsafe {
258            try_ffi!(sys::cusolverDnXsyevjGetResidual(
259                ctx.as_raw(),
260                self.as_raw(),
261                &raw mut residual,
262            ))?;
263        }
264        Ok(residual)
265    }
266
267    /// Returns the number of executed `syevj` or `sygvj` sweeps.
268    /// This accessor does not support `syevj_batched`.
269    /// Calling this after `syevj_batched` returns [`Status::NotSupported`].
270    ///
271    /// # Errors
272    ///
273    /// Returns an error if the info handle was used with `syevj_batched`,
274    /// which does not report a sweep count.
275    pub fn executed_sweeps(&self, ctx: &Context) -> Result<i32> {
276        ctx.bind()?;
277
278        let mut sweeps = 0;
279        unsafe {
280            try_ffi!(sys::cusolverDnXsyevjGetSweeps(
281                ctx.as_raw(),
282                self.as_raw(),
283                &raw mut sweeps,
284            ))?;
285        }
286        Ok(sweeps)
287    }
288
289    pub fn as_raw(&self) -> sys::syevjInfo_t {
290        self.handle
291    }
292}
293
294impl Drop for SyevjInfo {
295    fn drop(&mut self) {
296        unsafe {
297            if let Err(err) = try_ffi!(sys::cusolverDnDestroySyevjInfo(self.handle)) {
298                #[cfg(debug_assertions)]
299                eprintln!("failed to destroy cusolver syevj info: {err}");
300            }
301        }
302    }
303}
304
305pub fn xsyevd_buffer_size<TA: DataTypeLike, TW: DataTypeLike>(
306    ctx: &Context,
307    params: &Params,
308    mode: EigenMode,
309    fill_mode: FillMode,
310    n: usize,
311    a: &DeviceMemory<TA>,
312    lda: usize,
313    w: &DeviceMemory<TW>,
314) -> Result<WorkspaceSizes> {
315    xsyevd_raw_buffer_size(
316        ctx,
317        params,
318        mode,
319        fill_mode,
320        n,
321        TA::data_type(),
322        a,
323        lda,
324        TW::data_type(),
325        w,
326        TA::data_type(),
327    )
328}
329
330pub fn xsyevd<TA: DataTypeLike, TW: DataTypeLike>(
331    ctx: &Context,
332    params: &Params,
333    mode: EigenMode,
334    fill_mode: FillMode,
335    n: usize,
336    a: &mut DeviceMemory<TA>,
337    lda: usize,
338    w: &mut DeviceMemory<TW>,
339    workspace: ByteWorkspaceMut<'_>,
340    dev_info: &mut DeviceMemory<i32>,
341) -> Result<()> {
342    xsyevd_raw(
343        ctx,
344        params,
345        mode,
346        fill_mode,
347        n,
348        TA::data_type(),
349        a,
350        lda,
351        TW::data_type(),
352        w,
353        TA::data_type(),
354        workspace,
355        dev_info,
356    )
357}
358
359pub fn xsyevdx_buffer_size<TA: DataTypeLike, TR: Copy + Default, TW: DataTypeLike>(
360    ctx: &Context,
361    params: &Params,
362    mode: EigenMode,
363    fill_mode: FillMode,
364    selection: EigenSelection<TR>,
365    n: usize,
366    a: &DeviceMemory<TA>,
367    lda: usize,
368    w: &DeviceMemory<TW>,
369) -> Result<SelectionWorkspaceSizes> {
370    let (range, value_range, index_range) = selection_parts(selection);
371    xsyevdx_raw_buffer_size(
372        ctx,
373        params,
374        mode,
375        range,
376        fill_mode,
377        n,
378        TA::data_type(),
379        a,
380        lda,
381        value_range,
382        index_range,
383        TW::data_type(),
384        w,
385        TA::data_type(),
386    )
387}
388
389pub fn xsyevdx<TA: DataTypeLike, TR: Copy + Default, TW: DataTypeLike>(
390    ctx: &Context,
391    params: &Params,
392    mode: EigenMode,
393    fill_mode: FillMode,
394    selection: EigenSelection<TR>,
395    n: usize,
396    a: &mut DeviceMemory<TA>,
397    lda: usize,
398    w: &mut DeviceMemory<TW>,
399    workspace: ByteWorkspaceMut<'_>,
400    dev_info: &mut DeviceMemory<i32>,
401) -> Result<usize> {
402    let (range, value_range, index_range) = selection_parts(selection);
403    xsyevdx_raw(
404        ctx,
405        params,
406        mode,
407        range,
408        fill_mode,
409        n,
410        TA::data_type(),
411        a,
412        lda,
413        value_range,
414        index_range,
415        TW::data_type(),
416        w,
417        TA::data_type(),
418        workspace,
419        dev_info,
420    )
421}
422
423pub fn xsyev_batched_buffer_size<TA: DataTypeLike, TW: DataTypeLike>(
424    ctx: &Context,
425    params: &Params,
426    mode: EigenMode,
427    fill_mode: FillMode,
428    n: usize,
429    a: MatrixRef<'_, TA>,
430    w: &DeviceMemory<TW>,
431    batch_count: usize,
432) -> Result<WorkspaceSizes> {
433    ctx.bind()?;
434    validate_xsyev_batched_buffers(
435        n,
436        a.data.byte_len(),
437        a.leading_dimension,
438        TA::data_type(),
439        w.byte_len(),
440        TW::data_type(),
441        batch_count,
442    )?;
443    let mut device_bytes = 0;
444    let mut host_bytes = 0;
445    unsafe {
446        try_ffi!(sys::cusolverDnXsyevBatched_bufferSize(
447            ctx.as_raw(),
448            params.as_raw(),
449            mode.into(),
450            fill_mode.into(),
451            to_i64(n, "n")?,
452            TA::data_type().into(),
453            a.data.as_ptr().cast(),
454            to_i64(a.leading_dimension, "lda")?,
455            TW::data_type().into(),
456            w.as_ptr().cast(),
457            TA::data_type().into(),
458            &raw mut device_bytes,
459            &raw mut host_bytes,
460            to_i64(batch_count, "batch_count")?,
461        ))?;
462    }
463    Ok(WorkspaceSizes::new(
464        device_bytes as usize,
465        host_bytes as usize,
466    ))
467}
468
469/// Computes eigenvalues and eigenvectors of a sequence of symmetric
470/// (Hermitian) $n \times n$ matrices.
471///
472/// where $\Lambda\_j$ is a real $n \times n$ diagonal matrix. $V\_j$ is an $n \times n$ unitary matrix.
473/// The diagonal elements of $\Lambda\_j$ are the eigenvalues of $A\_j$ in ascending order.
474///
475/// `syevBatched` performs an eigendecomposition on each matrix.
476/// It requires all matrices to have the same size `n` and be packed contiguously.
477///
478/// Each matrix is column-major with leading dimension `lda`, so the formula for random access is $A\_{k}\operatorname{(i,j)} = {A\lbrack\ i\ +\ lda\cdot j\ +\ lda\cdot n\cdot k\rbrack}$.
479///
480/// `w` contains the eigenvalues of each matrix contiguously.
481///
482/// The formula for random access of `W` is $W\_{k}\operatorname{(j)} = {W\lbrack\ j\ +\ n\cdot k\rbrack}$.
483///
484/// Provide device and host workspace through `workspace`.
485/// Use [`xsyev_batched_buffer_size`] to determine the required sizes for
486/// `workspace.device` and `workspace.host`.
487///
488/// `dev_info` has one entry per batch item.
489/// If the call returns [`Status::InvalidValue`], `dev_info[0] == -i` indicates that the `i`th parameter is invalid.
490/// Otherwise, `dev_info[i] > 0` indicates that `syevBatched` did not converge on the `i`th matrix.
491///
492/// If `mode` is [`EigenMode::Vector`], $A\_{j}$ contains the orthonormal eigenvectors of the matrix $A\_{j}$.
493///
494/// The problem size is limited by
495/// `n * lda * batch_count <= INT32_MAX` primarily due to the current implementation
496/// constraints.
497///
498/// **Algorithms supported by [`xsyev_batched`]**
499///
500/// | Algorithm | Notes |
501/// | --- | --- |
502/// | [`AlgorithmMode::Default`](crate::types::AlgorithmMode::Default) | Default; may switch between algorithms for best performance. |
503/// | [`AlgorithmMode::Algorithm1`](crate::types::AlgorithmMode::Algorithm1) | Uses a single algorithm for consistent accuracy over all `n`. |
504///
505/// List of input arguments for [`xsyev_batched_buffer_size`] and [`xsyev_batched`]:
506///
507/// The generic operation has three data types: `data_type_a` is the data type
508/// of matrix `A`, `data_type_w` is the data type of `W`, and `compute_type` is
509/// the compute type of the operation.
510/// [`xsyev_batched`] only supports the following four combinations:
511///
512/// **Valid combination of data type and compute type**
513///
514/// | **data_type_a** | **data_type_w** | **compute_type** | **Meaning** |
515/// | --- | --- | --- | --- |
516/// | [`DataType::F32`] | [`DataType::F32`] | [`DataType::F32`] | `SSYEVBATCHED` |
517/// | [`DataType::F64`] | [`DataType::F64`] | [`DataType::F64`] | `DSYEVBATCHED` |
518/// | [`DataType::ComplexF32`] | [`DataType::F32`] | [`DataType::ComplexF32`] | `CSYEVBATCHED` |
519/// | [`DataType::ComplexF64`] | [`DataType::F64`] | [`DataType::ComplexF64`] | `ZSYEVBATCHED` |
520///
521/// # Errors
522///
523/// Returns an error if cuSOLVER has not been initialized, if the
524/// matrix dimensions, leading dimension, eigen mode, fill mode, or batch size
525/// are invalid, or if cuSOLVER reports an internal failure.
526pub fn xsyev_batched<TA: DataTypeLike, TW: DataTypeLike>(
527    ctx: &Context,
528    params: &Params,
529    mode: EigenMode,
530    fill_mode: FillMode,
531    n: usize,
532    a: MatrixMut<'_, TA>,
533    w: &mut DeviceMemory<TW>,
534    batch_count: usize,
535    workspace: ByteWorkspaceMut<'_>,
536    dev_info: &mut DeviceMemory<i32>,
537) -> Result<()> {
538    ctx.bind()?;
539    validate_xsyev_batched_buffers(
540        n,
541        a.data.byte_len(),
542        a.leading_dimension,
543        TA::data_type(),
544        w.byte_len(),
545        TW::data_type(),
546        batch_count,
547    )?;
548    require_info_buffer_len(dev_info, batch_count)?;
549    let workspace_sizes =
550        xsyev_batched_buffer_size(ctx, params, mode, fill_mode, n, a.as_ref(), w, batch_count)?;
551    require_workspace_bytes(workspace.device.byte_len(), workspace_sizes.device_bytes)?;
552    require_host_workspace(workspace.host.len(), workspace_sizes.host_bytes)?;
553    unsafe {
554        try_ffi!(sys::cusolverDnXsyevBatched(
555            ctx.as_raw(),
556            params.as_raw(),
557            mode.into(),
558            fill_mode.into(),
559            to_i64(n, "n")?,
560            TA::data_type().into(),
561            a.data.as_mut_ptr().cast(),
562            to_i64(a.leading_dimension, "lda")?,
563            TW::data_type().into(),
564            w.as_mut_ptr().cast(),
565            TA::data_type().into(),
566            workspace.device.as_mut_ptr().cast(),
567            workspace_sizes.device_bytes as _,
568            workspace.host.as_mut_ptr().cast(),
569            workspace_sizes.host_bytes as _,
570            dev_info.as_mut_ptr().cast(),
571            to_i64(batch_count, "batch_count")?,
572        ))?;
573    }
574    Ok(())
575}
576
577pub fn xgeev_buffer_size<TA: DataTypeLike, TW: DataTypeLike, TV: DataTypeLike>(
578    ctx: &Context,
579    params: &Params,
580    n: usize,
581    a: MatrixRef<'_, TA>,
582    eigenvalues: &DeviceMemory<TW>,
583    right_vectors: Option<MatrixRef<'_, TV>>,
584) -> Result<WorkspaceSizes> {
585    ctx.bind()?;
586    validate_xgeev_inputs(
587        n,
588        a.data.byte_len(),
589        a.leading_dimension,
590        TA::data_type(),
591        eigenvalues.byte_len(),
592        TW::data_type(),
593        matrix_ref_parts(right_vectors),
594        TV::data_type(),
595    )?;
596    let (vr_ptr, ldvr) = optional_xgeev_matrix_ptr(matrix_ref_parts(right_vectors))?;
597    let mut device_bytes = 0;
598    let mut host_bytes = 0;
599    unsafe {
600        try_ffi!(sys::cusolverDnXgeev_bufferSize(
601            ctx.as_raw(),
602            params.as_raw(),
603            EigenMode::NoVector.into(),
604            if right_vectors.is_some() {
605                EigenMode::Vector
606            } else {
607                EigenMode::NoVector
608            }
609            .into(),
610            to_i64(n, "n")?,
611            TA::data_type().into(),
612            a.data.as_ptr().cast(),
613            to_i64(a.leading_dimension, "lda")?,
614            TW::data_type().into(),
615            eigenvalues.as_ptr().cast(),
616            TA::data_type().into(),
617            ptr::null(),
618            1,
619            TV::data_type().into(),
620            vr_ptr.cast(),
621            ldvr,
622            TA::data_type().into(),
623            &raw mut device_bytes,
624            &raw mut host_bytes,
625        ))?;
626    }
627    Ok(WorkspaceSizes::new(
628        device_bytes as usize,
629        host_bytes as usize,
630    ))
631}
632
633/// Computes the eigenvalues and, optionally, the left and/or right eigenvectors
634/// of an n-by-n real non-symmetric or complex non-Hermitian matrix `A`.
635/// The right eigenvector `v(j)` of `A` satisfies
636///
637/// where `w(j)` is its eigenvalue.
638/// The left eigenvalue `u(j)` of `A` satisfies
639///
640/// where $u(j)^{H}$ denotes the conjugate-transpose of `u(j)`.
641///
642/// The computed eigenvectors are normalized to have Euclidean norm equal to 1 and largest component real.
643///
644/// If `A` is real-valued, there are two options to return the eigenvalues in `W`.
645/// The first option sets all data types to real-valued types.
646/// Then `W` holds `2*n` entries.
647/// The first n entries hold the real parts and the last n entries hold the imaginary parts.
648/// The LAPACK interface with separate arrays for the real parts `WR` and the
649/// imaginary parts `WI` can be recovered by setting pointers `WR = W` and
650/// `WI = W+n`.
651/// The second option uses a complex data type for `W`.
652/// Then `W` is n entries long; each real eigenvalue is stored as a complex number and for each complex conjugate pair, both eigenvalues are returned.
653/// The computation is still executed fully in real arithmetic.
654///
655/// Provide device and host workspace through `workspace`.
656/// Use [`xgeev_buffer_size`] to determine the required sizes for
657/// `workspace.device` and `workspace.host`.
658///
659/// If the reported `info` value is `-i`, the `i`th parameter is invalid.
660/// If `info == 0`, the QR algorithm converged; `W` contains the computed eigenvalues of `A`, and any requested left or right eigenvectors have been computed.
661/// If `info == i` with `i > 0`, the QR algorithm failed to compute all eigenvalues and no eigenvectors were computed.
662/// The elements `i + 1:n` of `W` contain eigenvalues that converged.
663///
664/// - `geev` only supports the computation of right eigenvectors.
665///   Therefore, `jobvl` must be [`EigenMode::NoVector`].
666///
667/// - `geev` uses balancing to improve the conditioning of the eigenvalues and eigenvectors.
668///
669/// - `geev` is a hybrid CPU-GPU algorithm.
670///   Best performance is attained with pinned host memory.
671///
672/// Currently, [`xgeev`] supports only the default algorithm.
673///
674/// **Table of algorithms supported by [`xgeev`]**
675///
676/// | Algorithm | Notes |
677/// | --- | --- |
678/// | [`AlgorithmMode::Default`](crate::types::AlgorithmMode::Default) | Default algorithm. |
679///
680/// List of input arguments for [`xgeev_buffer_size`] and [`xgeev`]:
681///
682/// The generic operation has five data types: `data_type_a` is the data type
683/// of matrix `A`, `data_type_w` is the data type of `W`, `data_type_vl` is the
684/// data type of matrix `VL`, `data_type_vr` is the data type of matrix `VR`,
685/// and `compute_type` is the compute type of the operation.
686/// [`xgeev`] only supports the following four combinations:
687///
688/// **Valid combination of data type and compute type**
689///
690/// | **data_type_a** | **data_type_w** | **data_type_vl** | **data_type_vr** | **compute_type** | **Meaning** |
691/// | --- | --- | --- | --- | --- | --- |
692/// | [`DataType::F32`] | [`DataType::F32`] | [`DataType::F32`] | [`DataType::F32`] | [`DataType::F32`] | `SGEEV` |
693/// | [`DataType::F32`] | [`DataType::ComplexF32`] | [`DataType::F32`] | [`DataType::F32`] | [`DataType::F32`] | 32F mixed real-complex |
694/// | [`DataType::F64`] | [`DataType::F64`] | [`DataType::F64`] | [`DataType::F64`] | [`DataType::F64`] | `DGEEV` |
695/// | [`DataType::F64`] | [`DataType::ComplexF64`] | [`DataType::F64`] | [`DataType::F64`] | [`DataType::F64`] | 64F mixed real-complex |
696/// | [`DataType::ComplexF32`] | [`DataType::ComplexF32`] | [`DataType::ComplexF32`] | [`DataType::ComplexF32`] | [`DataType::ComplexF32`] | `CGEEV` |
697/// | [`DataType::ComplexF64`] | [`DataType::ComplexF64`] | [`DataType::ComplexF64`] | [`DataType::ComplexF64`] | [`DataType::ComplexF64`] | `ZGEEV` |
698///
699/// # Errors
700///
701/// Returns an error if cuSOLVER has not been initialized, if the
702/// matrix dimensions, leading dimensions, or requested eigenvector modes are
703/// invalid, or if cuSOLVER reports an internal failure.
704pub fn xgeev<TA: DataTypeLike, TW: DataTypeLike, TV: DataTypeLike>(
705    ctx: &Context,
706    params: &Params,
707    n: usize,
708    a: MatrixMut<'_, TA>,
709    eigenvalues: &mut DeviceMemory<TW>,
710    right_vectors: Option<MatrixMut<'_, TV>>,
711    workspace: ByteWorkspaceMut<'_>,
712    dev_info: &mut DeviceMemory<i32>,
713) -> Result<()> {
714    ctx.bind()?;
715    validate_xgeev_inputs(
716        n,
717        a.data.byte_len(),
718        a.leading_dimension,
719        TA::data_type(),
720        eigenvalues.byte_len(),
721        TW::data_type(),
722        matrix_mut_ref_parts(right_vectors.as_ref()),
723        TV::data_type(),
724    )?;
725    require_info_buffer(dev_info)?;
726    let workspace_sizes = xgeev_buffer_size(
727        ctx,
728        params,
729        n,
730        a.as_ref(),
731        eigenvalues,
732        matrix_mut_ref_option(right_vectors.as_ref()),
733    )?;
734    require_workspace_bytes(workspace.device.byte_len(), workspace_sizes.device_bytes)?;
735    require_host_workspace(workspace.host.len(), workspace_sizes.host_bytes)?;
736    let (vr_ptr, ldvr) = optional_xgeev_matrix_mut_ptr(matrix_mut_parts(right_vectors))?;
737    unsafe {
738        try_ffi!(sys::cusolverDnXgeev(
739            ctx.as_raw(),
740            params.as_raw(),
741            EigenMode::NoVector.into(),
742            if vr_ptr.is_null() {
743                EigenMode::NoVector
744            } else {
745                EigenMode::Vector
746            }
747            .into(),
748            to_i64(n, "n")?,
749            TA::data_type().into(),
750            a.data.as_mut_ptr().cast(),
751            to_i64(a.leading_dimension, "lda")?,
752            TW::data_type().into(),
753            eigenvalues.as_mut_ptr().cast(),
754            TA::data_type().into(),
755            ptr::null_mut(),
756            1,
757            TV::data_type().into(),
758            vr_ptr.cast(),
759            ldvr,
760            TA::data_type().into(),
761            workspace.device.as_mut_ptr().cast(),
762            workspace_sizes.device_bytes as _,
763            workspace.host.as_mut_ptr().cast(),
764            workspace_sizes.host_bytes as _,
765            dev_info.as_mut_ptr().cast(),
766        ))?;
767    }
768    Ok(())
769}
770
771fn xsyevd_raw_buffer_size<TA, TW>(
772    ctx: &Context,
773    params: &Params,
774    mode: EigenMode,
775    fill_mode: FillMode,
776    n: usize,
777    a_type: DataType,
778    a: &DeviceMemory<TA>,
779    lda: usize,
780    w_type: DataType,
781    w: &DeviceMemory<TW>,
782    compute_type: DataType,
783) -> Result<WorkspaceSizes> {
784    ctx.bind()?;
785    validate_xsyevd_buffers(n, a.byte_len(), lda, a_type, w.byte_len(), w_type)?;
786    let mut device_bytes = 0;
787    let mut host_bytes = 0;
788    unsafe {
789        try_ffi!(sys::cusolverDnXsyevd_bufferSize(
790            ctx.as_raw(),
791            params.as_raw(),
792            mode.into(),
793            fill_mode.into(),
794            to_i64(n, "n")?,
795            a_type.into(),
796            a.as_ptr().cast(),
797            to_i64(lda, "lda")?,
798            w_type.into(),
799            w.as_ptr().cast(),
800            compute_type.into(),
801            &raw mut device_bytes,
802            &raw mut host_bytes,
803        ))?;
804    }
805    Ok(WorkspaceSizes::new(
806        device_bytes as usize,
807        host_bytes as usize,
808    ))
809}
810
811fn xsyevd_raw<TA, TW>(
812    ctx: &Context,
813    params: &Params,
814    mode: EigenMode,
815    fill_mode: FillMode,
816    n: usize,
817    a_type: DataType,
818    a: &mut DeviceMemory<TA>,
819    lda: usize,
820    w_type: DataType,
821    w: &mut DeviceMemory<TW>,
822    compute_type: DataType,
823    workspace: ByteWorkspaceMut<'_>,
824    dev_info: &mut DeviceMemory<i32>,
825) -> Result<()> {
826    ctx.bind()?;
827    validate_xsyevd_buffers(n, a.byte_len(), lda, a_type, w.byte_len(), w_type)?;
828    require_info_buffer(dev_info)?;
829    let workspace_sizes = xsyevd_raw_buffer_size(
830        ctx,
831        params,
832        mode,
833        fill_mode,
834        n,
835        a_type,
836        a,
837        lda,
838        w_type,
839        w,
840        compute_type,
841    )?;
842    require_workspace_bytes(workspace.device.byte_len(), workspace_sizes.device_bytes)?;
843    require_host_workspace(workspace.host.len(), workspace_sizes.host_bytes)?;
844    unsafe {
845        try_ffi!(sys::cusolverDnXsyevd(
846            ctx.as_raw(),
847            params.as_raw(),
848            mode.into(),
849            fill_mode.into(),
850            to_i64(n, "n")?,
851            a_type.into(),
852            a.as_mut_ptr().cast(),
853            to_i64(lda, "lda")?,
854            w_type.into(),
855            w.as_mut_ptr().cast(),
856            compute_type.into(),
857            workspace.device.as_mut_ptr().cast(),
858            workspace_sizes.device_bytes as _,
859            workspace.host.as_mut_ptr().cast(),
860            workspace_sizes.host_bytes as _,
861            dev_info.as_mut_ptr().cast(),
862        ))?;
863    }
864    Ok(())
865}
866
867fn xsyevdx_raw_buffer_size<TA, TR, TW>(
868    ctx: &Context,
869    params: &Params,
870    mode: EigenMode,
871    range: EigenRange,
872    fill_mode: FillMode,
873    n: usize,
874    a_type: DataType,
875    a: &DeviceMemory<TA>,
876    lda: usize,
877    value_range: Option<(TR, TR)>,
878    index_range: Option<(usize, usize)>,
879    w_type: DataType,
880    w: &DeviceMemory<TW>,
881    compute_type: DataType,
882) -> Result<SelectionWorkspaceSizes>
883where
884    TR: Copy + Default,
885{
886    ctx.bind()?;
887    validate_xsyevd_buffers(n, a.byte_len(), lda, a_type, w.byte_len(), w_type)?;
888    validate_xsyevdx_value_type::<TR>(w_type)?;
889    let (mut vl, mut vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
890    let mut meig = 0;
891    let mut device_bytes = 0;
892    let mut host_bytes = 0;
893    unsafe {
894        try_ffi!(sys::cusolverDnXsyevdx_bufferSize(
895            ctx.as_raw(),
896            params.as_raw(),
897            mode.into(),
898            range.into(),
899            fill_mode.into(),
900            to_i64(n, "n")?,
901            a_type.into(),
902            a.as_ptr().cast(),
903            to_i64(lda, "lda")?,
904            (&raw mut vl).cast(),
905            (&raw mut vu).cast(),
906            to_i64(il, "il")?,
907            to_i64(iu, "iu")?,
908            &raw mut meig,
909            w_type.into(),
910            w.as_ptr().cast(),
911            compute_type.into(),
912            &raw mut device_bytes,
913            &raw mut host_bytes,
914        ))?;
915    }
916    Ok(SelectionWorkspaceSizes::new(
917        to_usize(meig, "meig")?,
918        device_bytes as usize,
919        host_bytes as usize,
920    ))
921}
922
923fn xsyevdx_raw<TA, TR, TW>(
924    ctx: &Context,
925    params: &Params,
926    mode: EigenMode,
927    range: EigenRange,
928    fill_mode: FillMode,
929    n: usize,
930    a_type: DataType,
931    a: &mut DeviceMemory<TA>,
932    lda: usize,
933    value_range: Option<(TR, TR)>,
934    index_range: Option<(usize, usize)>,
935    w_type: DataType,
936    w: &mut DeviceMemory<TW>,
937    compute_type: DataType,
938    workspace: ByteWorkspaceMut<'_>,
939    dev_info: &mut DeviceMemory<i32>,
940) -> Result<usize>
941where
942    TR: Copy + Default,
943{
944    ctx.bind()?;
945    validate_xsyevd_buffers(n, a.byte_len(), lda, a_type, w.byte_len(), w_type)?;
946    validate_xsyevdx_value_type::<TR>(w_type)?;
947    require_info_buffer(dev_info)?;
948    let workspace_sizes = xsyevdx_raw_buffer_size(
949        ctx,
950        params,
951        mode,
952        range,
953        fill_mode,
954        n,
955        a_type,
956        a,
957        lda,
958        value_range,
959        index_range,
960        w_type,
961        w,
962        compute_type,
963    )?;
964    require_workspace_bytes(
965        workspace.device.byte_len(),
966        workspace_sizes.workspace.device_bytes,
967    )?;
968    require_host_workspace(workspace.host.len(), workspace_sizes.workspace.host_bytes)?;
969    let (mut vl, mut vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
970    let mut meig_raw = 0;
971    unsafe {
972        try_ffi!(sys::cusolverDnXsyevdx(
973            ctx.as_raw(),
974            params.as_raw(),
975            mode.into(),
976            range.into(),
977            fill_mode.into(),
978            to_i64(n, "n")?,
979            a_type.into(),
980            a.as_mut_ptr().cast(),
981            to_i64(lda, "lda")?,
982            (&raw mut vl).cast(),
983            (&raw mut vu).cast(),
984            to_i64(il, "il")?,
985            to_i64(iu, "iu")?,
986            &raw mut meig_raw,
987            w_type.into(),
988            w.as_mut_ptr().cast(),
989            compute_type.into(),
990            workspace.device.as_mut_ptr().cast(),
991            workspace_sizes.workspace.device_bytes as _,
992            workspace.host.as_mut_ptr().cast(),
993            workspace_sizes.workspace.host_bytes as _,
994            dev_info.as_mut_ptr().cast(),
995        ))?;
996    }
997    debug_assert_eq!(workspace_sizes.selection_size, to_usize(meig_raw, "meig")?);
998    Ok(workspace_sizes.selection_size)
999}
1000
1001pub fn ssyevd_buffer_size(
1002    ctx: &Context,
1003    mode: EigenMode,
1004    fill_mode: FillMode,
1005    n: usize,
1006    a: &DeviceMemory<f32>,
1007    lda: usize,
1008    w: &DeviceMemory<f32>,
1009) -> Result<usize> {
1010    ctx.bind()?;
1011    validate_syev_buffers(n, a.len(), lda, w.len())?;
1012    let mut lwork = 0;
1013    unsafe {
1014        try_ffi!(sys::cusolverDnSsyevd_bufferSize(
1015            ctx.as_raw(),
1016            mode.into(),
1017            fill_mode.into(),
1018            to_i32(n, "n")?,
1019            a.as_ptr().cast(),
1020            to_i32(lda, "lda")?,
1021            w.as_ptr().cast(),
1022            &raw mut lwork,
1023        ))?;
1024    }
1025    to_usize(lwork, "lwork")
1026}
1027
1028pub fn dsyevd_buffer_size(
1029    ctx: &Context,
1030    mode: EigenMode,
1031    fill_mode: FillMode,
1032    n: usize,
1033    a: &DeviceMemory<f64>,
1034    lda: usize,
1035    w: &DeviceMemory<f64>,
1036) -> Result<usize> {
1037    ctx.bind()?;
1038    validate_syev_buffers(n, a.len(), lda, w.len())?;
1039    let mut lwork = 0;
1040    unsafe {
1041        try_ffi!(sys::cusolverDnDsyevd_bufferSize(
1042            ctx.as_raw(),
1043            mode.into(),
1044            fill_mode.into(),
1045            to_i32(n, "n")?,
1046            a.as_ptr().cast(),
1047            to_i32(lda, "lda")?,
1048            w.as_ptr().cast(),
1049            &raw mut lwork,
1050        ))?;
1051    }
1052    to_usize(lwork, "lwork")
1053}
1054
1055pub fn cheevd_buffer_size(
1056    ctx: &Context,
1057    mode: EigenMode,
1058    fill_mode: FillMode,
1059    n: usize,
1060    a: &DeviceMemory<Complex32>,
1061    lda: usize,
1062    w: &DeviceMemory<f32>,
1063) -> Result<usize> {
1064    ctx.bind()?;
1065    validate_syev_buffers(n, a.len(), lda, w.len())?;
1066    let mut lwork = 0;
1067    unsafe {
1068        try_ffi!(sys::cusolverDnCheevd_bufferSize(
1069            ctx.as_raw(),
1070            mode.into(),
1071            fill_mode.into(),
1072            to_i32(n, "n")?,
1073            a.as_ptr().cast(),
1074            to_i32(lda, "lda")?,
1075            w.as_ptr().cast(),
1076            &raw mut lwork,
1077        ))?;
1078    }
1079    to_usize(lwork, "lwork")
1080}
1081
1082pub fn zheevd_buffer_size(
1083    ctx: &Context,
1084    mode: EigenMode,
1085    fill_mode: FillMode,
1086    n: usize,
1087    a: &DeviceMemory<Complex64>,
1088    lda: usize,
1089    w: &DeviceMemory<f64>,
1090) -> Result<usize> {
1091    ctx.bind()?;
1092    validate_syev_buffers(n, a.len(), lda, w.len())?;
1093    let mut lwork = 0;
1094    unsafe {
1095        try_ffi!(sys::cusolverDnZheevd_bufferSize(
1096            ctx.as_raw(),
1097            mode.into(),
1098            fill_mode.into(),
1099            to_i32(n, "n")?,
1100            a.as_ptr().cast(),
1101            to_i32(lda, "lda")?,
1102            w.as_ptr().cast(),
1103            &raw mut lwork,
1104        ))?;
1105    }
1106    to_usize(lwork, "lwork")
1107}
1108
1109/// Use the matching buffer-size helper to calculate the workspace length before
1110/// allocating workspace.
1111///
1112/// The S and D data types are real valued single and double precision, respectively.
1113///
1114/// The C and Z data types are complex valued single and double precision, respectively.
1115///
1116/// Computes eigenvalues and eigenvectors of a symmetric (Hermitian)
1117/// $n \times n$ matrix `A`.
1118/// The standard symmetric eigenvalue problem is $A Q = Q \Lambda$, where `Λ` is a real $n \times n$ diagonal matrix.
1119/// `V` is an $n \times n$ unitary matrix.
1120/// The diagonal elements of `Λ` are the eigenvalues of `A` in ascending order.
1121///
1122/// Provide workspace through `workspace`.
1123/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
1124/// The workspace size in bytes is `size_of::<T>() * lwork`.
1125///
1126/// If the reported `dev_info` value is `-i`, the `i`th parameter is invalid.
1127/// If `dev_info = i` (greater than zero), `i` off-diagonal elements of an intermediate tridiagonal form did not converge to zero.
1128///
1129/// If `mode` is [`EigenMode::Vector`], `A` contains the orthonormal eigenvectors of the matrix `A`.
1130/// The eigenvectors are computed by a divide and conquer algorithm.
1131///
1132/// # Errors
1133///
1134/// Returns an error if cuSOLVER has not been initialized, if the
1135/// matrix dimensions, leading dimension, eigen mode, or fill mode are invalid,
1136/// if the current GPU architecture is unsupported, or if cuSOLVER reports an
1137/// internal failure.
1138pub fn ssyevd(
1139    ctx: &Context,
1140    mode: EigenMode,
1141    fill_mode: FillMode,
1142    n: usize,
1143    a: &mut DeviceMemory<f32>,
1144    lda: usize,
1145    w: &mut DeviceMemory<f32>,
1146    workspace: &mut DeviceMemory<f32>,
1147    dev_info: &mut DeviceMemory<i32>,
1148) -> Result<()> {
1149    ctx.bind()?;
1150    validate_syev_buffers(n, a.len(), lda, w.len())?;
1151    require_info_buffer(dev_info)?;
1152    let lwork = ssyevd_buffer_size(ctx, mode, fill_mode, n, a, lda, w)?;
1153    require_workspace(workspace.len(), lwork)?;
1154    unsafe {
1155        try_ffi!(sys::cusolverDnSsyevd(
1156            ctx.as_raw(),
1157            mode.into(),
1158            fill_mode.into(),
1159            to_i32(n, "n")?,
1160            a.as_mut_ptr().cast(),
1161            to_i32(lda, "lda")?,
1162            w.as_mut_ptr().cast(),
1163            workspace.as_mut_ptr().cast(),
1164            to_i32(lwork, "lwork")?,
1165            dev_info.as_mut_ptr().cast(),
1166        ))?;
1167    }
1168    Ok(())
1169}
1170
1171/// Use the matching buffer-size helper to calculate the workspace length before
1172/// allocating workspace.
1173///
1174/// The S and D data types are real valued single and double precision, respectively.
1175///
1176/// The C and Z data types are complex valued single and double precision, respectively.
1177///
1178/// Computes eigenvalues and eigenvectors of a symmetric (Hermitian)
1179/// $n \times n$ matrix `A`.
1180/// The standard symmetric eigenvalue problem is $A Q = Q \Lambda$, where `Λ` is a real $n \times n$ diagonal matrix.
1181/// `V` is an $n \times n$ unitary matrix.
1182/// The diagonal elements of `Λ` are the eigenvalues of `A` in ascending order.
1183///
1184/// Provide workspace through `workspace`.
1185/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
1186/// The workspace size in bytes is `size_of::<T>() * lwork`.
1187///
1188/// If the reported `dev_info` value is `-i`, the `i`th parameter is invalid.
1189/// If `dev_info = i` (greater than zero), `i` off-diagonal elements of an intermediate tridiagonal form did not converge to zero.
1190///
1191/// If `mode` is [`EigenMode::Vector`], `A` contains the orthonormal eigenvectors of the matrix `A`.
1192/// The eigenvectors are computed by a divide and conquer algorithm.
1193///
1194/// # Errors
1195///
1196/// Returns an error if cuSOLVER has not been initialized, if the
1197/// matrix dimensions, leading dimension, eigen mode, or fill mode are invalid,
1198/// if the current GPU architecture is unsupported, or if cuSOLVER reports an
1199/// internal failure.
1200pub fn dsyevd(
1201    ctx: &Context,
1202    mode: EigenMode,
1203    fill_mode: FillMode,
1204    n: usize,
1205    a: &mut DeviceMemory<f64>,
1206    lda: usize,
1207    w: &mut DeviceMemory<f64>,
1208    workspace: &mut DeviceMemory<f64>,
1209    dev_info: &mut DeviceMemory<i32>,
1210) -> Result<()> {
1211    ctx.bind()?;
1212    validate_syev_buffers(n, a.len(), lda, w.len())?;
1213    require_info_buffer(dev_info)?;
1214    let lwork = dsyevd_buffer_size(ctx, mode, fill_mode, n, a, lda, w)?;
1215    require_workspace(workspace.len(), lwork)?;
1216    unsafe {
1217        try_ffi!(sys::cusolverDnDsyevd(
1218            ctx.as_raw(),
1219            mode.into(),
1220            fill_mode.into(),
1221            to_i32(n, "n")?,
1222            a.as_mut_ptr().cast(),
1223            to_i32(lda, "lda")?,
1224            w.as_mut_ptr().cast(),
1225            workspace.as_mut_ptr().cast(),
1226            to_i32(lwork, "lwork")?,
1227            dev_info.as_mut_ptr().cast(),
1228        ))?;
1229    }
1230    Ok(())
1231}
1232
1233pub fn cheevd(
1234    ctx: &Context,
1235    mode: EigenMode,
1236    fill_mode: FillMode,
1237    n: usize,
1238    a: &mut DeviceMemory<Complex32>,
1239    lda: usize,
1240    w: &mut DeviceMemory<f32>,
1241    workspace: &mut DeviceMemory<Complex32>,
1242    dev_info: &mut DeviceMemory<i32>,
1243) -> Result<()> {
1244    ctx.bind()?;
1245    validate_syev_buffers(n, a.len(), lda, w.len())?;
1246    require_info_buffer(dev_info)?;
1247    let lwork = cheevd_buffer_size(ctx, mode, fill_mode, n, a, lda, w)?;
1248    require_workspace(workspace.len(), lwork)?;
1249    unsafe {
1250        try_ffi!(sys::cusolverDnCheevd(
1251            ctx.as_raw(),
1252            mode.into(),
1253            fill_mode.into(),
1254            to_i32(n, "n")?,
1255            a.as_mut_ptr().cast(),
1256            to_i32(lda, "lda")?,
1257            w.as_mut_ptr().cast(),
1258            workspace.as_mut_ptr().cast(),
1259            to_i32(lwork, "lwork")?,
1260            dev_info.as_mut_ptr().cast(),
1261        ))?;
1262    }
1263    Ok(())
1264}
1265
1266pub fn zheevd(
1267    ctx: &Context,
1268    mode: EigenMode,
1269    fill_mode: FillMode,
1270    n: usize,
1271    a: &mut DeviceMemory<Complex64>,
1272    lda: usize,
1273    w: &mut DeviceMemory<f64>,
1274    workspace: &mut DeviceMemory<Complex64>,
1275    dev_info: &mut DeviceMemory<i32>,
1276) -> Result<()> {
1277    ctx.bind()?;
1278    validate_syev_buffers(n, a.len(), lda, w.len())?;
1279    require_info_buffer(dev_info)?;
1280    let lwork = zheevd_buffer_size(ctx, mode, fill_mode, n, a, lda, w)?;
1281    require_workspace(workspace.len(), lwork)?;
1282    unsafe {
1283        try_ffi!(sys::cusolverDnZheevd(
1284            ctx.as_raw(),
1285            mode.into(),
1286            fill_mode.into(),
1287            to_i32(n, "n")?,
1288            a.as_mut_ptr().cast(),
1289            to_i32(lda, "lda")?,
1290            w.as_mut_ptr().cast(),
1291            workspace.as_mut_ptr().cast(),
1292            to_i32(lwork, "lwork")?,
1293            dev_info.as_mut_ptr().cast(),
1294        ))?;
1295    }
1296    Ok(())
1297}
1298
1299pub fn ssyevj_buffer_size(
1300    ctx: &Context,
1301    mode: EigenMode,
1302    fill_mode: FillMode,
1303    n: usize,
1304    a: &DeviceMemory<f32>,
1305    lda: usize,
1306    w: &DeviceMemory<f32>,
1307    params: &SyevjInfo,
1308) -> Result<usize> {
1309    ctx.bind()?;
1310    validate_syev_buffers(n, a.len(), lda, w.len())?;
1311    let mut lwork = 0;
1312    unsafe {
1313        try_ffi!(sys::cusolverDnSsyevj_bufferSize(
1314            ctx.as_raw(),
1315            mode.into(),
1316            fill_mode.into(),
1317            to_i32(n, "n")?,
1318            a.as_ptr().cast(),
1319            to_i32(lda, "lda")?,
1320            w.as_ptr().cast(),
1321            &raw mut lwork,
1322            params.as_raw(),
1323        ))?;
1324    }
1325    to_usize(lwork, "lwork")
1326}
1327
1328pub fn dsyevj_buffer_size(
1329    ctx: &Context,
1330    mode: EigenMode,
1331    fill_mode: FillMode,
1332    n: usize,
1333    a: &DeviceMemory<f64>,
1334    lda: usize,
1335    w: &DeviceMemory<f64>,
1336    params: &SyevjInfo,
1337) -> Result<usize> {
1338    ctx.bind()?;
1339    validate_syev_buffers(n, a.len(), lda, w.len())?;
1340    let mut lwork = 0;
1341    unsafe {
1342        try_ffi!(sys::cusolverDnDsyevj_bufferSize(
1343            ctx.as_raw(),
1344            mode.into(),
1345            fill_mode.into(),
1346            to_i32(n, "n")?,
1347            a.as_ptr().cast(),
1348            to_i32(lda, "lda")?,
1349            w.as_ptr().cast(),
1350            &raw mut lwork,
1351            params.as_raw(),
1352        ))?;
1353    }
1354    to_usize(lwork, "lwork")
1355}
1356
1357pub fn cheevj_buffer_size(
1358    ctx: &Context,
1359    mode: EigenMode,
1360    fill_mode: FillMode,
1361    n: usize,
1362    a: &DeviceMemory<Complex32>,
1363    lda: usize,
1364    w: &DeviceMemory<f32>,
1365    params: &SyevjInfo,
1366) -> Result<usize> {
1367    ctx.bind()?;
1368    validate_syev_buffers(n, a.len(), lda, w.len())?;
1369    let mut lwork = 0;
1370    unsafe {
1371        try_ffi!(sys::cusolverDnCheevj_bufferSize(
1372            ctx.as_raw(),
1373            mode.into(),
1374            fill_mode.into(),
1375            to_i32(n, "n")?,
1376            a.as_ptr().cast(),
1377            to_i32(lda, "lda")?,
1378            w.as_ptr().cast(),
1379            &raw mut lwork,
1380            params.as_raw(),
1381        ))?;
1382    }
1383    to_usize(lwork, "lwork")
1384}
1385
1386pub fn zheevj_buffer_size(
1387    ctx: &Context,
1388    mode: EigenMode,
1389    fill_mode: FillMode,
1390    n: usize,
1391    a: &DeviceMemory<Complex64>,
1392    lda: usize,
1393    w: &DeviceMemory<f64>,
1394    params: &SyevjInfo,
1395) -> Result<usize> {
1396    ctx.bind()?;
1397    validate_syev_buffers(n, a.len(), lda, w.len())?;
1398    let mut lwork = 0;
1399    unsafe {
1400        try_ffi!(sys::cusolverDnZheevj_bufferSize(
1401            ctx.as_raw(),
1402            mode.into(),
1403            fill_mode.into(),
1404            to_i32(n, "n")?,
1405            a.as_ptr().cast(),
1406            to_i32(lda, "lda")?,
1407            w.as_ptr().cast(),
1408            &raw mut lwork,
1409            params.as_raw(),
1410        ))?;
1411    }
1412    to_usize(lwork, "lwork")
1413}
1414
1415/// Use the matching buffer-size helper to calculate the workspace length before
1416/// allocating workspace.
1417///
1418/// The S and D data types are real valued single and double precision, respectively.
1419///
1420/// The C and Z data types are complex valued single and double precision, respectively.
1421///
1422/// Computes eigenvalues and eigenvectors of a symmetric (Hermitian)
1423/// $n \times n$ matrix `A`.
1424/// The standard symmetric eigenvalue problem is $A Q = Q \Lambda$, where `Λ` is a real $n \times n$ diagonal matrix.
1425/// `Q` is an $n \times n$ unitary matrix.
1426/// The diagonal elements of `Λ` are the eigenvalues of `A` in ascending order.
1427///
1428/// `syevj` computes the same result as `syevd`, but uses the Jacobi method
1429/// instead of the QR algorithm.
1430/// The Jacobi method's parallelism gives GPUs better performance on small and
1431/// medium-size matrices.
1432/// `syevj` can also be configured to approximate results up to a chosen
1433/// accuracy.
1434///
1435/// `syevj` iteratively generates a sequence of unitary matrices that transform `A` toward $A = Q(W + E)Q^{H}$, where `W` is diagonal and `E` is symmetric with a zero diagonal.
1436///
1437/// During the iterations, the Frobenius norm of `E` decreases monotonically.
1438/// As `E` goes down to zero, `W` is the set of eigenvalues.
1439/// In practice, the Jacobi method stops when the off-diagonal residual is below the configured tolerance `eps`.
1440///
1441/// `syevj` has two parameters to control the accuracy.
1442/// The first parameter is the tolerance (`eps`).
1443/// The default value is machine accuracy, but [`SyevjInfo::set_tolerance`] can set an a priori tolerance.
1444/// The maximum-sweep parameter is the maximum number of sweeps, which controls the
1445/// number of Jacobi iterations.
1446/// The default value is 100, but [`SyevjInfo::set_max_sweeps`] can set a different bound.
1447/// Experiments show that 15 sweeps are typically enough to converge to machine
1448/// accuracy.
1449/// `syevj` stops when either the tolerance or the maximum number of sweeps is
1450/// reached.
1451///
1452/// The Jacobi method has quadratic convergence, so accuracy is not
1453/// proportional to the number of sweeps.
1454/// To target a specific accuracy, configure the tolerance.
1455///
1456/// After `syevj`, callers can query the residual with [`SyevjInfo::residual`] and the number of executed sweeps with [`SyevjInfo::executed_sweeps`].
1457/// However, the residual is the Frobenius norm of `E`, not the accuracy of each individual eigenvalue.
1458///
1459/// Provide workspace through `workspace`, just like [`ssyevd`] and [`dsyevd`].
1460/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
1461/// The workspace size in bytes is `size_of::<T>() * lwork`.
1462///
1463/// If the reported `info` value is `-i`, the `i`th parameter is invalid.
1464/// If `info == n + 1`, `syevj` did not converge within the given tolerance and maximum sweep count.
1465///
1466/// If the tolerance is too small, `syevj` may not converge.
1467/// Use a tolerance no smaller than machine accuracy.
1468///
1469/// If `mode` is [`EigenMode::Vector`], `A` contains the orthonormal eigenvectors `V`.
1470///
1471/// # Errors
1472///
1473/// Returns an error if cuSOLVER has not been initialized, if the
1474/// matrix dimensions, leading dimension, eigen mode, or fill mode are invalid,
1475/// or if cuSOLVER reports an internal failure.
1476pub fn ssyevj(
1477    ctx: &Context,
1478    mode: EigenMode,
1479    fill_mode: FillMode,
1480    n: usize,
1481    a: &mut DeviceMemory<f32>,
1482    lda: usize,
1483    w: &mut DeviceMemory<f32>,
1484    workspace: &mut DeviceMemory<f32>,
1485    dev_info: &mut DeviceMemory<i32>,
1486    params: &SyevjInfo,
1487) -> Result<()> {
1488    ctx.bind()?;
1489    validate_syev_buffers(n, a.len(), lda, w.len())?;
1490    require_info_buffer(dev_info)?;
1491    let lwork = ssyevj_buffer_size(ctx, mode, fill_mode, n, a, lda, w, params)?;
1492    require_workspace(workspace.len(), lwork)?;
1493    unsafe {
1494        try_ffi!(sys::cusolverDnSsyevj(
1495            ctx.as_raw(),
1496            mode.into(),
1497            fill_mode.into(),
1498            to_i32(n, "n")?,
1499            a.as_mut_ptr().cast(),
1500            to_i32(lda, "lda")?,
1501            w.as_mut_ptr().cast(),
1502            workspace.as_mut_ptr().cast(),
1503            to_i32(lwork, "lwork")?,
1504            dev_info.as_mut_ptr().cast(),
1505            params.as_raw(),
1506        ))?;
1507    }
1508    Ok(())
1509}
1510
1511/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
1512///
1513/// The S and D data types are real valued single and double precision, respectively.
1514///
1515/// The C and Z data types are complex valued single and double precision, respectively.
1516///
1517/// Computes eigenvalues and eigenvectors of a symmetric (Hermitian) $n \times n$ matrix `A`.
1518/// The standard symmetric eigenvalue problem is $A Q = Q \Lambda$, where `Λ` is a real $n \times n$ diagonal matrix.
1519/// `Q` is an $n \times n$ unitary matrix.
1520/// The diagonal elements of `Λ` are the eigenvalues of `A` in ascending order.
1521///
1522/// `syevj` computes the same symmetric eigenvalue problem as `syevd`.
1523/// The difference is that `syevd` uses a QR algorithm and `syevj` uses the Jacobi method.
1524/// The Jacobi method gives GPUs better parallelism on small and medium-size matrices.
1525/// Callers can configure `syevj` to target a chosen accuracy.
1526///
1527/// `syevj` iteratively generates a sequence of unitary matrices that transform `A` toward $A = Q(W + E)Q^{H}$, where `W` is diagonal and `E` is symmetric with a zero diagonal.
1528///
1529/// During the iterations, the Frobenius norm of `E` decreases monotonically.
1530/// As `E` goes down to zero, `W` is the set of eigenvalues.
1531/// In practice, the Jacobi method stops when the off-diagonal residual is below the configured tolerance `eps`.
1532///
1533/// `syevj` has two parameters to control the accuracy.
1534/// The first parameter is the tolerance (`eps`).
1535/// The default value is machine accuracy, but [`SyevjInfo::set_tolerance`] can set an a priori tolerance.
1536/// The maximum-sweep parameter is the maximum number of sweeps, which controls the number of Jacobi iterations.
1537/// The default value is 100, but [`SyevjInfo::set_max_sweeps`] can set a different bound.
1538/// Experiments show that 15 sweeps are enough to converge to machine accuracy.
1539/// `syevj` stops when either the tolerance or the maximum number of sweeps is reached.
1540///
1541/// The Jacobi method has quadratic convergence, so the accuracy is not proportional to the number of sweeps.
1542/// To guarantee a target accuracy, configure only the tolerance.
1543///
1544/// After `syevj`, callers can query the residual with [`SyevjInfo::residual`] and the number of executed sweeps with [`SyevjInfo::executed_sweeps`].
1545/// However, the residual is the Frobenius norm of `E`, not the accuracy of each individual eigenvalue.
1546///
1547/// Provide workspace through `workspace`, just like [`ssyevd`] and [`dsyevd`].
1548/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
1549/// The workspace size in bytes is `size_of::<T>() * lwork`.
1550///
1551/// If the reported `info` value is `-i`, the `i`th parameter is invalid.
1552/// If `info == n + 1`, `syevj` did not converge within the given tolerance and maximum sweep count.
1553///
1554/// If the tolerance is too small, `syevj` may not converge.
1555/// Use a tolerance no smaller than machine accuracy.
1556///
1557/// If `mode` is [`EigenMode::Vector`], `A` contains the orthonormal eigenvectors `V`.
1558///
1559/// # Errors
1560///
1561/// Returns an error if cuSOLVER has not been initialized, if the
1562/// matrix dimensions, leading dimension, eigen mode, or fill mode are invalid,
1563/// or if cuSOLVER reports an internal failure.
1564pub fn dsyevj(
1565    ctx: &Context,
1566    mode: EigenMode,
1567    fill_mode: FillMode,
1568    n: usize,
1569    a: &mut DeviceMemory<f64>,
1570    lda: usize,
1571    w: &mut DeviceMemory<f64>,
1572    workspace: &mut DeviceMemory<f64>,
1573    dev_info: &mut DeviceMemory<i32>,
1574    params: &SyevjInfo,
1575) -> Result<()> {
1576    ctx.bind()?;
1577    validate_syev_buffers(n, a.len(), lda, w.len())?;
1578    require_info_buffer(dev_info)?;
1579    let lwork = dsyevj_buffer_size(ctx, mode, fill_mode, n, a, lda, w, params)?;
1580    require_workspace(workspace.len(), lwork)?;
1581    unsafe {
1582        try_ffi!(sys::cusolverDnDsyevj(
1583            ctx.as_raw(),
1584            mode.into(),
1585            fill_mode.into(),
1586            to_i32(n, "n")?,
1587            a.as_mut_ptr().cast(),
1588            to_i32(lda, "lda")?,
1589            w.as_mut_ptr().cast(),
1590            workspace.as_mut_ptr().cast(),
1591            to_i32(lwork, "lwork")?,
1592            dev_info.as_mut_ptr().cast(),
1593            params.as_raw(),
1594        ))?;
1595    }
1596    Ok(())
1597}
1598
1599pub fn cheevj(
1600    ctx: &Context,
1601    mode: EigenMode,
1602    fill_mode: FillMode,
1603    n: usize,
1604    a: &mut DeviceMemory<Complex32>,
1605    lda: usize,
1606    w: &mut DeviceMemory<f32>,
1607    workspace: &mut DeviceMemory<Complex32>,
1608    dev_info: &mut DeviceMemory<i32>,
1609    params: &SyevjInfo,
1610) -> Result<()> {
1611    ctx.bind()?;
1612    validate_syev_buffers(n, a.len(), lda, w.len())?;
1613    require_info_buffer(dev_info)?;
1614    let lwork = cheevj_buffer_size(ctx, mode, fill_mode, n, a, lda, w, params)?;
1615    require_workspace(workspace.len(), lwork)?;
1616    unsafe {
1617        try_ffi!(sys::cusolverDnCheevj(
1618            ctx.as_raw(),
1619            mode.into(),
1620            fill_mode.into(),
1621            to_i32(n, "n")?,
1622            a.as_mut_ptr().cast(),
1623            to_i32(lda, "lda")?,
1624            w.as_mut_ptr().cast(),
1625            workspace.as_mut_ptr().cast(),
1626            to_i32(lwork, "lwork")?,
1627            dev_info.as_mut_ptr().cast(),
1628            params.as_raw(),
1629        ))?;
1630    }
1631    Ok(())
1632}
1633
1634pub fn zheevj(
1635    ctx: &Context,
1636    mode: EigenMode,
1637    fill_mode: FillMode,
1638    n: usize,
1639    a: &mut DeviceMemory<Complex64>,
1640    lda: usize,
1641    w: &mut DeviceMemory<f64>,
1642    workspace: &mut DeviceMemory<Complex64>,
1643    dev_info: &mut DeviceMemory<i32>,
1644    params: &SyevjInfo,
1645) -> Result<()> {
1646    ctx.bind()?;
1647    validate_syev_buffers(n, a.len(), lda, w.len())?;
1648    require_info_buffer(dev_info)?;
1649    let lwork = zheevj_buffer_size(ctx, mode, fill_mode, n, a, lda, w, params)?;
1650    require_workspace(workspace.len(), lwork)?;
1651    unsafe {
1652        try_ffi!(sys::cusolverDnZheevj(
1653            ctx.as_raw(),
1654            mode.into(),
1655            fill_mode.into(),
1656            to_i32(n, "n")?,
1657            a.as_mut_ptr().cast(),
1658            to_i32(lda, "lda")?,
1659            w.as_mut_ptr().cast(),
1660            workspace.as_mut_ptr().cast(),
1661            to_i32(lwork, "lwork")?,
1662            dev_info.as_mut_ptr().cast(),
1663            params.as_raw(),
1664        ))?;
1665    }
1666    Ok(())
1667}
1668
1669pub fn ssyevj_batched_buffer_size(
1670    ctx: &Context,
1671    mode: EigenMode,
1672    fill_mode: FillMode,
1673    n: usize,
1674    a: &DeviceMemory<f32>,
1675    lda: usize,
1676    w: &DeviceMemory<f32>,
1677    params: &SyevjInfo,
1678    batch_count: usize,
1679) -> Result<usize> {
1680    ctx.bind()?;
1681    validate_syevj_batched_buffers(n, a.len(), lda, w.len(), batch_count)?;
1682    let mut lwork = 0;
1683    unsafe {
1684        try_ffi!(sys::cusolverDnSsyevjBatched_bufferSize(
1685            ctx.as_raw(),
1686            mode.into(),
1687            fill_mode.into(),
1688            to_i32(n, "n")?,
1689            a.as_ptr().cast(),
1690            to_i32(lda, "lda")?,
1691            w.as_ptr().cast(),
1692            &raw mut lwork,
1693            params.as_raw(),
1694            to_i32(batch_count, "batch_count")?,
1695        ))?;
1696    }
1697    to_usize(lwork, "lwork")
1698}
1699
1700pub fn dsyevj_batched_buffer_size(
1701    ctx: &Context,
1702    mode: EigenMode,
1703    fill_mode: FillMode,
1704    n: usize,
1705    a: &DeviceMemory<f64>,
1706    lda: usize,
1707    w: &DeviceMemory<f64>,
1708    params: &SyevjInfo,
1709    batch_count: usize,
1710) -> Result<usize> {
1711    ctx.bind()?;
1712    validate_syevj_batched_buffers(n, a.len(), lda, w.len(), batch_count)?;
1713    let mut lwork = 0;
1714    unsafe {
1715        try_ffi!(sys::cusolverDnDsyevjBatched_bufferSize(
1716            ctx.as_raw(),
1717            mode.into(),
1718            fill_mode.into(),
1719            to_i32(n, "n")?,
1720            a.as_ptr().cast(),
1721            to_i32(lda, "lda")?,
1722            w.as_ptr().cast(),
1723            &raw mut lwork,
1724            params.as_raw(),
1725            to_i32(batch_count, "batch_count")?,
1726        ))?;
1727    }
1728    to_usize(lwork, "lwork")
1729}
1730
1731pub fn cheevj_batched_buffer_size(
1732    ctx: &Context,
1733    mode: EigenMode,
1734    fill_mode: FillMode,
1735    n: usize,
1736    a: &DeviceMemory<Complex32>,
1737    lda: usize,
1738    w: &DeviceMemory<f32>,
1739    params: &SyevjInfo,
1740    batch_count: usize,
1741) -> Result<usize> {
1742    ctx.bind()?;
1743    validate_syevj_batched_buffers(n, a.len(), lda, w.len(), batch_count)?;
1744    let mut lwork = 0;
1745    unsafe {
1746        try_ffi!(sys::cusolverDnCheevjBatched_bufferSize(
1747            ctx.as_raw(),
1748            mode.into(),
1749            fill_mode.into(),
1750            to_i32(n, "n")?,
1751            a.as_ptr().cast(),
1752            to_i32(lda, "lda")?,
1753            w.as_ptr().cast(),
1754            &raw mut lwork,
1755            params.as_raw(),
1756            to_i32(batch_count, "batch_count")?,
1757        ))?;
1758    }
1759    to_usize(lwork, "lwork")
1760}
1761
1762pub fn zheevj_batched_buffer_size(
1763    ctx: &Context,
1764    mode: EigenMode,
1765    fill_mode: FillMode,
1766    n: usize,
1767    a: &DeviceMemory<Complex64>,
1768    lda: usize,
1769    w: &DeviceMemory<f64>,
1770    params: &SyevjInfo,
1771    batch_count: usize,
1772) -> Result<usize> {
1773    ctx.bind()?;
1774    validate_syevj_batched_buffers(n, a.len(), lda, w.len(), batch_count)?;
1775    let mut lwork = 0;
1776    unsafe {
1777        try_ffi!(sys::cusolverDnZheevjBatched_bufferSize(
1778            ctx.as_raw(),
1779            mode.into(),
1780            fill_mode.into(),
1781            to_i32(n, "n")?,
1782            a.as_ptr().cast(),
1783            to_i32(lda, "lda")?,
1784            w.as_ptr().cast(),
1785            &raw mut lwork,
1786            params.as_raw(),
1787            to_i32(batch_count, "batch_count")?,
1788        ))?;
1789    }
1790    to_usize(lwork, "lwork")
1791}
1792
1793/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
1794///
1795/// The S and D data types are real valued single and double precision, respectively.
1796///
1797/// The C and Z data types are complex valued single and double precision, respectively.
1798///
1799/// Computes eigenvalues and eigenvectors of a sequence of symmetric (Hermitian) $n \times n$ matrices
1800///
1801/// where $\Lambda\_{j}$ is a real $n \times n$ diagonal matrix. $Q\_j$ is an $n \times n$ unitary matrix.
1802/// The diagonal elements of $\Lambda\_j$ are the eigenvalues of $A\_j$ in either ascending order or non-sorting order.
1803///
1804/// `syevj_batched` performs `syevj` on each matrix.
1805/// It requires that all matrices are of the same size `n` and are packed contiguously,
1806///
1807/// Each matrix is column-major with leading dimension `lda`, so the formula for random access is $A\_{k}\operatorname{(i,j)} = {A\lbrack\ i\ +\ lda\cdot j\ +\ lda\cdot n\cdot k\rbrack}$.
1808///
1809/// The `W` parameter also contains the eigenvalues of each matrix contiguously,
1810///
1811/// The formula for random access of `W` is $W\_{k}\operatorname{(j)} = {W\lbrack\ j\ +\ n\cdot k\rbrack}$.
1812///
1813/// Except for tolerance and maximum sweeps, `syevj_batched` can either sort the eigenvalues in ascending order (default) or choose as-is (without sorting) with [`SyevjInfo::set_sort_eigenvalues`].
1814/// If several tiny matrices are packed into diagonal blocks of one matrix, the non-sorting option can separate the spectra of those tiny matrices.
1815///
1816/// `syevj_batched` cannot report residual and executed sweeps through [`SyevjInfo::residual`] and [`SyevjInfo::executed_sweeps`].
1817/// Calling either accessor returns [`Status::NotSupported`].
1818/// Compute the residual explicitly when needed.
1819///
1820/// Provide workspace through `workspace`.
1821/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
1822/// The workspace size in bytes is `size_of::<T>() * lwork`.
1823///
1824/// `dev_info` has one entry per batch item.
1825/// If the call returns [`Status::InvalidValue`], `dev_info[0] == -i` indicates that the `i`th parameter is invalid.
1826/// Otherwise, `dev_info[i] == n + 1` indicates that `syevj_batched` did not converge on the `i`th matrix within the given tolerance and maximum sweep count.
1827///
1828/// If `mode` is [`EigenMode::Vector`], $A\_j$ contains the orthonormal eigenvectors $V\_j$.
1829///
1830/// # Errors
1831///
1832/// Returns an error if cuSOLVER has not been initialized, if the
1833/// matrix dimensions, leading dimension, eigen mode, fill mode, or batch size
1834/// are invalid, or if cuSOLVER reports an internal failure.
1835pub fn ssyevj_batched(
1836    ctx: &Context,
1837    mode: EigenMode,
1838    fill_mode: FillMode,
1839    n: usize,
1840    a: &mut DeviceMemory<f32>,
1841    lda: usize,
1842    w: &mut DeviceMemory<f32>,
1843    workspace: &mut DeviceMemory<f32>,
1844    dev_info: &mut DeviceMemory<i32>,
1845    params: &SyevjInfo,
1846    batch_count: usize,
1847) -> Result<()> {
1848    ctx.bind()?;
1849    validate_syevj_batched_buffers(n, a.len(), lda, w.len(), batch_count)?;
1850    require_info_buffer_len(dev_info, batch_count)?;
1851    let lwork =
1852        ssyevj_batched_buffer_size(ctx, mode, fill_mode, n, a, lda, w, params, batch_count)?;
1853    require_workspace(workspace.len(), lwork)?;
1854    unsafe {
1855        try_ffi!(sys::cusolverDnSsyevjBatched(
1856            ctx.as_raw(),
1857            mode.into(),
1858            fill_mode.into(),
1859            to_i32(n, "n")?,
1860            a.as_mut_ptr().cast(),
1861            to_i32(lda, "lda")?,
1862            w.as_mut_ptr().cast(),
1863            workspace.as_mut_ptr().cast(),
1864            to_i32(lwork, "lwork")?,
1865            dev_info.as_mut_ptr().cast(),
1866            params.as_raw(),
1867            to_i32(batch_count, "batch_count")?,
1868        ))?;
1869    }
1870    Ok(())
1871}
1872
1873/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
1874///
1875/// The S and D data types are real valued single and double precision, respectively.
1876///
1877/// The C and Z data types are complex valued single and double precision, respectively.
1878///
1879/// Computes eigenvalues and eigenvectors of a sequence of symmetric (Hermitian) $n \times n$ matrices
1880///
1881/// where $\Lambda\_{j}$ is a real $n \times n$ diagonal matrix. $Q\_j$ is an $n \times n$ unitary matrix.
1882/// The diagonal elements of $\Lambda\_j$ are the eigenvalues of $A\_j$ in either ascending order or non-sorting order.
1883///
1884/// `syevj_batched` performs `syevj` on each matrix.
1885/// It requires that all matrices are of the same size `n` and are packed contiguously,
1886///
1887/// Each matrix is column-major with leading dimension `lda`, so the formula for random access is $A\_{k}\operatorname{(i,j)} = {A\lbrack\ i\ +\ lda\cdot j\ +\ lda\cdot n\cdot k\rbrack}$.
1888///
1889/// The `W` parameter also contains the eigenvalues of each matrix contiguously,
1890///
1891/// The formula for random access of `W` is $W\_{k}\operatorname{(j)} = {W\lbrack\ j\ +\ n\cdot k\rbrack}$.
1892///
1893/// Except for tolerance and maximum sweeps, `syevj_batched` can either sort the eigenvalues in ascending order (default) or choose as-is (without sorting) with [`SyevjInfo::set_sort_eigenvalues`].
1894/// If several tiny matrices are packed into diagonal blocks of one matrix, the non-sorting option can separate the spectra of those tiny matrices.
1895///
1896/// `syevj_batched` cannot report residual and executed sweeps through [`SyevjInfo::residual`] and [`SyevjInfo::executed_sweeps`].
1897/// Calling either accessor returns [`Status::NotSupported`].
1898/// Compute the residual explicitly when needed.
1899///
1900/// Provide workspace through `workspace`.
1901/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
1902/// The workspace size in bytes is `size_of::<T>() * lwork`.
1903///
1904/// `dev_info` has one entry per batch item.
1905/// If the call returns [`Status::InvalidValue`], `dev_info[0] == -i` indicates that the `i`th parameter is invalid.
1906/// Otherwise, `dev_info[i] == n + 1` indicates that `syevj_batched` did not converge on the `i`th matrix within the given tolerance and maximum sweep count.
1907///
1908/// If `mode` is [`EigenMode::Vector`], $A\_j$ contains the orthonormal eigenvectors $V\_j$.
1909///
1910/// # Errors
1911///
1912/// Returns an error if cuSOLVER has not been initialized, if the
1913/// matrix dimensions, leading dimension, eigen mode, fill mode, or batch size
1914/// are invalid, or if cuSOLVER reports an internal failure.
1915pub fn dsyevj_batched(
1916    ctx: &Context,
1917    mode: EigenMode,
1918    fill_mode: FillMode,
1919    n: usize,
1920    a: &mut DeviceMemory<f64>,
1921    lda: usize,
1922    w: &mut DeviceMemory<f64>,
1923    workspace: &mut DeviceMemory<f64>,
1924    dev_info: &mut DeviceMemory<i32>,
1925    params: &SyevjInfo,
1926    batch_count: usize,
1927) -> Result<()> {
1928    ctx.bind()?;
1929    validate_syevj_batched_buffers(n, a.len(), lda, w.len(), batch_count)?;
1930    require_info_buffer_len(dev_info, batch_count)?;
1931    let lwork =
1932        dsyevj_batched_buffer_size(ctx, mode, fill_mode, n, a, lda, w, params, batch_count)?;
1933    require_workspace(workspace.len(), lwork)?;
1934    unsafe {
1935        try_ffi!(sys::cusolverDnDsyevjBatched(
1936            ctx.as_raw(),
1937            mode.into(),
1938            fill_mode.into(),
1939            to_i32(n, "n")?,
1940            a.as_mut_ptr().cast(),
1941            to_i32(lda, "lda")?,
1942            w.as_mut_ptr().cast(),
1943            workspace.as_mut_ptr().cast(),
1944            to_i32(lwork, "lwork")?,
1945            dev_info.as_mut_ptr().cast(),
1946            params.as_raw(),
1947            to_i32(batch_count, "batch_count")?,
1948        ))?;
1949    }
1950    Ok(())
1951}
1952
1953pub fn cheevj_batched(
1954    ctx: &Context,
1955    mode: EigenMode,
1956    fill_mode: FillMode,
1957    n: usize,
1958    a: &mut DeviceMemory<Complex32>,
1959    lda: usize,
1960    w: &mut DeviceMemory<f32>,
1961    workspace: &mut DeviceMemory<Complex32>,
1962    dev_info: &mut DeviceMemory<i32>,
1963    params: &SyevjInfo,
1964    batch_count: usize,
1965) -> Result<()> {
1966    ctx.bind()?;
1967    validate_syevj_batched_buffers(n, a.len(), lda, w.len(), batch_count)?;
1968    require_info_buffer_len(dev_info, batch_count)?;
1969    let lwork =
1970        cheevj_batched_buffer_size(ctx, mode, fill_mode, n, a, lda, w, params, batch_count)?;
1971    require_workspace(workspace.len(), lwork)?;
1972    unsafe {
1973        try_ffi!(sys::cusolverDnCheevjBatched(
1974            ctx.as_raw(),
1975            mode.into(),
1976            fill_mode.into(),
1977            to_i32(n, "n")?,
1978            a.as_mut_ptr().cast(),
1979            to_i32(lda, "lda")?,
1980            w.as_mut_ptr().cast(),
1981            workspace.as_mut_ptr().cast(),
1982            to_i32(lwork, "lwork")?,
1983            dev_info.as_mut_ptr().cast(),
1984            params.as_raw(),
1985            to_i32(batch_count, "batch_count")?,
1986        ))?;
1987    }
1988    Ok(())
1989}
1990
1991pub fn zheevj_batched(
1992    ctx: &Context,
1993    mode: EigenMode,
1994    fill_mode: FillMode,
1995    n: usize,
1996    a: &mut DeviceMemory<Complex64>,
1997    lda: usize,
1998    w: &mut DeviceMemory<f64>,
1999    workspace: &mut DeviceMemory<Complex64>,
2000    dev_info: &mut DeviceMemory<i32>,
2001    params: &SyevjInfo,
2002    batch_count: usize,
2003) -> Result<()> {
2004    ctx.bind()?;
2005    validate_syevj_batched_buffers(n, a.len(), lda, w.len(), batch_count)?;
2006    require_info_buffer_len(dev_info, batch_count)?;
2007    let lwork =
2008        zheevj_batched_buffer_size(ctx, mode, fill_mode, n, a, lda, w, params, batch_count)?;
2009    require_workspace(workspace.len(), lwork)?;
2010    unsafe {
2011        try_ffi!(sys::cusolverDnZheevjBatched(
2012            ctx.as_raw(),
2013            mode.into(),
2014            fill_mode.into(),
2015            to_i32(n, "n")?,
2016            a.as_mut_ptr().cast(),
2017            to_i32(lda, "lda")?,
2018            w.as_mut_ptr().cast(),
2019            workspace.as_mut_ptr().cast(),
2020            to_i32(lwork, "lwork")?,
2021            dev_info.as_mut_ptr().cast(),
2022            params.as_raw(),
2023            to_i32(batch_count, "batch_count")?,
2024        ))?;
2025    }
2026    Ok(())
2027}
2028
2029pub fn ssygvj_buffer_size(
2030    ctx: &Context,
2031    eig_type: EigenType,
2032    mode: EigenMode,
2033    fill_mode: FillMode,
2034    n: usize,
2035    a: &DeviceMemory<f32>,
2036    lda: usize,
2037    b: &DeviceMemory<f32>,
2038    ldb: usize,
2039    w: &DeviceMemory<f32>,
2040    params: &SyevjInfo,
2041) -> Result<usize> {
2042    ctx.bind()?;
2043    validate_sygvj_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2044    let mut lwork = 0;
2045    unsafe {
2046        try_ffi!(sys::cusolverDnSsygvj_bufferSize(
2047            ctx.as_raw(),
2048            eig_type.into(),
2049            mode.into(),
2050            fill_mode.into(),
2051            to_i32(n, "n")?,
2052            a.as_ptr().cast(),
2053            to_i32(lda, "lda")?,
2054            b.as_ptr().cast(),
2055            to_i32(ldb, "ldb")?,
2056            w.as_ptr().cast(),
2057            &raw mut lwork,
2058            params.as_raw(),
2059        ))?;
2060    }
2061    to_usize(lwork, "lwork")
2062}
2063
2064pub fn dsygvj_buffer_size(
2065    ctx: &Context,
2066    eig_type: EigenType,
2067    mode: EigenMode,
2068    fill_mode: FillMode,
2069    n: usize,
2070    a: &DeviceMemory<f64>,
2071    lda: usize,
2072    b: &DeviceMemory<f64>,
2073    ldb: usize,
2074    w: &DeviceMemory<f64>,
2075    params: &SyevjInfo,
2076) -> Result<usize> {
2077    ctx.bind()?;
2078    validate_sygvj_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2079    let mut lwork = 0;
2080    unsafe {
2081        try_ffi!(sys::cusolverDnDsygvj_bufferSize(
2082            ctx.as_raw(),
2083            eig_type.into(),
2084            mode.into(),
2085            fill_mode.into(),
2086            to_i32(n, "n")?,
2087            a.as_ptr().cast(),
2088            to_i32(lda, "lda")?,
2089            b.as_ptr().cast(),
2090            to_i32(ldb, "ldb")?,
2091            w.as_ptr().cast(),
2092            &raw mut lwork,
2093            params.as_raw(),
2094        ))?;
2095    }
2096    to_usize(lwork, "lwork")
2097}
2098
2099pub fn chegvj_buffer_size(
2100    ctx: &Context,
2101    eig_type: EigenType,
2102    mode: EigenMode,
2103    fill_mode: FillMode,
2104    n: usize,
2105    a: &DeviceMemory<Complex32>,
2106    lda: usize,
2107    b: &DeviceMemory<Complex32>,
2108    ldb: usize,
2109    w: &DeviceMemory<f32>,
2110    params: &SyevjInfo,
2111) -> Result<usize> {
2112    ctx.bind()?;
2113    validate_sygvj_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2114    let mut lwork = 0;
2115    unsafe {
2116        try_ffi!(sys::cusolverDnChegvj_bufferSize(
2117            ctx.as_raw(),
2118            eig_type.into(),
2119            mode.into(),
2120            fill_mode.into(),
2121            to_i32(n, "n")?,
2122            a.as_ptr().cast(),
2123            to_i32(lda, "lda")?,
2124            b.as_ptr().cast(),
2125            to_i32(ldb, "ldb")?,
2126            w.as_ptr().cast(),
2127            &raw mut lwork,
2128            params.as_raw(),
2129        ))?;
2130    }
2131    to_usize(lwork, "lwork")
2132}
2133
2134pub fn zhegvj_buffer_size(
2135    ctx: &Context,
2136    eig_type: EigenType,
2137    mode: EigenMode,
2138    fill_mode: FillMode,
2139    n: usize,
2140    a: &DeviceMemory<Complex64>,
2141    lda: usize,
2142    b: &DeviceMemory<Complex64>,
2143    ldb: usize,
2144    w: &DeviceMemory<f64>,
2145    params: &SyevjInfo,
2146) -> Result<usize> {
2147    ctx.bind()?;
2148    validate_sygvj_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2149    let mut lwork = 0;
2150    unsafe {
2151        try_ffi!(sys::cusolverDnZhegvj_bufferSize(
2152            ctx.as_raw(),
2153            eig_type.into(),
2154            mode.into(),
2155            fill_mode.into(),
2156            to_i32(n, "n")?,
2157            a.as_ptr().cast(),
2158            to_i32(lda, "lda")?,
2159            b.as_ptr().cast(),
2160            to_i32(ldb, "ldb")?,
2161            w.as_ptr().cast(),
2162            &raw mut lwork,
2163            params.as_raw(),
2164        ))?;
2165    }
2166    to_usize(lwork, "lwork")
2167}
2168
2169/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
2170///
2171/// The S and D data types are real valued single and double precision, respectively.
2172///
2173/// The C and Z data types are complex valued single and double precision, respectively.
2174///
2175/// Computes eigenvalues and eigenvectors of a symmetric (Hermitian) $n \times n$ matrix-pair (`A`,`B`).
2176/// The generalized symmetric-definite eigenvalue problem depends on
2177/// [`EigenType`].
2178///
2179/// where the matrix `B` is positive definite.
2180/// `Λ` is a real $n \times n$ diagonal matrix.
2181/// The diagonal elements of `Λ` are the eigenvalues of (`A`, `B`) in ascending order.
2182/// `V` is an $n \times n$ orthogonal matrix.
2183/// The eigenvectors are normalized according to the selected generalized
2184/// eigenvalue problem.
2185///
2186/// This computes the same generalized symmetric eigenvalue problem as `sygvd`, but `sygvj` uses `syevj` where `sygvd` uses `syevd`.
2187/// Therefore, `sygvj` inherits properties of `syevj`; use [`SyevjInfo::set_tolerance`] and [`SyevjInfo::set_max_sweeps`] to configure tolerance and maximum sweeps.
2188///
2189/// However the meaning of residual is different from `syevj`.
2190/// `sygvj` first computes Cholesky factorization of matrix `B`,
2191///
2192/// transform the problem to standard eigenvalue problem, then calls `syevj`.
2193///
2194/// For example, the standard eigenvalue problem of type I is
2195///
2196/// where matrix `M` is symmetric
2197///
2198/// The residual is the result of `syevj` on matrix `M`, not `A`.
2199///
2200/// Provide workspace through `workspace`.
2201/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
2202/// The workspace size in bytes is `size_of::<T>() * lwork`.
2203///
2204/// If the reported `info` value is `-i`, the `i`th parameter is invalid.
2205/// If `info == i` with `0 < i <= n`, `B` is not positive definite, the factorization of `B` could not be completed, and no eigenvalues or eigenvectors were computed.
2206/// If `info == n + 1`, `syevj` did not converge within the given tolerance and maximum sweep count.
2207/// In this case, the eigenvalues and eigenvectors are still computed because non-convergence comes from the configured tolerance or maximum sweep count.
2208///
2209/// If `mode` is [`EigenMode::Vector`], `A` contains the orthogonal eigenvectors `V`.
2210///
2211/// # Errors
2212///
2213/// Returns an error if cuSOLVER has not been initialized, if the
2214/// matrix dimensions, leading dimensions, generalized eigenproblem type, eigen
2215/// mode, or fill mode are invalid, or if cuSOLVER reports an internal failure.
2216pub fn ssygvj(
2217    ctx: &Context,
2218    eig_type: EigenType,
2219    mode: EigenMode,
2220    fill_mode: FillMode,
2221    n: usize,
2222    a: &mut DeviceMemory<f32>,
2223    lda: usize,
2224    b: &mut DeviceMemory<f32>,
2225    ldb: usize,
2226    w: &mut DeviceMemory<f32>,
2227    workspace: &mut DeviceMemory<f32>,
2228    dev_info: &mut DeviceMemory<i32>,
2229    params: &SyevjInfo,
2230) -> Result<()> {
2231    ctx.bind()?;
2232    validate_sygvj_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2233    require_info_buffer(dev_info)?;
2234    let lwork = ssygvj_buffer_size(ctx, eig_type, mode, fill_mode, n, a, lda, b, ldb, w, params)?;
2235    require_workspace(workspace.len(), lwork)?;
2236    unsafe {
2237        try_ffi!(sys::cusolverDnSsygvj(
2238            ctx.as_raw(),
2239            eig_type.into(),
2240            mode.into(),
2241            fill_mode.into(),
2242            to_i32(n, "n")?,
2243            a.as_mut_ptr().cast(),
2244            to_i32(lda, "lda")?,
2245            b.as_mut_ptr().cast(),
2246            to_i32(ldb, "ldb")?,
2247            w.as_mut_ptr().cast(),
2248            workspace.as_mut_ptr().cast(),
2249            to_i32(lwork, "lwork")?,
2250            dev_info.as_mut_ptr().cast(),
2251            params.as_raw(),
2252        ))?;
2253    }
2254    Ok(())
2255}
2256
2257/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
2258///
2259/// The S and D data types are real valued single and double precision, respectively.
2260///
2261/// The C and Z data types are complex valued single and double precision, respectively.
2262///
2263/// Computes eigenvalues and eigenvectors of a symmetric (Hermitian) $n \times n$ matrix-pair (`A`,`B`).
2264/// The generalized symmetric-definite eigenvalue problem depends on
2265/// [`EigenType`].
2266///
2267/// where the matrix `B` is positive definite.
2268/// `Λ` is a real $n \times n$ diagonal matrix.
2269/// The diagonal elements of `Λ` are the eigenvalues of (`A`, `B`) in ascending order.
2270/// `V` is an $n \times n$ orthogonal matrix.
2271/// The eigenvectors are normalized according to the selected generalized
2272/// eigenvalue problem.
2273///
2274/// This computes the same generalized symmetric eigenvalue problem as `sygvd`, but `sygvj` uses `syevj` where `sygvd` uses `syevd`.
2275/// Therefore, `sygvj` inherits properties of `syevj`; use [`SyevjInfo::set_tolerance`] and [`SyevjInfo::set_max_sweeps`] to configure tolerance and maximum sweeps.
2276///
2277/// However the meaning of residual is different from `syevj`.
2278/// `sygvj` first computes Cholesky factorization of matrix `B`,
2279///
2280/// transform the problem to standard eigenvalue problem, then calls `syevj`.
2281///
2282/// For example, the standard eigenvalue problem of type I is
2283///
2284/// where matrix `M` is symmetric
2285///
2286/// The residual is the result of `syevj` on matrix `M`, not `A`.
2287///
2288/// Provide workspace through `workspace`.
2289/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
2290/// The workspace size in bytes is `size_of::<T>() * lwork`.
2291///
2292/// If the reported `info` value is `-i`, the `i`th parameter is invalid.
2293/// If `info == i` with `0 < i <= n`, `B` is not positive definite, the factorization of `B` could not be completed, and no eigenvalues or eigenvectors were computed.
2294/// If `info == n + 1`, `syevj` did not converge within the given tolerance and maximum sweep count.
2295/// In this case, the eigenvalues and eigenvectors are still computed because non-convergence comes from the configured tolerance or maximum sweep count.
2296///
2297/// If `mode` is [`EigenMode::Vector`], `A` contains the orthogonal eigenvectors `V`.
2298///
2299/// # Errors
2300///
2301/// Returns an error if cuSOLVER has not been initialized, if the
2302/// matrix dimensions, leading dimensions, generalized eigenproblem type, eigen
2303/// mode, or fill mode are invalid, or if cuSOLVER reports an internal failure.
2304pub fn dsygvj(
2305    ctx: &Context,
2306    eig_type: EigenType,
2307    mode: EigenMode,
2308    fill_mode: FillMode,
2309    n: usize,
2310    a: &mut DeviceMemory<f64>,
2311    lda: usize,
2312    b: &mut DeviceMemory<f64>,
2313    ldb: usize,
2314    w: &mut DeviceMemory<f64>,
2315    workspace: &mut DeviceMemory<f64>,
2316    dev_info: &mut DeviceMemory<i32>,
2317    params: &SyevjInfo,
2318) -> Result<()> {
2319    ctx.bind()?;
2320    validate_sygvj_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2321    require_info_buffer(dev_info)?;
2322    let lwork = dsygvj_buffer_size(ctx, eig_type, mode, fill_mode, n, a, lda, b, ldb, w, params)?;
2323    require_workspace(workspace.len(), lwork)?;
2324    unsafe {
2325        try_ffi!(sys::cusolverDnDsygvj(
2326            ctx.as_raw(),
2327            eig_type.into(),
2328            mode.into(),
2329            fill_mode.into(),
2330            to_i32(n, "n")?,
2331            a.as_mut_ptr().cast(),
2332            to_i32(lda, "lda")?,
2333            b.as_mut_ptr().cast(),
2334            to_i32(ldb, "ldb")?,
2335            w.as_mut_ptr().cast(),
2336            workspace.as_mut_ptr().cast(),
2337            to_i32(lwork, "lwork")?,
2338            dev_info.as_mut_ptr().cast(),
2339            params.as_raw(),
2340        ))?;
2341    }
2342    Ok(())
2343}
2344
2345pub fn chegvj(
2346    ctx: &Context,
2347    eig_type: EigenType,
2348    mode: EigenMode,
2349    fill_mode: FillMode,
2350    n: usize,
2351    a: &mut DeviceMemory<Complex32>,
2352    lda: usize,
2353    b: &mut DeviceMemory<Complex32>,
2354    ldb: usize,
2355    w: &mut DeviceMemory<f32>,
2356    workspace: &mut DeviceMemory<Complex32>,
2357    dev_info: &mut DeviceMemory<i32>,
2358    params: &SyevjInfo,
2359) -> Result<()> {
2360    ctx.bind()?;
2361    validate_sygvj_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2362    require_info_buffer(dev_info)?;
2363    let lwork = chegvj_buffer_size(ctx, eig_type, mode, fill_mode, n, a, lda, b, ldb, w, params)?;
2364    require_workspace(workspace.len(), lwork)?;
2365    unsafe {
2366        try_ffi!(sys::cusolverDnChegvj(
2367            ctx.as_raw(),
2368            eig_type.into(),
2369            mode.into(),
2370            fill_mode.into(),
2371            to_i32(n, "n")?,
2372            a.as_mut_ptr().cast(),
2373            to_i32(lda, "lda")?,
2374            b.as_mut_ptr().cast(),
2375            to_i32(ldb, "ldb")?,
2376            w.as_mut_ptr().cast(),
2377            workspace.as_mut_ptr().cast(),
2378            to_i32(lwork, "lwork")?,
2379            dev_info.as_mut_ptr().cast(),
2380            params.as_raw(),
2381        ))?;
2382    }
2383    Ok(())
2384}
2385
2386pub fn zhegvj(
2387    ctx: &Context,
2388    eig_type: EigenType,
2389    mode: EigenMode,
2390    fill_mode: FillMode,
2391    n: usize,
2392    a: &mut DeviceMemory<Complex64>,
2393    lda: usize,
2394    b: &mut DeviceMemory<Complex64>,
2395    ldb: usize,
2396    w: &mut DeviceMemory<f64>,
2397    workspace: &mut DeviceMemory<Complex64>,
2398    dev_info: &mut DeviceMemory<i32>,
2399    params: &SyevjInfo,
2400) -> Result<()> {
2401    ctx.bind()?;
2402    validate_sygvj_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2403    require_info_buffer(dev_info)?;
2404    let lwork = zhegvj_buffer_size(ctx, eig_type, mode, fill_mode, n, a, lda, b, ldb, w, params)?;
2405    require_workspace(workspace.len(), lwork)?;
2406    unsafe {
2407        try_ffi!(sys::cusolverDnZhegvj(
2408            ctx.as_raw(),
2409            eig_type.into(),
2410            mode.into(),
2411            fill_mode.into(),
2412            to_i32(n, "n")?,
2413            a.as_mut_ptr().cast(),
2414            to_i32(lda, "lda")?,
2415            b.as_mut_ptr().cast(),
2416            to_i32(ldb, "ldb")?,
2417            w.as_mut_ptr().cast(),
2418            workspace.as_mut_ptr().cast(),
2419            to_i32(lwork, "lwork")?,
2420            dev_info.as_mut_ptr().cast(),
2421            params.as_raw(),
2422        ))?;
2423    }
2424    Ok(())
2425}
2426
2427pub fn ssygvd_buffer_size(
2428    ctx: &Context,
2429    eig_type: EigenType,
2430    mode: EigenMode,
2431    fill_mode: FillMode,
2432    n: usize,
2433    a: &DeviceMemory<f32>,
2434    lda: usize,
2435    b: &DeviceMemory<f32>,
2436    ldb: usize,
2437    w: &DeviceMemory<f32>,
2438) -> Result<usize> {
2439    ctx.bind()?;
2440    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2441    let mut lwork = 0;
2442    unsafe {
2443        try_ffi!(sys::cusolverDnSsygvd_bufferSize(
2444            ctx.as_raw(),
2445            eig_type.into(),
2446            mode.into(),
2447            fill_mode.into(),
2448            to_i32(n, "n")?,
2449            a.as_ptr().cast(),
2450            to_i32(lda, "lda")?,
2451            b.as_ptr().cast(),
2452            to_i32(ldb, "ldb")?,
2453            w.as_ptr().cast(),
2454            &raw mut lwork,
2455        ))?;
2456    }
2457    to_usize(lwork, "lwork")
2458}
2459
2460pub fn dsygvd_buffer_size(
2461    ctx: &Context,
2462    eig_type: EigenType,
2463    mode: EigenMode,
2464    fill_mode: FillMode,
2465    n: usize,
2466    a: &DeviceMemory<f64>,
2467    lda: usize,
2468    b: &DeviceMemory<f64>,
2469    ldb: usize,
2470    w: &DeviceMemory<f64>,
2471) -> Result<usize> {
2472    ctx.bind()?;
2473    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2474    let mut lwork = 0;
2475    unsafe {
2476        try_ffi!(sys::cusolverDnDsygvd_bufferSize(
2477            ctx.as_raw(),
2478            eig_type.into(),
2479            mode.into(),
2480            fill_mode.into(),
2481            to_i32(n, "n")?,
2482            a.as_ptr().cast(),
2483            to_i32(lda, "lda")?,
2484            b.as_ptr().cast(),
2485            to_i32(ldb, "ldb")?,
2486            w.as_ptr().cast(),
2487            &raw mut lwork,
2488        ))?;
2489    }
2490    to_usize(lwork, "lwork")
2491}
2492
2493pub fn chegvd_buffer_size(
2494    ctx: &Context,
2495    eig_type: EigenType,
2496    mode: EigenMode,
2497    fill_mode: FillMode,
2498    n: usize,
2499    a: &DeviceMemory<Complex32>,
2500    lda: usize,
2501    b: &DeviceMemory<Complex32>,
2502    ldb: usize,
2503    w: &DeviceMemory<f32>,
2504) -> Result<usize> {
2505    ctx.bind()?;
2506    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2507    let mut lwork = 0;
2508    unsafe {
2509        try_ffi!(sys::cusolverDnChegvd_bufferSize(
2510            ctx.as_raw(),
2511            eig_type.into(),
2512            mode.into(),
2513            fill_mode.into(),
2514            to_i32(n, "n")?,
2515            a.as_ptr().cast(),
2516            to_i32(lda, "lda")?,
2517            b.as_ptr().cast(),
2518            to_i32(ldb, "ldb")?,
2519            w.as_ptr().cast(),
2520            &raw mut lwork,
2521        ))?;
2522    }
2523    to_usize(lwork, "lwork")
2524}
2525
2526pub fn zhegvd_buffer_size(
2527    ctx: &Context,
2528    eig_type: EigenType,
2529    mode: EigenMode,
2530    fill_mode: FillMode,
2531    n: usize,
2532    a: &DeviceMemory<Complex64>,
2533    lda: usize,
2534    b: &DeviceMemory<Complex64>,
2535    ldb: usize,
2536    w: &DeviceMemory<f64>,
2537) -> Result<usize> {
2538    ctx.bind()?;
2539    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2540    let mut lwork = 0;
2541    unsafe {
2542        try_ffi!(sys::cusolverDnZhegvd_bufferSize(
2543            ctx.as_raw(),
2544            eig_type.into(),
2545            mode.into(),
2546            fill_mode.into(),
2547            to_i32(n, "n")?,
2548            a.as_ptr().cast(),
2549            to_i32(lda, "lda")?,
2550            b.as_ptr().cast(),
2551            to_i32(ldb, "ldb")?,
2552            w.as_ptr().cast(),
2553            &raw mut lwork,
2554        ))?;
2555    }
2556    to_usize(lwork, "lwork")
2557}
2558
2559/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
2560///
2561/// The S and D data types are real valued single and double precision, respectively.
2562///
2563/// The C and Z data types are complex valued single and double precision, respectively.
2564///
2565/// Computes eigenvalues and eigenvectors of a symmetric (Hermitian) $n \times n$ matrix-pair (`A`,`B`).
2566/// The generalized symmetric-definite eigenvalue problem depends on
2567/// [`EigenType`].
2568///
2569/// where the matrix `B` is positive definite.
2570/// `Λ` is a real $n \times n$ diagonal matrix.
2571/// The diagonal elements of `Λ` are the eigenvalues of (`A`, `B`) in ascending order.
2572/// `V` is an $n \times n$ orthogonal matrix.
2573/// The eigenvectors are normalized according to the selected generalized
2574/// eigenvalue problem.
2575///
2576/// Provide workspace through `workspace`.
2577/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
2578/// The workspace size in bytes is `size_of::<T>() * lwork`.
2579///
2580/// If the reported `dev_info` value is `-i`, the `i`th parameter is invalid.
2581/// If `dev_info = i` with `0 < i <= n` and `mode` is
2582/// [`EigenMode::NoVector`], `i` off-diagonal elements of an intermediate
2583/// tridiagonal form did not converge to zero.
2584/// If `dev_info = N + i` with `i > 0`, then the leading minor of order `i` of
2585/// `B` is not positive definite.
2586/// The factorization of `B` could not be completed and no eigenvalues or eigenvectors were computed.
2587///
2588/// If `mode` is [`EigenMode::Vector`], `A` contains the orthogonal eigenvectors of the matrix `A`.
2589/// The eigenvectors are computed by divide and conquer algorithm.
2590///
2591/// # Errors
2592///
2593/// Returns an error if cuSOLVER has not been initialized, if the
2594/// matrix dimensions, leading dimensions, generalized eigenproblem type, eigen
2595/// mode, or fill mode are invalid, if the current GPU architecture is
2596/// unsupported, or if cuSOLVER reports an internal failure.
2597pub fn ssygvd(
2598    ctx: &Context,
2599    eig_type: EigenType,
2600    mode: EigenMode,
2601    fill_mode: FillMode,
2602    n: usize,
2603    a: &mut DeviceMemory<f32>,
2604    lda: usize,
2605    b: &mut DeviceMemory<f32>,
2606    ldb: usize,
2607    w: &mut DeviceMemory<f32>,
2608    workspace: &mut DeviceMemory<f32>,
2609    dev_info: &mut DeviceMemory<i32>,
2610) -> Result<()> {
2611    ctx.bind()?;
2612    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2613    require_info_buffer(dev_info)?;
2614    let lwork = ssygvd_buffer_size(ctx, eig_type, mode, fill_mode, n, a, lda, b, ldb, w)?;
2615    require_workspace(workspace.len(), lwork)?;
2616    unsafe {
2617        try_ffi!(sys::cusolverDnSsygvd(
2618            ctx.as_raw(),
2619            eig_type.into(),
2620            mode.into(),
2621            fill_mode.into(),
2622            to_i32(n, "n")?,
2623            a.as_mut_ptr().cast(),
2624            to_i32(lda, "lda")?,
2625            b.as_mut_ptr().cast(),
2626            to_i32(ldb, "ldb")?,
2627            w.as_mut_ptr().cast(),
2628            workspace.as_mut_ptr().cast(),
2629            to_i32(lwork, "lwork")?,
2630            dev_info.as_mut_ptr().cast(),
2631        ))?;
2632    }
2633    Ok(())
2634}
2635
2636/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
2637///
2638/// The S and D data types are real valued single and double precision, respectively.
2639///
2640/// The C and Z data types are complex valued single and double precision, respectively.
2641///
2642/// Computes eigenvalues and eigenvectors of a symmetric (Hermitian) $n \times n$ matrix-pair (`A`,`B`).
2643/// The generalized symmetric-definite eigenvalue problem depends on
2644/// [`EigenType`].
2645///
2646/// where the matrix `B` is positive definite.
2647/// `Λ` is a real $n \times n$ diagonal matrix.
2648/// The diagonal elements of `Λ` are the eigenvalues of (`A`, `B`) in ascending order.
2649/// `V` is an $n \times n$ orthogonal matrix.
2650/// The eigenvectors are normalized according to the selected generalized
2651/// eigenvalue problem.
2652///
2653/// Provide workspace through `workspace`.
2654/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
2655/// The workspace size in bytes is `size_of::<T>() * lwork`.
2656///
2657/// If the reported `dev_info` value is `-i`, the `i`th parameter is invalid.
2658/// If `dev_info = i` with `0 < i <= n` and `mode` is
2659/// [`EigenMode::NoVector`], `i` off-diagonal elements of an intermediate
2660/// tridiagonal form did not converge to zero.
2661/// If `dev_info = N + i` with `i > 0`, then the leading minor of order `i` of
2662/// `B` is not positive definite.
2663/// The factorization of `B` could not be completed and no eigenvalues or eigenvectors were computed.
2664///
2665/// If `mode` is [`EigenMode::Vector`], `A` contains the orthogonal eigenvectors of the matrix `A`.
2666/// The eigenvectors are computed by divide and conquer algorithm.
2667///
2668/// # Errors
2669///
2670/// Returns an error if cuSOLVER has not been initialized, if the
2671/// matrix dimensions, leading dimensions, generalized eigenproblem type, eigen
2672/// mode, or fill mode are invalid, if the current GPU architecture is
2673/// unsupported, or if cuSOLVER reports an internal failure.
2674pub fn dsygvd(
2675    ctx: &Context,
2676    eig_type: EigenType,
2677    mode: EigenMode,
2678    fill_mode: FillMode,
2679    n: usize,
2680    a: &mut DeviceMemory<f64>,
2681    lda: usize,
2682    b: &mut DeviceMemory<f64>,
2683    ldb: usize,
2684    w: &mut DeviceMemory<f64>,
2685    workspace: &mut DeviceMemory<f64>,
2686    dev_info: &mut DeviceMemory<i32>,
2687) -> Result<()> {
2688    ctx.bind()?;
2689    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2690    require_info_buffer(dev_info)?;
2691    let lwork = dsygvd_buffer_size(ctx, eig_type, mode, fill_mode, n, a, lda, b, ldb, w)?;
2692    require_workspace(workspace.len(), lwork)?;
2693    unsafe {
2694        try_ffi!(sys::cusolverDnDsygvd(
2695            ctx.as_raw(),
2696            eig_type.into(),
2697            mode.into(),
2698            fill_mode.into(),
2699            to_i32(n, "n")?,
2700            a.as_mut_ptr().cast(),
2701            to_i32(lda, "lda")?,
2702            b.as_mut_ptr().cast(),
2703            to_i32(ldb, "ldb")?,
2704            w.as_mut_ptr().cast(),
2705            workspace.as_mut_ptr().cast(),
2706            to_i32(lwork, "lwork")?,
2707            dev_info.as_mut_ptr().cast(),
2708        ))?;
2709    }
2710    Ok(())
2711}
2712
2713pub fn chegvd(
2714    ctx: &Context,
2715    eig_type: EigenType,
2716    mode: EigenMode,
2717    fill_mode: FillMode,
2718    n: usize,
2719    a: &mut DeviceMemory<Complex32>,
2720    lda: usize,
2721    b: &mut DeviceMemory<Complex32>,
2722    ldb: usize,
2723    w: &mut DeviceMemory<f32>,
2724    workspace: &mut DeviceMemory<Complex32>,
2725    dev_info: &mut DeviceMemory<i32>,
2726) -> Result<()> {
2727    ctx.bind()?;
2728    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2729    require_info_buffer(dev_info)?;
2730    let lwork = chegvd_buffer_size(ctx, eig_type, mode, fill_mode, n, a, lda, b, ldb, w)?;
2731    require_workspace(workspace.len(), lwork)?;
2732    unsafe {
2733        try_ffi!(sys::cusolverDnChegvd(
2734            ctx.as_raw(),
2735            eig_type.into(),
2736            mode.into(),
2737            fill_mode.into(),
2738            to_i32(n, "n")?,
2739            a.as_mut_ptr().cast(),
2740            to_i32(lda, "lda")?,
2741            b.as_mut_ptr().cast(),
2742            to_i32(ldb, "ldb")?,
2743            w.as_mut_ptr().cast(),
2744            workspace.as_mut_ptr().cast(),
2745            to_i32(lwork, "lwork")?,
2746            dev_info.as_mut_ptr().cast(),
2747        ))?;
2748    }
2749    Ok(())
2750}
2751
2752pub fn zhegvd(
2753    ctx: &Context,
2754    eig_type: EigenType,
2755    mode: EigenMode,
2756    fill_mode: FillMode,
2757    n: usize,
2758    a: &mut DeviceMemory<Complex64>,
2759    lda: usize,
2760    b: &mut DeviceMemory<Complex64>,
2761    ldb: usize,
2762    w: &mut DeviceMemory<f64>,
2763    workspace: &mut DeviceMemory<Complex64>,
2764    dev_info: &mut DeviceMemory<i32>,
2765) -> Result<()> {
2766    ctx.bind()?;
2767    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2768    require_info_buffer(dev_info)?;
2769    let lwork = zhegvd_buffer_size(ctx, eig_type, mode, fill_mode, n, a, lda, b, ldb, w)?;
2770    require_workspace(workspace.len(), lwork)?;
2771    unsafe {
2772        try_ffi!(sys::cusolverDnZhegvd(
2773            ctx.as_raw(),
2774            eig_type.into(),
2775            mode.into(),
2776            fill_mode.into(),
2777            to_i32(n, "n")?,
2778            a.as_mut_ptr().cast(),
2779            to_i32(lda, "lda")?,
2780            b.as_mut_ptr().cast(),
2781            to_i32(ldb, "ldb")?,
2782            w.as_mut_ptr().cast(),
2783            workspace.as_mut_ptr().cast(),
2784            to_i32(lwork, "lwork")?,
2785            dev_info.as_mut_ptr().cast(),
2786        ))?;
2787    }
2788    Ok(())
2789}
2790
2791pub fn ssygvdx_selected_buffer_size(
2792    ctx: &Context,
2793    eig_type: EigenType,
2794    mode: EigenMode,
2795    fill_mode: FillMode,
2796    selection: EigenSelection<f32>,
2797    n: usize,
2798    a: &DeviceMemory<f32>,
2799    lda: usize,
2800    b: &DeviceMemory<f32>,
2801    ldb: usize,
2802    w: &DeviceMemory<f32>,
2803) -> Result<(usize, usize)> {
2804    ctx.bind()?;
2805    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2806    let (range, value_range, index_range) = selection_parts(selection);
2807    let (vl, vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
2808    let mut meig = 0;
2809    let mut lwork = 0;
2810    unsafe {
2811        try_ffi!(sys::cusolverDnSsygvdx_bufferSize(
2812            ctx.as_raw(),
2813            eig_type.into(),
2814            mode.into(),
2815            range.into(),
2816            fill_mode.into(),
2817            to_i32(n, "n")?,
2818            a.as_ptr().cast(),
2819            to_i32(lda, "lda")?,
2820            b.as_ptr().cast(),
2821            to_i32(ldb, "ldb")?,
2822            vl,
2823            vu,
2824            to_i32(il, "il")?,
2825            to_i32(iu, "iu")?,
2826            &raw mut meig,
2827            w.as_ptr().cast(),
2828            &raw mut lwork,
2829        ))?;
2830    }
2831    Ok((to_usize(meig, "meig")?, to_usize(lwork, "lwork")?))
2832}
2833
2834pub fn dsygvdx_selected_buffer_size(
2835    ctx: &Context,
2836    eig_type: EigenType,
2837    mode: EigenMode,
2838    fill_mode: FillMode,
2839    selection: EigenSelection<f64>,
2840    n: usize,
2841    a: &DeviceMemory<f64>,
2842    lda: usize,
2843    b: &DeviceMemory<f64>,
2844    ldb: usize,
2845    w: &DeviceMemory<f64>,
2846) -> Result<(usize, usize)> {
2847    ctx.bind()?;
2848    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2849    let (range, value_range, index_range) = selection_parts(selection);
2850    let (vl, vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
2851    let mut meig = 0;
2852    let mut lwork = 0;
2853    unsafe {
2854        try_ffi!(sys::cusolverDnDsygvdx_bufferSize(
2855            ctx.as_raw(),
2856            eig_type.into(),
2857            mode.into(),
2858            range.into(),
2859            fill_mode.into(),
2860            to_i32(n, "n")?,
2861            a.as_ptr().cast(),
2862            to_i32(lda, "lda")?,
2863            b.as_ptr().cast(),
2864            to_i32(ldb, "ldb")?,
2865            vl,
2866            vu,
2867            to_i32(il, "il")?,
2868            to_i32(iu, "iu")?,
2869            &raw mut meig,
2870            w.as_ptr().cast(),
2871            &raw mut lwork,
2872        ))?;
2873    }
2874    Ok((to_usize(meig, "meig")?, to_usize(lwork, "lwork")?))
2875}
2876
2877pub fn chegvdx_selected_buffer_size(
2878    ctx: &Context,
2879    eig_type: EigenType,
2880    mode: EigenMode,
2881    fill_mode: FillMode,
2882    selection: EigenSelection<f32>,
2883    n: usize,
2884    a: &DeviceMemory<Complex32>,
2885    lda: usize,
2886    b: &DeviceMemory<Complex32>,
2887    ldb: usize,
2888    w: &DeviceMemory<f32>,
2889) -> Result<(usize, usize)> {
2890    ctx.bind()?;
2891    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2892    let (range, value_range, index_range) = selection_parts(selection);
2893    let (vl, vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
2894    let mut meig = 0;
2895    let mut lwork = 0;
2896    unsafe {
2897        try_ffi!(sys::cusolverDnChegvdx_bufferSize(
2898            ctx.as_raw(),
2899            eig_type.into(),
2900            mode.into(),
2901            range.into(),
2902            fill_mode.into(),
2903            to_i32(n, "n")?,
2904            a.as_ptr().cast(),
2905            to_i32(lda, "lda")?,
2906            b.as_ptr().cast(),
2907            to_i32(ldb, "ldb")?,
2908            vl,
2909            vu,
2910            to_i32(il, "il")?,
2911            to_i32(iu, "iu")?,
2912            &raw mut meig,
2913            w.as_ptr().cast(),
2914            &raw mut lwork,
2915        ))?;
2916    }
2917    Ok((to_usize(meig, "meig")?, to_usize(lwork, "lwork")?))
2918}
2919
2920pub fn zhegvdx_selected_buffer_size(
2921    ctx: &Context,
2922    eig_type: EigenType,
2923    mode: EigenMode,
2924    fill_mode: FillMode,
2925    selection: EigenSelection<f64>,
2926    n: usize,
2927    a: &DeviceMemory<Complex64>,
2928    lda: usize,
2929    b: &DeviceMemory<Complex64>,
2930    ldb: usize,
2931    w: &DeviceMemory<f64>,
2932) -> Result<(usize, usize)> {
2933    ctx.bind()?;
2934    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
2935    let (range, value_range, index_range) = selection_parts(selection);
2936    let (vl, vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
2937    let mut meig = 0;
2938    let mut lwork = 0;
2939    unsafe {
2940        try_ffi!(sys::cusolverDnZhegvdx_bufferSize(
2941            ctx.as_raw(),
2942            eig_type.into(),
2943            mode.into(),
2944            range.into(),
2945            fill_mode.into(),
2946            to_i32(n, "n")?,
2947            a.as_ptr().cast(),
2948            to_i32(lda, "lda")?,
2949            b.as_ptr().cast(),
2950            to_i32(ldb, "ldb")?,
2951            vl,
2952            vu,
2953            to_i32(il, "il")?,
2954            to_i32(iu, "iu")?,
2955            &raw mut meig,
2956            w.as_ptr().cast(),
2957            &raw mut lwork,
2958        ))?;
2959    }
2960    Ok((to_usize(meig, "meig")?, to_usize(lwork, "lwork")?))
2961}
2962
2963/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
2964///
2965/// The S and D data types are real valued single and double precision, respectively.
2966///
2967/// The C and Z data types are complex valued single and double precision, respectively.
2968///
2969/// Computes all or selection of the eigenvalues and optionally eigenvectors of a symmetric (Hermitian) $n \times n$ matrix-pair (`A`,`B`).
2970/// The generalized symmetric-definite eigenvalue problem is selected by
2971/// `eig_type`, with positive-definite matrix `B`.
2972/// `Λ` is a real $n \times {h\_meig}$ diagonal matrix.
2973/// The diagonal elements of `Λ` are the eigenvalues of (`A`, `B`) in ascending order.
2974/// `V` is an $n \times {h\_meig}$ orthogonal matrix.
2975/// `h_meig` is the number of eigenvalues/eigenvectors computed by the operation, `h_meig` is equal to `n` when the whole spectrum (for example, `range` = [`EigenRange::All`]) is requested.
2976/// For [`EigenType::Type1`] and [`EigenType::Type2`], eigenvectors are
2977/// normalized so that $V^{T} B V = I$. For [`EigenType::Type3`], they are
2978/// normalized so that $V^{T} B^{-1} V = I$.
2979///
2980/// Provide workspace through `workspace`.
2981/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
2982/// The workspace size in bytes is `size_of::<T>() * lwork`.
2983///
2984/// If the reported `dev_info` value is `-i`, the `i`th parameter is invalid.
2985/// If `dev_info = i` with `0 < i <= n` and `mode` is
2986/// [`EigenMode::NoVector`], `i` off-diagonal elements of an intermediate
2987/// tridiagonal form did not converge to zero.
2988/// If `dev_info = n + i` with `i > 0`, then the leading minor of order `i` of
2989/// `B` is not positive definite.
2990/// The factorization of `B` could not be completed and no eigenvalues or eigenvectors were computed.
2991///
2992/// If `mode` is [`EigenMode::Vector`], `A` contains the orthogonal eigenvectors
2993/// of the matrix `A`.
2994/// The eigenvectors are computed by divide and conquer algorithm.
2995///
2996/// # Errors
2997///
2998/// Returns an error if cuSOLVER has not been initialized, if the
2999/// matrix dimensions, leading dimensions, generalized eigenproblem type, eigen
3000/// selection, eigen mode, or fill mode are invalid, if the current GPU
3001/// architecture is unsupported, or if cuSOLVER reports an internal failure.
3002pub fn ssygvdx_selected(
3003    ctx: &Context,
3004    eig_type: EigenType,
3005    mode: EigenMode,
3006    fill_mode: FillMode,
3007    selection: EigenSelection<f32>,
3008    n: usize,
3009    a: &mut DeviceMemory<f32>,
3010    lda: usize,
3011    b: &mut DeviceMemory<f32>,
3012    ldb: usize,
3013    w: &mut DeviceMemory<f32>,
3014    workspace: &mut DeviceMemory<f32>,
3015    dev_info: &mut DeviceMemory<i32>,
3016) -> Result<usize> {
3017    ctx.bind()?;
3018    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
3019    require_info_buffer(dev_info)?;
3020    let (range, value_range, index_range) = selection_parts(selection);
3021    let (meig, lwork) = ssygvdx_selected_buffer_size(
3022        ctx, eig_type, mode, fill_mode, selection, n, a, lda, b, ldb, w,
3023    )?;
3024    require_workspace(workspace.len(), lwork)?;
3025    let (vl, vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
3026    let mut meig_raw = 0;
3027    unsafe {
3028        try_ffi!(sys::cusolverDnSsygvdx(
3029            ctx.as_raw(),
3030            eig_type.into(),
3031            mode.into(),
3032            range.into(),
3033            fill_mode.into(),
3034            to_i32(n, "n")?,
3035            a.as_mut_ptr().cast(),
3036            to_i32(lda, "lda")?,
3037            b.as_mut_ptr().cast(),
3038            to_i32(ldb, "ldb")?,
3039            vl,
3040            vu,
3041            to_i32(il, "il")?,
3042            to_i32(iu, "iu")?,
3043            &raw mut meig_raw,
3044            w.as_mut_ptr().cast(),
3045            workspace.as_mut_ptr().cast(),
3046            to_i32(lwork, "lwork")?,
3047            dev_info.as_mut_ptr().cast(),
3048        ))?;
3049    }
3050    debug_assert_eq!(meig, to_usize(meig_raw, "meig")?);
3051    Ok(meig)
3052}
3053
3054/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
3055///
3056/// The S and D data types are real valued single and double precision, respectively.
3057///
3058/// The C and Z data types are complex valued single and double precision, respectively.
3059///
3060/// Computes all or selection of the eigenvalues and optionally eigenvectors of a symmetric (Hermitian) $n \times n$ matrix-pair (`A`,`B`).
3061/// The generalized symmetric-definite eigenvalue problem is selected by
3062/// `eig_type`, with positive-definite matrix `B`.
3063/// `Λ` is a real $n \times {h\_meig}$ diagonal matrix.
3064/// The diagonal elements of `Λ` are the eigenvalues of (`A`, `B`) in ascending order.
3065/// `V` is an $n \times {h\_meig}$ orthogonal matrix.
3066/// `h_meig` is the number of eigenvalues/eigenvectors computed by the operation, `h_meig` is equal to `n` when the whole spectrum (for example, `range` = [`EigenRange::All`]) is requested.
3067/// For [`EigenType::Type1`] and [`EigenType::Type2`], eigenvectors are
3068/// normalized so that $V^{T} B V = I$. For [`EigenType::Type3`], they are
3069/// normalized so that $V^{T} B^{-1} V = I$.
3070///
3071/// Provide workspace through `workspace`.
3072/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
3073/// The workspace size in bytes is `size_of::<T>() * lwork`.
3074///
3075/// If the reported `dev_info` value is `-i`, the `i`th parameter is invalid.
3076/// If `dev_info = i` with `0 < i <= n` and `mode` is
3077/// [`EigenMode::NoVector`], `i` off-diagonal elements of an intermediate
3078/// tridiagonal form did not converge to zero.
3079/// If `dev_info = n + i` with `i > 0`, then the leading minor of order `i` of
3080/// `B` is not positive definite.
3081/// The factorization of `B` could not be completed and no eigenvalues or eigenvectors were computed.
3082///
3083/// If `mode` is [`EigenMode::Vector`], `A` contains the orthogonal eigenvectors
3084/// of the matrix `A`.
3085/// The eigenvectors are computed by divide and conquer algorithm.
3086///
3087/// # Errors
3088///
3089/// Returns an error if cuSOLVER has not been initialized, if the
3090/// matrix dimensions, leading dimensions, generalized eigenproblem type, eigen
3091/// selection, eigen mode, or fill mode are invalid, if the current GPU
3092/// architecture is unsupported, or if cuSOLVER reports an internal failure.
3093pub fn dsygvdx_selected(
3094    ctx: &Context,
3095    eig_type: EigenType,
3096    mode: EigenMode,
3097    fill_mode: FillMode,
3098    selection: EigenSelection<f64>,
3099    n: usize,
3100    a: &mut DeviceMemory<f64>,
3101    lda: usize,
3102    b: &mut DeviceMemory<f64>,
3103    ldb: usize,
3104    w: &mut DeviceMemory<f64>,
3105    workspace: &mut DeviceMemory<f64>,
3106    dev_info: &mut DeviceMemory<i32>,
3107) -> Result<usize> {
3108    ctx.bind()?;
3109    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
3110    require_info_buffer(dev_info)?;
3111    let (range, value_range, index_range) = selection_parts(selection);
3112    let (meig, lwork) = dsygvdx_selected_buffer_size(
3113        ctx, eig_type, mode, fill_mode, selection, n, a, lda, b, ldb, w,
3114    )?;
3115    require_workspace(workspace.len(), lwork)?;
3116    let (vl, vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
3117    let mut meig_raw = 0;
3118    unsafe {
3119        try_ffi!(sys::cusolverDnDsygvdx(
3120            ctx.as_raw(),
3121            eig_type.into(),
3122            mode.into(),
3123            range.into(),
3124            fill_mode.into(),
3125            to_i32(n, "n")?,
3126            a.as_mut_ptr().cast(),
3127            to_i32(lda, "lda")?,
3128            b.as_mut_ptr().cast(),
3129            to_i32(ldb, "ldb")?,
3130            vl,
3131            vu,
3132            to_i32(il, "il")?,
3133            to_i32(iu, "iu")?,
3134            &raw mut meig_raw,
3135            w.as_mut_ptr().cast(),
3136            workspace.as_mut_ptr().cast(),
3137            to_i32(lwork, "lwork")?,
3138            dev_info.as_mut_ptr().cast(),
3139        ))?;
3140    }
3141    debug_assert_eq!(meig, to_usize(meig_raw, "meig")?);
3142    Ok(meig)
3143}
3144
3145pub fn chegvdx_selected(
3146    ctx: &Context,
3147    eig_type: EigenType,
3148    mode: EigenMode,
3149    fill_mode: FillMode,
3150    selection: EigenSelection<f32>,
3151    n: usize,
3152    a: &mut DeviceMemory<Complex32>,
3153    lda: usize,
3154    b: &mut DeviceMemory<Complex32>,
3155    ldb: usize,
3156    w: &mut DeviceMemory<f32>,
3157    workspace: &mut DeviceMemory<Complex32>,
3158    dev_info: &mut DeviceMemory<i32>,
3159) -> Result<usize> {
3160    ctx.bind()?;
3161    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
3162    require_info_buffer(dev_info)?;
3163    let (range, value_range, index_range) = selection_parts(selection);
3164    let (meig, lwork) = chegvdx_selected_buffer_size(
3165        ctx, eig_type, mode, fill_mode, selection, n, a, lda, b, ldb, w,
3166    )?;
3167    require_workspace(workspace.len(), lwork)?;
3168    let (vl, vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
3169    let mut meig_raw = 0;
3170    unsafe {
3171        try_ffi!(sys::cusolverDnChegvdx(
3172            ctx.as_raw(),
3173            eig_type.into(),
3174            mode.into(),
3175            range.into(),
3176            fill_mode.into(),
3177            to_i32(n, "n")?,
3178            a.as_mut_ptr().cast(),
3179            to_i32(lda, "lda")?,
3180            b.as_mut_ptr().cast(),
3181            to_i32(ldb, "ldb")?,
3182            vl,
3183            vu,
3184            to_i32(il, "il")?,
3185            to_i32(iu, "iu")?,
3186            &raw mut meig_raw,
3187            w.as_mut_ptr().cast(),
3188            workspace.as_mut_ptr().cast(),
3189            to_i32(lwork, "lwork")?,
3190            dev_info.as_mut_ptr().cast(),
3191        ))?;
3192    }
3193    debug_assert_eq!(meig, to_usize(meig_raw, "meig")?);
3194    Ok(meig)
3195}
3196
3197pub fn zhegvdx_selected(
3198    ctx: &Context,
3199    eig_type: EigenType,
3200    mode: EigenMode,
3201    fill_mode: FillMode,
3202    selection: EigenSelection<f64>,
3203    n: usize,
3204    a: &mut DeviceMemory<Complex64>,
3205    lda: usize,
3206    b: &mut DeviceMemory<Complex64>,
3207    ldb: usize,
3208    w: &mut DeviceMemory<f64>,
3209    workspace: &mut DeviceMemory<Complex64>,
3210    dev_info: &mut DeviceMemory<i32>,
3211) -> Result<usize> {
3212    ctx.bind()?;
3213    validate_sygvd_buffers(n, a.len(), lda, b.len(), ldb, w.len())?;
3214    require_info_buffer(dev_info)?;
3215    let (range, value_range, index_range) = selection_parts(selection);
3216    let (meig, lwork) = zhegvdx_selected_buffer_size(
3217        ctx, eig_type, mode, fill_mode, selection, n, a, lda, b, ldb, w,
3218    )?;
3219    require_workspace(workspace.len(), lwork)?;
3220    let (vl, vu, il, iu) = validate_xsyevdx_range(range, n, value_range, index_range)?;
3221    let mut meig_raw = 0;
3222    unsafe {
3223        try_ffi!(sys::cusolverDnZhegvdx(
3224            ctx.as_raw(),
3225            eig_type.into(),
3226            mode.into(),
3227            range.into(),
3228            fill_mode.into(),
3229            to_i32(n, "n")?,
3230            a.as_mut_ptr().cast(),
3231            to_i32(lda, "lda")?,
3232            b.as_mut_ptr().cast(),
3233            to_i32(ldb, "ldb")?,
3234            vl,
3235            vu,
3236            to_i32(il, "il")?,
3237            to_i32(iu, "iu")?,
3238            &raw mut meig_raw,
3239            w.as_mut_ptr().cast(),
3240            workspace.as_mut_ptr().cast(),
3241            to_i32(lwork, "lwork")?,
3242            dev_info.as_mut_ptr().cast(),
3243        ))?;
3244    }
3245    debug_assert_eq!(meig, to_usize(meig_raw, "meig")?);
3246    Ok(meig)
3247}
3248
3249pub fn ssyevdx_buffer_size(
3250    ctx: &Context,
3251    mode: EigenMode,
3252    range: EigenRange,
3253    fill_mode: FillMode,
3254    n: usize,
3255    a: &DeviceMemory<f32>,
3256    lda: usize,
3257    value_range: (f32, f32),
3258    index_range: (usize, usize),
3259    w: &DeviceMemory<f32>,
3260) -> Result<(usize, usize)> {
3261    ctx.bind()?;
3262    validate_syev_buffers(n, a.len(), lda, w.len())?;
3263    let mut meig = 0;
3264    let mut lwork = 0;
3265    unsafe {
3266        try_ffi!(sys::cusolverDnSsyevdx_bufferSize(
3267            ctx.as_raw(),
3268            mode.into(),
3269            range.into(),
3270            fill_mode.into(),
3271            to_i32(n, "n")?,
3272            a.as_ptr().cast(),
3273            to_i32(lda, "lda")?,
3274            value_range.0,
3275            value_range.1,
3276            to_i32(index_range.0, "il")?,
3277            to_i32(index_range.1, "iu")?,
3278            &raw mut meig,
3279            w.as_ptr().cast(),
3280            &raw mut lwork,
3281        ))?;
3282    }
3283    Ok((to_usize(meig, "meig")?, to_usize(lwork, "lwork")?))
3284}
3285
3286pub fn dsyevdx_buffer_size(
3287    ctx: &Context,
3288    mode: EigenMode,
3289    range: EigenRange,
3290    fill_mode: FillMode,
3291    n: usize,
3292    a: &DeviceMemory<f64>,
3293    lda: usize,
3294    value_range: (f64, f64),
3295    index_range: (usize, usize),
3296    w: &DeviceMemory<f64>,
3297) -> Result<(usize, usize)> {
3298    ctx.bind()?;
3299    validate_syev_buffers(n, a.len(), lda, w.len())?;
3300    let mut meig = 0;
3301    let mut lwork = 0;
3302    unsafe {
3303        try_ffi!(sys::cusolverDnDsyevdx_bufferSize(
3304            ctx.as_raw(),
3305            mode.into(),
3306            range.into(),
3307            fill_mode.into(),
3308            to_i32(n, "n")?,
3309            a.as_ptr().cast(),
3310            to_i32(lda, "lda")?,
3311            value_range.0,
3312            value_range.1,
3313            to_i32(index_range.0, "il")?,
3314            to_i32(index_range.1, "iu")?,
3315            &raw mut meig,
3316            w.as_ptr().cast(),
3317            &raw mut lwork,
3318        ))?;
3319    }
3320    Ok((to_usize(meig, "meig")?, to_usize(lwork, "lwork")?))
3321}
3322
3323pub fn cheevdx_buffer_size(
3324    ctx: &Context,
3325    mode: EigenMode,
3326    range: EigenRange,
3327    fill_mode: FillMode,
3328    n: usize,
3329    a: &DeviceMemory<Complex32>,
3330    lda: usize,
3331    value_range: (f32, f32),
3332    index_range: (usize, usize),
3333    w: &DeviceMemory<f32>,
3334) -> Result<(usize, usize)> {
3335    ctx.bind()?;
3336    validate_syev_buffers(n, a.len(), lda, w.len())?;
3337    let mut meig = 0;
3338    let mut lwork = 0;
3339    unsafe {
3340        try_ffi!(sys::cusolverDnCheevdx_bufferSize(
3341            ctx.as_raw(),
3342            mode.into(),
3343            range.into(),
3344            fill_mode.into(),
3345            to_i32(n, "n")?,
3346            a.as_ptr().cast(),
3347            to_i32(lda, "lda")?,
3348            value_range.0,
3349            value_range.1,
3350            to_i32(index_range.0, "il")?,
3351            to_i32(index_range.1, "iu")?,
3352            &raw mut meig,
3353            w.as_ptr().cast(),
3354            &raw mut lwork,
3355        ))?;
3356    }
3357    Ok((to_usize(meig, "meig")?, to_usize(lwork, "lwork")?))
3358}
3359
3360pub fn zheevdx_buffer_size(
3361    ctx: &Context,
3362    mode: EigenMode,
3363    range: EigenRange,
3364    fill_mode: FillMode,
3365    n: usize,
3366    a: &DeviceMemory<Complex64>,
3367    lda: usize,
3368    value_range: (f64, f64),
3369    index_range: (usize, usize),
3370    w: &DeviceMemory<f64>,
3371) -> Result<(usize, usize)> {
3372    ctx.bind()?;
3373    validate_syev_buffers(n, a.len(), lda, w.len())?;
3374    let mut meig = 0;
3375    let mut lwork = 0;
3376    unsafe {
3377        try_ffi!(sys::cusolverDnZheevdx_bufferSize(
3378            ctx.as_raw(),
3379            mode.into(),
3380            range.into(),
3381            fill_mode.into(),
3382            to_i32(n, "n")?,
3383            a.as_ptr().cast(),
3384            to_i32(lda, "lda")?,
3385            value_range.0,
3386            value_range.1,
3387            to_i32(index_range.0, "il")?,
3388            to_i32(index_range.1, "iu")?,
3389            &raw mut meig,
3390            w.as_ptr().cast(),
3391            &raw mut lwork,
3392        ))?;
3393    }
3394    Ok((to_usize(meig, "meig")?, to_usize(lwork, "lwork")?))
3395}
3396
3397/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
3398///
3399/// The S and D data types are real valued single and double precision, respectively.
3400///
3401/// The C and Z data types are complex valued single and double precision, respectively.
3402///
3403/// Computes all or selection of the eigenvalues and optionally eigenvectors of a symmetric (Hermitian) $n \times n$ matrix `A`.
3404/// The standard symmetric eigenvalue problem is $A V = V \Lambda$, where `Λ`
3405/// is a real `n×h_meig` diagonal matrix.
3406/// `V` is an `n×h_meig` unitary matrix.
3407/// `h_meig` is the number of eigenvalues/eigenvectors computed by the operation, `h_meig` is equal to `n` when the whole spectrum is requested, for example with `range` = [`EigenRange::All`].
3408/// The diagonal elements of `Λ` are the eigenvalues of `A` in ascending order.
3409///
3410/// Provide workspace through `workspace`.
3411/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
3412/// The workspace size in bytes is `size_of::<T>() * lwork`.
3413///
3414/// If the reported `dev_info` value is `-i`, the `i`th parameter is invalid.
3415/// If `dev_info = i` (greater than zero), `i` off-diagonal elements of an intermediate tridiagonal form did not converge to zero.
3416///
3417/// If `mode` is [`EigenMode::Vector`], `A` contains the orthonormal
3418/// eigenvectors of the matrix `A`.
3419/// The eigenvectors are computed by a divide and conquer algorithm.
3420///
3421/// # Errors
3422///
3423/// Returns an error if cuSOLVER has not been initialized, if the
3424/// matrix dimensions, leading dimension, eigen selection, eigen mode, or fill
3425/// mode are invalid, if the current GPU architecture is unsupported, or if
3426/// cuSOLVER reports an internal failure.
3427pub fn ssyevdx(
3428    ctx: &Context,
3429    mode: EigenMode,
3430    range: EigenRange,
3431    fill_mode: FillMode,
3432    n: usize,
3433    a: &mut DeviceMemory<f32>,
3434    lda: usize,
3435    value_range: (f32, f32),
3436    index_range: (usize, usize),
3437    w: &mut DeviceMemory<f32>,
3438    workspace: &mut DeviceMemory<f32>,
3439    dev_info: &mut DeviceMemory<i32>,
3440) -> Result<usize> {
3441    ctx.bind()?;
3442    validate_syev_buffers(n, a.len(), lda, w.len())?;
3443    require_info_buffer(dev_info)?;
3444    let (meig, lwork) = ssyevdx_buffer_size(
3445        ctx,
3446        mode,
3447        range,
3448        fill_mode,
3449        n,
3450        a,
3451        lda,
3452        value_range,
3453        index_range,
3454        w,
3455    )?;
3456    require_workspace(workspace.len(), lwork)?;
3457    let mut meig_raw = 0;
3458    unsafe {
3459        try_ffi!(sys::cusolverDnSsyevdx(
3460            ctx.as_raw(),
3461            mode.into(),
3462            range.into(),
3463            fill_mode.into(),
3464            to_i32(n, "n")?,
3465            a.as_mut_ptr().cast(),
3466            to_i32(lda, "lda")?,
3467            value_range.0,
3468            value_range.1,
3469            to_i32(index_range.0, "il")?,
3470            to_i32(index_range.1, "iu")?,
3471            &raw mut meig_raw,
3472            w.as_mut_ptr().cast(),
3473            workspace.as_mut_ptr().cast(),
3474            to_i32(lwork, "lwork")?,
3475            dev_info.as_mut_ptr().cast(),
3476        ))?;
3477    }
3478    debug_assert_eq!(meig, to_usize(meig_raw, "meig")?);
3479    Ok(meig)
3480}
3481
3482/// Use the matching buffer-size helper to calculate the sizes needed for pre-allocated workspace.
3483///
3484/// The S and D data types are real valued single and double precision, respectively.
3485///
3486/// The C and Z data types are complex valued single and double precision, respectively.
3487///
3488/// Computes all or selection of the eigenvalues and optionally eigenvectors of a symmetric (Hermitian) $n \times n$ matrix `A`.
3489/// The standard symmetric eigenvalue problem is $A V = V \Lambda$, where `Λ`
3490/// is a real `n×h_meig` diagonal matrix.
3491/// `V` is an `n×h_meig` unitary matrix.
3492/// `h_meig` is the number of eigenvalues/eigenvectors computed by the operation, `h_meig` is equal to `n` when the whole spectrum is requested, for example with `range` = [`EigenRange::All`].
3493/// The diagonal elements of `Λ` are the eigenvalues of `A` in ascending order.
3494///
3495/// Provide workspace through `workspace`.
3496/// Use the corresponding `*_buffer_size` helper to query the required workspace length.
3497/// The workspace size in bytes is `size_of::<T>() * lwork`.
3498///
3499/// If the reported `dev_info` value is `-i`, the `i`th parameter is invalid.
3500/// If `dev_info = i` (greater than zero), `i` off-diagonal elements of an intermediate tridiagonal form did not converge to zero.
3501///
3502/// If `mode` is [`EigenMode::Vector`], `A` contains the orthonormal
3503/// eigenvectors of the matrix `A`.
3504/// The eigenvectors are computed by a divide and conquer algorithm.
3505///
3506/// # Errors
3507///
3508/// Returns an error if cuSOLVER has not been initialized, if the
3509/// matrix dimensions, leading dimension, eigen selection, eigen mode, or fill
3510/// mode are invalid, if the current GPU architecture is unsupported, or if
3511/// cuSOLVER reports an internal failure.
3512pub fn dsyevdx(
3513    ctx: &Context,
3514    mode: EigenMode,
3515    range: EigenRange,
3516    fill_mode: FillMode,
3517    n: usize,
3518    a: &mut DeviceMemory<f64>,
3519    lda: usize,
3520    value_range: (f64, f64),
3521    index_range: (usize, usize),
3522    w: &mut DeviceMemory<f64>,
3523    workspace: &mut DeviceMemory<f64>,
3524    dev_info: &mut DeviceMemory<i32>,
3525) -> Result<usize> {
3526    ctx.bind()?;
3527    validate_syev_buffers(n, a.len(), lda, w.len())?;
3528    require_info_buffer(dev_info)?;
3529    let (meig, lwork) = dsyevdx_buffer_size(
3530        ctx,
3531        mode,
3532        range,
3533        fill_mode,
3534        n,
3535        a,
3536        lda,
3537        value_range,
3538        index_range,
3539        w,
3540    )?;
3541    require_workspace(workspace.len(), lwork)?;
3542    let mut meig_raw = 0;
3543    unsafe {
3544        try_ffi!(sys::cusolverDnDsyevdx(
3545            ctx.as_raw(),
3546            mode.into(),
3547            range.into(),
3548            fill_mode.into(),
3549            to_i32(n, "n")?,
3550            a.as_mut_ptr().cast(),
3551            to_i32(lda, "lda")?,
3552            value_range.0,
3553            value_range.1,
3554            to_i32(index_range.0, "il")?,
3555            to_i32(index_range.1, "iu")?,
3556            &raw mut meig_raw,
3557            w.as_mut_ptr().cast(),
3558            workspace.as_mut_ptr().cast(),
3559            to_i32(lwork, "lwork")?,
3560            dev_info.as_mut_ptr().cast(),
3561        ))?;
3562    }
3563    debug_assert_eq!(meig, to_usize(meig_raw, "meig")?);
3564    Ok(meig)
3565}
3566
3567pub fn cheevdx(
3568    ctx: &Context,
3569    mode: EigenMode,
3570    range: EigenRange,
3571    fill_mode: FillMode,
3572    n: usize,
3573    a: &mut DeviceMemory<Complex32>,
3574    lda: usize,
3575    value_range: (f32, f32),
3576    index_range: (usize, usize),
3577    w: &mut DeviceMemory<f32>,
3578    workspace: &mut DeviceMemory<Complex32>,
3579    dev_info: &mut DeviceMemory<i32>,
3580) -> Result<usize> {
3581    ctx.bind()?;
3582    validate_syev_buffers(n, a.len(), lda, w.len())?;
3583    require_info_buffer(dev_info)?;
3584    let (meig, lwork) = cheevdx_buffer_size(
3585        ctx,
3586        mode,
3587        range,
3588        fill_mode,
3589        n,
3590        a,
3591        lda,
3592        value_range,
3593        index_range,
3594        w,
3595    )?;
3596    require_workspace(workspace.len(), lwork)?;
3597    let mut meig_raw = 0;
3598    unsafe {
3599        try_ffi!(sys::cusolverDnCheevdx(
3600            ctx.as_raw(),
3601            mode.into(),
3602            range.into(),
3603            fill_mode.into(),
3604            to_i32(n, "n")?,
3605            a.as_mut_ptr().cast(),
3606            to_i32(lda, "lda")?,
3607            value_range.0,
3608            value_range.1,
3609            to_i32(index_range.0, "il")?,
3610            to_i32(index_range.1, "iu")?,
3611            &raw mut meig_raw,
3612            w.as_mut_ptr().cast(),
3613            workspace.as_mut_ptr().cast(),
3614            to_i32(lwork, "lwork")?,
3615            dev_info.as_mut_ptr().cast(),
3616        ))?;
3617    }
3618    debug_assert_eq!(meig, to_usize(meig_raw, "meig")?);
3619    Ok(meig)
3620}
3621
3622pub fn zheevdx(
3623    ctx: &Context,
3624    mode: EigenMode,
3625    range: EigenRange,
3626    fill_mode: FillMode,
3627    n: usize,
3628    a: &mut DeviceMemory<Complex64>,
3629    lda: usize,
3630    value_range: (f64, f64),
3631    index_range: (usize, usize),
3632    w: &mut DeviceMemory<f64>,
3633    workspace: &mut DeviceMemory<Complex64>,
3634    dev_info: &mut DeviceMemory<i32>,
3635) -> Result<usize> {
3636    ctx.bind()?;
3637    validate_syev_buffers(n, a.len(), lda, w.len())?;
3638    require_info_buffer(dev_info)?;
3639    let (meig, lwork) = zheevdx_buffer_size(
3640        ctx,
3641        mode,
3642        range,
3643        fill_mode,
3644        n,
3645        a,
3646        lda,
3647        value_range,
3648        index_range,
3649        w,
3650    )?;
3651    require_workspace(workspace.len(), lwork)?;
3652    let mut meig_raw = 0;
3653    unsafe {
3654        try_ffi!(sys::cusolverDnZheevdx(
3655            ctx.as_raw(),
3656            mode.into(),
3657            range.into(),
3658            fill_mode.into(),
3659            to_i32(n, "n")?,
3660            a.as_mut_ptr().cast(),
3661            to_i32(lda, "lda")?,
3662            value_range.0,
3663            value_range.1,
3664            to_i32(index_range.0, "il")?,
3665            to_i32(index_range.1, "iu")?,
3666            &raw mut meig_raw,
3667            w.as_mut_ptr().cast(),
3668            workspace.as_mut_ptr().cast(),
3669            to_i32(lwork, "lwork")?,
3670            dev_info.as_mut_ptr().cast(),
3671        ))?;
3672    }
3673    debug_assert_eq!(meig, to_usize(meig_raw, "meig")?);
3674    Ok(meig)
3675}
3676
3677fn validate_syev_buffers(n: usize, a_len: usize, lda: usize, w_len: usize) -> Result<()> {
3678    validate_square_matrix(n, a_len, lda)?;
3679    if w_len < n {
3680        return Err(Error::InvalidVectorShape);
3681    }
3682    Ok(())
3683}
3684
3685fn validate_square_matrix(n: usize, len: usize, lda: usize) -> Result<()> {
3686    validate_matrix(n, n, len, lda)
3687}
3688
3689fn validate_matrix(rows: usize, cols: usize, len: usize, lda: usize) -> Result<()> {
3690    if rows == 0 || cols == 0 {
3691        return Err(Error::InvalidMatrixShape);
3692    }
3693    if lda < rows {
3694        return Err(Error::InvalidLeadingDimension);
3695    }
3696    let required = lda.checked_mul(cols).ok_or(Error::InvalidMatrixShape)?;
3697    if len < required {
3698        return Err(Error::InvalidMatrixShape);
3699    }
3700    Ok(())
3701}
3702
3703fn require_workspace(actual: usize, required: usize) -> Result<()> {
3704    if actual < required {
3705        return Err(Error::InsufficientWorkspaceSize { required, actual });
3706    }
3707    Ok(())
3708}
3709
3710fn require_workspace_bytes(actual: usize, required: usize) -> Result<()> {
3711    if actual < required {
3712        return Err(Error::InsufficientWorkspaceSize { required, actual });
3713    }
3714    Ok(())
3715}
3716
3717fn require_host_workspace(actual: usize, required: usize) -> Result<()> {
3718    if actual < required {
3719        return Err(Error::InsufficientWorkspaceSize { required, actual });
3720    }
3721    Ok(())
3722}
3723
3724fn require_info_buffer(dev_info: &DeviceMemory<i32>) -> Result<()> {
3725    if dev_info.is_empty() {
3726        return Err(Error::InvalidVectorShape);
3727    }
3728    Ok(())
3729}
3730
3731fn validate_xsyevd_buffers(
3732    n: usize,
3733    a_bytes: usize,
3734    lda: usize,
3735    a_type: DataType,
3736    w_bytes: usize,
3737    w_type: DataType,
3738) -> Result<()> {
3739    validate_x_matrix(n, n, a_bytes, lda, a_type)?;
3740    validate_x_vector(n, w_bytes, w_type)?;
3741    Ok(())
3742}
3743
3744fn validate_xsyev_batched_buffers(
3745    n: usize,
3746    a_bytes: usize,
3747    lda: usize,
3748    a_type: DataType,
3749    w_bytes: usize,
3750    w_type: DataType,
3751    batch_count: usize,
3752) -> Result<()> {
3753    if batch_count == 0 {
3754        return Err(Error::InvalidMatrixShape);
3755    }
3756    validate_x_matrix(
3757        n,
3758        n.checked_mul(batch_count)
3759            .ok_or(Error::InvalidMatrixShape)?,
3760        a_bytes,
3761        lda,
3762        a_type,
3763    )?;
3764    validate_x_vector(
3765        n.checked_mul(batch_count)
3766            .ok_or(Error::InvalidVectorShape)?,
3767        w_bytes,
3768        w_type,
3769    )?;
3770    let problem_size = n
3771        .checked_mul(lda)
3772        .and_then(|value| value.checked_mul(batch_count))
3773        .ok_or(Error::InvalidMatrixShape)?;
3774    if problem_size > i32::MAX as usize {
3775        return Err(Error::OutOfRange {
3776            name: "batched problem size".into(),
3777        });
3778    }
3779    Ok(())
3780}
3781
3782fn validate_syevj_batched_buffers(
3783    n: usize,
3784    a_len: usize,
3785    lda: usize,
3786    w_len: usize,
3787    batch_count: usize,
3788) -> Result<()> {
3789    if batch_count == 0 {
3790        return Err(Error::InvalidMatrixShape);
3791    }
3792    let matrix_cols = n
3793        .checked_mul(batch_count)
3794        .ok_or(Error::InvalidMatrixShape)?;
3795    validate_matrix(n, matrix_cols, a_len, lda)?;
3796    let eigenvalues = n
3797        .checked_mul(batch_count)
3798        .ok_or(Error::InvalidVectorShape)?;
3799    if w_len < eigenvalues {
3800        return Err(Error::InvalidVectorShape);
3801    }
3802    Ok(())
3803}
3804
3805fn validate_sygvj_buffers(
3806    n: usize,
3807    a_len: usize,
3808    lda: usize,
3809    b_len: usize,
3810    ldb: usize,
3811    w_len: usize,
3812) -> Result<()> {
3813    validate_square_matrix(n, a_len, lda)?;
3814    validate_square_matrix(n, b_len, ldb)?;
3815    if w_len < n {
3816        return Err(Error::InvalidVectorShape);
3817    }
3818    Ok(())
3819}
3820
3821fn validate_sygvd_buffers(
3822    n: usize,
3823    a_len: usize,
3824    lda: usize,
3825    b_len: usize,
3826    ldb: usize,
3827    w_len: usize,
3828) -> Result<()> {
3829    validate_square_matrix(n, a_len, lda)?;
3830    validate_square_matrix(n, b_len, ldb)?;
3831    if w_len < n {
3832        return Err(Error::InvalidVectorShape);
3833    }
3834    Ok(())
3835}
3836
3837fn validate_x_matrix(
3838    rows: usize,
3839    cols: usize,
3840    bytes: usize,
3841    lda: usize,
3842    data_type: DataType,
3843) -> Result<()> {
3844    if rows == 0 || cols == 0 {
3845        return Err(Error::InvalidMatrixShape);
3846    }
3847    if lda < rows {
3848        return Err(Error::InvalidLeadingDimension);
3849    }
3850    let required = lda
3851        .checked_mul(cols)
3852        .and_then(|count| count.checked_mul(data_type.size_in_bytes()))
3853        .ok_or(Error::InvalidMatrixShape)?;
3854    if bytes < required {
3855        return Err(Error::InvalidMatrixShape);
3856    }
3857    Ok(())
3858}
3859
3860fn validate_x_vector(len: usize, bytes: usize, data_type: DataType) -> Result<()> {
3861    let required = len
3862        .checked_mul(data_type.size_in_bytes())
3863        .ok_or(Error::InvalidVectorShape)?;
3864    if bytes < required {
3865        return Err(Error::InvalidVectorShape);
3866    }
3867    Ok(())
3868}
3869
3870fn validate_xsyevdx_value_type<T>(w_type: DataType) -> Result<()> {
3871    if size_of::<T>() != w_type.size_in_bytes() {
3872        return Err(Error::InvalidEigRange);
3873    }
3874    Ok(())
3875}
3876
3877type EigenSelectionParts<T> = (EigenRange, Option<(T, T)>, Option<(usize, usize)>);
3878
3879fn validate_xsyevdx_range<T>(
3880    range: EigenRange,
3881    n: usize,
3882    value_range: Option<(T, T)>,
3883    index_range: Option<(usize, usize)>,
3884) -> Result<(T, T, usize, usize)>
3885where
3886    T: Default,
3887{
3888    match range {
3889        EigenRange::All => {
3890            if value_range.is_some() || index_range.is_some() {
3891                return Err(Error::InvalidEigRange);
3892            }
3893            Ok((T::default(), T::default(), 0, 0))
3894        }
3895        EigenRange::Value => {
3896            if index_range.is_some() {
3897                return Err(Error::InvalidEigRange);
3898            }
3899            let Some((vl, vu)) = value_range else {
3900                return Err(Error::InvalidEigRange);
3901            };
3902            Ok((vl, vu, 0, 0))
3903        }
3904        EigenRange::Index => {
3905            if value_range.is_some() {
3906                return Err(Error::InvalidEigRange);
3907            }
3908            let Some((il, iu)) = index_range else {
3909                return Err(Error::InvalidEigRange);
3910            };
3911            if il == 0 || iu == 0 || il > iu || iu > n {
3912                return Err(Error::InvalidEigRange);
3913            }
3914            Ok((T::default(), T::default(), il, iu))
3915        }
3916    }
3917}
3918
3919fn selection_parts<T>(selection: EigenSelection<T>) -> EigenSelectionParts<T> {
3920    match selection {
3921        EigenSelection::All => (EigenRange::All, None, None),
3922        EigenSelection::ByValue { lower, upper } => (EigenRange::Value, Some((lower, upper)), None),
3923        EigenSelection::ByIndex { start, end } => (EigenRange::Index, None, Some((start, end))),
3924    }
3925}
3926
3927fn matrix_ref_parts<T>(matrix: Option<MatrixRef<'_, T>>) -> Option<(&DeviceMemory<T>, usize)> {
3928    matrix.map(|matrix| (matrix.data, matrix.leading_dimension))
3929}
3930
3931fn matrix_mut_parts<T>(matrix: Option<MatrixMut<'_, T>>) -> Option<(&mut DeviceMemory<T>, usize)> {
3932    matrix.map(|matrix| (matrix.data, matrix.leading_dimension))
3933}
3934
3935fn matrix_mut_ref_parts<'a, T>(
3936    matrix: Option<&'a MatrixMut<'a, T>>,
3937) -> Option<(&'a DeviceMemory<T>, usize)> {
3938    matrix.map(|matrix| (&*matrix.data, matrix.leading_dimension))
3939}
3940
3941fn matrix_mut_ref_option<'a, T>(matrix: Option<&'a MatrixMut<'a, T>>) -> Option<MatrixRef<'a, T>> {
3942    matrix.map(MatrixMut::as_ref)
3943}
3944
3945fn validate_xgeev_inputs<TV>(
3946    n: usize,
3947    a_bytes: usize,
3948    lda: usize,
3949    a_type: DataType,
3950    w_bytes: usize,
3951    w_type: DataType,
3952    right_vectors: Option<(&DeviceMemory<TV>, usize)>,
3953    vr_type: DataType,
3954) -> Result<()> {
3955    validate_x_matrix(n, n, a_bytes, lda, a_type)?;
3956    validate_xgeev_eigenvalues(n, w_bytes, a_type, w_type)?;
3957    if let Some((vr, ldvr)) = right_vectors {
3958        validate_x_matrix(n, n, vr.byte_len(), ldvr, vr_type)?;
3959    }
3960    Ok(())
3961}
3962
3963fn validate_xgeev_eigenvalues(
3964    n: usize,
3965    w_bytes: usize,
3966    a_type: DataType,
3967    w_type: DataType,
3968) -> Result<()> {
3969    let expected_len = match (a_type, w_type) {
3970        (DataType::F32, DataType::F32) | (DataType::F64, DataType::F64) => {
3971            n.checked_mul(2).ok_or(Error::InvalidVectorShape)?
3972        }
3973        (DataType::F32, DataType::ComplexF32)
3974        | (DataType::F64, DataType::ComplexF64)
3975        | (DataType::ComplexF32, DataType::ComplexF32)
3976        | (DataType::ComplexF64, DataType::ComplexF64) => n,
3977        _ => return Err(Error::InvalidVectorShape),
3978    };
3979    validate_x_vector(expected_len, w_bytes, w_type)
3980}
3981
3982fn optional_xgeev_matrix_ptr<T: DataTypeLike>(
3983    matrix: Option<(&DeviceMemory<T>, usize)>,
3984) -> Result<(*const T, i64)> {
3985    match matrix {
3986        Some((matrix, ld)) => Ok((matrix.as_ptr().cast(), to_i64(ld, "ldvr")?)),
3987        None => Ok((ptr::null(), 1)),
3988    }
3989}
3990
3991fn optional_xgeev_matrix_mut_ptr<T: DataTypeLike>(
3992    matrix: Option<(&mut DeviceMemory<T>, usize)>,
3993) -> Result<(*mut T, i64)> {
3994    match matrix {
3995        Some((matrix, ld)) => Ok((matrix.as_mut_ptr().cast(), to_i64(ld, "ldvr")?)),
3996        None => Ok((ptr::null_mut(), 1)),
3997    }
3998}
3999
4000fn require_info_buffer_len(dev_info: &DeviceMemory<i32>, required: usize) -> Result<()> {
4001    if dev_info.len() < required {
4002        return Err(Error::InvalidVectorShape);
4003    }
4004    Ok(())
4005}
4006
4007#[cfg(test)]
4008mod tests {
4009    use super::{EigenSelection, selection_parts};
4010    use crate::types::EigenRange;
4011
4012    #[test]
4013    fn eigen_selection_all_maps_cleanly() {
4014        let (range, values, indices) = selection_parts::<f32>(EigenSelection::All);
4015        assert_eq!(range, EigenRange::All);
4016        assert_eq!(values, None);
4017        assert_eq!(indices, None);
4018    }
4019
4020    #[test]
4021    fn eigen_selection_index_maps_cleanly() {
4022        let (range, values, indices) =
4023            selection_parts::<f64>(EigenSelection::ByIndex { start: 2, end: 5 });
4024        assert_eq!(range, EigenRange::Index);
4025        assert_eq!(values, None);
4026        assert_eq!(indices, Some((2, 5)));
4027    }
4028}
singe_cusolver/eigen.rs

singe_cusolver/
eigen.rs
