Code covered by the BSD License  

Highlights from
Resistor Calculator

image thumbnail
from Resistor Calculator by Ligong Han
An equivalent resistance calculator for resistor network

RNet
function RNet
set(gcf,'color',[0    0.2850    0.3985],'menubar','none',...%[.3 .81 .66]
    'numbertitle','off','name','Circuit Calculator');
rvalue = uicontrol('style','edit','position',[30,50,80,20],...
    'background','white');
rok = uicontrol('style','toggle','string','UPDATE',...
    'position',[30,25,80,20]);
nok = uicontrol('style','toggle','string','NODE',...
    'position',[450,50,80,20]);
calc = uicontrol('style','toggle','string','CALC',...
    'position',[450,25,80,20]);
rdisp = uicontrol('style','text','position',[140,25,280,45]);
noksqare = [8.70,0.12;10.48,0.12;10.48,0.65;8.70,0.65];
drawnow
axis([0 10 0 10]);
axis off
cm = [0 0 0];
pos = [];%pos(k,:) is position of node labeled k
nn = 0;
ne = 0;
stt = true;%true: on, ready to insert; false: off, waiting for the 2nd node
hold on
while 1%~get(nok,'value')
    %display current operation
    if stt
        set(rdisp,'string',['Insert a new edge: The ',order_disp(ne+1),...
            ' edge. First node.']);
    else
        set(rdisp,'string',['Insert a new edge: The ',order_disp(ne),...
            ' edge. Second node.']);
    end
    [x,y] = ginput(1);
    if get(nok,'value')||((x-8.70)*(x-10.48)<0&&(y-0.12)*(y-0.65)<0)
        break;
    end
    %search
    %newnode=true;
    k = 1;
    while k <= nn
        if ((x-pos(k,1))^2+(y-pos(k,2))^2)<=0.25^2
            k = -k;%newnode=false;
            break;
        end
        k = k+1;
    end
    %insert an edge
    hold on
    if k > 0%newnode
        nn = nn+1;
        plot(x,y,'ro');
        text(x+.15,y,['N.',num2str(k)]);
        pos = [pos;x,y];
        k = -k;
    end
    if stt
        %1st node
        ne = ne+1;
        cm(ne,1) = -k;
        %stt=false;
    else
        %2nd node
        cm(ne,2) = -k;
        %R value
        set(rdisp,'string','Input resistance value and click UPDATE.');
        while ~get(rok,'value')
            %waiting
            pause(0.01);
        end
        r = str2double(get(rvalue,'string'));
        %set(rvalue,'string',[]);%% clear the text
        cm(ne,3) = r;
        x_1 = pos(cm(ne,1),1);
        y_1 = pos(cm(ne,1),2);
        plot([x_1,x],[y_1,y],'y');
        ang = atan((y-y_1)/(x-x_1)); a = .2;
        t = [a,pi-a,pi+a,2*pi-a]+ang;
        x_m = (x+x_1)/2;
        y_m = (y+y_1)/2;
        fill(.72*cos(t)+x_m,.72*sin(t)+y_m,'b');
        text(x_m,y_m,[num2str(r),' Ohm']);
        set(rok,'value',false);
        %stt=true;
    end
    stt = ~stt;
    pause(0.01);
end
%input terminal
t = zeros(1,2);
tdisp = 'AB';
for k = 1:2
    set(rdisp,'string',['Input Terminal ',tdisp(k),'.']);
    [x,y] = ginput(1);
    for kk = 1:nn
        if ((x-pos(kk,1))^2+(y-pos(kk,2))^2) <= 0.25^2
            break;
        end
    end
    t(k) = kk;
    plot(pos(kk,1),pos(kk,2),'ks');
    text(x+.12,y-.425,['Terminal ',tdisp(k)]);
end
set(rdisp,'string','Click CALC to calculate equivalent resistance.');
while ~get(calc,'value')
    pause(0.01);
end
%clac
figure('name','Calculation Report')%
netMatrix = cm;
terminalA = t(1);
terminalB = t(2);
% Req =- 1;
edgeNum = size(netMatrix,1);
node_temp = netMatrix(1:2*edgeNum);
nodeLabel = node_temp(1);
num = 2;
for k = 2:2*edgeNum%create vector nodeInd
    if(~size(find(~(nodeLabel-node_temp(k))),2))
        nodeLabel(num) = node_temp(k);
        num = num+1;
    end
