A numerical implementation of two-phase flows of Newtonian fluids with a non-linear viscoelastic interface is validated and applied to the case of uniaxial extension of a drop in a matrix fluid. Second-order convergence in space and time is obtained with a Lagrangian-based interface tracking finite element method. The flow problem is analysed using dimensionless groups based on the relative magnitudes of the viscoelastic interfacial extra stress, the interfacial tension and the viscous stress of the bulk fluids. After fitting the intrinsic viscoelastic stress-strain behaviour of interfaces in shear to experimental results from the literature, the influence of interfacial rheology on the drop shape and interfacial stress is investigated. The drop shape is not significantly influenced by the viscoelastic properties of the interface if the interfacial viscoelastic stress, interfacial tension and bulk viscous stress are of the same order of magnitude. However, the interface develops distinct stress profiles for varying interfacial viscoelastic properties. For relatively large viscoelastic interface stress compared to interfacial tension and bulk viscous stress, simulations become unstable. Eventually, the interfacial viscoelastic stress exceeds the interfacial tension and bulk viscous stress, and possibly buckling of the interface occurs, caused by compressive stresses at the drop tip.

We present a fully-coupled numerical scheme for computing steady and time-dependent viscoelastic flows. The scheme relies on the contravariant deformation tensor formulation and uses a Newton-Raphson iteration to solve the non-linear system of equations. The contravariant reformulation allows for the computation and implementation of the analytical Jacobian relatively easily, especially compared to other reformulations such as the log-conformation. The contravariant deformation tensor rotates in steady state shearing flows, which is solved hereby "resetting"it as a pre-processing step in the numerical scheme, rather than a post-processing step. We use the finite element method with standard stabilization techniques (SUPG and DEVSS-G) for the spatial discretization. The numerical scheme is tested in three viscoelastic flow problems which are studied in terms of stability and accuracy: planar Couette flow, 2D flow around a cylinder and 3D flow around a sphere. For all problems, quadratic convergence is observed in both the difference between iterations and the residuals during the Newton-Raphson procedure. Moreover, we observe that the residuals are several orders smaller than the difference between iterations. A distinct advantage of the numerical scheme presented here, is that it significantly relaxes the requirement on the time-step size in time-dependent problems, as compared to explicit or semi-implicit methods. Moreover, steady states can be efficiently computed if the initial guess in the Newton-Raphson iteration is close enough to the solution.

A framework is presented for the formulation of a class of continuum constitutive models for sharp interfaces with non-linear viscoelastic behaviour due to a considerable isotropic interfacial microstructure. For the formulation of a thermodynamically consistent elastoviscoplastic interface constitutive model we adapt an approach successful in describing the behaviour of bulk polymer glasses. The model has a clear separation between dilatation and shear, and is used to predict phenomena related to the plasticity of interfaces observed in the experimental literature, which is relevant for many applications. Stress-strain predictions in standard interfacial rheological flows, i.e. shear and dilatation, are investigated numerically. A predominantly elastic response is obtained at small deformations, with a transition to primarily plastic flow at high stress levels. In interfacial shear flow, strain softening and eventually a plastic plateau occur upon further deformation beyond the yield point. The yield stress and strain and (the relative strength of) the stress overshoot in interfacial shear flow are shown to be controlled by two dimensionless groups of parameters in the model. In interfacial dilatation, the model predicts elastoviscoplastic behaviour with a stress maximum and a decreasing stress without a plateau at even larger deformations. These phenomena are studied for various choices for the parameters in the model.

