Table of Contents

Circuit Solver

I'll write a description at some point.

Code

#!/usr/bin/perl
use strict;
use warnings;
 
use Math::Matrix;
use List::Util qw(max);
use Data::Dumper;
 
# Merges two matrices into one bigger matrix, putting them side-by-side.
# Copied from: http://stackoverflow.com/questions/608406/matrices-of-matrices-within-perl/608549#608549
sub hmerge(\@\@;\@\@\@\@\@\@) {
    my @ret;
    for my $i (0 .. max map $#$_, @_) {
        push @ret, [map @{$$_[$i]}, @_];
    }
    @ret;
}
 
# Returns an array of arrays filled with the same value. This can be passed to Math::Matrix.
sub filled_matrix
{
    my $value = shift;
    my $m = shift;
    my $n = shift;
    my @arrays;
 
    for (my $i = 0; $i < $n; ++$i)
    {
    my @temp;
    for (my $j = 0; $j < $m; ++$j)
    {
        push(@temp,$value);
    }
    push(@arrays, \@temp);
    }
 
    return \@arrays;
}
 
sub zeros
{
    return filled_matrix(0,shift,shift);
}
 
my $w;
my @components;
my $nodes = 0;
 
print "Reading components...\n";
 
while(<>)
{
    # Something along the lines of:
    # <R|V><index> <source> <drain> <value>
    if (/^(.\d+)\s+(\d+)\s+(\d+)\s+(\d+)/)
    {
    push(@components, [$1, $2, $3, $4]);
 
    # Tracks the largest node index.
    $nodes = $2 if ($2 > $nodes);
    $nodes = $3 if ($3 > $nodes);
    }
}
 
# Creates a matrix of all of the branch connections.
# The souce would be considered the positive end of a component,
# which really only matters for the voltage source.
# The ground node (node 0) is ignored, since it won't affect the
# calculations.
print "Making A...\n";
my @Aarray;
for (my $i = 1; $i < $nodes+1; ++$i)
{
    my @nodearray;
    foreach my $c (@components)
    {
    if    ($c->[1] == $i) { push(@nodearray, 1); } # Source
    elsif ($c->[2] == $i) { push(@nodearray,-1); } # Drain
    else                  { push(@nodearray, 0); } # No connection
    }
    push(@Aarray, \@nodearray);
}
my $A = Math::Matrix->new(@Aarray);
 
# $A is an M x N matrix. (Width x Height)
my ($An, $Am) = $A->size;
 
# $M always seems to be the identity matrix, but the paper mentions it specifically,
# so it might be best to give it its own variable.
my $M = Math::Matrix->eye($Am);
 
# $N and $us are all 0s, except:
# -  $N is the negative of the resistor values along the diagonal
# -  $us is has any voltage source values in the respective locations
my $N = Math::Matrix->new(@{zeros($Am,$Am)});
my $us = Math::Matrix->new(@{zeros(1,$Am)});
for (my $i = 0; $i < @components; ++$i)
{
    if ($components[$i][0] =~ /^V/)
    {
    $us->[$i][0] = $components[$i][3];
    }
    else
    {
    $N->[$i][$i] = -$components[$i][3];
    }
}
 
# $T is a 3x3 matrix composed of:
# /  0   0   A \
# | -A`  I   0 |
# \  0   M   N /
# Where I is the identity matrix, and 0 is a zero-filled matrix.
print "Making T...\n";
my $T = Math::Matrix->new(hmerge(@{zeros($An,$An)},         @{zeros($Am,$An)},        @{$A}),
              hmerge(@{$A->transpose->negative},@{Math::Matrix->eye($Am)},@{zeros($Am,$Am)}),
              hmerge(@{zeros($An,$Am)},         @{$M},                    @{$N}));
 
# $u is essentially $us extended to have the same number of rows as $T.
my $u = Math::Matrix->new(@{zeros(1,$An)},@{zeros(1,$Am)},@{$us});
 
print "Making w...\n";
# $w is the vector of node and branch voltages and branch currents. (This is what we're solving for.)
$w = $T->concat($u);
 
print "Solving...\n";
 
my @solution = @{$w = $w->solve};
 
# Print the equivalent resistance.
print abs($solution[$nodes][0]/$solution[$nodes+@components][0]) ."\n";

To Do

  1. Use a sparse matrix library. This will likely make things faster and use less memory!
  2. Support AC circuits
    • When a capacitor or inductor is seen, take the real value of the impedance as resistance.
    • Need to do a frequency sweep
  3. Support transistor circuits
    • Need to find/make a matrix library that will allow adding a variable instead of just a value (dependent voltage sources)
 
circuit_solver.txt · Last modified: 2010/03/13 05:14 by Ryan Fox
 
Except where otherwise noted, content on this wiki is licensed under the following license:CC Attribution-Share Alike 3.0 Unported
Recent changes RSS feed Donate Powered by PHP Valid XHTML 1.0 Valid CSS Driven by DokuWiki