n = @(x, z, T) 1./(exp(x./T).*z + 1);
f = @(x) (3/2)*sqrt(x);
eqchem = @(z, T) integral(@(x) f(x)*n(x, z, T), 0, inf) - 1;
At this point the Mathematica code provided breaks down. It uses z without z being a parameter. We could add it as a parameter to zSolution, an assumption that it should have been written as
zSolution[T_, z_] := ....
but then we need to pass in two parameters in the reference to zSolution that appears in points. We can see that there is explicit range of T values being iterated over in the Table construction, but there is no match for z there. Is the table built to be symbolic in z? Or is the {0., 1.} intended to be a z range, or is that {0., 1.} intended to be a final (or leading) point in the table being constructed?
Perhaps I am misunderstanding the meaning of the zSolution code. The /. operator looks like you might have intended a rule replacement, but the right side of the ./ does not appear to be in the form of a Mathematica rule. Is the intention that zSolution should be the z such that z is a root of eqchem near 0.01 ?