The dynamics of a cell suspended in a Newtonian liquid subjected to a pressure-driven flow at non-negligible inertia in cylindrical and square cross section microfluidic channels is studied through three-dimensional arbitrary Lagrangian–Eulerian finite-element numerical simulations. The cell is modeled through the neo-Hookean hyper-elastic constitutive equation, which can describe biological particles undergoing moderate deformations. The cell-to-channel relative dimension is fixed to 0.2, whereas the Reynolds number Re, measuring the relative importance of liquid inertial and viscous forces, and the elastic capillary number Cae, measuring the relative importance of liquid viscous stress and solid elastic stress, are varied by several orders of magnitude. In a cylindrical tube, the cell migrates transversally to the flow direction until reaching a radial equilibrium position depending on Re and Cae. Given Re, the softer the cell (i.e., the larger Cae) the closer its equilibrium position to the tube axis, thus allowing for the separation of healthy and diseased cells which have similar dimensions but different mechanical properties. In a channel with a square cross section, a much more complex dynamics is found. Depending on Re and Cae, the cell can either migrate to the channel centerline, to the closest median of the channel cross section (thus, four equilibrium positions can be identified due to symmetry), to the closest diagonal (again, four equilibrium positions), or to an intermediate position in between the median and the diagonal (eight equilibrium positions).

The blend morphology model developed by Wong et al., based on Peters et al., is used to investigate the development of the disperse polymer blend morphology in twin-screw extruder flow. First, the model is written in a point-wise form suitable for using in conjunction with particle tracking. Particle tracking methods are used to generate trajectories along the flow field. Macroscopic droplet populations are placed along these trajectories and the velocity gradient tensor is extracted and applied on the point-wise blend morphology model. Very large morphology differences arise between trajectories that pass through the middle gap and those that do not. In the global distribution of (macroscopically averaged, monodisperse) droplet sizes, two distinct peaks appear due to these different trajectories. Given enough number of screw rotations, a droplet population can reach almost every position in the twin-screw extruder and travel along both types of particle trajectories. The effect of varying the gap size is that the largest droplets are unaffected, but the smallest droplets are smaller for a smaller gap size due to the higher maximum shear rate. The effect of varying the viscosity ratio on the global droplet size distribution is found to be nonlinear and is strongly determined by the Grace curve. The effect on polydisperse droplet populations is found to be that trajectories that do not pass through the gap evolve toward a single peak, whereas trajectories that do pass through the gap lead to a split into two peaks that ultimately rejoin as one peak. It is concluded that the initial position of a population in the twin-screw extruder has a very large effect on the developing transient blend morphology, though future work should be done on the importance of the initial position on the steady-state blend morphology after a very large number of screw rotations.

External gear pumps are used in fluid transport systems because of their tight clearances and accurate flow control. These tight clearances are a challenge for numerical studies in terms of spatial discretization. In earlier work, the flow of a viscous fluid in an external gear pump is computed using the finite element method (FEM). An element size based on the respective distance between boundaries is proposed. In this study, results based on the earlier work are compared to extrusion experiments of EPDM. The aim of this study is not only to validate the numerical simulations, but also to determine what material characteristics need to be taken into account for an accurate output prediction of the external gear pump. Especially the introduction of shear-thinning behavior results in an improvement of the amplitude of the pressure difference fluctuation. Taking into account compressibility, alters the torque fluctuation in such a way that it mimics the experiments. Unfortunately, the fluctuation in torque still has a too high amplitude. Eventually, simulations are performed including shear-thinning behavior, a temperature- and pressure-dependent viscosity, and compressibility. The effect of measuring the material behavior using oscillatory or shear experiments is shown. Furthermore, the simulations are applied to a second EPDM. Finally, different processing conditions are tested. For the simulations, only qualitative agreement is found, possibly as a result of the no slip boundary condition.

Profile extrusion is a production technique to make polymer profiles and hollow tubes. At a certain distance after the die exit the polymer extrudate is guided into a cooling chamber, where the polymer is cooled extensively to set the part dimensions. Uniform cooling in the cooling chamber is required to avoid warpage and twisting of the extrudate due to differential shrinkage. The effect of extrudate cooling on the shape of the extrudate and the residual stresses in the final extrudate product are numerically studied. Two different cases are studied: a 2D axisymmetric problem of a viscoelastic fluid exiting an annular die and a 3D problem of a hollow square shaped die. The extrudate is extensively cooled from a distance Lc$L_{\text{c}}$ from the die exit. Extrudate cooling is modeled by subjecting the extrudate to heat loss via convection. For LdieLc$z>L_{\text{c}}$ the Nusselt number is varied to mimic different cooling speeds. The Weissenberg number Wi, the Nusselt number and Lc$L_{\text{c}}$ are varied in this work. Results show that for high values of Wi a significant skin-effect of high residual stresses appear near the extrudate edges for high Nusselt numbers and small values of Lc$L_{\text{c}}$. It is also shown that for high Weissenberg numbers the shape of the extrudate is significantly altered due to extrudate swell, when the extrudate is cooled slower. For fast cooling and small values of Lc$L_{\text{c}}$, the polymer does not have time to relax and the elastic stresses are "frozen in." This leads to less extrudate swell but large residual stresses near the extrudate edges.

