Circuit Solver - V

A little progress. At least the ODE is solvable now. I used the strategy of row operations to get rid of the elements with low (or zero) time constants from as many rows as possible. This is the code (click on "view raw" below the code box to see it in a new window):



So for a circuit:

The matrices are transformed as follows:

Original matrices are:
------------------------------------------------------------------------------------------------------------------
e=
0.0  0.0  0.0  0.0
0.0  0.2  0.0  0.1
0.0  0.0  0.0001  0.0001
0.0  0.1  0.0001  0.2001

a=
10000000.0  0.0  10000000.0  10000000.0
0.0  20.0  0.0  10.0
10000000.0  0.0  10000010.1  10000000.1
10000000.0  10.0  10000000.1  10000020.1

b=
1.0
0.0
0.0
0.0
------------------------------------------------------------------------------------------------------------------

These are manipulated to:
------------------------------------------------------------------------------------------------------------------
e=
0.0  0.0  0.0  0.0
0.0  0.2  0.0  0.1
0.0  0.0  0.0001  0.0001
0.0  0.1  0.0001  0.2001

a=
10000000.0  0.0  10000000.0  10000000.0
0.0  20.0  0.0  10.0
0.0  0.0  10.0999999996  0.0999999996275
0.0  10.0  0.0999999996275  20.0999999996

b=
1.0
0.0
-1.0
-1.0
------------------------------------------------------------------------------------------------------------------

With this the ODE solves with a time step of 10.0e-6 seconds. However, the stiff equation still is not solved satifsfactorily because of the entire row of large resistances. The way, loop current i1 would be calculated would be:

i1=(-10000000.0*i3-10000000.0*i4)/10000000.0

This would cause i1 to be of the same order as the other loop currents as this reduces to:

i1=-i3-i4

Which is wrong.

So this makes me wonder, why are the off-diagonal terms present in that case? Why not just do:

i1=u/10000000.0

Seems like that might be the solution. Because all I need is the sum total of all the resistance between the input signal in the loop to give me the approximate current. But this is then eventually ONLY the approximate current.

Circuit Solver - IV

Spent the weekend reading on different integration techniques and got only more confused.

So need to investigate a little more. Take the simple circuit below:

These are the parameters and the loops found:
--------------------------------------------------------------------------------------------------------------------
Resistor is  R4 = 10.000000  located at  9I
Inductor is  L4 =0.100000  located at  9K
Ammeter is  A5  located at  9J  with positive polarity towards 9K
Ammeter is  A2  located at  3H  with positive polarity towards 4H
Resistor is  R3 = 10.000000  located at  6I
Ammeter is  A4  located at  6J  with positive polarity towards 6K
Inductor is  L3 =0.100000  located at  6K
Ammeter is  A1  located at  1B  with positive polarity towards 1C
Resistor is  R1 = 0.100000  located at  1E
Resistor is  C1 = 10.000000  located at  4H
Inductor is  L1 =0.000100  located at  1F
Resistor is  R2 = 10.000000  located at  1I
Inductor is  L2 =0.100000  located at  1K
Ammeter is  A3  located at  1J  with positive polarity towards 1K
Resistor is  V1 = 100000.000000  located at  4D
Voltage Source is  V1 of 120.000000 V(peak), 60.000000 Hz(frequency) and 0.000000 (degrees phase shift)  located at  4A  with positive polarity towards 3A
**************************************************
Number of nodes 5
Number of branches 8
Number of loops 4
**************************************************
1D 2D 3D 4D 5D 6D 6C 6B 6A 5A 4A 3A 2A 1A 1B 1C 1D
6H 7H 8H 9H 9I 9J 9K 9L 8L 7L 6L 6K 6J 6I 6H
1H 1G 1F 1E 1D 2D 3D 4D 5D 6D 6E 6F 6G 6H 5H 4H 3H 2H 1H
1H 1G 1F 1E 1D 2D 3D 4D 5D 6D 6E 6F 6G 6H 7H 8H 9H 9I 9J 9K 9L 8L 7L 6L 5L 4L 3L 2L 1L 1K 1J 1I 1H
**************************************************
--------------------------------------------------------------------------------------------------------------------

Which is correct. However, these are the system matrices that result:
--------------------------------------------------------------------------------------------------------------------
a=
100000.0  0.0  100000.0  100000.0
0.0  20.0  0.0  10.0
100000.0  0.0  100010.1  100000.1
100000.0  10.0  100000.1  100020.1