end
%k node_temp num
% clear node_temp num;
nodeNum = size(nodeLabel,2);
value = inf*ones(nodeNum);
for k = 1:edgeNum%create matrix value
    edge_temp = netMatrix(k,:);
    p = find(nodeLabel == edge_temp(1));
    q = find(nodeLabel == edge_temp(2));
    R = edge_temp(3);
    value(p,q) = R;
    value(q,p) = R;
end
%k edge_temp p q R
% clear edge_temp p q R;
%draw graph
AInd = find(nodeLabel == terminalA);
BInd = find(nodeLabel == terminalB);
num_current = 0;
edgeLabel_original = zeros(nodeNum);
for p = 1:nodeNum%create original edgeLabel
    for q = p+1:nodeNum
        if value(p,q)<inf
            num_current = num_current+1;
            edgeLabel_original(p,q) = num_current;
            edgeLabel_original(q,p) = num_current;
        end
    end
end
%num_current p q
% clear num_current;
visited = false*ones(1,nodeNum);
queue = AInd;
visited(AInd) = true;
added = [];
BFSTree_original = [];
nodeList_original = [];
while ~isempty(queue)%create original BFSTree
    indv = queue(1);
    nodeList_original(end+1) = indv;
    queue(1) = [];
    tempVector = edgeLabel_original(indv,:);
    att = find(tempVector);
    for k = 1:size(att,2)
        if visited(att(k)) == false
            added(end+1) = att(k);
            visited(att(k)) = true;
            BFSTree_original(end+1) = edgeLabel_original(indv,att(k));
        end
    end
    queue = [queue,added];
    added = [];
end
%visited queue added indv att tempVector k
% clear visited queue added indv att tempVector;
%ShowingGraph
l = 5;
ang = 2*pi/nodeNum;
x_original = zeros(1,nodeNum); x_original(1) = 5;
y_original = zeros(1,nodeNum); y_original(1) = 0;
for k = 1:nodeNum-1
    x_original(k+1) = x_original(k)+l*cos(k*ang);
    y_original(k+1) = y_original(k)+l*sin(k*ang);
end
subplot(2,2,1);%Resistor Network
%axis off%%%
hold on
for k = 1:nodeNum
    text(x_original(k),y_original(k),['Node ',num2str(nodeLabel(k))]);
    plot(x_original(k),y_original(k),'o');
end
for k = 1:edgeNum
    edge_temp = netMatrix(k,:);
    p = find(nodeLabel == edge_temp(1));
    q = find(nodeLabel == edge_temp(2));
    R = edge_temp(3);
    plot([x_original(p),x_original(q)],[y_original(p),y_original(q)]);
    xr = (x_original(p)+2*x_original(q))/3;
    yr = (y_original(p)+2*y_original(q))/3;
    text(xr,yr,[num2str(R),' Ohm','  No.',num2str(edgeLabel_original(p,q))]);
end
plot(x_original(AInd),y_original(AInd),'rsquare');
plot(x_original(BInd),y_original(BInd),'rsquare');
title('Original Resistor Network');
hold off
subplot(2,2,2);%BFSTree
%axis off%%%
hold on
for k = 1:nodeNum
    text(x_original(k),y_original(k),['Node ',num2str(nodeLabel(k))]);
    plot(x_original(k),y_original(k),'o');
end
plot(x_original(AInd),y_original(AInd),'rsquare');
plot(x_original(BInd),y_original(BInd),'rsquare');
edgeInd_temp = triu(edgeLabel_original);
for k = 1:nodeNum-1
    indEdge = BFSTree_original(k);
    [p,q] = find(edgeInd_temp == indEdge);
    plot([x_original(p),x_original(q)],[y_original(p),y_original(q)],'r');