In an extrusion process, an external gear pump can be used to control the flow rate of the system. When extruding polymers, the viscosity is quite high, resulting in negligible inertia and thus laminar flow. The external gear pump contains two gears, one driven by a motor and one driven by means of contact with the other gear. In our previous work, the flow of a viscous fluid through an external gear pump was studied using the finite element method. Local mesh refinement was applied based on the respective distance between boundaries. Furthermore, the rotation of both gears was imposed. In this work, the rotation of one gear is imposed, whereas the other gear is freely rotating. However, the minimum distance between the gears is limited to a minimum value. When this value is reached, contact is assumed and also the rotation of second gear is imposed. A reversion of the torque on this gear results in a release of contact. In this manner, a quasi driver/driven situation is created in the numerical simulations. It is observed that contact is released periodically, and thus cannot be assumed present continuously, as is often prescribed. Non-Newtonian material properties, such as shear thinning and the pressure dependence of the density or the viscosity, alter how long contact is released during a tooth rotation.

In this paper we study the numerical stability of four positive (semi-)definite formulations for viscoelastic fluid flows: the LCR formulation of Fattal and Kupferman (2004), the Cholesky-log formulation of Vaithianathan and Collins (2003), the b-symmetric formulation of Balci et al. (2011) and the CDT formulation of Carrozza et al. (2019). The numerical stability is tested by evaluating the critical time step in four different flow problems: the flow around a cylinder confined between two walls, planar extrudate swell, an axisymmetric swell problem and the extension of the axisymmetric swell problem to three dimensions. Although LCR is a good performer throughout, the CDT formulation turns out to be the most stable/robust one. The numerical stability of the Cholesky-log and the b-symmetric formulation turns out to be much lower in at least one benchmark. This supports our hypothesis that a vector-like formulation increases numerical stability at least as much as the LCR formulation and that positive definiteness is not the main contribution to numerical stability.

The extrusion of highly filled elastomers is widely used in the automotive industry. In this paper, we numerically study the effect of thixotropy on 2D planar extrudate swell for constant and fluctuating flow rates, as well as the effect of thixotropy on the swell behavior of a 3D rectangular extrudate for a constant flowrate. To this end, we used the Finite Element Method. The state of the network structure in the material is described using a kinetic equation for a structure parameter. Rate and stress-controlled models for this kinetic equation are compared. The effect of thixotropy on extrudate swell is studied by varying the damage and recovery parameters in these models. It was found that thixotropy in general decreases extrudate swell. The stress-controlled approach always predicts a larger swell ratio compared to the rate-controlled approach for the Weissenberg numbers studied in this work. When the damage parameter in the models is increased, a less viscous fluid layer appears near the die wall, which decreases the swell ratio to a value lower than the Newtonian swell ratio. Upon further increasing the damage parameter, the high viscosity core layer becomes very small, leading to an increase in the swell ratio compared to smaller damage parameters, approaching the Newtonian value. The existence of a low-viscosity outer layer and a high-viscosity core in the die have a pronounced effect on the swell ratio for thixotropic fluids.

The blend morphology model developed by Wong et al. ( Rheologica Acta , 2019), based on Peters et al. ( J Rheol 45(3):659–689, 2001), is used to investigate the development of the polydispersity of the disperse polymer blend morphology in complex flow. First, the model is extended with additional morphological states. The extended model is tested for simple shear flow, where it is found that the droplet size distribution does not simply scale with the shear rate, because this scaling does not hold for coalescing droplets. Subsequently, the model is applied to Poiseuille flow, showing formation of distinct layers, which occurs in realistic pressure-driven flows. Finally, the model is applied on an eccentric cylinder flow, where histograms are made of the average droplet size throughout the domain. It is observed that outer cylinder rotation results in narrow distributions where the small droplets are relatively large, whereas inner cylinder rotation results in broad distributions where the small droplets are significantly smaller than in the case of outer cylinder rotation. Eccentricity seems to only have a minor effect if the maximum shear rate is held constant. The flow profile and history in combination with the maximum shear rate strongly determine how the polydisperse droplet size distribution develops.

