andes.core package

Submodules

andes.core.block module

class andes.core.block.Block(name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, namespace: str = 'local')[source]

Bases: object

Base class for control blocks.

Blocks are meant to be instantiated as Model attributes to provide pre-defined equation sets. Subclasses must overload the __init__ method to take custom inputs. Subclasses of Block must overload the define method to provide initialization and equation strings. Exported variables, services and blocks must be constructed into a dictionary self.vars at the end of the constructor.

Blocks can be nested. A block can have blocks but itself as attributes and therefore reuse equations. When a block has sub-blocks, the outer block must be constructed with a``name``.

Nested block works in the following way: the parent block modifies the sub-block's name attribute by prepending the parent block's name at the construction phase. The parent block then exports the sub-block as a whole. When the parent Model class picks up the block, it will recursively import the variables in the block and the sub-blocks correctly. See the example section for details.

Parameters:
name : str, optional

Block name

tex_name : str, optional

Block LaTeX name

info : str, optional

Block description.

namespace : str, local or parent

Namespace of the exported elements. If 'local', the block name will be prepended by the parent. If 'parent', the original element name will be used when exporting.

Warning

It is a good practice to avoid more than one level of nesting, to avoid multi-underscore variable names.

Examples

Example for two-level nested blocks. Suppose we have the following hierarchy

SomeModel  instance M
   |
LeadLag A  exports (x, y)
   |
Lag B      exports (x, y)

SomeModel instance M contains an instance of LeadLag block named A, which contains an instance of a Lag block named B. Both A and B exports two variables x and y.

In the code of Model, the following code is used to instantiate LeadLag

class SomeModel:
    def __init__(...)
        ...
        self.A = LeadLag(name='A',
                         u=self.foo1,
                         T1=self.foo2,
                         T2=self.foo3)

To use Lag in the LeadLag code, the following lines are found in the constructor of LeadLag

class LeadLag:
    def __init__(name, ...)
        ...
        self.B = Lag(u=self.y, K=self.K, T=self.T)
        self.vars = {..., 'A': self.A}

The __setattr__ magic of LeadLag takes over the construction and assigns A_B to B.name, given A's name provided at run time. self.A is exported with the internal name A at the end.

Again, the LeadLag instance name (A in this example) MUST be provided in SomeModel's constructor for the name prepending to work correctly. If there is more than one level of nesting, other than the leaf-level block, all parent blocks' names must be provided at instantiation.

