Documentation of the FEMTISE module.

FEMTISE.aprox_dirac_deltaMethod
aprox_dirac_delta(x,params)

Aim

  • Approximation of Dirac delta function as rectangular function

Arguments

  • x: independent variable (is a FE coordinate object)
  • params::Tuple{Float64,Float64,Int,Float64}: parameters
    • x₀::Float64: specific point to centre dirac function
    • δnorm::Float64: norm value to obtain normalized delta function
    • component::Int: component of x coordinate
    • Δx::Float64: thickness of rectangular function

Output

  • result::Vector: integral value

Example

x=0:0.1:10;
params=(5.0,1.0,1,0.1);
result=aprox_dirac_delta(x,params)
source
FEMTISE.bilineal_formsMethod
bilineal_forms(p,q,r,dΩ)

Aim

  • Returns bilineals forms (a(u,v) and b(u,v)) for eigenvalues 1D or 2D (equal masses)

Arguments

  • p::Function: kinetic energy function from Sturm Liouville_Formulation.
  • q::Function: potential energy function from Sturm Liouville_Formulation.
  • r::Function: weight or density function from Sturm Liouville_Formulation.
  • dΩ::Gridap.CellData.GenericMeasure: differential FE domain
source
FEMTISE.build_input_dataMethod
build_input_data(full_path_name)

Aim

  • Build input data from a file

Arguments

  • full_path_name::String: full path name to input data file
source
FEMTISE.collect_result_dataMethod
collect_result_data(switch_input_file,full_path_name)

Aim

  • Collect result data

Arguments

  • switch_input_file::Bool: switch input file
  • full_path_name::String: full path name
source
FEMTISE.coord_first_momentMethod
coord_first_moment(psi,TrialSpace,Omega,dOmega,x_component)

Aim

  • Compute first moment of coordinate (specific component)

Arguments

  • psi::Vector{CellField}: array of FE wave functions
  • TrialSpace::FESpace: Trial finite element space
  • Omega: Finite element domain
  • dOmega::Gridap.CellData.GenericMeasure: integration domain
  • x_component::Int: component of specific coordinate variable

Output

  • expval::Array{Float64}: array of expectation values of specific coordinate
source
FEMTISE.coord_second_momentMethod
coord_second_moment(psi,TrialSpace,Omega,dOmega,x_component)

Aim

  • Compute first moment of coordinate (specific component) or position variance

Arguments

  • psi::Vector{CellField}: array of FE wave functions
  • TrialSpace::FESpace: Trial finite element space
  • Omega: Finite element domain
  • dOmega::Gridap.CellData.GenericMeasure: integration domain
  • x_component::Int: component of specific coordinate variable

Output

  • var::Array{Float64}: array of variance values of specific coordinate
source
FEMTISE.coord_standar_deviationMethod
coord_standar_deviation(expval,var)

Aim

  • Compute standar deviation of coordinate (specific component)

Arguments

  • expval::Array{Float64}: array of expectation values of specific coordinate
  • var::Array{Float64}: array of variance values of specific coordinate

Output

  • sigma::Array{Float64}: array of standar deviation of specific coordinate

Example

expval=[1.0,2.0,3.0];
var=[0.1,0.2,0.3];
sigma=coord_standar_deviation(expval,var)
source
FEMTISE.coord_standar_deviationMethod
coord_standar_deviation(psi,TrialSpace,Omega,dOmega,x_component)

Aim

  • Compute standar deviation of coordinate (specific component)

Arguments

  • psi::Vector{CellField}: array of FE wave functions
  • TrialSpace::FESpace: Trial finite element space
  • Omega: Finite element domain
  • dOmega::Gridap.CellData.GenericMeasure: integration domain
  • x_component::Int: component of specific coordinate variable

Output

  • sigma::Array{Float64}: array of standar deviation of specific coordinate
source
FEMTISE.default_solver_eigen_problemMethod
default_solver_eigen_problem(params)

Aim

  • Function to resolve bidimensonal eigen problem over cartesian grid

Arguments

  • params::Params2D: parameters fo 2D potential
  • different_masses::Tuple: keyword to specify if we want to simulate two particles with diferent masses
source
FEMTISE.densityMethod
density(phi)

Aim

  • Compute probability density from set of specific wave functions (FE object)

Arguments

  • phi::Vector{CellField}: array of eigenstates (array of FE objects)

