API of the TensND library

Basis(v::AbstractMatrix{T}, ::Val{:cov})
Basis{dim, T<:Number}()
Basis(θ::T<:Number, ϕ::T<:Number, ψ::T<:Number)

Basis built from a square matrix v where columns correspond either to

• primal vectors ie eᵢ=v[:,i] if var=:cov as by default
• dual vectors ie eⁱ=v[:,i] if var=:cont.

Basis without any argument refers to the canonical basis (CanonicalBasis) in Rᵈⁱᵐ (by default dim=3 and T=Sym)

Basis can also be built from Euler angles (RotatedBasis) θ in 2D and (θ, ϕ, ψ) in 3D

The attributes of this object can be obtained by

• vecbasis(ℬ, :cov): square matrix defining the primal basis eᵢ=e[:,i]
• vecbasis(ℬ, :cont): square matrix defining the dual basis eⁱ=E[:,i]
• metric(ℬ, :cov): square matrix defining the covariant components of the metric tensor gᵢⱼ=eᵢ⋅eⱼ=g[i,j]
• metric(ℬ, :cont): square matrix defining the contravariant components of the metric tensor gⁱʲ=eⁱ⋅eʲ=gⁱʲ[i,j]

Examples

julia> v = Sym[1 0 0; 0 1 0; 0 1 1] ; ℬ = Basis(v)
Basis{3, Sym}
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 1  0  0
 0  1  0
 0  1  1
# dual basis: 3×3 Tensor{2, 3, Sym, 9}:
 1  0   0
 0  1  -1
 0  0   1
# covariant metric tensor: 3×3 SymmetricTensor{2, 3, Sym, 6}:
 1  0  0
 0  2  1
 0  1  1
# contravariant metric tensor: 3×3 SymmetricTensor{2, 3, Sym, 6}:
 1   0   0
 0   1  -1
 0  -1   2

julia> θ, ϕ, ψ = symbols("θ, ϕ, ψ", real = true) ; ℬʳ = Basis(θ, ϕ, ψ) ; display(vecbasis(ℬʳ, :cov))
3×3 Tensor{2, 3, Sym, 9}:
 -sin(ψ)⋅sin(ϕ) + cos(θ)⋅cos(ψ)⋅cos(ϕ)  -sin(ψ)⋅cos(θ)⋅cos(ϕ) - sin(ϕ)⋅cos(ψ)  sin(θ)⋅cos(ϕ)
  sin(ψ)⋅cos(ϕ) + sin(ϕ)⋅cos(θ)⋅cos(ψ)  -sin(ψ)⋅sin(ϕ)⋅cos(θ) + cos(ψ)⋅cos(ϕ)  sin(θ)⋅sin(ϕ)
                        -sin(θ)⋅cos(ψ)                          sin(θ)⋅sin(ψ)         cos(θ)

CanonicalBasis{dim, T}

Canonical basis of dimension dim (default: 3) and type T (default: Sym)

The attributes of this object can be obtained by

• vecbasis(ℬ, :cov): square matrix defining the primal basis eᵢ=e[:,i]=δᵢⱼ
• vecbasis(ℬ, :cont): square matrix defining the dual basis eⁱ=E[:,i]=δᵢⱼ
• metric(ℬ, :cov): square matrix defining the covariant components of the metric tensor gᵢⱼ=eᵢ⋅eⱼ=g[i,j]=δᵢⱼ
• metric(ℬ, :cont): square matrix defining the contravariant components of the metric tensor gⁱʲ=eⁱ⋅eʲ=gⁱʲ[i,j]=δᵢⱼ

Examples

julia> ℬ = CanonicalBasis()
CanonicalBasis{3, Sym}
# basis: 3×3 TensND.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1
# dual basis: 3×3 TensND.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1
# covariant metric tensor: 3×3 TensND.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1
# contravariant metric tensor: 3×3 TensND.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1

julia> ℬ₂ = CanonicalBasis{2, Float64}()
CanonicalBasis{2, Float64}
# basis: 2×2 TensND.LazyIdentity{2, Float64}:
 1.0  0.0
 0.0  1.0
# dual basis: 2×2 TensND.LazyIdentity{2, Float64}:
 1.0  0.0
 0.0  1.0
# covariant metric tensor: 2×2 TensND.LazyIdentity{2, Float64}:
 1.0  0.0
 0.0  1.0
# contravariant metric tensor: 2×2 TensND.LazyIdentity{2, Float64}:
 1.0  0.0
 0.0  1.0

CoorSystemNum{dim,T<:Real} <: AbstractCoorSystem{dim,T}

Numerical coordinate system for automatic-differentiation-based differential operators.

Stores three point-wise functions:

• χ_func(x) : Lamé coefficients at point x
• R_func(x) : rotation matrix (columns = unit vectors of the normalized basis) at x
• Γ_func(x) : Christoffel symbols at x → Array{T,3} with Γ[i,j,k] = Γᵏᵢⱼ

All three functions accept an AbstractVector of coordinates and are called at each evaluation point. Predefined constructors are available:

coorsys_cartesian_num(T=Float64)
coorsys_polar_num(T=Float64)
coorsys_cylindrical_num(T=Float64)
coorsys_spherical_num(T=Float64)

A generic constructor from an OM::Function (position vector) is also provided:

CoorSystemNum(OM::Function, dim::Integer, T::Type=Float64)

Unlike CoorSystemSym, the normalized basis, natural basis vectors, Lamé coefficients, and Christoffel symbols all depend on the evaluation point x₀ and are accessed via point-wise accessors normalized_basis(CS, x₀), natural_basis(CS, x₀), etc.

CoorSystemNum(OM_func::Function, dim::Integer, T::Type=Float64)

Generic constructor: build a CoorSystemNum from a position-vector function. OM_func(x) must return a Vector of Cartesian components given a coordinate vector x of length dim. Christoffel symbols are computed by nested automatic differentiation of the metric tensor g = J'J where J = ∂OM/∂x.

RotatedBasis(θ::T<:Number, ϕ::T<:Number, ψ::T<:Number)
RotatedBasis(θ::T<:Number)

Orthonormal basis of dimension dim (default: 3) and type T (default: Sym) built from Euler angles θ in 2D and (θ, ϕ, ψ) in 3D

Examples

julia> θ, ϕ, ψ = symbols("θ, ϕ, ψ", real = true) ; ℬʳ = RotatedBasis(θ, ϕ, ψ) ; display(vecbasis(ℬʳ, :cov))
3×3 Tensor{2, 3, Sym, 9}:
 -sin(ψ)⋅sin(ϕ) + cos(θ)⋅cos(ψ)⋅cos(ϕ)  -sin(ψ)⋅cos(θ)⋅cos(ϕ) - sin(ϕ)⋅cos(ψ)  sin(θ)⋅cos(ϕ)
  sin(ψ)⋅cos(ϕ) + sin(ϕ)⋅cos(θ)⋅cos(ψ)  -sin(ψ)⋅sin(ϕ)⋅cos(θ) + cos(ψ)⋅cos(ϕ)  sin(θ)⋅sin(ϕ)
                        -sin(θ)⋅cos(ψ)                          sin(θ)⋅sin(ψ)         cos(θ)

CoorSystemSym(OM::AbstractTens{1,dim,Sym},coords::NTuple{dim,Sym},bnorm::AbstractBasis{dim,Sym},χᵢ::NTuple{dim},
              tmp_coords::NTuple = (),params::NTuple = ();rules::Dict = Dict(),tmp_var::Dict = Dict(),to_coords::Dict = Dict()) where {dim}
CoorSystemSym(OM::AbstractTens{1,dim,Sym},coords::NTuple{dim,Sym},
              tmp_coords::NTuple = (),params::NTuple = ();rules::Dict = Dict(),tmp_var::Dict = Dict(),to_coords::Dict = Dict()) where {dim}

Define a new coordinate system either from

1. the position vector OM, the coordinates coords, the basis of unit vectors (𝐞ᵢ) bnorm and the Lamé coefficients χᵢ

   In this case the natural basis is formed by the vectors 𝐚ᵢ = χᵢ 𝐞ᵢ directly calculated from the input data.

2. or the position vector OM and the coordinates coords

   In this case the natural basis is formed by the vectors 𝐚ᵢ = ∂ᵢOM i.e. by the derivative of the position vector with respect to the iᵗʰ coordinate

Optional parameters can be provided:

• tmp_coords contains temporary variables depending on coordinates (in order to allow symbolic simplifications)
• params contains possible parameters involved in OM
• rules contains a Dict with substitution rules to facilitate the simplification of formulas
• tmp_var contains a Dict with substitution of coordinates by temporary variables
• to_coords indicates how to eliminate the temporary variables to come back to the actual coordinates before derivation for Examples

Examples

julia> ϕ, p = symbols("ϕ p", real = true) ;

julia> p̄, q, q̄, c = symbols("p̄ q q̄ c", positive = true) ;

julia> coords = (ϕ, p, q) ; tmp_coords = (p̄, q̄) ; params = (c,) ;

julia> OM = Tens(c * [p̄ * q̄ * cos(ϕ), p̄ * q̄ * sin(ϕ), p * q]) ;

julia> Spheroidal = CoorSystemSym(OM, coords, tmp_coords, params; tmp_var = Dict(1-p^2 => p̄^2, q^2-1 => q̄^2), to_coords = Dict(p̄ => √(1-p^2), q̄ => √(q^2-1))) ;

Tens{order,dim,T,A<:AbstractArray}

Tensor type of any order defined by

• a multidata of components (of any type heriting from AbstractArray, e.g. Tensor or SymmetricTensor)
• a basis of AbstractBasis type
• a tuple of variances (covariant :cov or contravariant :cont) of length equal to the order of the tensor

Examples

julia> ℬ = Basis(Sym[1 0 0; 0 1 0; 0 1 1]) ;

julia> T = Tens(metric(ℬ,:cov),ℬ,(:cov,:cov))
Tens{2, 3, Sym, SymmetricTensor{2, 3, Sym, 6}}
# data: 3×3 SymmetricTensor{2, 3, Sym, 6}:
 1  0  0
 0  2  1
 0  1  1
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 1  0  0
 0  1  0
 0  1  1
