In this talk, a coupling interface method (CIM) is proposed for solving complex interface problems, especially for Poisson-Boltzmann Equation. The coefficients, the source terms, and the solutions may be discontinuous or singular across the interfaces. The method uses adaptive-order strategy and is extended to high dimensions through a dimension-by-dimension approach. To connect information from each dimension, a coupled equation for the principal derivatives is derived through the jump conditions in each coordinate direction. The cross derivatives are approximated by one-side interpolation. The method is easy to implement and flexible to integrate with other approaches. We compare our method with some existing methods. Numerical tests are carefully performed to show the efficiency, robustness and accuracy of our method.