Skip to content

Commit 029c57a

Browse files
phmbressanphmbressan
authored
Merge pull request #330 from Space-Enterprise-at-Berkeley/enh/function-updates
Added Piecewise class and other new Function methods
2 parents a9f9fab + 66bd7e6 commit 029c57a

File tree

1 file changed

+255
-5
lines changed

1 file changed

+255
-5
lines changed

‎rocketpy/Function.py

Lines changed: 255 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -8,7 +8,7 @@
88

99
import matplotlib.pyplot as plt
1010
import numpy as np
11-
from scipy import integrate, linalg
11+
from scipy import integrate, linalg, optimize
1212

1313

1414
class Function:
@@ -2024,9 +2024,259 @@ def integral(self, a, b, numerical=False):
20242024
ans, _ = integrate.quad(self, a, b, epsabs=0.1, limit=10000)
20252025
return ans
20262026

2027-
# Not implemented
20282027
def differentiate(self, x, dx=1e-6):
2028+
"""Evaluates the derivative of a 1-D Function at x with a step size of dx.
2029+
2030+
Parameters
2031+
----------
2032+
x : float
2033+
Point at which to evaluate the derivative.
2034+
dx : float
2035+
Step size to use in the derivative calculation. Default is 1e-6.
2036+
2037+
2038+
Returns
2039+
-------
2040+
ans : float
2041+
Evaluated derivative.
2042+
"""
20292043
return (self.getValue(x + dx) - self.getValue(x - dx)) / (2 * dx)
2030-
# h = (10)**-300
2031-
# z = x + h*1j
2032-
# return self(z).imag/h
2044+
2045+
def derivativeFunction(self):
2046+
"""Returns a Function object which gives the derivative of the Function object.
2047+
2048+
Returns
2049+
-------
2050+
result : Function
2051+
A Function object which gives the derivative of self.
2052+
"""
2053+
# Check if Function object source is array
2054+
if isinstance(self.source, np.ndarray):
2055+
# Operate on grid values
2056+
Ys = np.diff(self.source[:, 1]) / np.diff(self.source[:, 0])
2057+
Xs = self.source[:-1, 0] + np.diff(self.source[:, 0]) / 2
2058+
source = np.concatenate(([Xs], [Ys])).transpose()
2059+
2060+
# Retrieve inputs, outputs and interpolation
2061+
inputs = self.__inputs__[:]
2062+
outputs = "d(" + self.__outputs__[0] + ")/d(" + inputs[0] + ")"
2063+
outputs = "(" + outputs + ")"
2064+
interpolation = "linear"
2065+
2066+
# Create new Function object
2067+
return Function(source, inputs, outputs, interpolation)
2068+
else:
2069+
return Function(lambda x: self.differentiate(x))
2070+
2071+
def integralFunction(self, lower=None, upper=None, datapoints=100):
2072+
"""Returns a Function object representing the integral of the Function object.
2073+
2074+
Parameters
2075+
----------
2076+
lower : scalar, optional
2077+
The lower limit of the interval in which the function is to be
2078+
plotted. If the Function is given by a dataset, the default
2079+
value is the start of the dataset.
2080+
upper : scalar, optional
2081+
The upper limit of the interval in which the function is to be
2082+
plotted. If the Function is given by a dataset, the default
2083+
value is the end of the dataset.
2084+
datapoints : int, optional
2085+
The number of points in which the integral will be evaluated for
2086+
plotting it, which draws lines between each evaluated point.
2087+
The default value is 100.
2088+
2089+
Returns
2090+
-------
2091+
result : Function
2092+
The integral of the Function object.
2093+
"""
2094+
# Check if lower and upper are given
2095+
if lower is None:
2096+
if isinstance(self.source, np.ndarray):
2097+
lower = self.source[0, 0]
2098+
else:
2099+
raise ValueError("Lower limit must be given if source is a function.")
2100+
if upper is None:
2101+
if isinstance(self.source, np.ndarray):
2102+
upper = self.source[-1, 0]
2103+
else:
2104+
raise ValueError("Upper limit must be given if source is a function.")
2105+
2106+
# Create a new Function object
2107+
xData = np.linspace(lower, upper, datapoints)
2108+
yData = np.zeros(datapoints)
2109+
for i in range(datapoints):
2110+
yData[i] = self.integral(lower, xData[i])
2111+
return Function(np.concatenate(([xData], [yData])).transpose(),
2112+
inputs=self.__inputs__,
2113+
outputs=[o + " Integral" for o in self.__outputs__])
2114+
2115+
def inverseFunction(self, lower=None, upper=None, datapoints=100):
2116+
"""
2117+
Returns the inverse of the Function. The inverse function of F is a function that undoes the operation of F. The
2118+
inverse of F exists if and only if F is bijective. Makes the domain the range and the range the domain.
2119+
2120+
Parameters
2121+
----------
2122+
lower : float
2123+
Lower limit of the new domain. Only required if the Function's source is a callable instead of a list of points.
2124+
upper : float
2125+
Upper limit of the new domain. Only required if the Function's source is a callable instead of a list of points.
2126+
2127+
Returns
2128+
-------
2129+
result : Function
2130+
A Function whose domain and range have been inverted.
2131+
"""
2132+
if isinstance(self.source, np.ndarray):
2133+
# Swap the columns
2134+
source = np.concatenate(([self.source[:, 1]], [self.source[:, 0]])).transpose()
2135+
2136+
return Function(
2137+
source,
2138+
inputs=self.__outputs__,
2139+
outputs=self.__inputs__,
2140+
interpolation=self.__interpolation__
2141+
)
2142+
else:
2143+
return Function(
2144+
lambda x: self.findOptimalInput(x),
2145+
inputs=self.__outputs__,
2146+
outputs=self.__inputs__,
2147+
interpolation=self.__interpolation__
2148+
)
2149+
2150+
def findOptimalInput(self, val):
2151+
"""
2152+
Finds the optimal input for a given output.
2153+
2154+
Parameters
2155+
----------
2156+
val : float
2157+
The value of the output.
2158+
2159+
Returns
2160+
-------
2161+
result : ndarray
2162+
The value of the input which gives the output closest to val.
2163+
"""
2164+
return optimize.fmin(lambda x: np.abs(self.getValue(x) - val), 0, ftol=1e-6, disp=False)
2165+
2166+
def compose(self, func, lower=None, upper=None, datapoints=100):
2167+
"""
2168+
Returns a Function object which is the result of inputing a function into a function
2169+
(i.e. f(g(x))). The domain will become the domain of the input function and the range
2170+
will become the range of the original function.
2171+
2172+
Parameters
2173+
----------
2174+
func : Function
2175+
The function to be inputed into the function.
2176+
lower : float
2177+
Lower limit of the new domain.
2178+
upper : float
2179+
Upper limit of the new domain.
2180+
2181+
Returns
2182+
-------
2183+
result : Function
2184+
The result of inputing the function into the function.
2185+
"""
2186+
# Check if the input is a function
2187+
if not isinstance(func, Function):
2188+
raise TypeError("Input must be a Function object.")
2189+
2190+
# Checks to make sure lower bound is given
2191+
# If not it will start at the higher of the two lower bounds
2192+
if lower is None:
2193+
if isinstance(self.source, np.ndarray):
2194+
lower = self.source[0, 0]
2195+
if isinstance(func.source, np.ndarray):
2196+
lower = func.source[0, 0] if lower is None else max(lower, func.source[0, 0])
2197+
if lower is None:
2198+
raise ValueError("If Functions.source is a <class 'function'>, must provide bounds")
2199+
2200+
# Checks to make sure upper bound is given
2201+
# If not it will end at the lower of the two upper bounds
2202+
if upper is None:
2203+
if isinstance(self.source, np.ndarray):
2204+
upper = self.source[-1, 0]
2205+
if isinstance(func.source, np.ndarray):
2206+
upper = func.source[-1, 0] if upper is None else min(upper, func.source[-1, 0])
2207+
if upper is None:
2208+
raise ValueError("If Functions.source is a <class 'function'>, must provide bounds")
2209+
2210+
# Create a new Function object
2211+
xData = np.linspace(lower, upper, datapoints)
2212+
yData = np.zeros(datapoints)
2213+
for i in range(datapoints):
2214+
yData[i] = self.getValue(func.getValue(xData[i]))
2215+
return Function(np.concatenate(([xData], [yData])).T,
2216+
inputs=func.__inputs__,
2217+
outputs=self.__outputs__,
2218+
interpolation=self.__interpolation__,
2219+
extrapolation=self.__extrapolation__)
2220+
2221+
2222+
class PiecewiseFunction(Function):
2223+
def __new__(cls, source, inputs=["Scalar"], outputs=["Scalar"], interpolation="akima", extrapolation=None, datapoints=50):
2224+
"""
2225+
Creates a piecewise function from a dictionary of functions. The keys of the dictionary
2226+
must be tuples that represent the domain of the function. The domains must be disjoint.
2227+
The piecewise function will be evaluated at datapoints points to create Function object.
2228+
2229+
Parameters
2230+
----------
2231+
source: dictionary
2232+
A dictionary of Function objects, where the keys are the domains.
2233+
inputs : list
2234+
A list of strings that represent the inputs of the function.
2235+
outputs: list
2236+
A list of strings that represent the outputs of the function.
2237+
interpolation: str
2238+
The type of interpolation to use. The default value is 'akima'.
2239+
extrapolation: str
2240+
The type of extrapolation to use. The default value is None.
2241+
datapoints: int
2242+
The number of points in which the piecewise function will be
2243+
evaluated to create a base function. The default value is 100.
2244+
"""
2245+
# Check if source is a dictionary
2246+
if not isinstance(source, dict):
2247+
raise TypeError("source must be a dictionary")
2248+
2249+
# Check if all keys are tuples
2250+
for key in source.keys():
2251+
if not isinstance(key, tuple):
2252+
raise TypeError("keys of source must be tuples")
2253+
2254+
# Check if all domains are disjoint
2255+
for key1 in source.keys():
2256+
for key2 in source.keys():
2257+
if key1 != key2:
2258+
if key1[0] < key2[1] and key1[1] > key2[0]:
2259+
raise ValueError("domains must be disjoint")
2260+
2261+
# Crate Function
2262+
def calcOutput(func, inputs):
2263+
o = np.zeros(len(inputs))
2264+
for j in range(len(inputs)):
2265+
o[j] = func.getValue(inputs[j])
2266+
return o
2267+
2268+
inputData = []
2269+
outputData = []
2270+
for key in sorted(source.keys()):
2271+
i = np.linspace(key[0], key[1], datapoints)
2272+
i = i[~np.in1d(i, inputData)]
2273+
inputData = np.concatenate((inputData, i))
2274+
2275+
f = Function(source[key])
2276+
outputData = np.concatenate((outputData, calcOutput(f, i)))
2277+
2278+
return Function(np.concatenate(([inputData], [outputData])).T,
2279+
inputs=inputs,
2280+
outputs=outputs,
2281+
interpolation=interpolation,
2282+
extrapolation=extrapolation)

0 commit comments

Comments
 (0)