VISCOSITY_MODEL

Specifies a viscosity model for the flow equation.

Type

AcuSolve Command

Syntax

VISCOSITY_MODEL("name") {parameters...}

Qualifier

User-given name.

Parameters

type (enumerated) [=none]

Type of the viscosity.

constant or const: Constant viscosity. Requires viscosity.
ramped or ramp: Ramp the viscosity to a constant value. Requires viscosity.
power_law: Non-Newtonian power law viscosity model. Requires power_law_viscosity, power_law_time_constant, power_law_index and power_law_lower_strain_rate.
bingham: Non-Newtonian Bingham viscosity model. Requires bingham_viscosity, bingham_yield_stress, bingham_stress_growth_exponent, bingham_time_constant, bingham_index and bingham_infinite_shear_viscosity.
carreau: Non-Newtonian Carreau-Yasuda viscosity model. Requires carreau_zero_shear_viscosity, carreau_infinite_shear_viscosity, carreau_time_constant, carreau_index and carreau_transition_index.
gas_kinetic: Gas kinetic viscosity model for gases only. Requires GAS_KINETIC_MODEL.
sutherland: Sutherland temperature dependent viscosity model. Requires sutherland_reference_viscosity, sutherland_reference_temperature, sutherland_constant.
piecewise_linear or linear: Piecewise linear curve fit. Requires curve_fit_values and curve_fit_variable.
cubic_spline or spline: Cubic spline curve fit. Requires curve_fit_values and curve_fit_variable.
user_function or user: User-defined function. Requires user_function, user_values and user_strings.

viscosity or visc (real) >=0 [=0]

Constant value of the dynamic viscosity. Used with constant and ramped types.

power_law_viscosity (real) >=0 [=0]

Power law model reference viscosity. Used with power_law type.

power_law_time_constant (real) >=0 [=0]

Power law model time constant. Used with power_law type.

power_law_index (real) [=1]

Power law model power index. Used with power_law type.

power_law_lower_strain_rate (real) >=0 [=0]

Power law model lower limit of strain rate. Used with power_law type.

bingham_viscosity (real) >0 [=1]

Bingham model reference viscosity. Used with bingham type.

bingham_yield_stress (real) >=0 [=0]

Bingham model yield stress. Used with bingham type.

bingham_stress_growth_exponent (real) >0 [=500]

Bingham model unyielded stress region growth rate exponent. Used with bingham type.

bingham_time_constant (real) >0 [=1]

Bingham model power law time constant. Used with bingham type.

bingham_index (real) [=1]

Bingham model power law index. Used with bingham type.

bingham_infinite_shear_viscosity (real) [=0]

Bingham model minimum viscosity. Used with bingham type.

carreau_zero_shear_viscosity (real) >=0 [=0]

Carreau model viscosity at zero shear rate. Used with carreau type.

carreau_infinite_shear_viscosity (real) >=0 [=0]

Carreau model viscosity at infinite shear rate. Used with carreau type.

carreau_time_constant (real) >=0 [=0]

Carreau model power law time constant. Used with carreau type.

carreau_index (real) [=1]

Carreau model power law index. Used with carreau type.

carreau_transition_index (real) >0 [=2]

Carreau model transition index. Used with carreau type.

carreau_temperature_dependence_type (enumerated) [=none]

Types to add the temperature dependence terms to the Carreau Yasuda viscosity material to account for shear heating due to viscous dissipation. Used with carreau type when temperature = advective_diffusive under the EQUATION command and viscous_heating = on in ELEMENT_SET.

Arrhenius: Arrhenius model. Requires activation energy.
William_Landel_Ferry: William Landel Ferry (WLF) model. Requires William_Landel_Ferry_coefficient_C1, William_Landel_Ferry_coefficient_C2, glass_transition_temperature.
exponential: Exponential model. Requires exponential_reference_temperature and exponential_beta_coefficient.

activation_energy (real) >0 [=16628]

