We develop a numerical scheme for nonlinear Poisson-Boltzmann equation. First, we regularize the solution of PBE to remove the singularity. We introduce the discontinuous bubble function to treat the nonhomogeneous jump conditions of the regularized solution. Next, starting with an initial guess, we apply linearization to treat the nonlinearity. Then, we discretize the discontinuous bubble and the bilinear form of PBE. Finally, we solve the discretized linear problem by IFEM. This process is repeated by updating the previous approximation. We carry out numerical experiments. We observe optimal convergence rate for all examples.