# var: (:cov, :cov)

julia> components(T,(:cont,:cov),b)
3×3 Matrix{Sym}:
 1  0  0
 0  1  0
 0  0  1

TensOrtho{T} <: AbstractTens{4,3,T}

Orthotropic 4th-order tensor with material frame (e₁,e₂,e₃) and 9 independent elastic constants (C₁₁,C₂₂,C₃₃,C₁₂,C₁₃,C₂₃,C₄₄,C₅₅,C₆₆) where C₄₄=C₂₃₂₃, C₅₅=C₁₃₁₃, C₆₆=C₁₂₁₂:

ℂ = C₁₁P₁⊗P₁ + C₂₂P₂⊗P₂ + C₃₃P₃⊗P₃
  + C₁₂(P₁⊗P₂+P₂⊗P₁) + C₁₃(P₁⊗P₃+P₃⊗P₁) + C₂₃(P₂⊗P₃+P₃⊗P₂)
  + 2C₄₄(P₂⊠ˢP₃) + 2C₅₅(P₁⊠ˢP₃) + 2C₆₆(P₁⊠ˢP₂)

with Pₘ = eₘ⊗eₘ. The Kelvin-Mandel matrix in the material frame is block-diagonal:

[[C₁₁,C₁₂,C₁₃, 0,   0,   0  ],
 [C₁₂,C₂₂,C₂₃, 0,   0,   0  ],
 [C₁₃,C₂₃,C₃₃, 0,   0,   0  ],
 [ 0,  0,  0,  2C₄₄, 0,   0  ],
 [ 0,  0,  0,   0,  2C₅₅, 0  ],
 [ 0,  0,  0,   0,   0,  2C₆₆]]

TensOrtho(KMmat::AbstractMatrix, frame)

Build a TensOrtho from a 6×6 Kelvin-Mandel matrix expressed in the material frame. The matrix must have the block-diagonal orthotropic structure: upper-left 3×3 for normal stresses and lower-right 3×3 diagonal for shear.

TensOrtho(C11,C22,C33,C12,C13,C23,C44,C55,C66, frame)

Orthotropic tensor from the 9 elastic constants in the material frame frame.

TensTI{order,T,N} <: AbstractTens{order,3,T}

Transversely isotropic tensor of order order (always dim=3) with symmetry axis n, parametrised like TensISO{order,dim,T,N}.

Order 2 (N=2): data = (a, b) → 𝐀 = a·nT + b·nₙ where nₙ = n⊗n, nT = 𝟏 − nₙ.

• a: transverse coefficient (plane ⊥ n)
• b: axial coefficient (along n)
• When a = b: isotropic, equivalent to TensISO{2,3,T}(a)

Order 4 is handled separately by TensWalpole{T,N} (Walpole basis with 5 or 6 coefficients).

Constructor

TensTI{2}(a, b, n) → TensTI{2,T,2}

Construct a TI 2nd-order tensor a·nT + b·nₙ with symmetry axis n. The axis n can be a Vector, NTuple{3}, or any AbstractTens of order 1.

Examples

julia> n = [0., 0., 1.];

julia> A = TensTI{2}(5.0, 8.0, n);

julia> getarray(A)
3×3 Matrix{Float64}:
 5.0  0.0  0.0
 0.0  5.0  0.0
 0.0  0.0  8.0

julia> tr(A)
18.0

julia> isISO(A)
false

julia> inv(A).data
(0.2, 0.125)

julia> B = TensTI{2}(5.0, 5.0, n); isISO(B)
true

TensWalpole{T,N} <: AbstractTens{4,3,T}

Transversely isotropic 4th-order tensor stored in the Walpole basis {W₁,…,W₆} with symmetry axis n (assumed unit vector):

L = ℓ₁W₁ + ℓ₂W₂ + ℓ₃W₃ + ℓ₄W₄ + ℓ₅W₅ + ℓ₆W₆

where (nₙ = n⊗n, nT = 𝟏 − nₙ):

N=5 (major-symmetric, ℓ₃=ℓ₄): data=(ℓ₁,ℓ₂,ℓ₃,ℓ₅,ℓ₆). N=6 (general): data=(ℓ₁,ℓ₂,ℓ₃,ℓ₄,ℓ₅,ℓ₆).

Synthetic notation: L ≡ ([[ℓ₁,ℓ₃],[ℓ₄,ℓ₂]], ℓ₅, ℓ₆).

• Double contraction: (L⊡M)_mat = L_mat × M_mat, (L⊡M)₅ = ℓ₅m₅, (L⊡M)₆ = ℓ₆m₆
• Inverse: (L⁻¹)_mat = (L_mat)⁻¹, 1/ℓ₅, 1/ℓ₆

TensWalpole(ℓ₁,ℓ₂,ℓ₃,ℓ₅,ℓ₆, n) → TensWalpole{T,5}

Major-symmetric Walpole tensor (ℓ₃ = ℓ₄), 5 independent scalars, with axis n.

TensWalpole(ℓ₁,ℓ₂,ℓ₃,ℓ₄,ℓ₅,ℓ₆, n) → TensWalpole{T,6}

General (not necessarily major-symmetric) Walpole tensor with axis n.

inv(t::TensOrtho) → TensOrtho

Inverse via the KM matrix in the material frame (block-diagonal, efficiently invertible).

inv(t::TensTI{2,T,2}) → TensTI{2,T,2}

Inverse: (a·nT + b·nₙ)⁻¹ = (1/a)·nT + (1/b)·nₙ.

Examples

julia> A = TensTI{2}(5.0, 8.0, [0., 0., 1.]);

julia> inv(A).data
(0.2, 0.125)

inv(t::TensWalpole{T,5}) → TensWalpole{T,5}
inv(t::TensWalpole{T,6}) → TensWalpole{T,6}

Inverse via the 2×2 Walpole matrix and scalar inverses for ℓ₅, ℓ₆.

dot(t1::AbstractTens{order1,dim}, t2::AbstractTens{order2,dim})

Define a contracted product between two tensors

a ⋅ b = aⁱbⱼ

normalize(ℬ::AbstractBasis, var = cov)

Build a basis after normalization of column vectors of input matrix v where columns define either

• primal vectors ie eᵢ=v[:,i]/norm(v[:,i]) if var = :cov as by default
• dual vector ie eⁱ=v[:,i]/norm(v[:,i]) if var = :cont.

tr(t::TensTI{2}) → scalar

Trace: tr(a·nT + b·nₙ) = 2a + b.

DIV(t::AbstractTens{order,dim,Sym},CS::CoorSystemSym{dim}) where {order,dim,T<:Number}

Calculate the divergence of T with respect to the coordinate system CS

DIV(f::Function, CS::CoorSystemNum) -> Function

Return a function x₀ -> divergence as a scalar (if f is a vector field) or AbstractTens (if f is a tensor field of order ≥ 2).

f may return a plain Array or an AbstractTens (components are extracted via Array).

DIV(T::AbstractTens{order,dim,Sym},SM::SubManifoldSym{dim}) where {order,dim}

Calculate the divergence of T with respect to the coordinate system SM

GRAD(t::Union{t,AbstractTens{order,dim,T}},CS::CoorSystemSym{dim}) where {order,dim,T<:Number}

Calculate the gradient of t with respect to the coordinate system CS

GRAD(f::Function, CS::CoorSystemNum) -> Function

Return a function x₀ -> gradient as an AbstractTens in the normalized basis at x₀.

• If f returns a scalar, the gradient is a rank-1 AbstractTens.
• If f returns a rank-n AbstractTens or array, the gradient is a rank-n+1 AbstractTens.

GRAD(T::Union{Sym,AbstractTens{order,dim,Sym}},SM::SubManifoldSym{dim}) where {order,dim}

Calculate the gradient of T with respect to the coordinate system SM

HESS(t::Union{Sym,AbstractTens{order,dim,Sym}},CS::CoorSystemSym{dim}) where {order,dim,T<:Number}

Calculate the Hessian of T with respect to the coordinate system CS

HESS(f::Function, CS::CoorSystemNum) -> Function

Return a function x₀ -> Hessian as an AbstractTens{2} in the normalized basis. Uses the raw (plain-array) double gradient internally.

HESS(T::Union{Sym,AbstractTens{order,dim,Sym}},SM::SubManifoldSym{dim}) where {order,dim}

Calculate the Hessian of T with respect to the coordinate system SM

ISO(::Val{dim} = Val(3), ::Val{T} = Val(Sym))

Return the three fourth-order isotropic tensors 𝕀, 𝕁, 𝕂

Examples

julia> 𝕀, 𝕁, 𝕂 = ISO() ;

KM(t::TensOrtho)

Returns the 6×6 Kelvin-Mandel matrix in the canonical frame. Use KM_material(t) for the block-diagonal form in the material frame.

KM(t::TensWalpole)

Kelvin-Mandel (6×6) matrix of the Walpole tensor.

KM(t::AbstractTens{order,dim}; kwargs...)
KM(t::AbstractTens{order,dim}, var::NTuple{order,Symbol}, b::AbstractBasis{dim}; kwargs...)

Write the components of a second or fourth order tensor in Kelvin-Mandel notation

Examples

julia> σ = Tens(SymmetricTensor{2,3}((i, j) -> symbols("σ$i$j", real = true))) ;

julia> KM(σ)
6-element Vector{Sym}:
         σ11
         σ22
         σ33
      √2⋅σ32
      √2⋅σ31
      √2⋅σ21

julia> C = Tens(SymmetricTensor{4,3}((i, j, k, l) -> symbols("C$i$j$k$l", real = true))) ;

