Version 1 of a small operating-point circuit simulator

Updated 2008-05-28 23:52:24 by rvb

A small operating point circuit simulator in incr tcl rvb

The following simple circuit simulator class determines the operating point of a small circuit composed of linear and/or nonlinear elements. It uses element classes which are listed following the SCKT class. The available elements are voltage source, current source, resistor, diode, MOSFET, and controlled voltage and current sources. The simulator uses a modified nodal formulation, and a Newton algorithm.

The simulator can be used to determine the individual contributions of different modeling effects on a circuit's operating point. For example, bulk-bias sensitivity, channel-length shortening, and short-channel effects in a cascode current mirror (example 3 below) all throw off the simple hand-calculations using a simple MOSFET square-law model. To investigate these individual contributions, one simply changes the model parameters controlling the effect to a neutral value.

The element classes can be used to produce Current-Voltage characteristics, or other parameters, such as transconductance, versus operating point.

Following the class definitions for the elements is a Matrix class which is used in the solver. The Matrix class can "spill" its internal representation out to the la package representation, so the routines in that package can also be used for the solver.

Three utility routines follow the Matrix class definition.

Finally, there are three examples of the simple circuit operating point simulator.

References:

  1. Computer Methods for Circuit Analysis and Design, J. Vlach, K. Singhal, 1994, Van Nostrand Reinhold, NY.
  2. An Algorithm for Modified Nodal Analysis [L1 ]
  3. Linear Algebra Package la

Simple Operating Point Circuit Simulator Class

 ########################################################################
 # CLASS: SCKT.dcd
 # PURPOSE: small circuit simulator (operating point)
 # AUTHOR:  rvb
 # ----------------------------------------------------------------------
 # NOTES:
 ########################################################################
 class SCKT {
   public variable showlevel 0 ;# display results: 0 none, 1 no wait, 2 wait
   public variable nitl 500    ;# maximum number of iterations
   public variable vlim 1.0    ;# maximum dv update
   public variable ilim 1e-3   ;# maximum di update
   public variable vtol 1e-4   ;# voltage tolerance
   public variable itol 1e-8   ;# current tolerance
   public variable temp 25.0   ;# temperature
   public variable freq 0.0    ;# frequency
   private variable _elements  ;# element names
   private variable _nodes     ;# node names
   private variable _vsrcs     ;# vsrc names
   private variable _nnodes    ;# number of nodes
   private variable _nvsrcs    ;# number of vsrcs
   private variable _nmna      ;# dimension of MNA array
   private variable _A         ;# conductance matrix
   private variable _Y         ;# current matrix
   private variable _X         ;# dV, dI matrix
   private variable _V         ;# solution vector
   #=====================================================================
   # array SELEM_TYPES of SCKT element types, key = first letter of inst
   #=====================================================================
   common SELEM_TYPES
   array set SELEM_TYPES {
     i SISRC v SVSRC r SRES d SDIO m SMOS e SVCVS f SCCCS g SVCCS h SCCVS
   }
   #=====================================================================
   # METHOD: constructor
   # PURPOSE: read in circuit lines
   # NOTES :
   #  * line continuations not supported
   #  * i v r e f g h don't support PAR=VAL instance format
   #  * d m    require model and PAR=VAL instance format
   #  * each element is created anew: fails if any element already exists
   #=====================================================================
   constructor {cktlines args} {
     set _A  [code [Matrix #auto]]
     set _Y  [code [Matrix #auto]]
     set _X  [code [Matrix #auto]]
     set _V  [code [Matrix #auto]]
     set _elements {}
     set _nodes {}
     set _vsrcs {}
     set   Nodes(0) 1
     unset Nodes(0)
     foreach line [split $cktlines "\n"] {
       set line [string trim [string tolower $line]]
       if {[string length $line] == 0} {continue}
       set ctype [string range $line 0 0]
       if {[info exists SELEM_TYPES($ctype)]} {
         set type   [set SELEM_TYPES($ctype)]
         set name   [lshift line]
         set inodes {}
         set ipars  {}
         set model  ""
         foreach inode [set ::${type}::nodelist] {
           lappend inodes -$inode [lshift line]
         }
         switch -- $ctype {
           i {lappend ipars -I [lshift line]}
           v {lappend ipars -V [lshift line]}
           r {lappend ipars -R [lshift line]}
           e {lappend ipars -A [lshift line]}
           f {lappend ipars -VC [lshift line] -A [lshift line]}
           g {lappend ipars -A [lshift line]}
           h {lappend ipars -VC [lshift line] -A [lshift line]}
           d - m {
             set model  [lshift line]
             regsub -all -- "=" $line " " line
             foreach {par val} $line {
               lappend ipars -[string toupper $par] $val
             }
           }
         }
         lappend _elements $name
         set Element($name-type)    $type
         set Element($name-inodes)  $inodes
         set Element($name-ipars)   $ipars
         set Element($name-model)   $model
         foreach {node_config node} $inodes {
           if {$node ne "0"} {
             set Node($node) 1
           }
         }
         if {$type eq "SVSRC" || $type eq "SVCVS" || $type eq "SCCVS"} {
           lappend _vsrcs $name
         }
       } else {
         switch -regexp -- $line {
           {^[*#]} {
             set comment [string trim [string range $line 1 end]]
             puts "$comment"
           }
           {^.model} {
             set model [lshift line]
             set model [lshift line]
             regsub -all -- "=" $line " " line
             set mpars {}
             foreach {par val} $line {
               lappend mpars -[string toupper $par] $val
             }
             set Model($model) $mpars
           }
         }
       }
     }
     set _nodes [lsort [array names Node]]
     foreach name $_elements {
       set type   $Element($name-type)
       set inodes $Element($name-inodes)
       set ipars  $Element($name-ipars)
       set model  $Element($name-model)
       set mpars {}
       if {$model ne ""} {
         set mpars $Model($model)
       }
       eval $type $name $inodes $ipars $mpars
     }
     set _nnodes [llength $_nodes]
     set _nvsrcs [llength $_vsrcs]
     set _nmna [expr $_nnodes + $_nvsrcs]
     $_V zero $_nmna 1
     eval configure $args
   }
   #=====================================================================
   # METHOD: destructor
   # PURPOSE: delete constituents
   #=====================================================================
   destructor {
     delete object _A
     delete object _Y
     delete object _X
     delete object _V
     foreach element $_elements {
       delete object $element
     }
   }
   #=====================================================================
   # METHOD: node_index
   # PURPOSE: find array index for a node
   #=====================================================================
   method node_index {node} {
     if {$node eq "0"} {
       return -1
     }
     set nindex [lsearch $_nodes $node]
     if {$nindex < 0} {
       error "node \"$node\" not found"
     }
     return $nindex
   }
   #=====================================================================
   # METHOD: vsrc_index
   # PURPOSE: find array index for a vsrc
   #=====================================================================
   method vsrc_index {vsrc} {
     set vindex [lsearch $_vsrcs $vsrc]
     if {$vindex < 0} {
       error "vsrc \"$vsrc\" not found"
     }
     return [expr $vindex + $_nnodes]
   }
   #=====================================================================
   # METHOD: node_voltage
   # PURPOSE: find node voltage in solution vector
   #=====================================================================
   method node_voltage {node} {
     set index [node_index $node]
     if {$index < 0} {return 0.0}
     return [$_V get $index 0]
   }
   #=====================================================================
   # METHOD: vsrc_current
   # PURPOSE: find vsrc current in solution vector
   #=====================================================================
   method vsrc_current {vsrc} {
     set index [vsrc_index $vsrc]
     if {$index < 0} {return 0.0}
     return [$_V get $index 0]
   }
   #=====================================================================
   # METHOD: stamp
   # NOTES:
   #   [0] = [J][V] + [A B] [dV]
   #   [0]   [K][I] + [C 0] [dI]
   #=====================================================================
   method stamp {} {
     foreach element $_elements {
       foreach item [$element stamp [code $this]] {
         if {[regexp {@(.):} $item match block]} {continue}
         lvset {index value} [split $item :]
         switch -- $block {
           G {
             lvset {n1 n2} [split $index ,]
             set i1 [node_index $n1]
             set i2 [node_index $n2]
             if {$i1 < 0 || $i2 < 0} {continue}
             $_A +put $i1 $i2 $value
           }
           B {
             lvset {n1 v2} [split $index ,]
             set i1 [node_index $n1]
             set i2 [vsrc_index $v2]
             if {$i1 < 0} {continue}
             $_A +put $i1 $i2 $value
           }
           C {
             lvset {v1 n2} [split $index ,]
             set i1 [vsrc_index $v1]
             set i2 [node_index $n2]
             if {$i2 < 0} {continue}
             $_A +put $i1 $i2 $value
           }
           J {
             set i1 [node_index $index]
             if {$i1 < 0} {continue}
             $_Y +put $i1 0 $value
           }
           K {
             set i1 [vsrc_index $index]
             $_Y +put $i1 0 $value
           }
         }
       }
     }
   }
   #=====================================================================
   # METHOD: solve
   # PURPOSE: solve circuit operating point
   #=====================================================================
   method solve {} {
     puts [separator %]
     puts "% solving [namespace tail $this]"
     puts [separator %]
     #-------------------------------------------------------------------
     # initialize solution vector
     #-------------------------------------------------------------------
     $_V zero $_nmna 1
     set conv 0
     set iter 0
     while {1} {
       if {[incr iter] > $nitl} {
         puts " exceeded iteration limit"
         break
       }
       #-----------------------------------------------------------------
       # initialize, stamp conductance and current arrays, solve for delta
       #-----------------------------------------------------------------
       $_A zero $_nmna $_nmna
       $_Y zero $_nmna 1
       stamp
       $_X solve $_A $_Y
       #-----------------------------------------------------------------
       # develop and apply update
       #-----------------------------------------------------------------
       set dvmax 0.0
       set dimax 0.0
       for {set i 0} {$i < $_nmna} {incr i} {
         set q [expr abs([$_X get $i 0])]
         if {$i < $_nnodes} {
           if {$q > $dvmax} {
             set dvmax $q
           }
         } else {
           if {$q > $dimax} {
             set dimax $q
           }
         }
       }
       if {$dvmax > $vlim} {
         set ratio [expr $vlim/$dvmax]
       } else {
         set ratio 1.0
       }
       if {$dimax > $ilim} {
         set r [expr $ilim/$dimax]
         if {$r < $ratio} {set ratio $r}
       }
       set sumdv 0.0
       set sumdi 0.0
       for {set i 0} {$i < $_nmna} {incr i} {
         if {$i < $_nnodes} {
           set dv [$_X get $i 0]
           set dv [expr -$dv*$ratio]
           set sumdv [expr $sumdv + abs($dv)]
           $_V +put $i 0 $dv
         } else {
           set di [$_X get $i 0]
           set di [expr -$di*$ratio]
           set sumdi [expr $sumdi + abs($di)]
           $_V +put $i 0 $di
         }
       }
       if {$showlevel > 0} {
         showiter "interation = $iter"
       }
       #-----------------------------------------------------------------
       # convergence test
       #-----------------------------------------------------------------
       if {$sumdv < $vtol && $sumdi < $itol} {
         set conv 1
         puts " converged"
         break
       }
     }
     showiter "solution in $iter iterations" 1
   }
   #=====================================================================
   # METHOD: showiter
   # PURPOSE: display results at each iteration
   #=====================================================================
   method showiter {message {final 0}} {
     puts [separator =]
     puts "% $message"
     puts [separator =]
     if {$showlevel > 0} {
       puts "% conductance"
       $_A show
       puts "% current"
       $_Y show
       puts "% dV"
       $_X show
       puts "% solution"
     }
     foreach node $_nodes {
       puts "  V($node) = [node_voltage $node]"
     }
     foreach vsrc $_vsrcs {
       puts "  I($vsrc) = [vsrc_current $vsrc]"
     }
     if {!$final && $showlevel == 2} {
       puts "\n enter \"n\" to quit, or <enter> to continue"
       if {[string tolower [gets stdin]] eq "n"} {
         exit
       }
     }
   }
   #=====================================================================
   # METHOD: showelem
   # PURPOSE: display operating point results for all elements
   #=====================================================================
   method showelem {} {
     foreach element $_elements {
       $element show
     }
   }
 }

Simple Operating Point Circuit Simulator Element Classes

 #=======================================================================
 # CLASS: SELEM
 # PURPOSE: base class for SCKT elements
 #=======================================================================
 class SELEM {
   protected variable _Info

   method get {what} {
     if {[info exists _Info($what)]} {
       return $_Info($what)
     } else {
       puts " can't get parameter \$what\" "
       return ""
     }
   }
   method save {what value} {
     set _Info($what) $value
   }
   method show {} {
     set line {}
     set name [namespace tail $this]
     foreach what [lsort [array names _Info]] {
       lappend line $what=$_Info($what)
     }
     puts "  element \"$name\": [join $line]"
   }
 }


 #=======================================================================
 # CLASS: SISRC
 # PURPOSE: independent current source
 #=======================================================================
 class SISRC {
   inherit SELEM
   public variable POS 0
   public variable NEG 0
   common nodelist {POS NEG}
   public variable I 0.0

   constructor {args} {
     eval configure $args
   }

   method stamp {ckt} {
     set p_i $I
     set m_i [expr -$p_i]
     return "
       @J: $POS:$p_i $NEG:$m_i
     "
   }
 }


 #=======================================================================
 # CLASS: SVSRC
 # PURPOSE: independent voltage source
 #=======================================================================
 class SVSRC {
   inherit SELEM
   public variable POS 0
   public variable NEG 0
   common nodelist {POS NEG}
   public variable V 0.0

   constructor {args} {
     eval configure $args
   }

   method stamp {ckt} {
     set name [namespace tail $this]
     set vpos [$ckt node_voltage $POS]
     set vneg [$ckt node_voltage $NEG]
     set p_i  [$ckt vsrc_current $name]
     set m_i  [expr -$p_i]
     set dvx  [expr $vpos - $vneg - $V]
     return "
       @J: $POS:$p_i $NEG:$m_i
       @K: $name:$dvx
       @B: $POS,$name:1 $NEG,$name:-1
       @C: $name,$POS:1 $name,$NEG:-1
     "
   }
 }


 #=======================================================================
 # CLASS: SRES
 # PURPOSE: linear resistor
 #=======================================================================
 class SRES {
   inherit SELEM
   public variable POS 0
   public variable NEG 0
   common nodelist {POS NEG}
   public variable R 0.0 {
     if {$R < 1e-12 || $R > 1e12} {
       error "R is out of range"
     }
   }

   constructor {args} {
     eval configure $args
   }

   method evaluate {vpos vneg} {
     set  ir [expr ($vpos - $vneg)/$R]
     save ir $ir
     save gr [expr 1.0/$R]
     return $ir
   }

   method stamp {ckt} {
     set vpos [$ckt node_voltage $POS]
     set vneg [$ckt node_voltage $NEG]
     evaluate $vpos $vneg
     set p_ir [get ir]
     set p_gr [get gr]
     set m_ir [expr -$p_ir]
     set m_gr [expr -$p_gr]
     return "
       @J: $POS:$p_ir $NEG:$m_ir
       @G: $POS,$POS:$p_gr $POS,$NEG:$m_gr
           $NEG,$POS:$m_gr $NEG,$NEG:$p_gr
     "
   }
 }


 #=======================================================================
 # CLASS: SDIO
 # PURPOSE: diode
 #=======================================================================
 class SDIO {
   inherit SELEM
   public variable POS 0
   public variable NEG 0
   common nodelist {POS NEG}
   public variable A 1.0 {
     if {$A < 0.01 || $A > 1e4} {
       error "A is out of range"
     }
   }
   public variable IS 1e-15 {
     if {$IS < 1e-20 || $IS > 1e-3} {
       error "IS is out of range"
     }
   }
   public variable vdelta 0.01 {
     if {$vdelta < 1e-8 || $vdelta > 1e2} {
       error "vdelta is out of range"
     }
   }

   constructor {args} {
     eval configure $args
   }

   method evaluate {vpos vneg} {
     set VT 0.026
     set Id [expr $A*$IS*exp(($vpos - $vneg)/$VT)]
     save id $Id
     return $Id
   }

   method evaluate_conductances {vpos vneg} {
     set idd [evaluate [expr $vpos + $vdelta] $vneg]
     set id0 [evaluate $vpos $vneg]
     set gd  [expr ($idd - $id0)/$vdelta]
     save gd $gd
   }

   method stamp {ckt} {
     set vpos [$ckt node_voltage $POS]
     set vneg [$ckt node_voltage $NEG]
     evaluate_conductances $vpos $vneg
     set p_id [get id]
     set p_gd [get gd]
     set m_id [expr -$p_id]
     set m_gd [expr -$p_gd]
     return "
       @J: $POS:$p_id $NEG:$m_id
       @G: $POS,$POS:$p_gd $POS,$NEG:$m_gd
           $NEG,$POS:$m_gd $NEG,$NEG:$p_gd
     "
   }
 }


 #=======================================================================
 # CLASS: SMOS
 # PURPOSE: MOSFET
 #=======================================================================
 class SMOS {
   inherit SELEM
   public variable D 0
   public variable G 0
   public variable S 0
   public variable B 0
   common nodelist {D G S B}
   public variable CHAN 1.0 {
     if {$CHAN != 1.0 && $CHAN != -1.0} {
       error "CHAN is out of range"
     }
   }
   public variable W 1.0 {
     if {$W < 0.01 || $W > 1e4} {
       error "W is out of range"
     }
   }
   public variable L 1.0 {
     if {$L < 0.01 || $L > 1e4} {
       error "L is out of range"
     }
   }
   public variable KP 200.0 {
     if {$KP < 1.0 || $KP > 1e4} {
       error "KP is out of range"
     }
   }
   public variable VTH 0.5 {
     if {$VTH < 0.2 || $VTH > 2.0} {
       error "VTH is out of range"
     }
   }
   public variable PHI 0.7 {
     if {$PHI < 0.5 || $PHI > 1.2} {
       error "PHI is out of range"
     }
   }
   public variable GAMMA 0.7 {
     if {$GAMMA < 0.0 || $GAMMA > 2.0} {
       error "GAMMA is out of range"
     }
   }
   public variable LAMBDA 0.01 {
     if {$LAMBDA < 0.0 || $LAMBDA > 0.5} {
       error "LAMBDA is out of range"
     }
   }
   public variable THETAS 0.01 {
     if {$THETAS < 0.0 || $THETAS > 5.0} {
       error "THETAS is out of range"
     }
   }
   public variable THETAC 0.10 {
     if {$THETAC < 0.0 || $THETAC > 5.0} {
       error "THETAC is out of range"
     }
   }
   public variable ETA 0.0 {
     if {$ETA < 0.0 || $ETA > 0.2} {
       error "ETA is out of range"
     }
   }
   public variable SST 80.0 {
     if {$SST < 20.0 || $SST > 200.0} {
       error "SST is out of range"
     }
   }
   public variable vdelta 0.01 {
     if {$vdelta < 1e-8 || $vdelta > 1e2} {
       error "vdelta is out of range"
     }
   }

   constructor {args} {
     eval configure $args
   }

   method evaluate {vd vg vs vb} {
     set Vbs [expr $CHAN*($vb-$vs)]
     set Vgs [expr $CHAN*($vg-$vs)]
     set Vds [expr $CHAN*($vd-$vs)]
     set pol $CHAN
     if {$Vds < 0.0} {
       set pol [expr -$CHAN]
       set Vbs [expr $CHAN*($vb-$vd)]
       set Vds [expr $CHAN*($vs-$vd)]
       set Vgs [expr $CHAN*($vg-$vd)]
     }
     set arg  [expr $PHI-$Vbs]
     set arg  [expr ($arg<0.0)?0.0:$arg]
     set rsp  [expr sqrt($arg)]
     set abc  [expr 1.0+$GAMMA/(2.0*$rsp)]
     set Vsth [expr $SST/1151.0]
     set Vthe [expr $VTH+$GAMMA*($rsp-sqrt($PHI))-$ETA*$Vds]
     set Vgse [expr $Vgs-$Vthe]
     set Vgsl [expr $Vsth*log(1+exp($Vgse/$Vsth))]
     set Vdss [expr $Vgsl/($abc+$THETAC*$Vgsl/2.0)]
     set Vdsl [expr ($Vds < $Vdss)?$Vds : $Vdss]
     set advd [expr abs($Vds-$Vdsl)]
     set avdl [expr abs($Vdsl)]
     set vmob [expr 1.0+$THETAS*($Vgsl+2.0*$GAMMA*$rsp)]
     set hmob [expr 1.0+$THETAC*$avdl]
     set cls  [expr 1.0+$LAMBDA*$advd]
     set beta [expr $W/$L*$KP*1e-6]
     set Ids  [expr $pol*$beta*$cls/($vmob*$hmob)*($Vgsl-$abc*$avdl/2.0)*$Vdsl]
     save id   $Ids
     save vth  $Vthe
     save vsat $Vdss
     return $Ids
   }

   method evaluate_conductances {vd vg vs vb} {
     set idg  [evaluate $vd [expr $vg + $vdelta] $vs $vb]
     set idd  [evaluate [expr $vd + $vdelta] $vg $vs $vb]
     set idb  [evaluate $vd $vg $vs [expr $vb + $vdelta]]
     set id0  [evaluate $vd $vg $vs $vb]
     save gm  [expr ($idg - $id0)/$vdelta]
     save gds [expr ($idd - $id0)/$vdelta]
     save gmb [expr ($idb - $id0)/$vdelta]
   }

   method stamp {ckt} {
     set vd [$ckt node_voltage $D]
     set vg [$ckt node_voltage $G]
     set vs [$ckt node_voltage $S]
     set vb [$ckt node_voltage $B]
     evaluate_conductances $vd $vg $vs $vb
     set p_id  [get id ]
     set p_gm  [get gm ]
     set p_gds [get gds]
     set p_gmb [get gmb]
     set m_id  [expr -$p_id ]
     set m_gm  [expr -$p_gm ]
     set m_gds [expr -$p_gds]
     set m_gmb [expr -$p_gmb]
     return "
       @J: $D:$p_id $S:$m_id
       @G: $D,$G:$p_gm  $D,$S:$m_gm  $S,$G:$m_gm  $S,$S:$p_gm
           $D,$D:$p_gds $D,$S:$m_gds $S,$D:$m_gds $S,$S:$p_gds
           $D,$B:$p_gmb $D,$S:$m_gmb $S,$B:$m_gmb $S,$S:$p_gmb
     "
   }
 }
 #=======================================================================
 # CLASS: SVCCS
 # PURPOSE: voltage-controlled current source
 #=======================================================================
 class SVCCS {
   inherit SELEM
   public variable POS 0
   public variable NEG 0
   public variable CPOS 0
   public variable CNEG 0
   common nodelist {POS NEG CPOS CNEG}
   public variable A 1.0

   constructor {args} {
     eval configure $args
   }

   method stamp {ckt} {
     set vcpos [$ckt node_voltage $CPOS]
     set vcneg [$ckt node_voltage $CNEG]
     set p_a  $A
     set m_a  [expr -$p_a]
     set p_i  [expr $A*($vcpos - $vcneg)]
     set m_i  [expr -$p_i]
     return "
       @J: $POS:$p_i $NEG:$m_i
       @G: $POS,$CPOS:$p_a $POS,$CNEG:$m_a $NEG,$CPOS:$m_a $NEG,$CNEG:$p_a
     "
   }
 }


 #=======================================================================
 # CLASS: SVCVS
 # PURPOSE: voltage-controlled voltage source
 #=======================================================================
 class SVCVS {
   inherit SELEM
   public variable POS 0
   public variable NEG 0
   public variable CPOS 0
   public variable CNEG 0
   common nodelist {POS NEG CPOS CNEG}
   public variable A 1.0

   constructor {args} {
     eval configure $args
   }

   method stamp {ckt} {
     set name [namespace tail $this]
     set vpos  [$ckt node_voltage $POS]
     set vneg  [$ckt node_voltage $NEG]
     set vcpos [$ckt node_voltage $CPOS]
     set vcneg [$ckt node_voltage $CNEG]
     set p_i   [$ckt vsrc_current $name]
     set m_i   [expr -$p_i]
     set p_a   $A
     set m_a   [expr -$A]
     set dvx   [expr $vpos - $vneg - $A*($vcpos - $vcneg)]
     return "
       @J: $POS:$p_i $NEG:$m_i
       @K: $name:$dvx
       @B: $POS,$name:1 $NEG,$name:-1
       @C: $name,$POS:1 $name,$NEG:-1
           $name,$CPOS:$m_a $name,$CNEG:$p_a
     "
   }
 }

Matrix class

 ########################################################################
 # Class: Matrix
 # PURPOSE: matrix operations
 # AUTHOR:  rvb
 # ----------------------------------------------------------------------
 # NOTES:
 #  * baseindex = 0 or 1 for matrix indices from 0,0 or 1,1
 ########################################################################
 class Matrix {
   public variable baseindex 0
   private variable _nrow
   private variable _ncol
   private variable _mxdata
   private variable _maxindex
   #=====================================================================
   # constructor
   #=====================================================================
   constructor {args} {
     eval configure $args
   }
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   # access methods
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   #=====================================================================
   # METHOD: nrows
   # RETURN: number of matrix rows
   #=====================================================================
   method nrows {} {
     return $_nrow
   }
   #=====================================================================
   # METHOD: ncols
   # RETURN: number of matrix columns
   #=====================================================================
   method ncols {} {
     return $_ncol
   }
   #=====================================================================
   # METHOD: spill
   # PURPOSE: spill matrix representation for copying
   # RETURN: LA package-compatible list
   # NOTES:
   #   * LA package representation:
   #     {2 <nrows> <ncols> <data>}
   #=====================================================================
   method spill {} {
     return [concat 2 $_nrow $_ncol $_mxdata]
   }
   #=====================================================================
   # METHOD: fill
   # PURPOSE: fill matrix representation for copying
   # ARGUMENTS:
   #   % LA package-compatible list
   # NOTES:
   #   * LA package representation:
   #     {2 <nrows> <ncols> <data>}
   #=====================================================================
   method fill {mxrep} {
     set type [lindex $mxrep 0]
     set nrow [lindex $mxrep 1]
     set ncol [lindex $mxrep 2]
     _dimension $nrow $ncol
     set _mxdata [lrange $mxrep 3 end]
   }
   #=====================================================================
   # METHOD: get
   # PURPOSE: get a matrix data value at index i, j
   #=====================================================================
   method get {i j} {
     set k [_locate $i $j]
     return [lindex $_mxdata $k]
   }
   #=====================================================================
   # METHOD: put
   # PURPOSE: put a matrix data value into index i, j
   #=====================================================================
   method put {i j value} {
     set k [_locate $i $j]
     set _mxdata [lreplace $_mxdata $k $k $value]
   }
   #=====================================================================
   # METHOD: +put
   # PURPOSE: increment a matrix data value at index i, j by value
   #=====================================================================
   method +put {i j value} {
     put $i $j [expr $value + [get $i $j]]
   }
   #=====================================================================
   # METHOD: *put
   # PURPOSE: scale a matrix data value at index i, j by value
   #=====================================================================
   method *put {i j value} {
     put $i $j [expr $value * [get $i $j]]
   }
   #=====================================================================
   # METHOD: show
   # PURPOSE: display matrix values
   #=====================================================================
   method show {} {
     set nr [nrows]
     set nc [ncols]
     puts "[namespace tail $this] ($nr x $nc) ="
     for {set i 0} {$i < $nr} {incr i} {
       set line {}
       for {set j 0} {$j < $nc} {incr j} {
         lappend line [format "%9.3g" [get $i $j]]
       }
       puts "  \[ [join $line] \]"
     }
   }
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   # basic
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   #=====================================================================
   # METHOD: _dimension (private)
   # PURPOSE: set matrix size and initialize matrix data to 0.0
   #=====================================================================
   private method _dimension {nrow ncol} {
     set _nrow $nrow
     set _ncol $ncol
     set _maxindex [expr $nrow*$ncol - 1]
     set _mxdata {}
     for {set i 0} {$i <= $_maxindex} {incr i} {
       lappend _mxdata 0.0
     }
   }
   #=====================================================================
   # METHOD: _locate (private)
   # PURPOSE: find index into matrix data array using matrix indices
   # ARGUMENTS:
   #   % i j
   #       matrix indices for row and column
   # NOTES:
   #   * matrix indices are
   #       0 to (nrow - 1) for baseindex = 0, or
   #       1 to nrow       for baseindex = 1
   #=====================================================================
   private method _locate {i j} {
     set i [expr int($i)]
     set j [expr int($j)]
     set index [expr ($i - $baseindex)*$_ncol + ($j - $baseindex)]
     if {$index > $_maxindex} {
       fatal "_locate" \
         "indices \"$i,$j\" are out of range" \
         "index = $index maxindex = $_maxindex"
     }
     return $index
   }
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   # specific matrices
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   #=====================================================================
   # METHOD: zero
   # PURPOSE: set array to zero matrix (all zero elements)
   #=====================================================================
   method zero {nrow ncol} {
     _dimension $nrow $ncol
   }
   #=====================================================================
   # METHOD: identity
   # PURPOSE: set array to identity matrix
   #=====================================================================
   method identity {nrow ncol} {
     if {$ncol != $nrow} {
       fatal "identity" \
         "nrow = \"$nrow\" ncol = \"$ncol\"" \
         "can't make identity of non-square matrix"
     }
     _dimension $nrow $ncol
     for {set i 0} {$i < $nrow} {incr i} {
       put $i $i 1.0
     }
   }
   #=====================================================================
   # METHOD: :=
   # PURPOSE: set matrix with {0 2 3; 1 2 3} spec, or to another mx
   # ARGUMENTS:
   #   % matrixspec
   #       list of rows separated by newline or ; 
   #     or: another matrix name
   # NOTES:
   #   * skips empty rows
   #=====================================================================
   method := {matrixspec} {
     if {[itcl::is object -class Matrix $matrixspec]} {
       set nar [$matrixspec nrows]
       set nac [$matrixspec ncols]
       zero $nar $nac
       fill [$matrixspec spill]
       return
     }
     set nrow 0
     set ncol 0
     set values {}
     foreach row [split $matrixspec ";\n"] {
       set ncoli [llength $row]
       if {$ncoli > 0} {
         incr nrow
         if {$nrow == 1} {
           set ncol $ncoli
         }
         if {$ncoli != $ncol} {
           set irow [expr $nrow + $baseindex - 1]
           fatal ":=" \
             "number of values in row \"$irow\"" \
             "  is different from previous rows"
         }
         set values [concat $values $row]
       }
     }
     _dimension $nrow $ncol
     set _mxdata $values
   }
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   # matrix operations
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   #=====================================================================
   # METHOD: +
   # PURPOSE: add one matrix to another
   #=====================================================================
   method + {A B} {
     set nar [$A nrows]
     set nac [$A ncols]
     set nbr [$B nrows]
     set nbc [$B ncols]
     if {($nar != $nbr) || ($nac != $nbc)} {
       fatal "+" \
         "matrix \"$A\" has \"$nar\" rows  and \"$nac\" cols" \
         "matrix \"$B\" has \"$nbr\" rows  and \"$nbc\" cols" \
         "cannot add $A to $B"
     }
     set bia [$A cget -baseindex]
     set bib [$B cget -baseindex]
     $A configure -baseindex 0
     $B configure -baseindex 0
     set mx [Matrix #auto]
     $mx zero $nar $nac
     for {set i 0} {$i < $nar} {incr i} {
       for {set j 0} {$j < $nac} {incr j} {
         $mx put $i $j [expr [$A get $i $j] + [$B get $i $j]]
       }
     }
     := $mx
     delete object $mx
     $A configure -baseindex $bia
     $B configure -baseindex $bib
   }
   #=====================================================================
   # METHOD: *
   # PURPOSE: multiply one matrix by another
   #=====================================================================
   method * {A B} {
     set nar [$A nrows]
     set nac [$A ncols]
     set nbr [$B nrows]
     set nbc [$B ncols]
     if {$nac != $nbr} {
       fatal "*" \
         "matrix \"$A\" has \"$nac\" columns" \
         "matrix \"$B\" has \"$nbr\" rows" \
         "cannot multiply $A with $B"
     }
     set bia [$A cget -baseindex]
     set bib [$B cget -baseindex]
     $A configure -baseindex 0
     $B configure -baseindex 0
     set mx [Matrix #auto]
     $mx zero $nar $nbc
     for {set i 0} {$i < $nar} {incr i} {
       for {set k 0} {$k < $nbc} {incr k} {
         set c 0.0
         for {set j 0} {$j < $nac} {incr j} {
           set c [expr $c + [$A get $i $j] * [$B get $j $k]]
         }
         $mx put $i $k $c
       }
     }
     := $mx
     delete object $mx
     $A configure -baseindex $bia
     $B configure -baseindex $bib
   }
   #=====================================================================
   # METHOD: scale
   # PURPOSE: scale matrix by value
   #=====================================================================
   method scale {value} {
     for {set i 0} {$i < [nrows]} {incr i} {
       for {set j 0} {$j < [ncols]} {incr j} {
         *put $i $j $value
       }
     }
   }
   #=====================================================================
   # METHOD: transpose
   # PURPOSE: transpose a matrix
   #=====================================================================
   method transpose {A} {
     set nar [$A nrows]
     set nac [$A ncols]
     set bia [$A cget -baseindex]
     $A configure -baseindex 0
     set mx [Matrix #auto]
     $mx zero $nac $nar
     for {set i 0} {$i < $nar} {incr i} {
       for {set j 0} {$j < $nac} {incr j} {
         $mx put $j $i [$A get $i $j]
       }
     }
     := $mx
     delete object $mx
     $A configure -baseindex $bia
   }
   #=====================================================================
   # METHOD: invert
   # PURPOSE: invert a matrix and return the determinant
   #=====================================================================
   method invert {A} {
     set nar [$A nrows]
     set nac [$A ncols]
     if {$nar != $nac} {
       fatal "invert" \
         "matrix \"$A\" has \"$nar\" rows  and \"$nac\" cols" \
         "can't invert non-square matrix"
     }
     set bia [$A cget -baseindex]
     $A configure -baseindex 0
     set mx [Matrix #auto]
     set nuc [expr $nar*2]
     $mx zero $nar $nuc
     set n $nar
     #-------------------------------------------------------------------
     # fill working array with A:I
     #-------------------------------------------------------------------
     for {set i 0} {$i < $n} {incr i} {
       for {set j 0} {$j < $n} {incr j} {
         $mx put $i $j [$A get $i $j]
       }
       $mx put $i [expr $i+$n] 1.0
     }
     #-------------------------------------------------------------------
     # invert:
     #   1. find biggest element in column for pivot
     #   2. interchange rows with pivot row
     #-------------------------------------------------------------------
     set det 1.0
     for {set i 0} {$i < $n} {incr i} {
       set pivot 0.0
       set ipivot $i
       for {set j $i} {$j < $n} {incr j} {
         set q [$mx get $j $i]
         if {abs($q) > abs($pivot)} {
           set pivot $q
           set ipivot $j
         }
       }
       set det [expr $det*$pivot]
       if {$ipivot != $i} {
         set det [expr -$det]
       }
       if {$det == 0.0} {
         fatal "invert" "\"$A\" is a singular matrix"
       }
       for {set j $i} {$j < $nuc} {incr j} {
         set d [$mx get $ipivot $j]
         if {$ipivot != $i} {
           $mx put $ipivot $j [$mx get $i $j]
         }
         $mx put $i $j [expr 1.0*$d/$pivot]
       }
       for {set j 0} {$j < $n} {incr j} {
         set d [$mx get $j $i]
         if {($j != $i) && ($d != 0.0)} {
           for {set k $i} {$k < $nuc} {incr k} {
             $mx put $j $k [expr [$mx get $j $k] - $d*[$mx get $i $k]]
           }
         }
       }
     }
     #-------------------------------------------------------------------
     # extract inverse from working array
     #-------------------------------------------------------------------
     zero $nar $nar
     for {set i 0} {$i < $n} {incr i} {
       for {set j 0} {$j < $n} {incr j} {
         put $i $j [$mx get $i [expr $j+$n]]
       }
     }
     delete object $mx
     $A configure -baseindex $bia
     return $det
   }
   #=====================================================================
   # METHOD: solve
   # PURPOSE: solve linear system B=A*X
   #=====================================================================
   method solve {A B} {
     set nar [$A nrows]
     set nac [$A ncols]
     set nbr [$B nrows]
     set nbc [$B ncols]
     if {$nar != $nac} {
       fatal "invert" \
         "matrix \"$A\" has \"$nar\" rows  and \"$nac\" cols" \
         "cannot solve $B = $A * X"
     }
     if {$nac != $nbr} {
       fatal "solve" \
         "matrix \"$A\" has \"$nar\" rows  and \"$nac\" cols" \
         "matrix \"$B\" has \"$nbr\" rows  and \"$nbc\" cols" \
         "cannot solve $B = $A * X"
     }
     set bia [$A cget -baseindex]
     set bib [$B cget -baseindex]
     $A configure -baseindex 0
     $B configure -baseindex 0
     set mx [Matrix #auto]
     $mx zero $nar $nac
     zero $nbr $nbc
     #-------------------------------------------------------------------
     # decomposition
     #-------------------------------------------------------------------
     $mx := $A
     set nm1 [expr $nar-1]
     set ipivots {}
     for {set k 0} {$k < $nm1} {incr k} {
       set kp1 [expr $k+1]
       set ipivot $k    
       set apivot [expr abs([$mx get $ipivot $k])]
       for {set i $kp1} {$i < $nar} {incr i} {
         set q [expr abs([$mx get $i $k])]
         if {$q > $apivot} {
           set ipivot $i
           set apivot $q
         }
       }
       lappend ipivots $ipivot
       set t [$mx get $ipivot $k]
       if {$ipivot != $k} {
         $mx put $ipivot $k [$mx get $k $k]
         $mx put $k $k $t
       }
       if {$t == 0.0} {
         fatal "solve" "\"$A\" is a singular matrix"
       }
       for {set i $kp1} {$i < $nar} {incr i} {
         $mx put $i $k [expr -1.0*[$mx get $i $k]/$t]
       }
       for {set j $kp1} {$j < $nar} {incr j} {
         set t [$mx get $ipivot $j]
         if {$ipivot != $k} {
           $mx put $ipivot $j [$mx get $k $j]
           $mx put $k $j $t
         }
         for {set i $kp1} {$i < $nar} {incr i} {
           $mx +put $i $j [expr [$mx get $i $k]*$t]
         }
       }
     }
     #-------------------------------------------------------------------
     # back-substitution:
     #-------------------------------------------------------------------
     := $B
     for {set bcol 0} {$bcol < $nbc} {incr bcol} {
       for {set k 0} {$k < $nm1} {incr k} {
         set kp1 [expr $k + 1]
         set ipivot [lindex $ipivots $k]
         set t [get $ipivot $bcol]
         if {$k != $ipivot} {
           put $ipivot $bcol [get $k $bcol]
           put $k $bcol $t
         }
         for {set i $kp1} {$i < $nar} {incr i} {
           +put $i $bcol [expr [$mx get $i $k]*$t]
         }
       }
       for {set kb 0} {$kb < $nm1} {incr kb} {
         set k [expr $nm1 - $kb]
         set a [$mx get $k $k]
         set b [get $k $bcol]
         set t [expr 1.0*$b/$a]
         put $k $bcol $t
         for {set i 0} {$i < $k} {incr i} {
           +put $i $bcol [expr -[$mx get $i $k]*$t]
         }
       }
       put 0 $bcol [expr 1.0*[get 0 $bcol] / [$mx get 0 0]]
     }
     delete object $mx
     $A configure -baseindex $bia
     $B configure -baseindex $bib
   }
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   # utility
   #%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   #=====================================================================
   # METHOD: Matrix::warning
   # PURPOSE : print warning message
   #=====================================================================
   method warning {method args} {
     set str [join $args " "]
     set name [namespace tail $this]
     puts "Matrix::$method $name warning: $str\n"
   }
   #=====================================================================
   # METHOD: Matrix::fatal
   # PURPOSE : print error message, then die
   #=====================================================================
   method fatal {method args} {
     set sep [separator % 72]
     set str [join $args "\n  "]
     set name [namespace tail $this]
     puts "\n$sep\nMatrix::$method $name error:\n  $str\n$sep\n"
     exit -1
   }
 }

Utility routines used in the class definitions

 #=======================================================================
 # PROC    : lshift
 # PURPOSE : shift list and return first element
 # AUTHOR  : rvb
 # ----------------------------------------------------------------------
 # ARGUMENTS :
 #   % inputlist
 #       List to be shifted.
 # RESULTS :
 #   * Sets inputlist to 2nd to last elements of original inputlist
 #   * Returns first element in inputlist
 # EXAMPLE-CALL :
 #{
 #  while {[llength $argv] > 0} {
 #    set arg [lshift argv]
 #    switch -- $arg {
 #      -lib  {set lib [lshift argv]}
 #      -show {set show 1}
 #      default {lappend tests $arg}
 #    }
 #  }
 #}
 #=======================================================================
 proc lshift {inputlist} {
   upvar $inputlist argv
   set arg  [lindex $argv 0]
   set argv [lrange $argv 1 end]
   return $arg
 }

 #=======================================================================
 # PROC    : lvset
 # PURPOSE : set list of variables in first list to values in second list
 # AUTHOR  : rvb
 # ----------------------------------------------------------------------
 # ARGUMENTS :
 #   % varlist
 #       List of variables to be assigned values
 #   % valuelist
 #       List of values to be assigned
 # RESULTS :
 #   * Sets variables in varlist to values in vallist
 #   * Error if lengths of the two lists aren't the same
 # EXAMPLE-CALL :
 #{
 #  lvset {a b c} {1.2 "test 1" 3}
 #}
 #=======================================================================
 proc lvset {varlist vallist} {
   if {[llength $varlist] != [llength $vallist]} {
     error "number of values is different from number of variables"
   }
   foreach var $varlist val $vallist {
     uplevel "set $var {$val}"
   }
 }

 #=============================================================================
 # PROC    : separator
 # PURPOSE : make a string of characters for a separator
 # AUTHOR  : rvb
 # ---------------------------------------------------------------------------
 # ARGUMENTS :
 #   % char
 #       character to be repeated nchars times
 #   % nchars
 #       number of characters in separator
 #       (default is 78)
 # RESULTS :
 #   * Returns separator string, composed of nchars characters.
 #=============================================================================
 proc separator {char {nchars 78}} {
   set s "";for {set i 0} {$i<$nchars} {incr i} {set s $s$char}
   return $s
 }

Example 1: voltage source and resistor

 SCKT ckt {
   r1 1 0 2.0
   v1 1 0 1.2
 } -showlevel 0 -vlim 10 -ilim 10
 ckt solve

Results from Example 1:

  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  % solving ckt
  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   converged
  ==========================================================================
  % solution in 2 iterations
  ==========================================================================
    V(1) = 1.2
    I(v1) = -0.6

Example 2: voltage source and diode

 SCKT ckt {
   v1 1 0 0.6
   d1 1 0 ndio
   .model ndio IS=1e-12
 }
 ckt solve

Results from Example 2:

  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  % solving ckt
  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   converged
  ==========================================================================
  % solution in 13 iterations
  ==========================================================================
    V(1) = 0.6
    I(v1) = -0.0105239881794353

Example 3: cascode current mirror

 SCKT ckt {
   mdd A M 0 0 nsmos W=4.0 L=1.0
   mxd M X A A nsmos W=4.0 L=1.0
   mxx X X 0 0 nsmos W=1.0 L=1.0
   mmz B M 0 0 nsmos W=4.0 L=1.0
   mxz Z X B B nsmos W=4.0 L=1.0
   vz  Z 0 1.2
   im  0 M 50e-6
   ix  0 X 50e-6
   .model nsmos  KP=250 VTH=0.4 GAMMA=0 LAMBDA=0.01
 } -showlevel 0 -vlim 10
 ckt solve
 ckt showelem

Results from Example 3:

  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  % solving ckt
  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   converged
  ==========================================================================
  % solution in 10 iterations
  ==========================================================================
    V(a) = 0.33227721920471
    V(b) = 0.333063880432942
    V(m) = 0.721044469558005
    V(x) = 1.05322837111299
    V(z) = 1.2
    I(vz) = -5.00003932690847e-05
    element "mdd": gds=4.99922454479687e-07 gm=0.0003066863882724 gmb=0.0 id=5.00000622114379e-05 vsat=0.316633336264721 vth=0.4
    element "mxd": gds=4.99665374840399e-07 gm=0.00030678184656523 gmb=0.0 id=5.00026251046766e-05 vsat=0.316551820324572 vth=0.4
    element "mxx": gds=4.97906443390473e-07 gm=0.00014855907832262 gmb=0.0 id=5.00000912727292e-05 vsat=0.632573783883089 vth=0.4
    element "mmz": gds=4.99922454489851e-07 gm=0.00030668883219441 gmb=0.0 id=5.00004599916131e-05 vsat=0.316633336264721 vth=0.4
    element "mxz": gds=4.97261814260259e-07 gm=0.00030751525675164 gmb=0.0 id=5.00002454620565e-05 vsat=0.315788673803137 vth=0.4
    element "vz": 
    element "im": 
    element "ix": 

See also Playing with circuits