end
title('Original BFSTree');
hold off
% clear p q R edge_temp xr yr edgeInd_temp indEdge k;
% clear x_original y_original nodeList_original edgeLabel_original BFSTree_original;
%-------------------------%
%simplify
k = 1;
while k <= nodeNum;%scan
    AInd = find(nodeLabel == terminalA);
    BInd = find(nodeLabel == terminalB);
    if (k == AInd) || (k == BInd)
        k = k+1;
        continue;
    end
    degree = size(find(value(k,:) < inf),2);
    if degree == 2
        att = find(value(k,:)<inf);
        Rpq = sum(value(k,att));
        R_original = value(att(1),att(2));
        if R_original < inf
            Rpq = Rpq*R_original/(Rpq+R_original);
            edgeNum = edgeNum-1;
        end
        edgeNum = edgeNum-1;
        value(att(1),att(2)) = Rpq;
        value(att(2),att(1)) = Rpq;
        value(k,:) = [];
        value(:,k) = [];
        nodeLabel(k) = [];
        nodeNum = nodeNum-1;
        k = 1;
        continue;
    end
    k = k+1;
end
%k AInd BInd degree att Rpq R_original
% clear k AInd BInd degree att Rpq R_original;
ind = 1;
while ind<=nodeNum%simplify
    AInd=find(nodeLabel==terminalA);
    BInd=find(nodeLabel==terminalB);
    if (ind==AInd)||(ind==BInd)
        ind=ind+1;
        continue;
    end
    degree=size(find(value(ind,:)<inf),2);
    if degree==3%Y2Delta
        att=find(value(ind,:)<inf);
        R=value(ind,att);
        numerator=R(1)*R(2)+R(2)*R(3)+R(3)*R(1);
        Rpq=numerator/R(3);
        Rqr=numerator/R(1);
        Rrp=numerator/R(2);
        R_original=value(att(1),att(2));%serR
        if R_original<inf
            Rpq=Rpq*R_original/(Rpq+R_original);
            edgeNum=edgeNum-1;
        end
        edgeNum=edgeNum+1;
        value(att(1),att(2))=Rpq;
        value(att(2),att(1))=Rpq;
        R_original=value(att(2),att(3));%setR
        if R_original<inf
            Rqr=Rqr*R_original/(Rqr+R_original);
            edgeNum=edgeNum-1;
        end
        edgeNum=edgeNum+1;
        value(att(2),att(3))=Rqr;
        value(att(3),att(2))=Rqr;
        R_original=value(att(3),att(1));%setR
        if R_original<inf
            Rrp=Rrp*R_original/(Rrp+R_original);
            edgeNum=edgeNum-1;
        end
        edgeNum=edgeNum+1;
        value(att(3),att(1))=Rrp;
        value(att(1),att(3))=Rrp;
        value(ind,:)=[];
        value(:,ind)=[];
        nodeLabel(ind)=[];
        nodeNum=nodeNum-1;
        edgeNum=edgeNum-3;
        k=1;
        while k<=nodeNum;%scan
            AInd=find(nodeLabel==terminalA);
            BInd=find(nodeLabel==terminalB);
            if (k==AInd)||(k==BInd)
                k=k+1;
                continue;
            end
            degree=size(find(value(k,:)<inf),2);
            if degree==2
                att=find(value(k,:)<inf);
                Rpq=sum(value(k,att));
                R_original=value(att(1),att(2));
                if R_original<inf
                    Rpq=Rpq*R_original/(Rpq+R_original);
                    edgeNum=edgeNum-1;
                end
                edgeNum=edgeNum-1;
                value(att(1),att(2))=Rpq;
                value(att(2),att(1))=Rpq;
                value(k,:)=[];
                value(:,k)=[];
                nodeLabel(k)=[];
                nodeNum=nodeNum-1;
                k=1;
                continue;
            end
            k=k+1;
        end
        ind=1;
        continue;
    end
    ind=ind+1;
end
%ind AInd BInd degree att R numerator Rpq Rqr Rrp R_original k
% clear degree att R numerator Rpq Rqr Rrp R_original k ind
%-------------------------%
AInd=find(nodeLabel==terminalA);
BInd=find(nodeLabel==terminalB);
%draw graph
num_current=0;
edgeLabel=zeros(nodeNum);
for p=1:nodeNum%create edgeLabel
    for q=p+1:nodeNum
        if value(p,q)<inf
            num_current=num_current+1;
            edgeLabel(p,q)=num_current;
            edgeLabel(q,p)=num_current;
        end
    end
end
%num_current p q
% clear num_current p q;
visited=false*ones(1,nodeNum);
queue=AInd;
visited(AInd)=true;
added=[];
BFSTree=[];
nodeList=[];
while ~isempty(queue)%create BFSTree
    indv=queue(1);
    nodeList(end+1)=indv;
    queue(1)=[];
    tempVector=edgeLabel(indv,:);
    att=find(tempVector);
    for k=1:size(att,2)
        if visited(att(k))==false
            added(end+1)=att(k);
            visited(att(k))=true;
            BFSTree(end+1)=edgeLabel(indv,att(k));
        end
    end
    queue=[queue,added];
    added=[];
