Gradient

Description

This treatment computes the gradient of variables located at nodes. It can also compute the curl, divergence, Qcriterion, Lambda2 and Lambdaci of vectors since they are closely related to the velocity gradient.

Two computing methods are available.

The first method (vtk = False) is based on the Green-Ostrogradski theorem to compute the gradient and the volume at the cell centers. i.e. \(\displaystyle V = \frac{1}{3} \oint_S C.n \, dS\) with V the volume, C one point on the surface S, and n the outward normal. It handles both structured and unstructured grids.

The second method (vtk = True) is based on the VTK library to compute the gradients at nodes. It handles only unstructured grids. Note that the antares.treatment.TreatmentUnstructure.TreatmentUnstructure may be used for structured grids.

Parameters

  • base: Base

    The input base on which the gradient will be computed. This base is modified in-place.

  • coordinates: list(str)

    The coordinates used for gradient.

  • variables: list(str), default= []

    The variables for which gradient will be computed. If not given, gradient will be computed on all variables except the coordinates.

  • vtk: bool, default= False

    If False, use the method based on Green-Ostrogradski theorem. if True, use the VTK library.

  • curl: list(list(str, list(str))), default= []

    The list of list of the name of the curl to compute and the list of the vector variables from which the curl will be computed.

    • example: [ [‘vorticity’, [‘U’, ‘V’, ‘W’] ] ]

    The variables vorticity_<coordinate_name> and vorticity_modulus will be computed. The vector variables must be ordered as the coordinates.

  • divergence: list(list(str, list(str))), default= []

    The list of list of the name of the divergence to compute and the list of the vector variables from which the divergence will be computed.

    • example: [ [‘dilatation’, [‘U’, ‘V’, ‘W’] ] ]

    The variable dilatation will be computed. The vector variables must be ordered as the coordinates.

  • Qcriterion: list(str), default= []

    The list of velocity variables to compute the variable Qcrit.

    • example: [‘U’, ‘V’, ‘W’]

    The vector variables must be ordered as the coordinates.

  • Lambda2: list(str), default= []

    The list of velocity variables to compute the variable Lambda2.

    • example: [‘U’, ‘V’, ‘W’]

    The vector variables must be ordered as the coordinates.

  • L2_eig: bool, default= False

    If True, the eigenvectors are saved in the variables: Lambda2_V1y, Lambda2_V1z, Lambda2_V2x, Lambda2_V2y, Lambda2_V2z

  • Lambdaci: list(str), default= []

    The list of velocity variables to compute the variable Lambdaci.

    • example: [‘U’, ‘V’, ‘W’]

    The vector variables must be ordered as the coordinates.

  • Lci_eig: bool, default= False

    If True, the eigenvectors are saved in the variables: Lambdaci_Vx, Lambdaci_Vy, Lambdaci_Vz

Preconditions

Zones can be either structured or unstructured without vtk. Zones must be unstructured with vtk.

Zones can contain multiple instants and a shared instant.

Gradient computation with the first method has only been implemented in 3D, and 2D with triangles.

Coordinates and variables must be available at nodes.

The computations of curl, divergence, Qcriterion, Lambda2 and Lambdaci assume 3D cartesian coordinates.

The computations of Lambda2 and Lambdaci may not work with vtk = True.

Postconditions

The input base is modified in-place.

If the input variables are shared, then the output gradients are shared.

The computed gradient is either at cells or at nodes depending on the method.

If curl, divergence, Qcriterion, Lambda2 and Lambdaci of vectors are computed, intermediate computed gradients are kept in the base.

The computed gradient are added to the input base with the names ‘grad_<var_name>_<coordinate_name>’.

The volume variable is added to the input base with the name (‘cell_volume’, ‘cell’).

The curl variables are added to the input base with the names ‘<curl_name>_<coordinate_name>’ and ‘<curl_name>_modulus’

The divergence \(\nabla \cdot \vec{U} =\frac{\partial U_{x}}{\partial x} + \frac{\partial U_{y}}{\partial y} + \frac{\partial U_{z}}{\partial z}\) is added to the input base with the names ‘<div_name>’.

