Both numerical simulation and theoretical analysis of seawater intrusion in coastal regions are of great theoretical importance in environmental sciences. The mathematical model can be described as a coupled system of three dimensional nonlinear partial differential equations with initial-boundary value problems. In this paper, according to the actual conditions of molecular and three-dimensional characteristic of the problem, we construct the characteristic finite element alternating-direction schemes which can be divided into three continuous one-dimensional problems. By making use of tensor product algorithm, and priori estimation theory and techniques, the optimal order estimates in H1 norm are derived for the error in the approximate solution. Both numerical simulation and theoretical analysis of seawater intrusion in coastal regions are of great theoretical importance in environmental sciences. The mathematical model can be described as a coupled system of three dimensional nonlinear partial differential equations with initial-boundary value problems. In this paper, according to make the use of tensor product algorithm, and priori estimation theory and techniques, the optimal order estimates in H1 norm are derived for the error in the approximate solution.