Output

  • rho::Vector{CellField}: Array of probability densities (array of FE objects)
source
FEMTISE.eigen_problemMethod
eigen_problem(weakform_k,weakform_m,test,trial; <keyword arguments>)

Aim

  • Define eigen problem as an input to solve function where we compute eigen problem by Arpack eigs function.

Arguments

  • weakform_k::Function: forma bilineal lado izquierdo de la formulación débil
  • weakform_m::Function: forma bilineal lado derecho de la formulación débil
  • test::FESpace: espacio de prueba, puede ser MultiFieldFESpace
  • trial::FESpace: espacio de solución, puede ser MultiFieldFESpace
  • nev::Int=10: número de autovalores requeridos
  • tol::Float64=10e-6: relative tolerance for convergence of Ritz values
  • maxiter::Integer=100: maximum number of iterations
  • explicittransform::Symbol=:none: shift and invert should be explicitly invoked in julia code
    • =:auto:
    • =:shiftinvert:
  • sigma::Float64=1.0: the level shift used in inverse iteration.
  • which::Symbol=:LM: eigenvalues of largest magnitude (default)
    • =:SM: eigenvalues of smallest magnitude
    • =:LR: eigenvalues of largest real part
    • =:SR: eigenvalues of smallest real part
    • =:LI: eigenvalues of largest imaginary part (nonsymmetric or complex matrix only)
    • =:SI: eigenvalues of smallest imaginary part (nonsymmetric or complex matrix only)
    • =:BE: compute half of the eigenvalues from each end of the spectrum, biased in favor of the high end. (real symmetric matrix only)

Returns

  • EigenProblem: problem definition
source
FEMTISE.eigen_values_and_eigen_vectorsMethod
eigen_values_and_eigen_vectors(p₁,p₂,q,r,dΩ,USpace,VSpace; <keyword arguments>)

Aim

  • Compute matrix eigen problem and return eigenvalues and eigenvectors

Arguments

  • p₁,p₂,q,r::Function: functions from Sturm Liouville formulation
  • dΩ::Gridap.CellData.GenericMeasure: measure of FE grid
  • USpace::FESpace: trial FE Space
  • VSpace::FESpace: test FE Space
  • params::Tuple=(nev,tol,maxiter,explicittransform,sigma): params to Arpack eigs function.
    • nev::Integer=10: quantity of eigendata to calculate
    • tol::Float64=10e-9: relative tolerance for convergence of Ritz values
    • maxiter::Integer=100: maximum number of iterations
    • explicittransform::Symbol=:none: shift and invert should be explicitly invoked in julia code
    • sigma::Float64=1.0: the level shift used in inverse iteration.

Example

using Gridap
using FEMTISE
using LinearAlgebra

# Define functions
p₁(x)=1.0
p₂(x)=0.0
q(x)=0.0
r(x)=1.0
# Define measure
Ω=1.0:0.1:10.0
dΩ=Gridap.CellData.CellData(Ω)
# Define FE Spaces
USpace=fe_space(1,1,Ω)
VSpace=fe_space(1,1,Ω)
# Compute eigenvalues and eigenvectors
ϵ,ϕ=eigen_values_and_eigen_vectors(p₁,p₂,q,r,dΩ,USpace,VSpace)
source
FEMTISE.eigen_values_and_eigen_vectorsMethod
eigen_values_and_eigen_vectors(p,q,r,dΩ,USpace,VSpace; <keyword arguments>)

Aim

  • Compute matrix eigen problem and return eigenvalues and eigenvectors

Arguments

  • p,q,r::Function: functions from Sturm Liouville formulation
  • dΩ::Gridap.CellData.GenericMeasure: measure of FE grid
  • USpace::FESpace: trial FE Space
  • VSpace::FESpace: test FE Space
  • params::Tuple=(nev,tol,maxiter,explicittransform,sigma): params to Arpack eigs function.
    • nev::Integer=10: quantity of eigendata to calculate
    • tol::Float64=10e-9: relative tolerance for convergence of Ritz values
    • maxiter::Integer=100: maximum number of iterations
    • explicittransform::Symbol=:none: shift and invert should be explicitly invoked in julia code
    • sigma::Float64=1.0: the level shift used in inverse iteration.

Example

using Gridap
using FEMTISE
using LinearAlgebra

