#!/usr/bin/env python

# Homework 2, problem 1.  Conventions and notation follow the C
# version, hw2.1.c, except that we don't use the crazy Numerical
# Recipes numbering convention.
#
# Run with
#
#	hw2.1.py  N
#
# where N is the number of nodes in the network.

# The system consists of N fully connected nodes.  If we set V_1 = 0
# and ignore the equation for I_1, the problem reduces to an
# (N-1)x(N-1) matrix, for nodes 2 to N, which is matrix A below.
# Similarly the current vector I starts at I_2.  However, python wants
# to number starting at 0, hence the offsets...

import os
import sys

from numpy import zeros
from numpy.linalg import solve

def resistance(i, j):		# i and j are nodes (1:N)
    if i == j:
        return 0
    else:
        return min(i,j) + 2*max(i,j)

N = 5
if len(sys.argv) > 1:
  N = int(sys.argv[1])

n1 = N - 1

# Set up R and I.

R = zeros((n1,n1))
for i in range(n1):		# i = row, j = column
    for j in range(n1):		# 0:N-2 --> nodes 2:N
        if (i == j):
            R[i,j] = 0
            for k in range(N):	# k spans all nodes (0:N-1 --> nodes 1:N)
                if k != i+1:
                    R[i,j] -= 1.0/resistance(i+2, k+1)
        else:
            R[i,j] = 1.0/resistance(i+2, j+2)

I = zeros(n1)
I[0] = 1			# I[0] refers to node 2

# print R

# Solve R.V = I.

V = solve(R, I)
print "N =", N, "  R = -V2 =", -V[0]