|
6 | 6 | from numbers import Number |
7 | 7 | import numpy as np |
8 | 8 | from copy import deepcopy |
| 9 | +import itertools |
9 | 10 | from warnings import warn |
10 | 11 | from collections import abc, namedtuple |
11 | 12 | from .sim_result import SimResult, LazySimResult |
@@ -419,7 +420,7 @@ def __next_state(self, x, u, dt) -> dict: |
419 | 420 | """ |
420 | 421 |
|
421 | 422 | # Calculate next state and add process noise |
422 | | - next_state = self.apply_process_noise(self.next_state(x, u, dt)) |
| 423 | + next_state = self.apply_process_noise(self.next_state(x, u, dt), dt) |
423 | 424 |
|
424 | 425 | # Apply Limits |
425 | 426 | return self.apply_limits(next_state) |
@@ -744,9 +745,10 @@ def simulate_to_threshold(self, future_loading_eqn, first_output = None, thresho |
744 | 745 | raise ProgModelInputException("'dt' must be a number or function, was a {}".format(type(config['dt']))) |
745 | 746 | if isinstance(config['dt'], Number) and config['dt'] < 0: |
746 | 747 | raise ProgModelInputException("'dt' must be positive, was {}".format(config['dt'])) |
747 | | - if not isinstance(config['save_freq'], Number): |
| 748 | + if not isinstance(config['save_freq'], Number) and not isinstance(config['save_freq'], tuple): |
748 | 749 | raise ProgModelInputException("'save_freq' must be a number, was a {}".format(type(config['save_freq']))) |
749 | | - if config['save_freq'] <= 0: |
| 750 | + if (isinstance(config['save_freq'], Number) and config['save_freq'] <= 0) or \ |
| 751 | + (isinstance(config['save_freq'], tuple) and config['save_freq'][1] <= 0): |
750 | 752 | raise ProgModelInputException("'save_freq' must be positive, was {}".format(config['save_freq'])) |
751 | 753 | if not isinstance(config['save_pts'], abc.Iterable): |
752 | 754 | raise ProgModelInputException("'save_pts' must be list or array, was a {}".format(type(config['save_pts']))) |
@@ -795,9 +797,20 @@ def check_thresholds(thresholds_met): |
795 | 797 | saved_states = [] |
796 | 798 | saved_outputs = [] |
797 | 799 | saved_event_states = [] |
798 | | - save_freq = config['save_freq'] |
799 | 800 | horizon = t+config['horizon'] |
800 | | - next_save = t+save_freq |
| 801 | + if isinstance(config['save_freq'], tuple): |
| 802 | + # Tuple used to specify start and frequency |
| 803 | + t_step = config['save_freq'][1] |
| 804 | + # Use starting time or the next multiple |
| 805 | + t_start = config['save_freq'][0] |
| 806 | + start = max(t_start, t - (t-t_start)%t_step) |
| 807 | + iterator = itertools.count(start, t_step) |
| 808 | + else: |
| 809 | + # Otherwise - start is t0 |
| 810 | + t_step = config['save_freq'] |
| 811 | + iterator = itertools.count(t, t_step) |
| 812 | + next(iterator) # Skip current time |
| 813 | + next_save = next(iterator) |
801 | 814 | save_pt_index = 0 |
802 | 815 | save_pts = config['save_pts'] |
803 | 816 | save_pts.append(1e99) # Add last endpoint |
@@ -839,13 +852,15 @@ def next_time(t, x): |
839 | 852 |
|
840 | 853 | while t < horizon: |
841 | 854 | dt = next_time(t, x) |
842 | | - t = t + dt |
| 855 | + t = t + dt/2 |
| 856 | + # Use state at midpoint of step to best represent the load during the duration of the step |
843 | 857 | u = future_loading_eqn(t, x) |
| 858 | + t = t + dt/2 |
844 | 859 | x = next_state(x, u, dt) |
845 | 860 |
|
846 | 861 | # Save if at appropriate time |
847 | 862 | if (t >= next_save): |
848 | | - next_save += save_freq |
| 863 | + next_save = next(iterator) |
849 | 864 | update_all() |
850 | 865 | if (t >= save_pts[save_pt_index]): |
851 | 866 | save_pt_index += 1 |
|
0 commit comments