julia> KM(C)
6×6 Matrix{Sym}:
         C₁₁₁₁     C₁₁₂₂     C₁₁₃₃  √2⋅C₁₁₃₂  √2⋅C₁₁₃₁  √2⋅C₁₁₂₁
         C₂₂₁₁     C₂₂₂₂     C₂₂₃₃  √2⋅C₂₂₃₂  √2⋅C₂₂₃₁  √2⋅C₂₂₂₁
         C₃₃₁₁     C₃₃₂₂     C₃₃₃₃  √2⋅C₃₃₃₂  √2⋅C₃₃₃₁  √2⋅C₃₃₂₁
      √2⋅C₃₂₁₁  √2⋅C₃₂₂₂  √2⋅C₃₂₃₃   2⋅C₃₂₃₂   2⋅C₃₂₃₁   2⋅C₃₂₂₁
      √2⋅C₃₁₁₁  √2⋅C₃₁₂₂  √2⋅C₃₁₃₃   2⋅C₃₁₃₂   2⋅C₃₁₃₁   2⋅C₃₁₂₁
      √2⋅C₂₁₁₁  √2⋅C₂₁₂₂  √2⋅C₂₁₃₃   2⋅C₂₁₃₂   2⋅C₂₁₃₁   2⋅C₂₁₂₁

KM(v::AllIsotropic{dim}; kwargs...)

Kelvin-Mandel vector or matrix representation

KM_material(t::TensOrtho)

Returns the 6×6 Kelvin-Mandel matrix in the material frame (block-diagonal).

LAPLACE(t::Union{Sym,AbstractTens{order,dim,Sym}},CS::CoorSystemSym{dim}) where {order,dim,T<:Number}

Calculate the Laplace operator of T with respect to the coordinate system CS

LAPLACE(f::Function, CS::CoorSystemNum) -> Function

Return a function x₀ -> Laplacian (scalar). Uses the raw (plain-array) gradient and divergence internally, so it is fully compatible with ForwardDiff differentiation.

LAPLACE(T::Union{Sym,AbstractTens{order,dim,Sym}},SM::SubManifoldSym{dim}) where {order,dim}

Calculate the Laplace operator of T with respect to the coordinate system SM

LeviCivita(T::Type{<:Number} = Sym)

Builds an Array{T,3} of Levi-Civita Symbol ϵᵢⱼₖ = (i-j) (j-k) (k-i) / 2

Examples

julia> ε = LeviCivita(Sym)
3×3×3 Array{Sym, 3}:
[:, :, 1] =
 0   0  0
 0   0  1
 0  -1  0

[:, :, 2] =
 0  0  -1
 0  0   0
 1  0   0

[:, :, 3] =
  0  1  0
 -1  0  0
  0  0  0

SYMGRAD(t::Union{T,AbstractTens{order,dim,T}},CS::CoorSystemSym{dim}) where {order,dim,T<:Number}

Calculate the symmetrized gradient of T with respect to the coordinate system CS

SYMGRAD(f::Function, CS::CoorSystemNum) -> Function

Return a function x₀ -> symmetric gradient as an AbstractTens{2} in the normalized basis at x₀. Applies to vector-valued functions f.

SYMGRAD(T::Union{Sym,AbstractTens{order,dim,Sym}},SM::SubManifoldSym{dim}) where {order,dim}

Calculate the symmetrized gradient of T with respect to the coordinate system SM

Walpole(n)           → (W₁,W₂,W₃,W₄,W₅,W₆)
Walpole(n; sym=true) → (W₁ˢ,W₂ˢ,W₃ˢ,W₄ˢ,W₅ˢ) where W₃ˢ = W₃+W₄

_KM_of_array(A::AbstractArray{T,4}) → Matrix{T}

Compute the 6×6 Kelvin-Mandel matrix of a 3×3×3×3 array.

_KM_of_array(A::AbstractArray{T,2}) → Matrix{T}

Compute the 6-vector (or just return the 3×3 matrix) Kelvin-Mandel for a 2nd-order tensor.

_KM_rotation(θ, ϕ, ψ) → SMatrix{6,6}

6×6 Kelvin-Mandel (Bond) rotation matrix from Euler angles (θ, ϕ, ψ). Compatible with ForwardDiff.

Transforms a 6×6 KM matrix as: C' = Q' · C · Q where Q = _KM_rotation(...).

Uses the standard Bond matrix construction from R = _rot3_raw(θ, ϕ, ψ).

Examples

julia> Q = _KM_rotation(0.3, 0.5, 0.0);

julia> size(Q)
(6, 6)

_angles_from_n(n) → (θ, ϕ)

Extract spherical angles (θ, ϕ) from a unit vector n:

• θ = atan(√(n₁² + n₂²), n₃) (polar angle from z-axis)
• ϕ = atan(n₂, n₁) (azimuthal angle)

Examples

julia> _angles_from_n((0.0, 0.0, 1.0))
(0.0, 0.0)

julia> _angles_from_n((1.0, 0.0, 0.0))
(1.5707963267948966, 0.0)

Apply the Christoffel correction to a rank-order tensor in natural contravariant components. t_arr : array of natural contravariant components (size dim^order) Γᵢ : Γ[i,:,:] slice (dim×dim matrix, entry [k,j] = Γᵏᵢₖ ... = Γᵏᵢⱼ)

_build_ORTHO_KM(C₁₁, C₂₂, C₃₃, C₁₂, C₁₃, C₂₃, C₄₄, C₅₅, C₆₆) → SMatrix{6,6}

Build a 6×6 KM matrix in the material frame from 9 orthotropic parameters.

Examples

julia> B = _build_ORTHO_KM(10.0, 8.0, 12.0, 3.0, 2.5, 1.5, 2.0, 3.0, 3.5);

julia> B[4,4]  # = 2*C₄₄ = 4
4.0

_build_TI_KM(ℓ₁, ℓ₂, ℓ₃, ℓ₅, ℓ₆) → SMatrix{6,6}

Build a 6×6 KM matrix (in the frame where e₃ is the TI axis) from 5 major-symmetric Walpole coefficients.

Examples

julia> B = _build_TI_KM(12.0, 13.0, sqrt(2)*2.5, 7.0, 4.0);

julia> B[1,1]  # = (ℓ₂ + ℓ₅)/2 = (13+7)/2 = 10
10.0

Compute Christoffel symbols from position vector function OM via nested AD. OM_func(x) returns a vector of Cartesian coordinates.

Compute Christoffel symbols for an orthogonal system from the Jacobian of the Lamé-coefficient function.

Convention: Γ[i,j,k] = Γᵏᵢⱼ.

_n_from_angles(θ, ϕ) → NTuple{3}

Build a unit vector from spherical angles: n = (sinθ·cosϕ, sinθ·sinϕ, cosθ).

Examples

julia> _n_from_angles(0.0, 0.0)
(0.0, 0.0, 1.0)

_project_ORTHO_KM(C::AbstractMatrix) → NTuple{9}

Extract 9 orthotropic parameters from a 6×6 KM matrix in the material frame. Returns (C₁₁, C₂₂, C₃₃, C₁₂, C₁₃, C₂₃, C₄₄, C₅₅, C₆₆).

For non-symmetric KM matrices (non-major-symmetric tensors), the off-diagonal entries are averaged: Cᵢⱼ = (C[i,j] + C[j,i]) / 2.

Examples

julia> C = Float64[10 3 2.5 0 0 0; 3 8 1.5 0 0 0; 2.5 1.5 12 0 0 0;
                   0 0 0 4 0 0; 0 0 0 0 6 0; 0 0 0 0 0 7];

julia> _project_ORTHO_KM(C)
(10.0, 8.0, 12.0, 3.0, 2.5, 1.5, 2.0, 3.0, 3.5)

_project_TI_KM(C::AbstractMatrix) → (ℓ₁, ℓ₂, ℓ₃, ℓ₅, ℓ₆)

Project a 6×6 KM matrix (assumed to be in a frame where e₃ is the TI symmetry axis) onto the major-symmetric TI subspace. Returns the 5 Walpole coefficients.

Formulas (from ECHOES tensor_ti.h):

• c = (C[1,1] + C[2,2]) / 2
• d = (C[1,2] + C[2,1]) / 2
• ℓ₁ = C[3,3]
• ℓ₂ = c + d
• ℓ₃ = (C[1,3] + C[2,3] + C[3,1] + C[3,2]) / (2√2)
• ℓ₅ = (c − d + C[6,6]) / 2
• ℓ₆ = (C[4,4] + C[5,5]) / 2

Examples

julia> C = Float64[10 3 2.5 0 0 0; 3 10 2.5 0 0 0; 2.5 2.5 12 0 0 0;
                   0 0 0 4 0 0; 0 0 0 0 4 0; 0 0 0 0 0 7];

julia> _project_TI_KM(C)
(12.0, 13.0, 3.5355339059327378, 7.0, 4.0)

_rot3_raw(θ, ϕ, ψ) → SMatrix{3,3}

Rotation matrix (RotZYZ convention: R = Rz(ϕ)·Ry(θ)·Rz(ψ)) built from explicit cos/sin, compatible with ForwardDiff Dual numbers.

The third column of the matrix is the direction (sinθ cosϕ, sinθ sinϕ, cosθ).

Examples

julia> R = _rot3_raw(0.3, 0.5, 0.0)
3×3 SMatrix{3, 3, Float64, 9} with indices SOneTo(3)×SOneTo(3):
 ...

Return an array of shape (dim, dim, ...) (order times) where entry [i₁,…,iₙ] = χᵢ₁ ⋯ χᵢₙ.

angles(M::AbstractMatrix{T})

Determine the Euler angles corresponding to the input matrix supposed to be a rotation matrix or at least a similarity

Examples

julia> θ, ϕ, ψ = symbols("θ, ϕ, ψ", real = true) ; ℬʳ = RotatedBasis(θ, ϕ, ψ) ; display(vecbasis(ℬʳ, :cov))
3×3 Tensor{2, 3, Sym, 9}:
 -sin(ψ)⋅sin(ϕ) + cos(θ)⋅cos(ψ)⋅cos(ϕ)  -sin(ψ)⋅cos(θ)⋅cos(ϕ) - sin(ϕ)⋅cos(ψ)  sin(θ)⋅cos(ϕ)
  sin(ψ)⋅cos(ϕ) + sin(ϕ)⋅cos(θ)⋅cos(ψ)  -sin(ψ)⋅sin(ϕ)⋅cos(θ) + cos(ψ)⋅cos(ϕ)  sin(θ)⋅sin(ϕ)
                        -sin(θ)⋅cos(ψ)                          sin(θ)⋅sin(ψ)         cos(θ)