e=
0.0  0.0  0.0  0.0
0.0  0.2  0.0  0.1
0.0  0.0  0.0001  0.0001
0.0  0.1  0.0001  0.2001

b=
1.0
0.0
0.0
0.0
--------------------------------------------------------------------------------------------------------------------

I chose resistor V1 as a large resistor. And that has caused all the problems. It appears in three out of four loops. And the entire ODE becomes stiff. One glaring element is in row 4, column 2 of matrix a. This is 10.0 while all other elements are 100000 or greater. This can be one way of reducing the system.

The resistor V1 should appear in at least one loop. Suppose it appears in a loop and we neglect the dynamics of that loop and consider a static equation. This may be OK because if you have a ridiculously low time constant, why solve it unless you want to. So can we reduce the time constants of the other equations to values that we can integrate peacefully? In the above example, resistor V1 appears quite conveniently in a static equation. So we can get rid of it from the other equations.

Performing row operations:

R3=R3-R1
R4=R4-R1

will get rid of these small time constants.

So this makes we wonder, if I need to do a preprocessing to isolate the time constants as far as possible. And if that needs to be done, what should be logic?

Circuit Solver - III

With the SourceForge (SF from now on) project up and going, I will dedicate this blog to description rather than posting code.

One major issue. Choosing the simulation time step. The simulation time step can be chosen larger for speed but choosing a time step too large will result in instability. The concept is as follows.

Each loop has a time constant that is set by the ratio of total inductance to total resistance (L/R) in the loop. In order to simulate the system accurately, you need to choose a time step that is smaller than the smaller L/R ratio in the system. This has the potential of creating a few problems because if you choose a small parasitic inductance in series and a large resistance in parallel, a loop containing these will have a ridiculously low time constants that will need a simulation time step of nanoseconds.

So this means, maybe I need to choose the loops with the largest inductances instead. Or inform the user that the lowest time constant is "x" and a decent value of simulation time step is "y".

Circuit Solver - II

Since the number of programs is increasing and will continue to do so, I guess it is time to create a SourceForge project. Here is the link:
http://sourceforge.net/projects/pythonpowerelec/?source=directory
The programs can be found in the "Files" section.

I will keep posting program updates on the blog anyway.

Circuit Solver

I finally got a basic circuit solver going! What has been done so far:

1. Only voltage source, resistor and inductor have been included. Capacitor has not yet been tested.
2. Ammeters can be introduced in the circuit wherever needed.

So here is the code (click "view raw" below the code box to see the code in a new window):


The main comments are:

1. The loop currents are calculated at every simulation time step. Currents in every branch and every element will not be calculated unless necessary. For example, at a later stage, saturation in the inductor or heating in the resistor can be included. So then, current will be calculated and used.

2. If you want to know the currents, you connect an ammeter. Voltmeters need to be coded. Ammeter was simple as I need to manipulate the loop currents. But voltmeters, I need to think how to use them. If I connect a voltmeter across two nodes, it will be seen as a branch. So it will appear in loop equations. This means, to make sure the currents don't change a large resistance will have to be chosen. However, a simpler method will be to look for the voltmeter tag and exclude that branch from the system matrix. I can then back calculate later.

A sample circuit:

Circuit Parameters - VII

A few changes to the code for taking circuit parameters.

1. I forgot completely "jump" labels. They do not need an object, but the code should ignore them while reading the cells or they will generate "Component not found" errors.
2. Take the circuit layour file name as input from the user and generate the parameters files uniquely for the circuit.

The next stage is to generate the loop matrix for the circuit. The loops have been identified. However, every loop will interact with other loops. So, this has been coded.

Next stage is to read the parameters of every component and generate the ODE for solving the circuit KVL laws.

Anyway, here is the code (click on "view raw" below the code box for code in a new window)

Circuit Parameters - VI

There is a couple of major design flaw in the previous approach.

1. The component names are automatically generated. It is always better to let the user name the components as they may correspond to PCB designs etc.
2. Any change in the circuit topology, even an extension of a branch would be seen a new circuit with default values. Better to check if any of the components in the new circuit have parameters in the parameter file and start with those parameters.

So, a change in the code. The components will now be names followed by an underscore ("_") and then the name of the component. Example: Resistor_R1, voltagesource_v1. Case is not important in the component names though it is in the tags.

The unique identifier continues to be cell position. However, for a particular component type, the same tag can't be used. For example: resistor_r1 and resistor_r1 would be illegal. However, two different components could have the same tag. For example: Resistor_br and inductor_br are OK. If the user want to name the components according to function, this might be a good option.

