Commit ea58dd96 authored by Sabine Auer's avatar Sabine Auer
Browse files

microgrid test case added

parent 5ea3e694
......@@ -92,9 +92,8 @@ def define_network_rhs(node_list, Y):
length = len(Y)
total_length = 2*length
coupling_sp = sps.csr_matrix(Y)
def network_rhs(y, t):
coupling_sp = sps.csr_matrix(Y)
dydt = np.empty(total_length + length)
v = y[:total_length].view(np.complex128)
......@@ -198,6 +197,7 @@ if __name__ == "__main__":
from scipy.optimize import root
result = root(root_rhs, ic)
if result.success:
print(rhs(result.x, 0.))
......
from __future__ import division, print_function, absolute_import, unicode_literals
import numpy as np
import baobap as bao
from numba import njit
import os
from complex_current_and_nodes import SwingEquationNode,define_network_rhs,define_network_rhs_codegen, define_root_rhs
default_dir = os.path.join("simulation_data", "microgrid")
def load_PYPSA(filename):
data = np.load(filename)
input_power = data['input_power']
Y = data['Y']
system_size = len(input_power)
node_list = [SwingEquationNode(3, infeed=input_power[i], H=0.1, damping=0.01) for i in range(system_size)]
return node_list,Y
def define_gen_rc(brp, rhs):
system_size = brp.system_dimension
def generate_run_conditions(batch, run):
ic = np.ones(3 * system_size, dtype=np.float64) * 0.5
root_rhs = define_root_rhs(2 * system_size, rhs)
from scipy.optimize import root
result = root(root_rhs, ic)
if result.success:
ic = result.x
ic[2 * system_size + batch+1] += 0.1
else:
print("failed")
return ic, ()
return generate_run_conditions
# This function observes the results of the simulation run and returns its observations in a dict
# noinspection PyUnusedLocal
def lorenz_ob(time_series, rc):
return{"max_values": np.max(time_series, axis=0)}
def main(sim_dir=default_dir, create_test_data=True, run_test=True,flag_baobab=True):
node_list, Y = load_PYPSA('microgrid_testcase.npz')
rhs = define_network_rhs(node_list, Y)
# nrhs = define_network_rhs_codegen(node_list, Y)
if flag_baobab:
# Run PYPSD code with BAOBAB
np.random.seed(0)
result_file, batch_dir = bao.prep_dir(sim_dir)
times = np.linspace(0, 10, 100)
brp = bao.BatchRunParameters(number_of_batches=10, simulations_per_batch=10, system_dimension=100)
ob = bao.combine_obs(bao.run_condition_observer, lorenz_ob)
bao.run_experiment(result_file, brp, None, define_gen_rc(brp,rhs), ob, times, figures=(True, 1), rhs=rhs)
res = bao.load_all_fields_from_results(result_file)
if bao.am_master:
if create_test_data:
for key in res.keys():
np.save(os.path.join(sim_dir, key), res[key])
if run_test:
for key in res.keys():
filename = key + ".npy"
with open(os.path.join(sim_dir, filename), mode='rb') as f:
arr = np.load(f)
assert np.allclose(arr, res[key])
else:
# Run Test PYPSD without BAOBAB
times = np.linspace(0., 10., 100)
from scipy.integrate import odeint
brp = bao.BatchRunParameters(number_of_batches=10, simulations_per_batch=10, system_dimension=100)
generate_run_conditions = define_gen_rc(brp,rhs)
rc = generate_run_conditions(0,0)
#node_list[0].infeed += 0.1
rhs = define_network_rhs(node_list, Y)
states = odeint(rhs, rc[0], times)
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
plt.figure()
plt.plot(times, states[:, 2*brp.system_dimension:])
fig = plt.figure()
ax = fig.gca(projection='3d')
ax.plot(states[:, 0], states[:, 1], times)
ax.plot(states[:, 2], states[:, 3], times)
plt.show()
if __name__ == "__main__":
main()
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment