Code covered by the BSD License

### Highlights fromCircuit Analysis Toolbox

from Circuit Analysis Toolbox by John McDermid
The circuit analysis toolbox allows you to perform an AC analysis of a circuit.

tableau
classdef tableau
% Creates an object for the simplified tableau format.
% The tableau matrix is system of equations contrived from two matrix
% relations. The first is the relationship between elemental values
% described by the simple equation ki*Ib-kv*Vb=s. Ib is the current
% through component and Vb is the voltage across the component. The
% second is A*IB=0 where A is the reduced incidence matrix and IB is a
% vector of branch currents. A*IB=0 is a matrix statement of Kirkoff's
% current law for all nodes.
%
% Values of ki, kv, and s are picked to match each component type.
%       Resistor
%           ki=R,  kv=1,  s=0,   making R=Vb/Ib
%       Inductor
%           ki=j*2*pi*f*L,  kv=1,  s=0,  making Vb/Ib=j*2*pi*f*L=j*XL
%       Capacitor
%           ki=1,  kv=j*2*pi*f*C,  s=0  making Vb/Ib=1/(j*2*pi*f*C)=j*XC
%       Independent Current Source
%           ki=1,  kv=0,  s=I   making Ib=I
%       Independent Voltage Source
%           ki=0,  kv=1,  s=V   making -Vb=V (note that the current in
%                                             the voltage source flows in the
%                                             opposite direction to the
%                                             voltage accross it and
%                                             hence the minus sign)
%       Current Dependent Current Source
%           ki=[1 -Beta],  kv=[0 0],  s=0, making [1 -Beta)*[Ib
%           Ib_cntrl]'=0
%                                       or Ib=Beta*Ib_cntrl
%       Voltage Dependent Current Source
%           ki=1,  kv=gm,  s=0  making Ib-gm*Vb_cntrl
%                                     or Ib=gm*Vb_cntrl
%       Voltage Dependent Voltage Source
%           ki=[0 0],  kv=[1 Av],  s=0,  -[1 -Av]*[Vb VB_cntrl]'=0
%                                            or -Vb=Av*Vb_cntrl
%                                             (note that the current in
%                                             the voltage source flows in the
%                                             opposite direction to the
%                                             voltage accross it and
%                                             hence the minus sign)
%
% Kv and Ki (constructed from kv and ki) are essentially diagonal
% matricies except for the dependent sources.
%
% The tableau matrix is constructed in the following way:
%       tableau=[Kv -Ki*A';A Z] where Z is a matrix of zeros and of proper
%                                dimension to make tableau a square matrix
% The S column vector is constructed from the individual s and 0
%
% PROPERTIES:
%         ckt                Structure of circuit properties
%         typesAvailable     A list of the types available in this circuit
%         typeIndx           An index to the list of types indexed by type name
%         desigIndx          An index to the designators in A, Ki, and Kv by designator name
%         nodeIndx           An index to nodes in Aa by designator name
%         VnNodeIndx         An index to nodes in A by designator name
%         TmIndx             An index to all elements in the tableau matrix
%                              including designators and nodes by name
%         desigList          The list of designators in matrix order
%         nodeList           The list of nodes in matrix order
%         branchIndx         Index to all branch currents in tableau matrix
%         nodeVoltIndx       Index to all node voltages in the tableau matrix
%         Aa                 The complete incidence matrix
%         accInAa            A logical vector indicating which nodes in Aa have access
%         accInA             A logical vector indicating which nodes in A have access
%         A                  The reduced incidence matrix with the reference node removed (Sparse)
%         Ki                 A sparse matrix holding impedance values
%         Kv                 A sparse matrix holding conductance values
%         Z                  A sparse array of zeros
%         S                  A sparse vector of source voltages and currents
%         deltaKi            A sparse matrix containing a "1" at each value location in Ki
%         deltaKv            A sparse matrix containing a "1" at each value location in Kv
%         deltaS             A sparse vector containing a "1" at each value location in S
%
%
% EXAMPLE:
%        tableauObject = tableau(circuitObject)

% AUTHOR:
%       John McDermid
% CREATED:
%       October 19, 2009

properties (SetAccess = protected)
% tableau Properties:
ckt               % Structure of circuit properties
typesAvailable    % A list of the types availble in this circuit
typeIndx          % An index to the list of types indexed by type name
desigIndx         % An index to the designators in A, Ki, and Kv by designator name
nodeIndx          % An index to nodes in Aa by designator name
VnNodeIndx        % An index to nodes in A by designator name
TmIndx            % An index to all elements in the tableau matrix
%   including designators and nodes by name
desigList         % The list of designators in matrix order
nodeList          % The list of nodes in matrix order
branchIndx        % Index to all branch currents in tableau matrix
nodeVoltIndx      % Index to all node voltages in the tableau matrix
accInAa           % A logical vector indicating which nodes in Aa have access
accInA            % A logical vector indicating which nodes in A have access
Aa                % The complete incidence matrix (Sparse)
A                 % The reduced incidence matrix with the reference node removed (Sparse)
Ki                % A sparse matrix holding impedance values
Kv                % A sparse matrix holding conductance values
Z                 % A sparse array of zeros
S                 % A sparse vector of source voltages and currents
deltaKi           % A sparse matrix containing a "1" at each value location in Ki
deltaKv           % A sparse matrix containing a "1" at each value location in Kv
deltaS            % A sparse vector containing a "1" at each value location in S
end

methods
function obj=tableau(cirObj)
% Constructor for a tableau object.
% Requires a circuit object as input.

% Find the index of types in the cirDesc structure so that
% components can be grouped according to type and place them in
% a container so that the indexof types can be retrieved by
% type name
uniqueTypes=unique({cirObj.cirDesc.type});
M=numel(uniqueTypes);            % Number of unique component types in the circuit
typeIndexCirObj{M}=[];           % Pre-allocate for speed
N=numel({cirObj.cirDesc.type});  % Number of components in the circuit
for J=1:M                        % Assemble the component index for each unique type
row=[];
for I=1:N
if cirObj.cirDesc(I).type==uniqueTypes{J}
row=[row I];   %#ok<AGROW>
end
end
typeIndexCirObj{J}=row;
end
types=containers.Map(uniqueTypes,typeIndexCirObj);
% Build the ckt property (a structure) from the unique types
for I=1:M
typeName=char(uniqueTypes(I));
obj.ckt.(typeName).designator={cirObj.cirDesc(types(typeName)).designator};
obj.ckt.(typeName).fromNode={cirObj.cirDesc(types(typeName)).fromNode};
obj.ckt.(typeName).toNode={cirObj.cirDesc(types(typeName)).toNode};
obj.ckt.(typeName).value=[cirObj.cirDesc(types(typeName)).value];
obj.ckt.(typeName).plusTol=[cirObj.cirDesc(types(typeName)).plusTol];
obj.ckt.(typeName).negTol=[cirObj.cirDesc(types(typeName)).negTol];
obj.ckt.(typeName).cntrlBranch={cirObj.cirDesc(types(typeName)).cntrlBranch};
end
% Find the order of unique types that matches the order of the type names
%  available
typeNamesAvailable={'V' 'I' 'R' 'L' 'C' 'IDI' 'VDI' 'VDV'};
currentOrder(M)=0;
% Created an order for the unique types frm the order in type
% names available
for I=1:M
currentOrder(I)=find(strcmp(uniqueTypes(I),typeNamesAvailable));
end
[~,typeOrder]=sort(currentOrder);
% Make a container for designators that indexes the branch
% dimension in the A matrix. Construct it in the order defined
% by "obj.typesAvailable".
obj.typesAvailable=uniqueTypes(typeOrder);
% Makes a container for designators and a container for types
% and store them in the object
typeIndex=containers.Map;
desigIndex=containers.Map;
count=0;
for I=obj.typesAvailable
indexMap=[];
for J=1:numel(obj.ckt.(I{:}).designator)
count=count+1;
indexMap=[indexMap count]; %#ok<AGROW>
keyName=char(obj.ckt.(I{:}).designator(J));
desig_List(count)={keyName}; %#ok<AGROW>
desigIndex(keyName)=count;
end
typeIndex(char(I{:}))=indexMap;
end
obj.typeIndx=typeIndex;
obj.desigIndx=desigIndex;
obj.desigList=desig_List;
% Make a container for the index of nodes in the Aa matrix
nNodes=numel(cirObj.nodeList);
nodeIndex=containers.Map;
for I=1:nNodes
node=cirObj.nodeList(I).list;
nodeIndex(node)=I;
end
obj.nodeIndx=nodeIndex;
% Make a logical vector for accessible nodes in Aa
obj.accInAa=[cirObj.nodeList.accessibleFlag];
% Make a logical vector for for the accessible nodes in A
obj.accInA=obj.accInAa(~[cirObj.nodeList.refNodeIndex]);
% Create the incidence matrix. The unreduced incidence
% matrix is created first and then the reference (ground) node is
% removed to generate the reduced incidence matrix
A_unReduced(obj.nodeIndx.Count,obj.desigIndx.Count)=0;
for I=obj.typesAvailable
typeName=char(I);
N=numel(obj.ckt.(typeName).fromNode);
%Define the "FROM" node value
for J=1:N
indexName=char(obj.ckt.(typeName).fromNode(J));
rowIndexA=nodeIndex(indexName);
indexDesig=char(obj.ckt.(typeName).designator(J));
colIndexA=obj.desigIndx(indexDesig);
A_unReduced(rowIndexA,colIndexA)=1;
end
end
for I=obj.typesAvailable
typeName=char(I);
N=numel(obj.ckt.(typeName).toNode);
% Define the "TO" node value
for J=1:N
indexName=char(obj.ckt.(typeName).toNode(J));
rowIndexA=nodeIndex(indexName);
indexDesig=char(obj.ckt.(typeName).designator(J));
colIndexA=obj.desigIndx(indexDesig);
A_unReduced(rowIndexA,colIndexA)=-1;
end
end
obj.Aa=sparse(A_unReduced);
% Create the reduced A matrix by removing the row defined by
% the reference node
A_unReduced(nodeIndex(cirObj.refNode),:)=[];
obj.A=sparse(A_unReduced);
% Create the index for the reduced incidence matrix
VnNodeIndex=containers.Map;
nCount=0;
for I=1:nNodes
node=cirObj.nodeList(I).list;
if logical(cirObj.nodeList(I).refNodeIndex)==0
nCount=nCount+1;
VnNodeIndex(node)=nCount;
node_List(nCount)={node};  %#ok<AGROW>
end
end
obj.VnNodeIndx=VnNodeIndex;
obj.nodeList=node_List;
% Construct Ki, Kv, S, and Z matrices
Ki_(obj.desigIndx.Count,obj.desigIndx.Count)=0;
Kv_(obj.desigIndx.Count,obj.desigIndx.Count)=0;
deltaKi_=Ki_;
deltaKv_=Kv_;
S_(double(obj.desigIndx.Count)+double(obj.VnNodeIndx.Count))=0;
S_=S_';         % Make "S" a column Vector
deltaS_=S_;
for I=obj.typesAvailable
typeName=char(I);
Indx=typeIndex(typeName);
switch typeName
case 'V'
S_(Indx)=obj.ckt.V.value;
deltaS_(Indx)=ones(numel(Indx),1);
Kv_(Indx,Indx)=eye(numel(Indx));
case 'I'
S_(Indx)=obj.ckt.I.value;
deltaS_(Indx)=ones(numel(Indx),1);
Ki_(Indx)=eye(numel(Indx));
case 'R'
Ki_(Indx,Indx)=diag(obj.ckt.R.value);
deltaKi_(Indx,Indx)=eye(numel(Indx));
Kv_(Indx,Indx)=eye(numel(Indx));
case 'L'
Ki_(Indx,Indx)=diag(1j*obj.ckt.L.value);
deltaKi_(Indx,Indx)=eye(numel(Indx));
Kv_(Indx,Indx)=eye(numel(Indx));
case 'C'
Ki_(Indx,Indx)=eye(numel(Indx));
Kv_(Indx,Indx)=diag(1j*obj.ckt.C.value);
deltaKv_(Indx,Indx)=eye(numel(Indx));
case 'IDI'
Ki_(Indx,Indx)=eye(numel(Indx));
for J=1:numel(Indx)
dependentBranch=obj.desigIndx(char(obj.ckt.IDI.cntrlBranch(J)));
independentBranch=Indx(J);
Ki_(independentBranch,dependentBranch)=-obj.ckt.IDI.value(J);
deltaKi_(independentBranch,dependentBranch)=1;
end
case 'VDI'
Ki_(Indx,Indx)=eye(numel(Indx));
for J=1:numel(Indx)
dependentBranch=obj.desigIndx(char(obj.ckt.VDI.cntrlBranch(J)));
independentBranch=Indx(J);
Kv_(independentBranch,dependentBranch)=obj.ckt.VDI.value(J);
deltaKv_(independentBranch,dependentBranch)=1;
end
case 'VDV'
Kv_(Indx,Indx)=eye(numel(Indx));
for J=1:numel(Indx)
dependentBranch=obj.desigIndx(char(obj.ckt.VDV.cntrlBranch(J)));
independentBranch=Indx(J);
Kv_(independentBranch,dependentBranch)=obj.ckt.VDV.value(J);
deltaKv_(independentBranch,dependentBranch)=1;
end
end
end
obj.S=sparse(S_);
obj.Ki=sparse(Ki_);
obj.Kv=sparse(Kv_);
obj.deltaS=sparse(deltaS_);
obj.deltaKi=sparse(deltaKi_);
obj.deltaKv=sparse(deltaKv_);
obj.Z=sparse(zeros(obj.VnNodeIndx.Count,obj.VnNodeIndx.Count));
% Create an index for the completed tableau_matrix
Tindex=containers.Map;
for I=keys(obj.desigIndx)
charI=char(I);
Tindex(charI)=obj.desigIndx(charI);
obj.branchIndx=1:double(obj.desigIndx.Count);
end
count=0;
for I=keys(obj.VnNodeIndx)
count=count+1;
charI=char(I);
val=obj.VnNodeIndx(charI)+double(obj.desigIndx.Count);
obj.nodeVoltIndx(count)=val;
Tindex(charI)=val;
end
obj.TmIndx=Tindex;
end

function matrixT=tableau_matrix(T,f,varyKi,varyKv)
% Creates a tableau matrix from a tableau object at the
% specified frequency. Matrix T is a tableau object. Parameter
% f is the frequency used to create the tableau matrix.
% Variable "varyKi" and "varyKv" are parameters supplied by the
% function "monteCarlo" and provide a random permutation of
% parameters within the defined tolerance level for each
% component in Kv and Ki. Parameters T and f must be supplied.
% If varyKi and varyKv are not supplied, the tableau matrix is
% created from the nominal values of the components. If they
% are supplied, the tableau matrix will be created from a
% random distribution of component values within the stated
% component tolerance.
% Example:
%      matrixT=tableau_matrix(T,f)
%
%       or optionally
%
%      matrixT=tableau_matrix(T,f,varyKi,varyKv)
%

% Find all capacitors and set the appropriate parts of Kv to
% the desired frequency.
if any(strcmp('C',T.typesAvailable))
Indx=T.typeIndx('C');
T.Kv(Indx,Indx)=diag(1j*2*pi*f*T.ckt.C.value);
end
% Find all the inductors and set the appropriate parts of Ki to
% the desired frequency.
if any(strcmp('L',T.typesAvailable))
Indx=T.typeIndx('L');
T.Ki(Indx,Indx)=diag(1j*2*pi*f*T.ckt.L.value);
end
% If varyKi and varyKv are not supplied, use the nominal values
% of the components.
if nargin==2
matrixT=[T.Ki -T.Kv*T.A';T.A T.Z];
else
matrixT=[T.Ki.*varyKi -T.Kv.*varyKv*T.A';T.A T.Z];
end
end

function data=bodeplot(tableauObj,name,startFreq,stopFreq,nPoints)

% BODEPLOT creates a bode plot for the node or branch specified
% by name. The parameter "name" is a string specifying a branch
% or node name. The frequency range is specified by the start
% and stop frequencies. The number of points in the frequency
% range may optionally be specified. (Default is 50) The
% resulting plot is not docked in the desk top.
%
% Example:
% bodeplot(tableauObj,name,startFreq,stopFreq)
%
%  or optionally
%
% bodeplot(tableauObj,name,startFreq,stopFreq,nPoints)
%

% Check for optional arguments, define the frequency range, and
% predefine the output array IVmat for speed.
if nargin==4
nPoints=50;
elseif nargin~=5
error('bodeplot requires 4 or 5 parameters!')
end
freq=logspace(log10(startFreq),log10(stopFreq),nPoints);
IVmat(numel(tableauObj.S),nPoints)=0;
% Invert the matrix for each of the specified frequencies
for n=1:nPoints
T=tableau_matrix(tableauObj,freq(n));
IVmat(:,n)=T\tableauObj.S;
end
% Create the bode plot as two subplots and label each plot
indx=tableauObj.TmIndx(name);
magnitude=abs(IVmat(indx,:));
phase=180/pi*(angle(IVmat(indx,:)));
% Make numeric data available in the workspace before plotting
% if requested
if nargout==1
data=[freq;magnitude;phase]';
end
h=figure;
set(h,'WindowStyle','normal')
subplot(2,1,1);loglog(freq,magnitude);
xlabel('Frequency (Hz)');
if ismember(tableauObj.TmIndx(name),tableauObj.branchIndx)
ylabel(['Magnitude of Current ' name ' (Amps)']);
title(['Magnitude of the Current Through ' name]);
else
ylabel(['Voltage Magnitude At ' name ' (Volts)']);
title(['Magnitude of the Voltage At Node ' name]);
end
subplot(2,1,2);semilogx(freq,phase);
xlabel('Frequency (Hz)');
if ismember(tableauObj.TmIndx(name),tableauObj.branchIndx)
ylabel(['Phase of Current ' name ' (Degrees)']);
title(['Phase of the Current Through ' name]);
else
ylabel(['Voltage Phase At ' name ' (Degrees)']);
title(['Phase of the Voltage At Node ' name]);
end
set(h,...
'Position',[560 162 560 786],...
'Name',['Bode Plot At ' name],...
'NumberTitle','off');
end

function mCarloBodeplot(tableauObj,name,startFreq,stopFreq,nCkts,nPoints)
% MCARLOBODEPLOT creates a bode plot for the node or branch
% specified by string "name", one line on the graph for each
% circuit simulated. The frequency range is specified by the
% start and stop frequencies. The number of circuits simulated
% is specified by "nCkts". The optional variable "nPoints"
% (Default is 50) specifies the number of points in the
% frequency sweep. The resulting plot is not docked in the desk
% top.
%
% Example:
%      mCarloBodeplot(tableauObj,name,startFreq,stopFreq,nCkts)
%
%        or optionally
%
%      mCarloBodeplot(tableauObj,name,startFreq,stopFreq,nCkts,nPoints)
%

% Check for optional arguments, define the frequency range, and
% initialize IVmat, magnitude, and phase for speed.
if nargin==5
nPoints=50;
elseif nargin~=6
error('mCarloBodeplot requires 5 or 6 parameters!')
end
freq=logspace(log10(startFreq),log10(stopFreq),nPoints);
IVmat(numel(tableauObj.S),nPoints)=0;
magnitude(nCkts,nPoints)=0;
phase(nCkts,nPoints)=0;
% Create all the circuit simulations
indx=tableauObj.TmIndx(name);
for m=1:nCkts
% Create the matricies to vary component parameters
[varyKi,varyKv,varyS]=monteCarlo(tableauObj);
% Invert the matrix for each of the specified frequencies
for n=1:nPoints
T=tableau_matrix(tableauObj,freq(n),varyKi,varyKv);
IVmat(:,n)=T\(tableauObj.S.*varyS);
end
% Create the bode plot as two subplots and label each plot
magnitude(m,:)=abs(IVmat(indx,:));
phase(m,:)=180/pi*(angle(IVmat(indx,:)));
end
% Create the gain and phase plots for all circuits
h=figure;
set(h,'WindowStyle','normal')
subplot(2,1,1);loglog(freq,magnitude','color','blue');
xlabel('Frequency (Hz)');
if ismember(tableauObj.TmIndx(name),tableauObj.branchIndx)
ylabel(['Magnitude of Current ' name ' (Amps)']);
title(['Magnitude of the Current Through ' name]);
else
ylabel(['Voltage Magnitude At ' name ' (Volts)']);
title(['Magnitude of the Voltage At Node ' name]);
end
subplot(2,1,2);semilogx(freq,phase','color','blue');
xlabel('Frequency (Hz)');
if ismember(tableauObj.TmIndx(name),tableauObj.branchIndx)
ylabel(['Phase of Current ' name ' (Degrees)']);
title(['Phase of the Current Through ' name]);
else
ylabel(['Voltage Phase At ' name ' (Degrees)']);
title(['Phase of the Voltage At Node ' name]);
end
set(h,...
'Position',[560 162 560 786],...
'Name',['Bode Plot At ' name],...
'NumberTitle','off');
end

function [varyKi,varyKv,varyS]=monteCarlo(obj,Kinit,Sinit)
% Constructs a matrix that will randomly vary S, Ki, and Kv
% within the defined tolerance range of pTol and nTol. The
% input is a tableau object. The random variation is a uniform
% distribution over the stated range. To vary S (which is a
% sparse matrix) it is multiplied by 'varyS' which is also a
% sparse matrix. The multiplication is S.*varyS which is an
% element-by-element multiplication. Since both arrays are the
% same size and sparse, only the non-zero elements are actually
% multiplied. Each time the monteCarlo function is called a new
% set of multiplying values is created.
%
% Kinit and Sinit are optional parameters and should be
% initialized as sparse matricies (value of zero) to prevent
% the overhead of initializing the matricies every time.
%
% Example:
%      [varyKi,varyKv,varyS]=monteCarlo(obj)
%
%        or optionally
%
%      [varyKi,varyKv,varyS]=monteCarlo(obj,Kinit,Sinit)
%

% Initialize S_, Ki_, and Kv_ for performance
switch nargin
case 1
Ki_(obj.desigIndx.Count,obj.desigIndx.Count)=0;
Kv_(obj.desigIndx.Count,obj.desigIndx.Count)=0;
S_(double(obj.desigIndx.Count)+double(obj.VnNodeIndx.Count),1)=0;
case 3
Ki_=Kinit;
Kv_=Kinit;
S_=Sinit;
otherwise
error('monteCarlo requires 1 or 3 arguments')
end
for I=obj.typesAvailable
typeName=char(I);
Indx=obj.typeIndx(typeName);
switch typeName
case 'V'
multVal=range(obj.ckt.V.negTol,obj.ckt.V.plusTol);
S_(Indx)=multVal;
Kv_(Indx,Indx)=eye(numel(Indx));
case 'I'
multVal=range(obj.ckt.I.negTol,obj.ckt.I.plusTol);
S_(Indx)=multVal;
Ki_(Indx)=eye(numel(Indx));
case 'R'
multVal=range(obj.ckt.R.negTol,obj.ckt.R.plusTol);
Ki_(Indx,Indx)=diag(multVal);
Kv_(Indx,Indx)=eye(numel(Indx));
case 'L'
multVal=range(obj.ckt.L.negTol,obj.ckt.L.plusTol);
Ki_(Indx,Indx)=diag(multVal);
Kv_(Indx,Indx)=eye(numel(Indx));
case 'C'
Ki_(Indx,Indx)=eye(numel(Indx));
multVal=range(obj.ckt.C.negTol,obj.ckt.C.plusTol);
Kv_(Indx,Indx)=diag(multVal);
case 'IDI'
Ki_(Indx,Indx)=eye(numel(Indx));
for J=1:numel(Indx)
dependentBranch=obj.desigIndx(char(obj.ckt.IDI.cntrlBranch(J)));
independentBranch=Indx(J);
multVal=range(obj.ckt.IDI.negTol,obj.ckt.IDI.plusTol);
Ki_(independentBranch,dependentBranch)=multVal(J);
end
case 'VDI'
Ki_(Indx,Indx)=eye(numel(Indx));
for J=1:numel(Indx)
dependentBranch=obj.desigIndx(char(obj.ckt.VDI.cntrlBranch(J)));
independentBranch=Indx(J);
multVal=range(obj.ckt.VDI.negTol,obj.ckt.VDI.plusTol);
Kv_(independentBranch,dependentBranch)=multVal(J);
end
case 'VDV'
Kv_(Indx,Indx)=eye(numel(Indx));
for J=1:numel(Indx)
dependentBranch=obj.desigIndx(char(obj.ckt.VDV.cntrlBranch(J)));
independentBranch=Indx(J);
multVal=range(obj.ckt.VDV.negTol,obj.ckt.VDV.plusTol);
Kv_(independentBranch,dependentBranch)=multVal(J);
end
end
end
% Make the matricies sparse
varyS=sparse(S_);
varyKi=sparse(Ki_);
varyKv=sparse(Kv_);
% Function to create a random multiplier around unity that
% covers the tolerance range for the component.
function multVal=range(nTol,pTol)
n=numel(nTol);
multVal=1+(-nTol+(pTol+nTol).*rand(1,n))/100;
end
end

function deltaT_IV=deltaTableau_matrix(T,f,IVn)
% A circuit with its components at their nominal values can be
% represented in tableau format as
%
%      Tn*IVn = S
%
% where Tn is the tableau matrix, IVn is the vector of branch
% currents and node voltages, and S is the vector of
% independent sources. The variation of branch currents and
% node voltages can be explored by letting the components vary
% from their nominal value. (Here we will assume that the
% independent sources remain constant since source variation
% would simply be a linear response to the stimulus.) In this
% situation, the preturbed circuit represented in tableau
% format is
%
%      (Tn+deltaT)*(IVn+deltaIV)=S
%
% When the first equation is subtracted from the second
% equation, the result is
%
%      (Tn+deltaT)*(IVn+deltaIV)-Tn*IVn=0
%
% This equation can be simplified to
%
%      Tn*deltaIV=-deltaT*(IVn+deltaIV)
%
% Since the tableau matrix has an inverse, then
%
%      deltaIV=-inv(Tn)*deltaT*(IVn+deltaIV)
%
% The tableau system of equations has a property not shared by
% all methods of formulating circuit equations; a component
% value appears in one and only one row of the tableau format
% at a single location within the row. Each location in a row
% of deltaT has a number representing the change-in-value for
% that component. For each row, the change in component value
% can be extracted from deltaT and replaced by the appropriate
% value of (IVn+delatIV) because there is only one matrix
% element in a row. The equation becomes
%
%      deltaIV=-inv(Tn)*deltaT_IV*deltaComp
%
% where deltaT_IV is the matrix where the change in component
% value has been replaced by (IVn+delatIV) and deltaComp is the
% vector of changes in component values.
%
% In this form, the matrix
%
%       -inv(Tn)*deltaT_IV
%
% maps the change-in-component-value space to the
% change-in-current-and-voltage space.
%
% Example:
% deltaT_IV=deltaTableau_matrix(T,f,IVn)
%

% Extract deltaKv and deltaKi to a local variable, preserving
% the property value
delta_Kv=T.deltaKv;
delta_Ki=T.deltaKi;
% Find all capacitors and set the appropriate parts of
% T.deltaKv to the desired frequency.
if any(strcmp('C',T.typesAvailable))
Indx=T.typeIndx('C');
delta_Kv(Indx,Indx)=1j*2*pi*f*T.deltaKv(Indx,Indx);
end
% Find all the inductors and set the appropriate parts of
% T.deltaKi to the desired frequency.
if any(strcmp('L',T.typesAvailable))
Indx=T.typeIndx('L');
delta_Ki(Indx,Indx)=1j*2*pi*f*T.deltaKi(Indx,Indx);
end
% Create the beginning of deltaT_IV. Frequency and imaginary
%  components are present, (IV+deltaIV) are not yet added.
[m,n]=size(T.A);
Zero=sparse(m,n);
deltaT_IV=[delta_Ki -delta_Kv*T.A';Zero T.Z];
% Multiply each non-zero location by the appropriate value in
% IVn
L=length(deltaT_IV);
for I=1:L
col=find(deltaT_IV(I,:));
deltaT_IV(I,col)=IVn(I)*deltaT_IV(I,col);
end
end

function deltaComp=deltaComp_vector(T)
% This function builds a structure that contains three vectors.
% These vectors are deltaComp.nTol, deltaComp.pTol, and
% deltaComp.compVal. The vectors deltaComp.nTol and
% deltaComp.pTol are the fractional limits
% (tolerance_percent/100) of the components. The
% deltaComp.compVal vector is the actual component value. The
% vectors do not contain values for independent sources.
% Example:
%      deltaComp=deltaComp_vector(T)

% the values.
L=T.TmIndx.Count;
nTol(L,1)=0;
pTol(L,1)=0;
compVal(L,1)=0;
desig{L,1}='';
% Build the various elements from the types available
for I=T.typesAvailable
typeName=char(I);
Indx=T.typeIndx(typeName);
switch typeName
case 'R'
nTol(Indx)=T.ckt.R.negTol;
pTol(Indx)=T.ckt.R.plusTol;
compVal(Indx)=T.ckt.R.value;
desig(Indx)=T.ckt.R.designator;
case 'L'
nTol(Indx)=T.ckt.L.negTol;
pTol(Indx)=T.ckt.L.plusTol;
compVal(Indx)=T.ckt.L.value;
desig(Indx)=T.ckt.L.designator;
case 'C'
nTol(Indx)=T.ckt.C.negTol;
pTol(Indx)=T.ckt.C.plusTol;
compVal(Indx)=T.ckt.C.value;
desig(Indx)=T.ckt.C.designator;
case 'IDI'
for J=1:numel(Indx)
nTol(Indx(J))=T.ckt.IDI.negTol(J);
pTol(Indx(J))=T.ckt.IDI.plusTol(J);
compVal(Indx(J))=T.ckt.IDI.value(J);
desig(Indx(J))=T.ckt.IDI.designator(J);
end
case 'VDI'
for J=1:numel(Indx)
nTol(Indx(J))=T.ckt.VDI.negTol(J);
pTol(Indx(J))=T.ckt.VDI.plusTol(J);
compVal(Indx(J))=T.ckt.VDI.value(J);
desig(Indx(J))=T.ckt.VDI.designator(J);
end
case 'VDV'
for J=1:numel(Indx)
nTol(Indx(J))=T.ckt.VDV.negTol(J);
pTol(Indx(J))=T.ckt.VDV.plusTol(J);
compVal(Indx(J))=T.ckt.VDV.value(J);
desig(Indx(J))=T.ckt.VDV.designator(J);
end
end
end
% Assign the elements of deltaComp to a structure
deltaComp.nTol=nTol/100;
deltaComp.pTol=pTol/100;
deltaComp.compVal=compVal;
deltaComp.desig=desig;
end
end
end