fastsim_core/
drive_cycle.rs

1pub mod maneuvers;
2pub mod manipulation_utils;
3
4use crate::drive_cycle::manipulation_utils::{
5    speed_for_constant_jerk, ConstantJerkTrajectory, CycleCache,
6};
7use crate::imports::*;
8use crate::prelude::*;
9use fastsim_2::cycle::RustCycle as Cycle2;
10use std::cmp;
11
12#[serde_api]
13#[derive(Clone, Debug, Deserialize, Serialize, PartialEq, Default)]
14#[non_exhaustive]
15#[serde(deny_unknown_fields)]
16#[cfg_attr(feature = "pyo3", pyclass(module = "fastsim", subclass, eq))]
17/// Container
18pub struct Cycle {
19    /// Name of cycle (can be left empty)
20    #[serde(default, skip_serializing_if = "String::is_empty")]
21    pub name: String,
22    /// inital elevation
23    pub init_elev: Option<si::Length>,
24    /// simulation time
25    pub time: Vec<si::Time>,
26    /// prescribed speed
27    #[serde(alias = "speed_mps")]
28    pub speed: Vec<si::Velocity>,
29    // TODO: consider trapezoidal integration scheme
30    /// calculated prescribed distance based on RHS integral of time and speed
31    #[serde(default, skip_serializing_if = "Vec::is_empty")]
32    pub dist: Vec<si::Length>,
33    /// road grade (expressed as a decimal, not percent)
34    #[serde(default, skip_serializing_if = "Vec::is_empty")]
35    pub grade: Vec<si::Ratio>,
36    // TODO: consider trapezoidal integration scheme
37    // TODO: @mokeefe, please check out how elevation is handled
38    /// calculated prescribed elevation based on RHS integral distance and grade
39    #[serde(default, skip_serializing_if = "Vec::is_empty")]
40    pub elev: Vec<si::Length>,
41    /// road charging/discharing capacity
42    #[serde(default, skip_serializing_if = "Vec::is_empty")]
43    pub pwr_max_chrg: Vec<si::Power>,
44    /// ambient air temperature w.r.t. to time (rather than spatial position)
45    #[serde(default, skip_serializing_if = "Vec::is_empty")]
46    pub temp_amb_air: Vec<si::Temperature>,
47    /// solar heat load w.r.t. to time (rather than spatial position)
48    #[serde(default, skip_serializing_if = "Vec::is_empty")]
49    pub pwr_solar_load: Vec<si::Power>,
50    // TODO: add provision for optional time-varying aux load
51    /// grade interpolator
52    #[serde(default, skip_serializing_if = "Option::is_none")]
53    pub grade_interp: Option<InterpolatorEnumOwned<f64>>,
54    /// elevation interpolator
55    #[serde(default, skip_serializing_if = "Option::is_none")]
56    pub elev_interp: Option<InterpolatorEnumOwned<f64>>,
57}
58
59#[pyo3_api]
60impl Cycle {
61    #[pyo3(name = "len")]
62    /// return the length of the cycle
63    fn len_py(&self) -> PyResult<usize> {
64        Ok(self.len_checked()?)
65    }
66
67    #[pyo3(name = "to_microtrips", signature=(stop_speed_m_per_s=None))]
68    /// convert cycle to a list of microtrips.
69    /// If stop speed is specified, it signifies the speed at or below which
70    /// a vehicle should be considered as stopped. This can be useful when
71    /// processing real-world data.
72    fn to_microtrips_py(&self, stop_speed_m_per_s: Option<f64>) -> PyResult<Vec<Cycle>> {
73        let stop_speed = stop_speed_m_per_s.map(|v| v * uc::MPS);
74        Ok(self.to_microtrips(stop_speed))
75    }
76
77    #[pyo3(name = "extend_time", signature=(absolute_time_s=None, time_fraction=None))]
78    /// extend cycle with idle time.
79    /// This is useful when a cycle's duration needs to be extended.
80    /// - absolute_time_s: optional time to extend the cycle
81    /// - time_fraction: optional fraction of cycle's duration to add to cycle.
82    ///
83    /// NOTE: if both absolute time and time fraction are specified, they
84    /// both add to extend the cycle. For example, if we have a 100 s cycle
85    /// and specify an absolute_time_s of 10 and time_fraction of 0.5, the
86    /// resulting cycle will have a duration of 160 s = 100.0 + (10 + 100.0 * 0.5)
87    fn extend_time_py(
88        &mut self,
89        absolute_time_s: Option<f64>,
90        time_fraction: Option<f64>,
91    ) -> PyResult<Cycle> {
92        let absolute_time = absolute_time_s.map(|t| t * uc::S);
93        let time_fraction = time_fraction.map(|f| f * uc::R);
94        Ok(self.extend_time(absolute_time, time_fraction))
95    }
96
97    #[pyo3(name = "dt_at_i")]
98    /// time step duration at step i.
99    pub fn dt_at_i_py(&self, i: usize) -> PyResult<f64> {
100        let i = std::cmp::max(1, i);
101        let dt = if i < self.time.len() {
102            self.time[i].get::<si::second>() - self.time[i - 1].get::<si::second>()
103        } else {
104            0.0
105        };
106        Ok(dt)
107    }
108
109    #[pyo3(name = "ending_idle_time_s")]
110    /// calculate and return the ending "idle" time of a cycle.
111    /// "Idle" time is defined as the amount of contiguous time
112    /// at the end of a cycle where the vehicle is not moving.
113    pub fn ending_idle_time_py(&self) -> PyResult<f64> {
114        let dt_end_idle = self.ending_idle_time();
115        Ok(dt_end_idle.get::<si::second>())
116    }
117
118    #[pyo3(name = "trim_ending_idle", signature=(idle_to_keep_s=None))]
119    /// trim ending "idle" time from a cycle.
120    /// The "idle" time is the time the vehicle is not moving.
121    /// - idle_to_keep_s: the amount of time to keep
122    ///
123    /// NOTE: if idle_to_keep_s is specified, the ending idle duration
124    /// will be UP TO this idle_to_keep_s amount but could be less if
125    /// there is insufficient idle time.
126    pub fn trim_ending_idle_py(&self, idle_to_keep_s: Option<f64>) -> PyResult<Cycle> {
127        let idle_to_keep = idle_to_keep_s.map(|idle| idle * uc::S);
128        Ok(self.trim_ending_idle(idle_to_keep))
129    }
130
131    #[pyo3(name = "average_speed_m_per_s", signature=(while_moving=None))]
132    /// calculate and return the average speed of the cycle in (m/s).
133    /// - while_moving: if specified and true, calculate the speed only
134    ///   while the vehicle is moving. Otherwise, calculate the average speed
135    ///   including stopped time.
136    pub fn average_speed_py(&self, while_moving: Option<bool>) -> PyResult<f64> {
137        let while_moving = while_moving.unwrap_or(false);
138        let vavg = self.average_speed(while_moving);
139        Ok(vavg.get::<si::meter_per_second>())
140    }
141
142    #[pyo3(name = "average_step_speeds_m_per_s")]
143    /// calculate and return the average speeds per time-step in (m/s).
144    pub fn average_step_speeds_py(&self) -> PyResult<Vec<f64>> {
145        Ok(self
146            .average_step_speeds()
147            .iter()
148            .map(|v| v.get::<si::meter_per_second>())
149            .collect())
150    }
151
152    #[pyo3(name = "average_step_speed_in_m_per_s_at")]
153    /// calculate the average step speed at the given step in (m/s).
154    pub fn average_step_speed_at_py(&self, i: usize) -> PyResult<f64> {
155        Ok(self.average_step_speed_at(i).get::<si::meter_per_second>())
156    }
157
158    #[pyo3(name = "resample")]
159    /// create a new cycle with the values resampled to the given time-step
160    /// duration.
161    pub fn resample_py(&self, time_step_s: f64) -> PyResult<Cycle> {
162        let time_step = time_step_s.max(0.01) * uc::S;
163        Ok(self.resample(time_step))
164    }
165}
166
167lazy_static! {
168    pub static ref ELEV_DEFAULT: si::Length = 400. * uc::FT;
169}
170
171impl Init for Cycle {
172    /// Sets `self.dist` and `self.elev`
173    /// # Assumptions
174    /// - if `init_elev.is_none()`, then defaults to [static@ELEV_DEFAULT]
175    fn init(&mut self) -> Result<(), Error> {
176        let _ = self
177            .len_checked()
178            .map_err(|err| Error::InitError(format_dbg!(err)))?;
179
180        if !self.temp_amb_air.is_empty() {
181            if self.temp_amb_air.len() != self.time.len() {
182                return Err(Error::InitError(format_dbg!()));
183            }
184        } else {
185            self.temp_amb_air = vec![*TE_STD_AIR; self.time.len()];
186        }
187
188        // calculate distance from RHS integral of speed and time
189        self.dist = {
190            self.time
191                .diff()
192                .iter()
193                .zip(&self.speed)
194                .scan(0. * uc::M, |dist, (dt, speed)| {
195                    *dist += *dt * *speed;
196                    Some(*dist)
197                })
198                .collect()
199        };
200
201        // populate grade if not provided
202        if self.grade.is_empty() {
203            self.grade = vec![
204                si::Ratio::ZERO;
205                self.len_checked()
206                    .map_err(|err| Error::InitError(format_dbg!(err)))?
207            ]
208        };
209        // calculate elevation from RHS integral of grade and distance
210        self.init_elev = self.init_elev.or_else(|| Some(*ELEV_DEFAULT));
211        self.elev = self
212            .grade
213            .iter()
214            .zip(&self.dist)
215            .scan(
216                // already guaranteed to be `Some`
217                self.init_elev.unwrap(),
218                |elev, (grade, dist)| {
219                    // TODO: Kyle, check this
220                    *elev += *dist * *grade;
221                    Some(*elev)
222                },
223            )
224            .collect();
225        let g0 = if !self.grade.is_empty() {
226            self.grade[0]
227        } else {
228            0.0 * uc::R
229        };
230        if self.grade.iter().all(|&g| g != g0) {
231            self.grade_interp = Some(
232                InterpolatorEnum::new_1d(
233                    self.dist.iter().map(|x| x.get::<si::meter>()).collect(),
234                    self.grade.iter().map(|y| y.get::<si::ratio>()).collect(),
235                    strategy::Linear,
236                    Extrapolate::Error,
237                )
238                .map_err(|e| Error::NinterpError(e.to_string()))?,
239            );
240
241            self.elev_interp = Some(
242                InterpolatorEnum::new_1d(
243                    self.dist.iter().map(|x| x.get::<si::meter>()).collect(),
244                    self.elev.iter().map(|y| y.get::<si::meter>()).collect(),
245                    strategy::Linear,
246                    Extrapolate::Error,
247                )
248                .map_err(|e| Error::NinterpError(e.to_string()))?,
249            );
250        } else {
251            self.grade_interp = Some(InterpolatorEnum::new_0d(g0.get::<si::ratio>()));
252            self.elev_interp = Some(InterpolatorEnum::new_0d(
253                self.init_elev.unwrap().get::<si::meter>(),
254            ));
255        }
256
257        Ok(())
258    }
259}
260
261impl SerdeAPI for Cycle {
262    const ACCEPTED_BYTE_FORMATS: &'static [&'static str] = &[
263        #[cfg(feature = "csv")]
264        "csv",
265        #[cfg(feature = "json")]
266        "json",
267        #[cfg(feature = "msgpack")]
268        "msgpack",
269        #[cfg(feature = "toml")]
270        "toml",
271        #[cfg(feature = "yaml")]
272        "yaml",
273    ];
274    const ACCEPTED_STR_FORMATS: &'static [&'static str] = &[
275        #[cfg(feature = "csv")]
276        "csv",
277        #[cfg(feature = "json")]
278        "json",
279        #[cfg(feature = "toml")]
280        "toml",
281        #[cfg(feature = "yaml")]
282        "yaml",
283    ];
284    #[cfg(feature = "resources")]
285    const RESOURCES_SUBDIR: &'static str = "cycles";
286
287    /// Write (serialize) an object into anything that implements [`std::io::Write`]
288    ///
289    /// # Arguments:
290    ///
291    /// * `wtr` - The writer into which to write object data
292    /// * `format` - The target format, any of those listed in [`ACCEPTED_BYTE_FORMATS`](`SerdeAPI::ACCEPTED_BYTE_FORMATS`)
293    ///
294    fn to_writer<W: std::io::Write>(&self, mut wtr: W, format: &str) -> Result<(), Error> {
295        match format.trim_start_matches('.').to_lowercase().as_str() {
296            #[cfg(feature = "csv")]
297            "csv" => {
298                let mut wtr = csv::Writer::from_writer(wtr);
299                for i in 0..self
300                    .len_checked()
301                    .map_err(|err| Error::SerdeError(format_dbg!(err)))?
302                {
303                    wtr.serialize(CycleElement {
304                        // unchecked indexing should be ok because of `self.len()`
305                        time: self.time[i],
306                        speed: self.speed[i],
307                        grade: if !self.grade.is_empty() {
308                            Some(self.grade[i])
309                        } else {
310                            None
311                        },
312                        pwr_max_charge: if !self.pwr_max_chrg.is_empty() {
313                            Some(self.pwr_max_chrg[i])
314                        } else {
315                            None
316                        },
317                        temp_amb_air: if !self.temp_amb_air.is_empty() {
318                            Some(self.temp_amb_air[i])
319                        } else {
320                            None
321                        },
322                        pwr_solar_load: if !self.pwr_solar_load.is_empty() {
323                            Some(self.pwr_solar_load[i])
324                        } else {
325                            None
326                        },
327                    })
328                    .map_err(|err| Error::SerdeError(format_dbg!(err)))?;
329                }
330                wtr.flush()
331                    .map_err(|err| Error::SerdeError(format_dbg!(err)))?
332            }
333            #[cfg(feature = "json")]
334            "json" => serde_json::to_writer(wtr, self)
335                .map_err(|err| Error::SerdeError(format_dbg!(err)))?,
336            #[cfg(feature = "toml")]
337            "toml" => {
338                let toml_string = self
339                    .to_toml()
340                    .map_err(|err| Error::SerdeError(format_dbg!(err)))?;
341                wtr.write_all(toml_string.as_bytes())
342                    .map_err(|err| Error::SerdeError(format_dbg!(err)))?;
343            }
344            #[cfg(feature = "yaml")]
345            "yaml" | "yml" => serde_yaml::to_writer(wtr, self)
346                .map_err(|err| Error::SerdeError(format_dbg!(err)))?,
347            _ => Err(Error::SerdeError(format!(
348                "Unsupported format {format:?}, must be one of {:?}",
349                Self::ACCEPTED_BYTE_FORMATS,
350            )))?,
351        }
352        Ok(())
353    }
354
355    /// Deserialize an object from anything that implements [`std::io::Read`]
356    ///
357    /// # Arguments:
358    ///
359    /// * `rdr` - The reader from which to read object data
360    /// * `format` - The source format, any of those listed in [`ACCEPTED_BYTE_FORMATS`](`SerdeAPI::ACCEPTED_BYTE_FORMATS`)
361    ///
362    fn from_reader<R: std::io::Read>(
363        rdr: &mut R,
364        format: &str,
365        skip_init: bool,
366    ) -> Result<Self, Error> {
367        let mut deserialized: Self =
368            match format.trim_start_matches('.').to_lowercase().as_str() {
369                #[cfg(feature = "csv")]
370                "csv" => {
371                    // Create empty cycle to be populated
372                    let mut cyc = Self::default();
373                    let mut rdr = csv::Reader::from_reader(rdr);
374                    for result in rdr.deserialize() {
375                        cyc.push(result.map_err(|err| Error::SerdeError(format_dbg!(err)))?)
376                            .map_err(|err| Error::SerdeError(format!("{err}")))?;
377                    }
378                    cyc
379                }
380                #[cfg(feature = "json")]
381                "json" => serde_json::from_reader(rdr)
382                    .map_err(|err| Error::SerdeError(format!("{err}")))?,
383                #[cfg(feature = "toml")]
384                "toml" => {
385                    let mut buf = String::new();
386                    rdr.read_to_string(&mut buf)
387                        .map_err(|err| Error::SerdeError(format_dbg!(err)))?;
388                    Self::from_toml(buf, skip_init)
389                        .map_err(|err| Error::SerdeError(format_dbg!(err)))?
390                }
391                #[cfg(feature = "yaml")]
392                "yaml" | "yml" => serde_yaml::from_reader(rdr)
393                    .map_err(|err| Error::SerdeError(format_dbg!(err)))?,
394                _ => {
395                    return Err(Error::SerdeError(format!(
396                        "Unsupported format {format:?}, must be one of {:?}",
397                        Self::ACCEPTED_BYTE_FORMATS
398                    )))
399                }
400            };
401        if !skip_init {
402            deserialized.init()?;
403        }
404        Ok(deserialized)
405    }
406
407    /// Write (serialize) an object into a string
408    ///
409    /// # Arguments:
410    ///
411    /// * `format` - The target format, any of those listed in [`ACCEPTED_STR_FORMATS`](`SerdeAPI::ACCEPTED_STR_FORMATS`)
412    ///
413    fn to_str(&self, format: &str) -> anyhow::Result<String> {
414        match format.trim_start_matches('.').to_lowercase().as_str() {
415            #[cfg(feature = "csv")]
416            "csv" => self.to_csv(),
417            #[cfg(feature = "json")]
418            "json" => self.to_json(),
419            #[cfg(feature = "toml")]
420            "toml" => self.to_toml(),
421            #[cfg(feature = "yaml")]
422            "yaml" | "yml" => self.to_yaml(),
423            _ => bail!(
424                "Unsupported format {format:?}, must be one of {:?}",
425                Self::ACCEPTED_STR_FORMATS
426            ),
427        }
428    }
429
430    /// Read (deserialize) an object from a string
431    ///
432    /// # Arguments:
433    ///
434    /// * `contents` - The string containing the object data
435    /// * `format` - The source format, any of those listed in [`ACCEPTED_STR_FORMATS`](`SerdeAPI::ACCEPTED_STR_FORMATS`)
436    ///
437    fn from_str<S: AsRef<str>>(contents: S, format: &str, skip_init: bool) -> anyhow::Result<Self> {
438        Ok(
439            match format.trim_start_matches('.').to_lowercase().as_str() {
440                #[cfg(feature = "csv")]
441                "csv" => Self::from_csv(contents, skip_init)?,
442                #[cfg(feature = "json")]
443                "json" => Self::from_json(contents, skip_init)?,
444                #[cfg(feature = "toml")]
445                "toml" => Self::from_toml(contents, skip_init)?,
446                #[cfg(feature = "yaml")]
447                "yaml" | "yml" => Self::from_yaml(contents, skip_init)?,
448                _ => bail!(
449                    "Unsupported format {format:?}, must be one of {:?}",
450                    Self::ACCEPTED_STR_FORMATS
451                ),
452            },
453        )
454    }
455}
456
457impl Cycle {
458    /// rust-internal time steps at i
459    pub fn dt_at_i(&self, i: usize) -> anyhow::Result<si::Time> {
460        Ok(*self.time.get(i).with_context(|| format_dbg!())?
461            - *self.time.get(i - 1).with_context(|| format_dbg!())?)
462    }
463
464    /// return the length of the cycle
465    pub fn len_checked(&self) -> anyhow::Result<usize> {
466        ensure!(
467            self.time.len() == self.speed.len(),
468            format!(
469                "{}\n`time` and `speed` fields do not have same `len()`",
470                format_dbg!()
471            )
472        );
473        ensure!(
474            self.dist.is_empty() || self.time.len() == self.dist.len(),
475            format!(
476                "{}\n`time` and `dist` fields do not have same `len()`",
477                format_dbg!()
478            )
479        );
480        ensure!(
481            self.grade.is_empty() || self.time.len() == self.grade.len(),
482            format!(
483                "{}\n`time` and `grade` fields do not have same `len()`",
484                format_dbg!()
485            )
486        );
487        ensure!(
488            self.elev.is_empty() || self.grade.len() == self.elev.len(),
489            format!(
490                "{}\n`grade` and `elev` fields do not have same `len()`",
491                format_dbg!()
492            )
493        );
494        ensure!(
495            self.pwr_max_chrg.is_empty() || self.time.len() == self.pwr_max_chrg.len(),
496            format!(
497                "{}\n`time` and `pwr_max_chrg` fields do not have same `len()`",
498                format_dbg!()
499            )
500        );
501        ensure!(
502            self.temp_amb_air.is_empty() || self.time.len() == self.temp_amb_air.len(),
503            format!(
504                "{}\n`time` and `temp_amb_air` fields do not have same `len()`",
505                format_dbg!()
506            )
507        );
508        Ok(self.time.len())
509    }
510
511    /// return true if the cycle is empty, else false
512    pub fn is_empty(&self) -> anyhow::Result<bool> {
513        Ok(self.len_checked().with_context(|| format_dbg!())? == 0)
514    }
515
516    /// append the given cycle element
517    pub fn push(&mut self, element: CycleElement) -> anyhow::Result<()> {
518        // TODO: maybe automate generation of this function as derive macro
519        // TODO: maybe automate `ensure!` that all vec fields are same length before returning result
520        // TODO: make sure all fields are being updated as appropriate
521        self.time.push(element.time);
522        self.speed.push(element.speed);
523        match element.grade {
524            Some(grade) => self.grade.push(grade),
525            None => self.grade.push(si::Ratio::ZERO),
526        }
527        match element.pwr_max_charge {
528            Some(pwr_max_chrg) => self.pwr_max_chrg.push(pwr_max_chrg),
529            None => self.pwr_max_chrg.push(si::Power::ZERO),
530        }
531        match element.temp_amb_air {
532            Some(temp_amb_air) => self.temp_amb_air.push(temp_amb_air),
533            None => self.temp_amb_air.push(*TE_STD_AIR),
534        }
535        match element.pwr_solar_load {
536            Some(pwr_solar_load) => self.pwr_solar_load.push(pwr_solar_load),
537            None => self.pwr_solar_load.push(si::Power::ZERO),
538        }
539        Ok(())
540    }
541
542    /// extend the cycle by a vector of elements
543    pub fn extend(&mut self, vec: Vec<CycleElement>) -> anyhow::Result<()> {
544        self.time.extend(vec.iter().map(|x| x.time).clone());
545        todo!();
546        // self.time.extend(vec.iter().map(|x| x.time).clone());
547        // match (&mut self.grade, vec.grade) {
548        //     (Some(grade_mut), Some(grade)) => grade_mut.push(grade),
549        //     (None, Some(_)) => {
550        //         bail!("Element and Cycle `grade` fields must both be `Some` or `None`")
551        //     }
552        //     (Some(_), None) => {
553        //         bail!("Element and Cycle `grade` fields must both be `Some` or `None`")
554        //     }
555        //     _ => {}
556        // }
557        // match (&mut self.pwr_max_chrg, vec.pwr_max_charge) {
558        //     (Some(pwr_max_chrg_mut), Some(pwr_max_chrg)) => pwr_max_chrg_mut.push(pwr_max_chrg),
559        //     (None, Some(_)) => {
560        //         bail!("Element and Cycle `pwr_max_chrg` fields must both be `Some` or `None`")
561        //     }
562        //     (Some(_), None) => {
563        //         bail!("Element and Cycle `pwr_max_chrg` fields must both be `Some` or `None`")
564        //     }
565        //     _ => {}
566        // }
567        // self.speed.push(vec.speed);
568        // Ok(())
569    }
570
571    /// trim the cycle to the given start_idx and end_idx.
572    ///
573    /// NOTE: ending cycle will include start_idx but NOT end_idx
574    pub fn trim(&mut self, start_idx: Option<usize>, end_idx: Option<usize>) -> anyhow::Result<()> {
575        let start_idx = start_idx.unwrap_or_default();
576        let len = self.len_checked().with_context(|| format_dbg!())?;
577        let end_idx = end_idx.unwrap_or(len);
578        ensure!(end_idx <= len, format_dbg!(end_idx <= len));
579
580        self.time = self.time[start_idx..end_idx].to_vec();
581        self.speed = self.speed[start_idx..end_idx].to_vec();
582        Ok(())
583    }
584
585    /// Write (serialize) cycle to a CSV string
586    #[cfg(feature = "csv")]
587    pub fn to_csv(&self) -> anyhow::Result<String> {
588        let mut buf = Vec::with_capacity(self.len_checked().with_context(|| format_dbg!())?);
589        self.to_writer(&mut buf, "csv")?;
590        Ok(String::from_utf8(buf)?)
591    }
592
593    /// Read (deserialize) an object from a CSV string
594    ///
595    /// # Arguments
596    ///
597    /// * `json_str` - JSON-formatted string to deserialize from
598    ///
599    #[cfg(feature = "csv")]
600    fn from_csv<S: AsRef<str>>(csv_str: S, skip_init: bool) -> anyhow::Result<Self> {
601        let mut csv_de = Self::from_reader(&mut csv_str.as_ref().as_bytes(), "csv", skip_init)?;
602        if !skip_init {
603            csv_de.init()?;
604        }
605        Ok(csv_de)
606    }
607
608    pub fn to_fastsim2(&self) -> anyhow::Result<Cycle2> {
609        let cyc2 = Cycle2 {
610            name: self.name.clone(),
611            time_s: self.time.iter().map(|t| t.get::<si::second>()).collect(),
612            mps: self
613                .speed
614                .iter()
615                .map(|s| s.get::<si::meter_per_second>())
616                .collect(),
617            grade: self.grade.iter().map(|g| g.get::<si::ratio>()).collect(),
618            orphaned: false,
619            road_type: vec![0.; self.len_checked().with_context(|| format_dbg!())?].into(),
620        };
621
622        Ok(cyc2)
623    }
624
625    /// convert cycle to a vector of CycleElement
626    pub fn to_elements(&self) -> Vec<CycleElement> {
627        let mut result = Vec::with_capacity(self.time.len());
628        for idx in 0..self.time.len() {
629            let element = CycleElement {
630                time: self.time[idx],
631                speed: self.speed[idx],
632                grade: if self.grade.is_empty() {
633                    None
634                } else {
635                    Some(self.grade[idx])
636                },
637                pwr_max_charge: if self.pwr_max_chrg.is_empty() {
638                    None
639                } else {
640                    Some(self.pwr_max_chrg[idx])
641                },
642                temp_amb_air: if self.temp_amb_air.is_empty() {
643                    None
644                } else {
645                    Some(self.temp_amb_air[idx])
646                },
647                pwr_solar_load: if self.pwr_solar_load.is_empty() {
648                    None
649                } else {
650                    Some(self.pwr_solar_load[idx])
651                },
652            };
653            result.push(element);
654        }
655        result
656    }
657
658    /// Convert cycle into a vector of "microtrips".
659    /// A microtrip is a start to a subsequent stop plus any idle time.
660    /// - stop_speed: the speed at or below which vehicle is considered "stopped"
661    ///
662    /// RETURN: vector of cycles with each cycle being a "microtrip".
663    pub fn to_microtrips(&self, stop_speed: Option<si::Velocity>) -> Vec<Cycle> {
664        let stop_speed = stop_speed.unwrap_or(1e-6 * uc::MPS);
665        let mut microtrips = Vec::new();
666        let mut current = Cycle {
667            name: self.name.clone(),
668            init_elev: self.init_elev,
669            time: vec![],
670            speed: vec![],
671            dist: vec![],
672            grade: vec![],
673            elev: vec![],
674            pwr_max_chrg: vec![],
675            temp_amb_air: vec![],
676            pwr_solar_load: vec![],
677            grade_interp: self.grade_interp.clone(),
678            elev_interp: self.elev_interp.clone(),
679        };
680        let elements = self.to_elements();
681        let mut moving: bool = false;
682        for element in &elements {
683            if element.speed > stop_speed && !moving && current.time.len() > 1 {
684                current.init().unwrap();
685                let last_idx = current.time.len() - 1;
686                let last_time = current.time[last_idx];
687                let last_speed = current.speed[last_idx];
688                let last_grade = if last_idx >= current.grade.len() {
689                    None
690                } else {
691                    Some(current.grade[last_idx])
692                };
693                let last_elevation = if last_idx >= current.elev.len() {
694                    None
695                } else {
696                    Some(current.elev[last_idx])
697                };
698                let last_temperature = if last_idx >= current.temp_amb_air.len() {
699                    None
700                } else {
701                    Some(current.temp_amb_air[last_idx])
702                };
703                let last_solar_load = if last_idx >= current.pwr_solar_load.len() {
704                    None
705                } else {
706                    Some(current.pwr_solar_load[last_idx])
707                };
708                let last_charge_power = if last_idx >= current.pwr_max_chrg.len() {
709                    None
710                } else {
711                    Some(current.pwr_max_chrg[last_idx])
712                };
713                current.time = current.time.iter().map(|t| *t - current.time[0]).collect();
714                microtrips.push(current.clone());
715                current = Cycle {
716                    name: self.name.clone(),
717                    init_elev: last_elevation,
718                    time: vec![last_time],
719                    speed: vec![last_speed],
720                    dist: vec![],
721                    grade: if let Some(g) = last_grade {
722                        vec![g]
723                    } else {
724                        vec![]
725                    },
726                    elev: vec![],
727                    pwr_max_chrg: if let Some(p) = last_charge_power {
728                        vec![p]
729                    } else {
730                        vec![]
731                    },
732                    temp_amb_air: if let Some(temp) = last_temperature {
733                        vec![temp]
734                    } else {
735                        vec![]
736                    },
737                    pwr_solar_load: if let Some(p) = last_solar_load {
738                        vec![p]
739                    } else {
740                        vec![]
741                    },
742                    grade_interp: self.grade_interp.clone(),
743                    elev_interp: self.elev_interp.clone(),
744                };
745            }
746            current
747                .push(element.clone())
748                .expect("Push shouldn't have an error path");
749            moving = element.speed > stop_speed;
750        }
751        if current.time.len() > 1 {
752            current.time = current.time.iter().map(|t| *t - current.time[0]).collect();
753            current.init().unwrap();
754            microtrips.push(current.clone());
755        }
756        microtrips
757    }
758
759    /// Determine average speed of cycle.
760    /// -- while_moving: if true, only takes average while moving
761    ///
762    /// RETURN: average speed
763    pub fn average_speed(&self, while_moving: bool) -> si::Velocity {
764        let mut d = si::Length::ZERO;
765        let mut t = si::Time::ZERO;
766        for idx in 1..self.speed.len() {
767            let dt = self.time[idx] - self.time[idx - 1];
768            let vavg = 0.5 * (self.speed[idx] + self.speed[idx - 1]);
769            let dd = vavg * dt;
770            let no_move = (dd.get::<si::meter>().ceil() as i32) == 0;
771            d += dd;
772            t += if while_moving && no_move {
773                si::Time::ZERO
774            } else {
775                dt
776            };
777        }
778        if t > si::Time::ZERO {
779            d / t
780        } else {
781            si::Velocity::ZERO
782        }
783    }
784
785    /// Return the average step speeds of the cycle as vector of velicities.
786    /// NOTE: the average speed from sample i-1 to i will appear as entry i.
787    /// RETURN: vector of velocities representing average step speeds.
788    pub fn average_step_speeds(&self) -> Vec<si::Velocity> {
789        let mut result = Vec::with_capacity(self.time.len());
790        result.push(0.0 * uc::MPS);
791        for i in 1..self.time.len() {
792            result.push(0.5 * (self.speed[i] + self.speed[i - 1]));
793        }
794        result
795    }
796
797    /// Calculate the average step speed at step i
798    /// (i.e., from sample point i-1 to i)
799    pub fn average_step_speed_at(&self, i: usize) -> si::Velocity {
800        if i >= self.speed.len() {
801            return 0.0 * uc::MPS;
802        }
803        0.5 * (self.speed[i] + self.speed[i - 1])
804    }
805
806    /// The distances traveled over each step using trapezoidal
807    /// integration.
808    pub fn trapz_step_distances(&self) -> Vec<si::Length> {
809        let mut result = Vec::with_capacity(self.time.len());
810        result.push(0.0 * uc::M);
811        for i in 1..self.time.len() {
812            let step_time = self.time[i] - self.time[i - 1];
813            let average_speed = 0.5 * (self.speed[i] + self.speed[i - 1]);
814            result.push(step_time * average_speed);
815        }
816        result
817    }
818
819    /// The elevation climb each step using trapezoidal integration.
820    pub fn trapz_step_elevations(&self) -> Vec<si::Length> {
821        let mut result = Vec::with_capacity(self.time.len());
822        result.push(0.0 * uc::M);
823        for i in 1..self.time.len() {
824            let step_time = self.time[i].get::<si::second>() - self.time[i - 1].get::<si::second>();
825            let average_speed = 0.5
826                * (self.speed[i].get::<si::meter_per_second>()
827                    + self.speed[i - 1].get::<si::meter_per_second>());
828            let step_dist = step_time * average_speed;
829            let gr = self.grade[i].get::<si::ratio>();
830            let dh = gr.atan().cos() * step_dist * gr;
831            result.push(dh * uc::M);
832        }
833        result
834    }
835
836    /// The distance traveled from start to the beginning of step i
837    /// (i.e., distance traveled up to sample point i-1)
838    pub fn trapz_step_start_distance(&self, step: usize) -> si::Length {
839        let mut distance = 0.0 * uc::M;
840        let step_max = cmp::min(step, self.time.len());
841        for i in 1..step_max {
842            let step_time = self.time[i] - self.time[i - 1];
843            let average_speed = 0.5 * (self.speed[i] + self.speed[i - 1]);
844            distance += step_time * average_speed;
845        }
846        distance
847    }
848
849    /// The distance traveled during the given step
850    /// (i.e., distance from sample point i-1 to i for step i)
851    pub fn trapz_distance_for_step(&self, step: usize) -> si::Length {
852        let average_speed = self.average_step_speed_at(step);
853        let elapsed_time = self.time[step] - self.time[step - 1];
854        average_speed * elapsed_time
855    }
856
857    /// Calculate the distance from step i_start to the start of step i_end
858    /// (i.e., distance from sample point i_start - 1 to i_end - 1)
859    pub fn trapz_distance_over_range(&self, step0: usize, step1: usize) -> si::Length {
860        let distances = self.trapz_step_distances();
861        let last_i = cmp::max(distances.len() - 1, 0);
862        let i_start = cmp::min(step0, last_i);
863        let i_end = cmp::min(step1, last_i);
864        let mut distance = 0.0 * uc::M;
865        for d in &distances[cmp::min(i_start, i_end)..cmp::max(i_start, i_end)] {
866            distance += *d;
867        }
868        distance
869    }
870
871    /// Calculate the time in a cycle spent moving
872    /// - stopped_speed_m_per_s: the speed above which we are considered to be moving
873    ///
874    /// RETURN: the time spent moving in seconds
875    pub fn time_spent_moving(&self, stopped_speed: Option<si::Velocity>) -> si::Time {
876        let stop_speed = stopped_speed.unwrap_or(0.0 * uc::MPS);
877        let mut result = 0.0 * uc::S;
878        for i in 1..self.time.len() {
879            let step_time = self.time[i] - self.time[i - 1];
880            if self.speed[i] > stop_speed || self.speed[i - 1] > stop_speed {
881                result += step_time;
882            }
883        }
884        result
885    }
886
887    /// Create distance and target speeds by microtrip.
888    /// Splits cycle into microtrips and returns a list of
889    /// 2-tuples of:
890    /// (distance from start in meters, target speed in m/s)
891    /// The distance is measured to the start of the microtrip.
892    ///
893    /// # Parameters
894    ///
895    /// * `blend_factor`: from 0.0 to 1.0
896    ///    - if 0.0, use the average speed of the microtrip
897    ///    - if 1.0, use the average speed while moving (i.e., no stopped time)
898    ///    - otherwise, something in between
899    /// * `min_target_speed`: the minimum target speed allowed
900    ///
901    /// # Result
902    ///
903    /// List of 2-tuple of (distance from start, target speed).
904    /// A tuple represents the distance from start of the start
905    /// of the given microtrip and its target speed.
906    ///
907    /// # Notes
908    ///
909    /// * target speed per microtrip is not allowed to be
910    ///   below the `min_target_speed`
911    pub fn distance_and_target_speeds_by_microtrip(
912        &self,
913        stop_speed: Option<si::Velocity>,
914        blend_factor: f64,
915        min_target_speed: si::Velocity,
916    ) -> Vec<(si::Length, si::Velocity)> {
917        let blend_factor = blend_factor.clamp(0.0, 1.0);
918        let mut result = Vec::new();
919        let microtrips = self.to_microtrips(stop_speed);
920        let mut distance_at_start = 0.0 * uc::M;
921        let t0 = 0.0 * uc::S;
922        let v0 = 0.0 * uc::MPS;
923        let d0 = 0.0 * uc::M;
924        for mt in microtrips {
925            let distance = mt
926                .trapz_step_distances()
927                .iter()
928                .fold(0.0 * uc::M, |total, dist| total + *dist);
929            let last_index = cmp::max(mt.time.len() - 1, 0);
930            let end_time = mt.time[last_index];
931            let start_time = mt.time[0];
932            let total_time = end_time - start_time;
933            let moving_time = mt.time_spent_moving(stop_speed);
934            let average_speed = if total_time > t0 {
935                distance / total_time
936            } else {
937                v0
938            };
939            let moving_average_speed = if moving_time > t0 {
940                distance / moving_time
941            } else {
942                v0
943            };
944            let target_speed =
945                blend_factor * (moving_average_speed - average_speed) + average_speed;
946            let target_speed = if target_speed > min_target_speed {
947                target_speed
948            } else {
949                min_target_speed
950            };
951            if distance > d0 {
952                result.push((distance_at_start, target_speed));
953                distance_at_start += distance;
954            }
955        }
956        result
957    }
958
959    /// Add idle time to Cycle.
960    /// By "idle" time, we mean "stopped" time (i.e., vehicle not moving).
961    pub fn extend_time(
962        &self,
963        absolute_time: Option<si::Time>,
964        time_fraction: Option<si::Ratio>,
965    ) -> Cycle {
966        let absolute_time = absolute_time.unwrap_or(0.0 * uc::S);
967        let time_fraction = time_fraction.unwrap_or(0.0 * uc::R);
968        let mut ts = self.time.clone();
969        let mut vs = self.speed.clone();
970        let mut gs = self.grade.clone();
971        let mut ps = self.pwr_max_chrg.clone();
972        let mut temps = self.temp_amb_air.clone();
973        let mut ss = self.pwr_solar_load.clone();
974        let t_end = *ts.last().unwrap();
975        let extra_time_s = (absolute_time.get::<si::second>()
976            + time_fraction.get::<si::ratio>() * t_end.get::<si::second>())
977        .round() as i32;
978        if extra_time_s == 0 {
979            return self.clone();
980        }
981        let dt = 1.0 * uc::S;
982        let dt_s = dt.get::<si::second>();
983        let mut idx = 1;
984        loop {
985            let dt_extra_s = dt_s * idx as f64;
986            if dt_extra_s > extra_time_s as f64 {
987                break;
988            }
989            ts.push(t_end + dt_extra_s * uc::S);
990            vs.push(0.0 * uc::MPS);
991            if !gs.is_empty() {
992                gs.push(0.0 * uc::R);
993            }
994            if !ps.is_empty() {
995                ps.push(*ps.last().unwrap());
996            }
997            if !temps.is_empty() {
998                temps.push(*temps.last().unwrap());
999            }
1000            if !ss.is_empty() {
1001                ss.push(*ss.last().unwrap());
1002            }
1003            idx += 1;
1004        }
1005        let mut cyc = Cycle {
1006            name: self.name.clone(),
1007            init_elev: self.init_elev,
1008            time: ts,
1009            speed: vs,
1010            dist: vec![],
1011            grade: gs,
1012            elev: vec![],
1013            pwr_max_chrg: vec![],
1014            grade_interp: self.grade_interp.clone(),
1015            elev_interp: self.elev_interp.clone(),
1016            temp_amb_air: temps,
1017            pwr_solar_load: ss,
1018        };
1019        cyc.init().unwrap();
1020        cyc
1021    }
1022
1023    /// Create a cache object for faster computations on Cycle.
1024    pub fn build_cache(&self) -> CycleCache {
1025        CycleCache::new(self)
1026    }
1027
1028    /// Returns the average grade over the given range of distances.
1029    /// - distance_start: the distance at start of evaluation area
1030    /// - delta_distance: distance traveled from distance_start
1031    /// - cache: optional CycleCache which can save computation time
1032    ///
1033    /// RETURN: average grade (rise over run) for the given range.
1034    ///
1035    /// NOTE: grade is assumed to be constant from just after the
1036    /// previous sample point until the current sample point (inclusive).
1037    /// That is, grade[i] applies from distance, d, of (d[i - 1], d[i]]
1038    pub fn average_grade_over_range(
1039        &self,
1040        distance_start: si::Length,
1041        delta_distance: si::Length,
1042        cache: Option<&CycleCache>,
1043    ) -> si::Ratio {
1044        let tol = 1e-6;
1045        match &cache {
1046            Some(cc) => {
1047                let dd_m = delta_distance.get::<si::meter>();
1048                if cc.grade_all_zero {
1049                    0.0 * uc::R
1050                } else if dd_m <= tol {
1051                    let dist_m = distance_start.get::<si::meter>();
1052                    cc.interp_grade(dist_m) * uc::R
1053                } else {
1054                    let dist0_m = distance_start.get::<si::meter>();
1055                    let dist1_m = dist0_m + dd_m;
1056                    let e0 = cc.interp_elevation(dist0_m);
1057                    let e1 = cc.interp_elevation(dist1_m);
1058                    ((e1 - e0) / dd_m).asin().tan() * uc::R
1059                }
1060            }
1061            None => {
1062                let zero_grade = 0.0 * uc::R;
1063                let grade_all_zero = {
1064                    let mut all0 = true;
1065                    for idx in 0..self.grade.len() {
1066                        if self.grade[idx] != zero_grade {
1067                            all0 = false;
1068                            break;
1069                        }
1070                    }
1071                    all0
1072                };
1073                if grade_all_zero {
1074                    0.0 * uc::R
1075                } else {
1076                    let delta_dists_m: Vec<f64> = self
1077                        .trapz_step_distances()
1078                        .iter()
1079                        .map(|dd| dd.get::<si::meter>())
1080                        .collect();
1081                    let trapz_distances_m = {
1082                        let mut d = 0.0;
1083                        let mut result = Vec::with_capacity(delta_dists_m.len());
1084                        for dd in &delta_dists_m {
1085                            d += *dd;
1086                            result.push(d);
1087                        }
1088                        result
1089                    };
1090                    let dist0_m = distance_start.get::<si::meter>();
1091                    let dd_m = delta_distance.get::<si::meter>();
1092                    let dist1_m = dist0_m + dd_m;
1093                    if dd_m < tol {
1094                        if dist0_m < trapz_distances_m[0] {
1095                            return self.grade[0];
1096                        }
1097                        let max_idx = self.grade.len() - 1;
1098                        if dist0_m > trapz_distances_m[max_idx] {
1099                            return self.grade[max_idx];
1100                        }
1101                        for idx in 1..self.time.len() {
1102                            if dist0_m > trapz_distances_m[idx - 1]
1103                                && dist0_m <= trapz_distances_m[idx]
1104                            {
1105                                return self.grade[idx];
1106                            }
1107                        }
1108                        self.grade[max_idx]
1109                    } else {
1110                        // NOTE: we use the following instead of delta_elev_m
1111                        // as it uses more precise trapezoidal diatance and
1112                        // elevation at sample points. This also uses the
1113                        // fully accurate trig functions in case we have large
1114                        // slope angles. This level of rigor may be overkill.
1115                        let trapz_elevations_m = {
1116                            let delta_elevs_m: Vec<f64> = self
1117                                .grade
1118                                .iter()
1119                                .zip(delta_dists_m)
1120                                .map(|(g, dd)| {
1121                                    let gr = g.get::<si::ratio>();
1122                                    gr.atan().cos() * dd * gr
1123                                })
1124                                .collect();
1125                            let mut result = Vec::with_capacity(delta_elevs_m.len());
1126                            let mut elev_m = 0.0;
1127                            for de in &delta_elevs_m {
1128                                elev_m += *de;
1129                                result.push(elev_m);
1130                            }
1131                            result
1132                        };
1133                        let interp: InterpolatorEnum<ndarray::OwnedRepr<f64>> =
1134                            InterpolatorEnum::new_1d(
1135                                trapz_distances_m.clone().into(),
1136                                trapz_elevations_m.clone().into(),
1137                                strategy::Linear,
1138                                Extrapolate::Clamp,
1139                            )
1140                            .unwrap();
1141                        let e0_m = interp.interpolate(&[dist0_m]).unwrap();
1142                        let e1_m = interp.interpolate(&[dist1_m]).unwrap();
1143                        ((e1_m - e0_m) / dd_m).asin().tan() * uc::R
1144                    }
1145                }
1146            }
1147        }
1148    }
1149
1150    /// Calculate the distance to next stop from `distance`.
1151    /// - distance: the distance to calculate distance-to-stop from
1152    ///
1153    /// RETURN: returns the distance to the next stop from `distance`
1154    ///
1155    /// NOTE: distance may be negative if we're beyond the last stop
1156    pub fn calc_distance_to_next_stop_from(
1157        &self,
1158        distance: si::Length,
1159        cache: Option<&CycleCache>,
1160    ) -> si::Length {
1161        let tol = 1e-6;
1162        let distance_m = distance.get::<si::meter>();
1163        match cache {
1164            Some(cc) => {
1165                for (&d_m, &v) in cc.trapz_distances_m.iter().zip(self.speed.iter()) {
1166                    let v_mps = v.get::<si::meter_per_second>();
1167                    if (v_mps < tol) && (d_m > (distance_m + tol)) {
1168                        return (d_m - distance_m) * uc::M;
1169                    }
1170                }
1171                (*cc.trapz_distances_m.last().unwrap_or(&0.0) * uc::M) - distance
1172            }
1173            None => {
1174                let ds_m = {
1175                    let mut result = Vec::with_capacity(self.time.len());
1176                    let mut d_m = 0.0;
1177                    for dd in self.trapz_step_distances() {
1178                        let dd_m = dd.get::<si::meter>();
1179                        d_m += dd_m;
1180                        result.push(d_m);
1181                    }
1182                    result
1183                };
1184                for (&d_m, &v) in ds_m.iter().zip(self.speed.iter()) {
1185                    let v_mps = v.get::<si::meter_per_second>();
1186                    if (v_mps < tol) && (d_m > (distance_m + tol)) {
1187                        return (d_m - distance_m) * uc::M;
1188                    }
1189                }
1190                *ds_m.last().unwrap_or(&0.0) * uc::M
1191            }
1192        }
1193    }
1194
1195    /// Modify the cycle using the given constant-jerk trajectory.
1196    /// - i: the index into the cycle to initiate modification
1197    ///   NOTE: THIS point is modified as trajectory is calculated as
1198    ///   starting at i-1
1199    /// - n: the number of steps ahead
1200    /// - jerk: the jerk (deriviative of acceleration with time)
1201    /// - accel0: the starting accelartion
1202    ///
1203    /// NOTE:
1204    /// - modifies the cycle in-place. Purpose is to allow hitting
1205    ///   a rendezvous point in time/speed in the future.
1206    /// - CAUTION: not robust against variable duration time-steps
1207    ///
1208    /// RETURN: the final modified speed
1209    pub fn modify_by_const_jerk_trajectory(
1210        &mut self,
1211        i: usize,
1212        n: usize,
1213        jerk: si::Jerk,
1214        accel0: si::Acceleration,
1215    ) -> si::Velocity {
1216        if n == 0 {
1217            return si::Velocity::ZERO;
1218        }
1219        let jerk_m_per_s3 = jerk.get::<si::meter_per_second_cubed>();
1220        let accel0_m_per_s2 = accel0.get::<si::meter_per_second_squared>();
1221        let num_samples = self.speed.len();
1222        if i >= num_samples {
1223            if num_samples > 0 {
1224                return self.speed[num_samples - 1];
1225            }
1226            return si::Velocity::ZERO;
1227        }
1228        let v0 = self.speed[i - 1].get::<si::meter_per_second>();
1229        let dt = (self.time[i] - self.time[i - 1]).get::<si::second>();
1230        let mut v = v0;
1231        for ni in 1..(n + 1) {
1232            let idx_to_set = (i - 1) + ni;
1233            if idx_to_set >= num_samples {
1234                break;
1235            }
1236            v = speed_for_constant_jerk(ni, v0, accel0_m_per_s2, jerk_m_per_s3, dt);
1237            self.speed[idx_to_set] = v.max(0.0) * uc::MPS;
1238        }
1239        self.init().unwrap();
1240        v * uc::MPS
1241    }
1242
1243    /// Modify cycle to add a braking trajectory that would cover the same
1244    /// distance as the given constant brake deceleration.
1245    /// - brake_accel: the brake acceleration (m/s2); must be negative
1246    /// - i: index where to initiate the stop trectory; start of the step
1247    /// - desired_distance_to_stop: the desired distance to stop within. If
1248    ///   not provided, it is calculated based on the braking deceleration.
1249    ///
1250    /// RETURN: (final speed of modified trajectory, number of steps to complete)
1251    /// - the final speed should be zero ideally
1252    /// - the number of time-steps required to complete the braking maneuver
1253    ///
1254    /// NOTE:
1255    /// - modifies the cycle in-place.
1256    pub fn modify_with_braking_trajectory(
1257        &mut self,
1258        brake_accel: si::Acceleration,
1259        i: usize,
1260        desired_distance_to_stop: Option<si::Length>,
1261    ) -> (si::Velocity, usize) {
1262        let brake_accel = if brake_accel > si::Acceleration::ZERO {
1263            -brake_accel
1264        } else {
1265            brake_accel
1266        };
1267        assert!(brake_accel < si::Acceleration::ZERO);
1268        if i >= self.time.len() {
1269            return (*self.speed.last().unwrap(), 0);
1270        }
1271        let i = if i < 1 { 1 } else { i };
1272        let v0 = self.speed[i - 1].get::<si::meter_per_second>();
1273        let dt = (self.time[i] - self.time[i - 1]).get::<si::second>();
1274        let brake_accel_m_per_s2 = brake_accel.get::<si::meter_per_second_squared>();
1275        // distance-to-stop (m)
1276        let dts_m = match desired_distance_to_stop {
1277            Some(value) => {
1278                let result = value.get::<si::meter>();
1279                if result > 0.0 {
1280                    result
1281                } else {
1282                    -0.5 * v0 * v0 / brake_accel_m_per_s2
1283                }
1284            }
1285            None => -0.5 * v0 * v0 / brake_accel_m_per_s2,
1286        };
1287        if dts_m <= 0.0 {
1288            return (v0 * uc::MPS, 0);
1289        }
1290        // time-to-stop (s)
1291        let tts_s = -v0 / brake_accel_m_per_s2;
1292        // number of steps to stop
1293        let n = (tts_s / dt).round() as usize;
1294        let n = if n < 2 { 2 } else { n }; // need at least 2 steps
1295        let traj =
1296            ConstantJerkTrajectory::from_speed_and_distance_targets(n, 0.0, v0, dts_m, 0.0, dt);
1297        let v_final = self.modify_by_const_jerk_trajectory(
1298            i,
1299            n,
1300            traj.jerk_m_per_s3 * uc::MPS3,
1301            traj.acceleration_m_per_s2 * uc::MPS2,
1302        );
1303        (v_final, n)
1304    }
1305
1306    /// Report the stopped time (i.e., idle) at the end of a cycle.
1307    ///
1308    /// RESULT: time vehicle is at zero speed at cycle end
1309    pub fn ending_idle_time(&self) -> si::Time {
1310        let mut result = si::Time::ZERO;
1311        let vzero = si::Velocity::ZERO;
1312        for idx in (1..self.time.len()).rev() {
1313            let v0 = self.speed[idx - 1];
1314            let v1 = self.speed[idx];
1315            if v0 != vzero || v1 != vzero {
1316                break;
1317            } else {
1318                let dt = self.time[idx] - self.time[idx - 1];
1319                result += dt;
1320            }
1321        }
1322        result
1323    }
1324
1325    /// Remove idel time at end of cycle except for the optionally
1326    /// specified duration.
1327    /// - idle_to_keep: optional duration of idle time to keep. Default is 0 s
1328    ///
1329    /// RESULT: a new cycle with idle time trimmed.
1330    pub fn trim_ending_idle(&self, idle_to_keep: Option<si::Time>) -> Cycle {
1331        let idle_to_keep = idle_to_keep.unwrap_or(si::Time::ZERO).max(si::Time::ZERO);
1332        let vzero = si::Velocity::ZERO;
1333        let mut idle_start_idx = 0;
1334        for idx in (1..self.time.len()).rev() {
1335            let v0 = self.speed[idx - 1];
1336            let v1 = self.speed[idx];
1337            if v0 != vzero || v1 != vzero {
1338                idle_start_idx = idx + 1;
1339                break;
1340            }
1341        }
1342        if idle_start_idx >= self.time.len() {
1343            return self.clone();
1344        }
1345        let end_idx = if idle_to_keep == si::Time::ZERO {
1346            idle_start_idx
1347        } else {
1348            let mut dt_idle = si::Time::ZERO;
1349            let mut idx_drop = idle_start_idx;
1350            for idx in idle_start_idx..self.time.len() {
1351                let dt = self.time[idx] - self.time[idx - 1];
1352                dt_idle += dt;
1353                if dt_idle > idle_to_keep {
1354                    idx_drop = idx;
1355                    break;
1356                }
1357            }
1358            idx_drop
1359        };
1360        let mut cyc = Cycle {
1361            name: self.name.clone(),
1362            time: self.time[0..end_idx].to_vec(),
1363            speed: self.speed[0..end_idx].to_vec(),
1364            init_elev: self.init_elev,
1365            grade: if self.grade.is_empty() {
1366                vec![]
1367            } else {
1368                self.grade[0..end_idx].to_vec()
1369            },
1370            dist: vec![],
1371            elev: vec![],
1372            pwr_max_chrg: if self.pwr_max_chrg.is_empty() {
1373                vec![]
1374            } else {
1375                self.pwr_max_chrg[0..end_idx].to_vec()
1376            },
1377            temp_amb_air: if self.temp_amb_air.is_empty() {
1378                vec![]
1379            } else {
1380                self.temp_amb_air[0..end_idx].to_vec()
1381            },
1382            pwr_solar_load: if self.pwr_solar_load.is_empty() {
1383                vec![]
1384            } else {
1385                self.pwr_solar_load[0..end_idx].to_vec()
1386            },
1387            grade_interp: None,
1388            elev_interp: None,
1389        };
1390        cyc.init().unwrap();
1391        cyc
1392    }
1393
1394    /// Resample cycle to a lower or higher frequency.
1395    /// - dt: the new step duration.
1396    ///
1397    /// RETURN: cycle
1398    /// NOTE: a value of dt <= 0 s implies to just clone the current cycle
1399    /// "as is"
1400    pub fn resample(&self, dt: si::Time) -> Cycle {
1401        if dt <= si::Time::ZERO {
1402            return self.clone();
1403        }
1404        let mut t = si::Time::ZERO;
1405        let speed_interp: InterpolatorEnum<OwnedRepr<f64>> = InterpolatorEnum::new_1d(
1406            self.time.iter().map(|x| x.get::<si::second>()).collect(),
1407            self.speed
1408                .iter()
1409                .map(|y| y.get::<si::meter_per_second>())
1410                .collect(),
1411            strategy::Linear,
1412            Extrapolate::Clamp,
1413        )
1414        .unwrap();
1415        let grade_interp: InterpolatorEnum<OwnedRepr<f64>> = InterpolatorEnum::new_1d(
1416            self.time.iter().map(|x| x.get::<si::second>()).collect(),
1417            self.grade.iter().map(|y| y.get::<si::ratio>()).collect(),
1418            strategy::RightNearest,
1419            Extrapolate::Clamp,
1420        )
1421        .unwrap();
1422        let temp_interp: Option<InterpolatorEnum<OwnedRepr<f64>>> =
1423            if self.temp_amb_air.len() == self.time.len() {
1424                Some(
1425                    InterpolatorEnum::new_1d(
1426                        self.time.iter().map(|t| t.get::<si::second>()).collect(),
1427                        self.temp_amb_air
1428                            .iter()
1429                            .map(|temp| temp.get::<si::kelvin_abs>())
1430                            .collect(),
1431                        strategy::Linear,
1432                        Extrapolate::Clamp,
1433                    )
1434                    .unwrap(),
1435                )
1436            } else {
1437                None
1438            };
1439        let solar_interp: Option<InterpolatorEnum<OwnedRepr<f64>>> =
1440            if self.pwr_solar_load.len() == self.time.len() {
1441                Some(
1442                    InterpolatorEnum::new_1d(
1443                        self.time.iter().map(|t| t.get::<si::second>()).collect(),
1444                        self.pwr_solar_load
1445                            .iter()
1446                            .map(|p| p.get::<si::kilowatt>())
1447                            .collect(),
1448                        strategy::Linear,
1449                        Extrapolate::Clamp,
1450                    )
1451                    .unwrap(),
1452                )
1453            } else {
1454                None
1455            };
1456        let chg_pwr_interp: Option<InterpolatorEnum<OwnedRepr<f64>>> =
1457            if self.pwr_max_chrg.len() == self.time.len() {
1458                Some(
1459                    InterpolatorEnum::new_1d(
1460                        self.time.iter().map(|t| t.get::<si::second>()).collect(),
1461                        self.pwr_max_chrg
1462                            .iter()
1463                            .map(|p| p.get::<si::kilowatt>())
1464                            .collect(),
1465                        strategy::Linear,
1466                        Extrapolate::Clamp,
1467                    )
1468                    .unwrap(),
1469                )
1470            } else {
1471                None
1472            };
1473        let mut ts = vec![];
1474        let mut vs = vec![];
1475        let mut gs = vec![];
1476        let mut pwr_chg = vec![];
1477        let mut temps = vec![];
1478        let mut solars = vec![];
1479        while t <= self.time[self.time.len() - 1] {
1480            ts.push(t);
1481            let t0 = t.get::<si::second>();
1482            let v = speed_interp.interpolate(&[t0]).unwrap();
1483            vs.push(v * uc::MPS);
1484            let g = grade_interp.interpolate(&[t0]).unwrap();
1485            gs.push(g * uc::R);
1486            if let Some(ref interp) = chg_pwr_interp {
1487                let pchg = interp.interpolate(&[t0]).unwrap();
1488                pwr_chg.push(pchg * uc::KW);
1489            }
1490            if let Some(ref interp) = temp_interp {
1491                let temp = interp.interpolate(&[t0]).unwrap();
1492                temps.push(temp * uc::KELVIN);
1493            }
1494            if let Some(ref interp) = solar_interp {
1495                let solar = interp.interpolate(&[t0]).unwrap();
1496                solars.push(solar * uc::KW);
1497            }
1498            t += dt;
1499        }
1500
1501        let mut cyc = Cycle {
1502            name: self.name.clone(),
1503            init_elev: self.init_elev,
1504            time: ts,
1505            speed: vs,
1506            dist: vec![],
1507            grade: gs,
1508            elev: vec![],
1509            pwr_max_chrg: pwr_chg,
1510            temp_amb_air: temps,
1511            pwr_solar_load: solars,
1512            grade_interp: None,
1513            elev_interp: None,
1514        };
1515        cyc.init().unwrap();
1516        cyc
1517    }
1518}
1519
1520impl TryFrom<CycleBuilder> for Cycle {
1521    type Error = anyhow::Error;
1522    fn try_from(value: CycleBuilder) -> anyhow::Result<Self, Self::Error> {
1523        let mut cyc = Self {
1524            name: value.name,
1525            init_elev: None,
1526            time: value.time,
1527            speed: value.speed,
1528            dist: Default::default(),
1529            grade: Default::default(),
1530            elev: Default::default(),
1531            pwr_max_chrg: Default::default(),
1532            temp_amb_air: Default::default(),
1533            pwr_solar_load: Default::default(),
1534            grade_interp: None,
1535            elev_interp: Default::default(),
1536        };
1537        cyc.init()?;
1538        Ok(cyc)
1539    }
1540}
1541
1542/// Trait for CycleBuilder and Cycle to support builder pattern
1543pub trait CBTrait {
1544    /// Return cycle with `grade`
1545    fn with_grade(&mut self, grade: Vec<si::Ratio>) -> anyhow::Result<Cycle>;
1546
1547    /// Return cycle with `temp_amb_air`
1548    fn with_temp_amb_air(&mut self, temp_amb_air: Vec<si::Temperature>) -> anyhow::Result<Cycle>;
1549
1550    // TODO: add more of these builder helpers
1551}
1552
1553impl CBTrait for Cycle {
1554    fn with_grade(&mut self, grade: Vec<si::Ratio>) -> anyhow::Result<Cycle> {
1555        ensure!(
1556            self.len_checked().with_context(|| format_dbg!())? == grade.len(),
1557            format!(
1558                "{}\n`self.len()`: `{}\n`grade.len()`",
1559                self.len_checked().with_context(|| format_dbg!())?,
1560                grade.len()
1561            )
1562        );
1563        self.grade = grade;
1564        Ok(self.clone())
1565    }
1566
1567    fn with_temp_amb_air(&mut self, temp_amb_air: Vec<si::Temperature>) -> anyhow::Result<Cycle> {
1568        ensure!(
1569            self.len_checked().with_context(|| format_dbg!())? == temp_amb_air.len(),
1570            format!(
1571                "{}\n`self.len()`: `{}\n`temp_amb_air.len()`",
1572                self.len_checked().with_context(|| format_dbg!())?,
1573                temp_amb_air.len()
1574            )
1575        );
1576        self.temp_amb_air = temp_amb_air;
1577        Ok(self.clone())
1578    }
1579}
1580
1581#[serde_api]
1582#[derive(Default, Debug, Serialize, Deserialize, PartialEq, Clone)]
1583#[non_exhaustive]
1584/// Simple cycle to be converted into [Cycle] with appropriate defaults
1585pub struct CycleBuilder {
1586    /// Name of cycle (can be left empty)
1587    #[serde(default, skip_serializing_if = "String::is_empty")]
1588    pub name: String,
1589    /// simulation time
1590    pub time: Vec<si::Time>,
1591    /// prescribed speed
1592    pub speed: Vec<si::Velocity>,
1593}
1594
1595impl CBTrait for CycleBuilder {
1596    fn with_grade(&mut self, grade: Vec<si::Ratio>) -> anyhow::Result<Cycle> {
1597        let mut cyc: Cycle = self.clone().try_into().with_context(|| format_dbg!())?;
1598        cyc.grade = grade;
1599        Ok(cyc)
1600    }
1601
1602    fn with_temp_amb_air(&mut self, temp_amb_air: Vec<si::Temperature>) -> anyhow::Result<Cycle> {
1603        let mut cyc: Cycle = self.clone().try_into().with_context(|| format_dbg!())?;
1604        cyc.temp_amb_air = temp_amb_air;
1605        Ok(cyc)
1606    }
1607}
1608
1609#[serde_api]
1610#[derive(Default, Debug, Serialize, Deserialize, PartialEq, Clone)]
1611#[non_exhaustive]
1612#[serde(deny_unknown_fields)]
1613#[cfg_attr(feature = "pyo3", pyclass(module = "fastsim", subclass, eq))]
1614/// Element of `Cycle`.  Used for vec-like operations.
1615pub struct CycleElement {
1616    /// simulation time \[s\]
1617    #[serde(alias = "cycSecs")]
1618    pub time: si::Time,
1619    /// simulation power \[W\]
1620    #[serde(alias = "speed_mps", alias = "cycMps")]
1621    pub speed: si::Velocity,
1622    // `dist` is not included here because it is derived in `Init::init`
1623    /// road grade
1624    #[serde(alias = "cycGrade")]
1625    pub grade: Option<si::Ratio>,
1626    // `elev` is not included here because it is derived in `Init::init`
1627    /// road charging/discharing capacity
1628    pub pwr_max_charge: Option<si::Power>,
1629    // TODO: make sure all fields in cycle are represented here, as appropriate
1630    /// ambient air temperature w.r.t. to time (rather than spatial position)
1631    pub temp_amb_air: Option<si::Temperature>,
1632    /// solar heat load w.r.t. to time (rather than spatial position)
1633    pub pwr_solar_load: Option<si::Power>,
1634}
1635
1636impl SerdeAPI for CycleElement {}
1637impl Init for CycleElement {}
1638
1639#[pyo3_api]
1640impl CycleElement {}
1641
1642#[cfg(test)]
1643mod tests {
1644    use super::{manipulation_utils::ConstantJerkTrajectory, *};
1645    /// Build, initialize, and return 2-element cycle
1646    fn mock_cyc_len_2() -> Cycle {
1647        let mut cyc = Cycle {
1648            name: String::new(),
1649            init_elev: None,
1650            time: (0..=2).map(|x| (x as f64) * uc::S).collect(),
1651            speed: (0..=2).map(|x| (x as f64) * uc::MPS).collect(),
1652            dist: vec![],
1653            grade: (0..=2).map(|x| (x as f64 * uc::R) / 100.).collect(),
1654            elev: vec![],
1655            pwr_max_chrg: vec![],
1656            grade_interp: Default::default(),
1657            elev_interp: Default::default(),
1658            temp_amb_air: Default::default(),
1659            pwr_solar_load: Default::default(),
1660        };
1661        cyc.init().unwrap();
1662        cyc
1663    }
1664
1665    fn make_two_triangles_cycle() -> Cycle {
1666        let mut cyc = Cycle {
1667            name: String::from("Two Triangles"),
1668            init_elev: Some(0.0 * uc::M),
1669            time: vec![
1670                0.0 * uc::S,
1671                10.0 * uc::S,
1672                20.0 * uc::S,
1673                30.0 * uc::S,
1674                40.0 * uc::S,
1675                50.0 * uc::S,
1676            ],
1677            speed: vec![
1678                0.0 * uc::MPS,
1679                4.0 * uc::MPS,
1680                0.0 * uc::MPS,
1681                0.0 * uc::MPS,
1682                5.0 * uc::MPS,
1683                0.0 * uc::MPS,
1684            ],
1685            dist: vec![],
1686            grade: vec![
1687                0.0 * uc::R,
1688                0.0 * uc::R,
1689                0.0 * uc::R,
1690                0.0 * uc::R,
1691                0.01 * uc::R,
1692                0.01 * uc::R,
1693            ],
1694            elev: vec![],
1695            pwr_max_chrg: vec![],
1696            grade_interp: Default::default(),
1697            elev_interp: Default::default(),
1698            temp_amb_air: Default::default(),
1699            pwr_solar_load: Default::default(),
1700        };
1701        cyc.init().unwrap();
1702        cyc
1703    }
1704
1705    #[test]
1706    fn test_init() {
1707        let cyc = mock_cyc_len_2();
1708        assert_eq!(
1709            cyc.dist,
1710            [0., 1., 3.] // meters
1711                .iter()
1712                .map(|x| *x * uc::M)
1713                .collect::<Vec<si::Length>>()
1714        );
1715        assert_eq!(
1716            cyc.elev,
1717            [121.92, 121.93, 121.99000000000001] // meters
1718                .iter()
1719                .map(|x| *x * uc::M)
1720                .collect::<Vec<si::Length>>()
1721        );
1722    }
1723
1724    #[test]
1725    fn test_to_elements() {
1726        let cyc = mock_cyc_len_2();
1727        let elements = cyc.to_elements();
1728        assert_eq!(elements.len(), 3);
1729        assert_eq!(elements[0].time, 0.0 * uc::S);
1730        assert_eq!(elements[2].time, cyc.time[2]);
1731        assert_eq!(elements[2].speed, cyc.speed[2]);
1732        assert_eq!(elements[2].grade.unwrap(), 0.02 * uc::R);
1733        assert!(elements[2].pwr_max_charge.is_none());
1734        assert_eq!(elements[2].temp_amb_air.unwrap(), *TE_STD_AIR);
1735        assert!(elements[2].pwr_solar_load.is_none());
1736    }
1737
1738    #[test]
1739    fn test_to_microtrips() {
1740        let cyc = make_two_triangles_cycle();
1741        let actual = cyc.to_microtrips(Some(0.01 * uc::MPH));
1742        assert_eq!(actual.len(), 2);
1743        let cyc0 = &actual[0];
1744        assert_eq!(
1745            cyc0.time,
1746            vec![0.0 * uc::S, 10.0 * uc::S, 20.0 * uc::S, 30.0 * uc::S]
1747        );
1748        assert_eq!(
1749            cyc0.speed,
1750            vec![0.0 * uc::MPS, 4.0 * uc::MPS, 0.0 * uc::MPS, 0.0 * uc::MPS]
1751        );
1752        assert_eq!(
1753            cyc0.grade,
1754            vec![0.0 * uc::R, 0.0 * uc::R, 0.0 * uc::R, 0.0 * uc::R]
1755        );
1756        let cyc1 = &actual[1];
1757        assert_eq!(cyc1.time, vec![0.0 * uc::S, 10.0 * uc::S, 20.0 * uc::S]);
1758        assert_eq!(
1759            cyc1.speed,
1760            vec![0.0 * uc::MPS, 5.0 * uc::MPS, 0.0 * uc::MPS]
1761        );
1762        assert_eq!(cyc1.grade, vec![0.0 * uc::R, 0.01 * uc::R, 0.01 * uc::R]);
1763    }
1764
1765    #[test]
1766    fn test_distance_and_target_speeds_by_microtrip() {
1767        let cyc = make_two_triangles_cycle();
1768        let expected = [
1769            (0.0 * uc::M, (40.0 / 20.0) * uc::MPS),
1770            (40.0 * uc::M, (50.0 / 20.0) * uc::MPS),
1771        ];
1772        let actual = cyc.distance_and_target_speeds_by_microtrip(None, 1.0, 0.0 * uc::MPS);
1773        assert_eq!(actual.len(), expected.len());
1774        for i in 0..expected.len() {
1775            assert_eq!(actual[i].0, expected[i].0);
1776            assert_eq!(actual[i].1, expected[i].1);
1777        }
1778        let expected = [
1779            (0.0 * uc::M, (40.0 / 30.0) * uc::MPS),
1780            (40.0 * uc::M, (50.0 / 20.0) * uc::MPS),
1781        ];
1782        let actual = cyc.distance_and_target_speeds_by_microtrip(None, 0.0, 0.0 * uc::MPS);
1783        assert_eq!(actual.len(), expected.len());
1784        for i in 0..expected.len() {
1785            assert_eq!(actual[i].0, expected[i].0);
1786            assert_eq!(actual[i].1, expected[i].1);
1787        }
1788    }
1789
1790    #[test]
1791    fn test_extending_cycle_time() {
1792        let cyc = make_two_triangles_cycle();
1793        let expected = {
1794            let mut c = Cycle {
1795                name: String::from("Two Triangles"),
1796                init_elev: Some(0.0 * uc::M),
1797                time: vec![
1798                    0.0 * uc::S,
1799                    10.0 * uc::S,
1800                    20.0 * uc::S,
1801                    30.0 * uc::S,
1802                    40.0 * uc::S,
1803                    50.0 * uc::S,
1804                    51.0 * uc::S,
1805                    52.0 * uc::S,
1806                    53.0 * uc::S,
1807                    54.0 * uc::S,
1808                    55.0 * uc::S,
1809                    56.0 * uc::S,
1810                    57.0 * uc::S,
1811                    58.0 * uc::S,
1812                ],
1813                speed: vec![
1814                    0.0 * uc::MPS,
1815                    4.0 * uc::MPS,
1816                    0.0 * uc::MPS,
1817                    0.0 * uc::MPS,
1818                    5.0 * uc::MPS,
1819                    0.0 * uc::MPS,
1820                    0.0 * uc::MPS,
1821                    0.0 * uc::MPS,
1822                    0.0 * uc::MPS,
1823                    0.0 * uc::MPS,
1824                    0.0 * uc::MPS,
1825                    0.0 * uc::MPS,
1826                    0.0 * uc::MPS,
1827                    0.0 * uc::MPS,
1828                ],
1829                dist: vec![],
1830                grade: vec![
1831                    0.0 * uc::R,
1832                    0.0 * uc::R,
1833                    0.0 * uc::R,
1834                    0.0 * uc::R,
1835                    0.01 * uc::R,
1836                    0.01 * uc::R,
1837                    0.0 * uc::R,
1838                    0.0 * uc::R,
1839                    0.0 * uc::R,
1840                    0.0 * uc::R,
1841                    0.0 * uc::R,
1842                    0.0 * uc::R,
1843                    0.0 * uc::R,
1844                    0.0 * uc::R,
1845                ],
1846                elev: vec![],
1847                pwr_max_chrg: vec![],
1848                grade_interp: Default::default(),
1849                elev_interp: Default::default(),
1850                temp_amb_air: Default::default(),
1851                pwr_solar_load: Default::default(),
1852            };
1853            c.init().unwrap();
1854            c
1855        };
1856        let absolute_time = Some(3.0 * uc::S);
1857        let time_fraction = Some(0.10 * uc::R);
1858        // extend by 3 s and 10% of existing time (i.e., 5 s)
1859        // = extend by 8 s
1860        let actual = cyc.extend_time(absolute_time, time_fraction);
1861        assert_eq!(actual, expected);
1862    }
1863
1864    /// Round the given number n to the given number of digits
1865    /// - n: the number to round
1866    /// - digits: the digits to round or defaults to 2; if not positive,
1867    fn round(n: f64, digits: Option<i32>) -> f64 {
1868        let digits = digits.unwrap_or(2);
1869        let digits = if digits < 0 { 0 } else { digits };
1870        let multiplier = 10.0_f64.powi(digits);
1871        (n * multiplier).round() / multiplier
1872    }
1873
1874    #[test]
1875    fn cycle_step_distances_are_as_expected() {
1876        let c = make_two_triangles_cycle();
1877        let expected = [
1878            0.0 * uc::M,
1879            20.0 * uc::M,
1880            20.0 * uc::M,
1881            0.0 * uc::M,
1882            25.0 * uc::M,
1883            25.0 * uc::M,
1884        ];
1885        let actual = c.trapz_step_distances();
1886        assert_eq!(actual.len(), expected.len());
1887        for i in 0..expected.len() {
1888            assert_eq!(actual[i], expected[i], "differ at step {i}");
1889        }
1890    }
1891
1892    #[test]
1893    fn cycle_elevations_are_as_expected() {
1894        let c = make_two_triangles_cycle();
1895        let dh = 0.01_f64.atan().cos() * 25.0_f64 * 0.01_f64;
1896        let expected = [
1897            0.0 * uc::M,
1898            0.0 * uc::M,
1899            0.0 * uc::M,
1900            0.0 * uc::M,
1901            dh * uc::M,
1902            dh * uc::M,
1903        ];
1904        let actual = c.trapz_step_elevations();
1905        assert_eq!(actual.len(), expected.len());
1906        for i in 0..expected.len() {
1907            assert_eq!(actual[i], expected[i], "differ at step {i}");
1908        }
1909    }
1910
1911    #[test]
1912    fn cycle_cache_yields_same_results() {
1913        let c = make_two_triangles_cycle();
1914        let cache = c.build_cache();
1915        let dist_m = 0.0;
1916        let e0_expected = 0.0;
1917        let e0_actual = cache.interp_elevation(dist_m);
1918        assert_eq!(e0_actual, e0_expected);
1919        let dist_m = 65.0;
1920        let e1_expected = 0.01_f64.atan().cos() * 25.0_f64 * 0.01_f64;
1921        let e1_actual = cache.interp_elevation(dist_m);
1922        assert_eq!(e1_actual, e1_expected);
1923    }
1924
1925    #[test]
1926    fn average_grade_over_range_is_correct() {
1927        let c = make_two_triangles_cycle();
1928        let cache = c.build_cache();
1929        let d0 = 40.0 * uc::M;
1930        let dd = 50.0 * uc::M;
1931        let expected0 = 0.01 * uc::R;
1932        let actual00 = c.average_grade_over_range(d0, dd, None);
1933        let actual00 = round(actual00.get::<si::ratio>(), Some(6)) * uc::R;
1934        assert_eq!(actual00, expected0);
1935        let actual01 = c.average_grade_over_range(d0, dd, Some(&cache));
1936        let actual01 = round(actual01.get::<si::ratio>(), Some(6)) * uc::R;
1937        assert_eq!(actual01, expected0);
1938    }
1939
1940    #[test]
1941    fn distance_to_next_stop_is_correct() {
1942        let c = make_two_triangles_cycle();
1943        let cache = c.build_cache();
1944        let d = 20.0 * uc::M;
1945        let expected = 20.0 * uc::M;
1946        let actual = c.calc_distance_to_next_stop_from(d, None);
1947        assert_eq!(actual, expected);
1948        let actual = c.calc_distance_to_next_stop_from(d, Some(&cache));
1949        assert_eq!(actual, expected);
1950        let d = 65.0 * uc::M;
1951        let expected = 25.0 * uc::M;
1952        let actual = c.calc_distance_to_next_stop_from(d, None);
1953        assert_eq!(actual, expected);
1954        let actual = c.calc_distance_to_next_stop_from(d, Some(&cache));
1955        assert_eq!(actual, expected);
1956        let d = 0.0 * uc::M;
1957        let expected = 40.0 * uc::M;
1958        let actual = c.calc_distance_to_next_stop_from(d, None);
1959        assert_eq!(actual, expected);
1960        let actual = c.calc_distance_to_next_stop_from(d, Some(&cache));
1961        assert_eq!(actual, expected);
1962    }
1963
1964    #[test]
1965    fn modifying_a_cycle_with_trajectory() {
1966        let c0 = make_two_triangles_cycle();
1967        let mut c = c0.clone();
1968        let n = 3;
1969        let d0 = 20.0; // units: m
1970        let v0 = 4.0; // units: m/s
1971        let dr = 65.0; // units: m
1972        let vr = 5.0; // units: m/s
1973        let dt = 10.0; // units: s
1974        let traj = ConstantJerkTrajectory::from_speed_and_distance_targets(n, d0, v0, dr, vr, dt);
1975        c.modify_by_const_jerk_trajectory(
1976            2,
1977            n,
1978            traj.jerk_m_per_s3 * uc::MPS3,
1979            traj.acceleration_m_per_s2 * uc::MPS2,
1980        );
1981        let expected = {
1982            let mut cyc = Cycle {
1983                name: String::from("Two Triangles"),
1984                init_elev: Some(0.0 * uc::M),
1985                time: vec![
1986                    0.0 * uc::S,
1987                    10.0 * uc::S,
1988                    20.0 * uc::S,
1989                    30.0 * uc::S,
1990                    40.0 * uc::S,
1991                    50.0 * uc::S,
1992                ],
1993                speed: vec![
1994                    0.0 * uc::MPS,
1995                    4.0 * uc::MPS,
1996                    traj.speed_at_step(1) * uc::MPS,
1997                    traj.speed_at_step(2) * uc::MPS,
1998                    5.0 * uc::MPS,
1999                    0.0 * uc::MPS,
2000                ],
2001                dist: vec![],
2002                grade: vec![
2003                    0.0 * uc::R,
2004                    0.0 * uc::R,
2005                    0.0 * uc::R,
2006                    0.0 * uc::R,
2007                    0.01 * uc::R,
2008                    0.01 * uc::R,
2009                ],
2010                elev: vec![],
2011                pwr_max_chrg: vec![],
2012                grade_interp: Default::default(),
2013                elev_interp: Default::default(),
2014                temp_amb_air: Default::default(),
2015                pwr_solar_load: Default::default(),
2016            };
2017            cyc.init().expect("initializaiton should not throw");
2018            cyc
2019        };
2020        assert_eq!(c.time.len(), expected.time.len());
2021        assert_eq!(c.speed.len(), expected.speed.len());
2022        assert_eq!(c.dist.len(), expected.dist.len());
2023        assert_eq!(c.grade.len(), expected.grade.len());
2024        for idx in 0..c.speed.len() {
2025            assert_eq!(c.time[idx], expected.time[idx]);
2026            assert_eq!(c.speed[idx], expected.speed[idx]);
2027            assert_eq!(c.dist[idx], expected.dist[idx]);
2028            assert_eq!(c.grade[idx], expected.grade[idx]);
2029        }
2030    }
2031
2032    #[test]
2033    pub fn modify_with_braking_trajectory() {
2034        let mut actual = {
2035            let mut cyc = Cycle {
2036                name: String::from("Test"),
2037                init_elev: Some(0.0 * uc::M),
2038                time: vec![
2039                    0.0 * uc::S,
2040                    1.0 * uc::S,
2041                    2.0 * uc::S,
2042                    3.0 * uc::S,
2043                    4.0 * uc::S,
2044                    5.0 * uc::S,
2045                ],
2046                speed: vec![
2047                    0.0 * uc::MPS,
2048                    4.0 * uc::MPS,
2049                    4.0 * uc::MPS,
2050                    1.0 * uc::MPS,
2051                    1.0 * uc::MPS,
2052                    0.0 * uc::MPS,
2053                ],
2054                dist: vec![],
2055                grade: vec![],
2056                elev: vec![],
2057                pwr_max_chrg: vec![],
2058                grade_interp: Default::default(),
2059                elev_interp: Default::default(),
2060                temp_amb_air: Default::default(),
2061                pwr_solar_load: Default::default(),
2062            };
2063            cyc.init().expect("initializaiton should not throw");
2064            cyc
2065        };
2066        let precision = Some(6);
2067        let (v_end, n_steps) =
2068            actual.modify_with_braking_trajectory((-4.0 / 3.0) * uc::MPS2, 3, Some(4.0 * uc::M));
2069        let v_end = round(v_end.get::<si::meter_per_second>(), precision);
2070        assert_eq!(v_end, 0.0);
2071        assert_eq!(n_steps, 3);
2072        let expected = {
2073            let n = 3;
2074            let d0 = 0.0;
2075            let v0 = 4.0;
2076            let dr = 4.0;
2077            let vr = 0.0;
2078            let dt = 1.0;
2079            let traj =
2080                ConstantJerkTrajectory::from_speed_and_distance_targets(n, d0, v0, dr, vr, dt);
2081            let mut cyc = Cycle {
2082                name: String::from("Test"),
2083                init_elev: Some(0.0 * uc::M),
2084                time: vec![
2085                    0.0 * uc::S,
2086                    1.0 * uc::S,
2087                    2.0 * uc::S,
2088                    3.0 * uc::S,
2089                    4.0 * uc::S,
2090                    5.0 * uc::S,
2091                ],
2092                speed: vec![
2093                    0.0 * uc::MPS,
2094                    4.0 * uc::MPS,
2095                    4.0 * uc::MPS,
2096                    traj.speed_at_step(1) * uc::MPS,
2097                    traj.speed_at_step(2) * uc::MPS,
2098                    traj.speed_at_step(3) * uc::MPS,
2099                ],
2100                dist: vec![],
2101                grade: vec![],
2102                elev: vec![],
2103                pwr_max_chrg: vec![],
2104                grade_interp: Default::default(),
2105                elev_interp: Default::default(),
2106                temp_amb_air: Default::default(),
2107                pwr_solar_load: Default::default(),
2108            };
2109            cyc.init().expect("initializaiton should not throw");
2110            cyc
2111        };
2112        assert_eq!(actual.time.len(), expected.time.len());
2113        for i in 0..actual.time.len() {
2114            let at = round(actual.time[i].get::<si::second>(), precision);
2115            let et = round(expected.time[i].get::<si::second>(), precision);
2116            let av = round(actual.speed[i].get::<si::meter_per_second>(), precision);
2117            let ev = round(expected.speed[i].get::<si::meter_per_second>(), precision);
2118            let ad = round(actual.dist[i].get::<si::meter>(), precision);
2119            let ed = round(expected.dist[i].get::<si::meter>(), precision);
2120            assert_eq!(at, et, "time@t={et}&i={i}");
2121            assert_eq!(av, ev, "speed@t={et}&i={i}");
2122            assert_eq!(ad, ed, "dist@t={et}&i={i}");
2123        }
2124    }
2125
2126    #[test]
2127    pub fn test_trim() {
2128        let c = make_two_triangles_cycle();
2129        let cyc = c.extend_time(Some(10.0 * uc::S), None);
2130        let dt_idle = cyc.ending_idle_time();
2131        assert_eq!(dt_idle, 10.0 * uc::S);
2132        // NOTE: extend_time adds time by 1.0 s increments so 10 points
2133        assert_eq!(cyc.time.len(), c.time.len() + 10);
2134        assert_eq!(*cyc.time.iter().last().unwrap(), 60.0 * uc::S);
2135        let cyc_trimmed = cyc.trim_ending_idle(None);
2136        assert_eq!(cyc_trimmed.time.len(), c.time.len());
2137    }
2138    type StructWithResources = Cycle;
2139
2140    #[test]
2141    fn test_resources() {
2142        let resource_list = StructWithResources::list_resources().unwrap();
2143        assert!(!resource_list.is_empty());
2144
2145        // verify that resources can all load
2146        for resource in resource_list {
2147            StructWithResources::from_resource(resource.clone(), false)
2148                .with_context(|| format_dbg!(resource))
2149                .unwrap();
2150        }
2151    }
2152
2153    #[test]
2154    fn test_resample() {
2155        let cyc0 = {
2156            let mut c = Cycle {
2157                name: String::from("a test"),
2158                time: vec![0.0 * uc::S, 10.0 * uc::S, 20.0 * uc::S],
2159                speed: vec![0.0 * uc::MPS, 10.0 * uc::MPS, 0.0 * uc::MPS],
2160                grade: vec![0.01 * uc::R, 0.01 * uc::R, -0.01 * uc::R],
2161                init_elev: None,
2162                dist: vec![],
2163                elev: vec![],
2164                pwr_max_chrg: vec![],
2165                temp_amb_air: vec![],
2166                pwr_solar_load: vec![],
2167                grade_interp: None,
2168                elev_interp: None,
2169            };
2170            c.init().unwrap();
2171            c
2172        };
2173        let cyc1 = cyc0.resample(1.0 * uc::S);
2174        assert_eq!(21, cyc1.time.len());
2175        assert_eq!(
2176            cyc1.time[cyc1.time.len() - 1],
2177            cyc0.time[cyc0.time.len() - 1]
2178        );
2179        assert_eq!(cyc1.time[0], cyc0.time[0]);
2180        assert_eq!(cyc1.time[0], 0.0 * uc::S);
2181        assert_eq!(cyc1.time[5], 5.0 * uc::S);
2182        assert_eq!(cyc1.speed[5], 5.0 * uc::MPS);
2183        assert_eq!(cyc1.grade[5], 0.01 * uc::R);
2184        assert_eq!(cyc1.time[10], 10.0 * uc::S);
2185        assert_eq!(cyc1.speed[10], 10.0 * uc::MPS);
2186        assert_eq!(cyc1.grade[10], 0.01 * uc::R);
2187        assert_eq!(cyc1.time[11], 11.0 * uc::S);
2188        assert_eq!(cyc1.speed[11], 9.0 * uc::MPS);
2189        assert_eq!(cyc1.grade[11], -0.01 * uc::R);
2190        assert_eq!(cyc1.time[20], 20.0 * uc::S);
2191        assert_eq!(cyc1.speed[20], 0.0 * uc::MPS);
2192        assert_eq!(cyc1.grade[20], -0.01 * uc::R);
2193    }
2194}
2195
2196lazy_static! {
2197    pub static ref CYC_ACCEL: Cycle = Cycle::try_from(CycleBuilder {
2198        name: String::from("accel test"),
2199        time: (0..300)
2200            .map(|t| (t as f64) * uc::S)
2201            .collect::<Vec<si::Time>>(),
2202        speed: vec![90.0 * uc::MPH; 300],
2203    })
2204    .unwrap();
2205}