julia> angles(ℬʳ)
(θ = θ, ϕ = ϕ, ψ = ψ)

argTI(t::TensWalpole) → (C₁₁₁₁, C₁₁₂₂, C₁₁₃₃, C₃₃₃₃, C₂₃₂₃)

Extract the 5 independent TI components from a Walpole tensor, directly from the stored coefficients (no array materialisation).

Inverse of tensTI:

• C₃₃₃₃ = ℓ₁
• C₁₁₁₁ = (ℓ₂ + ℓ₅)/2
• C₁₁₂₂ = (ℓ₂ − ℓ₅)/2
• C₁₁₃₃ = ℓ₃/√2
• C₂₃₂₃ = ℓ₆/2

See also argTI_eng.

argTI_Hoenig(𝕊::TensWalpole) → (E, ν₁, ν₂, H, Γ)

Extract the 5 Hoenig parameters from a TI compliance tensor.

See also tensTI_Hoenig, argTI_eng.

argTI_eng(𝕊::TensWalpole) → (E₁, E₃, ν₁₂, ν₃₁, G₃₁)

Extract engineering constants from a TI compliance tensor.

See also tensTI_eng, argTI.

best_sym_tens(t, n_or_frame; proj=(:ISO, :TI, :ORTHO), ε=1e-6)

Find the best symmetry of tensor t with a fixed symmetry axis n (for TI) or material frame frame (for ORTHO). No rotation optimization is performed.

• n_or_frame: a vector (axis for TI) or OrthonormalBasis{3} (frame for ORTHO). For ISO projection the extra argument is ignored.

Returns (projected, d, drel, sym).

Examples

julia> n = [0., 0., 1.];

julia> C = tensTI(10., 3., 2.5, 12., 2., n);

julia> _, _, drel, sym = best_sym_tens(C, n);

julia> sym == :TI && drel < 1e-12
true

best_sym_tens(t; proj=(:ISO, :TI, :ORTHO), ε=1e-6)

Find the best (most restrictive) symmetry of tensor t by trying each symmetry in proj (from most symmetric to least) and accepting the first whose relative projection error is below ε.

Requires NLopt for rotation-optimized TI and ORTHO projections (i.e. when no axis/frame is provided). For fixed-basis usage, see the variant with n_or_frame.

Returns (projected, d, drel, sym) where sym ∈ {:ISO, :TI, :ORTHO, :ANISO}.

Examples

julia> C = tensTI(10., 3., 2.5, 12., 2., [0., 0., 1.]);

julia> _, _, _, sym = best_sym_tens(C);

julia> sym
:ISO  # or :TI depending on the tensor

best_sym_tens(t, n_or_frame; proj=(:ISO, :TI, :ORTHO), ε=1e-6)

Find the best symmetry of tensor t with a fixed symmetry axis n (for TI) or material frame frame (for ORTHO). No rotation optimization is performed.

• n_or_frame: a vector (axis for TI) or OrthonormalBasis{3} (frame for ORTHO). For ISO projection the extra argument is ignored.

Returns (projected, d, drel, sym).

Examples

julia> n = [0., 0., 1.];

julia> C = tensTI(10., 3., 2.5, 12., 2., n);

julia> _, _, drel, sym = best_sym_tens(C, n);

julia> sym == :TI && drel < 1e-12
true

best_sym_tens(t; proj=(:ISO, :TI, :ORTHO), ε=1e-6)

Find the best (most restrictive) symmetry of tensor t by trying each symmetry in proj (from most symmetric to least) and accepting the first whose relative projection error is below ε.

Requires NLopt for rotation-optimized TI and ORTHO projections (i.e. when no axis/frame is provided). For fixed-basis usage, see the variant with n_or_frame.

Returns (projected, d, drel, sym) where sym ∈ {:ISO, :TI, :ORTHO, :ANISO}.

Examples

julia> C = tensTI(10., 3., 2.5, 12., 2., [0., 0., 1.]);

julia> _, _, _, sym = best_sym_tens(C);

julia> sym
:ISO  # or :TI depending on the tensor

change_tens(t::AbstractTens{order,dim,T},ℬ::AbstractBasis{dim},var::NTuple{order,Symbol})
change_tens(t::AbstractTens{order,dim,T},ℬ::AbstractBasis{dim})
change_tens(t::AbstractTens{order,dim,T},var::NTuple{order,Symbol})

Rewrite the same tensor with components corresponding to new variances and/or to a new basis

julia> ℬ = Basis(Sym[0 1 1; 1 0 1; 1 1 0]) ;

julia> TV = Tens(Tensor{1,3}(i->symbols("v$i",real=true)))
TensND.TensCanonical{1, 3, Sym, Vec{3, Sym}}
# data: 3-element Vec{3, Sym}:
 v₁
 v₂
 v₃
# basis: 3×3 TensND.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1
# var: (:cont,)

julia> factor.(components(TV, ℬ, (:cont,)))
3-element Vector{Sym}:
 -(v1 - v2 - v3)/2
  (v1 - v2 + v3)/2
  (v1 + v2 - v3)/2

julia> ℬ₀ = Basis(Sym[0 1 1; 1 0 1; 1 1 1]) ;

julia> TV0 = change_tens(TV, ℬ₀)
Tens{1, 3, Sym, Vec{3, Sym}}
# data: 3-element Vec{3, Sym}:
     -v₁ + v₃
     -v₂ + v₃
 v₁ + v₂ - v₃
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 0  1  1
 1  0  1
 1  1  1
# var: (:cont,)

change_tens_canon(t::AbstractTens{order,dim,T},var::NTuple{order,Symbol})

Rewrite the same tensor with components corresponding to the canonical basis

julia> ℬ = Basis(Sym[0 1 1; 1 0 1; 1 1 0]) ;

julia> TV = Tens(Tensor{1,3}(i->symbols("v$i",real=true)), ℬ)
Tens{1, 3, Sym, Vec{3, Sym}}
# data: 3-element Vec{3, Sym}:
 v₁
 v₂
 v₃
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 0  1  1
 1  0  1
 1  1  1
# var: (:cont,)

julia> TV0 = change_tens_canon(TV)
TensND.TensCanonical{1, 3, Sym, Vec{3, Sym}}
# data: 3-element Vec{3, Sym}:
      v₂ + v₃
      v₁ + v₃
 v₁ + v₂ + v₃
# basis: 3×3 TensND.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1
# var: (:cont,)

components(t::AbstractTens{order,dim,T},ℬ::AbstractBasis{dim},var::NTuple{order,Symbol})
components(t::AbstractTens{order,dim,T},ℬ::AbstractBasis{dim})
components(t::AbstractTens{order,dim,T},var::NTuple{order,Symbol})

Extract the components of a tensor for new variances and/or in a new basis

Examples

julia> ℬ = Basis(Sym[0 1 1; 1 0 1; 1 1 0]) ;

julia> TV = Tens(Tensor{1,3}(i->symbols("v$i",real=true)))
TensND.TensCanonical{1, 3, Sym, Vec{3, Sym}}
# data: 3-element Vec{3, Sym}:
 v₁
 v₂
 v₃
# basis: 3×3 TensND.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1
# var: (:cont,)

julia> factor.(components(TV, ℬ, (:cont,)))
3-element Vector{Sym}:
 -(v1 - v2 - v3)/2
  (v1 - v2 + v3)/2
  (v1 + v2 - v3)/2

julia> components(TV, ℬ, (:cov,))
3-element Vector{Sym}:
 v₂ + v₃
 v₁ + v₃
 v₁ + v₂

julia> simplify.(components(TV, normalize(ℬ), (:cov,)))
3-element Vector{Sym}:
 sqrt(2)*(v2 + v3)/2
 sqrt(2)*(v1 + v3)/2
 sqrt(2)*(v1 + v2)/2

julia> TT = Tens(Tensor{2,3}((i,j)->symbols("t$i$j",real=true)))
TensND.TensCanonical{2, 3, Sym, Tensor{2, 3, Sym, 9}}
# data: 3×3 Tensor{2, 3, Sym, 9}:
 t₁₁  t₁₂  t₁₃
 t₂₁  t₂₂  t₂₃
 t₃₁  t₃₂  t₃₃
# basis: 3×3 TensND.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1
# var: (:cont, :cont)

julia> components(TT, ℬ, (:cov,:cov))
3×3 Matrix{Sym}:
 t₂₂ + t₂₃ + t₃₂ + t₃₃  t₂₁ + t₂₃ + t₃₁ + t₃₃  t₂₁ + t₂₂ + t₃₁ + t₃₂
 t₁₂ + t₁₃ + t₃₂ + t₃₃  t₁₁ + t₁₃ + t₃₁ + t₃₃  t₁₁ + t₁₂ + t₃₁ + t₃₂
 t₁₂ + t₁₃ + t₂₂ + t₂₃  t₁₁ + t₁₃ + t₂₁ + t₂₃  t₁₁ + t₁₂ + t₂₁ + t₂₂

julia> factor.(components(TT, ℬ, (:cont,:cov)))
3×3 Matrix{Sym}:
 -(t12 + t13 - t22 - t23 - t32 - t33)/2  …  -(t11 + t12 - t21 - t22 - t31 - t32)/2
  (t12 + t13 - t22 - t23 + t32 + t33)/2      (t11 + t12 - t21 - t22 + t31 + t32)/2
  (t12 + t13 + t22 + t23 - t32 - t33)/2      (t11 + t12 + t21 + t22 - t31 - t32)/2

components_canon(t::AbstractTens)

Extract the components of a tensor in the canonical basis

contract(t::AbstractTens{order,dim}, i::Integer, j::Integer)

Calculate the tensor obtained after contraction with respect to the indices i and j

