function b=readrhs(n)

fid = fopen('rhs');
b = fread(fid, [n,n], 'double');
fclose(fid);
