Buckling of elastic columns in the presence of rigid obstacles is formulated in terms of an eigenvalue complementarity problem. The resolution of the eigenvalue complementarity problem for the computation of bifurcation loads and instability modes by a semismooth Newton method is proposed. The unilateral constraints materialized by the obstacles are rigorously satisfied without penalty-type techniques. The algorithm is able to compute all the lower order modes with a relatively small grid of initial values.