coorsys_cartesian(coords = symbols("x y z", real = true))

Return the cartesian coordinate system

Examples

julia> Cartesian = coorsys_cartesian() ; 𝐗 = getcoords(Cartesian) ; 𝐄 = unitvec(Cartesian) ; ℬ = getbasis(Cartesian)

julia> 𝛔 = Tens(SymmetricTensor{2,3}((i, j) -> SymFunction("σ$i$j", real = true)(𝐗...))) ;

julia> DIV(𝛔, CScar)
Tens.TensCanonical{1, 3, Sym, Vec{3, Sym}}
# data: 3-element Vec{3, Sym}:
 Derivative(σ11(x, y, z), x) + Derivative(σ21(x, y, z), y) + Derivative(σ31(x, y, z), z)
 Derivative(σ21(x, y, z), x) + Derivative(σ22(x, y, z), y) + Derivative(σ32(x, y, z), z)
 Derivative(σ31(x, y, z), x) + Derivative(σ32(x, y, z), y) + Derivative(σ33(x, y, z), z)
# basis: 3×3 Tens.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1
# var: (:cont,)

coorsys_cartesian_num(T=Float64)

Numerical Cartesian coordinate system in 3D. Coordinates: (x, y, z).

coorsys_cylindrical(coords = (symbols("r", positive = true), symbols("θ", real = true), symbols("z", real = true)); canonical = false)

Return the cylindrical coordinate system

Examples

julia> Cylindrical = coorsys_cylindrical() ; rθz = getcoords(Cylindrical) ; 𝐞ʳ, 𝐞ᶿ, 𝐞ᶻ = unitvec(Cylindrical) ; ℬᶜ = getbasis(Cylindrical)

julia> 𝐯 = Tens(Vec{3}(i -> SymFunction("v$(rθz[i])", real = true)(rθz...)), ℬᶜ) ;

julia> DIV(𝐯, Cylindrical)
                                                  ∂
                                    vr(r, θ, z) + ──(vθ(r, θ, z))
∂                 ∂                               ∂θ
──(vr(r, θ, z)) + ──(vz(r, θ, z)) + ─────────────────────────────
∂r                ∂z                              r

coorsys_cylindrical_num(T=Float64)

Numerical cylindrical coordinate system. Coordinates: (r, θ, z). Lamé coefficients: (1, r, 1).

coorsys_polar(coords = (symbols("r", positive = true), symbols("θ", real = true)); canonical = false)

Return the polar coordinate system

Examples

julia> Polar = coorsys_polar() ; r, θ = getcoords(Polar) ; 𝐞ʳ, 𝐞ᶿ = unitvec(Polar) ; ℬᵖ = getbasis(Polar)

julia> f = SymFunction("f", real = true)(r, θ) ;

julia> LAPLACE(f, Polar)
                               2
                              ∂
                             ───(f(r, θ))
                               2
               ∂             ∂θ
  2            ──(f(r, θ)) + ────────────
 ∂             ∂r                 r
───(f(r, θ)) + ──────────────────────────
  2                        r
∂r

coorsys_polar_num(T=Float64)

Numerical polar coordinate system. Coordinates: (r, θ). Lamé coefficients: (1, r).

coorsys_spherical(coords = (symbols("θ", real = true), symbols("ϕ", real = true), symbols("r", positive = true)); canonical = false)

Return the spherical coordinate system

Examples

julia> Spherical = coorsys_spherical() ; θ, ϕ, r = getcoords(Spherical) ; 𝐞ᶿ, 𝐞ᵠ, 𝐞ʳ = unitvec(Spherical) ; ℬˢ = getbasis(Spherical)

julia> for σⁱʲ ∈ ("σʳʳ", "σᶿᶿ", "σᵠᵠ") @eval $(Symbol(σⁱʲ)) = SymFunction($σⁱʲ, real = true)($r) end ;

julia> 𝛔 = σʳʳ * 𝐞ʳ ⊗ 𝐞ʳ + σᶿᶿ * 𝐞ᶿ ⊗ 𝐞ᶿ + σᵠᵠ * 𝐞ᵠ ⊗ 𝐞ᵠ ;

julia> div𝛔 = DIV(𝛔, Spherical)
Tens.TensRotated{1, 3, Sym, Vec{3, Sym}}
# data: 3-element Vec{3, Sym}:
                              (-σᵠᵠ(r) + σᶿᶿ(r))*cos(θ)/(r*sin(θ))
                                                                 0
 Derivative(σʳʳ(r), r) + (σʳʳ(r) - σᵠᵠ(r))/r + (σʳʳ(r) - σᶿᶿ(r))/r
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 cos(θ)⋅cos(ϕ)  -sin(ϕ)  sin(θ)⋅cos(ϕ)
 sin(ϕ)⋅cos(θ)   cos(ϕ)  sin(θ)⋅sin(ϕ)
       -sin(θ)        0         cos(θ)
# var: (:cont,)

julia> div𝛔 ⋅ 𝐞ʳ
d            σʳʳ(r) - σᵠᵠ(r)   σʳʳ(r) - σᶿᶿ(r)
──(σʳʳ(r)) + ─────────────── + ───────────────
dr                  r                 r

coorsys_spherical_num(T=Float64)

Numerical spherical coordinate system. Coordinates: (θ, ϕ, r) — same convention as coorsys_spherical(). Lamé coefficients: (r, r·sin(θ), 1).

coorsys_spheroidal(coords = (symbols("ϕ", real = true),symbols("p", real = true),symbols("q", positive = true),),
                        c = symbols("c", positive = true),tmp_coords = (symbols("p̄ q̄", positive = true)...,),)

Return the spheroidal coordinate system

Examples

julia> Spheroidal = coorsys_spheroidal() ; OM = getOM(Spheroidal)
Tens.TensCanonical{1, 3, Sym, Vec{3, Sym}}
# data: 3-element Vec{3, Sym}:
 c⋅p̄⋅q̄⋅cos(ϕ)
 c⋅p̄⋅q̄⋅sin(ϕ)
          c⋅p⋅q
# basis: 3×3 Tens.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1
# var: (:cont,)

julia> LAPLACE(OM[1]^2, Spheroidal)
2

fromISO(A::TensISO{4,3}, n) → TensWalpole{T,5}

Convert an isotropic 4th-order tensor αJ + βK into its Walpole representation.

Formulas: ℓ₁=(α+2β)/3, ℓ₂=(2α+β)/3 (note: dim=3 → these are (3k,2μ) related), ℓ₃=ℓ₄=√2(α−β)/3, ℓ₅=ℓ₆=β. Here α = data[1] and β = data[2] in TensISO (coefficients of J and K).

get_ℓ(t::TensWalpole) → NTuple{6}

Always returns a 6-tuple (ℓ₁,ℓ₂,ℓ₃,ℓ₄,ℓ₅,ℓ₆). For N=5 (symmetric), ℓ₃ = ℓ₄ is stored once so data[3] is duplicated.

getarray(t::TensOrtho{T}) → Array{T,4}

Compute the 3×3×3×3 component array in the canonical frame.

getarray(t::TensTI{2,T,2}) → Array{T,2}

Compute the 3×3 component array: a*(δᵢⱼ − nᵢnⱼ) + b*nᵢnⱼ.

Examples

julia> A = TensTI{2}(5.0, 8.0, [0., 0., 1.]);

julia> getarray(A)
3×3 Matrix{Float64}:
 5.0  0.0  0.0
 0.0  5.0  0.0
 0.0  0.0  8.0

getarray(t::TensWalpole{T}) → Array{T,4}

Compute the 3×3×3×3 component array from the Walpole coefficients and axis.

getaxis(t::TensWalpole) → NTuple{3}

Returns the symmetry axis as a 3-tuple.

init_cartesian(coords = symbols("x y z", real = true))

Returns the coordinates, unit vectors and basis of the cartesian basis

Examples

julia> coords, vectors, ℬ = init_cartesian() ; x, y, z = coords ; 𝐞₁, 𝐞₂, 𝐞₃ = vectors ;

init_cylindrical(coords = (symbols("r", positive = true), symbols("θ z", real = true)...); canonical = false)

Returns the coordinates, base vectors and basis of the cylindrical basis

Examples

julia> coords, vectors, ℬᶜ = init_cylindrical() ; r, θ, z = coords ; 𝐞ʳ, 𝐞ᶿ, 𝐞ᶻ = vectors ;

init_polar(coords = (symbols("r θ", real = true)); canonical = false)

Returns the coordinates, base vectors and basis of the polar basis

Examples

julia> coords, vectors, ℬᵖ = init_polar() ; r, θ = coords ; 𝐞ʳ, 𝐞ᶿ = vectors ;

init_rotated(coords = symbols("θ ϕ ψ", real = true); canonical = false)

Return the angles, base vectors and basis of the rotated basis. Note that here the coordinates are angles and do not represent a valid parametrization of ℝ³

Examples

julia> angles, vectors, ℬʳ = init_rotated() ; θ, ϕ, ψ = angles ; 𝐞ᶿ, 𝐞ᵠ, 𝐞ʳ = vectors ;

init_spherical(coords = (symbols("θ ϕ", real = true)..., symbols("r", positive = true)); canonical = false)

Return the coordinates, base vectors and basis of the spherical basis. Take care that the order of the 3 vectors is 𝐞ᶿ, 𝐞ᵠ, 𝐞ʳ so that the basis coincides with the canonical one when the angles are null and in consistency the coordinates are ordered as θ, ϕ, r.

Examples

julia> coords, vectors, ℬˢ = init_spherical() ; θ, ϕ, r = coords ; 𝐞ᶿ, 𝐞ᵠ, 𝐞ʳ  = vectors ;

invKM(v::AbstractVecOrMat; kwargs...)

Define a tensor from a Kelvin-Mandel vector or matrix representation

isTI(A)

Return true if A is a TensWalpole, indicating transverse isotropy.

isorthogonal(ℬ::AbstractBasis)

Check whether the basis ℬ is orthogonal

isorthonormal(ℬ::AbstractBasis)

