Diode model - loops and manipulations

For any circuit, the number of independent loops is Branches-Nodes+1. The loop finder function returns excess loops for all circuits except maybe the simplest circuits. Initially, I had removed the excess loops and assumed they were linear combinations of the previous loops. But when the circuit has a ladder structure, it may so happen that some of the essential loops appear right at the end and get thrown out.

This is what I found when I was trying to figure out why the three-phase diode bridge rectifier was not working. The only way to ensure that essential loops are not deleted is to let all of them be to begin with. With row operations, some loops will be eliminated. As and how loops are eliminated, they should be completely deleted and the system should be reduced.

So, first the loop finder function (click on "View Raw" below the code box to see the code in a new window):

	def find_loop(br_map, nd_list, lp_list, lp_iter, elem, lp_count, lp_limit):
	"""Find the loops from info on branches and
	nodes. The starting point is the first branch in br_map.
	The loops found need not be independent loops."""
	no_nodes=len(br_map)
	# First branch
	start_row=elem[0]
	start_col=elem[1]
	# Move right from that element
	# This is the first element
	# In a general sense, the direction is horiz
	loop_dir="horiz"
	# The termination condition is
	# that there should not be any element
	# in the nd_list. The nodes are deleted
	# as a completed loop contains them.
	# This is to ensure that all the nodes
	# are included in the loops found.
	# To ensure that parallel loops between
	# a few pair of nodes, do not cause
	# loops to be left out, additionally,
	# it is checked whether
	# Loops < Branches - Nodes + 1
	# while (nd_list or lp_count<lp_limit):
	while (lp_iter):
	# Will be executed if we are moving horizontally
	if (loop_dir == "horiz"):
	lp_count, lp_iter=loop_horiz(br_map, nd_list, lp_list, lp_iter, elem, lp_count)
	# Change direction to vertical
	loop_dir="vert"
	# Will be executed if we are moving vertically
	if (loop_dir == "vert"):
	lp_count, lp_iter=loop_vert(br_map, nd_list, lp_list, lp_iter, elem, lp_count)
	# Change direction to horizontal
	loop_dir="horiz"
	return lp_count

view raw loop_finder.py hosted with ❤ by GitHub

The difference is in the terminating condition:
if loop_iter:
As long as loop_iter is searching, let it search. So let it add as many valid loops as possible.
The number of excess loops can be pretty huge (x6).

Next comes the main circuit_solver.py. In this, another matrix has been conceived called the system_loop_map. This is to indicate which branches are stiff so as to eliminate stiff branches from as many loops as possible. The code for this has been put together in one block (click on "View Raw" below the code box to see the code in a new window):

	# Stiff ratio indicates whether a branch
	# is stiff.
	stiff_ratio=[]
	for c1 in range(len(branch_params)):
	stiff_ratio.append("no")
	max_res=abs(branch_params[0][-1][0][0])
	for c1 in range(1, len(branch_params)):
	if branch_params[c1][-1][0][0]>max_res:
	max_res=branch_params[c1][-1][0][0]
	# Calculates the time constants of the loops
	# from the diagonal elements as L/R ratios.
	for c1 in range(len(branch_params)):
	if branch_params[c1][-1][0][0]:
	if abs(branch_params[c1][-1][0][1]/branch_params[c1][-1][0][0])<0.1*dt:
	if branch_params[c1][-1][0][0]/max_res > dt:
	stiff_ratio[c1]="yes"
	system_loop_map=[]
	for c1 in range(len(system_loops_copy)):
	br_vector=[]
	for c2 in range(len(branch_params)):
	br_vector.append("no")
	system_loop_map.append(br_vector)
	for c1 in range(len(system_loops_copy)):
	for c3 in range(len(branch_params)):
	for c2 in range(len(system_loops_copy[c1][c1])):
	if branch_params[c3][:-1]==system_loops_copy[c1][c1][c2][:-1]:
	if stiff_ratio[c3]=="yes":
	system_loop_map[c1][c3]="stiff"
	else:
	system_loop_map[c1][c3]="yes"

view raw system_loop_map.py hosted with ❤ by GitHub

