Add custom control structures

[mdbartos]

Taking pumps as an example, analogous functions for the following are needed:

operator_flow_coefficients

pipedream/pipedream_solver/nsuperlink.py

Lines 978 to 1009 in 0fc0159

	def pump_flow_coefficients(self, u=None):
	"""
	Compute pump flow coefficients: alpha_up, beta_up, chi_up,
	alpha_dp, beta_dp, chi_dp.
	"""
	# Import instance variables
	H_j = self.H_j # Head at superjunction j
	_z_inv_j = self._z_inv_j # Invert elevation at superjunction j
	_J_up = self._J_up # Index of superjunction upstream of pump p
	_J_dp = self._J_dp # Index of superjunction downstream of pump p
	_z_p = self._z_p # Offset of pump inlet above upstream invert elevation
	_dHp_max = self._dHp_max # Maximum pump head difference
	_dHp_min = self._dHp_min # Minimum pump head difference
	_a_p = self._a_p # Pump curve parameter `a`
	_b_p = self._b_p # Pump curve parameter `b`
	_c_p = self._c_p # Pump curve parameter `c`
	_Qp = self._Qp # Current flow rate through pump p
	_alpha_p = self._alpha_p # Pump flow coefficient alpha_p
	_beta_p = self._beta_p # Pump flow coefficient beta_p
	_chi_p = self._chi_p # Pump flow coefficient chi_p
	# If no input signal, assume pump is closed
	if u is None:
	u = np.zeros(self.n_p, dtype=np.float64)
	# Check max/min head differences
	assert (_dHp_min <= _dHp_max).all()
	# Compute pump flow coefficients
	numba_pump_flow_coefficients(_alpha_p, _beta_p, _chi_p, H_j, _z_inv_j, _Qp, u,
	_z_p, _dHp_max, _dHp_min, _a_p, _b_p, _c_p, _J_up, _J_dp)
	# Export instance variables
	self._alpha_p = _alpha_p
	self._beta_p = _beta_p
	self._chi_p = _chi_p

solve_operator_flows

pipedream/pipedream_solver/nsuperlink.py

Lines 1468 to 1490 in 0fc0159

	def solve_pump_flows(self, u=None):
	"""
	Solve for pump discharges given superjunction heads at time t + dt.
	"""
	# Import instance variables
	H_j = self.H_j # Head at superjunction j
	_z_inv_j = self._z_inv_j # Invert elevation of superjunction j
	_J_up = self._J_up # Index of superjunction upstream of pump p
	_J_dp = self._J_dp # Index of superjunction downstream of pump p
	_z_p = self._z_p # Offset of pump inlet above upstream invert
	_dHp_max = self._dHp_max # Maximum pump head difference
	_dHp_min = self._dHp_min # Minimum pump head difference
	_a_p = self._a_p # Pump curve parameter `a`
	_b_p = self._b_p # Pump curve parameter `b`
	_c_p = self._c_p # Pump curve parameter `c`
	_Qp = self._Qp # Current flow rate through pump p
	# If no input signal, assume pump is closed
	if u is None:
	u = np.zeros(self.n_p, dtype=np.float64)
	# Compute pump flows
	_Qp_next = numba_solve_pump_flows(H_j, u, _z_inv_j, _z_p, _dHp_max,
	_dHp_min, _a_p, _b_p, _c_p, _J_up, _J_dp)
	self._Qp = _Qp_next

In sparse_matrix_equations, analogous lines for the following are needed:

pipedream/pipedream_solver/nsuperlink.py

Lines 1067 to 1078 in 0fc0159

	if n_p:
	P = self.P
	_J_up = self._J_up # Index of superjunction upstream of pump p
	_J_dp = self._J_dp # Index of superjunction downstream of pump p
	_alpha_p = self._alpha_p # Pump flow coefficient
	_beta_p = self._beta_p # Pump flow coefficient
	_chi_p = self._chi_p # Pump flow coefficient
	_alpha_upm = self._alpha_upm # Summation of pump flow coefficients
	_beta_dpl = self._beta_dpl # Summation of pump flow coefficients
	_chi_upl = self._chi_upl # Summation of pump flow coefficients
	_chi_dpm = self._chi_dpm # Summation of pump flow coefficients
	_P_diag = self._P_diag # Diagonal elements of matrix P

pipedream/pipedream_solver/nsuperlink.py

Lines 1145 to 1159 in 0fc0159

	if n_p:
	_alpha_up = _alpha_p
	_alpha_dp = _alpha_p
	_beta_up = _beta_p
	_beta_dp = _beta_p
	_chi_up = _chi_p
	_chi_dp = _chi_p
	_P_diag.fill(0)
	numba_clear_off_diagonals(P, bc, _J_up, _J_dp, n_p)
	# Set diagonal
	numba_create_OWP_matrix(P, _P_diag, bc, _J_up, _J_dp, _alpha_up,
	_alpha_dp, _beta_up, _beta_dp, M, n_p)
	# Set right-hand side
	numba_add_at(D, _J_up, -_chi_up)
	numba_add_at(D, _J_dp, _chi_dp)

In solve_sparse_matrix, a new matrix must be added to the sum:

pipedream/pipedream_solver/nsuperlink.py

Line 1198 in 0fc0159

l = A + O + W + P

And same goes for solve_banded_matrix:

pipedream/pipedream_solver/nsuperlink.py

Line 1243 in 0fc0159

l = A + O + W + P