Check whether the basis ℬ is orthonormal

metric(ℬ::AbstractBasis, var = :cov)

Return the covariant (if var = :cov) or contravariant (if var = :cont) metric matrix

natural_basis(CS::CoorSystemNum, x₀) → AbstractBasis

Return the natural (non-normalized) basis at point x₀ as an OrthogonalBasis (i.e. a RotatedBasis scaled by the Lamé coefficients).

natvec(CS::CoorSystemNum, x₀, i, var=:cov) → AbstractTens

Return the i-th natural basis vector at point x₀, either covariant (var=:cov, proportional to χᵢ) or contravariant (var=:cont, proportional to 1/χᵢ).

normalized_basis(CS::CoorSystemNum, x₀) → AbstractBasis

Return the normalized (orthonormal) basis at point x₀ as a RotatedBasis or CanonicalBasis (when the system is Cartesian at that point).

otimesul(t1::AbstractTens{order1,dim}, t2::AbstractTens{order2,dim})

Define a special tensor product between two tensors of at least second order

(𝐚 ⊠ˢ 𝐛) ⊡ 𝐩 = (𝐚 ⊠ 𝐛) ⊡ (𝐩 + ᵗ𝐩)/2 = 1/2(aⁱᵏbʲˡ+aⁱˡbʲᵏ) pₖₗ eᵢ⊗eⱼ

proj_tens(sym::Symbol, A::AbstractArray, arg) → (projected, d, drel)

Convenience dispatch: proj_tens(:TI, A, n) calls proj_tens(Val(:TI), A, n).

proj_tens(::Val{:ORTHO}, A::AbstractArray{T,4}, frame::OrthonormalBasis{3}) → (TensOrtho{T}, d, drel)

Project a 4th-order tensor A (3×3×3×3) onto the orthotropic subspace with fixed material frame frame. Returns a TensOrtho{T}.

Examples

julia> frame = CanonicalBasis{3,Float64}();

julia> t = TensOrtho(10., 8., 12., 3., 2.5, 1.5, 2., 3., 3.5, frame);

julia> B, d, drel = proj_tens(:ORTHO, getarray(t), frame);

julia> d < 1e-12
true

proj_tens(::Val{:ORTHO}, A::AbstractArray{T,4}) where {T<:AbstractFloat}

Find the best orthotropic approximation of A by optimizing over all possible material frames. Requires the NLopt package: using NLopt.

proj_tens(::Val{:ORTHO}, A::AbstractArray{T,2}, frame::OrthonormalBasis{3}) → (Array{T,2}, d, drel)

Project a 2nd-order tensor A (3×3) onto the orthotropic subspace with fixed material frame frame. The projection is diag(M₁₁, M₂₂, M₃₃) in the material frame.

Examples

julia> frame = CanonicalBasis{3,Float64}();

julia> A = [5. 1 2; 1 8 3; 2 3 12];

julia> B, d, drel = proj_tens(:ORTHO, A, frame);

julia> B ≈ diagm([5., 8., 12.])
true

proj_tens(::Val{:ORTHO}, A::AbstractArray{T,2}) where {T<:AbstractFloat}

Find the best orthotropic approximation of a 2nd-order tensor A by optimizing the material frame. Requires the NLopt package: using NLopt.

proj_tens(::Val{:TI}, A::AbstractArray{T,4}, n) → (TensWalpole{T,5}, d, drel)

Project a 4th-order tensor A (3×3×3×3) onto the transversely isotropic subspace with fixed symmetry axis n. Returns a major-symmetric TensWalpole{T,5}.

The projection minimises the Frobenius distance ‖B − A‖ over all TI tensors B with axis n.

Returns a 3-tuple (B, d, drel):

• B: the projected TensWalpole{T,5}
• d: absolute Frobenius distance ‖B − A‖
• drel: relative distance d / ‖A‖

Examples

julia> n = [0., 0., 1.];

julia> C = tensTI(10., 3., 2.5, 12., 2., n);

julia> B, d, drel = proj_tens(:TI, getarray(C), n);

julia> d < 1e-12
true

julia> argTI(B) == argTI(C)
true

proj_tens(::Val{:TI}, A::AbstractArray{T,4}) where {T<:AbstractFloat}

Find the best TI approximation of A by optimizing over all possible symmetry axes. Requires the NLopt package: using NLopt.

See also proj_tens(::Val{:TI}, A, n) for fixed-axis projection.

proj_tens(::Val{:TI}, A::AbstractArray{T,2}, n) → (TensTI{2,T,2}, d, drel)

Project a 2nd-order tensor A (3×3) onto the transversely isotropic subspace with fixed symmetry axis n. Returns a TensTI{2,T,2}.

In the rotated frame where e₃ = n:

• a = (M[1,1] + M[2,2]) / 2 (transverse)
• b = M[3,3] (axial)

Examples

julia> n = [0., 0., 1.];

julia> A = [5. 0 0; 0 5 0; 0 0 8];

julia> B, d, drel = proj_tens(:TI, A, n);

julia> B.data
(5.0, 8.0)

proj_tens(::Val{:TI}, A::AbstractArray{T,2}) where {T<:AbstractFloat}

Find the best TI approximation of a 2nd-order tensor A by optimizing the symmetry axis. Requires the NLopt package: using NLopt.

qcontract(t1::AbstractTens{order1,dim}, t2::AbstractTens{order2,dim})

Define a quadruple contracted product between two tensors

𝔸 ⊙ 𝔹 = AᵢⱼₖₗBⁱʲᵏˡ

Examples

julia> 𝕀 = t𝕀(Sym) ; 𝕁 = t𝕁(Sym) ; 𝕂 = t𝕂(Sym) ;

julia> 𝕀 ⊙ 𝕀
6

julia> 𝕁 ⊙ 𝕀
1

julia> 𝕂 ⊙ 𝕀
5

julia> 𝕂 ⊙ 𝕁
0

rot2(θ)

Return a 2D rotation matrix with respect to the angle θ

Examples

julia> rot2(θ)
2×2 Tensor{2, 2, Sym, 4}:
 cos(θ)  -sin(θ)
 sin(θ)   cos(θ)

rot3(θ, ϕ = 0, ψ = 0)

Return a rotation matrix with respect to the 3 Euler angles θ, ϕ, ψ

Examples

julia> cθ, cϕ, cψ, sθ, sϕ, sψ = symbols("cθ cϕ cψ sθ sϕ sψ", real = true) ;

julia> d = Dict(cos(θ) => cθ, cos(ϕ) => cϕ, cos(ψ) => cψ, sin(θ) => sθ, sin(ϕ) => sϕ, sin(ψ) => sψ) ;

julia> subs.(rot3(θ, ϕ, ψ),d...)
3×3 StaticArrays.SMatrix{3, 3, Sym, 9} with indices SOneTo(3)×SOneTo(3):
 cθ⋅cψ⋅cϕ - sψ⋅sϕ  -cθ⋅cϕ⋅sψ - cψ⋅sϕ  cϕ⋅sθ
 cθ⋅cψ⋅sϕ + cϕ⋅sψ  -cθ⋅sψ⋅sϕ + cψ⋅cϕ  sθ⋅sϕ
           -cψ⋅sθ              sθ⋅sψ     cθ

rot6(θ, ϕ = 0, ψ = 0)

Return a rotation matrix with respect to the 3 Euler angles θ, ϕ, ψ

Examples

julia> cθ, cϕ, cψ, sθ, sϕ, sψ = symbols("cθ cϕ cψ sθ sϕ sψ", real = true) ;

julia> d = Dict(cos(θ) => cθ, cos(ϕ) => cϕ, cos(ψ) => cψ, sin(θ) => sθ, sin(ϕ) => sϕ, sin(ψ) => sψ) ;

julia> R = Tens(subs.(rot3(θ, ϕ, ψ),d...))
Tens.TensCanonical{2, 3, Sym, Tensor{2, 3, Sym, 9}}
# data: 3×3 Tensor{2, 3, Sym, 9}:
 cθ⋅cψ⋅cϕ - sψ⋅sϕ  -cθ⋅cϕ⋅sψ - cψ⋅sϕ  cϕ⋅sθ
 cθ⋅cψ⋅sϕ + cϕ⋅sψ  -cθ⋅sψ⋅sϕ + cψ⋅cϕ  sθ⋅sϕ
           -cψ⋅sθ              sθ⋅sψ     cθ
# var: (:cont, :cont)
# basis: 3×3 Tens.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1

julia> RR = R ⊠ˢ R
Tens.TensCanonical{4, 3, Sym, SymmetricTensor{4, 3, Sym, 36}}
# data: 6×6 Matrix{Sym}:
                          (cθ*cψ*cϕ - sψ*sϕ)^2                            (-cθ*cϕ*sψ - cψ*sϕ)^2           cϕ^2*sθ^2                      √2⋅cϕ⋅sθ⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)                     √2⋅cϕ⋅sθ⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)                                   √2⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)
                          (cθ*cψ*sϕ + cϕ*sψ)^2                            (-cθ*sψ*sϕ + cψ*cϕ)^2           sθ^2*sϕ^2                      √2⋅sθ⋅sϕ⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ)                     √2⋅sθ⋅sϕ⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)                                   √2⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ)
                                     cψ^2*sθ^2                                        sθ^2*sψ^2                cθ^2                                       √2⋅cθ⋅sθ⋅sψ                                    -√2⋅cθ⋅cψ⋅sθ                                                              -sqrt(2)*cψ*sθ^2*sψ
             -√2⋅cψ⋅sθ⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)                √2⋅sθ⋅sψ⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ)    √2⋅cθ⋅sθ⋅sϕ                    cθ*(-cθ*sψ*sϕ + cψ*cϕ) + sθ^2*sψ*sϕ                   cθ*(cθ*cψ*sϕ + cϕ*sψ) - cψ*sθ^2*sϕ                            -cψ⋅sθ⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ) + sθ⋅sψ⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)
             -√2⋅cψ⋅sθ⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)                √2⋅sθ⋅sψ⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)    √2⋅cθ⋅cϕ⋅sθ                    cθ*(-cθ*cϕ*sψ - cψ*sϕ) + cϕ*sθ^2*sψ                   cθ*(cθ*cψ*cϕ - sψ*sϕ) - cψ*cϕ*sθ^2                            -cψ⋅sθ⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ) + sθ⋅sψ⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)
 √2⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)  √2⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ)  sqrt(2)*cϕ*sθ^2*sϕ  cϕ⋅sθ⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ) + sθ⋅sϕ⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)  cϕ⋅sθ⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ) + sθ⋅sϕ⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)  (cθ*cψ*cϕ - sψ*sϕ)*(-cθ*sψ*sϕ + cψ*cϕ) + (cθ*cψ*sϕ + cϕ*sψ)*(-cθ*cϕ*sψ - cψ*sϕ)
