The numerical simulation of multicomponent compressible flow in porous media is an important research topic in reservoir modeling. Traditional semi-implicit methods for such problems are usually conditionally stable, suffer from large splitting errors, and may accompany with violations of the boundedness requirement of the numerical solution. In this study we reformulate the original multicomponent equations into a nonlinear complementarity problem and discretize it using a fully implicit finite element method. We solve the resultant nonsmooth nonlinear system of equations arising at each time step by a parallel, scalable, and nonlinearly preconditioned semismooth Newton algorithm, which is able to preserve the boundedness of the solution and meanwhile treats the possibly imbalanced nonlinearity of the system. Some numerical results are presented to demonstrate the robustness and efficiency of the proposed algorithm on the Tianhe-2 supercomputer for both standard benchmarks as well as realistic problems in highly heterogeneous media.