|
1 | | -use crate::imp_prelude::*; |
2 | | -use num_complex::Complex; |
3 | | -use rawpointer::PointerExt; |
4 | | -use std::mem; |
5 | | -use std::ptr::NonNull; |
6 | | - |
7 | 1 | mod impl_numeric; |
8 | | - |
9 | | -impl<T, S, D> ArrayBase<S, D> |
10 | | -where |
11 | | - S: Data<Elem = Complex<T>>, |
12 | | - D: Dimension, |
13 | | -{ |
14 | | - /// Returns views of the real and imaginary components of the elements. |
15 | | - /// |
16 | | - /// ``` |
17 | | - /// use ndarray::prelude::*; |
18 | | - /// use num_complex::{Complex, Complex64}; |
19 | | - /// |
20 | | - /// let arr = array![ |
21 | | - /// [Complex64::new(1., 2.), Complex64::new(3., 4.)], |
22 | | - /// [Complex64::new(5., 6.), Complex64::new(7., 8.)], |
23 | | - /// [Complex64::new(9., 10.), Complex64::new(11., 12.)], |
24 | | - /// ]; |
25 | | - /// let Complex { re, im } = arr.view_re_im(); |
26 | | - /// assert_eq!(re, array![[1., 3.], [5., 7.], [9., 11.]]); |
27 | | - /// assert_eq!(im, array![[2., 4.], [6., 8.], [10., 12.]]); |
28 | | - /// ``` |
29 | | - pub fn view_re_im(&self) -> Complex<ArrayView<'_, T, D>> { |
30 | | - debug_assert!(self.pointer_is_inbounds()); |
31 | | - |
32 | | - let dim = self.dim.clone(); |
33 | | - |
34 | | - // Double the strides. In the zero-sized element case and for axes of |
35 | | - // length <= 1, we leave the strides as-is to avoid possible overflow. |
36 | | - let mut strides = self.strides.clone(); |
37 | | - if mem::size_of::<T>() != 0 { |
38 | | - for ax in 0..strides.ndim() { |
39 | | - if dim[ax] > 1 { |
40 | | - strides[ax] *= 2; |
41 | | - } |
42 | | - } |
43 | | - } |
44 | | - |
45 | | - let ptr_re: NonNull<T> = self.ptr.cast(); |
46 | | - let ptr_im: NonNull<T> = if self.is_empty() { |
47 | | - // In the empty case, we can just reuse the existing pointer since |
48 | | - // it won't be dereferenced anyway. It is not safe to offset by one |
49 | | - // since the allocation may be empty. |
50 | | - ptr_re |
51 | | - } else { |
52 | | - // In the nonempty case, we can safely offset into the first |
53 | | - // (complex) element. |
54 | | - unsafe { ptr_re.add(1) } |
55 | | - }; |
56 | | - |
57 | | - // `Complex` is `repr(C)` with only fields `re: T` and `im: T`. So, the |
58 | | - // real components of the elements start at the same pointer, and the |
59 | | - // imaginary components start at the pointer offset by one, with |
60 | | - // exactly double the strides. The new, doubled strides still meet the |
61 | | - // overflow constraints: |
62 | | - // |
63 | | - // - For the zero-sized element case, the strides are unchanged in |
64 | | - // units of bytes and in units of the element type. |
65 | | - // |
66 | | - // - For the nonzero-sized element case: |
67 | | - // |
68 | | - // - In units of bytes, the strides are unchanged. |
69 | | - // |
70 | | - // - Since `Complex<T>` for nonzero `T` is always at least 2 bytes, |
71 | | - // and the original strides did not overflow in units of bytes, we |
72 | | - // know that the new doubled strides will not overflow in units of |
73 | | - // `T`. |
74 | | - unsafe { |
75 | | - Complex { |
76 | | - re: ArrayView::new(ptr_re, dim.clone(), strides.clone()), |
77 | | - im: ArrayView::new(ptr_im, dim, strides), |
78 | | - } |
79 | | - } |
80 | | - } |
81 | | - |
82 | | - /// Returns mutable views of the real and imaginary components of the elements. |
83 | | - /// |
84 | | - /// ``` |
85 | | - /// use ndarray::prelude::*; |
86 | | - /// use num_complex::{Complex, Complex64}; |
87 | | - /// |
88 | | - /// let mut arr = array![ |
89 | | - /// [Complex64::new(1., 2.), Complex64::new(3., 4.)], |
90 | | - /// [Complex64::new(5., 6.), Complex64::new(7., 8.)], |
91 | | - /// [Complex64::new(9., 10.), Complex64::new(11., 12.)], |
92 | | - /// ]; |
93 | | - /// |
94 | | - /// let Complex { mut re, mut im } = arr.view_mut_re_im(); |
95 | | - /// assert_eq!(re, array![[1., 3.], [5., 7.], [9., 11.]]); |
96 | | - /// assert_eq!(im, array![[2., 4.], [6., 8.], [10., 12.]]); |
97 | | - /// |
98 | | - /// re[[0, 1]] = 13.; |
99 | | - /// im[[2, 0]] = 14.; |
100 | | - /// |
101 | | - /// assert_eq!(arr[[0, 1]], Complex64::new(13., 4.)); |
102 | | - /// assert_eq!(arr[[2, 0]], Complex64::new(9., 14.)); |
103 | | - /// ``` |
104 | | - pub fn view_mut_re_im(&mut self) -> Complex<ArrayViewMut<'_, T, D>> |
105 | | - where |
106 | | - S: DataMut, |
107 | | - { |
108 | | - self.ensure_unique(); |
109 | | - |
110 | | - let dim = self.dim.clone(); |
111 | | - |
112 | | - // Double the strides. In the zero-sized element case and for axes of |
113 | | - // length <= 1, we leave the strides as-is to avoid possible overflow. |
114 | | - let mut strides = self.strides.clone(); |
115 | | - if mem::size_of::<T>() != 0 { |
116 | | - for ax in 0..strides.ndim() { |
117 | | - if dim[ax] > 1 { |
118 | | - strides[ax] *= 2; |
119 | | - } |
120 | | - } |
121 | | - } |
122 | | - |
123 | | - let ptr_re: NonNull<T> = self.ptr.cast(); |
124 | | - let ptr_im: NonNull<T> = if self.is_empty() { |
125 | | - // In the empty case, we can just reuse the existing pointer since |
126 | | - // it won't be dereferenced anyway. It is not safe to offset by one |
127 | | - // since the allocation may be empty. |
128 | | - ptr_re |
129 | | - } else { |
130 | | - // In the nonempty case, we can safely offset into the first |
131 | | - // (complex) element. |
132 | | - unsafe { ptr_re.add(1) } |
133 | | - }; |
134 | | - |
135 | | - // `Complex` is `repr(C)` with only fields `re: T` and `im: T`. So, the |
136 | | - // real components of the elements start at the same pointer, and the |
137 | | - // imaginary components start at the pointer offset by one, with |
138 | | - // exactly double the strides. The new, doubled strides still meet the |
139 | | - // overflow constraints: |
140 | | - // |
141 | | - // - For the zero-sized element case, the strides are unchanged in |
142 | | - // units of bytes and in units of the element type. |
143 | | - // |
144 | | - // - For the nonzero-sized element case: |
145 | | - // |
146 | | - // - In units of bytes, the strides are unchanged. |
147 | | - // |
148 | | - // - Since `Complex<T>` for nonzero `T` is always at least 2 bytes, |
149 | | - // and the original strides did not overflow in units of bytes, we |
150 | | - // know that the new doubled strides will not overflow in units of |
151 | | - // `T`. |
152 | | - unsafe { |
153 | | - Complex { |
154 | | - re: ArrayViewMut::new(ptr_re, dim.clone(), strides.clone()), |
155 | | - im: ArrayViewMut::new(ptr_im, dim, strides), |
156 | | - } |
157 | | - } |
158 | | - } |
159 | | -} |
0 commit comments