3.2. Running simulations

3.2.1. Initializing

We can initialize an instance of the simulator by importing pyramses and invoking the pyramses.sim class:

import pyramses
ram = pyramses.sim()

You can get more information on running simulations in Running a simulation. The available calls to the pyramses.sim class are fully documented in Full list of functions.

3.2.2. Start, Pause and progress

To run a simulation, we first need to setup a proper test case (see Building a scenario and loading data). Then we can start the simulation with pyramses.sim.execSim

ram.execSim(case) # start and run until the end of the simulation

If we want to pause the simulation at a particular point, we can use:

ram.execSim(case, 10.0) # start the simulation and pause at time t=10 seconds

To continue the simulation and pause at a new point, we use pyramses.sim.contSim:

ram.contSim(20.0) # continue simulation until t=20.0 seconds

To simulate until the end:

ram.contSim(ram.getInfTime()) # simulate until the end (time horizon was reached or an early stop happened due to system violations or collapse)

3.2.3. End simulation

To end the simulation without simulating to the end of the time horizon, we can use pyramses.sim.endSim:

ram.endSim() # end simulation and exit

3.2.4. Full list of functions

class pyramses.sim(custLibDir=None)[source]

Simulation class.

addDisturb(t_dist, disturb)[source]

Add a new disturbance at a specific time. Follows the same structure as the disturbances in the dst files.

Parameters
  • t_dist (float) – time of the disturbance

  • disturb (str) – description of disturbance

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 80.0) # simulate until 80 seconds and pause
>>> ram.addDisturb(100.000, 'CHGPRM DCTL 1-1041  Vsetpt -0.05 0') # Decrease the setpoint of the DCTL by 0.015 pu, at t=100 s
>>> ram.addDisturb(100.000, 'CHGPRM DCTL 2-1042  Vsetpt -0.05 0')
>>> ram.addDisturb(100.000, 'CHGPRM DCTL 3-1043  Vsetpt -0.05 0')
>>> ram.addDisturb(100.000, 'CHGPRM DCTL 4-1044  Vsetpt -0.05 0')
>>> ram.addDisturb(100.000, 'CHGPRM DCTL 5-1045  Vsetpt -0.05 0')
>>> ram.contSim(ram.getInfTime()) # continue the simulation
addObserv(string)[source]

Add an element to be observed

Parameters

string (str) – the element with proper format

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt") # command file without any observables
>>> ram.execSim(case, 0.0) # start
>>> traj_filenm = 'obs.trj'
>>> ram.initObserv(traj_filenm)
>>> string = 'BUS *' # monitor all buses
>>> ram.addObserv(string)
contSim(pause=None)[source]

Continue the simulation until the given time.

The pause argument is optional and it gives the time until the simulation will stop.

defineSS(ssID, filter1, filter2, filter3)[source]

Define a subsytem using three filters. The resulting list is an intersection of the filters.

Parameters
  • ssID (int) – Number of the SS

  • filter1 (list of str or str) – Voltage levels to be included

  • filter2 (list of str or str) – Zones (which zone/zones will be included)

  • filter3 (list of str or str) – Bus names to be included

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 0.0)
>>> ram.defineSS(1, ['735'], [], []) # SS 1 with all buses at 735 kV, no zones, no list of buses

Note

An empty filter means it is deactivated and discarded.

endSim()[source]

End the simulation

execSim(cmd, pause=None)[source]

Run a simulation.

Parameters
  • cmd (type(cfg)) – provides the case description

  • t_pause (float) – pause time (optional). If not given, the simulation will run until the end as described in the dst file.

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case) # start simulation

Note

If you have an existing command file, you can pass it to the simulator as pyramses.cfg(“cmd.txt”)

finalObserv()[source]

Finalize observable selection, allocate buffer, and write header of trajectory file

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt") # command file without any observables
>>> ram.execSim(case, 0.0) # start
>>> traj_filenm = 'obs.trj'
>>> ram.initObserv(traj_filenm)
>>> string = 'BUS *' # monitor all buses
>>> ram.addObserv(string)
>>> ram.finalObserv()
getAllCompNames(comp_type)[source]

