ABINIT/dotproduct [ Functions ]

I've benchmarked this to be speedier than the intrinsic dot_product even on
big vectors. The point is that less check is performed.
MG: FIXME: Well, optized blas1 is for sure better than what you wrote!
Now I dont' have time to update ref files

m_numeric_tools/bisect_int [ Functions ]

Given an array AA(1:N), and a value x, returns the index j such that AA(j)<=x<= AA(j + 1).
AA must be monotonic, either increasing or decreasing. j=0 or
j=N is returned to indicate that x is out of range.

m_numeric_tools/bisect_rdp [ Functions ]

Given an array AA(1:N), and a value x, returns the index j such that AA(j)<=x<= AA(j + 1).
AA must be monotonic, either increasing or decreasing. j=0 or
j=N is returned to indicate that x is out of range.

nlev=Number of "levels" in the continued fraction.
term_type=Type of the terminator.
0 --> No terminator.
-1 --> Assume constant coefficients for a_i and b_i for i>nlev with a_inf = a(nlev) and b_inf = b(nleb)
1 --> Same as above but a_inf and b_inf are obtained by averaging over the nlev values.
aa(nlev)=Set of a_i coefficients.
bb(nlev)=Set of b_i coefficients.
nz=Number of points on the z-mesh.
zpts(nz)=z-mesh.

m_numeric_tools/cspint [ Functions ]

The routine is given the value of a function F(X) at a set of
nodes XTAB, and estimates
Integral ( A <= X <= B ) F(X) DX
by computing the cubic natural spline S(X) that interpolates
F(X) at the nodes, and then computing
Integral ( A <= X <= B ) S(X) DX
exactly.
Other output from the program includes the definite integral
from X(1) to X(I) of S(X), and the coefficients necessary for
the user to evaluate the spline S(X) at any point.
Modified:
30 October 2000
Reference:
Philip Davis and Philip Rabinowitz,
Methods of Numerical Integration,
Blaisdell Publishing, 1967.
Parameters:
Input, real (dp) FTAB(NTAB), contains the tabulated values of
the function, FTAB(I) = F(XTAB(I)).
Input, real (dp) XTAB(NTAB), contains the points at which the
function was evaluated. The XTAB's must be distinct and
in ascending order.
Input, integer NTAB, the number of entries in FTAB and
XTAB. NTAB must be at least 3.
Input, real (dp) A, lower limit of integration.
Input, real (dp) B, upper limit of integration.
Output, real (dp) Y(3,NTAB), will contain the coefficients
of the interpolating natural spline over each subinterval.
For XTAB(I) <= X <= XTAB(I+1),
S(X) = FTAB(I) + Y(1,I)*(X-XTAB(I))
+ Y(2,I)*(X-XTAB(I))**2
+ Y(3,I)*(X-XTAB(I))**3
Output, real (dp) E(NTAB), E(I) = the definite integral from
XTAB(1) to XTAB(I) of S(X).
Workspace, real (dp) WORK(NTAB).
Output, real (dp) RESULT, the estimated value of the integral.

m_numeric_tools/findmin [ Functions ]

Compute the minimum of a function whose value and derivative are known at two points.
Also deduce different quantities at this predicted point, and at the two other points
It uses a quartic interpolation, with the supplementary
condition that the second derivative vanishes at one and
only one point. See Schlegel, J. Comp. Chem. 3, 214 (1982) [[cite:Schlegel1982]].
For this option, lambda_1 must be 1 (new point),
and lambda_2 must be 0 (old point).
Also, if the derivative at the new point is more negative
than the derivative at the old point, the predicted
point cannot correspond to a minimum, but will be lambda=2.5_dp,
if the energy of the second point is lower than the energy
of the first point.

INPUTS

etotal_1=first value of the function
etotal_2=second value of the function
dedv_1=first value of the derivative
dedv_2=second value of the derivative
lambda_1=first value of the argument
lambda_2=second value of the argument

OUTPUT

dedv_predict=predicted value of the derivative (usually zero,
except if choice=4, if it happens that a minimum cannot be located,
and a trial step is taken)
d2edv2_predict=predicted value of the second derivative (not if choice=4)
d2edv2_1=first value of the second derivative (not if choice=4)
d2edv2_2=second value of the second derivative (not if choice=4)
etotal_predict=predicted value of the function
lambda_predict=predicted value of the argument
status= 0 if everything went normally ;
1 if negative second derivative
2 if some other problem

chmout(n*n+n)=redefined matrix with strictly real diagonal elements.
May be same storage location as chmin.
ierr=0 if no problem, 1 if the imaginary part of some element
too large (at present, stop in this case).

m_numeric_tools/hermitianize_dpc [ Functions ]

uplo=String describing which part of the matrix has been calculated.
Only the first character is tested (no case sensitive). Possible values are:
"All"= Full matrix is supplied in input
"Upper"=Upper triangle is in input. Lower triangle is reconstructed by symmetry.
"Lower"=Lower triangle is in input. Upper triangle is reconstructed by symmetry.

m_numeric_tools/hermitianize_spc [ Functions ]

uplo=String describing which part of the matrix has been calculated.
Only the first character is tested (no case sensitive). Possible values are:
"All"= Full matrix is supplied in input
"Upper"=Upper triangle is in input. Lower triangle is reconstructed by symmetry.
"Lower"=Lower triangle is in input. Upper triangle is reconstructed by symmetry.

nomega=number of real frequencies
omega(nomega)= real frequencies
eps(nomega)= function on the frequency grid (both real and imaginary part)
real part can be used to check whether the K-K relation is satisfied or not
method=method used to perform the integration
0= naive integration
1=simpson rule
only_check= if /=0 the real part of eps is checked against the imaginary part,
a final report in written but the array eps is not modified
if ==0 the real part of eps is overwritten using the
results obtained using the Kramers-Kronig relation

m_numeric_tools/list2blocks [ Functions ]

Given a list of integers, find the number of contiguos groups of values.
and returns the set of indices that can be used to loop over these groups
Example list = [1,2,3,5,6] --> blocks = [[1,3], [4,5]]

INPUTS

list(:)=List of integers
OUTPUTS
nblocks=Number of blocks
blocks(2,nblocks)=
allocatable array in input
in output:
blocks(1,i) gives the start of the i-th block
blocks(2,i) gives the end of the i-th block

m_numeric_tools/llsfit_svd [ Functions ]

Given a set of N data points (x,y) with individual standard deviations sigma_i,
use chi-square minimization to determine the M coefficients, par, of a function that
depends linearly on nfuncs functions, i.e f(x) = \sum_i^{nfuncs} par_i * func_i(x).
Solve the fitting equations using singular value decomposition of the design matrix as in Eq 14.3.17
of Numerical Recipes. The program returns values for the M fit parameters par, and chi-square.
The user supplies a subroutine funcs(x,nfuncs) that returns the M basis functions evaluated at xx.

m_numeric_tools/midpoint_ [ Functions ]

This routine computes the n-th stage of refinement of an extended midpoint rule.

INPUTS

func(external)=the name of the function to be integrated
xmin,xmax=the limits of integration
nn=integer defining the refinement of the mesh, each call adds (2/3)*3n-1 additional
interior points between xmin ans xmax

OUTPUT

See SIDE EFFECTS

SIDE EFFECTS

quad=the integral at the n-th stage.

PARENTS

CHILDREN

NOTES

When called with nn=1, the routine returns as quad the crudest estimate of the integral
Subsequent calls with nn=2,3,... (in that sequential order) will improve the accuracy of quad by adding
(2/3)*3n-1 additional interior points. quad should not be modified between sequential calls.
Subroutine is defined as recursive to allow multi-dimensional integrations

m_numeric_tools/polyn_interp [ Functions ]

Given arrays xa and ya of length N, and given a value x, return a value y, and an error estimate dy.
If P(x) is the polynomial of degree N-1 such that P(xai)=yai, i=1,...,N, then the returned value y=P(x).

INPUTS

xa(:)=abscissas in ascending order
ya(:)=ordinates
x=the point where the set of data has to be interpolated

m_numeric_tools/print_arr1d_spc [ Functions ]

arr(:)=vector/matrix to be printed
mode_paral(optional)=parallel mode, DEFAULT is "COLL"
"COLL" if all procs are calling the routine with the same message to be written only once
"PERS" if the procs are calling the routine with different mesgs each to be written,
or if one proc is calling the routine
unit(optional)=the unit number of the file, DEFAULT=std_out
max_r,max_c(optional)=Max number of rows and columns to be printed
(DEFAULT is 9, output format assumes to be less that 99, but there might be
problems with wrtout if message size exceeds 500 thus max number of elements should be ~60)