The Q-criterion \(Q=0.5({||\Omega||}^2_F - {||S||}^2_F)\) is added to the input base with the name ‘Qcrit’. \(\nabla \vec{U} = S + \Omega\) is the velocity gradient tensor with its symmetric part, the strain rate tensor \(S= 0.5(\nabla \vec{U} + (\nabla \vec{U})^T)\), and its antisymmetric part, the vorticity tensor \(\Omega = 0.5(\nabla \vec{U} - (\nabla \vec{U})^T)\). The principal invariants of the rank-two tensor \(\nabla \vec{U}\) of dimension three are the solutions of the cubic characteristic polynomial \(det(\nabla \vec{U} - \lambda I) = 0\), or \(\lambda^3 + P \lambda^2 + Q \lambda +R = 0\). The Q-criterion is the second invariant, \(0.5( (tr(\nabla \vec{U}))^2 -tr((\nabla \vec{U})^2) )\) where \(tr\) is the trace. A vortex exists where \(Q>0\).

The \(\lambda_2\)-criterion is added to the input base with the name ‘Lambda2’. The three eigenvalues of \(\Omega^2 + S^2\) are computed. \(\lambda_2\) is defined as the second eigenvalue. The two first eigenvectors ‘Lambda2_V1x’, ‘Lambda2_V1y’, ‘Lambda2_V1z’, ‘Lambda2_V2x’, ‘Lambda2_V2y’ and ‘Lambda2_V2z’ are added to the input base if L2_eig is True.

  1. Jeong and F. Hussain. On the identification of a vortex. Journal of Fluid Mechanics, 285:69-94, (1995).

The Lambdaci criterion is added to the input base with the name ‘Lambdaci’. Its eigenvectors ‘Lambdaci_Vx’, ‘Lambdaci_Vy’ and ‘Lambdaci_Vz’ are added to the input base if Lci_eig is True.

  1. Zhou, R. J. Adrian, S. Balachandar, and T. M. Kendall. Mechanisms for generating coherent packets of hairpin vortices in channel flow. Journal of Fluid Mechanics, 387:353-396, (1999).

Example

The following example shows gradient computations.

import antares
myt = antares.Treatment('gradient')
myt['base'] = base
myt['coordinates'] = ['x', 'y', 'z']
myt['variables'] = ['Density', 'Pressure']
myt['vtk'] = False
myt['curl'] = [['vorticity', ['U', 'V', 'W']]]
myt['divergence'] = [['dilatation', ['U', 'V', 'W']]]
myt['Qcrit'] = ['U', 'V', 'W']
myt['Lambda2'] = ['U', 'V', 'W']
myt['L2_eig'] = True
myt['Lambdaci'] = ['U', 'V', 'W']
myt['Lci_eig'] = True
myt.execute()

Main functions

class antares.treatment.TreatmentGradient.TreatmentGradient
execute()

Execute the treatment.

Returns:

the base containing the initial variables and the newly computed gradient variables

Return type:

Base

Example

"""
This example shows how to compute the gradient of variables on a 3D base.
Note that the gradient vector computed for each variable is stored
at cell centers in the input base.
"""
import os
if not os.path.isdir('OUTPUT'):
    os.makedirs('OUTPUT')

from antares import Reader, Treatment, Writer

# ------------------
# Reading the files
# ------------------
reader = Reader('bin_tp')
reader['filename'] = os.path.join('..', 'data', 'ROTOR37', 'GENERIC', 'flow_<zone>_ite0.dat')
ini_base = reader.read()

# ---------------------
# Gradient computation
# ---------------------
treatment = Treatment('gradient')
treatment['base'] = ini_base
result = treatment.execute()

result.compute_cell_volume()
# -------------------
# Writing the result
# -------------------
writer = Writer('bin_tp')
writer['filename'] = os.path.join('OUTPUT', 'ex_gradients.plt')
writer['base'] = result
writer.dump()

"""
This example shows how to use the Treatment gradient on a poiseuille flow.
Gradients, vorticity, dilatation and stress tensors are computed.
"""
import os
import antares
import numpy as np

OUTPUT = 'OUTPUT'
if not os.path.isdir(OUTPUT):
    os.makedirs(OUTPUT)

# Define case
Re = 100000.                  # [-] Reynolds number
Mu = 8.9e-4                   # [Pa.s] Dynamic Viscosity
Rho = 997.                    # [Kg.m-3] Density
H = 0.1                       # [m] Height
L = 0.1                       # [m] Length
l = L*0.01                    # [m] thickness
Nx, Ny, Nz = 100, 100, 2      # [-] Number of nodes
P0 = 0.                       # [Pa] Ref Pressure

# Compute some values
U_mean = Re*Mu/H/Rho
U_max = 1.5*U_mean
Q = U_mean*H
dpdx = -Q/(H**3)*12.*Mu
Tau_w = -H/2.*dpdx

# Define geometry
xmin, ymin, zmin = 0., 0., 0.
xmax, ymax, zmax = L, H, l

x = np.linspace(xmin, xmax, Nx)
y = np.linspace(ymin, ymax, Ny)
z = np.linspace(zmin, zmax, Nz)

X, Y, Z = np.meshgrid(x, y, z)

# Define velocity and pressure field
U = -H**2/2./Mu*dpdx*(Y/H*(1-Y/H))
V = np.zeros_like(U)
W = np.zeros_like(V)
Pressure = P0 + X*dpdx

# Initialisation
base = antares.Base()
zone = antares.Zone()
instant = antares.Instant()
base['zone_0'] = zone
base['zone_0'].shared[('x', 'node')] = X
base['zone_0'].shared[('y', 'node')] = Y
base['zone_0'].shared[('z', 'node')] = Z
base['zone_0']['instant_0'] = instant
base['zone_0']['instant_0'][('Pressure', 'node')] = Pressure
base['zone_0']['instant_0'][('U', 'node')] = U
base['zone_0']['instant_0'][('V', 'node')] = V
base['zone_0']['instant_0'][('W', 'node')] = W

# Use vtk or not
use_vtk = False

# Unstructure base if vtk
if use_vtk:
    t = antares.Treatment('unstructure')
    t['base'] = base
    base = t.execute()
    loc = 'node'
else:
    loc = 'cell'
    base.node_to_cell(variables=['Pressure'])

# Compute gradients
treatment = antares.Treatment('Gradient')
treatment['base'] = base
treatment['coordinates'] = ['x', 'y', 'z']
treatment['variables'] = ['Pressure']
treatment['vtk'] = use_vtk
treatment['curl'] = [['vorticity', ['U', 'V', 'W']]]
treatment['divergence'] = [['dilatation', ['U', 'V', 'W']]]
treatment['Qcrit'] = ['U', 'V', 'W']
treatment['Lambda2'] = ['U', 'V', 'W']
treatment['L2_eig'] = True
treatment['Lambdaci'] = ['U', 'V', 'W']
treatment['Lci_eig'] = True
base = treatment.execute()

# Compute stress tensors
base.set_computer_model('internal')

base.set_formula('Mu={:f}'.format(Mu))

base.compute('tau_xx = -2. / 3. * Mu * dilatation + Mu * (2 * grad_U_x)', location=loc)
base.compute('tau_yy = -2. / 3. * Mu * dilatation + Mu * (2 * grad_V_y)', location=loc)
base.compute('tau_zz = -2. / 3. * Mu * dilatation + Mu * (2 * grad_W_z)', location=loc)
base.compute('tau_xy = Mu * (grad_V_x + grad_U_y)', location=loc)
base.compute('tau_xz = Mu * (grad_W_x + grad_U_z)', location=loc)
base.compute('tau_yz = Mu * (grad_W_y + grad_V_z)', location=loc)

base.compute('sigma_xx = Pressure + tau_xx', location=loc)
base.compute('sigma_yy = Pressure + tau_yy', location=loc)
base.compute('sigma_zz = Pressure + tau_zz', location=loc)
base.compute('sigma_xy = tau_xy', location=loc)
base.compute('sigma_xz = tau_xz', location=loc)
base.compute('sigma_yz = tau_yz', location=loc)

# Write output
writer = antares.Writer('bin_vtk')
writer['filename'] = os.path.join(OUTPUT, 'results')
writer['base'] = base
writer.dump()