Of the many processes for rapid prototyping collectively referred to as 3D printing, fused deposition modeling (FDM) or fused filament fabrication (FFF), is perhaps the best known. In FDM complex, three-dimensional objects are built by depositing filaments of hot polymers that fuse together when they solidify. The process is a versatile way for generating fairly complex geometries, but although the physical processes behind FDM are easily understood, the quality of the final product depends sensitively on exactly how the process is carried out. The current understanding of how to select the various control parameters has mostly been achieved by experimentation, but detailed numerical simulations should be able to provide a much more rigorous way to do so.
In our first paper on modeling FFF/FDM  we introduced a numerical method for fully resolved simulations of the deposition of a filament of hot polymer onto a substrata and its cooling. Viscosity was taken to be a function of temperature so the object was essentially rigid, once it had cooled down. The conservation equations for mass, momentum and energy are solved on a fixed, structured and staggered grid, discretizing a domain containing two or more fluids, and the interface between the different fluids is tracked using connected marker particles. The method was extended to FFF/FDM by making the computations of the diffusion terms for momentum and temperature implicit to handle realistic material parameters; the viscosity of the hot polymer was taken to be a function of temperature and shear rate; and a volume source was incorporated to model the nozzle. The method was implemented for a rectangular domain initially containing air and the dynamics of both the polymer and the air was simulated. The performance and accuracy of the method was tested in several ways using the deposition of two short filaments, one on top of the other, and it was found that for governing parameters similar to those encountered under realistic conditions a converged solution could be obtained using about thirty control volumes across thefilament diameter. A more complex geometry, consisting of a two layer in-filled rectangle, was also simulated to show the capability of the method. The computational setup was, however, simplified in several ways. The most significant of those were that the injection nozzle was modeled using a simple volume source, and that except for the change in viscosity with temperature no model was included for the material behavior of the polymer as it solidified. The method was extended in , where a more realistic model for the nozzle was implemented, using an immersed boundary to represent a closed cylinder, into which the polymer was injected, and an equation for the material deformation gradients was evolved in time, to enable the inclusion of elastic stresses in the solidified polymer. A material model accounting for the shrinkage of the polymer, as it solidified, was also added, thus making it possible to compute residual stresses. The accuracy of the method was investigated by convergence tests, as in , and the capabilities demonstrated by simulating the fabrication of a two-layer part of a few side-by-side filaments. For other studies see .
In those papers the behavior of the molten polymer was modeled using relatively simple relationships for the stresses. It is well known that the orientation of the polymer chains can have a significant impact of the properties of the final object and to predict those it is necessary to solve for the evolution of the conformation tensor as the polymer filament is extruded and laid down. In our most recent paper , the methodology was extended to include the evolution of the confirmation tensor and the viscoelastic stresses, using the Rolie-Poly model to include the effects of the polymer, The implementation was first tested on the extrusion of a polymer from a nozzle, where we could do a comparison with other results, and then used to simulate the deposition of a few filaments on top of each other to examine what the effect of including the viscoelastic stresses had on the shape of the filaments. Although including the conformation tensor allows us to add viscoelastic stresses, we have not included the full Rolie-Poly equation. Thus, we in effect assume that the Rouse time is infinite, and ignore convective constraint release. Similarly, we did not include entanglement losses.
The video shows a simulation of the deposition of a few short filaments from , where the injection nozzle its modeled using an immersed boundary method and the solid stresses are found as the polymer cools down and solidifies.
 H. Xia, J. Lu, S. Dabiri, and G. Tryggvason. Fully resolved numerical simulations of fused deposition modeling. part i–fluid flow. Rapid Prototyping Journal, 24:463–476, 2018.
 H. Xia, J. Lu, and G. Tryggvason. Fully resolved numerical simulations of fused deposition modeling. part ii—solidification, residual stresses, and modeling of the nozzle. Rapid Prototyping Journal, 2018.
 H. Xia, J. Lu and G. Tryggvason. Simulations of Fused Filament Fabrication using a Front Tracking Method. International Journal of Heat and Mass Transfer. 138 (2019) 1310–1319.
 H. Xia, J. Lu, and G. Tryggvason. A numerical study of the effect of viscoelastic stresses in fused filament fabrication. Computer Methods in Applied Mechanics and Engineering, 346:242–259, 2018.