# Define functions
p(x)=1.0
q(x)=0.0
r(x)=1.0
# Define measure
Ω=1.0:0.1:10.0
dΩ=Gridap.CellData.CellData(Ω)
# Define FE Spaces
USpace=fe_space(1,1,Ω)
VSpace=fe_space(1,1,Ω)
# Compute eigenvalues and eigenvectors
ϵ,ϕ=eigen_values_and_eigen_vectors(p,q,r,dΩ,USpace,VSpace)
source
FEMTISE.eigenstates_normalizationMethod

eigenstates_normalization(ϕ,dΩ)

Aim

  • Check eigenstates normalization.

Arguments

  • ϕ::Vector{CellField}: array of FE eigenstates
  • dΩ::Gridap.CellData.GenericMeasure: integration domain

Output

  • nom_vector::Vector{Float64}: array of norm value for specific eigenstate
source
FEMTISE.fe_spacesMethod
fe_spaces(model,reff,grid_type; <keyword arguments>)

Aim

  • Create finite element (FE) spaces (Trial and Test spaces).

Arguments

  • BC_type::String="FullDirichlet": the type of boundary condition.
  • TypeData::Type=ComplexF64: the type of data to define FE spaces.
  • conf_type::Symbol=:H1: the regularity of the interpolation at the boundaries of cells in the mesh. (e.g.:L2,:H1,:C0,:Hgrad,)
source
FEMTISE.input_dataMethod
input_data(data_file_name)

Aim

  • Definition of input data form input.dat file using specific type structures.

Arguments

  • data_file_name::String: name of input data file
source
FEMTISE.interpolation_eigenstates!Method
interpolation_eigenstates!(eigen_states,USpace)

Aim

  • Interpolate eigenstates

Arguments

  • eigen_states::Vector{CellField}: eigenstates
  • USpace::FESpace: FE space
source
FEMTISE.load_input_dataMethod
load_input_data(full_path_input_data)

Aim

  • Load input data from a file

Arguments

  • full_path_input_data::String: full path name to input data file
source
FEMTISE.make_boundary_conditionsMethod
make_boundary_conditions(grid_type,BC_type,TypeData;<keyword arguments>)

Aim

  • Create boundary conditions for a given grid and type of boundary conditions

Arguments

  • grid_type::String: type of grid
  • BC_type::String: type of boundary conditions
  • TypeData::Type: type of data
  • <keyword arguments>:
    • homogeneous::Bool=true: if boundary conditions are homogeneous
    • params::Tuple=(nothing,nothing): parameters for boundary conditions

Returns

  • BC_values::Array{TypeData,1}: values of boundary conditions
  • BC_tags::Array{String,1}: tags of boundary conditions

Example

BC_values,BC_tags=make_boundary_conditions("simple_line","FullDirichlet",Float64)
source
FEMTISE.make_modelMethod
make_model(grid_type,params)

Aim

  • Create a mesh model using Gmsh

Arguments

  • grid_type::String: type of grid to create
  • params::Tuple: parameters to create the grid

Returns

  • model::GmshDiscreteModel: mesh model

References

  • https://gitlab.onelab.info/gmsh/gmsh/blob/master/api/gmsh.jl
source
FEMTISE.measuresMethod
measures(model,degree,tags_boundary)

Aim

The triangulation and integration aproximated Lebesgue measures

Arguments

  • model: FE grid model.
  • degree::Integer: degree of quadrature rule to use in the cells of triangulation.
  • tags_boundary: tags of boundary conditions.
source
FEMTISE.norm_l2Method
norm_l2(𝜳,dΩ)

Aim

  • Compute de L2 norm for specific FE wave function

Arguments

  • 𝜳::CellField: specific FE wave function
  • dΩ::Gridap.CellData.GenericMeasure: integration domain

Output

  • orthogonality_vector::Real: L2 norm value
source
FEMTISE.orthogonality_checkMethod

orthogonality_check(ϕ,dΩ,TrialSpace)

Aim

  • Check eigenstates orthogonality. Each eigenstate is a multifield FE object.

Arguments

  • ϕ::Vector{CellField}: array of FE eigenstates
  • dΩ::Gridap.CellData.GenericMeasure: integration domain
  • TrialSpace::FESpace: Trial FE space

Output

  • orthogonality_vector::Vector{Float64}: array of inner product between differents eigenstates
source
FEMTISE.orthogonality_checkMethod

orthogonality_check(ϕ,dΩ)

