We present a multi-dimensional finite element numerical model (COSFLOW) to simulate 1D canal, 2D overland, and 3D subsurface flow in South Florida. The developed model was designed for use in both regional scale analyses of various water resources projects and for the detailed design of specific projects. The diffusion wave approach was used to determine water flows on overland and in each canal reach controlled by hydraulic structures at its two ends. The Richard's equation was solved to compute the subsurface flow in both saturated/unsaturated zones. The flow governing equations were discretized with the Gelerkin finite element method, where the interaction between surface and subsurface waters is handled numerically through a coupling process. The developed numerical code has been incorporated into a modified GMS (Groundwater Modeling System) graphical user environment to allow accurate construction of computational domains and efficient use of the code for the user. From the developed hydrogeologic conceptual model, the South Dade Model that covered the area from just west of L-67 Extension eastward to the coast in South Florida was calibrated with water elevation data from observation wells in 1995. In the South Dade Model, the canal system consisted of all major canals in South Dada County where each canal reach contained many canal elements whose lengths were about 600 meters. The 2D overland flow mesh consisted of 4,720 nodes, which was formed by the surface layer of nodes and elements in the 3D subsurface mesh. The 3D finite element mesh for the subsurface flow consisted of 37,760 nodes and 65,429 triangular prism elements which were distributed over 7 discretized layers. The model calibration started with a steady 3D subsurface flow simulation using consistent boundary conditions; followed by overland/subsurface flow simulations; and finally the fully coupled 1D/2D/3D mode. Partial calibration results are presented. Suggestions to improve computer efficiency are also provided for future works.