Learning Center

Two-Phase Flow Modeling Guidelines

---

When modeling fluid flow in COMSOL Multiphysics^® using any of the physics interfaces that fall under the Two-Phase Flow, Level Set or Two-Phase Flow, Phase Field application areas, there are several guidelines you can follow to ensure that your model converges and does so both in a reasonable time and to a reasonable result. Here, you will find guidance for solving these types of fluid flow problems.

The Two-Phase Flow, Level Set and Two-Phase Flow, Phase Field branches of physics interfaces in the Select Physics window.

Background

Fluid flow is governed by the Navier–Stokes equations, which solve for fluid velocity and pressure fields everywhere within the modeling domain. The fluid properties (viscosity and density) are either constant or vary relatively smoothly in space as a function of temperature, pressure, shear rate, etc. However, if we want to model two immiscible fluids, then the fluid properties will vary significantly across the interface between the two fluids. The surface tension effects can become important, as well as the effect of contact angles at wetted walls.

To model this in COMSOL Multiphysics^®, we can use either the level set or phase field method. Both of these methods introduce an additional scalar field (the level set or phase field function) to the modeling domain. These scalar fields vary smoothly between -1 and +1 everywhere and are used to define the fluid viscosity and density in the governing Navier–Stokes equations. The transition (where the fields vary from 0 to 1 for the level set implementation and -1 to +1 for phase field) is quite abrupt in space, thus giving a good resolution of the two phases. The governing equations for the level set and phase field methods are a type of convection–diffusion equation, with the advective velocity coming from the Navier–Stokes equations. Such equations are, however, quite numerically challenging to solve due to the combination of the significant advective term, abrupt transition in the field, and strong coupling to the Navier–Stokes equations.

The level set function varies smoothly and defines the interface between the fluids. The boundary is advected by the fluid flow (back arrows.)

Options and General Guidelines for Modeling Two-Phase Flow

If you want to model two-phase flow that does not involve any topology changes (that is, no droplet breakup, free surface formation, etc.), and if the free surface does not undergo significant shape changes, you may be able to use the Laminar Two Phase Flow, Moving Mesh physics interface. This formulation, if it can be used, will require less computational resources. Regardless of the physics interface used, there are general guidelines that can be used when modeling this type of flow. These are listed below and organized by the modeling workflow step.

Geometry

If at all possible, begin your modeling in 2D or 2D axisymmetry. Such models will be significantly faster to solve than 3D models and serve as a useful test for trying out different model settings that will be applicable for 3D models.

If reasonable, avoid sharp edges and corners at any boundaries where you expect the interface between the two fluids to traverse. Introduce a fillet in the geometry so that there is a smooth transition between boundaries.

Physics

All fluid inlets should be in regions where there is a single phase entering the modeling domain, and the phase entering at that inlet should match the initial phase inside the adjacent domain. You may need to alter the geometry slightly to achieve this.

For both the level set and phase field approaches, there is a Two-Phase Flow, Level Set and Two-Phase Flow, Phase Field feature under the Multiphysics node. Make sure to choose the appropriate materials from the drop-down lists for Fluid 1 and Fluid 2 in the settings.

Within the Laminar Flow physics interface, make sure to use the incompressible formulation. You should also make sure that all specified velocities or pressures at inlets ramp smoothly up from consistent initial conditions, which is described in more detail here. For models involving gravity, make sure to enable the Include gravity option under the settings for the Laminar Flow physics interface, instead of using a Volume Force condition. This is so that the hydrostatic pressure is automatically taken into account under the Initial Values node and under other pressure conditions.

The Laminar Flow physics interface with the recommended settings, stated above.

If modeling bubbles with surface tension, make sure to define the additional initial pressure inside of the bubble that balances the surface tension. A reasonable value for this initial pressure can be computed from the Young–Laplace equation.

When using either the Level Set or Phase Field physics interfaces, there is a Level Set Model or Phase Field Model feature, which is used to define the Parameter controlling interface thickness. If this value is too small, it can lead to numerical instabilities. If it is too large, the interface movement is not captured correctly. The parameter controlling the interface thickness is set proportional to the maximum element size: hmax.

The Settings window for the Phase Field Model feature.

When using the Phase Field physics interface, the Phase Field Model feature also defines the Mobility tuning parameter. Again, too small a value leads to numerical instabilities, and too large a value will not capture the interface movement correctly. A good initial estimate for the Mobility tuning parameter is:

2*umax/(3*sqrt(2)*sigma)*(hmax/epsilon)

where umax is the expected maximum velocity magnitude; sigma is the surface tension coefficient; hmax is the value of the parameter controlling the maximum element size; and epsilon is the value of the Parameter controlling interface thickness.

When using the Level Set physics interface, the Level Set Model feature also defines the Reinitialization parameter. If the reinitialization parameter is too small, it can lead to numerical instabilities. On the other hand, if it is too big, then the interface movement is not captured correctly. A good estimate for the reinitialization parameter is the maximum expected velocity magnitude.

The Settings window for the Level Set Model feature.

When using the level set approach, there is a Wetted Wall multiphysics coupling boundary feature, and when using the phase field approach, there is a Wetted Wall node that specifies a Contact angle. This contact angle should be in agreement with the actual initial angle of the two phases, as defined by the geometry and the Phase Field/Level Set > Initial Values feature. If the contact angle as defined by the initial values is different, then the specified Contact angle should be ramped up in a way similar to the velocity ramping (see: Solving time dependent CFD simulations).

In the settings for the Laminar Flow, Level Set, and Phase Field physics interfaces, there is a Discretization section where you can change the element order. The default discretization is P1+P1 for the Laminar Flow physics interface, and Linear for the Level Set and Phase Field physics interfaces. Higher-order discretization will be significantly more computationally intensive and is not usually recommended. Instead, you can refine the mesh.

Mesh

The mesh should have roughly equal size throughout the modeling domain, and the mesh elements should be as isotropic as possible; that is, the elements should not be overly stretched or compressed in one direction. The mesh size must be small enough to give a good resolution of the fluid interface between phases. A good rule of thumb is to begin with a global element size that is one-tenth of the smallest expected droplet size, and to refine the mesh from there.

Define a parameter hmax that is used within the Mesh > Size node's settings to define the Maximum element size in all domains. Once you have arrived at a converged solution on a particular mesh, repeat the analysis with a smaller mesh size. Once the solution does not change significantly with mesh refinement, the solution can be converged with respect to the mesh.

If it is known that one subdomain of the entire modeling domain will only ever contain the same fluid phase, then a coarser mesh can be used in that subdomain.

There is not a significant difference between structured and unstructured meshes for modeling these types of problems.

Study

The study always consists of two steps. First, a Phase initialization step initializes the level set or phase field variable such that it is smoothly varying everywhere. Next, a Time Dependent step simultaneously solves the Navier–Stokes equations and the level set or phase field equation.

Within the Time Dependent study step settings, there is a Tolerance field that defaults to Physics controlled. This can be changed to User controlled, and tighter solver relative tolerances should be investigated to confirm that the solution is converged with respect to relative tolerance. That is, once a solution is computed, repeat the solution with a finer solver tolerance. Once the solution does not change significantly with solver tolerance as well as mesh refinement, the solution can be considered converged. You should not study finer solver tolerances before addressing all of the other above guidelines.

The Settings window for the Time Dependent study step where the Tolerance is set to User controlled.