So essentially, for every loop there is minimal information about every branch in the circuit - if it exists and if it does, is the branch stiff.

	# The concept here is to minimize the number of times
	# a stiff branch appears in loops. To do so a sys_loop_map
	# has been constructed. This has three options "stiff", "yes"
	# and "no". Using row operations, this sys_loop_map is made
	# into a upper triangular matrix only with respect to the
	# stiff branches.
	for c1 in range(len(sys_loop_map)):
	# Check if loop c1 has a stiff branch
	is_loop_stiff="no"
	for c2 in range(len(branch_info)):
	if sys_loop_map[c1][2]=="stiff":
	is_loop_stiff="yes"
	# Check if there is another loop after c1 (>c1)
	# which has a stiff branch "sooner" than c1 has.
	# All this is mere nomenclature.
	# If so, exchange the loops for the triangularization.
	if is_loop_stiff=="yes":
	c2=0
	while c2<len(branch_info) and sys_loop_map[c1][c2]!="stiff":
	for c3 in range(c1+1, len(sys_loop_map)):
	if sys_loop_map[c3][c2]=="stiff":
	c2=len(branch_info)-1
	for c4 in range(len(branch_info)):
	sys_loop_map[c1][c4], sys_loop_map[c3][c4] = sys_loop_map[c3][c4], sys_loop_map[c1][c4]
	sys_loops[c1][c1], sys_loops[c3][c3] = sys_loops[c3][c3], sys_loops[c1][c1]
	c2=c2+1
	# Look for the first stiff branch
	c2=0
	while c2<len(branch_info) and sys_loop_map[c1][c2]!="stiff":
	c2=c2+1
	if c2<len(branch_info):
	# Eliminate that stiff branch from subsequent loops.
	for c3 in range(c1+1, len(sys_loop_map)):
	if sys_loop_map[c3][c2]=="stiff":
	# Loop manipulations may cause loops
	# to change lenghts rapidly. So the try/except
	# statmements are used to avoid exceeding indices.
	c4=len(sys_loops[c1][c1])-1
	while c4>=0:
	try:
	sys_loops[c1][c1][c4]
	except:
	pass
	else:
	c5=len(sys_loops[c3][c3])-1
	while c5>=0:
	try:
	sys_loops[c3][c3][c5]
	except:
	pass
	else:
	if sys_loops[c1][c1][c4][:-1]==sys_loops[c3][c3][c5][:-1]:
	if branch_info[c2][:-1]==sys_loops[c1][c1][c4][:-1]:
	if sys_loops[c1][c1][c4][-1]==sys_loops[c3][c3][c5][-1]:
	loop_manipulate(sys_loops, c3, c1, "diff")
	else:
	loop_manipulate(sys_loops, c3, c1, "add")
	c5=c5-1
	c4=c4-1
	# Update the sys_loop_map info
	for c4 in range(len(branch_info)):
	if sys_loop_map[c1][c4]=="yes" and sys_loop_map[c3][c4]=="yes":
	sys_loop_map[c3][c4]="no"
	elif sys_loop_map[c1][c4]=="yes" and sys_loop_map[c3][c4]=="no":
	sys_loop_map[c3][c4]="yes"
	elif sys_loop_map[c1][c4]=="stiff" and sys_loop_map[c3][c4]=="stiff":
	sys_loop_map[c3][c4]="no"
	elif sys_loop_map[c1][c4]=="stiff" and sys_loop_map[c3][c4]=="no":
	sys_loop_map[c3][c4]="stiff"
	# This is to make sure, that stiff loops
	# are connected to an input.
	for c1 in range(len(sys_loops)):
	is_loop_stiff="no"
	has_input="no"
	for c2 in range(len(sys_loops[c1][c1])):
	# Check if any of the branches are stiff
	for c3 in range(len(branch_info)):
	if branch_info[c3][:-1]==sys_loops[c1][c1][c2][:-1]:
	if stiff_info[c3]=="yes":
	is_loop_stiff="yes"
	# Check if any branch has a non-zero B matrix entry.
	for c4 in range(len(branch_info[c3][-1][1])):
	if branch_info[c3][-1][1][c4]:
	has_input="yes"
	# If loop is stiff and has no input
	if is_loop_stiff=="yes" and has_input=="no":
	c2=0
	flag_input="no"
	# Check all other loops whether they are stiff
	# and whether they have inputs.
	while flag_input=="no" and c2<len(sys_loops):
	if not c1==c2:
	is_loop_stiff="no"
	flag_input="no"
	for c3 in range(len(sys_loops[c2][c2])):
	for c4 in range(len(branch_info)):
	if branch_info[c4][:-1]==sys_loops[c2][c2][c3][:-1]:
	if stiff_info[c4]=="yes":
	is_loop_stiff="yes"
	if is_loop_stiff=="no":
	for c5 in range(len(branch_info[c4][-1][1])):
	if branch_info[c4][-1][1][c5]:
	flag_input="yes"
	c2=c2+1
	# Perform a row operation with a loop that is non-stiff
	# and has an input.
	if is_loop_stiff=="no" and flag_input=="yes":
	c2=c2-1
	loop_manipulate(sys_loops, c1, c2, "add")

