Package with mathematical functions
Information
This package contains functions for commonly used
mathematical operations. The functions are used in
the blocks
Buildings.Utilities.Math.
Extends from Modelica.Icons.VariantsPackage (Icon for package containing variants).
Package Content
Average of a vector
Information
This block outputs the average of the vector.
Inputs
Type | Name | Default | Description |
Integer | nin | | Number of inputs |
Real | u[nin] | | Input vector |
Outputs
Type | Name | Description |
Real | y | Result |
Modelica definition
function average "Average of a vector"
input Integer nin "Number of inputs";
input Real u[nin] "Input vector";
output Real y "Result";
algorithm
y := sum(u)/nin;
end average;
Bicubic function
Information
This function computes
y = a1
+ a2 x1 + a3 x12
+ a4 x2 + a5 x22
+ a6 x1 x2
+ a7 x1^3
+ a8 x2^3
+ a9 x12 x2
+ a10 x1 x22
Inputs
Type | Name | Default | Description |
Real | a[10] | | Coefficients |
Real | x1 | | Independent variable |
Real | x2 | | Independent variable |
Outputs
Type | Name | Description |
Real | y | Result |
Modelica definition
function bicubic "Bicubic function"
input Real a[10] "Coefficients";
input Real x1 "Independent variable";
input Real x2 "Independent variable";
output Real y "Result";
protected
Real x1Sq "= x1^2";
Real x2Sq "= x2^2";
algorithm
x1Sq :=x1*x1;
x2Sq :=x2*x2;
y := a[1] + a[2] * x1 + a[3] * x1^2
+ a[4] * x2 + a[5] * x2^2
+ a[6] * x1 * x2
+ a[7] * x1Sq * x1
+ a[8] * x2Sq * x2
+ a[9] * x1Sq * x2
+ a[10] * x1 * x2Sq;
end bicubic;
Biquadratic function
Information
This function computes
y = a1 + a2 x1
+ a3 x12
+ a4 x2 + a5 x22
+ a6 x1 x2
Inputs
Type | Name | Default | Description |
Real | a[6] | | Coefficients |
Real | x1 | | Independent variable |
Real | x2 | | Independent variable |
Outputs
Type | Name | Description |
Real | y | Result |
Modelica definition
function biquadratic "Biquadratic function"
input Real a[6] "Coefficients";
input Real x1 "Independent variable";
input Real x2 "Independent variable";
output Real y "Result";
algorithm
y :=a[1] + x1*(a[2] + a[3]*x1) + x2*(a[4]+ a[5]*x2) + a[6]*x1*x2;
end biquadratic;
Replicates Boolean signals
Information
This function replicates the boolean input signal to an array of nout
identical output signals.
Inputs
Type | Name | Default | Description |
Integer | nout | 1 | Number of outouts |
Boolean | u | | Boolean input signal |
Outputs
Type | Name | Description |
Boolean | y[nout] | Boolean output signals |
Modelica definition
function booleanReplicator "Replicates Boolean signals"
input Integer nout=1 "Number of outouts";
input Boolean u "Boolean input signal";
output Boolean y[nout] "Boolean output signals";
algorithm
y :=fill(u, nout);
end booleanReplicator;
Interpolate using a cubic Hermite spline with linear extrapolation
Information
For x1 < x < x2, this function interpolates
using cubic hermite spline. For x outside this interval, the function
linearly extrapolates.
For how to use this function, see
Buildings.Utilities.Math.Functions.Examples.CubicHermite.
Inputs
Type | Name | Default | Description |
Real | x | | Abscissa value |
Real | x1 | | Lower abscissa value |
Real | x2 | | Upper abscissa value |
Real | y1 | | Lower ordinate value |
Real | y2 | | Upper ordinate value |
Real | y1d | | Lower gradient |
Real | y2d | | Upper gradient |
Outputs
Type | Name | Description |
Real | y | Interpolated ordinate value |
Modelica definition
function cubicHermiteLinearExtrapolation
"Interpolate using a cubic Hermite spline with linear extrapolation"
input Real x "Abscissa value";
input Real x1 "Lower abscissa value";
input Real x2 "Upper abscissa value";
input Real y1 "Lower ordinate value";
input Real y2 "Upper ordinate value";
input Real y1d "Lower gradient";
input Real y2d "Upper gradient";
output Real y "Interpolated ordinate value";
algorithm
if (x > x1 and x < x2) then
y:=Modelica.Fluid.Utilities.cubicHermite(
x=x,
x1=x1,
x2=x2,
y1=y1,
y2=y2,
y1d=y1d,
y2d=y2d);
elseif x <= x1 then
// linear extrapolation
y:=y1 + (x - x1)*y1d;
else
y:=y2 + (x - x2)*y2d;
end if;
end cubicHermiteLinearExtrapolation;
Replicates integer signals
Information
This function replicates the integer input signal to an array of
nout
identical output signals.
Inputs
Type | Name | Default | Description |
Integer | nout | | Number of outputs |
Integer | u | | Integer input signal |
Outputs
Type | Name | Description |
Integer | y[nout] | Integer output signals |
Modelica definition
function integerReplicator "Replicates integer signals"
input Integer nout "Number of outputs";
input Integer u "Integer input signal";
output Integer y[nout] "Integer output signals";
algorithm
y :=fill(u, nout);
end integerReplicator;
Function that approximates 1/x by a twice continuously differentiable function
Information
Function that approximates y=1 ⁄ x
inside the interval -δ ≤ x ≤ δ.
The approximation is twice continuously differentiable with a bounded derivative on the whole
real line.
See the package Examples
for the graph.
Inputs
Type | Name | Default | Description |
Real | x | | Abscissa value |
Real | delta | | Abscissa value below which approximation occurs |
Outputs
Type | Name | Description |
Real | y | Function value |
Modelica definition
function inverseXRegularized
"Function that approximates 1/x by a twice continuously differentiable function"
input Real x "Abscissa value";
input Real delta(min=0) "Abscissa value below which approximation occurs";
output Real y "Function value";
protected
Real delta2 "Delta^2";
Real x2_d2 "=x^2/delta^2";
algorithm
if (abs(x) > delta) then
y := 1/x;
else
delta2 :=delta*delta;
x2_d2 := x*x/delta2;
y := x/delta2 + x*abs(x/delta2/delta*(2 - x2_d2*(3 - x2_d2)));
end if;
end inverseXRegularized;
Returns true if the argument is a monotonic sequence
Information
This function returns true
if its argument is
monotonic increasing or decreasing, and false
otherwise.
If strict=true
, then strict monotonicity is tested,
otherwise weak monotonicity is tested.
Inputs
Type | Name | Default | Description |
Real | x[:] | | Sequence to be tested |
Boolean | strict | false | Set to true to test for strict monotonicity |
Outputs
Type | Name | Description |
Boolean | monotonic | True if x is monotonic increasing or decreasing |
Modelica definition
function isMonotonic
"Returns true if the argument is a monotonic sequence"
input Real x[:] "Sequence to be tested";
input Boolean strict=false "Set to true to test for strict monotonicity";
output Boolean monotonic "True if x is monotonic increasing or decreasing";
protected
Integer n=size(x, 1) "Number of data points";
algorithm
if n == 1 then
monotonic := true;
else
monotonic := true;
if strict then
if (x[1] >= x[n]) then
for i in 1:n - 1 loop
if (not x[i] > x[i + 1]) then
monotonic := false;
end if;
end for;
else
for i in 1:n - 1 loop
if (not x[i] < x[i + 1]) then
monotonic := false;
end if;
end for;
end if;
else
// not strict
if (x[1] >= x[n]) then
for i in 1:n - 1 loop
if (not x[i] >= x[i + 1]) then
monotonic := false;
end if;
end for;
else
for i in 1:n - 1 loop
if (not x[i] <= x[i + 1]) then
monotonic := false;
end if;
end for;
end if;
end if;
// strict
end if;
end isMonotonic;
Polynomial function
Information
This function computes a polynomial of arbitrary order.
The polynomial has the form
y = a1 + a2 x + a3 x2 + ...
Inputs
Type | Name | Default | Description |
Real | x | | Independent variable |
Real | a[:] | | Coefficients |
Outputs
Type | Name | Description |
Real | y | Result |
Modelica definition
function polynomial
"Polynomial function"
annotation(derivative=
Buildings.Utilities.Math.Functions.BaseClasses.der_polynomial);
input Real x
"Independent variable";
input Real a[:]
"Coefficients";
output Real y
"Result";
protected
parameter Integer n =
size(a, 1)-1;
Real xp[n+1]
"Powers of x";
algorithm
xp[1] :=1;
for i
in 1:n
loop
xp[i+1] :=xp[i]*x;
end for;
y :=a*xp;
end polynomial;
Power function that is linearized below a user-defined threshold
Information
Function that approximates y=xn
where 0 < n so that
- the function is defined and monotone increasing for all x.
- dy/dx is bounded and continuous everywhere (for n < 1).
For x < x0, this function replaces
y=xn by
a linear function that is continuously differentiable everywhere.
A typical use of this function is to replace
T = T4(1/4) in a radiation balance to ensure that the
function is defined everywhere. This can help solving the initialization problem
when a solver may be far from a solution and hence T4 < 0.
See the package Examples
for the graph.
Inputs
Type | Name | Default | Description |
Real | x | | Abscissa value |
Real | n | | Exponent |
Real | x0 | | Abscissa value below which linearization occurs |
Outputs
Type | Name | Description |
Real | y | Function value |
Modelica definition
function powerLinearized
"Power function that is linearized below a user-defined threshold"
input Real x "Abscissa value";
input Real n "Exponent";
input Real x0 "Abscissa value below which linearization occurs";
output Real y "Function value";
algorithm
if x > x0 then
y := x^n;
else
y := x0^n * (1-n) + n * x0^(n-1) * x;
end if;
end powerLinearized;
Function that is quadratic in first argument and linear in second argument
Information
This function computes
y = a1 + a2 x1
+ a3 x12
+ (a4 + a5 x1
+ a6 x12) x2
Inputs
Type | Name | Default | Description |
Real | a[6] | | Coefficients |
Real | x1 | | Independent variable for quadratic part |
Real | x2 | | Independent variable for linear part |
Outputs
Type | Name | Description |
Real | y | Result |
Modelica definition
function quadraticLinear
"Function that is quadratic in first argument and linear in second argument"
input Real a[6] "Coefficients";
input Real x1 "Independent variable for quadratic part";
input Real x2 "Independent variable for linear part";
output Real y "Result";
protected
Real x1Sq;
algorithm
x1Sq :=x1*x1;
y :=a[1] + a[2]*x1 + a[3]*x1Sq + (a[4] + a[5]*x1 + a[6]*x1Sq)*x2;
end quadraticLinear;
Power function, regularized near zero, but nonzero value for x=0
Information
Function that approximates y=|x|n where n > 0
so that
- y(0) is not equal to zero.
- dy/dx is bounded and continuous everywhere.
This function replaces y=|x|n in the interval
-δ...+δ by a 4-th order polynomial that has the same
function value and the first and second derivative at x=± δ.
A typical use of this function is to replace the
function for the convective heat transfer
coefficient for forced or free convection that is of the form
h=c |dT|n for some constant c and exponent
0 ≤ n ≤ 1.
By using this function, the original function
that has an infinite derivative near zero and that takes on zero
at the origin is replaced by a function with a bounded derivative and
a non-zero value at the origin. Physically,
the region -δ...+δ may be interpreted as the region
where heat conduction dominates convection in the boundary layer.
See the package Examples
for the graph.
Inputs
Type | Name | Default | Description |
Real | x | | Abscissa value |
Real | n | | Exponent |
Real | delta | 0.01 | Abscissa value where transition occurs |
Outputs
Type | Name | Description |
Real | y | Function value |
Modelica definition
function regNonZeroPower
"Power function, regularized near zero, but nonzero value for x=0"
annotation(derivative=
BaseClasses.der_regNonZeroPower);
input Real x
"Abscissa value";
input Real n
"Exponent";
input Real delta = 0.01
"Abscissa value where transition occurs";
output Real y
"Function value";
protected
Real a1;
Real a3;
Real a5;
Real delta2;
Real x2;
Real y_d
"=y(delta)";
Real yP_d
"=dy(delta)/dx";
Real yPP_d
"=d^2y(delta)/dx^2";
algorithm
if abs(x) > delta
then
y :=
abs(x)^n;
else
delta2 :=delta*delta;
x2 :=x*x;
y_d :=delta^n;
yP_d :=n*delta^(n - 1);
yPP_d :=n*(n - 1)*delta^(n - 2);
a1 := -(yP_d/delta - yPP_d)/delta2/8;
a3 := (yPP_d - 12 * a1 * delta2)/2;
a5 := (y_d - delta2 * (a3 + delta2 * a1));
y := a5 + x2 * (a3 + x2 * a1);
assert(a5>0, "Delta is too small for this exponent.");
end if;
end regNonZeroPower;
Once continuously differentiable approximation to exp(-|x|) in interval |x| < delta
Information
Function to provide a once continuously differentiable approximation
to exp(- |x| )
in the interval |x| < δ for some positive δ
Inputs
Type | Name | Default | Description |
Real | x | | Input argument |
Real | delta | | Transition point where approximation occurs |
Outputs
Type | Name | Description |
Real | y | Output argument |
Modelica definition
function smoothExponential
"Once continuously differentiable approximation to exp(-|x|) in interval |x| < delta"
input Real x "Input argument";
input Real delta "Transition point where approximation occurs";
output Real y "Output argument";
protected
Real absX "Absolute value of x";
Real a2 "Coefficient for approximating function";
Real a3 "Coefficient for approximating function";
Real e;
Real d2;
Real x2;
algorithm
absX :=abs(x);
if absX > delta then
y := exp(-absX);
else
d2 := delta*delta;
e := exp(-delta);
a2 := (delta*e-4*(1-e))/2/d2;
a3 := (e-1-a2*d2)/d2/d2;
x2 := x*x;
y := 1+x2*(a2+x2*a3);
end if;
end smoothExponential;
Once continuously differentiable approximation to the Heaviside function
Information
Once Lipschitz continuously differentiable approximation to the
Heaviside(.,.)
function.
Inputs
Type | Name | Default | Description |
Real | x | | Argument |
Real | delta | | Parameter used for scaling |
Outputs
Type | Name | Description |
Real | y | Result |
Modelica definition
function smoothHeaviside
"Once continuously differentiable approximation to the Heaviside function"
input Real x "Argument";
input Real delta "Parameter used for scaling";
output Real y "Result";
algorithm
y := Buildings.Utilities.Math.Functions.spliceFunction(1, 0, x, delta);
end smoothHeaviside;
Once continuously differentiable approximation to the limit function
Information
Once continuously differentiable approximation to the limit(.,.)
function.
The output is bounded to be in [l, u].
Note that the limit need not be respected, such as illustrated in
Buildings.Utilities.Math.Examples.SmoothMin.
Inputs
Type | Name | Default | Description |
Real | x | | Variable |
Real | l | | Low limit |
Real | u | | Upper limit |
Real | deltaX | | Width of transition interval |
Outputs
Type | Name | Description |
Real | y | Result |
Modelica definition
function smoothLimit
"Once continuously differentiable approximation to the limit function"
input Real x "Variable";
input Real l "Low limit";
input Real u "Upper limit";
input Real deltaX "Width of transition interval";
output Real y "Result";
protected
Real cor;
algorithm
cor :=deltaX/10;
y := Buildings.Utilities.Math.Functions.smoothMax(x,l+deltaX,cor);
y := Buildings.Utilities.Math.Functions.smoothMin(y,u-deltaX,cor);
end smoothLimit;
Once continuously differentiable approximation to the maximum function
Information
Once continuously differentiable approximation to the max(.,.)
function.
Note that the maximum need not be respected, such as illustrated in
Buildings.Utilities.Math.Examples.SmoothMin.
Inputs
Type | Name | Default | Description |
Real | x1 | | First argument |
Real | x2 | | Second argument |
Real | deltaX | | Width of transition interval |
Outputs
Type | Name | Description |
Real | y | Result |
Modelica definition
function smoothMax
"Once continuously differentiable approximation to the maximum function"
input Real x1 "First argument";
input Real x2 "Second argument";
input Real deltaX "Width of transition interval";
output Real y "Result";
algorithm
y := Buildings.Utilities.Math.Functions.spliceFunction(
pos=x1, neg=x2, x=x1-x2, deltax=deltaX);
end smoothMax;
Once continuously differentiable approximation to the minimum function
Information
Once continuously differentiable approximation to the min(.,.)
function.
Note that the minimum need not be respected, such as illustrated in
Buildings.Utilities.Math.Examples.SmoothMin.
Inputs
Type | Name | Default | Description |
Real | x1 | | First argument |
Real | x2 | | Second argument |
Real | deltaX | | Width of transition interval |
Outputs
Type | Name | Description |
Real | y | Result |
Modelica definition
function smoothMin
"Once continuously differentiable approximation to the minimum function"
input Real x1 "First argument";
input Real x2 "Second argument";
input Real deltaX "Width of transition interval";
output Real y "Result";
algorithm
y := Buildings.Utilities.Math.Functions.spliceFunction(
pos=x1, neg=x2, x=x2-x1, deltax=deltaX);
end smoothMin;
Information
Function to provide a once continuously differentialbe transition between
to arguments.
The function is adapted from
Modelica.Media.Air.MoistAir.Utilities.spliceFunction and provided here
for easier accessability to model developers.
Inputs
Type | Name | Default | Description |
Real | pos | | Argument of x > 0 |
Real | neg | | Argument of x < 0 |
Real | x | | Independent value |
Real | deltax | | Half width of transition interval |
Outputs
Type | Name | Description |
Real | out | Smoothed value |
Modelica definition
function spliceFunction
annotation(derivative=
BaseClasses.der_spliceFunction);
input Real pos
"Argument of x > 0";
input Real neg
"Argument of x < 0";
input Real x
"Independent value";
input Real deltax
"Half width of transition interval";
output Real out
"Smoothed value";
protected
Real scaledX1;
Real y;
constant Real asin1 =
Modelica.Math.asin(1);
algorithm
scaledX1 := x/deltax;
if scaledX1 <= -0.999999999
then
out := neg;
elseif scaledX1 >= 0.999999999
then
out := pos;
else
y := (
Modelica.Math.tanh(
Modelica.Math.tan(scaledX1*asin1)) + 1)/2;
out := pos*y + (1 - y)*neg;
end if;
end spliceFunction;
Function to compute the derivatives for cubic hermite spline interpolation
Information
This function computes the derivatives at the support points xi
that can be used as input for evaluating a cubic hermite spline.
If ensureMonotonicity=true
, then the support points yi
need to be monotone increasing (or increasing), and the computed derivatives
di are such that the cubic hermite is monotone increasing (or decreasing).
The algorithm to ensure monotonicity is based on the method described in Fritsch and Carlson (1980) for
ρ = ρ2.
This function is typically used with
Buildings.Utilities.Math.Functions.cubicHermiteLinearExtrapolation
which is used to evaluate the cubic spline.
Because in many applications, the shape of the spline depends on parameters,
this function has been implemented in such a way that all derivatives can be
computed at once and then stored for use during the time stepping,
in which the above function may be called.
References
F.N. Fritsch and R.E. Carlson, Monotone piecewise cubic interpolation.
SIAM J. Numer. Anal., 17 (1980), pp. 238?246.
Inputs
Type | Name | Default | Description |
Real | x[:] | | Support point, strict monotone increasing |
Real | y[size(x, 1)] | | Function values at x |
Boolean | ensureMonotonicity | isMonotonic(y, strict=false) | Set to true to ensure monotonicity of the cubic hermite |
Outputs
Type | Name | Description |
Real | d[size(x, 1)] | Derivative at the support points |
Modelica definition
function splineDerivatives
"Function to compute the derivatives for cubic hermite spline interpolation"
input Real x[:] "Support point, strict monotone increasing";
input Real y[size(x, 1)] "Function values at x";
input Boolean ensureMonotonicity=isMonotonic(y, strict=false)
"Set to true to ensure monotonicity of the cubic hermite";
output Real d[size(x, 1)] "Derivative at the support points";
protected
Integer n=size(x, 1) "Number of data points";
Real delta[n - 1] "Slope of secant line between data points";
Real alpha "Coefficient to ensure monotonicity";
Real beta "Coefficient to ensure monotonicity";
Real tau "Coefficient to ensure monotonicity";
algorithm
if (n>1) then
assert(x[1] < x[n], "x must be strictly increasing.
Received x[1] = " + String(x[1]) + "
x[" + String(n) + "] = " + String(x[n]));
// Check data
assert(isMonotonic(x, strict=true),
"x-values must be strictly monontone increasing or decreasing.");
if ensureMonotonicity then
assert(isMonotonic(y, strict=false),
"If ensureMonotonicity=true, y-values must be monontone increasing or decreasing.");
end if;
end if;
// Compute derivatives at the support points
if n == 1 then
// only one data point
d[1] :=0;
elseif n == 2 then
// linear function
d[1] := (y[2] - y[1])/(x[2] - x[1]);
d[2] := d[1];
else
// Slopes of the secant lines between i and i+1
for i in 1:n - 1 loop
delta[i] := (y[i + 1] - y[i])/(x[i + 1] - x[i]);
end for;
// Initial values for tangents at the support points.
// End points use one-sided derivatives
d[1] := delta[1];
d[n] := delta[n - 1];
for i in 2:n - 1 loop
d[i] := (delta[i - 1] + delta[i])/2;
end for;
end if;
// Ensure monotonicity
if n > 2 and ensureMonotonicity then
for i in 1:n - 1 loop
if (abs(delta[i]) < Modelica.Constants.small) then
d[i] := 0;
d[i + 1] := 0;
else
alpha := d[i]/delta[i];
beta := d[i + 1]/delta[i];
// Constrain derivative to ensure monotonicity in this interval
if (alpha^2 + beta^2) > 9 then
tau := 3/(alpha^2 + beta^2)^(1/2);
d[i] := delta[i]*alpha*tau;
d[i + 1] := delta[i]*beta*tau;
end if;
end if;
end for;
end if;
end splineDerivatives;
Integration using the trapezoidal rule
Information
This function computes a definite integral using the trapezoidal rule.
Inputs
Type | Name | Default | Description |
Integer | N | | Number of integrand points |
Real | f[:] | | Integrands |
Real | deltaX | | Width of interval for Trapezoidal integration |
Outputs
Type | Name | Description |
Real | result | Result |
Modelica definition
function trapezoidalIntegration
"Integration using the trapezoidal rule"
input Integer N "Number of integrand points";
input Real[:] f "Integrands";
input Real deltaX "Width of interval for Trapezoidal integration";
output Real result "Result";
algorithm
assert(N >= 2, "N must be no less than 2.");
result := 0;
for i in 1:N loop
result := result + f[i];
end for;
result := 2*result;
result := result - f[1] - f[N];
result := result*deltaX/2;
end trapezoidalIntegration;
Automatically generated Mon Jul 13 14:30:45 2015.