I maintain a math hobby - magic square web page _http://www.knechtmagicsquare.paulscomputing.com/_ (http://www.knechtmagicsquare.paulscomputing.com/) A very bright fellow in England sent me a algorithm written in python to solve a problem I was working on. I can't seem to get it to work ... syntax problems.... The segment below gives the explanation and the code: If you are interested in the problem this solves ... see section 1 on my website "Topographical Model For the Magic Square" Sorry about the delay in replying; I've been having email problems.]
> I read a article by your wife on magic squares. I have been fooling > around with magic square models. I read your bio and see you did some work > with random walk stuff. I have been wanting a algorithm ... ( perhaps a > Greedy algorithm) that would calculate the amount of retained water for any > order of magic square. ( see topographical model setion I on web link > below) Hm, interesting problem. Let's see. (Please forgive me if some or all of what I'm about to say is totally familiar and/or obvious to you.) If I've understood you right, your problem is as follows. Call the number written in each cell its "height". Define the "water height" of a cell to be the smallest N such that there's a path from that cell to the outside of the square using only cells of height <= N. So, when you pour an infinite amount of water over the square and let it settle, each cell ends up with water up to its water height, if that's bigger than its height. Then the "capacity" of the square is the sum of (water height - height) over cells for which that's positive, and you want to calculate the capacity of a given square. So, the following seems like a pretty good algorithm. It's a bit like Dijkstra's shortest-path algorithm. It needs a data structure called a "priority queue", which means a container that stores objects along with numbers called their "priorities", and has efficient ways of doing (1) add an object to the queue and (2) pull out the object with highest priority. We'll keep track of an upper bound on the water height for each cell, and use the queue to store cells for which we might be able to reduce the upper bound. No, actually we'll use it to store cells that might enable us to reduce their neighbours' upper bounds. 1. For each cell, set bound(cell) = n^2. (This is a very boring upper bound on the water height of the cell.) 2. For each cell on the edge of the square: 2a. Set bound(cell) = height(cell). 2b. Put that cell into the queue with priority -bound(cell). (So lower cells have higher priority. Actually, many priority queue implementations already use the "smaller number = higher priority" convention.) 3. While the queue isn't empty: 3a. Pull out the highest-priority cell from the queue. 3b. For each neighbour of that cell: 3b1. Let b = max(bound(cell),height(neighbour)). 3b2. If if bound(neighbour) > b: 3b2a. Set bound(neighbour) = b. 3b2b. Add the neighbour to the queue with priority -b. When this is finished, for every cell bound(cell) should equal (not just be an upper bound for) the cell's water height, because: 1. bound(cell) is always an upper bound on the cell's water height, because the only way it gets decreased is when we find that the cell has a neighbour whose water height is known to be <= its new value. (So there's a path from the cell to the edge, using nothing higher than that value, beginning with that neighbour.) 2. For cells along the edge, bound(cell) equals the cell's water height, because for these cells height = water height, and right at the start we set bound = height. 3. Suppose that when the algorithm finishes bound(cell) is actually bigger than cell's water height w. Then there is a path from the cell to the edge using nothing of height > w. Follow that path until you find a cell for which bound <= w. (There must always be one, since bound = height for edge cells and the last cell on our path is an edge cell and has height <= w.) Then we've found two adjacent cells of height <= w, one of which has bound <= w and one of which as bound > w. But at some point the former must have been put into the queue with bound <= w; when it came out, its neighbour would have been examined and had its bound decreased to <= w. Contradiction. >From 1 and 3, we find that at the end of the algorithm the bounds equal the water heights. Now we can just sum max(bound-height,0) over all cells, and get the total capacity of the square. Here's some Python code that does it (rather inefficiently): def process(square): n = len(square) # square is a list of lists bounds = [n*[n*n] for i in range(n)] queue = [] # lousy implementation, just for algorithm testing # head of queue is last element # big numbers for high priority for i in range(n-1): bounds[i][0] = square[i][0] queue.append((-bounds[i][0],i,0)) bounds[n-1][i] = square[n-1][i] queue.append((-bounds[n-1][i],n-1,i)) bounds[i+1][n-1] = square[i+1][n-1] queue.append((-bounds[i+1][n-1],i+1,n-1)) bounds[0][i+1] = square[0][i+1] queue.append((-bounds[0][i+1],0,i+1)) while queue: queue.sort() (b,i,j) = queue.pop() for (di,dj) in [(-1,0),(1,0),(0,-1),(0,1)]: u=i+di; v=j+dj if 0<=u<n and 0<=v<n: b1 = max(square[u][v],-b) if b1 < bounds[u][v]: bounds[u][v] = b1 queue.append((-b1,u,v)) return bounds def capacity(square): bounds = process(square) n = len(square) s = 0 for i in range(n): for j in range(n): delta = bounds[i][j]-square[i][j] if delta>0: s += delta return s Testing on one of your examples: >>> capacity([[6,26,49,34,35,13,12],[5,48,1,28,30,36,27], [47,2,19,20,21,29,37],[46,17,10,25,7,32,38],[45,24,22,3,33,9,39], [11,44,31,23,8,40,18],[15,14,43,42,41,16,4]]) 320 (Your page says 321, but I think 320 is correct. The sum of the heights in the submerged region is 394, not 393.) How efficient is this? Well, it's not quite obvious to me that a cell can't be put into the queue more than once, but I think that should be rare. So let's suppose that typically each cell goes in once. The queue should never be bigger than the number of cells aside from complications that arise when a cell is inserted more than once; it certainly shouldn't be bigger than, say, 4 times the number of cells. And each priority queue operation should take no worse than some modest multiple of log(queue_size) basic operations. So the total time to process the whole square shouldn't be worse than some small multiple of n^2 log n. That seems like about the best you could reasonably expect. (I have another approach in mind that might reduce that to n^2 f(n) where f grows even more slowly when n is very large, but at the cost of greater complexity when n is not very large. I bet you'll only be doing this for n not very large.) Note: the idiotic "sort the list every time" priority queue implementation in the Python code above is emphatically *not* suitable for production use. Find a proper priority queue implementation in your language of choice. -- g If some one could help me get this working ... I would be grateful. I have loaded Python Ver 2.6 ... but when I paste this code in it keeps giving me errors. Thanks Craig Knecht **************A Good Credit Score is 700 or Above. See yours in just 2 easy steps! (http://pr.atwola.com/promoclk/100126575x1222585064x1201462784/aol?redir=http://www.freecreditreport.com/pm/default.aspx?sc=668072&hmpgID=62&bcd= JunestepsfooterNO62)
-- http://mail.python.org/mailman/listinfo/python-list