scirs2_core/array_protocol/
operations.rs

1// Copyright (c) 2025, `SciRS2` Team
2//
3// Licensed under the Apache License, Version 2.0
4// (LICENSE-APACHE or http://www.apache.org/licenses/LICENSE-2.0)
5//
6
7//! Common array operations implemented with the array protocol.
8//!
9//! This module provides implementations of common array operations using
10//! the array protocol. These operations can work with any array type that
11//! implements the ArrayProtocol trait, including GPU arrays, distributed arrays,
12//! and custom third-party array implementations.
13
14use std::any::{Any, TypeId};
15use std::collections::HashMap;
16
17use ::ndarray::{Ix1, Ix2, IxDyn};
18
19use crate::array_protocol::{
20    get_implementing_args, ArrayFunction, ArrayProtocol, NdarrayWrapper, NotImplemented,
21};
22use crate::error::CoreError;
23// Note: num_traits not needed for current implementation
24
25/// Error type for array operations.
26#[derive(Debug, thiserror::Error)]
27pub enum OperationError {
28    /// The operation is not implemented for the given array types.
29    #[error("Operation not implemented: {0}")]
30    NotImplemented(String),
31    /// The array shapes are incompatible for the operation.
32    #[error("Shape mismatch: {0}")]
33    ShapeMismatch(String),
34    /// The array types are incompatible for the operation.
35    #[error("Type mismatch: {0}")]
36    TypeMismatch(String),
37    /// Other error during operation.
38    #[error("Operation error: {0}")]
39    Other(String),
40}
41
42impl From<NotImplemented> for OperationError {
43    fn from(_: NotImplemented) -> Self {
44        Self::NotImplemented("Operation not implemented for these array types".to_string())
45    }
46}
47
48impl From<CoreError> for OperationError {
49    fn from(err: CoreError) -> Self {
50        Self::Other(err.to_string())
51    }
52}
53
54// Define array operations using the array protocol
55
56// Define a macro for implementing array operations
57#[macro_export]
58macro_rules! array_function_dispatch {
59    // For normal functions
60    (fn $name:ident($($arg:ident: $arg_ty:ty),*) -> Result<$ret:ty, $err:ty> $body:block, $funcname:expr) => {
61        pub fn $name($($arg: $arg_ty),*) -> Result<$ret, $err> $body
62    };
63
64    // For normal functions with trailing commas
65    (fn $name:ident($($arg:ident: $arg_ty:ty,)*) -> Result<$ret:ty, $err:ty> $body:block, $funcname:expr) => {
66        pub fn $name($($arg: $arg_ty),*) -> Result<$ret, $err> $body
67    };
68
69    // For generic functions
70    (fn $name:ident<$($type_param:ident $(: $type_bound:path)?),*>($($arg:ident: $arg_ty:ty),*) -> Result<$ret:ty, $err:ty> $body:block, $funcname:expr) => {
71        pub fn $name <$($type_param $(: $type_bound)?),*>($($arg: $arg_ty),*) -> Result<$ret, $err> $body
72    };
73
74    // For generic functions with trailing commas
75    (fn $name:ident<$($type_param:ident $(: $type_bound:path)?),*>($($arg:ident: $arg_ty:ty,)*) -> Result<$ret:ty, $err:ty> $body:block, $funcname:expr) => {
76        pub fn $name <$($type_param $(: $type_bound)?),*>($($arg: $arg_ty),*) -> Result<$ret, $err> $body
77    };
78}
79
80// Matrix multiplication operation
81array_function_dispatch!(
82    fn matmul(
83        a: &dyn ArrayProtocol,
84        b: &dyn ArrayProtocol,
85    ) -> Result<Box<dyn ArrayProtocol>, OperationError> {
86        // Get implementing args
87        let boxed_a = Box::new(a.box_clone());
88        let boxed_b = Box::new(b.box_clone());
89        let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a, boxed_b];
90        let implementing_args =
91            get_implementing_args("scirs2::array_protocol::operations::matmul", &boxed_args);
92        if implementing_args.is_empty() {
93            // Comprehensive fallback implementation for ndarray types
94
95            // f64 types with Ix2 dimension (static dimension size)
96            if let (Some(a_array), Some(b_array)) = (
97                a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>(),
98                b.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>(),
99            ) {
100                let a_array_owned = a_array.as_array().clone();
101                let b_array_owned = b_array.as_array().clone();
102                let (m, k) = a_array_owned.dim();
103                let (_, n) = b_array_owned.dim();
104                let mut result = crate::ndarray::Array2::<f64>::zeros((m, n));
105                for i in 0..m {
106                    for j in 0..n {
107                        let mut sum = 0.0;
108                        for l in 0..k {
109                            sum += a_array_owned[[i, l]] * b_array_owned[[l, j]];
110                        }
111                        result[[i, j]] = sum;
112                    }
113                }
114                return Ok(Box::new(NdarrayWrapper::new(result)));
115            }
116
117            // f64 types with IxDyn dimension (dynamic dimension size)
118            if let (Some(a_array), Some(b_array)) = (
119                a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>(),
120                b.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>(),
121            ) {
122                let a_array_owned = a_array.as_array().to_owned();
123                let b_array_owned = b_array.as_array().to_owned();
124                let a_dim = a_array_owned.shape();
125                let b_dim = b_array_owned.shape();
126                if a_dim.len() != 2 || b_dim.len() != 2 || a_dim[1] != b_dim[0] {
127                    return Err(OperationError::ShapeMismatch(format!(
128                        "Invalid shapes for matmul: {a_dim:?} and {b_dim:?}"
129                    )));
130                }
131                let (m, k) = (a_dim[0], a_dim[1]);
132                let n = b_dim[1];
133                let mut result = crate::ndarray::Array2::<f64>::zeros((m, n));
134                for i in 0..m {
135                    for j in 0..n {
136                        let mut sum = 0.0;
137                        for l in 0..k {
138                            sum += a_array_owned[[i, l]] * b_array_owned[[l, j]];
139                        }
140                        result[[i, j]] = sum;
141                    }
142                }
143                return Ok(Box::new(NdarrayWrapper::new(result)));
144            }
145
146            // f32 types with Ix2 dimension
147            if let (Some(a_array), Some(b_array)) = (
148                a.as_any().downcast_ref::<NdarrayWrapper<f32, Ix2>>(),
149                b.as_any().downcast_ref::<NdarrayWrapper<f32, Ix2>>(),
150            ) {
151                let a_array_owned = a_array.as_array().clone();
152                let b_array_owned = b_array.as_array().clone();
153                let (m, k) = a_array_owned.dim();
154                let (_, n) = b_array_owned.dim();
155                let mut result = crate::ndarray::Array2::<f32>::zeros((m, n));
156                for i in 0..m {
157                    for j in 0..n {
158                        let mut sum = 0.0;
159                        for l in 0..k {
160                            sum += a_array_owned[[i, l]] * b_array_owned[[l, j]];
161                        }
162                        result[[i, j]] = sum;
163                    }
164                }
165                return Ok(Box::new(NdarrayWrapper::new(result)));
166            }
167
168            // f32 types with IxDyn dimension
169            if let (Some(a_array), Some(b_array)) = (
170                a.as_any().downcast_ref::<NdarrayWrapper<f32, IxDyn>>(),
171                b.as_any().downcast_ref::<NdarrayWrapper<f32, IxDyn>>(),
172            ) {
173                let a_array_owned = a_array.as_array().to_owned();
174                let b_array_owned = b_array.as_array().to_owned();
175                let a_dim = a_array_owned.shape();
176                let b_dim = b_array_owned.shape();
177                if a_dim.len() != 2 || b_dim.len() != 2 || a_dim[1] != b_dim[0] {
178                    return Err(OperationError::ShapeMismatch(format!(
179                        "Invalid shapes for matmul: {a_dim:?} and {b_dim:?}"
180                    )));
181                }
182                let (m, k) = (a_dim[0], a_dim[1]);
183                let n = b_dim[1];
184                let mut result = crate::ndarray::Array2::<f32>::zeros((m, n));
185                for i in 0..m {
186                    for j in 0..n {
187                        let mut sum = 0.0;
188                        for l in 0..k {
189                            sum += a_array_owned[[i, l]] * b_array_owned[[l, j]];
190                        }
191                        result[[i, j]] = sum;
192                    }
193                }
194                return Ok(Box::new(NdarrayWrapper::new(result)));
195            }
196
197            return Err(OperationError::NotImplemented(
198                "matmul not implemented for these array types".to_string(),
199            ));
200        }
201
202        // Delegate to the implementation
203        let array_ref = implementing_args[0].1;
204
205        let result = array_ref.array_function(
206            &ArrayFunction::new("scirs2::array_protocol::operations::matmul"),
207            &[TypeId::of::<Box<dyn ArrayProtocol>>()],
208            &[Box::new(a.box_clone()), Box::new(b.box_clone())],
209            &HashMap::new(),
210        )?;
211
212        // Try to downcast the result
213        match result.downcast::<Box<dyn ArrayProtocol>>() {
214            Ok(array) => Ok(*array),
215            Err(_) => Err(OperationError::Other(
216                "Failed to downcast result to ArrayProtocol".to_string(),
217            )),
218        }
219    },
220    "scirs2::array_protocol::operations::matmul"
221);
222
223// Element-wise addition operation
224array_function_dispatch!(
225    fn add(
226        a: &dyn ArrayProtocol,
227        b: &dyn ArrayProtocol,
228    ) -> Result<Box<dyn ArrayProtocol>, OperationError> {
229        // Get implementing args
230        let boxed_a = Box::new(a.box_clone());
231        let boxed_b = Box::new(b.box_clone());
232        let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a, boxed_b];
233        let implementing_args =
234            get_implementing_args("scirs2::array_protocol::operations::add", &boxed_args);
235        if implementing_args.is_empty() {
236            // Comprehensive fallback implementation for ndarray types
237
238            // f64 types
239            if let (Some(a_array), Some(b_array)) = (
240                a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix1>>(),
241                b.as_any().downcast_ref::<NdarrayWrapper<f64, Ix1>>(),
242            ) {
243                let result = a_array.as_array() + b_array.as_array();
244                return Ok(Box::new(NdarrayWrapper::new(result)));
245            }
246            if let (Some(a_array), Some(b_array)) = (
247                a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>(),
248                b.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>(),
249            ) {
250                let result = a_array.as_array() + b_array.as_array();
251                return Ok(Box::new(NdarrayWrapper::new(result)));
252            }
253            if let (Some(a_array), Some(b_array)) = (
254                a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>(),
255                b.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>(),
256            ) {
257                let result = a_array.as_array() + b_array.as_array();
258                return Ok(Box::new(NdarrayWrapper::new(result)));
259            }
260
261            // f32 types
262            if let (Some(a_array), Some(b_array)) = (
263                a.as_any().downcast_ref::<NdarrayWrapper<f32, Ix1>>(),
264                b.as_any().downcast_ref::<NdarrayWrapper<f32, Ix1>>(),
265            ) {
266                let result = a_array.as_array() + b_array.as_array();
267                return Ok(Box::new(NdarrayWrapper::new(result)));
268            }
269            if let (Some(a_array), Some(b_array)) = (
270                a.as_any().downcast_ref::<NdarrayWrapper<f32, Ix2>>(),
271                b.as_any().downcast_ref::<NdarrayWrapper<f32, Ix2>>(),
272            ) {
273                let result = a_array.as_array() + b_array.as_array();
274                return Ok(Box::new(NdarrayWrapper::new(result)));
275            }
276            if let (Some(a_array), Some(b_array)) = (
277                a.as_any().downcast_ref::<NdarrayWrapper<f32, IxDyn>>(),
278                b.as_any().downcast_ref::<NdarrayWrapper<f32, IxDyn>>(),
279            ) {
280                let result = a_array.as_array() + b_array.as_array();
281                return Ok(Box::new(NdarrayWrapper::new(result)));
282            }
283
284            // i32 types
285            if let (Some(a_array), Some(b_array)) = (
286                a.as_any().downcast_ref::<NdarrayWrapper<i32, Ix1>>(),
287                b.as_any().downcast_ref::<NdarrayWrapper<i32, Ix1>>(),
288            ) {
289                let result = a_array.as_array() + b_array.as_array();
290                return Ok(Box::new(NdarrayWrapper::new(result)));
291            }
292            if let (Some(a_array), Some(b_array)) = (
293                a.as_any().downcast_ref::<NdarrayWrapper<i32, Ix2>>(),
294                b.as_any().downcast_ref::<NdarrayWrapper<i32, Ix2>>(),
295            ) {
296                let result = a_array.as_array() + b_array.as_array();
297                return Ok(Box::new(NdarrayWrapper::new(result)));
298            }
299            if let (Some(a_array), Some(b_array)) = (
300                a.as_any().downcast_ref::<NdarrayWrapper<i32, IxDyn>>(),
301                b.as_any().downcast_ref::<NdarrayWrapper<i32, IxDyn>>(),
302            ) {
303                let result = a_array.as_array() + b_array.as_array();
304                return Ok(Box::new(NdarrayWrapper::new(result)));
305            }
306
307            // i64 types
308            if let (Some(a_array), Some(b_array)) = (
309                a.as_any().downcast_ref::<NdarrayWrapper<i64, Ix1>>(),
310                b.as_any().downcast_ref::<NdarrayWrapper<i64, Ix1>>(),
311            ) {
312                let result = a_array.as_array() + b_array.as_array();
313                return Ok(Box::new(NdarrayWrapper::new(result)));
314            }
315            if let (Some(a_array), Some(b_array)) = (
316                a.as_any().downcast_ref::<NdarrayWrapper<i64, Ix2>>(),
317                b.as_any().downcast_ref::<NdarrayWrapper<i64, Ix2>>(),
318            ) {
319                let result = a_array.as_array() + b_array.as_array();
320                return Ok(Box::new(NdarrayWrapper::new(result)));
321            }
322            if let (Some(a_array), Some(b_array)) = (
323                a.as_any().downcast_ref::<NdarrayWrapper<i64, IxDyn>>(),
324                b.as_any().downcast_ref::<NdarrayWrapper<i64, IxDyn>>(),
325            ) {
326                let result = a_array.as_array() + b_array.as_array();
327                return Ok(Box::new(NdarrayWrapper::new(result)));
328            }
329
330            return Err(OperationError::NotImplemented(
331                "add not implemented for these array types".to_string(),
332            ));
333        }
334
335        // Delegate to the implementation
336        let array_ref = implementing_args[0].1;
337
338        let result = array_ref.array_function(
339            &ArrayFunction::new("scirs2::array_protocol::operations::add"),
340            &[TypeId::of::<Box<dyn ArrayProtocol>>()],
341            &[Box::new(a.box_clone()), Box::new(b.box_clone())],
342            &HashMap::new(),
343        )?;
344
345        // Try to downcast the result
346        match result.downcast::<Box<dyn ArrayProtocol>>() {
347            Ok(array) => Ok(*array),
348            Err(_) => Err(OperationError::Other(
349                "Failed to downcast result to ArrayProtocol".to_string(),
350            )),
351        }
352    },
353    "scirs2::array_protocol::operations::add"
354);
355
356// Element-wise subtraction operation
357array_function_dispatch!(
358    fn subtract(
359        a: &dyn ArrayProtocol,
360        b: &dyn ArrayProtocol,
361    ) -> Result<Box<dyn ArrayProtocol>, OperationError> {
362        // Get implementing args
363        let boxed_a = Box::new(a.box_clone());
364        let boxed_b = Box::new(b.box_clone());
365        let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a, boxed_b];
366        let implementing_args =
367            get_implementing_args("scirs2::array_protocol::operations::subtract", &boxed_args);
368        if implementing_args.is_empty() {
369            // Comprehensive fallback implementation for ndarray types
370
371            // f64 types
372            if let (Some(a_array), Some(b_array)) = (
373                a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix1>>(),
374                b.as_any().downcast_ref::<NdarrayWrapper<f64, Ix1>>(),
375            ) {
376                let result = a_array.as_array() - b_array.as_array();
377                return Ok(Box::new(NdarrayWrapper::new(result)));
378            }
379            if let (Some(a_array), Some(b_array)) = (
380                a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>(),
381                b.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>(),
382            ) {
383                let result = a_array.as_array() - b_array.as_array();
384                return Ok(Box::new(NdarrayWrapper::new(result)));
385            }
386            if let (Some(a_array), Some(b_array)) = (
387                a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>(),
388                b.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>(),
389            ) {
390                let result = a_array.as_array() - b_array.as_array();
391                return Ok(Box::new(NdarrayWrapper::new(result)));
392            }
393
394            // f32 types
395            if let (Some(a_array), Some(b_array)) = (
396                a.as_any().downcast_ref::<NdarrayWrapper<f32, Ix1>>(),
397                b.as_any().downcast_ref::<NdarrayWrapper<f32, Ix1>>(),
398            ) {
399                let result = a_array.as_array() - b_array.as_array();
400                return Ok(Box::new(NdarrayWrapper::new(result)));
401            }
402            if let (Some(a_array), Some(b_array)) = (
403                a.as_any().downcast_ref::<NdarrayWrapper<f32, Ix2>>(),
404                b.as_any().downcast_ref::<NdarrayWrapper<f32, Ix2>>(),
405            ) {
406                let result = a_array.as_array() - b_array.as_array();
407                return Ok(Box::new(NdarrayWrapper::new(result)));
408            }
409            if let (Some(a_array), Some(b_array)) = (
410                a.as_any().downcast_ref::<NdarrayWrapper<f32, IxDyn>>(),
411                b.as_any().downcast_ref::<NdarrayWrapper<f32, IxDyn>>(),
412            ) {
413                let result = a_array.as_array() - b_array.as_array();
414                return Ok(Box::new(NdarrayWrapper::new(result)));
415            }
416
417            // i32 types
418            if let (Some(a_array), Some(b_array)) = (
419                a.as_any().downcast_ref::<NdarrayWrapper<i32, Ix1>>(),
420                b.as_any().downcast_ref::<NdarrayWrapper<i32, Ix1>>(),
421            ) {
422                let result = a_array.as_array() - b_array.as_array();
423                return Ok(Box::new(NdarrayWrapper::new(result)));
424            }
425            if let (Some(a_array), Some(b_array)) = (
426                a.as_any().downcast_ref::<NdarrayWrapper<i32, Ix2>>(),
427                b.as_any().downcast_ref::<NdarrayWrapper<i32, Ix2>>(),
428            ) {
429                let result = a_array.as_array() - b_array.as_array();
430                return Ok(Box::new(NdarrayWrapper::new(result)));
431            }
432            if let (Some(a_array), Some(b_array)) = (
433                a.as_any().downcast_ref::<NdarrayWrapper<i32, IxDyn>>(),
434                b.as_any().downcast_ref::<NdarrayWrapper<i32, IxDyn>>(),
435            ) {
436                let result = a_array.as_array() - b_array.as_array();
437                return Ok(Box::new(NdarrayWrapper::new(result)));
438            }
439
440            // i64 types
441            if let (Some(a_array), Some(b_array)) = (
442                a.as_any().downcast_ref::<NdarrayWrapper<i64, Ix1>>(),
443                b.as_any().downcast_ref::<NdarrayWrapper<i64, Ix1>>(),
444            ) {
445                let result = a_array.as_array() - b_array.as_array();
446                return Ok(Box::new(NdarrayWrapper::new(result)));
447            }
448            if let (Some(a_array), Some(b_array)) = (
449                a.as_any().downcast_ref::<NdarrayWrapper<i64, Ix2>>(),
450                b.as_any().downcast_ref::<NdarrayWrapper<i64, Ix2>>(),
451            ) {
452                let result = a_array.as_array() - b_array.as_array();
453                return Ok(Box::new(NdarrayWrapper::new(result)));
454            }
455            if let (Some(a_array), Some(b_array)) = (
456                a.as_any().downcast_ref::<NdarrayWrapper<i64, IxDyn>>(),
457                b.as_any().downcast_ref::<NdarrayWrapper<i64, IxDyn>>(),
458            ) {
459                let result = a_array.as_array() - b_array.as_array();
460                return Ok(Box::new(NdarrayWrapper::new(result)));
461            }
462
463            return Err(OperationError::NotImplemented(
464                "subtract not implemented for these array types".to_string(),
465            ));
466        }
467
468        // Delegate to the implementation
469        let array_ref = implementing_args[0].1;
470
471        let result = array_ref.array_function(
472            &ArrayFunction::new("scirs2::array_protocol::operations::subtract"),
473            &[TypeId::of::<Box<dyn ArrayProtocol>>()],
474            &[Box::new(a.box_clone()), Box::new(b.box_clone())],
475            &HashMap::new(),
476        )?;
477
478        // Try to downcast the result
479        match result.downcast::<Box<dyn ArrayProtocol>>() {
480            Ok(array) => Ok(*array),
481            Err(_) => Err(OperationError::Other(
482                "Failed to downcast result to ArrayProtocol".to_string(),
483            )),
484        }
485    },
486    "scirs2::array_protocol::operations::subtract"
487);
488
489// Element-wise multiplication operation
490array_function_dispatch!(
491    fn multiply(
492        a: &dyn ArrayProtocol,
493        b: &dyn ArrayProtocol,
494    ) -> Result<Box<dyn ArrayProtocol>, OperationError> {
495        // Get implementing args
496        let boxed_a = Box::new(a.box_clone());
497        let boxed_b = Box::new(b.box_clone());
498        let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a, boxed_b];
499        let implementing_args =
500            get_implementing_args("scirs2::array_protocol::operations::multiply", &boxed_args);
501        if implementing_args.is_empty() {
502            // Comprehensive fallback implementation for ndarray types
503
504            // f64 types
505            if let (Some(a_array), Some(b_array)) = (
506                a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix1>>(),
507                b.as_any().downcast_ref::<NdarrayWrapper<f64, Ix1>>(),
508            ) {
509                let result = a_array.as_array() * b_array.as_array();
510                return Ok(Box::new(NdarrayWrapper::new(result)));
511            }
512            if let (Some(a_array), Some(b_array)) = (
513                a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>(),
514                b.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>(),
515            ) {
516                let result = a_array.as_array() * b_array.as_array();
517                return Ok(Box::new(NdarrayWrapper::new(result)));
518            }
519            if let (Some(a_array), Some(b_array)) = (
520                a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>(),
521                b.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>(),
522            ) {
523                let result = a_array.as_array() * b_array.as_array();
524                return Ok(Box::new(NdarrayWrapper::new(result)));
525            }
526
527            // f32 types
528            if let (Some(a_array), Some(b_array)) = (
529                a.as_any().downcast_ref::<NdarrayWrapper<f32, Ix1>>(),
530                b.as_any().downcast_ref::<NdarrayWrapper<f32, Ix1>>(),
531            ) {
532                let result = a_array.as_array() * b_array.as_array();
533                return Ok(Box::new(NdarrayWrapper::new(result)));
534            }
535            if let (Some(a_array), Some(b_array)) = (
536                a.as_any().downcast_ref::<NdarrayWrapper<f32, Ix2>>(),
537                b.as_any().downcast_ref::<NdarrayWrapper<f32, Ix2>>(),
538            ) {
539                let result = a_array.as_array() * b_array.as_array();
540                return Ok(Box::new(NdarrayWrapper::new(result)));
541            }
542            if let (Some(a_array), Some(b_array)) = (
543                a.as_any().downcast_ref::<NdarrayWrapper<f32, IxDyn>>(),
544                b.as_any().downcast_ref::<NdarrayWrapper<f32, IxDyn>>(),
545            ) {
546                let result = a_array.as_array() * b_array.as_array();
547                return Ok(Box::new(NdarrayWrapper::new(result)));
548            }
549
550            // i32 types
551            if let (Some(a_array), Some(b_array)) = (
552                a.as_any().downcast_ref::<NdarrayWrapper<i32, Ix1>>(),
553                b.as_any().downcast_ref::<NdarrayWrapper<i32, Ix1>>(),
554            ) {
555                let result = a_array.as_array() * b_array.as_array();
556                return Ok(Box::new(NdarrayWrapper::new(result)));
557            }
558            if let (Some(a_array), Some(b_array)) = (
559                a.as_any().downcast_ref::<NdarrayWrapper<i32, Ix2>>(),
560                b.as_any().downcast_ref::<NdarrayWrapper<i32, Ix2>>(),
561            ) {
562                let result = a_array.as_array() * b_array.as_array();
563                return Ok(Box::new(NdarrayWrapper::new(result)));
564            }
565            if let (Some(a_array), Some(b_array)) = (
566                a.as_any().downcast_ref::<NdarrayWrapper<i32, IxDyn>>(),
567                b.as_any().downcast_ref::<NdarrayWrapper<i32, IxDyn>>(),
568            ) {
569                let result = a_array.as_array() * b_array.as_array();
570                return Ok(Box::new(NdarrayWrapper::new(result)));
571            }
572
573            // i64 types
574            if let (Some(a_array), Some(b_array)) = (
575                a.as_any().downcast_ref::<NdarrayWrapper<i64, Ix1>>(),
576                b.as_any().downcast_ref::<NdarrayWrapper<i64, Ix1>>(),
577            ) {
578                let result = a_array.as_array() * b_array.as_array();
579                return Ok(Box::new(NdarrayWrapper::new(result)));
580            }
581            if let (Some(a_array), Some(b_array)) = (
582                a.as_any().downcast_ref::<NdarrayWrapper<i64, Ix2>>(),
583                b.as_any().downcast_ref::<NdarrayWrapper<i64, Ix2>>(),
584            ) {
585                let result = a_array.as_array() * b_array.as_array();
586                return Ok(Box::new(NdarrayWrapper::new(result)));
587            }
588            if let (Some(a_array), Some(b_array)) = (
589                a.as_any().downcast_ref::<NdarrayWrapper<i64, IxDyn>>(),
590                b.as_any().downcast_ref::<NdarrayWrapper<i64, IxDyn>>(),
591            ) {
592                let result = a_array.as_array() * b_array.as_array();
593                return Ok(Box::new(NdarrayWrapper::new(result)));
594            }
595
596            return Err(OperationError::NotImplemented(
597                "multiply not implemented for these array types".to_string(),
598            ));
599        }
600
601        // Delegate to the implementation
602        let array_ref = implementing_args[0].1;
603
604        let result = array_ref.array_function(
605            &ArrayFunction::new("scirs2::array_protocol::operations::multiply"),
606            &[TypeId::of::<Box<dyn ArrayProtocol>>()],
607            &[Box::new(a.box_clone()), Box::new(b.box_clone())],
608            &HashMap::new(),
609        )?;
610
611        // Try to downcast the result
612        match result.downcast::<Box<dyn ArrayProtocol>>() {
613            Ok(array) => Ok(*array),
614            Err(_) => Err(OperationError::Other(
615                "Failed to downcast result to ArrayProtocol".to_string(),
616            )),
617        }
618    },
619    "scirs2::array_protocol::operations::multiply"
620);
621
622// Reduction operation: sum
623array_function_dispatch!(
624    fn sum(a: &dyn ArrayProtocol, axis: Option<usize>) -> Result<Box<dyn Any>, OperationError> {
625        // Get implementing args
626        let boxed_a = Box::new(a.box_clone());
627        let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a];
628        let implementing_args =
629            get_implementing_args("scirs2::array_protocol::operations::sum", &boxed_args);
630        if implementing_args.is_empty() {
631            // Fallback implementation for ndarray types
632            // Try with Ix2 dimension first (most common case)
633            if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>() {
634                match axis {
635                    Some(ax) => {
636                        let result = a_array.as_array().sum_axis(crate::ndarray::Axis(ax));
637                        return Ok(Box::new(NdarrayWrapper::new(result)));
638                    }
639                    None => {
640                        let result = a_array.as_array().sum();
641                        return Ok(Box::new(result));
642                    }
643                }
644            }
645            // Try with IxDyn dimension (used in tests)
646            else if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>() {
647                match axis {
648                    Some(ax) => {
649                        let result = a_array.as_array().sum_axis(crate::ndarray::Axis(ax));
650                        return Ok(Box::new(NdarrayWrapper::new(result)));
651                    }
652                    None => {
653                        let result = a_array.as_array().sum();
654                        return Ok(Box::new(result));
655                    }
656                }
657            }
658            return Err(OperationError::NotImplemented(
659                "sum not implemented for this array type".to_string(),
660            ));
661        }
662
663        // Delegate to the implementation
664        let mut kwargs = HashMap::new();
665        if let Some(ax) = axis {
666            kwargs.insert("axis".to_string(), Box::new(ax) as Box<dyn Any>);
667        }
668
669        let array_ref = implementing_args[0].1;
670
671        let result = array_ref.array_function(
672            &ArrayFunction::new("scirs2::array_protocol::operations::sum"),
673            &[TypeId::of::<Box<dyn Any>>()],
674            &[Box::new(a.box_clone())],
675            &kwargs,
676        )?;
677
678        Ok(result)
679    },
680    "scirs2::array_protocol::operations::sum"
681);
682
683// Transpose operation
684array_function_dispatch!(
685    fn transpose(a: &dyn ArrayProtocol) -> Result<Box<dyn ArrayProtocol>, OperationError> {
686        // Get implementing args
687        let boxed_a = Box::new(a.box_clone());
688        let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a];
689        let implementing_args =
690            get_implementing_args("scirs2::array_protocol::operations::transpose", &boxed_args);
691        if implementing_args.is_empty() {
692            // Fallback implementation for ndarray types
693            // Try with Ix2 dimension first (most common case)
694            if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>() {
695                let result = a_array.as_array().t().to_owned();
696                return Ok(Box::new(NdarrayWrapper::new(result)));
697            }
698            // Try with IxDyn dimension (used in tests)
699            else if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>() {
700                // For dynamic dimension, we need to check if it's a 2D array
701                let a_dim = a_array.as_array().shape();
702                if a_dim.len() != 2 {
703                    return Err(OperationError::ShapeMismatch(format!(
704                        "Transpose requires a 2D array, got shape: {a_dim:?}"
705                    )));
706                }
707
708                // Create a transposed array
709                let (m, n) = (a_dim[0], a_dim[1]);
710                let mut result = crate::ndarray::Array2::<f64>::zeros((n, m));
711
712                // Fill the transposed array
713                for i in 0..m {
714                    for j in 0..n {
715                        result[[j, i]] = a_array.as_array()[[i, j]];
716                    }
717                }
718
719                return Ok(Box::new(NdarrayWrapper::new(result)));
720            }
721            return Err(OperationError::NotImplemented(
722                "transpose not implemented for this array type".to_string(),
723            ));
724        }
725
726        // Delegate to the implementation
727        let array_ref = implementing_args[0].1;
728
729        let result = array_ref.array_function(
730            &ArrayFunction::new("scirs2::array_protocol::operations::transpose"),
731            &[TypeId::of::<Box<dyn ArrayProtocol>>()],
732            &[Box::new(a.box_clone())],
733            &HashMap::new(),
734        )?;
735
736        // Try to downcast the result
737        match result.downcast::<Box<dyn ArrayProtocol>>() {
738            Ok(array) => Ok(*array),
739            Err(_) => Err(OperationError::Other(
740                "Failed to downcast result to ArrayProtocol".to_string(),
741            )),
742        }
743    },
744    "scirs2::array_protocol::operations::transpose"
745);
746
747// Element-wise application of a function implementation
748#[allow(dead_code)]
749pub fn apply_elementwise<F>(
750    a: &dyn ArrayProtocol,
751    f: F,
752) -> Result<Box<dyn ArrayProtocol>, OperationError>
753where
754    F: Fn(f64) -> f64 + 'static,
755{
756    // Get implementing args
757    let boxed_a = Box::new(a.box_clone());
758    let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a];
759    let implementing_args = get_implementing_args(
760        "scirs2::array_protocol::operations::apply_elementwise",
761        &boxed_args,
762    );
763    if implementing_args.is_empty() {
764        // Fallback implementation for ndarray types
765        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>() {
766            let result = a_array.as_array().mapv(f);
767            return Ok(Box::new(NdarrayWrapper::new(result)));
768        }
769        return Err(OperationError::NotImplemented(
770            "apply_elementwise not implemented for this array type".to_string(),
771        ));
772    }
773
774    // For this operation, we need to handle the function specially
775    // In a real implementation, we would need to serialize the function or use a predefined set
776    // Here we'll just use the fallback implementation for simplicity
777    if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>() {
778        let result = a_array.as_array().mapv(f);
779        Ok(Box::new(NdarrayWrapper::new(result)))
780    } else {
781        Err(OperationError::NotImplemented(
782            "apply_elementwise not implemented for this array type".to_string(),
783        ))
784    }
785}
786
787// Concatenate operation
788array_function_dispatch!(
789    fn concatenate(
790        arrays: &[&dyn ArrayProtocol],
791        axis: usize,
792    ) -> Result<Box<dyn ArrayProtocol>, OperationError> {
793        if arrays.is_empty() {
794            return Err(OperationError::Other(
795                "No arrays provided for concatenation".to_string(),
796            ));
797        }
798
799        // Convert each array to Box<dyn Any>
800        let boxed_arrays: Vec<Box<dyn Any>> = arrays
801            .iter()
802            .map(|&a| Box::new(a.box_clone()) as Box<dyn Any>)
803            .collect();
804
805        let implementing_args = get_implementing_args(
806            "scirs2::array_protocol::operations::concatenate",
807            &boxed_arrays,
808        );
809        if implementing_args.is_empty() {
810            // Fallback implementation for ndarray types
811            // For simplicity, we'll handle just the 2D f64 case
812            let mut ndarray_arrays = Vec::new();
813            for &array in arrays {
814                if let Some(a) = array.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>() {
815                    ndarray_arrays.push(a.as_array().view());
816                } else {
817                    return Err(OperationError::TypeMismatch(
818                        "All arrays must be NdarrayWrapper<f64, Ix2>".to_string(),
819                    ));
820                }
821            }
822
823            let result = match crate::ndarray::stack(crate::ndarray::Axis(axis), &ndarray_arrays) {
824                Ok(arr) => arr,
825                Err(e) => return Err(OperationError::Other(format!("Concatenation failed: {e}"))),
826            };
827
828            return Ok(Box::new(NdarrayWrapper::new(result)));
829        }
830
831        // Delegate to the implementation
832        let array_boxed_clones: Vec<Box<dyn Any>> = arrays
833            .iter()
834            .map(|&a| Box::new(a.box_clone()) as Box<dyn Any>)
835            .collect();
836
837        let mut kwargs = HashMap::new();
838        kwargs.insert(axis.to_string(), Box::new(axis) as Box<dyn Any>);
839
840        let array_ref = implementing_args[0].1;
841
842        let result = array_ref.array_function(
843            &ArrayFunction::new("scirs2::array_protocol::operations::concatenate"),
844            &[TypeId::of::<Box<dyn ArrayProtocol>>()],
845            &array_boxed_clones,
846            &kwargs,
847        )?;
848
849        // Try to downcast the result
850        match result.downcast::<Box<dyn ArrayProtocol>>() {
851            Ok(array) => Ok(*array),
852            Err(_) => Err(OperationError::Other(
853                "Failed to downcast result to ArrayProtocol".to_string(),
854            )),
855        }
856    },
857    "scirs2::array_protocol::operations::concatenate"
858);
859
860// Reshape operation
861array_function_dispatch!(
862    fn reshape(
863        a: &dyn ArrayProtocol,
864        shape: &[usize],
865    ) -> Result<Box<dyn ArrayProtocol>, OperationError> {
866        // Get implementing args
867        let boxed_a = Box::new(a.box_clone());
868        let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a];
869        let implementing_args =
870            get_implementing_args("scirs2::array_protocol::operations::reshape", &boxed_args);
871        if implementing_args.is_empty() {
872            // Fallback implementation for ndarray types
873            // Try with Ix2 dimension first (most common case)
874            if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>() {
875                let result = match a_array.as_array().clone().into_shape_with_order(shape) {
876                    Ok(arr) => arr,
877                    Err(e) => {
878                        return Err(OperationError::ShapeMismatch(format!(
879                            "Reshape failed: {e}"
880                        )))
881                    }
882                };
883                return Ok(Box::new(NdarrayWrapper::new(result)));
884            }
885            // Try with IxDyn dimension (used in tests)
886            else if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>() {
887                let result = match a_array.as_array().clone().into_shape_with_order(shape) {
888                    Ok(arr) => arr,
889                    Err(e) => {
890                        return Err(OperationError::ShapeMismatch(format!(
891                            "Reshape failed: {e}"
892                        )))
893                    }
894                };
895                return Ok(Box::new(NdarrayWrapper::new(result)));
896            }
897            return Err(OperationError::NotImplemented(
898                "reshape not implemented for this array type".to_string(),
899            ));
900        }
901
902        // Delegate to the implementation
903        let mut kwargs = HashMap::new();
904        kwargs.insert(
905            "shape".to_string(),
906            Box::new(shape.to_vec()) as Box<dyn Any>,
907        );
908
909        let array_ref = implementing_args[0].1;
910
911        let result = array_ref.array_function(
912            &ArrayFunction::new("scirs2::array_protocol::operations::reshape"),
913            &[TypeId::of::<Box<dyn ArrayProtocol>>()],
914            &[Box::new(a.box_clone())],
915            &kwargs,
916        )?;
917
918        // Try to downcast the result
919        match result.downcast::<Box<dyn ArrayProtocol>>() {
920            Ok(array) => Ok(*array),
921            Err(_) => Err(OperationError::Other(
922                "Failed to downcast result to ArrayProtocol".to_string(),
923            )),
924        }
925    },
926    "scirs2::array_protocol::operations::reshape"
927);
928
929// Linear algebra operations
930
931// Type alias for SVD return type
932type SVDResult = (
933    Box<dyn ArrayProtocol>,
934    Box<dyn ArrayProtocol>,
935    Box<dyn ArrayProtocol>,
936);
937
938// SVD decomposition operation
939array_function_dispatch!(
940    fn svd(a: &dyn ArrayProtocol) -> Result<SVDResult, OperationError> {
941        // Get implementing args
942        let boxed_a = Box::new(a.box_clone());
943        let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a];
944        let implementing_args =
945            get_implementing_args("scirs2::array_protocol::operations::svd", &boxed_args);
946        if implementing_args.is_empty() {
947            // Fallback implementation for ndarray types
948            if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>() {
949                // Real SVD via the OxiBLAS-backed LAPACK bindings (only
950                // available when the `linalg` feature is enabled, since that
951                // feature gates the optional oxiblas-* dependencies).
952                #[cfg(feature = "linalg")]
953                {
954                    let svd_result = crate::linalg::svd_ndarray(a_array.as_array())
955                        .map_err(|e| OperationError::Other(format!("SVD failed: {e}")))?;
956                    return Ok((
957                        Box::new(NdarrayWrapper::new(svd_result.u)),
958                        Box::new(NdarrayWrapper::new(svd_result.s)),
959                        Box::new(NdarrayWrapper::new(svd_result.vt)),
960                    ));
961                }
962                #[cfg(not(feature = "linalg"))]
963                {
964                    return Err(OperationError::NotImplemented(
965                        "svd requires the `linalg` feature (OxiBLAS-backed decomposition) to be enabled"
966                            .to_string(),
967                    ));
968                }
969            }
970            return Err(OperationError::NotImplemented(
971                "svd not implemented for this array type".to_string(),
972            ));
973        }
974
975        // Delegate to the implementation
976        let array_ref = implementing_args[0].1;
977
978        let result = array_ref.array_function(
979            &ArrayFunction::new("scirs2::array_protocol::operations::svd"),
980            &[TypeId::of::<(
981                Box<dyn ArrayProtocol>,
982                Box<dyn ArrayProtocol>,
983                Box<dyn ArrayProtocol>,
984            )>()],
985            &[Box::new(a.box_clone())],
986            &HashMap::new(),
987        )?;
988
989        // Try to downcast the result
990        match result.downcast::<(
991            Box<dyn ArrayProtocol>,
992            Box<dyn ArrayProtocol>,
993            Box<dyn ArrayProtocol>,
994        )>() {
995            Ok(tuple) => Ok(*tuple),
996            Err(_) => Err(OperationError::Other(
997                "Failed to downcast result to SVD tuple".to_string(),
998            )),
999        }
1000    },
1001    "scirs2::array_protocol::operations::svd"
1002);
1003
1004// Inverse operation
1005array_function_dispatch!(
1006    fn inverse(a: &dyn ArrayProtocol) -> Result<Box<dyn ArrayProtocol>, OperationError> {
1007        // Get implementing args
1008        let boxed_a = Box::new(a.box_clone());
1009        let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a];
1010        let implementing_args =
1011            get_implementing_args("scirs2::array_protocol::operations::inverse", &boxed_args);
1012        if implementing_args.is_empty() {
1013            // Fallback implementation for ndarray types
1014            if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>() {
1015                let (m, n) = a_array.as_array().dim();
1016                if m != n {
1017                    return Err(OperationError::ShapeMismatch(
1018                        "Matrix must be square for inversion".to_string(),
1019                    ));
1020                }
1021
1022                // Real matrix inversion via the OxiBLAS-backed LAPACK
1023                // bindings (only available when the `linalg` feature is
1024                // enabled, since that feature gates the optional oxiblas-*
1025                // dependencies).
1026                #[cfg(feature = "linalg")]
1027                {
1028                    let inv = crate::linalg::inv_ndarray(a_array.as_array()).map_err(|e| {
1029                        OperationError::Other(format!("Matrix inversion failed: {e}"))
1030                    })?;
1031                    return Ok(Box::new(NdarrayWrapper::new(inv)));
1032                }
1033                #[cfg(not(feature = "linalg"))]
1034                {
1035                    return Err(OperationError::NotImplemented(
1036                        "inverse requires the `linalg` feature (OxiBLAS-backed decomposition) to be enabled"
1037                            .to_string(),
1038                    ));
1039                }
1040            }
1041            return Err(OperationError::NotImplemented(
1042                "inverse not implemented for this array type".to_string(),
1043            ));
1044        }
1045
1046        // Delegate to the implementation
1047        let array_ref = implementing_args[0].1;
1048
1049        let result = array_ref.array_function(
1050            &ArrayFunction::new("scirs2::array_protocol::operations::inverse"),
1051            &[TypeId::of::<Box<dyn ArrayProtocol>>()],
1052            &[Box::new(a.box_clone())],
1053            &HashMap::new(),
1054        )?;
1055
1056        // Try to downcast the result
1057        match result.downcast::<Box<dyn ArrayProtocol>>() {
1058            Ok(array) => Ok(*array),
1059            Err(_) => Err(OperationError::Other(
1060                "Failed to downcast result to ArrayProtocol".to_string(),
1061            )),
1062        }
1063    },
1064    "scirs2::array_protocol::operations::inverse"
1065);
1066
1067// Scalar multiplication operation (implemented without macro due to generic constraints)
1068#[allow(dead_code)]
1069pub fn multiply_by_scalar_f64(
1070    a: &dyn ArrayProtocol,
1071    scalar: f64,
1072) -> Result<Box<dyn ArrayProtocol>, OperationError> {
1073    // Get implementing args
1074    let boxed_a = Box::new(a.box_clone());
1075    let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a];
1076    let implementing_args = get_implementing_args(
1077        "scirs2::array_protocol::operations::multiply_by_scalar_f64",
1078        &boxed_args,
1079    );
1080    if implementing_args.is_empty() {
1081        // Fallback implementation for ndarray types
1082        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix1>>() {
1083            let result = a_array.as_array() * scalar;
1084            return Ok(Box::new(NdarrayWrapper::new(result)));
1085        }
1086        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>() {
1087            let result = a_array.as_array() * scalar;
1088            return Ok(Box::new(NdarrayWrapper::new(result)));
1089        }
1090        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>() {
1091            let result = a_array.as_array() * scalar;
1092            return Ok(Box::new(NdarrayWrapper::new(result)));
1093        }
1094        return Err(OperationError::NotImplemented(
1095            "multiply_by_scalar not implemented for this array type".to_string(),
1096        ));
1097    }
1098
1099    // Delegate to the implementation
1100    let mut kwargs = HashMap::new();
1101    kwargs.insert(scalar.to_string(), Box::new(scalar) as Box<dyn Any>);
1102
1103    let array_ref = implementing_args[0].1;
1104
1105    let result = array_ref.array_function(
1106        &ArrayFunction::new("scirs2::array_protocol::operations::multiply_by_scalar_f64"),
1107        &[TypeId::of::<Box<dyn ArrayProtocol>>()],
1108        &[Box::new(a.box_clone())],
1109        &kwargs,
1110    )?;
1111
1112    // Try to downcast the result
1113    match result.downcast::<Box<dyn ArrayProtocol>>() {
1114        Ok(array) => Ok(*array),
1115        Err(_) => Err(OperationError::Other(
1116            "Failed to downcast result to ArrayProtocol".to_string(),
1117        )),
1118    }
1119}
1120
1121// Scalar multiplication for f32
1122#[allow(dead_code)]
1123pub fn multiply_by_scalar_f32(
1124    a: &dyn ArrayProtocol,
1125    scalar: f32,
1126) -> Result<Box<dyn ArrayProtocol>, OperationError> {
1127    // Get implementing args
1128    let boxed_a = Box::new(a.box_clone());
1129    let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a];
1130    let implementing_args = get_implementing_args(
1131        "scirs2::array_protocol::operations::multiply_by_scalar_f32",
1132        &boxed_args,
1133    );
1134    if implementing_args.is_empty() {
1135        // Fallback implementation for ndarray types
1136        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f32, Ix1>>() {
1137            let result = a_array.as_array() * scalar;
1138            return Ok(Box::new(NdarrayWrapper::new(result)));
1139        }
1140        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f32, Ix2>>() {
1141            let result = a_array.as_array() * scalar;
1142            return Ok(Box::new(NdarrayWrapper::new(result)));
1143        }
1144        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f32, IxDyn>>() {
1145            let result = a_array.as_array() * scalar;
1146            return Ok(Box::new(NdarrayWrapper::new(result)));
1147        }
1148        return Err(OperationError::NotImplemented(
1149            "multiply_by_scalar not implemented for this array type".to_string(),
1150        ));
1151    }
1152
1153    // Delegate to the implementation
1154    let mut kwargs = HashMap::new();
1155    kwargs.insert(scalar.to_string(), Box::new(scalar) as Box<dyn Any>);
1156
1157    let array_ref = implementing_args[0].1;
1158
1159    let result = array_ref.array_function(
1160        &ArrayFunction::new("scirs2::array_protocol::operations::multiply_by_scalar_f32"),
1161        &[TypeId::of::<Box<dyn ArrayProtocol>>()],
1162        &[Box::new(a.box_clone())],
1163        &kwargs,
1164    )?;
1165
1166    // Try to downcast the result
1167    match result.downcast::<Box<dyn ArrayProtocol>>() {
1168        Ok(array) => Ok(*array),
1169        Err(_) => Err(OperationError::Other(
1170            "Failed to downcast result to ArrayProtocol".to_string(),
1171        )),
1172    }
1173}
1174
1175// Scalar division for f64
1176#[allow(dead_code)]
1177pub fn divide_by_scalar_f64(
1178    a: &dyn ArrayProtocol,
1179    scalar: f64,
1180) -> Result<Box<dyn ArrayProtocol>, OperationError> {
1181    // Get implementing args
1182    let boxed_a = Box::new(a.box_clone());
1183    let boxed_args: Vec<Box<dyn Any>> = vec![boxed_a];
1184    let implementing_args = get_implementing_args(
1185        "scirs2::array_protocol::operations::divide_by_scalar_f64",
1186        &boxed_args,
1187    );
1188    if implementing_args.is_empty() {
1189        // Fallback implementation for ndarray types
1190        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix1>>() {
1191            let result = a_array.as_array() / scalar;
1192            return Ok(Box::new(NdarrayWrapper::new(result)));
1193        }
1194        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, Ix2>>() {
1195            let result = a_array.as_array() / scalar;
1196            return Ok(Box::new(NdarrayWrapper::new(result)));
1197        }
1198        if let Some(a_array) = a.as_any().downcast_ref::<NdarrayWrapper<f64, IxDyn>>() {
1199            let result = a_array.as_array() / scalar;
1200            return Ok(Box::new(NdarrayWrapper::new(result)));
1201        }
1202        return Err(OperationError::NotImplemented(
1203            "divide_by_scalar not implemented for this array type".to_string(),
1204        ));
1205    }
1206
1207    // Delegate to the implementation
1208    let mut kwargs = HashMap::new();
1209    kwargs.insert(scalar.to_string(), Box::new(scalar) as Box<dyn Any>);
1210
1211    let array_ref = implementing_args[0].1;
1212
1213    let result = array_ref.array_function(
1214        &ArrayFunction::new("scirs2::array_protocol::operations::divide_by_scalar_f64"),
1215        &[TypeId::of::<Box<dyn ArrayProtocol>>()],
1216        &[Box::new(a.box_clone())],
1217        &kwargs,
1218    )?;
1219
1220    // Try to downcast the result
1221    match result.downcast::<Box<dyn ArrayProtocol>>() {
1222        Ok(array) => Ok(*array),
1223        Err(_) => Err(OperationError::Other(
1224            "Failed to downcast result to ArrayProtocol".to_string(),
1225        )),
1226    }
1227}
1228
1229#[cfg(test)]
1230mod tests {
1231    use super::*;
1232    use crate::array_protocol::{self, NdarrayWrapper};
1233    use ::ndarray::{array, Array2};
1234
1235    #[test]
1236    fn test_operations_with_ndarray() {
1237        use ::ndarray::array;
1238
1239        // Initialize the array protocol system
1240        array_protocol::init();
1241
1242        // Create regular ndarrays
1243        let a = Array2::<f64>::eye(3);
1244        let b = Array2::<f64>::ones((3, 3));
1245
1246        // Wrap them in NdarrayWrapper
1247        let wrapped_a = NdarrayWrapper::new(a.clone());
1248        let wrapped_b = NdarrayWrapper::new(b.clone());
1249
1250        // Test matrix multiplication
1251        let matmul_result = matmul(&wrapped_a, &wrapped_b).expect("matmul should succeed");
1252        let matmul_array = matmul_result
1253            .as_any()
1254            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1255            .expect("matmul should return a NdarrayWrapper<f64, Ix2>");
1256        assert_eq!(matmul_array.as_array(), &a.dot(&b));
1257
1258        // Test addition
1259        let add_result = add(&wrapped_a, &wrapped_b).expect("add should succeed");
1260        let add_array = add_result
1261            .as_any()
1262            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1263            .expect("add should return a NdarrayWrapper<f64, Ix2>");
1264        assert_eq!(add_array.as_array(), &(a.clone() + b.clone()));
1265
1266        // Test multiplication
1267        let multiply_result = multiply(&wrapped_a, &wrapped_b).expect("multiply should succeed");
1268        let multiply_array = multiply_result
1269            .as_any()
1270            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1271            .expect("multiply should return a NdarrayWrapper<f64, Ix2>");
1272        assert_eq!(multiply_array.as_array(), &(a.clone() * b.clone()));
1273
1274        // Test sum
1275        let sum_result = sum(&wrapped_a, None).expect("sum should succeed");
1276        let sum_value = sum_result
1277            .downcast_ref::<f64>()
1278            .expect("sum should return an f64");
1279        assert_eq!(*sum_value, a.sum());
1280
1281        // Test transpose
1282        let transpose_result = transpose(&wrapped_a).expect("transpose should succeed");
1283        let transpose_array = transpose_result
1284            .as_any()
1285            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1286            .expect("transpose should return a NdarrayWrapper<f64, Ix2>");
1287        assert_eq!(transpose_array.as_array(), &a.t().to_owned());
1288
1289        // Test reshape. `reshape`'s dispatch passes `shape` as a `Vec<usize>`
1290        // (a runtime-length container), and `ndarray`'s `IntoDimension` maps
1291        // `Vec<usize>` to `IxDyn` regardless of how many elements it holds —
1292        // not to a fixed-size `Ix1`/`Ix2`/etc., even when the vec happens to
1293        // have exactly 1 element. So the result is a `NdarrayWrapper<f64,
1294        // IxDyn>`, not `NdarrayWrapper<f64, Ix1>`.
1295        let c = array![[1., 2., 3.], [4., 5., 6.]];
1296        let wrapped_c = NdarrayWrapper::new(c.clone());
1297        let reshape_result = reshape(&wrapped_c, &[6]).expect("reshape should succeed");
1298        let result_array = reshape_result
1299            .as_any()
1300            .downcast_ref::<NdarrayWrapper<f64, IxDyn>>()
1301            .expect("reshape should return a NdarrayWrapper<f64, IxDyn>");
1302        let expected = c
1303            .clone()
1304            .into_shape_with_order(6)
1305            .expect("Operation failed")
1306            .into_dyn();
1307        assert_eq!(result_array.as_array(), &expected);
1308    }
1309
1310    /// Regression test for the `get_implementing_args` dispatch bug: before
1311    /// the fix, `implementing_args` was never empty for `NdarrayWrapper`
1312    /// arguments (since `NdarrayWrapper` boxed as `Box<dyn ArrayProtocol>`
1313    /// always downcast successfully), so `subtract`/`multiply_by_scalar_f64`/
1314    /// `divide_by_scalar_f64` always delegated straight into
1315    /// `NdarrayWrapper::array_function` — which doesn't have a match arm for
1316    /// any of these operations — and so always returned
1317    /// `Err(OperationError::NotImplemented(..))` instead of ever reaching
1318    /// the correct hand-written fallback code below. Uses non-constant data
1319    /// so a fabricated all-same-value result couldn't pass by accident.
1320    #[test]
1321    fn test_subtract_and_scalar_ops_reach_fallback() {
1322        array_protocol::init();
1323
1324        let a = Array2::from_shape_fn((3, 3), |(i, j)| (i * 3 + j) as f64 + 1.0);
1325        let b = Array2::from_shape_fn((3, 3), |(i, j)| (i as f64) * 0.5 - (j as f64) * 0.25);
1326
1327        let wrapped_a = NdarrayWrapper::new(a.clone());
1328        let wrapped_b = NdarrayWrapper::new(b.clone());
1329
1330        let subtract_result = subtract(&wrapped_a, &wrapped_b).expect("subtract should succeed");
1331        let subtract_array = subtract_result
1332            .as_any()
1333            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1334            .expect("subtract should return a NdarrayWrapper<f64, Ix2>");
1335        assert_eq!(subtract_array.as_array(), &(a.clone() - b.clone()));
1336
1337        let scaled =
1338            multiply_by_scalar_f64(&wrapped_a, 2.5).expect("multiply_by_scalar_f64 should succeed");
1339        let scaled_array = scaled
1340            .as_any()
1341            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1342            .expect("multiply_by_scalar_f64 should return a NdarrayWrapper<f64, Ix2>");
1343        assert_eq!(scaled_array.as_array(), &(a.clone() * 2.5));
1344
1345        let divided =
1346            divide_by_scalar_f64(&wrapped_a, 4.0).expect("divide_by_scalar_f64 should succeed");
1347        let divided_array = divided
1348            .as_any()
1349            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1350            .expect("divide_by_scalar_f64 should return a NdarrayWrapper<f64, Ix2>");
1351        assert_eq!(divided_array.as_array(), &(a.clone() / 4.0));
1352    }
1353
1354    /// Regression test for the fake `svd`/`inverse` placeholders: before the
1355    /// fix, both always returned an identity/ones "decomposition" regardless
1356    /// of the input, and — until the `get_implementing_args` dispatch fix —
1357    /// were also unreachable dead code. Uses a non-symmetric, non-constant
1358    /// matrix so an identity/ones fabrication provably fails reconstruction.
1359    #[cfg(feature = "linalg")]
1360    #[test]
1361    fn test_svd_reconstructs_original_matrix() {
1362        array_protocol::init();
1363
1364        let a = Array2::from_shape_fn((4, 3), |(i, j)| ((i + 1) as f64) * 1.7 - (j as f64) * 0.9);
1365        let wrapped_a = NdarrayWrapper::new(a.clone());
1366
1367        let (u, s, vt) = svd(&wrapped_a).expect("svd should succeed");
1368        let u = u
1369            .as_any()
1370            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1371            .expect("U should be a NdarrayWrapper<f64, Ix2>")
1372            .as_array()
1373            .clone();
1374        let s = s
1375            .as_any()
1376            .downcast_ref::<NdarrayWrapper<f64, ::ndarray::Ix1>>()
1377            .expect("S should be a NdarrayWrapper<f64, Ix1>")
1378            .as_array()
1379            .clone();
1380        let vt = vt
1381            .as_any()
1382            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1383            .expect("Vt should be a NdarrayWrapper<f64, Ix2>")
1384            .as_array()
1385            .clone();
1386
1387        // A real SVD must NOT be the identity/ones placeholder: the old
1388        // placeholder returned `s` as all-ones regardless of input, which a
1389        // non-constant, non-orthogonal input matrix like `a` can never
1390        // actually produce.
1391        assert!(
1392            s.iter().any(|&v| (v - 1.0).abs() > 1e-6),
1393            "singular values look like the old all-ones placeholder: {s:?}"
1394        );
1395
1396        // Build Sigma with whatever shape U/Vt actually came back with
1397        // (implementations differ on full vs. economy SVD), then verify the
1398        // defining reconstruction property A = U * Sigma * V^T.
1399        let mut sigma = Array2::<f64>::zeros((u.ncols(), vt.nrows()));
1400        for i in 0..s.len() {
1401            sigma[[i, i]] = s[i];
1402        }
1403        let reconstructed = u.dot(&sigma).dot(&vt);
1404
1405        for i in 0..a.nrows() {
1406            for j in 0..a.ncols() {
1407                assert!(
1408                    (reconstructed[[i, j]] - a[[i, j]]).abs() < 1e-8,
1409                    "SVD reconstruction mismatch at ({i}, {j}): {reconstructed_val} vs {orig_val}",
1410                    reconstructed_val = reconstructed[[i, j]],
1411                    orig_val = a[[i, j]]
1412                );
1413            }
1414        }
1415    }
1416
1417    #[cfg(feature = "linalg")]
1418    #[test]
1419    fn test_inverse_is_a_real_inverse() {
1420        array_protocol::init();
1421
1422        // A non-symmetric, non-trivial invertible matrix.
1423        let a = array![[4.0, 3.0, 0.0], [2.0, 5.0, 1.0], [1.0, 0.0, 3.0]];
1424        let wrapped_a = NdarrayWrapper::new(a.clone());
1425
1426        let inv_result = inverse(&wrapped_a).expect("inverse should succeed");
1427        let inv_array = inv_result
1428            .as_any()
1429            .downcast_ref::<NdarrayWrapper<f64, Ix2>>()
1430            .expect("inverse should return a NdarrayWrapper<f64, Ix2>")
1431            .as_array()
1432            .clone();
1433
1434        // The old placeholder always returned the identity matrix, which
1435        // would trivially (and wrongly) satisfy `a.dot(&identity) == a`, not
1436        // `a.dot(&inv) == identity` — assert the real property instead.
1437        let product = a.dot(&inv_array);
1438        for i in 0..3 {
1439            for j in 0..3 {
1440                let expected = if i == j { 1.0 } else { 0.0 };
1441                assert!(
1442                    (product[[i, j]] - expected).abs() < 1e-8,
1443                    "A * inv(A) mismatch at ({i}, {j}): {actual}",
1444                    actual = product[[i, j]]
1445                );
1446            }
1447        }
1448    }
1449}
scirs2_core/array_protocol/operations.rs

scirs2_core/array_protocol/
operations.rs