# var: (:cont, :cont, :cont, :cont)
# basis: 3×3 Tens.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1

julia> R6 = invKM(subs.(KM(rot6(θ, ϕ, ψ)),d...))
Tens.TensCanonical{4, 3, Sym, SymmetricTensor{4, 3, Sym, 36}}
# data: 6×6 Matrix{Sym}:
                          (cθ*cψ*cϕ - sψ*sϕ)^2                            (-cθ*cϕ*sψ - cψ*sϕ)^2           cϕ^2*sθ^2                      √2⋅cϕ⋅sθ⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)                     √2⋅cϕ⋅sθ⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)                                   √2⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)
                          (cθ*cψ*sϕ + cϕ*sψ)^2                            (-cθ*sψ*sϕ + cψ*cϕ)^2           sθ^2*sϕ^2                      √2⋅sθ⋅sϕ⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ)                     √2⋅sθ⋅sϕ⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)                                   √2⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ)
                                     cψ^2*sθ^2                                        sθ^2*sψ^2                cθ^2                                       √2⋅cθ⋅sθ⋅sψ                                    -√2⋅cθ⋅cψ⋅sθ                                                              -sqrt(2)*cψ*sθ^2*sψ
             -√2⋅cψ⋅sθ⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)                √2⋅sθ⋅sψ⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ)    √2⋅cθ⋅sθ⋅sϕ                    cθ*(-cθ*sψ*sϕ + cψ*cϕ) + sθ^2*sψ*sϕ                   cθ*(cθ*cψ*sϕ + cϕ*sψ) - cψ*sθ^2*sϕ                            -cψ⋅sθ⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ) + sθ⋅sψ⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)
             -√2⋅cψ⋅sθ⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)                √2⋅sθ⋅sψ⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)    √2⋅cθ⋅cϕ⋅sθ                    cθ*(-cθ*cϕ*sψ - cψ*sϕ) + cϕ*sθ^2*sψ                   cθ*(cθ*cψ*cϕ - sψ*sϕ) - cψ*cϕ*sθ^2                            -cψ⋅sθ⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ) + sθ⋅sψ⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)
 √2⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ)  √2⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ)  sqrt(2)*cde Liv Lehn ϕ*sθ^2*sϕ  cϕ⋅sθ⋅(-cθ⋅sψ⋅sϕ + cψ⋅cϕ) + sθ⋅sϕ⋅(-cθ⋅cϕ⋅sψ - cψ⋅sϕ)  cϕ⋅sθ⋅(cθ⋅cψ⋅sϕ + cϕ⋅sψ) + sθ⋅sϕ⋅(cθ⋅cψ⋅cϕ - sψ⋅sϕ)  (cθ*cψ*cϕ - sψ*sϕ)*(-cθ*sψ*sϕ + cψ*cϕ) + (cθ*cψ*sϕ + cϕ*sψ)*(-cθ*cϕ*sψ - cψ*sϕ)
# var: (:cont, :cont, :cont, :cont)
# basis: 3×3 Tens.LazyIdentity{3, Sym}:
 1  0  0
 0  1  0
 0  0  1

julia> R6 == RR
true

sotimes(t1::AbstractTens{order1,dim}, t2::AbstractTens{order2,dim})

Define a symmetric tensor product between two tensors

(aⁱeᵢ) ⊗ˢ (bʲeⱼ) = 1/2(aⁱbʲ + aʲbⁱ) eᵢ⊗eⱼ

tensId2(::Val{dim}, ::Val{T}) where {dim,T<:Number}

Identity tensor of second order 𝟏ᵢⱼ = δᵢⱼ = 1 if i=j otherwise 0

Examples

julia> 𝟏 = t𝟏() ; KM(𝟏)
6-element Vector{Sym}:
 1
 1
 1
 0
 0
 0

julia> 𝟏.data
3×3 SymmetricTensor{2, 3, Sym, 6}:
 1  0  0
 0  1  0
 0  0  1

tensId4(::Val{dim} = Val(3), ::Val{T} = Val(Sym))

Symmetric identity tensor of fourth order 𝕀 = 𝟏 ⊠ˢ 𝟏 i.e. (𝕀)ᵢⱼₖₗ = (δᵢₖδⱼₗ+δᵢₗδⱼₖ)/2

Examples

julia> 𝕀 = t𝕀() ; KM(𝕀)
6×6 Matrix{Sym}:
 1  0  0  0  0  0
 0  1  0  0  0  0
 0  0  1  0  0  0
 0  0  0  1  0  0
 0  0  0  0  1  0
 0  0  0  0  0  1

tensJ4(::Val{dim} = Val(3), ::Val{T} = Val(Sym))

Spherical projector of fourth order 𝕁 = (𝟏 ⊗ 𝟏) / dim i.e. (𝕁)ᵢⱼₖₗ = δᵢⱼδₖₗ/dim

Examples

julia> 𝕁 = t𝕁() ; KM(𝕁)
6×6 Matrix{Sym}:
 1/3  1/3  1/3  0  0  0
 1/3  1/3  1/3  0  0  0
 1/3  1/3  1/3  0  0  0
   0    0    0  0  0  0
   0    0    0  0  0  0
   0    0    0  0  0  0

tensK4(::Val{dim} = Val(3), ::Val{T} = Val(Sym))

Deviatoric projector of fourth order 𝕂 = 𝕀 - 𝕁 i.e. (𝕂)ᵢⱼₖₗ = (δᵢₖδⱼₗ+δᵢₗδⱼₖ)/2 - δᵢⱼδₖₗ/dim

Examples

julia> 𝕂 = t𝕂() ; KM(𝕂)
6×6 Matrix{Sym}:
  2/3  -1/3  -1/3  0  0  0
 -1/3   2/3  -1/3  0  0  0
 -1/3  -1/3   2/3  0  0  0
    0     0     0  1  0  0
    0     0     0  0  1  0
    0     0     0  0  0  1

tensTI(C₁₁₁₁, C₁₁₂₂, C₁₁₃₃, C₃₃₃₃, C₂₃₂₃, n) → TensWalpole{T,5}

Construct a major-symmetric TI 4th-order tensor from its 5 independent components and symmetry axis n. Works for both stiffness and compliance tensors (the formula is the same).

Walpole coefficients:

• ℓ₁ = C₃₃₃₃
• ℓ₂ = C₁₁₁₁ + C₁₁₂₂
• ℓ₃ = √2 C₁₁₃₃
• ℓ₅ = C₁₁₁₁ − C₁₁₂₂
• ℓ₆ = 2 C₂₃₂₃

See also argTI, tensTI_eng.

tensTI_Hoenig(E, ν₁, ν₂, H, Γ, n) → TensWalpole{T,5}

Construct the TI compliance tensor from 5 Hoenig parameters (Hoenig, 1978) and symmetry axis n.

• E : transverse Young's modulus (= 1/S₁₁₁₁)
• ν₁ : in-plane Poisson's ratio (= −E S₁₁₂₂)
• ν₂ : axial-transverse Poisson's ratio (= −E S₁₁₃₃)
• H : axial-to-transverse modulus ratio (= 1/(E S₃₃₃₃))
• Γ : shear anisotropy parameter (= (1+ν₁)/(2 E S₂₃₂₃))

Compliance components:

• S₁₁₁₁ = 1/E
• S₁₁₂₂ = −ν₁/E
• S₁₁₃₃ = −ν₂/E
• S₃₃₃₃ = 1/(E H)
• S₂₃₂₃ = (1+ν₁)/(2 E Γ)

To obtain the stiffness tensor, invert the result: inv(tensTI_Hoenig(…)).

See also argTI_Hoenig, tensTI_eng, tensTI.

tensTI_eng(E₁, E₃, ν₁₂, ν₃₁, G₃₁, n) → TensWalpole{T,5}

Construct the TI compliance tensor from 5 engineering constants and symmetry axis n.

• E₁ : transverse Young's modulus (isotropic plane)
• E₃ : axial Young's modulus (symmetry axis)
• ν₁₂: in-plane Poisson's ratio
• ν₃₁: axial-transverse Poisson's ratio (ν₃₁/E₃ = ν₁₃/E₁)
• G₃₁: axial shear modulus

To obtain the stiffness tensor, invert the result: inv(tensTI_eng(…)).

See also argTI_eng, tensTI.

tensW1(n) → TensWalpole{T,6}   (W₁ = nₙ⊗nₙ, coeffs (1,0,0,0,0,0))

tensW2(n) → TensWalpole{T,6}   (W₂ = (nT⊗nT)/2, coeffs (0,1,0,0,0,0))

tensW3(n) → TensWalpole{T,6}   (W₃ = (nₙ⊗nT)/√2, coeffs (0,0,1,0,0,0))

tensW4(n) → TensWalpole{T,6}   (W₄ = (nT⊗nₙ)/√2, coeffs (0,0,0,1,0,0))

tensW5(n) → TensWalpole{T,6}   (W₅ = nT⊠ˢnT − (nT⊗nT)/2, coeffs (0,0,0,0,1,0))

tensW6(n) → TensWalpole{T,6}   (W₆ = nT⊠ˢnₙ + nₙ⊠ˢnT, coeffs (0,0,0,0,0,1))

unitvec(CS::CoorSystemNum, x₀, i) → AbstractTens