When A is picked up by SomeModel.__setattr__, B is captured from A's exports. Recursively, B's variables are exported, Recall that B.name is now A_B, following the naming rule (parent block's name + variable name), B's internal variables become A_B_x and A_B_y.

In this way, B's define() needs no modification since the naming rule is the same. For example, B's internal y is always {self.name}_y, although B has gotten a new name A_B.

class_name

Return the class name.

define()[source]

Function for setting the initialization and equation strings for internal variables. This method must be implemented by subclasses.

The equations should be written with the "final" variable names. Let's say the block instance is named blk (kept at self.name of the block), and an internal variable v is defined. The internal variable will be captured as blk_v by the parent model. Therefore, all equations should use {self.name}_v to represent variable v, where {self.name} is the name of the block at run time.

On the other hand, the names of externally provided parameters or variables are obtained by directly accessing the name attribute. For example, if self.T is a parameter provided through the block constructor, {self.T.name} should be used in the equation.

See also

PIController.define
Equations for the PI Controller block

Examples

An internal variable v has a trivial equation T = v, where T is a parameter provided to the block constructor.

In the model, one has

class SomeModel():
    def __init__(...)
        self.input = Algeb()
        self.T = Param()

        self.blk = ExampleBlock(u=self.input, T=self.T)

In the ExampleBlock function, the internal variable is defined in the constructor as

class ExampleBlock():
    def __init__(...):
        self.v = Algeb()
        self.vars = {'v', self.v}

In the define, the equation is provided as

def define(self):
    self.v.v_str = '{self.T.name}'
    self.v.e_str = '{self.T.name} - {self.name}_v'

In the parent model, v from the block will be captured as blk_v, and the equation will evaluate into

self.blk_v.v_str = 'T'
self.blk_v.e_str = 'T - blk_v'
static enforce_tex_name(fields)[source]

Enforce tex_name is not None

export()[source]

Method for exporting instances defined in this class in a dictionary. This method calls the define method first and returns self.vars.

Returns:
dict

Keys are the (last section of the) variable name, and the values are the attribute instance.

f_numeric(**kwargs)[source]

Function call to update differential equation values.

This function should modify the e value of block State and ExtState in place.

g_numeric(**kwargs)[source]

Function call to update algebraic equation values.

This function should modify the e value of block Algeb and ExtAlgeb in place.

j_numeric()[source]

This function stores the constant and variable jacobian information in corresponding lists.

Constant jacobians are stored by indices and values in, for example, ifxc, jfxc and vfxc. Value scalars or arrays are stored in vfxc.

Variable jacobians are stored by indices and functions. The function shall return the value of the corresponding jacobian elements.

j_reset()[source]

Helper function to clear the lists holding the numerical Jacobians.

This function should be only called once at the beginning of j_numeric in blocks.

class andes.core.block.DeadBand1(u, center, lower, upper, gain=1.0, enable=True, name=None, tex_name=None, info=None, namespace='local')[source]

Bases: andes.core.block.Block

Deadband type 1.

Parameters:
center

Default value when within the deadband. If the input is an error signal, center should be set to zero.

gain

Gain multiplied to DeadBand discrete block's output.

Notes

Block diagram

      |   /
______|__/___   -> Gain -> DeadBand1_y
   /  |
  /   |
define()[source]

Notes

Implemented equation:

\[0 = center + z_u * (u - upper) + z_l * (u - lower) - y\]
class andes.core.block.Gain(u, K, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Gain block.

     ┌───┐
u -> │ K │ -> y
     └───┘

Exports an algebraic output y.

define()[source]

Implemented equation and the initial condition are

\[\begin{split}y = K u \\ y^{(0)} = K u^{(0)}\end{split}\]
class andes.core.block.GainLimiter(u, K, R, lower, upper, no_lower=False, no_upper=False, sign_lower=1, sign_upper=1, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Gain followed by a limiter and another gain.

Exports the limited output y, unlimited output x, and HardLimiter lim.

     ┌─────┐         upper  ┌─────┐
     │     │        /¯¯¯¯¯  │     │
u -> │  K  │ -> x  /   ->   │  R  │ -> y
     │     │ _____/         │     │
     └─────┘ lower          └─────┘
Parameters:
u : str, BaseVar

Input variable, or an equation string for constructing an anonymous variable

K : str, BaseParam, BaseService

Initial gain for u before limiter

R : str, BaseParam, BaseService

Post limiter gain

define()[source]

TODO: write docstring

class andes.core.block.HVGate(u1, u2, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

High Value Gate. Outputs the maximum of two inputs.

      ┌─────────┐
u1 -> │ HV Gate │
      │         │ ->  y
u2 -> │  (MAX)  │
      └─────────┘
define()[source]

Implemented equations and initial conditions

\[0 = s_0^{sl} u_1 + s_1^{sl} u_2 - y y_0 = maximum(u_1, u_2)\]

Notes

In the implementation, one should not use

self.y.v_str = f'maximum({self.u1.name}, {self.u2.name})',

because SymPy processes this equation to {self.u1.name}. Not sure if this is a bug or intended.

class andes.core.block.Integrator(u, T, K, y0, check_init=True, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Integrator block.

     ┌──────┐
u -> │ K/sT │ -> y
     └──────┘

Exports a differential variable y.

The initial output needs to be specified through y0.

define()[source]

Implemented equation and the initial condition are

\[\begin{split}\dot{y} = K u \\ y^{(0)} = 0\end{split}\]
class andes.core.block.IntegratorAntiWindup(u, T, K, y0, lower, upper, name=None, tex_name=None, info=None, no_warn=False)[source]

Bases: andes.core.block.Block

Integrator block with anti-windup limiter.

           upper
          /¯¯¯¯¯
     ┌──────┐
u -> │ K/sT │ -> y
     └──────┘
   _____/
   lower

Exports a differential variable y and an AntiWindup lim. The initial output must be specified through y0.

define()[source]

Implemented equation and the initial condition are

\[\begin{split}\dot{y} = K u \\ y^{(0)} = 0\end{split}\]
class andes.core.block.LVGate(u1, u2, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Low Value Gate. Outputs the minimum of the two inputs.

      ┌─────────┐
u1 -> │ LV Gate |
      │         | ->  y
u2 -> │  (MIN)  |
      └─────────┘
define()[source]

Implemented equations and initial conditions

\[0 = s_0^{sl} u_1 + s_1^{sl} u_2 - y y_0 = minimum(u_1, u_2)\]

Notes

Same problem as HVGate as minimum does not sympify correctly.

class andes.core.block.Lag(u, T, K, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Lag (low pass filter) transfer function.

     ┌────────┐
     │    K   │
u -> │ ────── │ -> y
     │ 1 + sT │
     └────────┘

Exports one state variable y as the output.

Parameters:
K

Gain

T

Time constant

u

Input variable

define()[source]

Notes

Equations and initial values are

\[\begin{split}T \dot{y} &= (Ku - y) \\ y^{(0)} &= K u\end{split}\]
class andes.core.block.Lag2ndOrd(u, K, T1, T2, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Second order lag transfer function (low-pass filter)

     ┌──────────────────┐
     │         K        │
u -> │ ──────────────── │ -> y
     │ 1 + sT1 + s^2 T2 │
     └──────────────────┘

Exports one two state variables (x, y), where y is the output.

Parameters:
u

Input

K

Gain

T1

First order time constant

T2

Second order time constant

define()[source]

Notes

Implemented equations and initial values are

\[\begin{split}T_2 \dot{x} &= Ku - y - T_1 x \\ \dot{y} &= x \\ x^{(0)} &= 0 \\ y^{(0)} &= K u\end{split}\]
class andes.core.block.LagAWFreeze(u, T, K, lower, upper, freeze, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.LagAntiWindup

Lag with anti-windup limiter and state freeze.

The output y is a state variable.

define()[source]

Notes

Equations and initial values are

\[\begin{split}T \dot{y} &= (1 - freeze) (Ku - y) \\ y^{(0)} &= K u\end{split}\]
class andes.core.block.LagAntiWindup(u, T, K, lower, upper, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Lag (low pass filter) transfer function block with an anti-windup limiter.

             upper
           /¯¯¯¯¯¯
     ┌────────┐
     │    K   │
u -> │ ────── │ -> y
     │ 1 + sT │
     └────────┘
   ______/
   lower

Exports one state variable y as the output and one AntiWindup instance lim.

Parameters:
K

Gain

T

Time constant

u

Input variable

define()[source]

Notes

Equations and initial values are

\[\begin{split}T \dot{y} &= (Ku - y) \\ y^{(0)} &= K u\end{split}\]
class andes.core.block.LagAntiWindupRate(u, T, K, lower, upper, rate_lower, rate_upper, no_lower=False, no_upper=False, rate_no_lower=False, rate_no_upper=False, rate_lower_cond=None, rate_upper_cond=None, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Lag (low pass filter) transfer function block with a rate limiter and an anti-windup limiter.

             upper && rate_upper
           /¯¯¯¯¯¯
     ┌────────┐
     │    K   │
u -> │ ────── │ -> y
     │ 1 + sT │
     └────────┘
   ______/
   lower & rate_lower

Exports one state variable y as the output and one AntiWindupRate instance lim.

Parameters:
K

Gain

T

Time constant

u

Input variable

define()[source]

Notes

Equations and initial values are

\[\begin{split}T \dot{y} &= (Ku - y) \\ y^{(0)} &= K u\end{split}\]
class andes.core.block.LagFreeze(u, T, K, freeze, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Lag

Lag with a state freeze input.

define()[source]

Notes

Equations and initial values are

\[\begin{split}T \dot{y} &= (1 - freeze) * (Ku - y) \\ y^{(0)} &= K u\end{split}\]
class andes.core.block.LagRate(u, T, K, rate_lower, rate_upper, rate_no_lower=False, rate_no_upper=False, rate_lower_cond=None, rate_upper_cond=None, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Lag (low pass filter) transfer function block with a rate limiter and an anti-windup limiter.

            rate_upper
           /
     ┌────────┐
     │    K   │
u -> │ ────── │ -> y
     │ 1 + sT │
     └────────┘
         /
rate_lower

Exports one state variable y as the output and one AntiWindupRate instance lim.

Parameters:
K

Gain

T

Time constant

u

Input variable

define()[source]

Notes

Equations and initial values are

\[\begin{split}T \dot{y} &= (Ku - y) \\ y^{(0)} &= K u\end{split}\]
class andes.core.block.LeadLag(u, T1, T2, K=1, zero_out=True, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Lead-Lag transfer function block in series implementation

     ┌───────────┐
     │   1 + sT1 │
u -> │ K ─────── │ -> y
     │   1 + sT2 │
     └───────────┘

Exports two variables: internal state x and output algebraic variable y.

Parameters:
T1 : BaseParam

Time constant 1

T2 : BaseParam

Time constant 2

zero_out : bool

True to allow zeroing out lead-lag as a pass through (when T1=T2=0)

Notes

To allow zeroing out lead-lag as a pure gain, set zero_out to True.

define()[source]

Notes

Implemented equations and initial values

\[\begin{split}T_2 \dot{x'} &= (u - x') \\ T_2 y &= K T_1 (u - x') + K T_2 x' + E_2 \, , \text{where} \\ E_2 = & \left\{\begin{matrix} (y - K x') &\text{ if } T_1 = T_2 = 0 \& zero\_out=True \\ 0& \text{ otherwise } \end{matrix}\right. \\ x'^{(0)} & = u\\ y^{(0)} & = Ku\\\end{split}\]
class andes.core.block.LeadLag2ndOrd(u, T1, T2, T3, T4, zero_out=False, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Second-order lead-lag transfer function block

     ┌──────────────────┐
     │ 1 + sT3 + s^2 T4 │
u -> │ ──────────────── │ -> y
     │ 1 + sT1 + s^2 T2 │
     └──────────────────┘

Exports two internal states (x1 and x2) and output algebraic variable y.

# TODO: instead of implementing zero_out using LessThan and an additional term, consider correcting all parameters to 1 if all are 0.

define()[source]

Notes

Implemented equations and initial values are

\[\begin{split}T_2 \dot{x}_1 &= u - x_2 - T_1 x_1 \\ \dot{x}_2 &= x_1 \\ T_2 y &= T_2 x_2 + T_2 T_3 x_1 + T_4 (u - x_2 - T_1 x_1) + E_2 \, , \text{ where} \\ E_2 = & \left\{\begin{matrix} (y - x_2) &\text{ if } T_1 = T_2 = T_3 = T_4 = 0 \& zero\_out=True \\ 0& \text{ otherwise } \end{matrix}\right. \\ x_1^{(0)} &= 0 \\ x_2^{(0)} &= y^{(0)} = u\end{split}\]
class andes.core.block.LeadLagLimit(u, T1, T2, lower, upper, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Lead-Lag transfer function block with hard limiter (series implementation)

     ┌─────────┐          upper
     │ 1 + sT1 │         /¯¯¯¯¯
u -> │ ─────── │ -> ynl / -> y
     │ 1 + sT2 │  _____/
     └─────────┘  lower

Exports four variables: state x, output before hard limiter ynl, output y, and AntiWindup lim.

define()[source]

Notes

Implemented control block equations (without limiter) and initial values

\[\begin{split}T_2 \dot{x'} &= (u - x') \\ T_2 y &= T_1 (u - x') + T_2 x' \\ x'^{(0)} &= y^{(0)} = u\end{split}\]
class andes.core.block.LimiterGain(u, K, lower, upper, no_lower=False, no_upper=False, sign_lower=1, sign_upper=1, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Limiter followed by a gain.

Exports the limited output y, unlimited output x, and HardLimiter lim.

           upper ┌─────┐
          /¯¯¯¯¯ │     │
u  ->    /   ->  │  K  │ -> y
   _____/        │     │
   lower         └─────┘

Deprecated since version 1.5.0: LimiterGain will be removed in ANDES 1.5.0. it is replaced by GainLimiter because the latter supports pre- and post-gains.

define()[source]

Function for setting the initialization and equation strings for internal variables. This method must be implemented by subclasses.

The equations should be written with the "final" variable names. Let's say the block instance is named blk (kept at self.name of the block), and an internal variable v is defined. The internal variable will be captured as blk_v by the parent model. Therefore, all equations should use {self.name}_v to represent variable v, where {self.name} is the name of the block at run time.

On the other hand, the names of externally provided parameters or variables are obtained by directly accessing the name attribute. For example, if self.T is a parameter provided through the block constructor, {self.T.name} should be used in the equation.

See also

PIController.define
Equations for the PI Controller block

Examples

An internal variable v has a trivial equation T = v, where T is a parameter provided to the block constructor.

In the model, one has

class SomeModel():
    def __init__(...)
        self.input = Algeb()
        self.T = Param()

        self.blk = ExampleBlock(u=self.input, T=self.T)

In the ExampleBlock function, the internal variable is defined in the constructor as

class ExampleBlock():
    def __init__(...):
        self.v = Algeb()
        self.vars = {'v', self.v}

In the define, the equation is provided as

def define(self):
    self.v.v_str = '{self.T.name}'
    self.v.e_str = '{self.T.name} - {self.name}_v'

In the parent model, v from the block will be captured as blk_v, and the equation will evaluate into

self.blk_v.v_str = 'T'
self.blk_v.e_str = 'T - blk_v'
class andes.core.block.PIAWHardLimit(u, kp, ki, aw_lower, aw_upper, lower, upper, no_lower=False, no_upper=False, ref=0.0, x0=0.0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.PIController

PI controller with anti-windup limiter on the integrator and hard limit on the output.

Limits lower and upper are on the final output, and aw_lower aw_upper are on the integrator.

define()[source]

Define equations for the PI Controller.

Notes

One state variable xi and one algebraic variable y are added.

Equations implemented are

\[\begin{split}\dot{x_i} &= k_i * (u - ref) \\ y &= x_i + k_p * (u - ref)\end{split}\]
class andes.core.block.PIController(u, kp, ki, ref=0.0, x0=0.0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Proportional Integral Controller.

The controller takes an error signal as the input. It takes an optional ref signal, which will be subtracted from the input.

Parameters:
u : BaseVar

The input variable instance

kp : BaseParam

The proportional gain parameter instance

ki : [type]

The integral gain parameter instance

define()[source]

Define equations for the PI Controller.

Notes

One state variable xi and one algebraic variable y are added.

Equations implemented are

\[\begin{split}\dot{x_i} &= k_i * (u - ref) \\ y &= x_i + k_p * (u - ref)\end{split}\]
class andes.core.block.PIControllerNumeric(u, kp, ki, ref=0.0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

A PI Controller implemented with numerical function calls.

ref must not be a variable.

define()[source]

Skip the symbolic definition

f_numeric(**kwargs)[source]

Function call to update differential equation values.

This function should modify the e value of block State and ExtState in place.

g_numeric(**kwargs)[source]

Function call to update algebraic equation values.

This function should modify the e value of block Algeb and ExtAlgeb in place.

j_numeric()[source]

This function stores the constant and variable jacobian information in corresponding lists.

Constant jacobians are stored by indices and values in, for example, ifxc, jfxc and vfxc. Value scalars or arrays are stored in vfxc.

Variable jacobians are stored by indices and functions. The function shall return the value of the corresponding jacobian elements.

class andes.core.block.PIFreeze(u, kp, ki, freeze, ref=0.0, x0=0.0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.PIController

PI controller with state freeze.

Freezes state when the corresponding freeze == 1.

Notes

Tested in experimental.TestPITrackAW.PIFreeze.

define()[source]

Notes

One state variable xi and one algebraic variable y are added.

Equations implemented are

\[\begin{split}\dot{x_i} &= k_i * (u - ref) \\ y &= (1-freeze) * (x_i + k_p * (u - ref)) + freeze * y\end{split}\]
class andes.core.block.PITrackAW(u, kp, ki, ks, lower, upper, no_lower=False, no_upper=False, ref=0.0, x0=0.0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

PI with tracking anti-windup limiter

define()[source]

Function for setting the initialization and equation strings for internal variables. This method must be implemented by subclasses.

The equations should be written with the "final" variable names. Let's say the block instance is named blk (kept at self.name of the block), and an internal variable v is defined. The internal variable will be captured as blk_v by the parent model. Therefore, all equations should use {self.name}_v to represent variable v, where {self.name} is the name of the block at run time.

On the other hand, the names of externally provided parameters or variables are obtained by directly accessing the name attribute. For example, if self.T is a parameter provided through the block constructor, {self.T.name} should be used in the equation.

See also

PIController.define
Equations for the PI Controller block

Examples

An internal variable v has a trivial equation T = v, where T is a parameter provided to the block constructor.

In the model, one has

class SomeModel():
    def __init__(...)
        self.input = Algeb()
        self.T = Param()

        self.blk = ExampleBlock(u=self.input, T=self.T)

In the ExampleBlock function, the internal variable is defined in the constructor as

class ExampleBlock():
    def __init__(...):
        self.v = Algeb()
        self.vars = {'v', self.v}

In the define, the equation is provided as

def define(self):
    self.v.v_str = '{self.T.name}'
    self.v.e_str = '{self.T.name} - {self.name}_v'

In the parent model, v from the block will be captured as blk_v, and the equation will evaluate into

self.blk_v.v_str = 'T'
self.blk_v.e_str = 'T - blk_v'
class andes.core.block.PITrackAWFreeze(u, kp, ki, ks, lower, upper, freeze, no_lower=False, no_upper=False, ref=0.0, x0=0.0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.PITrackAW

PI controller with tracking anti-windup limiter and state freeze.

define()[source]

Function for setting the initialization and equation strings for internal variables. This method must be implemented by subclasses.

The equations should be written with the "final" variable names. Let's say the block instance is named blk (kept at self.name of the block), and an internal variable v is defined. The internal variable will be captured as blk_v by the parent model. Therefore, all equations should use {self.name}_v to represent variable v, where {self.name} is the name of the block at run time.

On the other hand, the names of externally provided parameters or variables are obtained by directly accessing the name attribute. For example, if self.T is a parameter provided through the block constructor, {self.T.name} should be used in the equation.

See also

PIController.define
Equations for the PI Controller block

Examples

An internal variable v has a trivial equation T = v, where T is a parameter provided to the block constructor.

In the model, one has

class SomeModel():
    def __init__(...)
        self.input = Algeb()
        self.T = Param()

        self.blk = ExampleBlock(u=self.input, T=self.T)

In the ExampleBlock function, the internal variable is defined in the constructor as

class ExampleBlock():
    def __init__(...):
        self.v = Algeb()
        self.vars = {'v', self.v}

In the define, the equation is provided as

def define(self):
    self.v.v_str = '{self.T.name}'
    self.v.e_str = '{self.T.name} - {self.name}_v'

In the parent model, v from the block will be captured as blk_v, and the equation will evaluate into

self.blk_v.v_str = 'T'
self.blk_v.e_str = 'T - blk_v'
class andes.core.block.Piecewise(u, points: Union[List[T], Tuple], funs: Union[List[T], Tuple], name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Piecewise block. Outputs an algebraic variable y.

This block takes a list of N points, [x0, x1, ...x_{n-1}] to define N+1 ranges, namely (-inf, x0), (x0, x1), ..., (x_{n-1}, +inf). and a list of N+1 function strings [fun0, ..., fun_n].

Inputs that fall within each range applies the corresponding function. The first range (-inf, x0) applies fun_0, and the last range (x_{n-1}, +inf) applies the last function fun_n.

Parameters:
points : list, tuple

A list of piecewise points. Need to be provided in the constructor function.

funs : list, tuple

A list of strings for the piecewise functions. Need to be provided in the overloaded define function.

define()[source]

Build the equation string for the piecewise equations.

self.funs needs to be provided with the function strings corresponding to each range.

class andes.core.block.Washout(u, T, K, name=None, tex_name=None, info=None)[source]

Bases: andes.core.block.Block

Washout filter (high pass) block.

     ┌────────┐
     │   sK   │
u -> │ ────── │ -> y
     │ 1 + sT │
     └────────┘

Exports state x (symbol x') and output algebraic variable y.

define()[source]

Notes

Equations and initial values:

\[\begin{split}T \dot{x'} &= (u - x') \\ T y &= K (u - x') \\ x'^{(0)} &= u \\ y^{(0)} &= 0\end{split}\]
class andes.core.block.WashoutOrLag(u, T, K, name=None, zero_out=True, tex_name=None, info=None)[source]

Bases: andes.core.block.Washout

Washout with the capability to convert to Lag when K = 0.

Can be enabled with zero_out. Need to provide name to construct.

Exports state x (symbol x'), output algebraic variable y, and a LessThan block LT.

Parameters:
zero_out : bool, optional

If True, sT will become 1, and the washout will become a low-pass filter. If False, functions as a regular Washout.

define()[source]

Notes

Equations and initial values:

\[\begin{split}T \dot{x'} &= (u - x') \\ T y = z_0 K (u - x') + z_1 T x \\ x'^{(0)} &= u \\ y^{(0)} &= 0\end{split}\]

where z_0 is a flag array for the greater-than-zero elements, and z_1 is that for the less-than or equal-to zero elements.

andes.core.discrete module

class andes.core.discrete.AntiWindup(u, lower, upper, enable=True, no_warn=False, no_lower=False, no_upper=False, sign_lower=1, sign_upper=1, name=None, tex_name=None, info=None, state=None)[source]

Bases: andes.core.discrete.Limiter

Anti-windup limiter.

Anti-windup limiter prevents the wind-up effect of a differential variable. The derivative of the differential variable is reset if it continues to increase in the same direction after exceeding the limits. During the derivative return, the limiter will be inactive

if x > xmax and x dot > 0: x = xmax and x dot = 0
if x < xmin and x dot < 0: x = xmin and x dot = 0

This class takes one more optional parameter for specifying the equation.

Parameters:
state : State, ExtState

A State (or ExtState) whose equation value will be checked and, when condition satisfies, will be reset by the anti-windup-limiter.

check_eq()[source]

Check the variables and equations and set the limiter flags. Reset differential equation values based on limiter flags.

Notes

The current implementation reallocates memory for self.x_set in each call. Consider improving for speed. (TODO)

check_var(*args, **kwargs)[source]

This function is empty. Defers check_var to check_eq.

class andes.core.discrete.AntiWindupRate(u, lower, upper, rate_lower, rate_upper, no_lower=False, no_upper=False, rate_no_lower=False, rate_no_upper=False, rate_lower_cond=None, rate_upper_cond=None, enable=True, name=None, tex_name=None, info=None)[source]

Bases: andes.core.discrete.AntiWindup, andes.core.discrete.RateLimiter

Anti-windup limiter with rate limits

check_eq()[source]

Check the variables and equations and set the limiter flags. Reset differential equation values based on limiter flags.

Notes

The current implementation reallocates memory for self.x_set in each call. Consider improving for speed. (TODO)

class andes.core.discrete.Average(u, mode='step', delay=0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.discrete.Delay

Compute the average of a BaseVar over a period of time or a number of samples.

check_var(dae_t, *args, **kwargs)[source]

This function is called in l_update_var before evaluating equations.

It should update internal flags only.

class andes.core.discrete.DeadBand(u, center, lower, upper, enable=True, equal=False, zu=0.0, zl=0.0, zi=0.0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.discrete.Limiter

The basic deadband type.

Parameters:
u : NumParam

The pre-deadband input variable

center : NumParam

Neutral value of the output

lower : NumParam

Lower bound

upper : NumParam

Upper bound

enable : bool

Enabled if True; Disabled and works as a pass-through if False.

Notes

Input changes within a deadband will incur no output changes. This component computes and exports three flags.

Three flags computed from the current input:
  • zl: True if the input is below the lower threshold
  • zi: True if the input is within the deadband
  • zu: True if is above the lower threshold

Initial condition:

All three flags are initialized to zero. All flags are updated during check_var when enabled. If the deadband component is not enabled, all of them will remain zero.

Examples

Exported deadband flags need to be used in the algebraic equation corresponding to the post-deadband variable. Assume the pre-deadband input variable is var_in and the post-deadband variable is var_out. First, define a deadband instance db in the model using

self.db = DeadBand(u=self.var_in, center=self.dbc,
                   lower=self.dbl, upper=self.dbu)

To implement a no-memory deadband whose output returns to center when the input is within the band, the equation for var can be written as

var_out.e_str = 'var_in * (1 - db_zi) + \
                 (dbc * db_zi) - var_out'
check_var(*args, **kwargs)[source]

Notes

Updates three flags: zi, zu, zl based on the following rules:

zu:
1 if u > upper; 0 otherwise.
zl:
1 if u < lower; 0 otherwise.
zi:
not(zu or zl);
class andes.core.discrete.DeadBandRT(u, center, lower, upper, enable=True)[source]

Bases: andes.core.discrete.DeadBand

Deadband with flags for directions of return.

Parameters:
u : NumParam

The pre-deadband input variable

center : NumParam

Neutral value of the output

lower : NumParam

Lower bound

upper : NumParam

Upper bound

enable : bool

Enabled if True; Disabled and works as a pass-through if False.

Notes

Input changes within a deadband will incur no output changes. This component computes and exports five flags. The additional two flags on top of DeadBand indicate the direction of return:

  • zur: True if the input is/has been within the deadband and was returned from the upper threshold
  • zlr: True if the input is/has been within the deadband and was returned from the lower threshold

Initial condition:

All five flags are initialized to zero. All flags are updated during check_var when enabled. If the deadband component is not enabled, all of them will remain zero.

Examples

To implement a deadband whose output is pegged at the nearest deadband bounds, the equation for var can be provided as

var_out.e_str = 'var_in * (1 - db_zi) + \
                 dbl * db_zlr + \
                 dbu * db_zur - var_out'
check_var(*args, **kwargs)[source]

Notes

Updates five flags: zi, zu, zl; zur, and zlr based on the following rules:

zu:
1 if u > upper; 0 otherwise.
zl:
1 if u < lower; 0 otherwise.
zi:
not(zu or zl);
zur:
  • set to 1 when (previous zu + present zi == 2)
  • hold when (previous zi == zi)
  • clear otherwise
zlr:
  • set to 1 when (previous zl + present zi == 2)
  • hold when (previous zi == zi)
  • clear otherwise
class andes.core.discrete.Delay(u, mode='step', delay=0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.discrete.Discrete

Delay class to memorize past variable values.

Delay allows to impose a predefined "delay" (in either steps or seconds) for an input variable. The amount of delay is a scalar and has to be fixed at model definition in the current implementation.

check_var(dae_t, *args, **kwargs)[source]

This function is called in l_update_var before evaluating equations.

It should update internal flags only.

list2array(n)[source]

Allocate memory for storage arrays.

class andes.core.discrete.Derivative(u, name=None, tex_name=None, info=None)[source]

Bases: andes.core.discrete.Delay

Compute the derivative of an algebraic variable using numerical differentiation.

check_var(dae_t, *args, **kwargs)[source]

This function is called in l_update_var before evaluating equations.

It should update internal flags only.

class andes.core.discrete.Discrete(name=None, tex_name=None, info=None, no_warn=False, min_iter=2, err_tol=0.01)[source]

Bases: object

Base discrete class.

Discrete classes export flag arrays (usually boolean) .

check_eq()[source]

This function is called in l_check_eq after updating equations.

It updates internal flags, set differential equations, and record pegged variables.

check_iter_err(niter=None, err=None)[source]

Check if the minimum iteration or maximum error is reached so that this discrete block should be enabled.

Only when both niter and err are given, (niter < min_iter) , and (err > err_tol) it will return False.

This logic will start checking the discrete states if called from an external solver that does not feed niter or err at each step.

Returns:
bool

True if it should be enabled, False otherwise

check_var(*args, **kwargs)[source]

This function is called in l_update_var before evaluating equations.

It should update internal flags only.

class_name
get_names()[source]

Available symbols from this class

get_tex_names()[source]

Return tex_names of exported flags.

TODO: Fix the bug described in the warning below.

Returns:
list

A list of tex_names for all exported flags.

Warning

If underscore _ appears in both flag tex_name and self.tex_name (for example, when this discrete is within a block), the exported tex_name will become invalid for SymPy. Variable name substitution will fail.

get_values()[source]
list2array(n)[source]
warn_init_limit()[source]

Warn if initialized at limits.

class andes.core.discrete.HardLimiter(u, lower, upper, enable=True, name=None, tex_name=None, info=None, min_iter: int = 2, err_tol: float = 0.01, no_lower=False, no_upper=False, sign_lower=1, sign_upper=1, equal=True, no_warn=False, zu=0.0, zl=0.0, zi=1.0)[source]

Bases: andes.core.discrete.Limiter

Hard limiter for algebraic or differential variable. This class is an alias of Limiter.

class andes.core.discrete.LessThan(u, bound, equal=False, enable=True, name=None, tex_name=None, info=None, cache=False, z0=0, z1=1)[source]

Bases: andes.core.discrete.Discrete

Less than (<) comparison function.

Exports two flags: z1 and z0. For elements satisfying the less-than condition, the corresponding z1 = 1. z0 is the element-wise negation of z1.

Notes

The default z0 and z1, if not enabled, can be set through the constructor.

check_var(*args, **kwargs)[source]

If enabled, set flags based on inputs. Use cached values if enabled.

class andes.core.discrete.Limiter(u, lower, upper, enable=True, name=None, tex_name=None, info=None, min_iter: int = 2, err_tol: float = 0.01, no_lower=False, no_upper=False, sign_lower=1, sign_upper=1, equal=True, no_warn=False, zu=0.0, zl=0.0, zi=1.0)[source]

Bases: andes.core.discrete.Discrete

Base limiter class.

This class compares values and sets limit values. Exported flags are zi, zl and zu.

Parameters:
u : BaseVar

Input Variable instance

lower : BaseParam

Parameter instance for the lower limit

upper : BaseParam

Parameter instance for the upper limit

no_lower : bool

True to only use the upper limit

no_upper : bool

True to only use the lower limit

sign_lower: 1 or -1

Sign to be multiplied to the lower limit

sign_upper: bool

Sign to be multiplied to the upper limit

equal : bool

True to include equal signs in comparison (>= or <=).

no_warn : bool

Disable initial limit warnings

zu : 0 or 1

Default value for zu if not enabled

zl : 0 or 1

Default value for zl if not enabled

zi : 0 or 1

Default value for zi if not enabled

Notes

If not enabled, the default flags are zu = zl = 0, zi = 1.

Attributes:
zl : array-like

Flags of elements violating the lower limit; A array of zeros and/or ones.

zi : array-like

Flags for within the limits

zu : array-like

Flags for violating the upper limit

check_var(*args, **kwargs)[source]

Evaluate the flags.

class andes.core.discrete.RateLimiter(u, lower, upper, enable=True, no_lower=False, no_upper=False, lower_cond=None, upper_cond=None, name=None, tex_name=None, info=None)[source]

Bases: andes.core.discrete.Discrete

Rate limiter for a differential variable.

RateLimiter does not export any variable. It directly modifies the differential equation value.

Warning

RateLimiter cannot be applied to a state variable that already undergoes an AntiWindup limiter. Use AntiWindupRate for a rate-limited anti-windup limiter.

Notes

RateLimiter inherits from Discrete to avoid internal naming conflicts with Limiter.

check_eq()[source]

This function is called in l_check_eq after updating equations.

It updates internal flags, set differential equations, and record pegged variables.

class andes.core.discrete.Sampling(u, interval=1.0, offset=0.0, name=None, tex_name=None, info=None)[source]

Bases: andes.core.discrete.Discrete

Sample an input variable repeatedly at a given time interval.

check_var(dae_t, *args, **kwargs)[source]

Check and update the output.

Notes

Present output stored in v. Output of the last step is stored in _last_v. Time for the last output is stored in _last_t.

Initially, store v and _last_v.

If time progresses and dae_t is a multiple of period, update _last_v and then v. Record _last_t.

If time does not progress, update v.

If time rewinds, restore _last_v to v.

list2array(n)[source]
class andes.core.discrete.Selector(*args, fun, tex_name=None, info=None)[source]

Bases: andes.core.discrete.Discrete

Selection between two variables using the provided reduce function.

The reduce function should take the given number of arguments. An example function is np.maximum.reduce which can be used to select the maximum.

Names are in s0, s1.

Warning

A potential bug when more than two inputs are provided, and values in different inputs are equal. Only two inputs are allowed.

Notes

A common pitfall is the 0-based indexing in the Selector flags. Note that exported flags start from 0. Namely, s0 corresponds to the first variable provided for the Selector constructor.

Examples

Example 1: select the largest value between v0 and v1 and put it into vmax.

After the definitions of v0 and v1, define the algebraic variable vmax for the largest value, and a selector vs

self.vmax = Algeb(v_str='maximum(v0, v1)',
                  tex_name='v_{max}',
                  e_str='vs_s0 * v0 + vs_s1 * v1 - vmax')

self.vs = Selector(self.v0, self.v1, fun=np.maximum.reduce)

The initial value of vmax is calculated by maximum(v0, v1), which is the element-wise maximum in SymPy and will be generated into np.maximum(v0, v1). The equation of vmax is to select the values based on vs_s0 and vs_s1.

check_var(*args, **kwargs)[source]

Set the i-th variable's flags to 1 if the return of the reduce function equals the i-th input.

class andes.core.discrete.ShuntAdjust(*, v, lower, upper, bsw, gsw, dt, u, enable=True, min_iter=2, err_tol=0.01, name=None, tex_name=None, info=None, no_warn=False)[source]

Bases: andes.core.discrete.Discrete

Class for adjusting switchable shunts.

Parameters:
v : BaseVar

Voltage measurement

lower : BaseParam

Lower voltage bound

upper : BaseParam

Upper voltage bound

bsw : SwBlock

SwBlock instance for susceptance

gsw : SwBlock

SwBlock instance for conductance

dt : NumParam

Delay time

u : NumParam

Connection status

min_iter : int

Minimum iteration number to enable shunt switching

err_tol : float

Minimum iteration tolerance to enable switching

check_var(dae_t, *args, niter=None, err=None, **kwargs)[source]

Check voltage and perform shunt switching.

Parameters:
niter : int or None

Current iteration step

class andes.core.discrete.SortedLimiter(u, lower, upper, n_select: int = 5, name=None, tex_name=None, enable=True, abs_violation=True, min_iter: int = 2, err_tol: float = 0.01, zu=0.0, zl=0.0, zi=1.0, ql=0.0, qu=0.0)[source]

Bases: andes.core.discrete.Limiter

A limiter that sorts inputs based on the absolute or relative amount of limit violations.

Parameters:
n_select : int

the number of violations to be flagged, for each of over-limit and under-limit cases. If n_select == 1, at most one over-limit and one under-limit inputs will be flagged. If n_select is zero, heuristics will be used.

abs_violation : bool

True to use the absolute violation. False if the relative violation abs(violation/limit) is used for sorting. Since most variables are in per unit, absolute violation is recommended.

calc_select()[source]

Set n_select automatically.

check_var(*args, niter=None, err=None, **kwargs)[source]

Check for the largest and smallest n_select elements.

list2array(n)[source]

Initialize maximum and minimum n_select based on input size.

class andes.core.discrete.Switcher(u, options: Union[list, Tuple], info: str = None, name: str = None, tex_name: str = None, cache=True)[source]

Bases: andes.core.discrete.Discrete

Switcher based on an input parameter.

The switch class takes one v-provider, compares the input with each value in the option list, and exports one flag array for each option. The flags are 0-indexed.

Exported flags are named with _s0, _s1, ..., with a total number of len(options). See the examples section.

Notes

Switches needs to be distinguished from Selector.

Switcher is for generating flags indicating option selection based on an input parameter. Selector is for generating flags at run time based on variable values and a selection function.

Examples

The IEEEST model takes an input for selecting the signal. Options are 1 through 6. One can construct

self.IC = NumParam(info='input code 1-6')  # input code
self.SW = Switcher(u=self.IC, options=[0, 1, 2, 3, 4, 5, 6])

If the IC values from the data file ends up being

self.IC.v = np.array([1, 2, 2, 4, 6])

Then, the exported flag arrays will be

{'IC_s0': np.array([0, 0, 0, 0, 0]),
 'IC_s1': np.array([1, 0, 0, 0, 0]),
 'IC_s2': np.array([0, 1, 1, 0, 0]),
 'IC_s3': np.array([0, 0, 0, 0, 0]),
 'IC_s4': np.array([0, 0, 0, 1, 0]),
 'IC_s5': np.array([0, 0, 0, 0, 0]),
 'IC_s6': np.array([0, 0, 0, 0, 1])
}

where IC_s0 is used for padding so that following flags align with the options.

check_var(*args, **kwargs)[source]

Set the switcher flags based on inputs. Uses cached flags if cache is set to True.

list2array(n)[source]

This forces to evaluate Switcher upon System setup

andes.core.model module

Base class for building ANDES models.

class andes.core.model.Model(system=None, config=None)[source]

Bases: object

Base class for power system DAE models.

After subclassing ModelData, subclass Model` to complete a DAE model. Subclasses of Model defines DAE variables, services, and other types of parameters, in the constructor __init__.

Notes

To modify parameters or services use set(), which writes directly to the given attribute, or alter(), which converts parameters to system base like that for input data.

Examples

Take the static PQ as an example, the subclass of Model, PQ, should looks like

class PQ(PQData, Model):
    def __init__(self, system, config):
        PQData.__init__(self)
        Model.__init__(self, system, config)

Since PQ is calling the base class constructors, it is meant to be the final class and not further derived. It inherits from PQData and Model and must call constructors in the order of PQData and Model. If the derived class of Model needs to be further derived, it should only derive from Model and use a name ending with Base. See andes.models.synchronous.GENBASE.

Next, in PQ.__init__, set proper flags to indicate the routines in which the model will be used

self.flags.update({'pflow': True})

Currently, flags pflow and tds are supported. Both are False by default, meaning the model is neither used in power flow nor time-domain simulation. A very common pitfall is forgetting to set the flag.

Next, the group name can be provided. A group is a collection of models with common parameters and variables. Devices idx of all models in the same group must be unique. To provide a group name, use

self.group = 'StaticLoad'

The group name must be an existing class name in andes.models.group. The model will be added to the specified group and subject to the variable and parameter policy of the group. If not provided with a group class name, the model will be placed in the Undefined group.

Next, additional configuration flags can be added. Configuration flags for models are load-time variables specifying the behavior of a model. It can be exported to an andes.rc file and automatically loaded when creating the System. Configuration flags can be used in equation strings, as long as they are numerical values. To add config flags, use

self.config.add(OrderedDict((('pq2z', 1), )))

It is recommended to use OrderedDict instead of dict, although the syntax is verbose. Note that booleans should be provided as integers (1, or 0), since True or False is interpreted as a string when loaded from the rc file and will cause an error.

Next, it's time for variables and equations! The PQ class does not have internal variables itself. It uses its bus parameter to fetch the corresponding a and v variables of buses. Equation wise, it imposes an active power and a reactive power load equation.

To define external variables from Bus, use

self.a = ExtAlgeb(model='Bus', src='a',
                  indexer=self.bus, tex_name=r'\theta')
self.v = ExtAlgeb(model='Bus', src='v',
                  indexer=self.bus, tex_name=r'V')

Refer to the subsection Variables for more details.

The simplest PQ model will impose constant P and Q, coded as

self.a.e_str = "u * p"
self.v.e_str = "u * q"

where the e_str attribute is the equation string attribute. u is the connectivity status. Any parameter, config, service or variables can be used in equation strings.

Three additional scalars can be used in equations: - dae_t for the current simulation time can be used if the model has flag tds. - sys_f for system frequency (from system.config.freq). - sys_mva for system base mva (from system.config.mva).

The above example is overly simplified. Our PQ model wants a feature to switch itself to a constant impedance if the voltage is out of the range (vmin, vmax). To implement this, we need to introduce a discrete component called Limiter, which yields three arrays of binary flags, zi, zl, and zu indicating in range, below lower limit, and above upper limit, respectively.

First, create an attribute vcmp as a Limiter instance

self.vcmp = Limiter(u=self.v, lower=self.vmin, upper=self.vmax,
                     enable=self.config.pq2z)

where self.config.pq2z is a flag to turn this feature on or off. After this line, we can use vcmp_zi, vcmp_zl, and vcmp_zu in other equation strings.

self.a.e_str = "u * (p0 * vcmp_zi + " \
               "p0 * vcmp_zl * (v ** 2 / vmin ** 2) + " \
               "p0 * vcmp_zu * (v ** 2 / vmax ** 2))"

self.v.e_str = "u * (q0 * vcmp_zi + " \
               "q0 * vcmp_zl * (v ** 2 / vmin ** 2) + "\
               "q0 * vcmp_zu * (v ** 2 / vmax ** 2))"

Note that PQ.a.e_str can use the three variables from vcmp even before defining PQ.vcmp, as long as PQ.vcmp is defined, because vcmp_zi is just a string literal in e_str.

The two equations above implements a piecewise power injection equation. It selects the original power demand if within range, and uses the calculated power when out of range.

Finally, to let ANDES pick up the model, the model name needs to be added to models/__init__.py. Follow the examples in the OrderedDict, where the key is the file name, and the value is the class name.

Attributes:
num_params : OrderedDict

{name: instance} of numerical parameters, including internal and external ones

a_reset()[source]

Reset addresses to empty and reset flags.address to False.

alter(src, idx, value)[source]

Alter input parameter or service values.

If operates on a parameter, the input should be in the same base as that in the input file. This function will convert the new value to system-base per unit.

Parameters:
src : str

The parameter name to alter

idx : str, float, int

The device to alter

value : float

The desired value

class_name

Return the class name

doc(max_width=78, export='plain')[source]

Retrieve model documentation as a string.

e_clear()[source]

Clear equation value arrays associated with all internal variables.

f_numeric(**kwargs)[source]

Custom fcall functions. Modify equations directly.

f_update()[source]

Evaluate differential equations.

Notes

In-place equations: added to the corresponding DAE array. Non-inplace equations: in-place set to internal array to overwrite old values (and avoid clearing).

g_numeric(**kwargs)[source]

Custom gcall functions. Modify equations directly.

g_update()[source]

Evaluate algebraic equations.

get(src: str, idx, attr: str = 'v', allow_none=False, default=0.0)[source]

Get the value of an attribute of a model property.

The return value is self.<src>.<attr>[idx]

Parameters:
src : str

Name of the model property

idx : str, int, float, array-like

Indices of the devices

attr : str, optional, default='v'

The attribute of the property to get. v for values, a for address, and e for equation value.

allow_none : bool

True to allow None values in the indexer

default : float

If allow_none is true, the default value to use for None indexer.

Returns:
array-like

self.<src>.<attr>[idx]

get_init_order()[source]

Get variable initialization order and send to logger.info.

get_inputs(refresh=False)[source]

Get an OrderedDict of the inputs to the numerical function calls.

Parameters:
refresh : bool

Refresh the values in the dictionary. This is only used when the memory address of arrays changed. After initialization, all array assignments are inplace. To avoid overhead, refresh should not be used after initialization.

Returns:
OrderedDict

The input name and value array pairs in an OrderedDict

Notes

dae.t is now a numpy.ndarray which has stable memory. There is no need to refresh dat_t in this version.

get_md5()[source]

Return the md5 hash of concatenated equation strings.

get_times()[source]

Get event switch_times from TimerParam.

Returns:
list

A list containing all switching times defined in TimerParams

idx2uid(idx)[source]

Convert idx to the 0-indexed unique index.

Parameters:
idx : array-like, numbers, or str

idx of devices

Returns:
list

A list containing the unique indices of the devices

init(routine)[source]

Numerical initialization of a model.

Initialization sequence: 1. Sequential initialization based on the order of definition 2. Use Newton-Krylov method for iterative initialization 3. Custom init

j_numeric(**kwargs)[source]

Custom numeric update functions.

This function should append indices to _ifx, _jfx, and append anonymous functions to _vfx. It is only called once by store_sparse_pattern.

j_update()[source]

Update Jacobian elements.

Values are stored to Model.triplets[jname], where jname is a jacobian name.

Returns:
None
l_check_eq()[source]

Call the check_eq method of discrete components to update equation-dependent flags.

This function should be called after equation updates. AntiWindup limiters use it to append pegged states to the x_set list.

Returns:
None
l_update_var(dae_t, *args, niter=None, err=None, **kwargs)[source]

Call the check_var method of discrete components to update the internal status flags.

The function is variable-dependent and should be called before updating equations.

Returns:
None
list2array()[source]

Convert all the value attributes v to NumPy arrays.

Value attribute arrays should remain in the same address afterwards. Namely, all assignments to value array should be operated in place (e.g., with [:]).

numba_jitify(parallel=False, cache=False)[source]

Optionally convert self.calls.f and self.calls.g to JIT compiled functions.

This function can be turned on by setting System.config.numba to 1.

Warning

This feature is experimental and does not guarantee a speed up. In fact, the program will likely end up slower due to compilation.

post_init_check()[source]

Post init checking. Warns if values of InitChecker is not True.

prepare(quick=False, pycode_path=None)[source]

Symbolic processing and code generation.

refresh_inputs()[source]

This is the helper function to refresh inputs.

The functions collects object references into OrderedDict self._input and self._input_z.

Returns:
None
refresh_inputs_arg()[source]

Refresh inputs for each function with individual argument list.

s_numeric(**kwargs)[source]

Custom service value functions. Modify Service.v directly.

s_numeric_var(**kwargs)[source]

Custom variable service value functions. Modify VarService.v directly.

This custom numerical function is evaluated at each step/iteration before equation update.

s_update()[source]

Update service equation values.

This function is only evaluated at initialization. Service values are updated sequentially. The v attribute of services will be assigned at a new memory address.

s_update_post()[source]

Update post-initialization services.

s_update_var()[source]

Update VarService.

set(src, idx, attr, value)[source]

Set the value of an attribute of a model property.

Performs self.<src>.<attr>[idx] = value.

Parameters:
src : str

Name of the model property

idx : str, int, float, array-like

Indices of the devices

attr : str, optional, default='v'

The internal attribute of the property to get. v for values, a for address, and e for equation value.

value : array-like

New values to be set

Returns:
bool

True when successful.

set_in_use()[source]

Set the in_use attribute. Called at the end of System.collect_ref.

This function is overloaded by models with BackRef to disable calls when no model is referencing. Models with no back references will have internal variable addresses assigned but external addresses being empty.

For internal equations that has external variables, the row indices will be non-zeros, while the col indices will be empty, which causes an error when updating Jacobians.

Setting self.in_use to False when len(back_ref_instance.v) == 0 avoids this error. See COI.

solve_iter(name, kwargs)[source]

Solve iterative initialization.

solve_iter_single(name, inputs, pos)[source]

Solve iterative initialization for one given device.

store_sparse_pattern()[source]

Store rows and columns of the non-zeros in the Jacobians for building the sparsity pattern.

This function converts the internal 0-indexed equation/variable address to the numerical addresses for the loaded system.

Calling sequence: For each Jacobian name, fx, fy, gx and gy, store by a) generated constant and variable Jacobians c) user-provided constant and variable Jacobians, d) user-provided block constant and variable Jacobians

Notes

If self.n == 0, skipping this function will avoid appending empty lists/arrays and non-empty values, which, as a combination, is not accepted by kvxopt.spmatrix.

switch_action(dae_t)[source]

Call the switch actions.

Parameters:
dae_t : float

Current simulation time

Returns:
None

Warning

Timer exported from blocks are supposed to work but have not been tested.

v_numeric(**kwargs)[source]

Custom variable initialization function.

class andes.core.model.ModelCache[source]

Bases: object

Class for caching the return value of callback functions.

Check ModelCache.__dict__.keys() for fields.

add_callback(name: str, callback)[source]

Add a cache attribute and a callback function for updating the attribute.

Parameters:
name : str

name of the cached function return value

callback : callable

callback function for updating the cached attribute

refresh(name=None)[source]

Refresh the cached values

Parameters:
name : str, list, optional

name or list of cached to refresh, by default None for refreshing all

class andes.core.model.ModelCall[source]

Bases: object

Class for storing generated function calls, Jacobian calls, and arguments.

append_ijv(j_full_name, ii, jj, vv)[source]
clear_ijv()[source]
zip_ijv(j_full_name)[source]

Return a zipped iterator for the rows, cols and vals for the specified matrix name.

class andes.core.model.ModelData(*args, three_params=True, **kwargs)[source]

Bases: object

Class for holding parameter data for a model.

This class is designed to hold the parameter data separately from model equations. Models should inherit this class to define the parameters from input files.

Inherit this class to create the specific class for holding input parameters for a new model. The recommended name for the derived class is the model name with Data. For example, data for GENCLS should be named GENCLSData.

Parameters should be defined in the __init__ function of the derived class.

Refer to andes.core.param for available parameter types.

Notes

Three default parameters are pre-defined in ModelData and will be inherited by all models. They are

In rare cases one does not want to define these three parameters, one can pass three_params=True to the constructor of ModelData.

Examples

If we want to build a class PQData (for static PQ load) with three parameters, Vn, p0 and q0, we can use the following

from andes.core.model import ModelData, Model
from andes.core.param import IdxParam, NumParam

class PQData(ModelData):
    super().__init__()
    self.Vn = NumParam(default=110,
                       info="AC voltage rating",
                       unit='kV', non_zero=True,
                       tex_name=r'V_n')
    self.p0 = NumParam(default=0,
                       info='active power load in system base',
                       tex_name=r'p_0', unit='p.u.')
    self.q0 = NumParam(default=0,
                       info='reactive power load in system base',
                       tex_name=r'q_0', unit='p.u.')

In this example, all the three parameters are defined as andes.core.param.NumParam. In the full PQData class, other types of parameters also exist. For example, to store the idx of owner, PQData uses

self.owner = IdxParam(model='Owner', info="owner idx")
Attributes:
cache

A cache instance for different views of the internal data.

flags : dict

Flags to control the routine and functions that get called. If the model is using user-defined numerical calls, set f_num, g_num and j_num properly.

add(**kwargs)[source]

Add a device (an instance) to this model.

Parameters:
kwargs

model parameters are collected into the kwargs dictionary

Warning

This function is not intended to be used directly. Use the add method from System so that the index can be registered correctly.

as_df(vin=False)[source]

Export all parameters as a pandas.DataFrame object. This function utilizes as_dict for preparing data.

Returns:
DataFrame

A dataframe containing all model data. An uid column is added.

vin : bool

If True, export all parameters from original input (vin).

as_dict(vin=False)[source]

Export all parameters as a dict.

Returns:
dict

a dict with the keys being the ModelData parameter names and the values being an array-like of data in the order of adding. An additional uid key is added with the value default to range(n).

find_idx(keys, values, allow_none=False, default=False)[source]

Find idx of devices whose values match the given pattern.

Parameters:
keys : str, array-like, Sized

A string or an array-like of strings containing the names of parameters for the search criteria

values : array, array of arrays, Sized

Values for the corresponding key to search for. If keys is a str, values should be an array of elements. If keys is a list, values should be an array of arrays, each corresponds to the key.

allow_none : bool, Sized

Allow key, value to be not found. Used by groups.

default : bool

Default idx to return if not found (missing)

Returns:
list

indices of devices

find_param(prop)[source]

Find params with the given property and return in an OrderedDict.

Parameters:
prop : str

Property name

Returns:
OrderedDict
update_from_df(df, vin=False)[source]

Update parameter values from a DataFrame.

Adding devices are not allowed.

andes.core.param module

class andes.core.param.BaseParam(default: Union[float, str, int, None] = None, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, mandatory: bool = False, export: bool = True, iconvert: Optional[Callable] = None, oconvert: Optional[Callable] = None)[source]

Bases: object

The base parameter class.

This class provides the basic data structure and interfaces for all types of parameters. Parameters are from input files and in general constant once initialized.

Subclasses should overload the n() method for the total count of elements in the value array.

Parameters:
default : str or float, optional

The default value of this parameter if None is provided

name : str, optional

Parameter name. If not provided, it will be automatically set to the attribute name defined in the owner model.

tex_name : str, optional

LaTeX-formatted parameter name. If not provided, tex_name will be assigned the same as name.

info : str, optional

Descriptive information of parameter

mandatory : bool

True if this parameter is mandatory

export : bool

True if the parameter will be exported when dumping data into files. True for most parameters. False for BackRef.

Warning

The most distinct feature of BaseParam, DataParam and IdxParam is that values are stored in a list without conversion to array. BaseParam, DataParam or IdxParam are not allowed in equations.

Attributes:
v : list

A list holding all the values. The BaseParam class does not convert the v attribute into NumPy arrays.

property : dict

A dict containing the truth values of the model properties.

add(value=None)[source]

Add a new parameter value (from a new device of the owner model) to the v list.

Parameters:
value : str or float, optional

Parameter value of the new element. If None, the default will be used.

Notes

If the value is math.nan, it will set to None.

class_name

Return the class name.

get_names()[source]

Return self.name in a list.

This is a helper function to provide the same API as blocks or discrete components.

Returns:
list

A list only containing the name of the parameter

get_property(property_name: str)[source]

Check the boolean value of the given property. If the property does not exist in the dictionary, False will be returned.

Parameters:
property_name : str

Property name

Returns:
The truth value of the property.
n

Return the count of elements in the value array.

set(pos, attr, value)[source]

Set attributes of the BaseParam class to new values at the given positions.

Parameters:
pos : int, list of integers

Positions in arrays where the values should be set

attr : 'v', 'vin'

Name of the attribute to be set

value : str, float or list of above

New values

set_all(attr, value)[source]

Set attributes of the BaseParam class to new values for all positions.

Parameters:
attr : 'v', 'vin'

Name of the attribute to be set

value : list of str, float or int

New values

class andes.core.param.DataParam(default: Union[float, str, int, None] = None, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, mandatory: bool = False, export: bool = True, iconvert: Optional[Callable] = None, oconvert: Optional[Callable] = None)[source]

Bases: andes.core.param.BaseParam

An alias of the BaseParam class.

This class is used for string parameters or non-computational numerical parameters. This class does not provide a to_array method. All input values will be stored in v as a list.

See also

andes.core.param.BaseParam
Base parameter class
class andes.core.param.ExtParam(model: str, src: str, indexer=None, vtype=<class 'float'>, allow_none=False, default=0.0, **kwargs)[source]

Bases: andes.core.param.NumParam

A parameter whose values are retrieved from an external model or group.

Parameters:
model : str

Name of the model or group providing the original parameter

src : str

The source parameter name

indexer : BaseParam

A parameter defined in the model defining this ExtParam instance. indexer.v should contain indices into model.src.v. If is None, the source parameter values will be fully copied. If model is a group name, the indexer cannot be None.

Attributes:
parent_model : Model

The parent model providing the original parameter.

add(value=None)[source]

ExtParam has an empty add method.

Update parameter values provided by external models. This needs to be called before pu conversion.

Parameters:
ext_model : Model, Group

Instance of the parent model or group, provided by the System calling this method.

restore()[source]

ExtParam has an empty restore method

to_array()[source]

Convert to array when d_type is not str

class andes.core.param.IdxParam(default: Union[float, str, int, None] = None, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, mandatory: bool = False, unique: bool = False, export: bool = True, model: Optional[str] = None, iconvert: Optional[Callable] = None, oconvert: Optional[Callable] = None)[source]

Bases: andes.core.param.BaseParam

An alias of BaseParam with an additional storage of the owner model name

This class is intended for storing idx into other models. It can be used in the future for data consistency check.

Notes

This will be useful when, for example, one connects two TGs to one SynGen.

Examples

A PQ model connected to Bus model will have the following code

class PQModel(...):
    def __init__(...):
        ...
        self.bus = IdxParam(model='Bus')
add(value=None)[source]

Add a new parameter value (from a new device of the owner model) to the v list.

Parameters:
value : str or float, optional

Parameter value of the new element. If None, the default will be used.

Notes

If the value is math.nan, it will set to None.

class andes.core.param.NumParam(default: Union[float, str, Callable, None] = None, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, vrange: Union[List[T], Tuple, None] = None, vtype: Optional[Type[CT_co]] = <class 'float'>, iconvert: Optional[Callable] = None, oconvert: Optional[Callable] = None, non_zero: bool = False, non_positive: bool = False, non_negative: bool = False, mandatory: bool = False, power: bool = False, ipower: bool = False, voltage: bool = False, current: bool = False, z: bool = False, y: bool = False, r: bool = False, g: bool = False, dc_voltage: bool = False, dc_current: bool = False, export: bool = True)[source]

Bases: andes.core.param.BaseParam

A computational numerical parameter.

Parameters defined using this class will have their v field converted to a NumPy array after adding.

The original input values will be copied to vin, and the system-base per-unit conversion coefficients (through multiplication) will be stored in pu_coeff.

Parameters:
default : str or float, optional

The default value of this parameter if no value is provided

name : str, optional

Name of this parameter. If not provided, name will be set to the attribute name of the owner model.

tex_name : str, optional

LaTeX-formatted parameter name. If not provided, tex_name will be assigned the same as name.

info : str, optional

A description of this parameter

mandatory : bool

True if this parameter is mandatory

unit : str, optional

Unit of the parameter

vrange : list, tuple, optional

Typical value range

vtype : type, optional

Type of the v field. The default is float.

Other Parameters:
 
Sn : str

Name of the parameter for the device base power.

Vn : str

Name of the parameter for the device base voltage.

non_zero : bool

True if this parameter must be non-zero. non_zero can be combined with non_positive or non_negative.

non_positive : bool

True if this parameter must be non-positive.

non_negative : bool

True if this parameter must be non-negative.

mandatory : bool

True if this parameter must not be None.

power : bool

True if this parameter is a power per-unit quantity under the device base.

iconvert : callable

Callable to convert input data from excel or others to the internal v field.

oconvert : callable

Callable to convert input data from internal type to a serializable type.

ipower : bool

True if this parameter is an inverse-power per-unit quantity under the device base.

voltage : bool

True if the parameter is a voltage pu quantity under the device base.

current : bool

True if the parameter is a current pu quantity under the device base.

z : bool

True if the parameter is an AC impedance pu quantity under the device base.

y : bool

True if the parameter is an AC admittance pu quantity under the device base.

r : bool

True if the parameter is a DC resistance pu quantity under the device base.

g : bool

True if the parameter is a DC conductance pu quantity under the device base.

dc_current : bool

True if the parameter is a DC current pu quantity under device base.

dc_voltage : bool

True if the parameter is a DC voltage pu quantity under device base.

add(value=None)[source]

Add a value to the parameter value list.

In addition to BaseParam.add, this method checks for non-zero property and reset to default if is zero.

See also

BaseParam.add
the add method of BaseParam
restore()[source]

Restore parameter to the original input by copying self.vin to self.v.

pu_coeff will not be overwritten.

set_pu_coeff(coeff)[source]

Store p.u. conversion coefficient into self.pu_coeff and calculate the system-base per unit with self.v = self.vin * self.pu_coeff.

This function must be called after self.to_array.

Parameters:
coeff : np.ndarray

An array with the pu conversion coefficients

to_array()[source]

Converts field v to the NumPy array type. to enable array-based calculation.

Must be called after adding all elements. Store a copy of original input values to field vin. Set pu_coeff to all ones.

Warning

After this call, add will not be allowed to avoid unexpected issues.

class andes.core.param.TimerParam(callback: Optional[Callable] = None, default: Union[float, str, Callable, None] = None, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, non_zero: bool = False, mandatory: bool = False, export: bool = True)[source]

Bases: andes.core.param.NumParam

A parameter whose values are event occurrence times during the simulation.

The constructor takes an additional Callable self.callback for the action of the event. TimerParam has a default value of -1, meaning deactivated.

Examples

A connectivity status toggler class Toggler takes a parameter t for the toggle time. Inside Toggler.__init__, one would have

self.t = TimerParam()

The Toggler class also needs to define a method for togging the connectivity status

def _u_switch(self, is_time: np.ndarray):
    action = False
    for i in range(self.n):
        if is_time[i] and (self.u.v[i] == 1):
            instance = self.system.__dict__[self.model.v[i]]
            # get the original status and flip the value
            u0 = instance.get(src='u', attr='v', idx=self.dev.v[i])
            instance.set(src='u',
                         attr='v',
                         idx=self.dev.v[i],
                         value=1-u0)
            action = True
    return action

Finally, in Toggler.__init__, assign the function as the callback for self.t

self.t.callback = self._u_switch
is_time(dae_t)[source]

Element-wise check if the DAE time is the same as the parameter value. The current implementation uses np.equal.

Parameters:
dae_t : float

Current simulation time

Returns:
np.ndarray

The array containing the truth value of if the DAE time is close to the parameter value.

Notes

The previous implementation with np.isclose with default rtol=1e-5 mistakes the immediate pre- and post-event time as in-event when simulation time is greater than 10.

andes.core.service module

class andes.core.service.ApplyFunc(u, func, name=None, tex_name=None, info=None, cache=True)[source]

Bases: andes.core.service.BaseService

Class for applying a numerical function on a parameter..

Parameters:
u

Input parameter

func

A condition function that returns True or False.

Warning

This class is not ready.

v
class andes.core.service.BackRef(**kwargs)[source]

Bases: andes.core.service.BaseService

A special type of reference collector.

BackRef is used for collecting device indices of other models referencing the parent model of the BackRef. The v``field will be a list of lists, each containing the `idx of other models referencing each device of the parent model.

BackRef can be passed as indexer for params and vars, or shape for NumReduce and NumRepeat. See examples for illustration.

See also

andes.core.service.NumReduce
A more complete example using BackRef to build the COI model

Examples

A Bus device has an IdxParam of area, storing the idx of area to which the bus device belongs. In Bus.__init__(), one has

self.area = IdxParam(model='Area')

Suppose Bus has the following data

idx area Vn
1 1 110
2 2 220
3 1 345
4 1 500

The Area model wants to collect the indices of Bus devices which points to the corresponding Area device. In Area.__init__, one defines

self.Bus = BackRef()

where the member attribute name Bus needs to match exactly model name that Area wants to collect idx for. Similarly, one can define self.ACTopology = BackRef() to collect devices in the ACTopology group that references Area.

The collection of idx happens in andes.system.System._collect_ref_param(). It has to be noted that the specific Area entry must exist to collect model idx-dx referencing it. For example, if Area has the following data

idx
1

Then, only Bus 1, 3, and 4 will be collected into self.Bus.v, namely, self.Bus.v == [ [1, 3, 4] ].

If Area has data

idx
1
2

Then, self.Bus.v will end up with [ [1, 3, 4], [2] ].

class andes.core.service.BaseService(name: str = None, tex_name: str = None, info: str = None, vtype: Type[CT_co] = None)[source]

Bases: object

Base class for Service.

Service is a v-provider type for holding internal and temporary values. Subclasses need to implement v as a member attribute or using a property decorator.

Parameters:
name : str

Instance name

Attributes:
owner : Model

The hosting/owner model instance

assign_memory(n)[source]

Assign memory for self.v and set the array to zero.

Parameters:
n : int

Number of elements of the value array. Provided by caller (Model.list2array).

class_name

Return the class name

get_names()[source]

Return name in a list

Returns:
list

A list only containing the name of the service variable

n

Return the count of values in self.v.

Needs to be overloaded if v of subclasses is not a 1-dimensional array.

Returns:
int

The count of elements in this variable

class andes.core.service.ConstService(v_str: Optional[str] = None, v_numeric: Optional[Callable] = None, vtype: Optional[type] = None, name: Optional[str] = None, tex_name=None, info=None)[source]

Bases: andes.core.service.BaseService

A type of Service that stays constant once initialized.

ConstService are usually constants calculated from parameters. They are only evaluated once in the initialization phase before variables are initialized. Therefore, uninitialized variables must not be used in v_str`.

Parameters:
name : str

Name of the ConstService

v_str : str

An equation string to calculate the variable value.

v_numeric : Callable, optional

A callable which returns the value of the ConstService

Attributes:
v : array-like or a scalar

ConstService value

class andes.core.service.CurrentSign(bus, bus1, bus2, name=None, tex_name=None, info=None)[source]

Bases: andes.core.service.ConstService

Service for computing the sign of the current flowing through a series device.

With a given line connecting bus1 and bus2, one can compute the current flow using (v1*exp(1j*a1) - v2*exp(1j*a2)) / (r + jx) whose value is the outflow on bus1.

CurrentSign can be used to compute the sign to be multiplied depending on the observing bus. For each value in bus, the sign will be +1 if it appears in bus1 or -1 otherwise.

bus1          bus2
 *------>>-----*
bus(+)        bus(-)
check(**kwargs)[source]
class andes.core.service.DataSelect(optional, fallback, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None)[source]

Bases: andes.core.service.BaseService

Class for selecting values for optional DataParam or NumParam.

This service is a v-provider that uses optional DataParam if available with a fallback.

DataParam will be tested for None, and NumParam will be tested with np.isnan().

Notes

An use case of DataSelect is remote bus. One can do

self.buss = DataSelect(option=self.busr, fallback=self.bus)

Then, pass self.buss instead of self.bus as indexer to retrieve voltages.

Another use case is to allow an optional turbine rating. One can do

self.Tn = NumParam(default=None)
self.Sg = ExtParam(...)
self.Sn = DataSelect(Tn, Sg)
v
class andes.core.service.DeviceFinder(u, link, idx_name, name=None, tex_name=None, info=None)[source]

Bases: andes.core.service.BaseService

Service for finding indices of optionally linked devices.

If not provided, DeviceFinder will add devices at the beginning of System.setup.

Examples

IEEEST stabilizer takes an optional busf (IdxParam) for specifying the connected BusFreq, which is needed for mode 6. To avoid reimplementing BusFreq within IEEEST, one can do

self.busfreq = DeviceFinder(self.busf, link=self.buss, idx_name='bus')

where self.busf is the optional input, self.buss is the bus indices that busf should measure, and idx_name is the name of a BusFreq parameter through which the measured bus indices are specified. For each None values in self.busf, a BusFreq is created to measure the corresponding bus in self.buss.

That is, BusFreq.[idx_name].v = [link]. DeviceFinder will find / create BusFreq devices so that the returned list of BusFreq indices are connected to self.buss, respectively.

find_or_add(system)[source]

Find or add devices.

Points self.u.v to the found or newly added devices.

Find devices one by one. Devices previously added in this function can be used later without duplication.

v
class andes.core.service.EventFlag(u, vtype: Optional[type] = None, name: Optional[str] = None, tex_name=None, info=None)[source]

Bases: andes.core.service.VarService

Service to flag events when the input value changes. The typical input is a v-provider with binary values.

Implemented by providing self.check(**kwargs) as v_numeric. EventFlag.v stores the values of the input variable in the most recent iteration/step.

After the evaluation of self.check(), self.v will be updated.

check(**kwargs)[source]

Check status and set event flags.

Input values are compared with values in the memory.

class andes.core.service.ExtService(model: str, src: str, indexer: Union[andes.core.param.BaseParam, andes.core.service.BaseService], attr: str = 'v', allow_none: bool = False, default=0, name: str = None, tex_name: str = None, vtype=None, info: str = None)[source]

Bases: andes.core.service.BaseService

Service constants whose value is from an external model or group.

Parameters:
src : str

Variable or parameter name in the source model or group

model : str

A model name or a group name

indexer : IdxParam or BaseParam

An "Indexer" instance whose v field contains the idx of devices in the model or group.

Examples

A synchronous generator needs to retrieve the p and q values from static generators for initialization. ExtService is used for this purpose.

In a synchronous generator, one can define the following to retrieve StaticGen.p as p0:

class GENCLSModel(Model):
    def __init__(...):
        ...
        self.p0 = ExtService(src='p',
                             model='StaticGen',
                             indexer=self.gen,
                             tex_name='P_0')

Method to be called by System for getting values from the external model or group.

Parameters:
ext_model

An instance of a model or group provided by System

class andes.core.service.ExtendedEvent(u, t_ext: Union[int, float, andes.core.param.BaseParam, andes.core.service.BaseService] = 0.0, trig: str = 'rise', enable=True, v_disabled=0, extend_only=False, vtype: Optional[type] = None, name: Optional[str] = None, tex_name=None, info=None)[source]

Bases: andes.core.service.VarService

Service for indicating an event for an extended, predefined period of time following the event disappearance.

The triggering of an event, whether the rise or fall edge, is specified through trig. For example, if trig = rise, the change of the input from 0 to 1 will be considered as an input, whereas the subsequent change back to 0 will be considered as the event end.

ExtendedEvent.v stores the flags whether the extended time has completed. Outputs will become 1 once the event starts and return to 0 when the extended time ends.

Parameters:
u : v-provider

Triggering signal where the values are 0 or 1.

trig : str in ("rise", "fall")

Triggering edge for the beginning of an event. rise by default.

enable : bool or v-provider

If disabled, the output will be v_disabled

extend_only : bool

Only output during the extended period, not the event period.

Warning

The performance of this class needs to be optimized.

assign_memory(n)[source]

Assign memory for internal data.

check(**kwargs)[source]

Check if an extended event is in place.

Supplied as a v_numeric to VarService.

class andes.core.service.FlagCondition(u, func, flag=1, name=None, tex_name=None, info=None, cache=True)[source]

Bases: andes.core.service.BaseService

Class for flagging values based on a condition function.

By default, values whose condition function output equal that equal to True/1 will be flagged as 1. 0 otherwise.

Parameters:
u

Input parameter

func

A condition function that returns True or False.

flag : 1 by default, only 0 or 1 is accepted.

The flag for the inputs whose condition output is True.

Warning

This class is not ready.

FlagCondition can only be applied to BaseParam with cache=True. Applying to Service will fail unless cache is False (at a performance cost).

v
class andes.core.service.FlagGreaterThan(u, value=0.0, flag=1, equal=False, name=None, tex_name=None, info=None, cache=True)[source]

Bases: andes.core.service.FlagCondition

Service for flagging parameters > or >= the given value element-wise.

Parameters that satisfy the comparison (u > or >= value) will flagged as flag (1 by default).

class andes.core.service.FlagLessThan(u, value=0.0, flag=1, equal=False, name=None, tex_name=None, info=None, cache=True)[source]

Bases: andes.core.service.FlagCondition

Service for flagging parameters < or <= the given value element-wise.

Parameters that satisfy the comparison (u < or <= value) will flagged as flag (1 by default).

class andes.core.service.FlagValue(u, value, flag=0, name=None, tex_name=None, info=None, cache=True)[source]

Bases: andes.core.service.BaseService

Class for flagging values that equal to the given value.

By default, values that equal to value will be flagged as 0. Non-matching values will be flagged as 1.

Parameters:
u

Input parameter

value

Value to flag. Can be None, string, or a number.

flag : 0 by default, only 0 or 1 is accepted.

The flag for the matched ones

Warning

FlagNotNone can only be applied to BaseParam with cache=True. Applying to Service will fail unless cache is False (at a performance cost).

v
class andes.core.service.IdxRepeat(u, ref, **kwargs)[source]

Bases: andes.core.service.OperationService

Helper class to repeat IdxParam.

This class has the same functionality as andes.core.service.NumRepeat but only operates on IdxParam, DataParam or NumParam.

v

Return values stored in self._v. May be overloaded by subclasses.

class andes.core.service.InitChecker(u, lower=None, upper=None, equal=None, not_equal=None, enable=True, error_out=False, **kwargs)[source]

Bases: andes.core.service.OperationService

Class for checking init values against known typical values.

Instances will be stored in Model.services_post and Model.services_icheck, which will be checked in Model.post_init_check() after initialization.

Parameters:
u

v-provider to be checked

lower : float, BaseParam, BaseVar, BaseService

lower bound

upper : float, BaseParam, BaseVar, BaseService

upper bound

equal : float, BaseParam, BaseVar, BaseService

values that the value from v_str should equal

not_equal : float, BaseParam, BaseVar, BaseService

values that should not equal

enable : bool

True to enable checking

Examples

Let's say generator excitation voltages are known to be in the range of 1.6 - 3.0 per unit. One can add the following instance to GENBase

self._vfc = InitChecker(u=self.vf,
                        info='vf range',
                        lower=1.8,
                        upper=3.0,
                        )

lower and upper can also take v-providers instead of float values.

One can also pass float values from Config to make it adjustable as in our implementation of GENBase._vfc.

check()[source]

Check the bounds and equality conditions.

class andes.core.service.NumReduce(u, ref: andes.core.service.BackRef, fun: Callable, name=None, tex_name=None, info=None, cache=True)[source]

Bases: andes.core.service.OperationService

A helper Service type which reduces a linearly stored 2-D ExtParam into 1-D Service.

NumReduce works with ExtParam whose v field is a list of lists. A reduce function which takes an array-like and returns a scalar need to be supplied. NumReduce calls the reduce function on each of the lists and return all the scalars in an array.

Parameters:
u : ExtParam

Input ExtParam whose v contains linearly stored 2-dimensional values

ref : BackRef

The BackRef whose 2-dimensional shapes are used for indexing

fun : Callable

The callable for converting a 1-D array-like to a scalar

Examples

Suppose one wants to calculate the mean value of the Vn in one Area. In the Area class, one defines

class AreaModel(...):
    def __init__(...):
        ...
        # backward reference from `Bus`
        self.Bus = BackRef()

        # collect the Vn in an 1-D array
        self.Vn = ExtParam(model='Bus',
            src='Vn',
            indexer=self.Bus)

        self.Vn_mean = NumReduce(u=self.Vn,
            fun=np.mean,
            ref=self.Bus)

Suppose we define two areas, 1 and 2, the Bus data looks like

idx area Vn
1 1 110
2 2 220
3 1 345
4 1 500

Then, self.Bus.v is a list of two lists [ [1, 3, 4], [2] ]. self.Vn.v will be retrieved and linearly stored as [110, 345, 500, 220]. Based on the shape from self.Bus, numpy.mean() will be called on [110, 345, 500] and [220] respectively. Thus, self.Vn_mean.v will become [318.33, 220].

v

Return the reduced values from the reduction function in an array

Returns:
The array self._v storing the reduced values
class andes.core.service.NumRepeat(u, ref, **kwargs)[source]

Bases: andes.core.service.OperationService

A helper Service type which repeats a v-provider's value based on the shape from a BackRef

Examples

NumRepeat was originally designed for computing the inertia-weighted average rotor speed (center of inertia speed). COI speed is computed with

\[\omega_{COI} = \frac{ \sum{M_i * \omega_i} } {\sum{M_i}}\]

The numerator can be calculated with a mix of BackRef, ExtParam and ExtState. The denominator needs to be calculated with NumReduce and Service Repeat. That is, use NumReduce to calculate the sum, and use NumRepeat to repeat the summed value for each device.

In the COI class, one would have

class COIModel(...):
    def __init__(...):
        ...
        self.SynGen = BackRef()
        self.SynGenIdx = RefFlatten(ref=self.SynGen)
        self.M = ExtParam(model='SynGen',
                          src='M',
                          indexer=self.SynGenIdx)

        self.wgen = ExtState(model='SynGen',
                             src='omega',
                             indexer=self.SynGenIdx)

        self.Mt = NumReduce(u=self.M,
                                 fun=np.sum,
                                 ref=self.SynGen)

        self.Mtr = NumRepeat(u=self.Mt,
                               ref=self.SynGen)

        self.pidx = IdxRepeat(u=self.idx,ref=self.SynGen)

Finally, one would define the center of inertia speed as

self.wcoi = Algeb(v_str='1', e_str='-wcoi')

self.wcoi_sub = ExtAlgeb(model='COI',
                         src='wcoi',
                         e_str='M * wgen / Mtr',
                         v_str='M / Mtr',
                         indexer=self.pidx,
                         )

It is very worth noting that the implementation uses a trick to separate the average weighted sum into n sub-equations, each calculating the \((M_i * \omega_i) / (\sum{M_i})\). Since all the variables are preserved in the sub-equation, the derivatives can be calculated correctly.

v

Return the values of the repeated values in a sequential 1-D array

Returns:
The array, self._v storing the repeated values
class andes.core.service.NumSelect(optional, fallback, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None)[source]

Bases: andes.core.service.OperationService

Class for selecting values for optional NumParam.

NumSelect works with internal and external parameters.

Notes

One use case is to allow an optional turbine rating. One can do

self.Tn = NumParam(default=None)
self.Sg = ExtParam(...)
self.Sn = DataSelect(Tn, Sg)
v

Return values stored in self._v. May be overloaded by subclasses.

class andes.core.service.OperationService(name=None, tex_name=None, info=None)[source]

Bases: andes.core.service.BaseService

Base class for a type of Service which performs specific operations. OperationService may not use the assign_memory from BaseService, because it can have a different size.

This class cannot be used by itself.

See also

NumReduce
Service for Reducing linearly stored 2-D services into 1-D
NumRepeat
Service for repeating 1-D NumParam/ v-array following a sub-pattern
IdxRepeat
Service for repeating 1-D IdxParam/ v-list following a sub-pattern
v

Return values stored in self._v. May be overloaded by subclasses.

class andes.core.service.ParamCalc(param1, param2, func, name=None, tex_name=None, info=None, cache=True)[source]

Bases: andes.core.service.BaseService

Parameter calculation service.

Useful to create parameters calculated instantly from existing ones.

v
class andes.core.service.PostInitService(v_str: Optional[str] = None, v_numeric: Optional[Callable] = None, vtype: Optional[type] = None, name: Optional[str] = None, tex_name=None, info=None)[source]

Bases: andes.core.service.ConstService

Constant service that gets stored once after init.

This service is useful when one need to store initialization values stored in variables.

Examples

In ESST3A model, the vf variable is initialized followed by other variables. One can store the initial vf into vf0 so that equation vf - vf0 = 0 will hold.

self.vref0 = PostInitService(info='Initial reference voltage input',
                             tex_name='V_{ref0}',
                             v_str='vref',
                             )

Since all ConstService are evaluated before equation evaluation, without using PostInitService, one will need to create lots of ConstService to store values in the initialization path towards vf0, in order to correctly initialize vf.

class andes.core.service.RandomService(func=<built-in method rand of numpy.random.mtrand.RandomState object>, **kwargs)[source]

Bases: andes.core.service.BaseService

A service type for generating random numbers.

Parameters:
name : str

Name

func : Callable

A callable for generating the random variable.

Warning

The value will be randomized every time it is accessed. Do not use it if the value needs to be stable for each simulation step.

v

This class has v wrapped by a property decorator.

Returns:
array-like

Randomly generated service variables

class andes.core.service.RefFlatten(ref, **kwargs)[source]

Bases: andes.core.service.OperationService

A service type for flattening andes.core.service.BackRef into a 1-D list.

Examples

This class is used when one wants to pass BackRef values as indexer.

andes.models.coi.COI collects referencing andes.models.group.SynGen with

self.SynGen = BackRef(info='SynGen idx lists', export=False)

After collecting BackRefs, self.SynGen.v will become a two-level list of indices, where the first level correspond to each COI and the second level correspond to generators of the COI.

Convert self.SynGen into 1-d as self.SynGenIdx, which can be passed as indexer for retrieving other parameters and variables

self.SynGenIdx = RefFlatten(ref=self.SynGen)

self.M = ExtParam(model='SynGen', src='M',
                  indexer=self.SynGenIdx, export=False,
                  )
v

Return values stored in self._v. May be overloaded by subclasses.

class andes.core.service.Replace(old_val, flt, new_val, name=None, tex_name=None, info=None, cache=True)[source]

Bases: andes.core.service.BaseService

Replace parameters with new values if the function returns True

v
class andes.core.service.SwBlock(*, init, ns, blocks, ext_sel=None, name=None, tex_name=None, info=None)[source]

Bases: andes.core.service.OperationService

Service type for switched shunt blocks.

adjust(amount)[source]

Adjust capacitor banks by an amount.

check_data()[source]

Check data consistency.

find_sel()[source]

Determine the initial shunt selection level.

set_v()[source]

Set values to _v based on sel.

v

Return values stored in self._v. May be overloaded by subclasses.

class andes.core.service.VarHold(u, hold, vtype=None, name=None, tex_name=None, info=None)[source]

Bases: andes.core.service.VarService

Service for holding the input when the hold signal is on.

Parameters:
hold : v-provider, binary

Hold signal array with length equal to the input. For elements that are 1, the corresponding inputs are held until the hold signal returns to 0.

check(**kwargs)[source]

Custom v_numeric function for checking the hold signal and calculating outputs.

class andes.core.service.VarService(v_str: Optional[str] = None, v_numeric: Optional[Callable] = None, vtype: Optional[type] = None, name: Optional[str] = None, tex_name=None, info=None)[source]

Bases: andes.core.service.ConstService

Variable service that gets updated in each step/loop as variables change.

This class is useful when one has non-differentiable algebraic equations, which make use of abs(), re and im. Instead of creating Algeb, one can put the equation in VarService, which will be updated before solving algebraic equations.

Warning

VarService is not solved with other algebraic equations, meaning that there is one step "delay" between the algebraic variables and VarService. Use an algebraic variable whenever possible.

Examples

In ESST3A model, the voltage and current sensors (vd + jvq), (Id + jIq) estimate the sensed VE using equation

\[VE = | K_{PC}*(v_d + 1j v_q) + 1j (K_I + K_{PC}*X_L)*(I_d + 1j I_q)|\]

One can use VarService to implement this equation

self.VE = VarService(
    tex_name='V_E',
    info='VE',
    v_str='Abs(KPC*(vd + 1j*vq) + 1j*(KI + KPC*XL)*(Id + 1j*Iq))',
    )

andes.core.solver module

Sparse solvers wrapper.

class andes.core.solver.CuPySolver[source]

Bases: andes.core.solver.SciPySolver

CuPy lsqr solver (GPU-based).

solve(A, b)[source]

Solve linear systems.

Parameters:
A : scipy.csc_matrix

Sparse N-by-N matrix

b : numpy.ndarray

Dense 1-dimensional array of size N

Returns:
np.ndarray

Solution x to Ax = b

class andes.core.solver.KLUSolver[source]

Bases: andes.core.solver.SuiteSparseSolver

KLU solver.

Requires package kvxoptklu.

linsolve(A, b)[source]

Solve linear equation set Ax = b and returns the solutions in a 1-D array.

This function performs both symbolic and numeric factorizations every time, and can be slower than Solver.solve.

Parameters:
A

Sparse matrix

b

RHS of the equation

Returns:
The solution in a 1-D np array.
class andes.core.solver.SciPySolver[source]

Bases: object

Base class for scipy family solvers.

clear()[source]
linsolve(A, b)[source]

Exactly same functionality as solve.

solve(A, b)[source]

Solve linear systems.

Parameters:
A : scipy.csc_matrix

Sparse N-by-N matrix

b : numpy.ndarray

Dense 1-dimensional array of size N

Returns:
np.ndarray

Solution x to Ax = b

to_csc(A)[source]

Convert A to scipy.sparse.csc_matrix.

Parameters:
A : kvxopt.spmatrix

Sparse N-by-N matrix

Returns:
scipy.sparse.csc_matrix

Converted csc_matrix

class andes.core.solver.Solver(sparselib='umfpack')[source]

Bases: object

Sparse matrix solver class.

This class wraps UMFPACK, KLU, SciPy and CuPy solvers to provide an unified interface for solving sparse linear equations Ax = b.

Provides methods solve, linsolve and clear.

clear()[source]

Remove all cached objects.

linsolve(A, b)[source]

Solve linear equations without caching facorization. Performs full factorization each call.

Parameters:
A : kvxopt.spmatrix

Sparse N-by-N matrix

b : kvxopt.matrix or numpy.ndarray

Dense N-by-1 matrix

Returns:
numpy.ndarray

Dense N-by-1 array

solve(A, b)[source]

Solve linear equations and cache factorizations if possible.

Parameters:
A : kvxopt.spmatrix

Sparse N-by-N matrix

b : kvxopt.matrix or numpy.ndarray

Dense N-by-1 matrix

Returns:
numpy.ndarray

Dense N-by-1 array

class andes.core.solver.SpSolve[source]

Bases: andes.core.solver.SciPySolver

scipy.sparse.linalg.spsolve Solver.

solve(A, b)[source]

Solve linear systems.

Parameters:
A : scipy.csc_matrix

Sparse N-by-N matrix

b : numpy.ndarray

Dense 1-dimensional array of size N

Returns:
np.ndarray

Solution x to Ax = b

class andes.core.solver.SuiteSparseSolver[source]

Bases: object

Base SuiteSparse solver interface.

Need to be derived by specific solvers such as UMFPACK or KLU.

clear()[source]

Remove all cached PyCapsule of C objects

linsolve(A, b)[source]

Solve linear equation set Ax = b and returns the solutions in a 1-D array.

This function performs both symbolic and numeric factorizations every time, and can be slower than Solver.solve.

Parameters:
A

Sparse matrix

b

RHS of the equation

Returns:
The solution in a 1-D np array.
solve(A, b)[source]

Solve linear system Ax = b using numeric factorization N and symbolic factorization F. Store the solution in b.

This function caches the symbolic factorization in self.F and is faster in general. Will attempt Solver.linsolve if the cached symbolic factorization is invalid.

Parameters:
A

Sparse matrix for the equation set coefficients.

F

The symbolic factorization of A or a matrix with the same non-zero shape as A.

N

Numeric factorization of A.

b

RHS of the equation.

Returns:
numpy.ndarray

The solution in a 1-D ndarray

class andes.core.solver.UMFPACKSolver[source]

Bases: andes.core.solver.SuiteSparseSolver

UMFPACK solver.

Utilizes kvxopt.umfpack for factorization.

linsolve(A, b)[source]

Solve linear equation set Ax = b and returns the solutions in a 1-D array.

This function performs both symbolic and numeric factorizations every time, and can be slower than Solver.solve.

Parameters:
A

Sparse matrix

b

RHS of the equation

Returns:
The solution in a 1-D np array.

andes.core.common module

class andes.core.common.Config(name, dct=None, **kwargs)[source]

Bases: object

A class for storing system, model and routine configurations.

add(dct=None, **kwargs)[source]

Add config fields from a dictionary or keyword args.

Existing configs will NOT be overwritten.

add_extra(dest, dct=None, **kwargs)[source]

Add extra contents for config.

Parameters:
dest : str

Destination string in _alt, _help or _tex.

dct : OrderedDict, dict

key: value pairs

as_dict(refresh=False)[source]

Return the config fields and values in an OrderedDict.

Values are cached in self._dict unless refreshed.

check()[source]

Check the validity of config values.

doc(max_width=78, export='plain', target=False, symbol=True)[source]
load(config)[source]

Load from a ConfigParser object, config.

tex_names
class andes.core.common.DummyValue(value)[source]

Bases: object

Class for converting a scalar value to a dummy parameter with name and tex_name fields.

A DummyValue object can be passed to Block, which utilizes the name field to dynamically generate equations.

Notes

Pass a numerical value to the constructor for most use cases, especially when passing as a v-provider.

class andes.core.common.Indicator[source]

Bases: sympy.core.expr.Expr

Indicator class for printing SymPy Relational.

Relational expressions in SymPy need to be wrapped by Indicator.

Examples

To compare dae_t with 0, one need to use Indicator(dae_t < 0)`.

default_assumptions = {}
class andes.core.common.JacTriplet[source]

Bases: object

Storage class for Jacobian triplet lists.

append_ijv(j_full_name, ii, jj, vv)[source]

Append triplets to the given sparse matrix triplets.

Parameters:
j_full_name : str

Full name of the sparse Jacobian. If is a constant Jacobian, append 'c' to the Jacobian name.

ii : array-like

Row indices

jj : array-like

Column indices

vv : array-like

Value indices

clear_ijv()[source]

Clear stored triplets for all sparse Jacobian matrices

ijv(j_full_name)[source]

Return triplet lists in a tuple in the order or (ii, jj, vv)

merge(triplet)[source]

Merge another triplet into this one.

zip_ijv(j_full_name)[source]

Return a zip iterator in the order of (ii, jj, vv)

class andes.core.common.ModelFlags(collate=False, pflow=False, tds=False, pflow_init=None, tds_init=None, series=False, nr_iter=False, f_num=False, g_num=False, j_num=False, s_num=False, sv_num=False)[source]

Bases: object

Model flags.

Parameters:
collate : bool

True: collate variables by device; False: by variable. Non-collate (continuous memory) has faster computation speed.

pflow : bool

True: called during power flow

tds : bool

True if called during tds; if is False, dae_t cannot be used

pflow_init : bool or None

True if initialize pflow; False otherwise; None default to pflow

tds_init : bool or None

True if initialize tds; False otherwise; None default to tds

series : bool

True if is series device

nr_iter : bool

True if is series device

f_num : bool

True if the model defines f_numeric

g_num : bool

True if the model defines g_numeric

j_num : bool

True if the model defines j_numeric

s_num : bool

True if the model defines s_numeric

sv_num : bool

True if the model defines s_numeric_var

jited : bool

True if numba JIT code is generated

update(dct)[source]
andes.core.common.dummify(param)[source]

Dummify scalar parameter and return a DummyValue object. Do nothing for BaseParam instances.

Parameters:
param : float, int, str, BaseParam

parameter object or scalar value

Returns:
DummyValue(param) if param is a scalar; param itself, otherwise.

andes.core.var module

class andes.core.var.Algeb(name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, v_str: Union[str, float, None] = None, v_iter: Optional[str] = None, e_str: Optional[str] = None, discrete: Optional[andes.core.discrete.Discrete] = None, v_setter: Optional[bool] = False, e_setter: Optional[bool] = False, addressable: Optional[bool] = True, export: Optional[bool] = True, diag_eps: Optional[float] = 0.0, deps: Optional[List[T]] = None)[source]

Bases: andes.core.var.BaseVar

Algebraic variable class, an alias of the BaseVar.

Attributes:
e_code : str

Equation code string, equals string literal g

v_code : str

Variable code string, equals string literal y

e_code = 'g'
v_code = 'y'
class andes.core.var.AliasAlgeb(var, **kwargs)[source]

Bases: andes.core.var.ExtAlgeb

Alias algebraic variable. Essentially ExtAlgeb that links to a a model's own variable.

AliasAlgeb is useful when the final output of a model is from a block, but the model must provide the final output in a pre-defined name. Using AliasAlgeb, A model can avoid adding an additional variable with a dummy equations.

Like ExtVar, labels of AliasAlgeb will not be saved in the final output. When plotting from file, one need to look up the original variable name.

class andes.core.var.AliasState(var, **kwargs)[source]

Bases: andes.core.var.ExtState

Alias state variable.

Refer to the docs of AliasAlgeb.

class andes.core.var.BaseVar(name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, v_str: Union[str, float, None] = None, v_iter: Optional[str] = None, e_str: Optional[str] = None, discrete: Optional[andes.core.discrete.Discrete] = None, v_setter: Optional[bool] = False, e_setter: Optional[bool] = False, addressable: Optional[bool] = True, export: Optional[bool] = True, diag_eps: Optional[float] = 0.0, deps: Optional[List[T]] = None)[source]

Bases: object

Base variable class.

Derived classes State and Algeb should be used to build model variables.

Parameters:
name : str, optional

Variable name

info : str, optional

Descriptive information

unit : str, optional

Unit

tex_name : str

LaTeX-formatted variable name. If is None, use name instead.

discrete : Discrete

Associated discrete component. Will call check_var on the discrete component.

Attributes:
a : array-like

variable address

v : array-like

local-storage of the variable value

e : array-like

local-storage of the corresponding equation value

e_str : str

the string/symbolic representation of the equation

class_name
get_names()[source]
reset()[source]

Reset the internal numpy arrays and flags.

set_address(addr: numpy.ndarray, contiguous=False)[source]

Set the address of internal variables.

Parameters:
addr : np.ndarray

The assigned address for this variable

contiguous : bool, optional

If the addresses are contiguous

set_arrays(dae, inplace=True, alloc=True)[source]

Set the equation and values arrays.

Parameters:
dae : DAE

Reference to System.dae

class andes.core.var.ExtAlgeb(model: str, src: str, indexer: Union[List[T], numpy.ndarray, andes.core.param.BaseParam, andes.core.service.BaseService, None] = None, allow_none: Optional[bool] = False, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, v_str: Union[str, float, None] = None, v_iter: Optional[str] = None, e_str: Optional[str] = None, v_setter: Optional[bool] = False, e_setter: Optional[bool] = False, addressable: Optional[bool] = True, export: Optional[bool] = True, diag_eps: Optional[float] = 0.0)[source]

Bases: andes.core.var.ExtVar

External algebraic variable type.

e_code = 'g'
v_code = 'y'
class andes.core.var.ExtState(model: str, src: str, indexer: Union[List[T], numpy.ndarray, andes.core.param.BaseParam, andes.core.service.BaseService, None] = None, allow_none: Optional[bool] = False, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, v_str: Union[str, float, None] = None, v_iter: Optional[str] = None, e_str: Optional[str] = None, v_setter: Optional[bool] = False, e_setter: Optional[bool] = False, addressable: Optional[bool] = True, export: Optional[bool] = True, diag_eps: Optional[float] = 0.0)[source]

Bases: andes.core.var.ExtVar

External state variable type.

Warning

ExtState is not allowed to set t_const, as it will conflict with the source State variable. In fact, one should not set e_str for ExtState.

e_code = 'f'
t_const = None
v_code = 'x'
class andes.core.var.ExtVar(model: str, src: str, indexer: Union[List[T], numpy.ndarray, andes.core.param.BaseParam, andes.core.service.BaseService, None] = None, allow_none: Optional[bool] = False, name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, v_str: Union[str, float, None] = None, v_iter: Optional[str] = None, e_str: Optional[str] = None, v_setter: Optional[bool] = False, e_setter: Optional[bool] = False, addressable: Optional[bool] = True, export: Optional[bool] = True, diag_eps: Optional[float] = 0.0)[source]

Bases: andes.core.var.BaseVar

Externally defined algebraic variable

This class is used to retrieve the addresses of externally- defined variable. The e value of the ExtVar will be added to the corresponding address in the DAE equation.

Parameters:
model : str

Name of the source model

src : str

Source variable name

indexer : BaseParam, BaseService

A parameter of the hosting model, used as indices into the source model and variable. If is None, the source variable address will be fully copied.

allow_none : bool

True to allow None in indexer

Attributes:
parent_model : Model

The parent model providing the original parameter.

uid : array-like

An array containing the absolute indices into the parent_instance values.

e_code : str

Equation code string; copied from the parent instance.

v_code : str

Variable code string; copied from the parent instance.

Update variable addresses provided by external models

This method sets attributes including parent_model, parent_instance, uid, a, n, e_code and v_code. It initializes the e and v to zero.

Parameters:
ext_model : Model

Instance of the parent model

Returns:
None

Warning

link_external does not check if the ExtVar type is the same as the original variable to reduce performance overhead. It will be a silent error (a dimension too small error from dae.build_pattern) if a model uses ExtAlgeb to access a State, or vice versa.

set_address(addr, contiguous=False)[source]

Assigns address for equation RHS.

set_arrays(dae, inplace=True, alloc=True)[source]

Stride of dae.i for the RHS of external variables.

class andes.core.var.State(name: Optional[str] = None, tex_name: Optional[str] = None, info: Optional[str] = None, unit: Optional[str] = None, v_str: Union[str, float, None] = None, v_iter: Optional[str] = None, e_str: Optional[str] = None, discrete: Optional[andes.core.discrete.Discrete] = None, t_const: Union[andes.core.param.BaseParam, andes.core.common.DummyValue, andes.core.service.BaseService, None] = None, check_init: Optional[bool] = True, v_setter: Optional[bool] = False, e_setter: Optional[bool] = False, addressable: Optional[bool] = True, export: Optional[bool] = True, diag_eps: Optional[float] = 0.0, deps: Optional[List[T]] = None)[source]

Bases: andes.core.var.BaseVar

Differential variable class, an alias of the BaseVar.

Parameters:
t_const : BaseParam, DummyValue

Left-hand time constant for the differential equation. Time constants will not be evaluated as part of the differential equation. They will be collected to array dae.Tf to multiply to the right-hand side dae.f.

check_init : bool

True to check if the equation right-hand-side is zero initially. Disabling the checking can be used for integrators when the initial input may not be zero.

Attributes:
e_code : str

Equation code string, equals string literal f

v_code : str

Variable code string, equals string literal x

e_code = 'f'
v_code = 'x'

Module contents

Import subpackage classes