A direct numerical simulation of the three-dimensional elektrokinetic instability near a charge-selective surface (electric membrane, electrode, or system of micro- or nanochannels) has been carried out and analyzed. A special finite-difference method has been used for the space discretization along with a semi-implicit 31/3-step Runge-Kutta scheme for the integration in time. The calculations employ parallel computing. Three characteristic patterns, which correspond to the overlimiting currents, are observed: (a) two-dimensional electroconvective rolls, (b) three-dimensional regular hexagonal structures, and (c) three-dimensional structures of spatiotemporal chaos, which are a combination of unsteady hexagons, quadrangles, and triangles. The transition from (b) to (c) is accompanied by the generation of interacting two-dimensional solitary pulses.