function out=lapsolve1(bot,right,top,left)

% session 1, exercise 3, solves laplace in a square
% - bot,right,top,left are row vectors of equal length
%

% test: the vectors must have equal lengths
botl=length(bot); rightl=length(right); 
topl=length(top); leftl=length(left);
if (topl ~= rightl | topl ~= botl | topl ~= leftl) 
  error('error');
end

% if test passed...
len=botl;        % number of points per original row/column        
m=(len-1);       % number of subintervals per row/column
n=m-1;           % number of unknowns per row/column
l=m+1;           % number of points per row/column

% formation of A-matrix, this is essentially problem 1 (&2)
% this could be written as a function taking n as an argument
maindiag=4*ones(n*n,1);
sidediag=-1*ones((n*n)-1,1);
for i=1:n-1,              % set every n:th to zero
  sidediag(i*n)=0;
end
fardiag=-1*ones((n-1)*n,1);
A=diag(maindiag,0)+diag(sidediag,-1)+diag(sidediag,1)+diag(fardiag,-n)+diag(fardiag,n);

% formation of right hand side
B=[bot(3:l-2)' zeros(n-2,n-2) top(3:l-2)'];
B=[left(2:l-1);B;right(2:l-1)];
B(1,1)=B(1,1)+bot(2);
B(1,n)=B(1,n)+top(2);
B(n,1)=B(n,1)+bot(l-1);
B(n,n)=B(n,n)+top(l-1);
rhs=B(:);

sol=inv(A)*rhs;

solmat=reshape(sol,n,n);

solmat=[bot(2:l-1)' solmat top(2:l-1)'];
solmat=[left; solmat; right];

solmat(1,1)=(solmat(1,2)+solmat(2,1))/2;  % averages to corners
solmat(1,l)=(solmat(1,l-1)+solmat(2,l))/2;
solmat(l,1)=(solmat(l,2)+solmat(l-1,1))/2;
solmat(l,l)=(solmat(l,l-1)+solmat(l-1,l))/2;

out=solmat;