Here is the code (click on "view raw" below the box to see the entire code):

Circuit Parameters - V

Couple of updates to the previous code:

1. Some elements will also have a polarity. In this case, the voltage source definitely will have one. The capacitor could also have but I'll consider an ac capacitor for now.

2. There will have to be provision to check if everytime the program is run, whether the circuit has changed. If the circuit is the same but the parameters need to be updated, then nw_params.csv should not be filled with default values as the user will have to update everything all over again or will have to remember to copy nw_params.csv to another backup file. A more elaborate method can be thought of later where the actual topology can be read - nodes, branches, loops etc and if that remains the same along with elements in every branch, then circuit can be considered the same. Because otherwise a small change to a branch such as lengthening with more "wire" elements will be seen as a new circuit. For that matter, even if major changes are to be made, only the changed elements can be considered. All this is user interface and will come later.


So anyway, here is the updated code (click on "view raw" below the code box for code going out of the window):

Circuit Paramaters - IV

Now that the basic concept of entering parameters is established, this is the code for it.

Depending on the components found from the circuit spreadsheet, another spreadsheet is created with rows corresponding to every device found. Each device will have default parameters. So the user will have to change those. After all changes, the user can save it as another CSV file.

Each row in the spreadsheet will have the first three elements as - Component type (Resistor, Inductor etc), Component Reference (R1, L2 etc), Cell Position.

As before, the Cell position can be used to access the component object from the component object dictionary. The component classes will have a function to take the parameters which will be column 4 onwards.

Here is the code (click on "view raw" to see the part going out of the box):

Circuit Parameters - III

The next step is to take a sample circuit with a voltage source, resistors, inductors and capacitors and to automate the naming of the elements.

The unique  identification all these elements have is their position in the spreadsheet - only one element can occupy a cell. So the first part is to write a small function that will translate the cell information from [row, column] form to the form that can be read of from the spreadsheet.

That function is here:


So the string that marks the cell position being an immutable element can be the key in a dictionary. A dictionary can be made of all elements found in the circuit and these can be accessed at any time by the position.

Now an identification is chosen, what will be the value of this dictionary item? Ideally, it would be great if it could be the object itself. Then every element will be an object uniquely identified by its position and with the same methods for transferring data to and from the main solver.

This program seems to do this (click on "View Raw" below the code box to see the code going out of the window):


In this file, network_reader is the Python code that contains the circuit interpreter done so far. Each class so far contains only the information of its index (or serial number) and the position in the spreadsheet.
--------------------------------------------------------------------------------------
class Resistor:
    def __init__(self, res_index, res_pos):
        self.res_number=res_index
        self.res_pos=res_pos
    def display(self):
        print "Resistor is ",
        print "R"+str(self.res_number),
        print " located at ",
        print self.res_pos
--------------------------------------------------------------------------------------

All the class references are added to the dictionary:
--------------------------------------------------------------------------------------
component_list={"Resistor":Resistor, "Inductor":Inductor, "Capacitor":Capacitor, "Voltage_Source":Voltage_Source}
--------------------------------------------------------------------------------------

This is what component list looks like:
--------------------------------------------------------------------------------------
>>> component_list
{'Voltage_Source': <class __main__.Voltage_Source at 0x02C9F810>, 'Resistor': <class __main__.Resistor at 0x02C9F768>, 'Inductor': <class __main__.Inductor at 0x02C9F7A0>, 'Capacitor': <class __main__.Capacitor at 0x02C9F7D8>}
--------------------------------------------------------------------------------------


The next step is to generate a dictionary of all the components found in the spreadsheet:
--------------------------------------------------------------------------------------
components_found={}
for c1 in range(len(tst_mat)):
    for c2 in range(len(tst_mat[0])):
        if tst_mat[c1][c2]:
            if tst_mat[c1][c2].lower()!="wire":
                if tst_mat[c1][c2] in component_list.keys():
                    if tst_mat[c1][c2] not in components_found:
                        components_found[tst_mat[c1][c2]]=[nw_rd.csv_element([c1, c2])]
                    else:
                        components_found[tst_mat[c1][c2]].append(nw_rd.csv_element([c1, c2]))
                else:
                    print "Error! Component at %s doesn't exist." %nw_rd.csv_element([c1, c2])
--------------------------------------------------------------------------------------

And this is what the dictionary looks like:
--------------------------------------------------------------------------------------
>>> components_found
{'Voltage_Source': ['4A'], 'Resistor': ['1B', '1F'], 'Inductor': ['1C', '1H'], 'Capacitor': ['4E']}
--------------------------------------------------------------------------------------