Aim

  • Check eigenstates orthogonality

Arguments

  • ϕ::Vector{CellField}: array of FE eigenstates
  • dΩ::Gridap.CellData.GenericMeasure: integration domain

Output

  • orthogonality_vector::Vector{Float64}: array of inner product between differents eigenstates
source
FEMTISE.read_binMethod
read_bin(fileName;<keyword arguments>)

Aim

  • Read binary file

Arguments

  • fileName::String: name of file data
  • matrix_data::Bool: optional boolean keyword to specify matrix or array data
  • c_dim::Int: optional column number of matrix data
source
FEMTISE.read_binMethod
read_bin(io,::Type{T},n,matrix_data,c_dim)

Aim

  • Speeds up the read binary file

Arguments

  • io::IO: in/output variable
  • ::Type{T}: data type
  • n::Int: total number of elements in array/matrix data
  • matrix_data::Bool: boolean keyword to specify matrix or array data
  • c_dim::Int: column number of matrix data
source
FEMTISE.reduced_densityMethod
reduced_density(ϕ,r,model)

Aim

  • Compute partial probability density from set of specific wave functions (FE object)

Arguments

  • phi::Vector{CellField}: array of eigenstates (array of FE objects)
  • r::Tuple{Vector{Float64},Vector{Float64}}:
    • x_vector::Vector{Float64}:: array values of first coordinate (first partial domain)
    • y_vector::Vector{Float64}:: array values of second coordinate (second partial domain)
  • model::CartesianDiscreteModel:: Cartesian discrete 2D model

Output

  • result::Tuple{Vector{Float64},Vector{Float64}}
    • reduced_function_DOF1::Vector{Float64}: partial probability density of first coordinate
    • reduced_function_DOF2::Vector{Float64}: partial probability density of second coordinate
source
FEMTISE.reduced_integrationMethod
reduced_integration(FE_function,r_vector,Ω,dΩ)

Aim

  • Integration over reduced coordinate (partial domain) of specific two dimensional FE function

Arguments

  • FE_function::Vector{CellField}: specific 2D FE function to integrate over partial domain.
  • r_vector::Tuple{Vector{Float64},Vector{Float64}}
    • x_vector::Vector{Float64}:: array values of first coordinate (first partial domain)
    • y_vector::Vector{Float64}:: array values of second coordinate (second partial domain)
  • Omega: Finite element domain
  • dΩ::Gridap.CellData.GenericMeasure: full integration domain

Output

  • result::Tuple{Vector{Float64},Vector{Float64}}
    • reduced_function_DOF1::Vector{Float64}: integration of 2D FE function over first partial domain.
    • reduced_function_DOF2::Vector{Float64}: integration of 2D FE function over second partial domain.
source
FEMTISE.reduced_time_indep_entropyMethod
reduced_time_indep_entropy(r,rho,renyi_factor)

Aim

  • Compute time independent Reduced Rényi entropy for specific 2D FE wave function and Rényi factor

Arguments

  • r::Tuple{Vector{Float64},Vector{Float64}}:
    • x::Vector{Float64}: array values of first coordinate
    • y::Vector{Float64}: array values of second coordinate
  • rho::Matrix{Float64}: matrix values of probability density (where rows are rho(x) with fixed y and columns are rho(y) with fixed x)
  • renyi_factor::Float64: Rényi factor (need to be a Real number)

Output

  • S::Tuple{Vector{Float64},Vector{Float64}}:
    • Sx_vector::Vector{Float64}: reduced entropy as integration over first coordinate.
    • Sy_vector::Vector{Float64}: reduced entropy as integration over second coordinate.
source
FEMTISE.reduced_time_indep_entropyMethod
reduced_time_indep_entropy(r,rho)

Aim

  • Compute time independent Reduced Shannon entropy for specific 2D FE wave function

Arguments

  • r::Tuple{Vector{Float64},Vector{Float64}}:
    • x::Vector{Float64}: array values of first coordinate
    • y::Vector{Float64}: array values of second coordinate
  • rho::Matrix{Float64}: matrix values of probability density (where rows are rho(x) with fixed y and columns are rho(y) with fixed x)

Output

  • S::Tuple{Vector{Float64},Vector{Float64}}:
    • Sx_vector::Vector{Float64}: reduced entropy as integration over first coordinate.
    • Sy_vector::Vector{Float64}: reduced entropy as integration over second coordinate.