Arrhenius model constant to define the temperature sensitivity. Used with carreau type and carreau_temperature_dependence_type = Arrhenius.

William_Landel_Ferry_coefficient_C1 (real) > 0 [=17.44]

William-Landel-Ferry (WLF) model constant. Used with carreau type and carreau_temperature_dependence_type = William_Landel_Ferry.

William_Landel_Ferry_coefficient_C2 (real) > 0 [=51.6]

William-Landel-Ferry (WLF) model constant. Used with carreau type and carreau_temperature_dependence_type = William_Landel_Ferry.

glass_transition_temperature (real) > 0 [=320]

Glass transition temperature used in the William-Landel-Ferry (WLF) model. Used with carreau type and carreau_temperature_dependence_type = William_Landel_Ferry.

exponential_reference_temperature (real) > 0 [=350]

Reference temperature defined in the exponential model. Used with carreau type and carreau_temperature_dependence_type = exponential.

exponential_beta_coefficient (real) > 0 [=0.005]

Exponential model constant. Used with carreau type and carreau_temperature_dependence_type = exponential.

sutherland_reference_viscosity (real) > = 0 [=0.01716]

Sutherland model reference viscosity at sutherland_reference_temperature. Used with sutherland type.

sutherland_reference_temperature (real) > 0 [=273]

Sutherland model reference temperature corresponding to sutherland_reference_viscosity. Used with sutherland type.

sutherland_constant (real) > 0 [=111]

Sutherland model constant. Used with sutherland type.

curve_fit_values or curve_values array [={0,0}]

A two-column array of independent-variable/viscosity data values. Used with piecewise_linear and cubic_spline types.

curve_fit_variable or curve_var (enumerated) [=temperature]

Independent variable of the curve fit. Used with piecewise_linear and cubic_spline types.

x_coordinate or xcrd: X-component of coordinates.
y_coordinate or ycrd: Y-component of coordinates.
z_coordinate or zcrd: Z-component of coordinates.
x_reference_coordinate or xrefcrd: X-component of reference coordinates.
y_reference_coordinate or yrefcrd: Y-component of reference coordinates.
z_reference_coordinate or zrefcrd: Z-component of reference coordinates.
pressure or pres: Pressure.
temperature or temp: Temperature.
species_1 or spec1: Species 1.
species_2 or spec2: Species 2.
species_3 or spec3: Species 3.
species_4 or spec4: Species 4.
species_5 or spec5: Species 5.
species_6 or spec6: Species 6.
species_7 or spec7: Species 7.
species_8 or spec8: Species 8.
species_9 or spec9: Species 9.
strain_rate_second_invariant or strain_i2: Second invariant of the strain rate tensor.

user_function or user (string) [no default]

Name of the user-defined function. Used with user_function type.

user_values (array) [={}]

Array of values to be passed to the user-defined function. Used with user_function type.

user_strings (list) [={}]

Array of strings to be passed to the user-defined function. Used with user_function type.

multiplier_function (string) [=none]

User-given name of the multiplier function for scaling the viscosity. If none, no scaling is performed.

Description

This command specifies a dynamic viscosity model for the flow (momentum) equation. This model is only applicable to fluid element sets.

VISCOSITY_MODEL commands are referenced by MATERIAL_MODEL commands, which in turn are referenced by ELEMENT_SET commands:

VISCOSITY_MODEL( "my viscosity model" ) {
   type                                = constant 
   viscosity                           = 1.781e-5
 }
 MATERIAL_MODEL( "my material model" ) {
   viscosity_model                     = "my viscosity model"
   ...
 }
 ELEMENT_SET( "fluid elements" ) {
   material_model                      = "my material model"
   medium                              = fluid
   ...
}

The viscous stress tensor

is given by:(1)

where $μ$ is the dynamic viscosity and $\nabla_{s} u$ is the symmetric part of the strain-rate tensor.