Get list of all components of type comp_type

Parameters

comp_type (str) – the component type (BUS, SYNC, INJ, DCTL, BRANCH, TWOP, SHUNT, LOAD)

Returns

list of component names

Return type

list of str

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 0) # start simulation paused
>>> ram.getAllCompNames("BUS")
['B1',
'B2',
'B3']
getBranchCur(branchName)[source]

Return the currents of a list of branches

Parameters

branchName (list of str) – the names of branches

Returns

list of branch currents. These are x-y components at the origin and extremity respectively (ix_orig, iy_orig, ix_extr, iy_extr)

Return type

list of floats

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 10.0) # simulate until 10 seconds and pause
>>> branchName = ['1011-1013','1012-1014','1021-1022']
>>> ram.getBranchCur(branchName)
getBranchPow(branchName)[source]

Return the active and reactive powers of a list of branches

Parameters

branchName (list of str) – the names of branches

Returns

list of branch powers. These are active and reactive power at the origin and extremity respectively (p_orig, q_orig, p_extr, q_extr)

Return type

list of floats

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 10.0) # simulate until 10 seconds and pause
>>> branchName = ['1011-1013','1012-1014','1021-1022']
>>> ram.getBranchPow(branchName)
getBusPha(busNames)[source]

Return the voltage phase of a list of buses

Parameters

busNames (list of str) – the names of buses

Returns

list of bus voltage phase

Return type

list of floats

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 10.0) # simulate until 10 seconds and pause
>>> buses = ['g1','g2','g3']
>>> ram.getBusPha(buses)
[0.0000000000000000,
 10.0615442362180327,
 11.064702686997689]
getBusVolt(busNames)[source]

Return the voltage magnitude of a list of buses

Parameters

busNames (list of str) – the names of buses

Returns

list of bus voltage magnitudes

Return type

list of floats

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 10.0) # simulate until 10 seconds and pause
>>> buses = ['g1','g2','g3']
>>> ram.getBusVolt(buses)
[1.0736851673414456,
 1.0615442362180327,
 1.064702686997689]
getCompName(comp_type, num)[source]

Return the name of the num:superscript:th component of type comp_type.

Parameters
  • comp_type (str) – the component type (BUS, SYNC, INJ, DCTL, BRANCH, TWOP, SHUNT, LOAD)

  • num (int) – The number of the component to be fetched

Returns

The component name

Return type

srt

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 0) # start simulation paused
>>> ram.getCompName("BUS",1)
'B1'
getEndSim()[source]

Check if the simulation has ended

Returns

0 -> simulation is still working, 1 -> simulation has ended

Return type

integer

getInfTime()[source]

Get the maximum representable double from the simulator ()

Returns

maximum representable double

Return type

float

getJac()[source]

Gets the Jacobian matrix written in files

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 0) # start simulation paused
>>> ram.getJac()
getLastErr()[source]

Return the last error message issued by RAMSES.

Return type

str

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 0) # start simulation paused
>>> ram.getLastErr()
'This is an error msg'
getObs(comp_type, comp_name, obs_name)[source]

Get the value of a named observable.

Parameters
  • comp_type (list of str) – the types of components (‘EXC’,’TOR’,’INJ’,’TWOP’,’DCTL’,’SYN’)

  • comp_name (list of str) – the names of components

  • obs_name (list of str) – the names of observables

Returns

list of observable values

Return type

list of floats

Note

For the synchronous generator (‘SYN’) the accepted obs_name values are - ‘P’: Active power (MW) - ‘Q’: Reactive power (MVAr) -‘Omega’: Machine speed (pu) - ‘S’: Apparent power (MVA) - ‘SNOM’: Nominal apparent power (MVA) - ‘PNOM’: Nominal active power (MW)

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 10.0) # simulate until 10 seconds and pause
>>> comp_type = ['INJ','EXC','TOR']
>>> comp_name = ['L_11','g2','g3']
>>> obs_name = ['P','vf','Pm']
>>> ram.getObs(comp_type,comp_name, obs_name)
[199.73704732259995,
1.3372945282585218,
0.816671075203357]
getPrm(comp_type, comp_name, prm_name)[source]

