Efficient matrix inverse/solving linear equation with syms

Hi, sorry for the bad title. I am trying to solve the following equation as λ varies over ℝ. My current implementation seems slower than I would expect.

n=100;
lambda = 0:0.1:10;
syms a;
OneVector = repmat(1,n,1);
%Definitions of A and B - these don't really matter
A = diag(repmat(2,n,1)) + diag(repmat(-1,n,1),1) + diag(repmat(-1,n,1),-1) + diag(-1, -n +1) + diag(-1, n-1);
B = diag(rand(n, 1));
M = A + a*B;
%The whole process until now takes under a second for n = 100
x = linsolve(M, OneVector); %This is the long step, as it should be
y = x.^{-1}; %The actual vector I am interested in
SurfaceData = double(subs(y, a, lambda));
surf(SurfaceData); 

One immediate question I had was if it is more efficient to subsitute the syms before inverting the elements of the vector or to do it after as I am currently doing. The other question is if there is a faster method for doing this computation. I realize that the linsolve uses an algorithm and that its not reasonable to expect too much from it, but scaling it up from to is more than times the processing time.