A constant type applies a spatially constant viscosity model, as in the above example. This is the simplest case of a Newtonian model for the stress tensor.

The power_law type defines a non-Newtonian viscosity model in which the dynamic-viscosity is a function of the velocity field:(2)

μ = μ_{0} {[η^{2} \max (I_{2}, S_{\min}^{2})]}^{(n - 1) / 2}

where

μ_{0}

is the reference viscosity, given by power_law_viscosity;

η

is the time constant, given by power_law_time_constant; n is the power index, given by power_law_index; S_min is the lower limit of the strain rate, given by power_law_lower_strain_rate; and I₂ is the second invariant of the strain rate tensor:(3)

I_{2} = 2 \nabla_{s} u : \nabla_{s} u

For example,

VISCOSITY_MODEL( "my viscosity model" ) {
    type                 = constant
    viscosity            = 1.781e-5
}
MATERIAL_MODEL( "my material model" ) {
    viscosity_model      = "my viscosity model"
    ...
}
ELEMENT_SET( "fluid elements" ) {
    material_model       = "my material model"
    medium               = fluid
    ...
}

The bingham type defines another non-Newtonian viscosity model in which Papanastasiou's unyielded stress region and power law modifications are added to the basic Bingham viscosity model:(4)

μ = μ_{\infty} + (μ_{0} - μ_{\infty}) {(η g)}^{n - 1} + \frac{τ_{y}}{g} (1 - e^{- M g})

where

μ_{0}

is the reference viscosity, given by bingham_viscosity;

η

is the time constant, given by bingham_time_constant; n is the power index, given by bingham_index;

τ_{y}

is the yield stress, given by bingham_yield_stress; M is the unyielded stress region growth rate exponent, given by bingham_stress_growth_exponent;

μ

_∞ is the viscosity at infinite shear rate (or equivalently the minimum viscosity), given by bingham_infinite_shear_viscosity; and g is the square root of the second invariant of the strain rate tensor:(5)

g = \sqrt{2 \nabla_{s} u : \nabla_{s} u}

This may be used as a basic Bingham viscosity model,

μ = μ_{0} + τ_{y} / g

, as follows:

VISCOSITY_MODEL( "Basic Bingham Law" ) {
    type                     = bingham
    bingham_viscosity        = 0.001
    bingham_yield_stress     = 10
}

Papanastasiou's unyielded stress region growth rate exponent modification limits the viscosity at low shear rates to

μ = μ_{0} + τ_{y} M

and provides for a smooth transition from the unyielded to the yielded stress regions. For example:

VISCOSITY_MODEL( "Bingham Law with Papanastasiou exponent" ) {
    type                           = bingham
    bingham_viscosity              = 0.001
    bingham_yield_stress           = 10
    bingham_stress_growth_exponent = 1000
}

To add a power law model modification to the viscosity model:

VISCOSITY_MODEL( "Bingham Law with Papanastasiou exponent and power law" ) {
    type                           = bingham
    bingham_viscosity              = 0.001
    bingham_yield_stress           = 10
    bingham_stress_growth_exponent = 1000
    bingham_time_constant          = 1
    bingham_index                  = 0.75
}

The carreau type defines a third non-Newtonian viscosity model:(6)

μ = μ_{\infty} + (μ_{0} - μ_{\infty}) {(1 + {(η g)}^{a})}^{(η - 1) / a}

where $μ$ ₀ is the viscosity at zero shear rate, given by carreau_zero_shear_viscosity; $μ$ _∞ is the viscosity at infinite shear rate, given by carreau_infinite_shear_viscosity; $η$ is the time constant, given by carreau_time_constant; n is the power index, given by carreau_index; a is the transition index, given by carreau_transition_index; and g is the square root of the second invariant of the strain rate tensor (see above). In this form this model is known as the Carreau-Yasuda model. If a = 2, then the classic Carreau model is recovered.

For example,