From this list, each component can be given a name as they have a specific index in a list.

Finally, to create a dictionary of all the components found. Here the unique identification will be the position in the spreadsheet:
--------------------------------------------------------------------------------------
component_objects={}
for items in components_found.keys():
    for c1 in range(len(components_found[items])):
        component_objects[components_found[items][c1]]=component_list[items](c1+1, components_found[items][c1])
--------------------------------------------------------------------------------------

A quick explanation to what I have done:

1. components_found[items] are the different types of components - voltage source, resistors etc.
2. components_found[items][c1] are the positions of each of these components. So these are the keys of the dictionary component_objects.
3. component_list[items] is the value of the dictionary item corresponding to the component type - and this is a class constructor. This contructor takes the values of the index and the position.
This is the final result:
--------------------------------------------------------------------------------------
>>>
Inductor is  L1  located at  1C
Resistor is  R1  located at  1B
Resistor is  R2  located at  1F
Inductor is  L2  located at  1H
Capacitor is  C1  located at  4E
Voltage Source is  V1  located at  4A
--------------------------------------------------------------------------------------

For the circuit below

Circuit Parameters - II

The next step is to write the code such that all elements are identified and parameters are obtained. In C++ the concept that I followed was that every element was an object - a resistor, an inverter etc. There was three levels of data transfer to these objects. The first was when they were initialized. The second was when the parameters and initial conditions they contained was transferred to the system matrices for solving the ODEs and third was when the results of the ODE was transferred back to the objects to update the values of voltage/current etc.

In C++, the objects were created manually. In Python, I would like the objects to be created automatically from the network layout in the CSV file.

So, testing how to do this first step. A little piece of code to see how automatically the objects can be created and stored. This code is too ridiculously simple to put on Git, so this is it:

-----------------------------------------------------------------------------------------------------------------
#! /usr/bin/env python
import sys, math

class Resistor:
    def __init__(self, res_index):
        self.res_number=res_index
    def display(self):
        print "Resistor is ",
        print "R"+str(self.res_number)

res_list=[]
for c1 in range(1, 5):
    res_name="R"+str(c1)
    res_name=Resistor(c1)
    res_list.append(res_name)

print res_list
for res in res_list:
    res.display()
-----------------------------------------------------------------------------------------------------------------

Here, I am creating a class "Resistor" which has only one data element - its index. So the resistors will be called R1, R2, etc.

I create a list of resistors res_list and create resistors R1 to R4. Also, objects are created with that tag and added to res_list.

Here is the output:
-----------------------------------------------------------------------------------------------------------------
>>>
[<__main__.Resistor instance at 0x02A1F4E0>, <__main__.Resistor instance at 0x02A1F508>, <__main__.Resistor instance at 0x02A1F530>, <__main__.Resistor instance at 0x02A1F558>]
Resistor is  R1
Resistor is  R2
Resistor is  R3
Resistor is  R4
-----------------------------------------------------------------------------------------------------------------

res_list is a list with instances of the class Resistor. But importantly, when I access the method display of the class Resistor by res.display(), it does give me the correct resistor labels.

So now let me put everything in the spreadsheet into objects.

Circuit Parameters

After loop identification, now comes the time to input the parameters of the circuit.

One very nice thing about most circuit simulators is the dialog box that opens up with each element with all the parameters. So essentially, I need to put these features into my circuit simulator:

1. The user should only indicate the element at a given position.
2. We need to then generate a list of all the elements found, names for them, their positions in the spreadsheet and finally their parameters with some to be defaults.
3. With the case of elements like diodes, IGBTs and others that have polarities or certain ways in which they are connected, the polarity can be specified in the separate dialog.

For multiport elements like machines, transformers, there will have to be more than one tag and these will have to be translated. Things would have been much more convenient with a GUI.

This will be my first approach:

1. Create an index of names - "Resistor", "Inductor", "Capacitor", "Voltage_Source" etc.
2. The network interpreter will divide the circuit into three types of elements - "wire", "jump" and elements such as the above.
3. Each time an element is found, it will be assigned a name that creates an object of it - like R1, R2, C1, etc.
4. The program will create another spreadsheet with all the element names, their positions and parameter values in separate columns.
5. The user has to modify these values and continue with the simulation.

For passive circuits with only resistors, inductors, capacitors and voltage sources, I don't see much of a problem. However, there would be a problem with machines for example, a three-phase synchronous generator. So the stator will have three connections and the rotor will have two. How do I map these connections on a spreadsheet?