Gear pumps are frequently used for the transport of high-viscosity fluids, for instance, for the extrusion of polymers. The flow rate of this extrusion process is regularly controlled by an external gear pump. In this work, the 2D flow of a viscous fluid through such an external gear pump is studied using the finite element method. Local mesh refinement based on the respective distance between moving boundaries is essential to capture the relatively narrow clearances in the pump. The gear pump works against the pressure driven flow, therefore its performance is strongly dependent on material and processing parameters. The flow of Newtonian and shear-thinning fluids through the external gear pump is studied for a range of processing conditions. Pump curves are obtained that display the volumetric efficiency against the Hersey number, which is defined as viscosity times rotation speed divided by pressure difference over the pump. Analysis of the residence time of the fluid in the pump, shows that vortices are present in the inflow channel causing material to remain in the pump for longer times.

In this paper, an experimental strategy is presented to characterize the rheological behavior of filled, uncured rubber compounds. Oscillatory shear experiments on a regular plate-plate rheometer are combined with a phenomenological thixotropy model to obtain model parameters that can be used to describe the steady shear behavior. We compare rate- and stress-controlled kinetic equations for a structure parameter that determines the deformation history-dependent spectrum and, thus, the dynamic thixotropic behavior of the material. We keep the models as simple as possible and the characterization straightforward to maximize applicability. The model can be implemented in a finite element framework as a tool to simulate realistic rubber processing. This will be the topic of another work, currently under preparation. In shaping processes, such as rubber- and polymer extrusion, with realistic processing conditions, the range of shear rates is far outside the range obtained during rheological characterization. Based on some motivated choices, we will present an approach to extend this range.

In this paper we propose a novel approach to solve the inverse problem of three-dimensional die design for extrudate swell, using a real-time active control scheme. To this end, we envisioned a feedback connection between the corner-line finite element method, used to predict the positions of the free surfaces of the extrudate, and the controller. The corner-line method allows for local mesh refinement and transient flow to be taken into account (Spanjaards et al., 2019). We show the validity of this method by showing optimization results for 2D axisymmetric extrusion flows of a viscoelastic fluid for different Weissenberg numbers. In 3D we first give a proof of concept by showing the results of a Newtonian fluid exiting dies with increasing complexity in shape. Finally, we show that this method is able to obtain the desired extrudate shape of extrudates of a viscoelastic fluid for different Weissenberg numbers and different amounts of shear-thinning.

A 2D hysteretic Discrete Element Method (DEM) model is developed for simulating the flow of food particles, specifically with multi-material 3D food printing processes in mind. Particles are modeled as arbitrarily shaped polygons due to the diverse nature of food powders, which can be highly irregular in shape. The developed hysteretic force model is applicable to both convex and concave polygonal particles. It is adjusted to use a proportional weighted maximum intersection area upon splitting of contact areas. This results in a continuous force trajectory, which would otherwise not be guaranteed. The model is validated with in literature reported packing ratios for ellipses. Simulated deposition of sugar-shaped particles shows that for dense packings the fraction of splitting interactions occurring can be up to 7%. Furthermore, the influence of particle shape on the coordination number and packing density is shown with a simulation of the deposition of sugar-like material.