VISCOSITY_MODEL( "Carreau Model" ) {
    type                                = carreau
    carreau_zero_shear_viscosity        = 681.8
    carreau_infinite_shear_viscosity    = 0.17
    carreau_time_constant               = 31.7
    carreau_index                       = 0.11
    carreau_transition_index            = 2
}

The sutherland type defines a temperature dependent viscosity model according to the following equation:(7)

μ = μ_{0} {(\frac{T}{T_{0}})}^{3 / 2} \frac{T_{0} + S}{T + S}

where $μ_{0}$ is the sutherland_reference_viscosity at the specified sutherland_reference_temperature ( $T_{0}$ ), $T$ is the model temperature and $S$ is the sutherland_constant.

For example,

VISCOSITY_MODEL( "Sutherland Model" ) {
type = sutherland
sutherland_reference_viscosity = 1.716e-5
sutherland_reference_temperature = 273.0
sutherland_constant = 111
}

A ramped viscosity model is the same as the constant one, except that the viscosity is ramped down from 1000 times the viscosity value at time step one to the viscosity value at time step ten. For example,

VISCOSITY_MODEL( "ramped viscosity model" ) {
    type        = ramped
    viscosity   = 1.781e-5
}

starts with a viscosity of 1.781 x 10^-8 at time step one and reduces the viscosity until time step 10, at which time a viscosity of 1.781 x 10^-5 is used.

Viscosity models of types piecewise_linear and cubic_spline may be used to define viscosity as a function of a single independent variable. For example,

VISCOSITY_MODEL( "curve fit viscosity model" ) {
    type                 = piecewise_linear
    curve_fit_values     = { 273, 1.72e-5 ;
                             323, 1.95e-5 ;
                             373, 2.17e-5 ;
                             423, 2.38e-5 ; }
    curve_fit_variable   = temperature
}

defines viscosity as a function of temperature. In this case, the problem must contain a temperature equation; see the EQUATION command. The curve_fit_values parameter is a two-column array corresponding to the independent variable and the viscosity values. The independent variable values must be in ascending order. The limit point values of the curve fit are used when curve_fit_variable falls outside of the curve fit limits.

The curve_fit_values data may be read from a file. For the above example, the curve fit values may be placed in a file, such as viscosity.fit:

273       1.72e-5
323       1.95e-5
373       2.17e-5
423       2.38e-5

and read by:

VISCOSITY_MODEL( "curve fit viscosity model" ) {
    type                = piecewise_linear
    curve_fit_values    = Read( "viscosity.fit" )
    curve_fit_variable  = temperature
}

A viscosity of type user_function may be used to model more complex behaviors; the see AcuSolve User-Defined Functions Manual for a detailed description of user-defined functions.

For example, consider a simplified version of the non-Newtonian power law viscosity model discussed above. The input command may be given by:

VISCOSITY_MODEL( "UDF viscosity model" ) {
    type           = user_function
    user_function  = "usrViscosityExample"
    user_values    = { 9000,   # reference viscosity
                         12.5,   # time constant
                          0.48 }  # power index
}

where the user-defined function "usrViscosityExample" may be implemented as follows:

#include "acusim.h"
#include "udf.h"
UDF_PROTOTYPE( usrViscosityExample ) ;                                    /* function prototype               */
Void usrViscosityExample (
    UdfHd         udfHd,                                                  /* Opaque handle for accessing data */
	Real*         outVec,                                                 /* Output vector                    */
	Integer       nItems,                                                 /* Number of elements               */
	Integer       vecDim                                                  /* = 1                              */
) {
    Integer       elem ;                                                  /* an element counter               */
    Real          exp ;                                                   /* exponent ( (index-1)/2 )         */
    Real          index ;                                                 /* power index                      */
    Real          invar2 ;                                                /* second invariant of strain rate  */
    Real          tConst ;                                                /* time constant                    */
    Real          tConst2 ;                                               /* scaling factor (tConst ^ 2)      */
    Real          visc0 ;                                                 /* reference viscosity              */
    Real*         gradV ;                                                 /* gradient of velocity             */
    Real*         u11 ;                                                   /* partial u1 / partial x1          */
    Real*         u12 ;                                                   /* partial u1 / partial x2          */
    Real*         u13 ;                                                   /* partial u1 / partial x3          */
    Real*         u21 ;                                                   /* partial u2 / partial x1          */
    Real*         u22 ;                                                   /* partial u2 / partial x2          */
    Real*         u23 ;                                                   /* partial u2 / partial x3          */
    Real*         u31 ;                                                   /* partial u3 / partial x1          */
    Real*         u32 ;                                                   /* partial u3 / partial x2          */
    Real*         u33 ;                                                   /* partial u3 / partial x3          */
    Real*         usrVals ;                                               /* user values                      */
udfCheckNumUsrVals( udfHd, 3 ) ;                                          /* check for error                  */
    usrVals      = udfGetUsrVals( udfHd ) ;                               /* get the user vals                */
    visc0        = usrVals[0] ;                                           /* reference visc.                  */
    tConst       = usrVals[1] ;                                           /* time constant                    */
    index        = usrVals[2] ;                                           /* index                            */
    tConst2      = tConst * tConst ;                                      /* for efficiency                   */
    exp          = (index - 1) / 2 ;                                      /* for efficiency                   */
    gradV        = udfGetElmData( udfHd, UDF_ELM_GRAD_VELOCITY ) ;        /* get the grad-vel.                */
    u11          = &gradV[0*nItems] ;                                     /* localize u11                     */
    u12          = &gradV[1*nItems] ;                                     /* localize u12                     */
    u13          = &gradV[2*nItems] ;                                     /* localize u13                     */
    u21          = &gradV[3*nItems] ;                                     /* localize u21                     */
    u22          = &gradV[4*nItems] ;                                     /* localize u22                     */
    u23          = &gradV[5*nItems] ;                                     /* localize u23                     */
    u31          = &gradV[6*nItems] ;                                     /* localize u31                     */
    u32          = &gradV[7*nItems] ;                                     /* localize u32                     */
    u33          = &gradV[8*nItems] ;                                     /* localize u33                     */
    for ( elem = 0 ; elem < nItems ; elem++ ) {
        invar2       = 2.
                     * ( u11[elem] * u11[elem]
                     + u22[elem] * u22[elem]
                     + u33[elem] * u33[elem]
                     + (u12[elem]+u21[elem]) * (u12[elem]+u21[elem])
                     + (u23[elem]+u32[elem]) * (u23[elem]+u32[elem])
                     + (u31[elem]+u13[elem]) * (u31[elem]+u13[elem]) ;

         outVec[elem] = visc0 ;
         if ( tConst2 * invar2 > ) {
             outVec[elem] = visc0 * pow( tConst2 * invar2, exp ) ;
         }
    }
} /* end of usrViscosityExample() */

The dimension of the returned viscosity vector, outVec, is the number of elements.

The multiplier_function parameter may be used to uniformly scale the viscosity values. The value of this parameter refers to the user-given name of a MULTIPLIER_FUNCTION command in the input file. For example, a ramped viscosity may be specified by:

VISCOSITY_MODEL( "similar to ramped viscosity model" ) {
    type                   = constant
    viscosity              = 1.781e-5
    multiplier_funtion     = "ramped"
}
MULTIPLIER_FUNCTION( "ramped" ) {
    type                   = piecewise_log_linear
    curve_fit_values       = { 1, 1000 ; 10, 1 }
    curve_fit_variable     = time_step
}

In areas of turbulent flow the molecular viscosity $μ$ is replaced by the total viscosity $μ$ + $μ$ _t in the definition of the viscous stress tensor. The turbulent viscosity $μ$ _t is given by the turbulence model in use; see the EQUATION command for details.