Return the i-th unit vector of the normalized basis at point x₀.

vecbasis(ℬ::AbstractBasis, var = :cov)

Return the primal (if var = :cov) or dual (if var = :cont) basis

∂(f::Function, i::Integer, CS::CoorSystemNum, x₀)

Covariant derivative of f with respect to the i-th coordinate, evaluated at x₀.

f must accept an AbstractVector of length dim and return either a scalar or a plain Array of physical (normalized-frame) components of the tensor field. Returns a plain scalar or plain Array (suitable for further ForwardDiff differentiation).

See also GRAD, DIV, LAPLACE, HESS, SYMGRAD for high-level operators that wrap results in AbstractTens.

∂(t::AbstractTens{order,dim,T,A},xᵢ::T)

Return the derivative of the tensor t with respect to the variable x_i

Examples

julia> (θ, ϕ, r), (𝐞ᶿ, 𝐞ᵠ, 𝐞ʳ), ℬˢ = init_spherical() ;

julia> ∂(𝐞ʳ, ϕ) == sin(θ) * 𝐞ᵠ
true

julia> ∂(𝐞ʳ ⊗ 𝐞ʳ,θ)
Tens.TensRotated{2, 3, Sym, SymmetricTensor{2, 3, Sym, 6}}
# data: 3×3 SymmetricTensor{2, 3, Sym, 6}:
 0  0  1
 0  0  0
 1  0  0
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 cos(θ)⋅cos(ϕ)  -sin(ϕ)  sin(θ)⋅cos(ϕ)
 sin(ϕ)⋅cos(θ)   cos(ϕ)  sin(θ)⋅sin(ϕ)
       -sin(θ)        0         cos(θ)
# var: (:cont, :cont)

𝐞(i::Integer, dim::Int = 3, T::Type{<:Number} = Sym)

Vector of the canonical basis

Examples

julia> 𝐞(1)
Tens{1, 3, Sym, Sym, Vec{3, Sym}, CanonicalBasis{3, Sym}}
# data: 3-element Vec{3, Sym}:
 1
 0
 0
# var: (:cont,)
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 1  0  0
 0  1  0
 0  0  1

𝐞ˢ(i::Integer, θ::T = zero(Sym), ϕ::T = zero(Sym), ψ::T = zero(Sym); canonical = false)

Vector of the basis rotated with the 3 Euler angles θ, ϕ, ψ (spherical if ψ=0)

Examples

julia> θ, ϕ, ψ = symbols("θ, ϕ, ψ", real = true) ;

Tens{1, 3, Sym, Sym, Vec{3, Sym}, RotatedBasis{3, Sym}}
# data: 3-element Vec{3, Sym}:
 1
 0
 0
# var: (:cont,)
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 -sin(ψ)⋅sin(ϕ) + cos(θ)⋅cos(ψ)⋅cos(ϕ)  -sin(ψ)⋅cos(θ)⋅cos(ϕ) - sin(ϕ)⋅cos(ψ)  sin(θ)⋅cos(ϕ)
  sin(ψ)⋅cos(ϕ) + sin(ϕ)⋅cos(θ)⋅cos(ψ)  -sin(ψ)⋅sin(ϕ)⋅cos(θ) + cos(ψ)⋅cos(ϕ)  sin(θ)⋅sin(ϕ)
                        -sin(θ)⋅cos(ψ)                          sin(θ)⋅sin(ψ)         cos(θ)

𝐞ᵖ(i::Integer, θ::T = zero(Sym); canonical = false)

Vector of the polar basis

Examples

julia> θ = symbols("θ", real = true) ;

julia> 𝐞ᵖ(1, θ)
Tens{1, 2, Sym, Sym, Vec{2, Sym}, RotatedBasis{2, Sym}}
# data: 2-element Vec{2, Sym}:
 1
 0
# var: (:cont,)
# basis: 2×2 Tensor{2, 2, Sym, 4}:
 cos(θ)  -sin(θ)
 sin(θ)   cos(θ)

𝐞ᶜ(i::Integer, θ::T = zero(Sym); canonical = false)

Vector of the cylindrical basis

Examples

julia> θ = symbols("θ", real = true) ;

julia> 𝐞ᶜ(1, θ)
Tens{1, 3, Sym, Sym, Vec{3, Sym}, RotatedBasis{3, Sym}}
# data: 3-element Vec{3, Sym}:
 1
 0
 0
# var: (:cont,)
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 cos(θ)  -sin(θ)  0
 sin(θ)   cos(θ)  0
      0        0  1

dcontract(A::TensWalpole, B::TensISO{4,3}) → TensWalpole{T,6}
dcontract(A::TensISO{4,3}, B::TensWalpole) → TensWalpole{T,6}

dcontract(A::TensWalpole, B::TensWalpole) → TensWalpole{T,6}

Product rule via 2×2 matrix product + scalar products for ℓ₅, ℓ₆. Always returns N=6 since the product of two symmetric tensors need not be symmetric.

dcontract(t1::AbstractTens{order1,dim}, t2::AbstractTens{order2,dim})

Define a double contracted product between two tensors

𝛔 ⊡ 𝛆 = σⁱʲεᵢⱼ 𝛔 = ℂ ⊡ 𝛆

Examples

julia> 𝛆 = Tens(SymmetricTensor{2,3}((i, j) -> symbols("ε$i$j", real = true))) ;

julia> k, μ = symbols("k μ", real =true) ;

julia> ℂ = 3k * t𝕁() + 2μ * t𝕂() ;

julia> 𝛔 = ℂ ⊡ 𝛆
Tens{2, 3, Sym, Sym, SymmetricTensor{2, 3, Sym, 6}, CanonicalBasis{3, Sym}}
# data: 3×3 SymmetricTensor{2, 3, Sym, 6}:
 ε11*(k + 4*μ/3) + ε22*(k - 2*μ/3) + ε33*(k - 2*μ/3)                                              2⋅ε21⋅μ                                              2⋅ε31⋅μ
                                             2⋅ε21⋅μ  ε11*(k - 2*μ/3) + ε22*(k + 4*μ/3) + ε33*(k - 2*μ/3)                                              2⋅ε32⋅μ
                                             2⋅ε31⋅μ                                              2⋅ε32⋅μ  ε11*(k - 2*μ/3) + ε22*(k - 2*μ/3) + ε33*(k + 4*μ/3)
# var: (:cont, :cont)
# basis: 3×3 Tensor{2, 3, Sym, 9}:
 1  0  0
 0  1  0
 0  0  1

dotdot(v1::AbstractTens{order1,dim}, S::AbstractTens{orderS,dim}, v2::AbstractTens{order2,dim})

Define a bilinear operator 𝐯₁⋅𝕊⋅𝐯₂

Examples

julia> n = Tens(Sym[0, 0, 1]) ;

julia> k, μ = symbols("k μ", real =true) ;

julia> ℂ = 3k * t𝕁() + 2μ * t𝕂() ;

julia> dotdot(n,ℂ,n) # Acoustic tensor
3×3 Tens{2, 3, Sym, Sym, Tensor{2, 3, Sym, 9}, CanonicalBasis{3, Sym}}:
 μ  0          0
 0  μ          0
 0  0  k + 4*μ/3

otimes(A::TensTI{2}, B::TensISO{2,3}) → TensWalpole{T,6}

Tensor product of a TI 2nd-order tensor with a 3D isotropic 2nd-order tensor. The isotropic tensor λ·𝟏 is treated as TensTI{2}(λ,λ,n) with the axis of A.

otimes(A::TensTI{2}, B::TensTI{2}) → TensWalpole{T,6}

Tensor product of two TI 2nd-order tensors with the same axis. Falls back to generic otimes if axes differ.

(a₁·nT + b₁·nₙ) ⊗ (a₂·nT + b₂·nₙ)
= b₁b₂·W₁ + 2a₁a₂·W₂ + √2·b₁a₂·W₃ + √2·a₁b₂·W₄

otimes(A::TensISO{2,3}, B::TensTI{2}) → TensWalpole{T,6}

Tensor product of a 3D isotropic 2nd-order tensor with a TI 2nd-order tensor. The isotropic tensor λ·𝟏 is treated as TensTI{2}(λ,λ,n) with the axis of B.

otimes(A::TensTI{2}) → TensWalpole{T,5}

Self tensor product of a TI 2nd-order tensor. The result is always major-symmetric (ℓ₃ = ℓ₄) and lives in the Walpole basis with N=5.

(a·nT + b·nₙ) ⊗ (a·nT + b·nₙ)
= b²W₁ + 2a²W₂ + √2·ab·(W₃+W₄)

otimes(t1::AbstractTens{order1,dim}, t2::AbstractTens{order2,dim})

Define a tensor product between two tensors

(aⁱeᵢ) ⊗ (bʲeⱼ) = aⁱbʲ eᵢ⊗eⱼ

otimesu(t1::AbstractTens{order1,dim}, t2::AbstractTens{order2,dim})

Define a special tensor product between two tensors of at least second order

(𝐚 ⊠ 𝐛) ⊡ 𝐩 = 𝐚⋅𝐩⋅𝐛 = aⁱᵏbʲˡpₖₗ eᵢ⊗eⱼ

@set_coorsys CS
@set_coorsys(CS)

Set a coordinate system in order to avoid precising it in differential operators

Examples

julia> Spherical = coorsys_spherical() ; θ, ϕ, r = getcoords(Spherical) ; 𝐞ᶿ, 𝐞ᵠ, 𝐞ʳ = unitvec(Spherical) ; vec = ("𝐞ᶿ", "𝐞ᵠ", "𝐞ʳ") ;

julia> @set_coorsys Spherical

julia> intrinsic(GRAD(𝐞ʳ),vec)
(1/r)𝐞ᶿ⊗𝐞ᶿ + (1/r)𝐞ᵠ⊗𝐞ᵠ

julia> intrinsic(DIV(𝐞ʳ ⊗ 𝐞ʳ),vec)
(2/r)𝐞ʳ

julia> LAPLACE(1/r)
0