Flows with free surfaces and free boundaries arise in many industrial and biological applications. Examples are coating, polymer processing, ink-jet printing, DNA arrays, spraying, deformation of blood cells, blood flow in arteries and capillaries, and flow in the deep pulmonary alveoli. Most of these flows have two distinguishing features: (1) the fluid is complex (microstructured ); thus, the stress includes a visco-elastic term which is important and sometimes dominant, and (2) the surface forces are comparable to the viscous and elastic ones. Inertia is usually immaterial in these flows, because the relevant length scales are well below a millimeter. The surface and viscoelastic forces give rise to large non-diagonal contributions in the momentum equation. Other non-diagonal terms come from the coupling of the shape of the free boundaries to the velocity field, and strong dependence of the microstructure evolution on the velocity gradient. Thus, fully-coupled algorithms for solving the steady as well as time-dependent equations of the flow are desirable. I will discuss developments in applying mesoscopic models of microstructured liquids to three-dimensional free surface flows. In such models, the liquid microstructure is captured by tensors obeying convection-diffusion-generation equations—e.g., the gyration tensor of ensembles of polymer molecules, or the shape tensor of droplets or blood cells. Mesoscopic non-equilibrium thermodynamics ties the elastic stress to velocity-gradient-dependent terms in the microstructure evolution. This dependence yields general theories accounting for disparate microstructural models that are compatible with macroscopic transport phenomena and thermodynamics. Such theories can be incorporated into general three-dimensional finite element codes based on fully coupled formulations. Combining Newton’s method with GMRES and a Sparse Approximate Inverse Preconditioner yields a robust and efficient method for computing three-dimensional flows on low-cost parallel clusters. I will show results on model flows of polymer solutions, and discuss developments and connections to fine-grain, microscopic models of complex fluids where microstructure is tracked by using stochastic differential equations.