#hw26.txt, 5/3/14, Anthony Zaleski

Help:=proc():
print(`ExpD2(M,N)`):
end:

#ExpD2(M,N): Inputs positive integers, and outputs a list of lists,L, such that
#L[i][j] (1 ≤ i ≤ M-1,1 ≤ j ≤ N-1 )
#is the expected number of unit steps left for a random walker who starts at [i,j] 
#and moves with equal probability to one of its immediate neighbors
#[i-1,j],[i+1,j],[i,j-1],[i,j+1]
#and ends the walk when reaching one of the edges i=0, i=M, j=0, j=N.

ExpD2:=proc(M,N) local L,i,j,S,ij,eq,var:
var:={seq(seq(L[i,j],j=0..N),i=0..M)}:
eq:={seq(L[0,j]=0,j=0..N),seq(L[M,j]=0,j=0..N),seq(L[i,0]=0,i=0..M),seq(L[i,N]=0,i=0..M)}:

if nops(eq)=nops(var) then RETURN([seq([seq(0,j=0..N)],i=0..M)]): fi:

eq:=eq union {seq(seq(L[i,j]=1+(L[i-1,j]+L[i+1,j]+L[i,j-1]+L[i,j+1])/4,j=1..N-1),i=1..M-1)}:
var:=solve(eq,var):

if var=NULL then RETURN(FAIL): fi:
subs(var,[seq([seq(L[i,j],j=1..N-1)],i=1..M-1)]):
end: