Documentation of the FEMTISE module.
FEMTISE.aprox_dirac_delta
— Methodaprox_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}
: parametersx₀::Float64
: specific point to centre dirac functionδnorm::Float64
: norm value to obtain normalized delta functioncomponent::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)
FEMTISE.bilineal_forms
— Methodbilineal_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
FEMTISE.build_input_data
— Methodbuild_input_data(full_path_name)
Aim
- Build input data from a file
Arguments
full_path_name::String
: full path name to input data file
FEMTISE.collect_result_data
— Methodcollect_result_data(switch_input_file,full_path_name)
Aim
- Collect result data
Arguments
switch_input_file::Bool
: switch input filefull_path_name::String
: full path name
FEMTISE.coord_first_moment
— Methodcoord_first_moment(psi,TrialSpace,Omega,dOmega,x_component)
Aim
- Compute first moment of coordinate (specific component)
Arguments
psi::Vector{CellField}
: array of FE wave functionsTrialSpace::FESpace
: Trial finite element spaceOmega
: Finite element domaindOmega::Gridap.CellData.GenericMeasure
: integration domainx_component::Int
: component of specific coordinate variable
Output
expval::Array{Float64}
: array of expectation values of specific coordinate
FEMTISE.coord_second_moment
— Methodcoord_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 functionsTrialSpace::FESpace
: Trial finite element spaceOmega
: Finite element domaindOmega::Gridap.CellData.GenericMeasure
: integration domainx_component::Int
: component of specific coordinate variable
Output
var::Array{Float64}
: array of variance values of specific coordinate
FEMTISE.coord_standar_deviation
— Methodcoord_standar_deviation(expval,var)
Aim
- Compute standar deviation of coordinate (specific component)
Arguments
expval::Array{Float64}
: array of expectation values of specific coordinatevar::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)
FEMTISE.coord_standar_deviation
— Methodcoord_standar_deviation(psi,TrialSpace,Omega,dOmega,x_component)
Aim
- Compute standar deviation of coordinate (specific component)
Arguments
psi::Vector{CellField}
: array of FE wave functionsTrialSpace::FESpace
: Trial finite element spaceOmega
: Finite element domaindOmega::Gridap.CellData.GenericMeasure
: integration domainx_component::Int
: component of specific coordinate variable
Output
sigma::Array{Float64}
: array of standar deviation of specific coordinate
FEMTISE.default_solver_eigen_problem
— Methoddefault_solver_eigen_problem(params)
Aim
- Function to resolve unidimensonal eigen problem
Arguments
params::Params1D
: parameters of 1D potential
FEMTISE.default_solver_eigen_problem
— Methoddefault_solver_eigen_problem(params)
Aim
- Function to resolve bidimensonal eigen problem over cartesian grid
Arguments
params::Params2D
: parameters fo 2D potentialdifferent_masses::Tuple
: keyword to specify if we want to simulate two particles with diferent masses
FEMTISE.density
— Methoddensity(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)
FEMTISE.eigen_problem
— Methodeigen_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ébilweakform_m::Function
: forma bilineal lado derecho de la formulación débiltest::FESpace
: espacio de prueba, puede ser MultiFieldFESpacetrial::FESpace
: espacio de solución, puede ser MultiFieldFESpacenev::Int=10
: número de autovalores requeridostol::Float64=10e-6
: relative tolerance for convergence of Ritz valuesmaxiter::Integer=100
: maximum number of iterationsexplicittransform::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
FEMTISE.eigen_values_and_eigen_vectors
— Methodeigen_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 formulationdΩ::Gridap.CellData.GenericMeasure
: measure of FE gridUSpace::FESpace
: trial FE SpaceVSpace::FESpace
: test FE Spaceparams::Tuple=(nev,tol,maxiter,explicittransform,sigma)
: params to Arpack eigs function.nev::Integer=10
: quantity of eigendata to calculatetol::Float64=10e-9
: relative tolerance for convergence of Ritz valuesmaxiter::Integer=100
: maximum number of iterationsexplicittransform::Symbol=:none
: shift and invert should be explicitly invoked in julia codesigma::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)
FEMTISE.eigen_values_and_eigen_vectors
— Methodeigen_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 formulationdΩ::Gridap.CellData.GenericMeasure
: measure of FE gridUSpace::FESpace
: trial FE SpaceVSpace::FESpace
: test FE Spaceparams::Tuple=(nev,tol,maxiter,explicittransform,sigma)
: params to Arpack eigs function.nev::Integer=10
: quantity of eigendata to calculatetol::Float64=10e-9
: relative tolerance for convergence of Ritz valuesmaxiter::Integer=100
: maximum number of iterationsexplicittransform::Symbol=:none
: shift and invert should be explicitly invoked in julia codesigma::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)
FEMTISE.eigenstates_normalization
— Methodeigenstates_normalization(ϕ,dΩ)
Aim
- Check eigenstates normalization.
Arguments
ϕ::Vector{CellField}
: array of FE eigenstatesdΩ::Gridap.CellData.GenericMeasure
: integration domain
Output
nom_vector::Vector{Float64}
: array of norm value for specific eigenstate
FEMTISE.fe_spaces
— Methodfe_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,)
FEMTISE.input_data
— Methodinput_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
FEMTISE.interpolation_eigenstates!
— Methodinterpolation_eigenstates!(eigen_states,USpace)
Aim
- Interpolate eigenstates
Arguments
eigen_states::Vector{CellField}
: eigenstatesUSpace::FESpace
: FE space
FEMTISE.load_input_data
— Methodload_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
FEMTISE.load_results
— Methodload_results(id)
Aim
- Load results from a file
Arguments
id::InputData1D
: input data
FEMTISE.load_results
— Methodload_results(id)
Aim
- Load results from a file
Arguments
id::InputData2D
: input data
FEMTISE.load_results
— Methodload_results(id)
Aim
- Load results from a file
Arguments
id::InputData
: input data
FEMTISE.make_boundary_conditions
— Methodmake_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 gridBC_type::String
: type of boundary conditionsTypeData::Type
: type of data<keyword arguments>
:homogeneous::Bool=true
: if boundary conditions are homogeneousparams::Tuple=(nothing,nothing)
: parameters for boundary conditions
Returns
BC_values::Array{TypeData,1}
: values of boundary conditionsBC_tags::Array{String,1}
: tags of boundary conditions
Example
BC_values,BC_tags=make_boundary_conditions("simple_line","FullDirichlet",Float64)
FEMTISE.make_model
— Methodmake_model(grid_type,params)
Aim
- Create a mesh model using Gmsh
Arguments
grid_type::String
: type of grid to createparams::Tuple
: parameters to create the grid
Returns
model::GmshDiscreteModel
: mesh model
References
- https://gitlab.onelab.info/gmsh/gmsh/blob/master/api/gmsh.jl
FEMTISE.measures
— Methodmeasures(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.
FEMTISE.norm_l2
— Methodnorm_l2(𝜳,dΩ)
Aim
- Compute de L2 norm for specific FE wave function
Arguments
𝜳::CellField
: specific FE wave functiondΩ::Gridap.CellData.GenericMeasure
: integration domain
Output
orthogonality_vector::Real
: L2 norm value
FEMTISE.orthogonality_check
— Methodorthogonality_check(ϕ,dΩ,TrialSpace)
Aim
- Check eigenstates orthogonality. Each eigenstate is a multifield FE object.
Arguments
ϕ::Vector{CellField}
: array of FE eigenstatesdΩ::Gridap.CellData.GenericMeasure
: integration domainTrialSpace::FESpace
: Trial FE space
Output
orthogonality_vector::Vector{Float64}
: array of inner product between differents eigenstates
FEMTISE.orthogonality_check
— Methodorthogonality_check(ϕ,dΩ)
Aim
- Check eigenstates orthogonality
Arguments
ϕ::Vector{CellField}
: array of FE eigenstatesdΩ::Gridap.CellData.GenericMeasure
: integration domain
Output
orthogonality_vector::Vector{Float64}
: array of inner product between differents eigenstates
FEMTISE.read_bin
— Methodread_bin(fileName;<keyword arguments>)
Aim
- Read binary file
Arguments
fileName::String
: name of file datamatrix_data::Bool
: optional boolean keyword to specify matrix or array datac_dim::Int
: optional column number of matrix data
FEMTISE.read_bin
— Methodread_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 typen::Int
: total number of elements in array/matrix datamatrix_data::Bool
: boolean keyword to specify matrix or array datac_dim::Int
: column number of matrix data
FEMTISE.reduced_density
— Methodreduced_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 coordinatereduced_function_DOF2::Vector{Float64}
: partial probability density of second coordinate
FEMTISE.reduced_integration
— Methodreduced_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 domaindΩ::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.
FEMTISE.reduced_time_indep_entropy
— Methodreduced_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 coordinatey::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.
FEMTISE.reduced_time_indep_entropy
— Methodreduced_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 coordinatey::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.
FEMTISE.solve
— Methodsolve(prob)
Aim
- Compute eigen problem by Arpack eigs function and return eigenvalues and eigenvectors.
Arguments
prob::EigenProblem
: problem deinition
Returns
eigenvalues::Vector{Float64}
: eigenvalueseigenvectors::Vector{CellField}
: eigenvectors
FEMTISE.space_coord
— Methodspace_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.
FEMTISE.state_population
— Methodstate_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 spacedOmega::Gridap.CellData.GenericMeasure
: integration domain
Output
P::Tuple{Array{Float64},Array{Float64}}
:P1::Array{Float64}
: population of first stateP2::Array{Float64}
: population of second state
FEMTISE.time_indep_diff_mutual_information
— Methodtime_indep_diff_mutual_information(Sx,Sy,Sxy)
Aim
- Compute time independent Differential Mutual Information
Arguments
Sx::Vector{Float64}
: reduced entropy associated with first coordinateSy::Vector{Float64}
: reduced entropy associated with second coordinateSxy::Vector{Float64}
: total entropy associated with both coordinate
Output
I::Vector{Float64}
: array of mutual information values
FEMTISE.time_indep_entropy
— Methodtime_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 functionTrialSpace::FESpace
: Trial FE spacedOmega::Gridap.CellData.GenericMeasure
: integration domainrenyi_factor::Real
: Rényi factor (have to be a Real number)
Output
S::Vector{Float64}
: Array of entropy values.
FEMTISE.time_indep_entropy
— Methodtime_indep_entropy(psi,TrialSpace,dOmega)
Aim
- Compute time independent Shannon entropy for specific FE wave function
Arguments
psi::Vector{CellField}
: specific FE wave functionTrialSpace::FESpace
: Trial FE spacedOmega::Gridap.CellData.GenericMeasure
: integration domain
Output
S::Vector{Float64}
: Array of entropy values.
FEMTISE.triangulation_repair
— Methodtriangulation_repair(model,grid_type)
Aim
- Repair triangulation
Arguments
model::Model
: FE grid modelgrid_type::String
: FE grid type (could be "simple_line" or "Cartesian2D")
FEMTISE.uniform_trapezoidal_integration_method
— Methoduniform_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)
FEMTISE.write_bin
— Methodwrite_bin(x,fileName;<keyword arguments>)
Aim
- Write binary file from matrix data
Arguments
x::Matrix{T}
: matrix datafileName::String
: name of file dataexisting_file::Bool=false
: optional boolean keyword to delete or not delete existing data
FEMTISE.write_bin
— Methodwrite_bin(x,fileName;<keyword arguments>)
Aim
- Write binary file from array data
Arguments
x::Array{T,1}
: array/vector datafileName::String
: name of file dataexisting_file::Bool=false
: optional boolean keyword to delete or not delete existing data