source
FEMTISE.solveMethod
solve(prob)

Aim

  • Compute eigen problem by Arpack eigs function and return eigenvalues and eigenvectors.

Arguments

  • prob::EigenProblem: problem deinition

Returns

  • eigenvalues::Vector{Float64}: eigenvalues
  • eigenvectors::Vector{CellField}: eigenvectors
source
FEMTISE.space_coordMethod
space_coord(dom,Δr,n;<keyword arguments>)

Aim

  • Returns coordinate vector (r) and discrete points (pts) for 1D or 2D spaces.
    • if dimension=="1D" ⇒ dom=(x₁,x₂); Δr=Δx; n=nx
    • if dimension=="2d" ⇒ dom=(x₁,x₂,y₁,y₂); Δr=(Δx,Δy); n=(nx,ny)

Arguments

  • dom::Tuple: FE cartesian domain.
  • Δr:: discretization of FE space.
  • n:: number of FE in each direction.
source
FEMTISE.state_populationMethod
state_population(psi,TrialSpace,dOmega)

Aim

  • Compute the population of specific quantum state vector (2D)

Arguments

  • psi::Vector{CellField}: multifield quantum state (2D)
  • TrialSpace::FESpace: Trial finite element space
  • dOmega::Gridap.CellData.GenericMeasure: integration domain

Output

  • P::Tuple{Array{Float64},Array{Float64}}:
    • P1::Array{Float64}: population of first state
    • P2::Array{Float64}: population of second state
source
FEMTISE.time_indep_diff_mutual_informationMethod
time_indep_diff_mutual_information(Sx,Sy,Sxy)

Aim

  • Compute time independent Differential Mutual Information

Arguments

  • Sx::Vector{Float64}: reduced entropy associated with first coordinate
  • Sy::Vector{Float64}: reduced entropy associated with second coordinate
  • Sxy::Vector{Float64}: total entropy associated with both coordinate

Output

  • I::Vector{Float64}: array of mutual information values
source
FEMTISE.time_indep_entropyMethod
time_indep_entropy(psi,TrialSpace,dOmega,renyi_factor)

Aim

  • Compute time independent Rényi entropy for specific FE wave function and Rényi factor

Arguments

  • psi::Vector{CellField}: specific FE wave function
  • TrialSpace::FESpace: Trial FE space
  • dOmega::Gridap.CellData.GenericMeasure: integration domain
  • renyi_factor::Real: Rényi factor (have to be a Real number)

Output

  • S::Vector{Float64}: Array of entropy values.
source
FEMTISE.time_indep_entropyMethod
time_indep_entropy(psi,TrialSpace,dOmega)

Aim

  • Compute time independent Shannon entropy for specific FE wave function

Arguments

  • psi::Vector{CellField}: specific FE wave function
  • TrialSpace::FESpace: Trial FE space
  • dOmega::Gridap.CellData.GenericMeasure: integration domain

Output

  • S::Vector{Float64}: Array of entropy values.
source
FEMTISE.triangulation_repairMethod
triangulation_repair(model,grid_type)

Aim

  • Repair triangulation

Arguments

  • model::Model: FE grid model
  • grid_type::String: FE grid type (could be "simple_line" or "Cartesian2D")
source
FEMTISE.uniform_trapezoidal_integration_methodMethod
uniform_trapezoidal_integration_method(x_vec,fx_vec)

Aim

  • Numeric integration function using trapezoidal method considering uniform discretization.

Arguments

  • x_vec::Vector: integration domain (array of values)
  • fx_vec::Vector: array function evaluated in integration domain

Output

  • result::Vector: integral value

Example

x=0:0.1:10;
fx=sin.(x);
result=uniform_trapezoidal_integration_method(x,fx)
source
FEMTISE.write_binMethod
write_bin(x,fileName;<keyword arguments>)

Aim

  • Write binary file from matrix data

Arguments

  • x::Matrix{T}: matrix data
  • fileName::String: name of file data
  • existing_file::Bool=false: optional boolean keyword to delete or not delete existing data
source
FEMTISE.write_binMethod
write_bin(x,fileName;<keyword arguments>)

Aim

  • Write binary file from array data

Arguments

  • x::Array{T,1}: array/vector data
  • fileName::String: name of file data
  • existing_file::Bool=false: optional boolean keyword to delete or not delete existing data
source