view raw stiffness_removal.py hosted with ❤ by GitHub

So essentially, I use system_loops_map to make the system upper triangular as far as stiff branches are concerned. The next block I am not sure if it is needed. Whether it is necessary to make sure a stiff loop is connected to the input. I'll test it and try to get rid of it. For some reason, it looks like an ugly code block.

The last part of to reduce the size of the system by getting rid of redundant loops (click on "View Raw" below the code box to see the code in a new window):

	# This is to recalculate the sys_loops matrix.
	# First the diagonal loops are recalculated.
	# Then the off-diagonal (interactions)
	readjust_sys_loops(sys_loops, branch_info, stiff_info)
	# Re-initialize the matrices A, B, and E
	matrix_a.zeros(len(sys_loops),len(sys_loops))
	matrix_e.zeros(len(sys_loops),len(sys_loops))
	matrix_b.zeros(len(sys_loops),matrix_b.columns)
	# Recalculate the matrices A, B and E
	for c1 in range(len(sys_loops)):
	for c2 in range(len(sys_loops)):
	for c3 in range(len(sys_loops[c1][c2])):
	for c4 in range(len(branch_info)):
	if sys_loops[c1][c2][c3][:-1]==branch_info[c4][:-1]:
	if c1==c2:
	matrix_a.data[c1][c2]+=branch_info[c4][-1][0][0]
	matrix_e.data[c1][c2]+=branch_info[c4][-1][0][1]
	if sys_loops[c1][c2][c3][-1]=="forward":
	for c5 in range(matrix_b.columns):
	matrix_b.data[c1][c5]+=branch_info[c4][-1][1][c5]
	else:
	for c5 in range(matrix_b.columns):
	matrix_b.data[c1][c5]-=branch_info[c4][-1][1][c5]
	else:
	if sys_loops[c1][c2][c3][-1]=="forward":
	matrix_a.data[c1][c2]+=branch_info[c4][-1][0][0]
	matrix_e.data[c1][c2]+=branch_info[c4][-1][0][1]
	else:
	matrix_a.data[c1][c2]-=branch_info[c4][-1][0][0]
	matrix_e.data[c1][c2]-=branch_info[c4][-1][0][1]
	# This part if to eliminate the empty (redundant) loops
	# based on zero rows in A, E and B matrices.
	# At the same time, A, B and E matrices are resized.
	for c1 in range(len(sys_loops)-1, -1, -1):
	empty_loop="yes"
	for c2 in range(matrix_a.columns):
	if matrix_a.data[c1][c2]:
	empty_loop="no"
	for c2 in range(matrix_e.columns):
	if matrix_e.data[c1][c2]:
	empty_loop="no"
	for c2 in range(matrix_b.columns):
	if matrix_b.data[c1][c2]:
	empty_loop="no"
	if empty_loop=="yes":
	# If a loop is empty, the column c1 is first deleted
	# and then the row c1 is deleted.
	for c2 in range(len(sys_loops)-1, -1, -1):
	for c3 in range(len(sys_loops[c2][c1])-1, -1, -1):
	del sys_loops[c2][c1][c3]
	del sys_loops[c2][c1]
	for c2 in range(len(sys_loops[c1])-1, -1, -1):
	del sys_loops[c1][c2]
	del sys_loops[c1]
	# Same for the matrices A and E
	for c2 in range(matrix_a.rows-1, -1, -1):
	del matrix_a.data[c2][c1]
	del matrix_e.data[c2][c1]
	matrix_a.columns-=1
	matrix_e.columns-=1
	for c2 in range(matrix_a.columns-1, -1, -1):
	del matrix_a.data[c1][c2]
	del matrix_e.data[c1][c2]
	del matrix_a.data[c1]
	del matrix_e.data[c1]
	matrix_a.rows-=1
	matrix_e.rows-=1
	for c2 in range(matrix_b.columns-1, -1, -1):
	del matrix_b.data[c1][c2]
	del matrix_b.data[c1]
	matrix_b.rows-=1
	# Make sure that the stiff loops have no inductances
	# so that they are treated as static equations.
	for c1 in range(len(sys_loops)):
	for c2 in range(len(sys_loops[c1][c1])):
	for c3 in range(len(branch_info)):
	if sys_loops[c1][c1][c2][:-1]==branch_info[c3][:-1]:
	if stiff_info[c3]=="yes":
	for c4 in range(matrix_e.columns):
	matrix_e.data[c1][c4]=0.0

view raw sys_red.py hosted with ❤ by GitHub