end
% clear k att tempVector queue added visited indv;
%ShowingGraph
x=zeros(1,nodeNum);x(1)=5;
y=zeros(1,nodeNum);y(1)=0;
for k=1:nodeNum-1
    x(k+1)=x(k)+l*cos(k*ang);
    y(k+1)=y(k)+l*sin(k*ang);
end
subplot(2,2,3);%Resistor Network
%axis off%%%
hold on
for k=1:nodeNum
    text(x(k),y(k),['Node ',num2str(nodeLabel(k))]);
    plot(x(k),y(k),'o');
end
value_temp=[];
for p=1:nodeNum
    for q=p+1:nodeNum
        if value(p,q)<inf
            value_temp=[value_temp;fix(p) fix(q) value(p,q)];
        end
    end
end
% clear p q;
for k=1:edgeNum
    p=value_temp(k,1);
    q=value_temp(k,2);
    R=value_temp(k,3);
    plot([x(p),x(q)],[y(p),y(q)]);
    xr=(x(p)+2*x(q))/3;
    yr=(y(p)+2*y(q))/3;
    text(xr,yr,[num2str(R),' Ohm','  No.',num2str(edgeLabel(p,q))]);
end
plot(x(AInd),y(AInd),'rsquare');
plot(x(BInd),y(BInd),'rsquare');
title('Simplified Resistor Network');
hold off
edgeInd_temp=triu(edgeLabel);
subplot(2,2,4);%BFSTree
%axis off%%%
hold on
for k=1:nodeNum
    text(x(k),y(k),['Node ',num2str(nodeLabel(k))]);
    plot(x(k),y(k),'o');
end
plot(x(AInd),y(AInd),'rsquare');
plot(x(BInd),y(BInd),'rsquare');
for k=1:nodeNum-1
    indEdge=BFSTree(k);
    [p,q]=find(edgeInd_temp==indEdge);
    plot([x(p),x(q)],[y(p),y(q)],'r');
end
title('Simplified BFSTree');
hold off
% clear p q R xr yr edgeInd_temp value_temp indEdge k;
% clear l ang x y;
%KVL&KCL Equations
linkBranch=[];
for p=1:nodeNum%create linkBranch
    for q=p+1:nodeNum
        if (edgeLabel(p,q)~=0)&&(size(find(BFSTree==edgeLabel(p,q)),2)==0)
            linkBranch(end+1)=edgeLabel(p,q);
        end
    end
end
loops={};
loop_temp=[];%set, not list or array
loopValues={};
loopValue=[];
edgeInd_temp=triu(edgeLabel);
for k=1:size(linkBranch,2)%create loops
    loopValue=inf*ones(nodeNum);
    loop_temp=[BFSTree,linkBranch(k)];
    for kk=1:nodeNum%create loopValue
        [p,q]=find(edgeInd_temp==loop_temp(kk));
        loopValue(p,q)=value(p,q);
        loopValue(q,p)=value(p,q);
    end
    %loopNode=nodeLabel;%index of loopNode ~= index of node; set, not list or array
    num=nodeNum;
    nodeMark=false*ones(1,nodeNum);
    kk=1;
    while kk<=num%delete degree==1, kk is node index
        if nodeMark(kk)==true
            kk=kk+1;
            continue;
        end
        degree=size(find(loopValue(kk,:)<inf),2);
        if degree==1
            loop_temp(loop_temp==edgeLabel(kk,loopValue(kk,:)<inf))=[];
            loopValue(kk,:)=inf;
            loopValue(:,kk)=inf;
            nodeMark(kk)=true;
            kk=1;
            continue;
        end
        kk=kk+1;
    end
    loops{k}=loop_temp;
    loopValues{k}=loopValue;
