fastsim_core/
utils.rs

1//! Module containing miscellaneous utility functions.
2
3#[cfg(feature = "default")]
4use curl::easy::Easy;
5#[cfg(feature = "default")]
6use directories::ProjectDirs;
7use itertools::Itertools;
8use lazy_static::lazy_static;
9use ndarray::*;
10use regex::Regex;
11use std::collections::HashSet;
12#[cfg(feature = "default")]
13use std::io::Write;
14
15use crate::imports::*;
16#[cfg(feature = "pyo3")]
17use crate::pyo3imports::*;
18
19#[cfg(test)]
20pub fn resources_path() -> PathBuf {
21    let pb = PathBuf::from("../../python/fastsim/resources");
22    assert!(pb.exists());
23    pb
24}
25
26/// Error message for when user attempts to set value in a nested struct.
27pub const NESTED_STRUCT_ERR: &str = "Setting field value on nested struct not allowed.
28Assign nested struct to own variable, run the `reset_orphaned` method, and then
29modify field value. Then set the nested struct back inside containing struct.";
30
31pub fn diff(x: &Array1<f64>) -> Array1<f64> {
32    concatenate(
33        Axis(0),
34        &[
35            array![0.0].view(),
36            (&x.slice(s![1..]) - &x.slice(s![..-1])).view(),
37        ],
38    )
39    .unwrap()
40}
41
42/// Returns a new array with a constant added starting at xs\[i\] to the end. Values prior to xs\[i\] are unchanged.
43pub fn add_from(xs: &Array1<f64>, i: usize, val: f64) -> Array1<f64> {
44    let mut ys = Array1::zeros(xs.len());
45    for idx in 0..xs.len() {
46        if idx >= i {
47            ys[idx] = xs[idx] + val;
48        } else {
49            ys[idx] = xs[idx];
50        }
51    }
52    ys
53}
54
55/// Return first index of `arr` greater than `cut`
56pub fn first_grtr(arr: &[f64], cut: f64) -> Option<usize> {
57    let len = arr.len();
58    if len == 0 {
59        return None;
60    }
61    Some(arr.iter().position(|&x| x > cut).unwrap_or(len - 1)) // unwrap_or allows for default if not found
62}
63
64/// Return first index of `arr` equal to `cut`
65pub fn first_eq(arr: &[f64], cut: f64) -> Option<usize> {
66    let len = arr.len();
67    if len == 0 {
68        return None;
69    }
70    Some(arr.iter().position(|&x| x == cut).unwrap_or(len - 1)) // unwrap_or allows for default if not found
71}
72
73/// return max of 2 f64
74pub fn max(a: f64, b: f64) -> f64 {
75    a.max(b)
76}
77
78/// return min of 2 f64
79pub fn min(a: f64, b: f64) -> f64 {
80    a.min(b)
81}
82
83/// return max of arr
84pub fn arrmax(arr: &[f64]) -> f64 {
85    arr.iter().copied().fold(f64::NAN, f64::max)
86}
87
88/// return min of arr
89pub fn arrmin(arr: &[f64]) -> f64 {
90    arr.iter().copied().fold(f64::NAN, f64::min)
91}
92
93/// return true if the array is all zeros
94pub fn ndarrallzeros(arr: &Array1<f64>) -> bool {
95    arr.iter().all(|x| *x == 0.0)
96}
97
98/// return cumsum of arr
99pub fn ndarrcumsum(arr: &Array1<f64>) -> Array1<f64> {
100    arr.iter()
101        .scan(0.0, |acc, &x| {
102            *acc += x;
103            Some(*acc)
104        })
105        .collect()
106}
107
108/// return the unique values of the array
109pub fn ndarrunique(arr: &Array1<f64>) -> Array1<f64> {
110    let mut set = HashSet::new();
111    let mut new_arr = Vec::new();
112    let x_min = arr.min().unwrap();
113    let x_max = arr.max().unwrap();
114    let dx = if x_max == x_min { 1.0 } else { x_max - x_min };
115    for &x in arr.iter() {
116        let y = (((x - x_min) / dx) * (usize::MAX as f64)) as usize;
117        if !set.contains(&y) {
118            new_arr.push(x);
119            set.insert(y);
120        }
121    }
122    Array1::from_vec(new_arr)
123}
124
125/// interpolation algorithm from <http://www.cplusplus.com/forum/general/216928/>
126/// Arguments:
127/// x : value at which to interpolate
128pub fn interpolate(
129    x: &f64,
130    x_data_in: &Array1<f64>,
131    y_data_in: &Array1<f64>,
132    extrapolate: bool,
133) -> anyhow::Result<f64> {
134    ensure!(x_data_in.len() == y_data_in.len());
135    let mut new_x_data = Vec::new();
136    let mut new_y_data = Vec::new();
137    let mut last_x = x_data_in[0];
138    for idx in 0..x_data_in.len() {
139        if idx == 0 || (idx > 0 && x_data_in[idx] > last_x) {
140            last_x = x_data_in[idx];
141            new_x_data.push(x_data_in[idx]);
142            new_y_data.push(y_data_in[idx]);
143        }
144    }
145    let x_data = Array1::from_vec(new_x_data);
146    let y_data = Array1::from_vec(new_y_data);
147    let size = x_data.len();
148
149    let mut i = 0;
150    if x >= &x_data[size - 2] {
151        i = size - 2;
152    } else {
153        while x > &x_data[i + 1] {
154            i += 1;
155        }
156    }
157    let xl = &x_data[i];
158    let mut yl = &y_data[i];
159    let xr = &x_data[i + 1];
160    let mut yr = &y_data[i + 1];
161    if !extrapolate {
162        if x < xl {
163            yr = yl;
164        }
165        if x > xr {
166            yl = yr;
167        }
168    }
169    let dydx = (yr - yl) / (xr - xl);
170    Ok(yl + dydx * (x - xl))
171}
172
173/// interpolation algorithm from <http://www.cplusplus.com/forum/general/216928/>
174/// Arguments:
175/// x : value at which to interpolate
176pub fn interpolate_vectors(
177    x: &f64,
178    x_data_in: &Vec<f64>,
179    y_data_in: &Vec<f64>,
180    extrapolate: bool,
181) -> f64 {
182    assert!(x_data_in.len() == y_data_in.len());
183    let mut new_x_data = Vec::new();
184    let mut new_y_data = Vec::new();
185    let mut last_x = x_data_in[0];
186    for idx in 0..x_data_in.len() {
187        if idx == 0 || (idx > 0 && x_data_in[idx] > last_x) {
188            last_x = x_data_in[idx];
189            new_x_data.push(x_data_in[idx]);
190            new_y_data.push(y_data_in[idx]);
191        }
192    }
193    let x_data = new_x_data;
194    let y_data = new_y_data;
195    let size = x_data.len();
196
197    let mut i = 0;
198    if x >= &x_data[size - 2] {
199        i = size - 2;
200    } else {
201        while x > &x_data[i + 1] {
202            i += 1;
203        }
204    }
205    let xl = &x_data[i];
206    let mut yl = &y_data[i];
207    let xr = &x_data[i + 1];
208    let mut yr = &y_data[i + 1];
209    if !extrapolate {
210        if x < xl {
211            yr = yl;
212        }
213        if x > xr {
214            yl = yr;
215        }
216    }
217    let dydx = (yr - yl) / (xr - xl);
218    yl + dydx * (x - xl)
219}
220
221/// Generate all permutations of indices for a given *N*-dimensional array shape
222///
223/// # Arguments
224/// * `shape` - Reference to shape of the *N*-dimensional array, as returned by `ndarray::ArrayBase::shape()`
225///
226/// # Returns
227/// A `Vec<Vec<usize>>` where each inner `Vec<usize>` is one permutation of indices
228///
229/// # Example
230/// ```rust
231/// use fastsim_core::utils::get_index_permutations;
232/// let shape = [3, 2, 2];
233/// assert_eq!(
234///     get_index_permutations(&shape),
235///     [
236///         [0, 0, 0],
237///         [0, 0, 1],
238///         [0, 1, 0],
239///         [0, 1, 1],
240///         [1, 0, 0],
241///         [1, 0, 1],
242///         [1, 1, 0],
243///         [1, 1, 1],
244///         [2, 0, 0],
245///         [2, 0, 1],
246///         [2, 1, 0],
247///         [2, 1, 1],
248///     ]
249/// );
250/// ```
251///
252pub fn get_index_permutations(shape: &[usize]) -> Vec<Vec<usize>> {
253    if shape.is_empty() {
254        return vec![vec![]];
255    }
256    shape
257        .iter()
258        .map(|&len| 0..len)
259        .multi_cartesian_product()
260        .collect()
261}
262
263/// Multilinear interpolation function, accepting any dimensionality *N*.
264///
265/// # Arguments
266/// * `point` - An *N*-length array representing the interpolation point coordinates in each dimension
267/// * `grid` - A grid containing the coordinates for each dimension,
268///   i.e. `[[0.0, 1.0], [-0.5, 1.5]]` indicates x<sub>0</sub> = 0.0, x<sub>1</sub> = 1.0, y<sub>0</sub> = -0.5, y<sub>1</sub> = 1.5
269/// * `values` - An *N*-dimensional [`ndarray::ArrayD`] containing the values at given grid coordinates
270///
271/// # Errors
272/// This function returns an [`InterpolationError`] if any of the validation checks from [`validate_inputs`] fail,
273/// or if any values surrounding supplied `point` are `NaN`.
274///
275/// # Examples
276/// ## 1D Example
277/// ```rust
278/// use ndarray::prelude::*;
279/// use fastsim_core::utils::multilinear;
280///
281/// let grid = [vec![0.0, 1.0, 4.0]];
282/// let values = array![0.0, 2.0, 4.45].into_dyn();
283///
284/// let point_a = [0.82];
285/// assert_eq!(multilinear(&point_a, &grid, &values).unwrap(), 1.64);
286/// let point_b = [2.98];
287/// assert_eq!(multilinear(&point_b, &grid, &values).unwrap(), 3.617);
288/// let point_c = [grid[0][2]]; // returns value at x2
289/// assert_eq!(multilinear(&point_c, &grid, &values).unwrap(), values[2]);
290/// ```
291///
292/// ## 2D Example
293/// ```rust
294/// use ndarray::prelude::*;
295/// use fastsim_core::utils::multilinear;
296///
297/// let grid = [
298///     vec![0.0, 1.0, 2.0], // x0, x1, x2
299///     vec![0.0, 1.0, 2.0], // y0, y1, y2
300/// ];
301/// let values = array![
302///     [0.0, 2.0, 1.9], // (x0, y0), (x0, y1), (x0, y2)
303///     [2.0, 4.0, 3.1], // (x1, y0), (x1, y1), (x1, y2)
304///     [5.0, 0.0, 1.4], // (x2, y0), (x2, y1), (x2, y2)
305/// ]
306/// .into_dyn();
307///
308/// let point_a = [0.5, 0.5];
309/// assert_eq!(multilinear(&point_a, &grid, &values).unwrap(), 2.0);
310/// let point_b = [1.52, 0.36];
311/// assert_eq!(multilinear(&point_b, &grid, &values).unwrap(), 2.9696);
312/// let point_c = [grid[0][2], grid[1][1]]; // returns value at (x2, y1)
313/// assert_eq!(
314///     multilinear(&point_c, &grid, &values).unwrap(),
315///     values[[2, 1]]
316/// );
317/// ```
318///
319/// ## 2D Example with non-uniform dimension sizes
320/// ```rust
321/// use ndarray::prelude::*;
322/// use fastsim_core::utils::multilinear;
323///
324/// let grid = [
325///     vec![0.0, 1.0, 2.0], // x0, x1, x2
326///     vec![0.0, 1.0], // y0, y1
327/// ];
328/// let values = array![
329///     [0.0, 2.0], // f(x0, y0), f(x0, y1)
330///     [2.0, 4.0], // f(x1, y0), f(x1, y1)
331///     [5.0, 0.0], // f(x2, y0), f(x2, y1)
332/// ]
333/// .into_dyn();
334///
335/// let point_a = [0.5, 0.5];
336/// assert_eq!(multilinear(&point_a, &grid, &values).unwrap(), 2.0);
337/// let point_b = [1.52, 0.36];
338/// assert_eq!(multilinear(&point_b, &grid, &values).unwrap(), 2.9696);
339/// let point_c = [grid[0][2], grid[1][1]]; // returns value at (x2, y1)
340/// assert_eq!(
341///     multilinear(&point_c, &grid, &values).unwrap(),
342///     values[[2, 1]]
343/// );
344/// ```
345///
346/// ## 3D Example
347/// ```rust
348/// use ndarray::prelude::*;
349/// use fastsim_core::utils::multilinear;
350///
351/// let grid = [
352///     vec![0.0, 1.0, 2.0], // x0, x1, x2
353///     vec![0.0, 1.0, 2.0], // y0, y1, y2
354///     vec![0.0, 1.0, 2.0], // z0, z1, z2
355/// ];
356/// let values = array![
357///     [
358///         [0.0, 1.5, 3.0], // (x0, y0, z0), (x0, y0, z1), (x0, y0, z2)
359///         [2.0, 0.5, 1.4], // (x0, y1, z0), (x0, y1, z1), (x0, y1, z2)
360///         [1.9, 5.3, 2.2], // (x0, y2, z0), (x0, y0, z1), (x0, y2, z2)
361///     ],
362///     [
363///         [2.0, 5.1, 1.1], // (x1, y0, z0), (x1, y0, z1), (x1, y0, z2)
364///         [4.0, 1.0, 0.5], // (x1, y1, z0), (x1, y1, z1), (x1, y1, z2)
365///         [3.1, 0.9, 1.2], // (x1, y2, z0), (x1, y2, z1), (x1, y2, z2)
366///     ],
367///     [
368///         [5.0, 0.2, 5.1], // (x2, y0, z0), (x2, y0, z1), (x2, y0, z2)
369///         [0.7, 0.1, 3.2], // (x2, y1, z0), (x2, y1, z1), (x2, y1, z2)
370///         [1.4, 1.1, 0.0], // (x2, y2, z0), (x2, y2, z1), (x2, y2, z2)
371///     ],
372/// ]
373/// .into_dyn();
374///
375/// let point_a = [0.5, 0.5, 0.5];
376/// assert_eq!(multilinear(&point_a, &grid, &values).unwrap(), 2.0125);
377/// let point_b = [1.52, 0.36, 0.5];
378/// assert_eq!(multilinear(&point_b, &grid, &values).unwrap(), 2.46272);
379/// let point_c = [grid[0][2], grid[1][1], grid[2][0]]; // returns value at (x2, y1, z0)
380/// assert_eq!(
381///     multilinear(&point_c, &grid, &values).unwrap(),
382///     values[[2, 1, 0]]
383/// );
384/// ```
385///
386pub fn multilinear(point: &[f64], grid: &[Vec<f64>], values: &ArrayD<f64>) -> anyhow::Result<f64> {
387    // Dimensionality
388    let mut n = values.ndim();
389
390    // Validate inputs
391    anyhow::ensure!(
392        point.len() == n,
393        "Length of supplied `point` must be same as `values` dimensionality: {point:?} is not {n}-dimensional",
394    );
395    anyhow::ensure!(
396        grid.len() == n,
397        "Length of supplied `grid` must be same as `values` dimensionality: {grid:?} is not {n}-dimensional",
398    );
399    for i in 0..n {
400        anyhow::ensure!(
401            grid[i].len() == values.shape()[i],
402            "Supplied `grid` and `values` are not compatible shapes: dimension {i}, lengths {} != {}",
403            grid[i].len(),
404            values.shape()[i]
405        );
406        anyhow::ensure!(
407            grid[i].windows(2).all(|w| w[0] < w[1]),
408            "Supplied `grid` coordinates must be sorted and non-repeating: dimension {i}, {:?}",
409            grid[i]
410        );
411        anyhow::ensure!(
412            grid[i][0] <= point[i] && point[i] <= *grid[i].last().unwrap(),
413            "Supplied `point` must be within `grid` for dimension {i}: point[{i}] = {:?}, grid[{i}] = {:?}",
414            point[i],
415            grid[i],
416        );
417    }
418
419    // Point can share up to N values of a grid point, which reduces the problem dimensionality
420    // i.e. the point shares one of three values of a 3-D grid point, then the interpolation becomes 2-D at that slice
421    // or   if the point shares two of three values of a 3-D grid point, then the interpolation becomes 1-D
422    let mut point = point.to_vec();
423    let mut grid = grid.to_vec();
424    let mut values_view = values.view();
425    for dim in (0..n).rev() {
426        // Range is reversed so that removal doesn't affect indexing
427        if let Some(pos) = grid[dim]
428            .iter()
429            .position(|&grid_point| grid_point == point[dim])
430        {
431            point.remove(dim);
432            grid.remove(dim);
433            values_view.index_axis_inplace(Axis(dim), pos);
434        }
435    }
436    if values_view.len() == 1 {
437        // Supplied point is coincident with a grid point, so just return the value
438        return Ok(*values_view.first().unwrap());
439    }
440    // Simplified dimensionality
441    n = values_view.ndim();
442
443    // Extract the lower and upper indices for each dimension,
444    // as well as the fraction of how far the supplied point is between the surrounding grid points
445    let mut lower_idxs = Vec::with_capacity(n);
446    let mut interp_diffs = Vec::with_capacity(n);
447    for dim in 0..n {
448        let lower_idx = grid[dim]
449            .windows(2)
450            .position(|w| w[0] < point[dim] && point[dim] < w[1])
451            .unwrap();
452        let interp_diff =
453            (point[dim] - grid[dim][lower_idx]) / (grid[dim][lower_idx + 1] - grid[dim][lower_idx]);
454        lower_idxs.push(lower_idx);
455        interp_diffs.push(interp_diff);
456    }
457    // `interp_vals` contains all values surrounding the point of interest, starting with shape (2, 2, ...) in N dimensions
458    // this gets mutated and reduces in dimension each iteration, filling with the next values to interpolate with
459    // this ends up as a 0-dimensional array containing only the final interpolated value
460    let mut interp_vals = values_view
461        .slice_each_axis(|ax| {
462            let lower = lower_idxs[ax.axis.0];
463            Slice::from(lower..=lower + 1)
464        })
465        .to_owned();
466    let mut index_permutations = get_index_permutations(interp_vals.shape());
467    // This loop interpolates in each dimension sequentially
468    // each outer loop iteration the dimensionality reduces by 1
469    // `interp_vals` ends up as a 0-dimensional array containing only the final interpolated value
470    for (dim, diff) in interp_diffs.iter().enumerate() {
471        let next_dim = n - 1 - dim;
472        let next_shape = vec![2; next_dim];
473        // Indeces used for saving results of this dimensions interpolation results
474        // assigned to `index_permutations` at end of loop to be used for indexing in next iteration
475        let next_idxs = get_index_permutations(&next_shape);
476        let mut intermediate_arr = Array::default(next_shape);
477        for i in 0..next_idxs.len() {
478            // `next_idxs` is always half the length of `index_permutations`
479            let l = index_permutations[i].as_slice();
480            let u = index_permutations[next_idxs.len() + i].as_slice();
481            if dim == 0 {
482                anyhow::ensure!(
483                    !interp_vals[l].is_nan() && !interp_vals[u].is_nan(),
484                    "Surrounding value(s) cannot be NaN:\npoint = {point:?},\ngrid = {grid:?},\nvalues = {values:?}"
485                );
486            }
487            // This calculation happens 2^(n-1) times in the first iteration of the outer loop,
488            // 2^(n-2) times in the second iteration, etc.
489            intermediate_arr[next_idxs[i].as_slice()] =
490                interp_vals[l] * (1.0 - diff) + interp_vals[u] * diff;
491        }
492        index_permutations = next_idxs;
493        interp_vals = intermediate_arr;
494    }
495
496    // return the only value contained within the 0-dimensional array
497    Ok(*interp_vals.first().unwrap())
498}
499
500lazy_static! {
501    static ref TIRE_CODE_REGEX: Regex = Regex::new(
502        r"(?i)[P|LT|ST|T]?((?:[0-9]{2,3}\.)?[0-9]+)/((?:[0-9]{1,2}\.)?[0-9]+) ?[B|D|R]?[x|\-| ]?((?:[0-9]{1,2}\.)?[0-9]+)[A|B|C|D|E|F|G|H|J|L|M|N]?"
503    ).unwrap();
504}
505
506/// Calculate tire radius (in meters) from an [ISO metric tire code](https://en.wikipedia.org/wiki/Tire_code#ISO_metric_tire_codes)
507///
508/// # Arguments
509/// * `tire_code` - A string containing a parsable ISO metric tire code
510///
511/// # Examples
512/// ## Example 1:
513///
514/// ```rust
515/// // Note the floating point imprecision in the result
516/// use fastsim_core::utils::tire_code_to_radius;
517/// let tire_code = "225/70Rx19.5G";
518/// assert_eq!(tire_code_to_radius(&tire_code).unwrap(), 0.40514999999999995);
519/// ```
520///
521/// ## Example 2:
522///
523/// ```rust
524/// // Either `&str`, `&String`, or `String` can be passed
525/// use fastsim_core::utils::tire_code_to_radius;
526/// let tire_code = String::from("P205/60R16");
527/// assert_eq!(tire_code_to_radius(tire_code).unwrap(), 0.3262);
528/// ```
529///
530pub fn tire_code_to_radius<S: AsRef<str>>(tire_code: S) -> anyhow::Result<f64> {
531    let tire_code = tire_code.as_ref();
532    let captures = TIRE_CODE_REGEX.captures(tire_code).with_context(|| {
533        format!(
534            "Regex pattern does not match for {:?}: {:?}",
535            tire_code,
536            TIRE_CODE_REGEX.as_str(),
537        )
538    })?;
539    let width_mm: f64 = captures[1].parse()?;
540    let aspect_ratio: f64 = captures[2].parse()?;
541    let rim_diameter_in: f64 = captures[3].parse()?;
542
543    let sidewall_height_mm = width_mm * aspect_ratio / 100.0;
544    let radius_mm = (rim_diameter_in * 25.4) / 2.0 + sidewall_height_mm;
545
546    Ok(radius_mm / 1000.0)
547}
548
549/// Assumes the parent directory exists. Assumes file doesn't exist (i.e., newly created) or that it will be truncated if it does.
550#[cfg(feature = "default")]
551pub fn download_file_from_url(url: &str, file_path: &Path) -> anyhow::Result<()> {
552    let mut handle = Easy::new();
553    handle.follow_location(true)?;
554    handle.url(url)?;
555    let mut buffer = Vec::new();
556    {
557        let mut transfer = handle.transfer();
558        transfer.write_function(|data| {
559            buffer.extend_from_slice(data);
560            Ok(data.len())
561        })?;
562        let result = transfer.perform();
563        if result.is_err() {
564            return Err(anyhow!("Could not download from {}", url));
565        }
566    }
567    println!("Downloaded data from {}; bytes: {}", url, buffer.len());
568    if buffer.is_empty() {
569        return Err(anyhow!("No data available from {url}"));
570    }
571    {
572        let mut file = match File::create(file_path) {
573            Err(why) => {
574                return Err(anyhow!(
575                    "couldn't open {}: {}",
576                    file_path.to_str().unwrap(),
577                    why
578                ))
579            }
580            Ok(file) => file,
581        };
582        file.write_all(buffer.as_slice())?;
583    }
584    Ok(())
585}
586
587/// Creates/gets an OS-specific data directory and returns the path.
588#[cfg(feature = "default")]
589pub fn create_project_subdir<P: AsRef<Path>>(subpath: P) -> anyhow::Result<PathBuf> {
590    let proj_dirs = ProjectDirs::from("gov", "NREL", "fastsim").ok_or_else(|| {
591        anyhow!("Could not build path to project directory: \"gov.NREL.fastsim\"")
592    })?;
593    let path = PathBuf::from(proj_dirs.config_dir()).join(subpath);
594    std::fs::create_dir_all(path.as_path())?;
595    Ok(path)
596}
597
598/// Returns the path to the OS-specific data directory, if it exists.
599#[cfg(feature = "default")]
600pub fn path_to_cache() -> anyhow::Result<PathBuf> {
601    let proj_dirs = ProjectDirs::from("gov", "NREL", "fastsim").ok_or_else(|| {
602        anyhow!("Could not build path to project directory: \"gov.NREL.fastsim\"")
603    })?;
604    Ok(PathBuf::from(proj_dirs.config_dir()))
605}
606
607/// Deletes FASTSim data directory, clearing its contents. If subpath is
608/// provided, will only delete the subdirectory pointed to by the subpath,
609/// rather than deleting the whole data directory. If the subpath is an empty
610/// string, deletes the entire FASTSim directory.
611/// USE WITH CAUTION, as this function deletes ALL objects stored in the FASTSim
612/// data directory or provided subdirectory.
613/// # Arguments
614/// - subpath: Subpath to a subdirectory within the FASTSim data directory. If
615///   an empty string, the function will delete the whole FASTSim data
616///   directory, clearing all its contents.
617/// Note: it is not possible to delete single files using this function, only
618/// directories. If a single file needs deleting, the path_to_cache() function
619/// can be used to find the FASTSim data directory location. The file can then
620/// be found and manually deleted.
621#[cfg(feature = "default")]
622pub fn clear_cache<P: AsRef<Path>>(subpath: P) -> anyhow::Result<()> {
623    let path = path_to_cache()?.join(subpath);
624    Ok(std::fs::remove_dir_all(path)?)
625}
626
627/// takes an object from a url and saves it in the FASTSim data directory in a
628/// rust_objects folder
629/// WARNING: if there is a file already in the data subdirectory with the same
630/// name, it will be replaced by the new file
631/// to save to a folder other than rust_objects, define constant CACHE_FOLDER to
632/// be the desired folder name
633/// # Arguments
634/// - url: url (either as a string or url type) to object
635/// - subpath: path to subdirectory within FASTSim data directory. Suggested
636/// paths are "vehicles" for a RustVehicle, "cycles" for a RustCycle, and
637/// "rust_objects" for other Rust objects.
638/// Note: In order for the file to be save in the proper format, the URL needs
639/// to be a URL pointing directly to a file, for example a raw github URL.
640#[cfg(feature = "default")]
641pub fn url_to_cache<S: AsRef<str>, P: AsRef<Path>>(url: S, subpath: P) -> anyhow::Result<()> {
642    let url = url::Url::parse(url.as_ref())?;
643    let file_name = url
644        .path_segments()
645        .and_then(|segments| segments.last())
646        .with_context(|| "Could not parse filename from URL: {url:?}")?;
647    let data_subdirectory = create_project_subdir(subpath)
648        .with_context(|| "Could not find or build Fastsim data subdirectory.")?;
649    let file_path = data_subdirectory.join(file_name);
650    download_file_from_url(url.as_ref(), &file_path)?;
651    Ok(())
652}
653
654#[cfg(feature = "pyo3")]
655pub mod array_wrappers {
656    use crate::proc_macros::add_pyo3_api;
657
658    use super::*;
659    /// Helper struct to allow Rust to return a Python class that will indicate to the user that it's a clone.
660    #[add_pyo3_api]
661    #[derive(Default, Serialize, Deserialize, Clone, PartialEq, Eq)]
662    pub struct Pyo3ArrayU32(Array1<u32>);
663    impl SerdeAPI for Pyo3ArrayU32 {}
664
665    /// Helper struct to allow Rust to return a Python class that will indicate to the user that it's a clone.
666    #[add_pyo3_api]
667    #[derive(Default, Serialize, Deserialize, Clone, PartialEq, Eq)]
668    pub struct Pyo3ArrayI32(Array1<i32>);
669    impl SerdeAPI for Pyo3ArrayI32 {}
670
671    /// Helper struct to allow Rust to return a Python class that will indicate to the user that it's a clone.
672    #[add_pyo3_api]
673    #[derive(Default, Serialize, Deserialize, Clone, PartialEq)]
674    pub struct Pyo3ArrayF64(Array1<f64>);
675    impl SerdeAPI for Pyo3ArrayF64 {}
676
677    /// Helper struct to allow Rust to return a Python class that will indicate to the user that it's a clone.
678    #[add_pyo3_api]
679    #[derive(Default, Serialize, Deserialize, Clone, PartialEq, Eq)]
680    pub struct Pyo3ArrayBool(Array1<bool>);
681    impl SerdeAPI for Pyo3ArrayBool {}
682
683    /// Helper struct to allow Rust to return a Python class that will indicate to the user that it's a clone.
684    #[add_pyo3_api]
685    #[derive(Default, Serialize, Deserialize, Clone, PartialEq)]
686    pub struct Pyo3VecF64(Vec<f64>);
687
688    impl SerdeAPI for Pyo3VecF64 {}
689}
690
691#[cfg(feature = "pyo3")]
692pub use array_wrappers::*;
693
694#[cfg(test)]
695mod tests {
696    use super::*;
697    use crate::vehicle_utils::NETWORK_TEST_DISABLE_ENV_VAR_NAME;
698    use std::env;
699
700    #[test]
701    fn test_diff() {
702        assert_eq!(diff(&Array1::range(0.0, 3.0, 1.0)), array![0.0, 1.0, 1.0]);
703    }
704
705    #[test]
706    fn test_that_first_eq_finds_the_right_index_when_one_exists() {
707        let xs = [0.0, 1.2, 3.3, 4.4, 6.6];
708        let idx = first_eq(&xs, 3.3).unwrap();
709        let expected_idx = 2;
710        assert_eq!(idx, expected_idx)
711    }
712
713    #[test]
714    fn test_that_first_eq_yields_last_index_when_nothing_found() {
715        let xs = [0.0, 1.2, 3.3, 4.4, 6.6];
716        let idx = first_eq(&xs, 7.0).unwrap();
717        let expected_idx = xs.len() - 1;
718        assert_eq!(idx, expected_idx)
719    }
720
721    #[test]
722    fn test_that_first_grtr_finds_the_right_index_when_one_exists() {
723        let xs = [0.0, 1.2, 3.3, 4.4, 6.6];
724        let idx = first_grtr(&xs, 3.0).unwrap();
725        let expected_idx = 2;
726        assert_eq!(idx, expected_idx)
727    }
728
729    #[test]
730    fn test_that_first_grtr_yields_last_index_when_nothing_found() {
731        let xs = [0.0, 1.2, 3.3, 4.4, 6.6];
732        let idx = first_grtr(&xs, 7.0).unwrap();
733        let expected_idx = xs.len() - 1;
734        assert_eq!(idx, expected_idx)
735    }
736
737    #[test]
738    fn test_ndarrcumsum_expected_output() {
739        let xs = Array1::from_vec(vec![0.0, 1.0, 2.0, 3.0]);
740        let expected_ys = Array1::from_vec(vec![0.0, 1.0, 3.0, 6.0]);
741        let ys = ndarrcumsum(&xs);
742        for (i, (ye, y)) in expected_ys.iter().zip(ys.iter()).enumerate() {
743            assert_eq!(ye, y, "unequal at {}", i);
744        }
745    }
746
747    #[test]
748    fn test_add_from_yields_expected_output() {
749        let xs = Array1::from_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0]);
750        let mut expected_ys = Array1::from_vec(vec![1.0, 2.0, 3.0, 4.0, 5.0]);
751        let mut actual_ys = add_from(&xs, 100, 1.0);
752        assert_eq!(expected_ys.len(), actual_ys.len());
753        assert_eq!(expected_ys, actual_ys);
754        expected_ys = Array1::from_vec(vec![1.0, 2.0, 4.0, 5.0, 6.0]);
755        actual_ys = add_from(&xs, 2, 1.0);
756        assert_eq!(expected_ys.len(), actual_ys.len());
757        assert_eq!(expected_ys, actual_ys);
758    }
759
760    #[test]
761    fn test_ndarrunique_works() {
762        let xs = Array1::from_vec(vec![0.0, 1.0, 1.0, 2.0, 10.0, 10.0, 11.0]);
763        let expected = Array1::from_vec(vec![0.0, 1.0, 2.0, 10.0, 11.0]);
764        let actual = ndarrunique(&xs);
765        assert_eq!(expected.len(), actual.len());
766        for (ex, act) in expected.iter().zip(actual.iter()) {
767            assert_eq!(ex, act);
768        }
769    }
770    // #[test]
771    // fn test_that_argmax_does_the_right_thing_on_an_empty_array(){
772    //     let xs = Array::from_vec(vec![]);
773    //     let idx = first_grtr(&xs);
774    //     // unclear what should happen here; np.argmax throws a ValueError in the case of an empty vector
775    //     // ... possibly we should return an Option type?
776    //     let expected_idx:Option<usize> = None;
777    //     assert_eq!(idx, expected_idx);
778    // }
779
780    #[test]
781    fn test_that_interpolation_works() {
782        let xs = Array1::from_vec(vec![0.0, 1.0, 2.0, 3.0, 4.0]);
783        let ys = Array1::from_vec(vec![0.0, 10.0, 20.0, 30.0, 40.0]);
784        let x = 0.5;
785        let y_lookup = interpolate(&x, &xs, &ys, false).unwrap();
786        let expected_y_lookup = 5.0;
787        assert_eq!(expected_y_lookup, y_lookup);
788        let y_lookup = interpolate_vectors(&x, &xs.to_vec(), &ys.to_vec(), false);
789        assert_eq!(expected_y_lookup, y_lookup);
790    }
791
792    #[test]
793    fn test_that_interpolation_works_for_irrational_number() {
794        let xs = Array1::from_vec(vec![0.0, 1.0, 2.0, 3.0, 4.0]);
795        let ys = Array1::from_vec(vec![0.0, 10.0, 20.0, 30.0, 40.0]);
796        let x = 1.0 / 3.0;
797        let y_lookup = interpolate(&x, &xs, &ys, false).unwrap();
798        let expected_y_lookup = 3.3333333333;
799        assert!((expected_y_lookup - y_lookup).abs() < 1e-6);
800        let y_lookup = interpolate_vectors(&x, &xs.to_vec(), &ys.to_vec(), false);
801        assert!((expected_y_lookup - y_lookup).abs() < 1e-6);
802    }
803
804    #[test]
805    fn test_interpolate_with_small_vectors() {
806        let xs = Array1::from_vec(vec![0.0, 1.0]);
807        let ys = Array1::from_vec(vec![0.0, 10.0]);
808        let x = 0.5;
809        let y_lookup = interpolate(&x, &xs, &ys, false).unwrap();
810        let expected_y_lookup = 5.0;
811        assert!((expected_y_lookup - y_lookup).abs() < 1e-6);
812        let y_lookup = interpolate_vectors(&x, &xs.to_vec(), &ys.to_vec(), false);
813        assert!((expected_y_lookup - y_lookup).abs() < 1e-6);
814    }
815
816    #[test]
817    fn test_interpolate_when_lookup_is_at_end() {
818        let xs = Array1::from_vec(vec![0.0, 1.0]);
819        let ys = Array1::from_vec(vec![0.0, 10.0]);
820        let x = 1.0;
821        let y_lookup = interpolate(&x, &xs, &ys, false).unwrap();
822        let expected_y_lookup = 10.0;
823        assert!((expected_y_lookup - y_lookup).abs() < 1e-6);
824        let y_lookup = interpolate_vectors(&x, &xs.to_vec(), &ys.to_vec(), false);
825        assert!((expected_y_lookup - y_lookup).abs() < 1e-6);
826    }
827
828    #[test]
829    fn test_interpolate_when_lookup_is_past_end_without_extrapolate() {
830        let xs = Array1::from_vec(vec![0.0, 1.0]);
831        let ys = Array1::from_vec(vec![0.0, 10.0]);
832        let x = 1.01;
833        let y_lookup = interpolate(&x, &xs, &ys, false).unwrap();
834        let expected_y_lookup = 10.0;
835        assert!((expected_y_lookup - y_lookup).abs() < 1e-6);
836        let y_lookup = interpolate_vectors(&x, &xs.to_vec(), &ys.to_vec(), false);
837        assert!((expected_y_lookup - y_lookup).abs() < 1e-6);
838    }
839
840    #[test]
841    fn test_interpolate_with_x_data_that_repeats() {
842        let xs = Array1::from_vec(vec![0.0, 1.0, 1.0]);
843        let ys = Array1::from_vec(vec![0.0, 10.0, 10.0]);
844        let x = 1.0;
845        let y_lookup = interpolate(&x, &xs, &ys, false).unwrap();
846        let expected_y_lookup = 10.0;
847        assert_eq!(expected_y_lookup, y_lookup);
848        let y_lookup = interpolate_vectors(&x, &xs.to_vec(), &ys.to_vec(), false);
849        assert_eq!(expected_y_lookup, y_lookup);
850    }
851
852    #[test]
853    fn test_interpolate_with_non_evenly_spaced_x_data() {
854        let xs = Array1::from_vec(vec![0.0, 10.0, 100.0, 1000.0]);
855        let ys = Array1::from_vec(vec![0.0, 1.0, 2.0, 3.0]);
856        let x = 55.0;
857        let y_lookup = interpolate(&x, &xs, &ys, false).unwrap();
858        let expected_y_lookup = 1.5;
859        assert_eq!(expected_y_lookup, y_lookup);
860        let y_lookup = interpolate_vectors(&x, &xs.to_vec(), &ys.to_vec(), false);
861        assert_eq!(expected_y_lookup, y_lookup);
862    }
863
864    #[test]
865    fn test_path_to_cache() {
866        let path = path_to_cache().unwrap();
867        println!("{:?}", path);
868    }
869
870    #[test]
871    fn test_clear_cache() {
872        if env::var(NETWORK_TEST_DISABLE_ENV_VAR_NAME).is_ok() {
873            println!("SKIPPING: test_clear_cache");
874            return;
875        }
876        let temp_sub_dir = tempfile::TempDir::new_in(create_project_subdir("").unwrap()).unwrap();
877        let sub_dir_path = temp_sub_dir.path().to_str().unwrap();
878        let still_exists_before = std::fs::metadata(sub_dir_path).is_ok();
879        assert_eq!(still_exists_before, true);
880        url_to_cache("https://raw.githubusercontent.com/NREL/fastsim-vehicles/main/assets/2022_Tesla_Model_Y_RWD_example.yaml", "").unwrap();
881        clear_cache(sub_dir_path).unwrap();
882        let still_exists = std::fs::metadata(sub_dir_path).is_ok();
883        assert_eq!(still_exists, false);
884        let path_to_vehicle = path_to_cache()
885            .unwrap()
886            .join("2022_Tesla_Model_Y_RWD_example.yaml");
887        let vehicle_still_exists = std::fs::metadata(&path_to_vehicle).is_ok();
888        assert_eq!(vehicle_still_exists, true);
889        std::fs::remove_file(path_to_vehicle).unwrap();
890    }
891}
fastsim_core/utils.rs

fastsim_core/
utils.rs
