I think what you might be missing is that FIND is returning an index, so the eta value you seek is eta(F12_p). But I didn't really try your approach. Here's how I might have done it.
etaspan = -500:0.001:500;
f = @(x) (x.^(1/2))./(1+exp(x-etaspan));
y = integral(f,0,500,'ArrayValued',true);
[ys,idx] = sort(y);
findeta = @(y)etaspan(interp1(ys,idx,y,'nearest','extrap'));
FIND uses a linear search. INTERP1 uses a binary search, which should be faster. Now if you didn't want to use the table, fzero might work.
f = @(x,eta)(x.^(1/2))./(1+exp(x - eta));
findeta = @(y)fzero(@(eta)integral(@(x)f(x,eta),0,500) - y,0);