end
% clear nodeMark num;
direction={};
dirct=[];%index of this list means edge label
loopNodeSeq={};
nodeSeq=[];%set, not list or array
for k=1:size(linkBranch,2)%create direction; p->q,next
    [p,q]=find(edgeInd_temp==loops{k}(1));
    dirct=zeros(1,edgeNum);
    dirct(loops{k}(1))=1;
    nodeSeq=[p q];
    for kk=2:(size(loops{k},2))
        att=find(loopValues{k}(q,:)<inf);
        if att(1)==p
            next=att(2);
        else
            next=att(1);
        end
        nodeSeq(end+1)=next;
        if next>q
            dirct(edgeLabel(q,next))=1;
        else
            dirct(edgeLabel(q,next))=-1;
        end
        p=q;
        q=next;
    end
    direction{k}=dirct;
    loopNodeSeq{k}=nodeSeq;
end
% clear loop_temp loopValue loopValues dirct nodeSeq p q next k kk att nodeList degree linkBranch;
KVL=zeros(edgeNum-nodeNum+1,edgeNum);
KCL=zeros(nodeNum-1,edgeNum);
R=zeros(1,edgeNum);
k=1;
for p=1:nodeNum%create vector R
    for q=p+1:nodeNum
        if value(p,q)<inf
            R(k)=value(p,q);
            k=k+1;
        end
    end
end
for k=1:edgeNum-nodeNum+1%KVL
    KVL(k,:)=direction{k}.*R;
end
% clear p q k direction R;
current=zeros(edgeNum,1);
sumCurrent=zeros(edgeNum,1);
for k=1:nodeNum-1%KCL; k is node index
    att=find(value(k,:)<inf);
    if k==AInd
        sumCurrent(k)=1;
    elseif k==BInd
        sumCurrent(k)=-1;
    end
    for kk=1:size(att,2)
        if att(kk)>k
            KCL(k,edgeLabel(k,att(kk)))=1;
        else
            KCL(k,edgeLabel(k,att(kk)))=-1;
        end
    end
end
% clear k kk att loops loopNodeSeq;
equation=[KCL;KVL];
% clear KVL KCL;
%solve equations
%EQUATION*CURRENT=SUMCURRENT
%CURRENT=EQUATION\SUMCURRENT
current=equation\sumCurrent;
% clear equation sumCurrent;
%calculate Req
treeValue=inf*ones(nodeNum);
num=nodeNum;
nodeMark=false*ones(1,nodeNum);
nodeMark(AInd)=true;
nodeMark(BInd)=true;
for k=1:nodeNum-1%create treeValue
    [p,q]=find(edgeInd_temp==BFSTree(k));
    treeValue(p,q)=value(p,q);
    treeValue(q,p)=value(p,q);
end
k=1;
%create path
while k<=num%delete degree==1 except terminals, k is node index, RNet do not change
    if nodeMark(k)==true
        k=k+1;
        continue;
    end
    degree=size(find(treeValue(k,:)<inf),2);
    if degree==1
        treeValue(k,:)=inf;
        treeValue(:,k)=inf;
        nodeMark(k)=true;
        k=1;
        continue;
    end
    k=k+1;
end
% clear value BFSTree degree nodeMark num;
p=AInd;
q=find(treeValue(p,:)<inf);
Req=treeValue(p,q)*current(edgeLabel(p,q));
while q~=BInd%p->q, search along path
    att=find(treeValue(q,:)<inf);
    if att(1)==p
        next=att(2);
    else
        next=att(1);
    end
    p=q;
    q=next;
    Req=Req+treeValue(p,q)*current(edgeLabel(p,q));
end
disp(['equivalent resistance: Req = ',num2str(abs(Req)),' Ohm']);
% clear att next AInd BInd current edgeInd_temp edgeLabel edgeNum k netMatrix nodeLabel nodeNum p q terminalA terminalB treeValue;
set(rdisp,'string',['Equivalent Resistance: Req = ',num2str(abs(Req)),' Ohm']);

    function orderstr = order_disp(num)
        numStr = num2str(num);
        % onesdigit = num-10*floor(num/10);
        % tensdigit = floor((num-100*floor(num/100))/10);
        if ((length(numStr)>1)&&(numStr(end-1)~='1')) || length(numStr)<=1
            switch numStr(end)
                case '1'
                    ordersuffix = 'st';
                case '2'
                    ordersuffix = 'nd';
                case '3'
                    ordersuffix = 'rd';
                otherwise
                    ordersuffix = 'th';
            end
        else
            ordersuffix = 'th';
        end
        orderstr = [numStr ordersuffix];
    end

end

Contact us