Fluctuating viscoelasticity for conformation-tensor-based models is studied at equilibrium, in simple-shear deformation, and in uniaxial extension. The models studied are the upper-convected Maxwell model, the FENE-P model with finite chain-extensibility, and the Giesekus model with anisotropic drag. Using numerical simulations, the models are compared in detail both with each other and with analytical predictions for the Maxwell model. At equilibrium, the models differ only marginally, both in terms of static and dynamic characteristics. When deformed, the average mechanical response of the Maxwell model is unaffected by the strength of thermal fluctuations, while the mechanical response of the FENE-P and Giesekus models show a slight decrease the stronger the fluctuations in simple shear, whereas the decrease in uniaxial extension is marginal. For all models, the standard deviation of the mechanical response increases with increasing strength of fluctuations, and the magnitude of the standard deviation relative to the average for given fluctuation strength generally decreases the stronger the deformation, this effect being stronger for uniaxial extension than for simple-shear deformation.

Selective laser sintering (SLS) is an important additive manufacturing method for polymers and is increasingly used for industrial applications. To improve the quality of printed products a thorough understanding of the process is necessary. Simulations are a convenient way to assess the sintering in detail. In literature several numerical studies can be found on the sintering process of two particles. Both 2D and 3D geometries are used and sometimes the assumption of axisymmetry is made. If the effect of the laser beam in temperature-dependent sintering is taken into account, axisymmetry cannot be assumed anymore and full 3D simulations are required. However, these simulations are complex and computationally expensive. In this work, we assess if 2D simulations are representative for the 3D case of the temperature-dependent sintering process of two viscous particles. We find that 2D simulations are a good alternative for 3D simulations for the sintering problem of two particles, if the correct geometrical scaling factors are taken into account.

In the selective laser sintering of polymers, the most widely used powders are based on polyamide 12 (PA12), which is a semi-crystalline polymer. Because the mechanical properties of the printed parts depend largely on the microstructure, knowledge on the crystalline architecture is important. We developed a numerical model based on the finite element method to solve the flow, temperature and crystallization kinetics of PA12 powder during sintering using two different geometries. Our results show that the temperature plays a crucial role in the crystallization kinetics and that simplified 0D calculations can be used to study the crystallization kinetics if the temperature behavior in time at a certain location is known. With our choice of initial and boundary conditions, we found primarily crystals of the α′-phase.

Complex fluid–fluid interfaces determine for a large part the macroscopic material properties of foams and emulsions that appear in applications such as food, materials processing and consumer care products. As a step towards predicting these properties, a 2D axisymmetric and a 3D finite element model have been developed for simulating the dynamics of a single Newtonian drop in a Newtonian matrix fluid with a rheologically complex sharp interface in between. Interfaces with constant interfacial tension and with viscous, elastic and viscoelastic extra interfacial stresses are considered. The model has been validated by means of the method of manufactured solutions and by comparison with results from other studies using different discretisation methods. Higher-order convergence in space and time is obtained, demonstrating correct implementation of our numerical methods. Benchmark solutions for a drop with a Kelvin–Voigt interface under simple shear flow have been provided. Compared to a viscous interface, the drop deformation becomes smaller and the drop becomes less oriented in the direction of flow if interfacial elasticity is added.

A 3D transient non-isothermal finite element code is developed to predict the extrudate shape of viscoelastic fluids emerging from an asymmetric keyhole shaped die. The corner-line method is used to model the movement of the free surfaces. The code is tested using two benchmark problems. First the corner-line method is tested using a trumpet shaped object in a 3D uniaxial extensional flow. Secondly, the implementation of the energy balance and the viscoelastic material behaviour is tested using a non-isothermal pipeflow. For both benchmark problems convergence was obtained, giving confidence that the 3D non-isothermal swell problem is correctly implemented. The influence of shear-thinning, elasticity and temperature on the shape of the extrudate is systematically studied. Results are shown for isothermal flows as well as for non-isothermal flows, with isothermal and non-isothermal die walls. Results for isothermal die walls show increasing extrudate swelling with increasing elasticity and that the swelling opposes extrudate bending. Shear-thinning on the other hand, opposes swelling, which initially promotes bending, but also flattens the asymmetric velocity profile, leading to less extrudate bending for high amounts of shear-thinning. Furthermore, extrudate bending was observed even for purely viscous, isothermal extrudates, suggesting that bending is caused by asymmetry in the viscous stresses. Extrudate swelling can be influenced by the wall temperature of the die and non-isothermal die walls can lead to a change in bending direction.