Parametric-based, lumped watershed models have been widely employed for integrated surface and groundwater modelling to calculate surface runoff on various temporal and spatial scales of hydrologic regimes. Physics-based, process-level, distributed models that have the design capability to cover multimedia and multi-processes and are applicable to various scales have been practically nonexistent until late 1990s. It has long been recognized that only such models have the potential to further the understanding of the fundamental factors that take place in nature hydrologic regimes; to give mechanistic predictions; and most importantly to be able to couple and interact with weather/climate models. However, there are severe limitations with these models that inhibit their use. These are, among other things, the ad hoc approaches of coupling between various media, the simplification of modelling overland and/or river flow, and the excessive demand of computational time. This paper presents the development of an integrated media (river/stream networks, overland regime, and subsurface media), integrated processes (evaporation, evapotranspiration, infiltration, recharges, and flows) watershed model to address these issues. Rigorous coupling strategies are described for interactions among overland regime, rivers/streams/canals networks, and subsurface media. The necessities to include various options in modelling surface runoff and river hydraulics are emphasized. The options of selecting characteristic wave directions for two-dimensional problems are stated. The implementation of high performance computing to increase the computational speed is discussed. Four examples are used to demonstrate the flexibility and efficiency of the model as applied to a theoretical benchmark scale, a parallel computing, and two project-level large scale problems - one in Taiwan and the other in Florida.