#include "dec.ocm"

proc prompt(value str[])=
	seq i=[1 for str[byte 0]]
		output ! str[byte i]
:
def N=20 :

var n:
var A[N*N], x[N], k[N], y[N] :

proc initialise=
	var c:
	seq
		prompt("n?*n")
		c:='*s'
		decin(input, n, c)

		prompt("A?*n")
		seq i= [0 for n]
			seq j= [0 for n]
				decin(input, A[(i*n)+j], c)

		prompt("x?*n")
		seq i= [0 for n]
			decin(input, x[i], c)

		prompt("k?*n")
		seq i= [0 for n]
			decin(input, k[i], c)				:
		
proc produce.xj(value j, chan south) =
	-- north row: source of x values
	while true
		south ! x[j]						:

proc consume.yi(value i, chan east) =
	-- west column: read y values
	east ? y[i]						:

proc offset(value ki, chan west) =
	-- east column: source of k offsets
	while true
		west ! ki						:

proc multiplier(value aij, chan north, south, west, east) =
	-- middle: responsible for a values
	var xj, aij.times.xj, yi :
	seq
		north ? xj
		while true
			seq
				par
					south ! xj
					aij.times.xj:= aij*xj
					east ? yi
				par
					west ! yi+aij.times.xj
					north ? xj			:

proc sink(chan north) =
	-- south row: sink for unused outputs
	while true
		north ? any						:

seq
	initialise

	chan north.south[(N+1)*N], east.west[N*(N+1)] :
	par
		par j= [0 for n]	-- producer of co-ordinates x[j]
			produce.xj(j, north.south[j])

		par			-- the matrix multiplier
			par i= [0 for n]
				offset(k[i], east.west[(n*n)+i])
			par i= [0 for n]
				par j= [0 for n]
					multiplier(A[(n*i)+j],
						north.south[(n*i)+j],
						north.south[(n*(i+1))+j],
						east.west[i+(n*j)],
						east.west[i+(n*(j+1))])
			par j= [0 for n]
				sink(north.south[(n*n)+j])

		seq
			par i= [0 for n]-- consumer of transformed co-ordinates
				consume.yi(i, east.west[i])

			seq i= [0 for n]
				seq
					output ! 'y'; '['
					decout(output, i, 0)
					output ! ']'; '='
					decout(output, y[i], 5)
					output ! '*n'
			exit(0)