## Friday, December 28, 2012

### Solving ODEs

One of the major components in developing the circuit simulator is solving the ordinary differential equation (ODE) “E dx/dt = Ax + bu”. Moreover, the matrices need to be determined based on the network topology, values of different components and finally in the case of power electronics, the state of a device. So to begin with, a function to solve the ODE.

It works with the basic resistor, inductor, voltage source:

-------R1----L1---------------R3----L3------
|                             |                                 |
|                             |                                 |
|                             R2                               |
V                            |                                 |
|                             L2                                |
|                              |                                 |
|                              |                                 |
|----------------------|-------------------------|

I hate drawing circuits. I draw so many of them while publishing papers, I can't get myself to draw for this blog. But I think this should give an idea of what I am simulating.

A few things to consider:

1. The code has been written with a consideration for simulation speed. For example, as far as possible, the arguments passed to functions are as references and are used as such. Copy contructors are not invoked and these matrices are changed in the function and reflect in the main program which calls the function.
The reason, if I were to make the functions black bloxes which accept inputs without modifying them and return outputs to be used in the larger algorithm, I would need to copy the matrices each time the function is called. This would mean, copy constructors and destructors called with each function call. A simulation would run for 1 second with a time step of 20 microseconds, which means hundreds of thousands of function calls with the same number of times the constuctors and destructors are called. This would give a normal desktop a severe bellyache.
So, by using references wherever possible, the simulation runs faster as there are minimum number of copy constructors.

2. The function call:
mat_ode(e, a, b, [curr_state_vec, next_state_vec], u, 20.0e-6, ode_solver)
Looks so neat but hides the fact that several vectors are grouped together and passed as a single object:

ode_solver=[ode_k1, ode_k2, ode_k3, ode_k4, ode_dbydt]
Here I must say that I like this feature of Python. You pass anything into a function and it is treated as an object. As long as you know what to do inside the function, you are OK. Moreover, here "ode_solver" is a mutable object - a list of lists (with lists)! So inside the function:

k1=runge_vectors[0]
Will help me unpack the different components.

3. This means I need to take a look a the matrix.py class definition and try to make sure constructors are avoided as far as possible.

4. The ODE solves assuming that E is non-singular. But a look at the circuit will tell you that if L2=L3=0, E will be non-singular but will still be a valid circuit with one variable having a d/dt. So the ODE solver will have to check for singularity and solve spearately as a special case.

Anyway, this my guess will be the last post in 2012.