option nolet
print "*****************************************************"
print "* elek 2:Electrical networks(resistance),version 1.0*"
print "* peter.vlasschaert@gmail.com ,28/05/2018           *"
print "*****************************************************"
print " modification from :John G.KEMENY and THOMAS E.KURTZ*"
print "*****************************************************"
print
!*************************************************
! Modification : orginal program with line numbers
! John G.KEMENY and THOMAS E.KURTZ
! BASIC PROGRAMMING 'electrical networks ',16.2
!*************************************************
! no direct connection between node ' n=0 and n=7"
!*************************************************
! filename :elek2electricalnetworksfinalversion1.tru
! algo     :Ohms law , Kirchoff's law
!**************************************************
! X(I,J)  : current flow from I to J (neg value)
! P(I)    : potential at I
! R(I,J)  : resistance between I and J
! C(I,J)  : conductance = 1/R(I,J) between I and J
!***************************************************
! Ohms Law       : Electrical Networks
!***************************************************
! X(i,j) = 1/R(i,j)*[P(i)-P(j)] = C(i,j)*[P(i)-P(j)]
!***************************************************
! Kirchoff's Law : Electrical Networks
!***************************************************
! => current law of K irchoff
! sum( X(i,j) = 0, j <> i
!***************************************************

!*******************************************
! function ' rounding fuction ,n1 decimals '
!*******************************************
def fn_round(num,n1) = int(num*10^n1+.5)/10^n1
!
! init values
n  = 7    ! number of nodes electrical network 
v  = 24   ! voltage ( at node 1),b(1) = v
n1 = 3    ! number of decimals
! 
! matrix setup
dim c(1,1),b(1)
mat redim c(n,n),b(n)
mat c = zer(n,n)
mat b = zer(n)
c(1,1) = 1
c(n,n) = 1
b(1)   = v
!
! do loop 'exit ' => i= -1
do 
read i,j,r
if i <> -1 then
   c(i,j) = 1/r
   c(j,i) = 1/r
else
   exit do
end if
loop
!
! data
data 1,2,100 ! first
data 2,7,50
data 2,3,25
data 1,4,16
data 3,4,20
data 3,6,40
data 4,5,60
data 5,6,60
data 6,7,40   ! last
data -1,0,0   ! close data
print " conductance matrix : c(i,j) "
mat print c
!
! sum conductance
for i=2 to n-1
  s=0
 for j=1 to n
   s=s+c(i,j)
 next j
 c(i,i) = -s ! otherwise neg. values
 c(1,i) = 0
 c(n,i) = 0
next i
!
! calculation potentials at the nodes
dim u(1,1),p(1)
mat redim u(n,n),p(n)
mat u = inv(c)
mat p = u*b
print " node voltage vector : p(i) "
mat print p
!
! output node voltage
print "node","voltage"
print
for i=1 to n
    print i,fn_round(p(i),n1)
next i
print
print
print " current flow (i->j) 'neg',(j->i) 'pos' "
print
s1 = 0
for i=2 to n-1
  print "node";i
  for j=1 to n
   print fn_round(c(i,j)*(p(i)-p(j)),n1),
  next j
   print
   print
   !***************************************
   ! total current j=n,voltage between
   ! node = n (ground) and voltage node = 1
   !***************************************
   s1 = s1+c(i,n)*(p(i)-p(n))
   print
   print
next i
print 
print " total       current -  ";fn_round(s1,n1)
print " resistance  of circuit ";fn_round(v/s1,n1)  
END