m_numeric_tools/quadrature [ Functions ]

Driver routine to perform quadratures in finite domains using different techniques.
The routine improves the resolution of the grid until a given accuracy is reached

INPUTS

func(external)=the name of the function to be integrated
xmin,xmax=the limits of integration
npts=Initial number of points, only for Gauss-Legendre. At each step this number is doubled
accuracy=fractional accuracy required
ntrial=Max number of attempts
qopt=integer flag defining the algorithm for the quadrature:
1 for Trapezoidal rule, closed, O(1/N^2)
2 for Simpson based on trapezoidal,closed, O(1/N^4)
3 for Midpoint rule, open, O(1/N^2)
4 for midpoint rule with cancellation of leading error, open, O(1/N^4)
5 for Romberg integration (closed form) and extrapolation for h-->0 (order 10 is hard-coded)
6 for Romberg integration with midpoint rule and extrapolation for h-->0 (order 10 is hard-coded)
7 for Gauss-Legendre

m_numeric_tools/remove_copies [ Functions ]

Given an initial set of elements, set_in, return the subset of inequivalent items
packed in the firt n_out positions of set_in. Use the logical function is_equal
to define whether two items are equivalent.

INPUTS

n_in=Initial number of elements.
is_equal=logical function used to discern if two items are equal.

OUTPUT

n_out=Number of inequivalent items found.

SIDE EFFECTS

set_in(3,n_in)=
In input the initial set of n_in elements
In output set_in(3,1:n_out) contains the inequivalent elements found.

NOTES

The routines only deals with arrays of 3D-vectors although generalizing the
algorithm to nD-space is straightforward.

m_numeric_tools/simpson_cplx [ Functions ]

npts=Number of points.
step=Step of the mesh.
ff(npts)=Values of the integrand.

OUTPUT

simpson_cplx=Integral of ff.

NOTES

If npts is odd, the integration is done with the extended Simpson's rule (error = O^(step^4))
If npts is even, the last 4 four points are integrated separately via Simpson's 3/8 rule. Error = O(step^5)
while the first npts-3 points are integrared with the extended Simpson's rule.

m_numeric_tools/symmetrize_dpc [ Functions ]

uplo=String describing which part of the matrix has been calculated.
Only the first character is tested (no case sensitive). Possible values are:
"All"= Full matrix is supplied in input
"Upper"=Upper triangle is in input. Lower triangle is reconstructed by symmetry.
"Lower"=Lower triangle is in input. Upper triangle is reconstructed by symmetry.

m_numeric_tools/symmetrize_spc [ Functions ]

uplo=String describing which part of the matrix has been calculated.
Only the first character is tested (no case sensitive). Possible values are:
"All"= Full matrix is supplied in input
"Upper"=Upper triangle is in input. Lower triangle is reconstructed by symmetry.
"Lower"=Lower triangle is in input. Upper triangle is reconstructed by symmetry.

m_numeric_tools/trapezoidal_ [ Functions ]

Compute the n-th stage of refinement of an extended trapezoidal rule
adding 2^(n-2) additional interior point in the finite range of integration

INPUTS

func(external)=the name of the function to be integrated
xmin,xmax=the limits of integration
nn=integer defining the refinement of the mesh, each call adds 2^(n-2) additional interior points

OUTPUT

See SIDE EFFECTS

SIDE EFFECTS

quad=the integral at the n-th stage.

PARENTS

CHILDREN

NOTES

When called with nn=1, the routine returns the crudest estimate of the integral
Subsequent calls with nn=2,3,... (in that sequential order) will improve the accuracy
by adding 2^(n-2) additional interior points. Note that quad should not be modified between sequential calls.
Subroutine is defined as recursive to allow multi-dimensional integrations

m_numeric_tools/uniformrandom [ Functions ]

Returns a uniform random deviate between 0.0 and 1.0.
Set seed to any value < 0 to initialize or reinitialize sequence.
Parameters are chosen from integer overflow=2**23 (conservative).
For some documentation, see Numerical Recipes, 1986, p196.