xref: /libCEED/rust/libceed/src/vector.rs (revision 7e7773b5b3fd61233915b5fe058a2730d1639517)
1 // Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights
3 // reserved. See files LICENSE and NOTICE for details.
4 //
5 // This file is part of CEED, a collection of benchmarks, miniapps, software
6 // libraries and APIs for efficient high-order finite element and spectral
7 // element discretizations for exascale applications. For more information and
8 // source code availability see http://github.com/ceed.
9 //
10 // The CEED research is supported by the Exascale Computing Project 17-SC-20-SC,
11 // a collaborative effort of two U.S. Department of Energy organizations (Office
12 // of Science and the National Nuclear Security Administration) responsible for
13 // the planning and preparation of a capable exascale ecosystem, including
14 // software, applications, hardware, advanced system engineering and early
15 // testbed platforms, in support of the nation's exascale computing imperative
16 
17 //! A Ceed Vector constitutes the main data structure and serves as input/output
18 //! for Ceed Operators.
19 
20 use std::{
21     ops::{Deref, DerefMut},
22     os::raw::c_char,
23 };
24 
25 use crate::prelude::*;
26 
27 // -----------------------------------------------------------------------------
28 // CeedVector option
29 // -----------------------------------------------------------------------------
30 #[derive(Debug, Clone, Copy)]
31 pub enum VectorOpt<'a> {
32     Some(&'a Vector<'a>),
33     Active,
34     None,
35 }
36 /// Construct a VectorOpt reference from a Vector reference
37 impl<'a> From<&'a Vector<'_>> for VectorOpt<'a> {
38     fn from(vec: &'a Vector) -> Self {
39         debug_assert!(vec.ptr != unsafe { bind_ceed::CEED_VECTOR_NONE });
40         debug_assert!(vec.ptr != unsafe { bind_ceed::CEED_VECTOR_ACTIVE });
41         Self::Some(vec)
42     }
43 }
44 impl<'a> VectorOpt<'a> {
45     /// Transform a Rust libCEED VectorOpt into C libCEED CeedVector
46     pub(crate) fn to_raw(self) -> bind_ceed::CeedVector {
47         match self {
48             Self::Some(vec) => vec.ptr,
49             Self::None => unsafe { bind_ceed::CEED_VECTOR_NONE },
50             Self::Active => unsafe { bind_ceed::CEED_VECTOR_ACTIVE },
51         }
52     }
53 }
54 
55 // -----------------------------------------------------------------------------
56 // CeedVector context wrapper
57 // -----------------------------------------------------------------------------
58 #[derive(Debug)]
59 pub struct Vector<'a> {
60     ceed: &'a crate::Ceed,
61     pub(crate) ptr: bind_ceed::CeedVector,
62 }
63 impl From<&'_ Vector<'_>> for bind_ceed::CeedVector {
64     fn from(vec: &Vector) -> Self {
65         vec.ptr
66     }
67 }
68 
69 // -----------------------------------------------------------------------------
70 // Destructor
71 // -----------------------------------------------------------------------------
72 impl<'a> Drop for Vector<'a> {
73     fn drop(&mut self) {
74         let not_none_and_active = self.ptr != unsafe { bind_ceed::CEED_VECTOR_NONE }
75             && self.ptr != unsafe { bind_ceed::CEED_VECTOR_ACTIVE };
76 
77         if not_none_and_active {
78             unsafe { bind_ceed::CeedVectorDestroy(&mut self.ptr) };
79         }
80     }
81 }
82 
83 // -----------------------------------------------------------------------------
84 // Display
85 // -----------------------------------------------------------------------------
86 impl<'a> fmt::Display for Vector<'a> {
87     /// View a Vector
88     ///
89     /// ```
90     /// # use libceed::prelude::*;
91     /// # let ceed = libceed::Ceed::default_init();
92     /// let vec = libceed::vector::Vector::from_slice(&ceed, &[1., 2., 3.]).unwrap();
93     /// assert_eq!(
94     ///     vec.to_string(),
95     ///     "CeedVector length 3
96     ///     1.00000000
97     ///     2.00000000
98     ///     3.00000000
99     /// "
100     /// )
101     /// ```
102     fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
103         let mut ptr = std::ptr::null_mut();
104         let mut sizeloc = crate::MAX_BUFFER_LENGTH;
105         let format = CString::new("%12.8f").expect("CString::new failed");
106         let format_c: *const c_char = format.into_raw();
107         let cstring = unsafe {
108             let file = bind_ceed::open_memstream(&mut ptr, &mut sizeloc);
109             bind_ceed::CeedVectorView(self.ptr, format_c, file);
110             bind_ceed::fclose(file);
111             CString::from_raw(ptr)
112         };
113         cstring.to_string_lossy().fmt(f)
114     }
115 }
116 
117 // -----------------------------------------------------------------------------
118 // Implementations
119 // -----------------------------------------------------------------------------
120 impl<'a> Vector<'a> {
121     // Constructors
122     pub fn create(ceed: &'a crate::Ceed, n: usize) -> crate::Result<Self> {
123         let n = i32::try_from(n).unwrap();
124         let mut ptr = std::ptr::null_mut();
125         let ierr = unsafe { bind_ceed::CeedVectorCreate(ceed.ptr, n, &mut ptr) };
126         ceed.check_error(ierr)?;
127         Ok(Self { ceed, ptr })
128     }
129 
130     pub(crate) fn from_raw(
131         ceed: &'a crate::Ceed,
132         ptr: bind_ceed::CeedVector,
133     ) -> crate::Result<Self> {
134         Ok(Self { ceed, ptr })
135     }
136 
137     /// Create a Vector from a slice
138     ///
139     /// # arguments
140     ///
141     /// * `slice` - values to initialize vector with
142     ///
143     /// ```
144     /// # use libceed::prelude::*;
145     /// # let ceed = libceed::Ceed::default_init();
146     /// let vec = vector::Vector::from_slice(&ceed, &[1., 2., 3.]).unwrap();
147     /// assert_eq!(vec.length(), 3, "Incorrect length from slice");
148     /// ```
149     pub fn from_slice(ceed: &'a crate::Ceed, v: &[crate::Scalar]) -> crate::Result<Self> {
150         let mut x = Self::create(ceed, v.len())?;
151         x.set_slice(v)?;
152         Ok(x)
153     }
154 
155     /// Create a Vector from a mutable array reference
156     ///
157     /// # arguments
158     ///
159     /// * `slice` - values to initialize vector with
160     ///
161     /// ```
162     /// # use libceed::prelude::*;
163     /// # let ceed = libceed::Ceed::default_init();
164     /// let mut rust_vec = vec![1., 2., 3.];
165     /// let vec = libceed::vector::Vector::from_array(&ceed, &mut rust_vec).unwrap();
166     ///
167     /// assert_eq!(vec.length(), 3, "Incorrect length from slice");
168     /// ```
169     pub fn from_array(ceed: &'a crate::Ceed, v: &mut [crate::Scalar]) -> crate::Result<Self> {
170         let x = Self::create(ceed, v.len())?;
171         let (host, user_pointer) = (
172             crate::MemType::Host as bind_ceed::CeedMemType,
173             crate::CopyMode::UsePointer as bind_ceed::CeedCopyMode,
174         );
175         let v = v.as_ptr() as *mut crate::Scalar;
176         let ierr = unsafe { bind_ceed::CeedVectorSetArray(x.ptr, host, user_pointer, v) };
177         ceed.check_error(ierr)?;
178         Ok(x)
179     }
180 
181     /// Returns the length of a CeedVector
182     ///
183     /// ```
184     /// # use libceed::prelude::*;
185     /// # let ceed = libceed::Ceed::default_init();
186     /// let vec = ceed.vector(10).unwrap();
187     ///
188     /// let n = vec.length();
189     /// assert_eq!(n, 10, "Incorrect length");
190     /// ```
191     pub fn length(&self) -> usize {
192         let mut n = 0;
193         unsafe { bind_ceed::CeedVectorGetLength(self.ptr, &mut n) };
194         usize::try_from(n).unwrap()
195     }
196 
197     /// Returns the length of a CeedVector
198     ///
199     /// ```
200     /// # use libceed::prelude::*;
201     /// # let ceed = libceed::Ceed::default_init();
202     /// let vec = ceed.vector(10).unwrap();
203     /// assert_eq!(vec.len(), 10, "Incorrect length");
204     /// ```
205     pub fn len(&self) -> usize {
206         self.length()
207     }
208 
209     /// Set the CeedVector to a constant value
210     ///
211     /// # arguments
212     ///
213     /// * `val` - Value to be used
214     ///
215     /// ```
216     /// # use libceed::prelude::*;
217     /// # let ceed = libceed::Ceed::default_init();
218     /// let len = 10;
219     /// let mut vec = ceed.vector(len).unwrap();
220     ///
221     /// let val = 42.0;
222     /// vec.set_value(val).unwrap();
223     ///
224     /// vec.view().iter().for_each(|v| {
225     ///     assert_eq!(*v, val, "Value not set correctly");
226     /// });
227     /// ```
228     pub fn set_value(&mut self, value: crate::Scalar) -> crate::Result<i32> {
229         let ierr = unsafe { bind_ceed::CeedVectorSetValue(self.ptr, value) };
230         self.ceed.check_error(ierr)
231     }
232 
233     /// Set values from a slice of the same length
234     ///
235     /// # arguments
236     ///
237     /// * `slice` - values to into self; length must match
238     ///
239     /// ```
240     /// # use libceed::prelude::*;
241     /// # let ceed = libceed::Ceed::default_init();
242     /// let mut vec = ceed.vector(4).unwrap();
243     /// vec.set_slice(&[10., 11., 12., 13.]).unwrap();
244     ///
245     /// vec.view().iter().enumerate().for_each(|(i, v)| {
246     ///     assert_eq!(*v, 10. + i as Scalar, "Slice not set correctly");
247     /// });
248     /// ```
249     pub fn set_slice(&mut self, slice: &[crate::Scalar]) -> crate::Result<i32> {
250         assert_eq!(self.length(), slice.len());
251         let (host, copy_mode) = (
252             crate::MemType::Host as bind_ceed::CeedMemType,
253             crate::CopyMode::CopyValues as bind_ceed::CeedCopyMode,
254         );
255         let ierr = unsafe {
256             bind_ceed::CeedVectorSetArray(
257                 self.ptr,
258                 host,
259                 copy_mode,
260                 slice.as_ptr() as *mut crate::Scalar,
261             )
262         };
263         self.ceed.check_error(ierr)
264     }
265 
266     /// Sync the CeedVector to a specified memtype
267     ///
268     /// # arguments
269     ///
270     /// * `mtype` - Memtype to be synced
271     ///
272     /// ```
273     /// # use libceed::prelude::*;
274     /// # let ceed = libceed::Ceed::default_init();
275     /// let len = 10;
276     /// let mut vec = ceed.vector(len).unwrap();
277     ///
278     /// let val = 42.0;
279     /// vec.set_value(val);
280     /// vec.sync(MemType::Host).unwrap();
281     ///
282     /// vec.view().iter().for_each(|v| {
283     ///     assert_eq!(*v, val, "Value not set correctly");
284     /// });
285     /// ```
286     pub fn sync(&self, mtype: crate::MemType) -> crate::Result<i32> {
287         let ierr =
288             unsafe { bind_ceed::CeedVectorSyncArray(self.ptr, mtype as bind_ceed::CeedMemType) };
289         self.ceed.check_error(ierr)
290     }
291 
292     /// Create an immutable view
293     ///
294     /// ```
295     /// # use libceed::prelude::*;
296     /// # let ceed = libceed::Ceed::default_init();
297     /// let vec = ceed.vector_from_slice(&[10., 11., 12., 13.]).unwrap();
298     ///
299     /// let v = vec.view();
300     /// assert_eq!(v[0..2], [10., 11.]);
301     ///
302     /// // It is valid to have multiple immutable views
303     /// let w = vec.view();
304     /// assert_eq!(v[1..], w[1..]);
305     /// ```
306     pub fn view(&self) -> VectorView {
307         VectorView::new(self)
308     }
309 
310     /// Create an mutable view
311     ///
312     /// ```
313     /// # use libceed::prelude::*;
314     /// # let ceed = libceed::Ceed::default_init();
315     /// let mut vec = ceed.vector_from_slice(&[10., 11., 12., 13.]).unwrap();
316     ///
317     /// {
318     ///     let mut v = vec.view_mut();
319     ///     v[2] = 9.;
320     /// }
321     ///
322     /// let w = vec.view();
323     /// assert_eq!(w[2], 9., "View did not mutate data");
324     /// ```
325     pub fn view_mut(&mut self) -> VectorViewMut {
326         VectorViewMut::new(self)
327     }
328 
329     /// Return the norm of a CeedVector
330     ///
331     /// # arguments
332     ///
333     /// * `ntype` - Norm type One, Two, or Max
334     ///
335     /// ```
336     /// # use libceed::prelude::*;
337     /// # let ceed = libceed::Ceed::default_init();
338     /// let vec = ceed.vector_from_slice(&[1., 2., 3., 4.]).unwrap();
339     ///
340     /// let max_norm = vec.norm(NormType::Max).unwrap();
341     /// assert_eq!(max_norm, 4.0, "Incorrect Max norm");
342     ///
343     /// let l1_norm = vec.norm(NormType::One).unwrap();
344     /// assert_eq!(l1_norm, 10., "Incorrect L1 norm");
345     ///
346     /// let l2_norm = vec.norm(NormType::Two).unwrap();
347     /// assert!((l2_norm - 5.477) < 1e-3, "Incorrect L2 norm");
348     /// ```
349     pub fn norm(&self, ntype: crate::NormType) -> crate::Result<crate::Scalar> {
350         let mut res: crate::Scalar = 0.0;
351         let ierr = unsafe {
352             bind_ceed::CeedVectorNorm(self.ptr, ntype as bind_ceed::CeedNormType, &mut res)
353         };
354         self.ceed.check_error(ierr)?;
355         Ok(res)
356     }
357 
358     /// Compute x = alpha x for a CeedVector
359     ///
360     /// # arguments
361     ///
362     /// * `alpha` - scaling factor
363     ///
364     /// ```
365     /// # use libceed::prelude::*;
366     /// # let ceed = libceed::Ceed::default_init();
367     /// let mut vec = ceed.vector_from_slice(&[0., 1., 2., 3., 4.]).unwrap();
368     ///
369     /// vec = vec.scale(-1.0).unwrap();
370     /// vec.view().iter().enumerate().for_each(|(i, &v)| {
371     ///     assert_eq!(v, -(i as Scalar), "Value not set correctly");
372     /// });
373     /// ```
374     #[allow(unused_mut)]
375     pub fn scale(mut self, alpha: crate::Scalar) -> crate::Result<Self> {
376         let ierr = unsafe { bind_ceed::CeedVectorScale(self.ptr, alpha) };
377         self.ceed.check_error(ierr)?;
378         Ok(self)
379     }
380 
381     /// Compute y = alpha x + y for a pair of CeedVectors
382     ///
383     /// # arguments
384     ///
385     /// * `alpha` - scaling factor
386     /// * `x`     - second vector, must be different than self
387     ///
388     /// ```
389     /// # use libceed::prelude::*;
390     /// # let ceed = libceed::Ceed::default_init();
391     /// let x = ceed.vector_from_slice(&[0., 1., 2., 3., 4.]).unwrap();
392     /// let mut y = ceed.vector_from_slice(&[0., 1., 2., 3., 4.]).unwrap();
393     ///
394     /// y = y.axpy(-0.5, &x).unwrap();
395     /// y.view().iter().enumerate().for_each(|(i, &v)| {
396     ///     assert_eq!(v, (i as Scalar) / 2.0, "Value not set correctly");
397     /// });
398     /// ```
399     #[allow(unused_mut)]
400     pub fn axpy(mut self, alpha: crate::Scalar, x: &crate::Vector) -> crate::Result<Self> {
401         let ierr = unsafe { bind_ceed::CeedVectorAXPY(self.ptr, alpha, x.ptr) };
402         self.ceed.check_error(ierr)?;
403         Ok(self)
404     }
405 
406     /// Compute the pointwise multiplication w = x .* y for CeedVectors
407     ///
408     /// # arguments
409     ///
410     /// * `x` - first vector for product
411     /// * `y` - second vector for product
412     ///
413     /// ```
414     /// # use libceed::prelude::*;
415     /// # let ceed = libceed::Ceed::default_init();
416     /// let mut w = ceed.vector_from_slice(&[0., 1., 2., 3., 4.]).unwrap();
417     /// let x = ceed.vector_from_slice(&[0., 1., 2., 3., 4.]).unwrap();
418     /// let y = ceed.vector_from_slice(&[0., 1., 2., 3., 4.]).unwrap();
419     ///
420     /// w = w.pointwise_mult(&x, &y).unwrap();
421     /// w.view().iter().enumerate().for_each(|(i, &v)| {
422     ///     assert_eq!(v, (i as Scalar).powf(2.0), "Value not set correctly");
423     /// });
424     /// ```
425     #[allow(unused_mut)]
426     pub fn pointwise_mult(mut self, x: &crate::Vector, y: &crate::Vector) -> crate::Result<Self> {
427         let ierr = unsafe { bind_ceed::CeedVectorPointwiseMult(self.ptr, x.ptr, y.ptr) };
428         self.ceed.check_error(ierr)?;
429         Ok(self)
430     }
431 
432     /// Compute the pointwise multiplication w = w .* x for CeedVectors
433     ///
434     /// # arguments
435     ///
436     /// * `x` - second vector for product
437     ///
438     /// ```
439     /// # use libceed::prelude::*;
440     /// # let ceed = libceed::Ceed::default_init();
441     /// let mut w = ceed.vector_from_slice(&[0., 1., 2., 3., 4.]).unwrap();
442     /// let x = ceed.vector_from_slice(&[0., 1., 2., 3., 4.]).unwrap();
443     ///
444     /// w = w.pointwise_scale(&x).unwrap();
445     /// w.view().iter().enumerate().for_each(|(i, &v)| {
446     ///     assert_eq!(v, (i as Scalar).powf(2.0), "Value not set correctly");
447     /// });
448     /// ```
449     #[allow(unused_mut)]
450     pub fn pointwise_scale(mut self, x: &crate::Vector) -> crate::Result<Self> {
451         let ierr = unsafe { bind_ceed::CeedVectorPointwiseMult(self.ptr, self.ptr, x.ptr) };
452         self.ceed.check_error(ierr)?;
453         Ok(self)
454     }
455 
456     /// Compute the pointwise multiplication w = w .* w for a CeedVector
457     ///
458     /// ```
459     /// # use libceed::prelude::*;
460     /// # let ceed = libceed::Ceed::default_init();
461     /// let mut w = ceed.vector_from_slice(&[0., 1., 2., 3., 4.]).unwrap();
462     ///
463     /// w = w.pointwise_square().unwrap();
464     /// w.view().iter().enumerate().for_each(|(i, &v)| {
465     ///     assert_eq!(v, (i as Scalar).powf(2.0), "Value not set correctly");
466     /// });
467     /// ```
468     #[allow(unused_mut)]
469     pub fn pointwise_square(mut self) -> crate::Result<Self> {
470         let ierr = unsafe { bind_ceed::CeedVectorPointwiseMult(self.ptr, self.ptr, self.ptr) };
471         self.ceed.check_error(ierr)?;
472         Ok(self)
473     }
474 }
475 
476 // -----------------------------------------------------------------------------
477 // Vector Viewer
478 // -----------------------------------------------------------------------------
479 /// A (host) view of a Vector with Deref to slice.  We can't make
480 /// Vector itself Deref to slice because we can't handle the drop to
481 /// call bind_ceed::CeedVectorRestoreArrayRead().
482 #[derive(Debug)]
483 pub struct VectorView<'a> {
484     vec: &'a Vector<'a>,
485     array: *const crate::Scalar,
486 }
487 
488 impl<'a> VectorView<'a> {
489     /// Construct a VectorView from a Vector reference
490     fn new(vec: &'a Vector) -> Self {
491         let mut array = std::ptr::null();
492         unsafe {
493             bind_ceed::CeedVectorGetArrayRead(
494                 vec.ptr,
495                 crate::MemType::Host as bind_ceed::CeedMemType,
496                 &mut array,
497             );
498         }
499         Self {
500             vec: vec,
501             array: array,
502         }
503     }
504 }
505 
506 // Destructor
507 impl<'a> Drop for VectorView<'a> {
508     fn drop(&mut self) {
509         unsafe {
510             bind_ceed::CeedVectorRestoreArrayRead(self.vec.ptr, &mut self.array);
511         }
512     }
513 }
514 
515 // Data access
516 impl<'a> Deref for VectorView<'a> {
517     type Target = [crate::Scalar];
518     fn deref(&self) -> &[crate::Scalar] {
519         unsafe { std::slice::from_raw_parts(self.array, self.vec.len()) }
520     }
521 }
522 
523 // Viewing
524 impl<'a> fmt::Display for VectorView<'a> {
525     fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
526         write!(f, "VectorView({:?})", self.deref())
527     }
528 }
529 
530 // -----------------------------------------------------------------------------
531 // Vector Viewer Mutable
532 // -----------------------------------------------------------------------------
533 /// A mutable (host) view of a Vector with Deref to slice.
534 #[derive(Debug)]
535 pub struct VectorViewMut<'a> {
536     vec: &'a Vector<'a>,
537     array: *mut crate::Scalar,
538 }
539 
540 impl<'a> VectorViewMut<'a> {
541     /// Construct a VectorViewMut from a Vector reference
542     fn new(vec: &'a mut Vector) -> Self {
543         let mut ptr = std::ptr::null_mut();
544         unsafe {
545             bind_ceed::CeedVectorGetArray(
546                 vec.ptr,
547                 crate::MemType::Host as bind_ceed::CeedMemType,
548                 &mut ptr,
549             );
550         }
551         Self {
552             vec: vec,
553             array: ptr,
554         }
555     }
556 }
557 
558 // Destructor
559 impl<'a> Drop for VectorViewMut<'a> {
560     fn drop(&mut self) {
561         unsafe {
562             bind_ceed::CeedVectorRestoreArray(self.vec.ptr, &mut self.array);
563         }
564     }
565 }
566 
567 // Data access
568 impl<'a> Deref for VectorViewMut<'a> {
569     type Target = [crate::Scalar];
570     fn deref(&self) -> &[crate::Scalar] {
571         unsafe { std::slice::from_raw_parts(self.array, self.vec.len()) }
572     }
573 }
574 
575 // Mutable data access
576 impl<'a> DerefMut for VectorViewMut<'a> {
577     fn deref_mut(&mut self) -> &mut [crate::Scalar] {
578         unsafe { std::slice::from_raw_parts_mut(self.array, self.vec.len()) }
579     }
580 }
581 
582 // Viewing
583 impl<'a> fmt::Display for VectorViewMut<'a> {
584     fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
585         write!(f, "VectorViewMut({:?})", self.deref())
586     }
587 }
588 
589 // -----------------------------------------------------------------------------
590