Get the value of a named parameter

Parameters
  • comp_type (list of str) – the types of components (EXC, TOR, INJ, DCTL, TWOP)

  • comp_name (list of str) – the names of components

  • prm_name (list of str) – the names of parameters

Returns

list of parameter values

Return type

list of floats

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 10.0) # simulate until 10 seconds and pause
>>> comp_type = ['EXC','EXC','EXC']
>>> comp_name = ['g1','g2','g3']
>>> prm_name = ['V0','V0','V0']
>>> ram.getPrm(comp_type,comp_name, prm_name)
[1.0736851673414456,
 1.0615442362180327,
 1.064702686997689]
getPrmNames(comp_type, comp_name)[source]

Get the named parameters of a model

Parameters
  • comp_type (list of str or str) – the types of components (EXC, TOR, INJ, DCTL, TWOP)

  • comp_name (list of str or str) – the names of component instances

Returns

list of parameter names

Return type

list of lists of strings

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 0.0) # initialize and wait
>>> comp_type = ['EXC','EXC','EXC']
>>> comp_name = ['g1','g2','g3'] # name of synchronous machines
>>> ram.getPrmNames(comp_type,comp_name)
getSS(ssID)[source]

Retrieve the buses of a subsytem.

Parameters

ssID (int) – Number of the SS

Returns

list of buses

Return type

list of str

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 0.0)
>>> ram.defineSS(1, ['735'], [], []) # SS 1 with all buses at 735 kV
>>> ram.getSS(1) # get list of buses in SS 1
getSimTime()[source]

Get the current simulation time

Returns

current simulation time in RAMSES.

Return type

float

getTrfoSS(ssID, location, in_service, rettype)[source]

Retrieve transformer information of subsystem after applying some filters

Parameters
  • ssID (int) – Number of the SS

  • location (int) – 1 – both buses inside SS, 2 - tie transformers, 3 – 1 & 2

  • in_service (int) – 1 – transformers in service, 2 - all transformers

  • rettype (str) – Type of response (NAME, From, To, Status, Tap, Currentf, Currentt, Pf, Qf, Pt, Qt).

Returns

list of transformer names

Return type

list of str

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case, 0.0)
>>> ram.defineSS(1, ['735'], [], []) # SS 1 with all buses at 735 kV
>>> ram.getTrfoSS(1,3,2,'Status')

Note

Tap is not implemented yet.

Todo

Implement Taps after discussing with Lester.

get_MDL_no()[source]

Unload external DLL file with user defined models. Should be in current directory or absolute path.

Returns

list of observable values

Return type

list of floats

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> ram = pyramses.load_MDL("MDLs.dll")
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case) # simulate
>>> ram = pyramses.get_MDL_no()
initObserv(traj_filenm)[source]

Initialize observable selection (structure and trajectory file)

Parameters

traj_filenm (str) – the file to save the trajectory at

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> case = pyramses.cfg("cmd.txt") # command file without any observables
>>> ram.execSim(case, 0.0) # start
>>> traj_filenm = 'obs.trj'
>>> ram.initObserv(traj_filenm)
load_MDL()[source]

Load external DLL file with user defined models. Should be in current directory or absolute path.

Parameters

MDLName (str) – path to file

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> ram = pyramses.load_MDL("MDLs.dll")
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case) # simulate
pauseSim(t_pause)[source]

Pause the simulation at the t_pause time

Parameters

t_pause (float) – pause time

unload_MDL()[source]

Unload external DLL file with user defined models. Should be in current directory or absolute path.

Parameters

MDLName (str) – path to file

Example

>>> import pyramses
>>> ram = pyramses.sim()
>>> ram = pyramses.load_MDL("MDLs.dll")
>>> case = pyramses.cfg("cmd.txt")
>>> ram.execSim(case) # simulate
>>> ram = pyramses.unload_MDL("MDLs.dll")