Psi4
Files | Functions
libqt: The Quantum-Trio Miscellaneous Library

Files

file  3d_array.cc
 Routines for 3d arrays.
 
file  blas_intfc.cc
 The PSI3 BLAS1 interface routines.
 
file  blas_intfc23.cc
 Interface to all BLAS routinesAutogenerated by Rob Parrish on 1/24/2011.
 
file  blas_intfc23_mangle.h
 The PSI3 BLAS2 and BLAS3 interface routines.
 
file  blas_intfc_mangle.h
 The PSI3 BLAS1 interface routines.
 
file  cc_excited.cc
 Determine if a wavefunction is a CC-excited wavefunction type.
 
file  david.cc
 In-core Davidson-Liu diagonalization of symm matrices.
 
file  dirprd_block.cc
 Take a direct product of two matrices.
 
file  dot_block.cc
 Take dot product of two block matrices.
 
file  fill_sym_matrix.cc
 Fill a symmetric matrix from a lower triangle.
 
file  invert.cc
 Invert a small matrix.
 
file  lapack_intfc.cc
 Interface to all LAPACK routinesAutogenerated by Rob Parrish on 1/23/2011.
 
file  lapack_intfc.h
 Interface to LAPACK routinesRollin A. King and T. Daniel Crawford August 2001 - January 2002.
 
file  lapack_intfc_mangle.h
 The PSI3 LAPACK interface routines.
 
file  mat_print.cc
 Print a matrix to a file in a formatted style.
 
file  newmm_rking.cc
 Matrix multiplication routine (deprecated)
 
file  normalize.cc
 Normalize a set of vectors.
 
file  pople.cc
 Pople's method for solving linear equations.
 
file  probabil.cc
 Contains some probability functions.
 
file  qt.h
 Header file for the Quantum Trio LibraryDavid Sherrill 1994.
 
file  ras_set.cc
 Obtain orbital space and reordering for CI/MCSCF wavefunctions.
 
file  reorder_qt.cc
 Obtain the QT orbital reordering array between Pitzer and correlated order.
 
file  schmidt.cc
 Gram-Schmidt orthogonalize a set of vectors.
 
file  libqt/schmidt_add.cc
 Gram-Schmidt orthogonalize vector and add to set.
 
file  slaterdset.h
 Header file for SlaterDetSetsEdward Valeev, June 2002.
 
file  solve_pep.cc
 Solve a 2x2 pseudo-eigenvalue problem.
 
file  timer.cc
 Obtain user and system timings for blocks of codeTIMER.CC: These functions allow one to obtain user and system timings for arbitrary blocks of code. If a code block is called repeatedly during the course of program execution, the timer functions will report the block's cumulative execution time and the number of calls. In addition, one may time multiple code blocks simultaneously, and even ``overlap'' timers. Timing data is written to the file "timer.dat" at the end of timer execution, i.e., when timer_done() is called.
 

Functions

double *** psi::init_3d_array (int p, int q, int r)
 
void psi::free_3d_array (double ***A, int p, int q)
 
void psi::C_DSWAP (size_t length, double *x, int inc_x, double *y, int inc_y)
 
void psi::C_DAXPY (size_t length, double a, double *x, int inc_x, double *y, int inc_y)
 
void psi::C_DCOPY (size_t length, double *x, int inc_x, double *y, int inc_y)
 
void psi::C_DSCAL (size_t length, double alpha, double *vec, int inc)
 
void psi::C_DROT (size_t length, double *x, int inc_x, double *y, int inc_y, double costheta, double sintheta)
 
double psi::C_DDOT (size_t length, double *x, int inc_x, double *y, int inc_y)
 
double psi::C_DNRM2 (size_t length, double *x, int inc_x)
 
double psi::C_DASUM (size_t length, double *x, int inc_x)
 
size_t psi::C_IDAMAX (size_t length, double *x, int inc_x)
 
int psi::cc_excited (const char *wfn)
 
int psi::cc_excited (std::string wfn)
 
int psi::david (double **A, int N, int M, double *eps, double **v, double cutoff, int print)
 
void psi::dirprd_block (double **A, double **B, int rows, int cols)
 
double psi::dot_block (double **A, double **B, int rows, int cols, double alpha)
 
void psi::fill_sym_matrix (double **A, int size)
 
double psi::invert_matrix (double **a, double **y, int N, std::string out)
 
int psi::C_DGEEV (int n, double **a, int lda, double *wr, double *wi, double **vl, int ldvl, double **vr, int ldvr, double *work, int lwork, int info)
 
int psi::C_DGESV (int n, int nrhs, double *a, int lda, int *ipiv, double *b, int ldb)
 
int psi::C_DGETRF (int m, int n, double *a, int lda, int *ipiv)
 
int psi::C_DPOTRF (char uplo, int n, double *a, int lda)
 
int psi::C_DGETRI (int n, double *a, int lda, int *ipiv, double *work, int lwork)
 
int psi::C_DPOTRI (char uplo, int n, double *a, int lda)
 
int psi::C_DPOTRS (char uplo, int n, int nrhs, double *a, int lda, double *b, int ldb)
 
int psi::C_DGESVD (char jobu, char jobvt, int m, int n, double *A, int lda, double *s, double *u, int ldu, double *vt, int ldvt, double *work, int lwork)
 
int psi::C_DSYEV (char jobz, char uplo, int n, double *a, int lda, double *w, double *work, int lwork)
 
int psi::mat_print (double **matrix, int rows, int cols, std::string out)
 
void psi::normalize (double **A, int rows, int cols)
 
int psi::pople (double **A, double *x, int dimen, int, double tolerance, std::string out, int print_lvl)
 
double psi::combinations (int n, int k)
 
double psi::factorial (int n)
 
int psi::ras_set3 (int nirreps, int nmo, int *orbspi, int *docc, int *socc, int *frdocc, int *fruocc, int *restrdocc, int *restruocc, int **ras_opi, int *core_guess, int *order, int ras_type, bool is_mcscf, Options &options)
 
void psi::reorder_qt (int *docc_in, int *socc_in, int *frozen_docc_in, int *frozen_uocc_in, int *order, int *orbs_per_irrep, int nirreps)
 
void psi::reorder_qt_uhf (int *docc, int *socc, int *frozen_docc, int *frozen_uocc, int *order_alpha, int *order_beta, int *orbspi, int nirreps)
 
void psi::schmidt (double **A, int rows, int cols, std::string out_fname)
 
int psi::schmidt_add (double **A, int rows, int cols, double *v)
 
void psi::solve_2x2_pep (double **H, double S, double *evals, double **evecs)
 
void psi::timer_init (void)
 
void psi::timer_done (void)
 
void psi::timer_on (const std::string &key)
 
void psi::timer_off (const std::string &key)
 
void psi::parallel_timer_on (const std::string &key, int thread_rank)
 
void psi::parallel_timer_off (const std::string &key, int thread_rank)
 

Detailed Description

Function Documentation

double psi::C_DASUM ( size_t  length,
double *  x,
int  inc_x 
)

This function returns the sum of the absolute value of this vector.

Parameters
lengthNumber of elements in x.
xA pointer to the beginning of the data in x. Must be of at least length (1+(N-1)*abs(inc_x).
inc_xhow many places to skip to get to next element in x
Returns
the sum of the absolute value
void psi::C_DAXPY ( size_t  length,
double  a,
double *  x,
int  inc_x,
double *  y,
int  inc_y 
)

This function performs y = a * x + y.

Steps every inc_x in x and every inc_y in y (normally both 1).

Parameters
lengthlength of arrays
ascalar a to multiply vector x
xvector x
inc_xhow many places to skip to get to next element in x
yvector y
inc_yhow many places to skip to get to next element in y
void psi::C_DCOPY ( size_t  length,
double *  x,
int  inc_x,
double *  y,
int  inc_y 
)

This function copies x into y.

Steps every inc_x in x and every inc_y in y (normally both 1).

Parameters
length= length of array
x= vector x
inc_x= how many places to skip to get to next element in x
y= vector y
inc_y= how many places to skip to get to next element in y
double psi::C_DDOT ( size_t  length,
double *  x,
int  inc_x,
double *  y,
int  inc_y 
)

This function returns the dot product of two vectors, x and y.

Parameters
lengthNumber of elements in x and y.
xA pointer to the beginning of the data in x. Must be of at least length (1+(N-1)*abs(inc_x).
inc_xhow many places to skip to get to next element in x
yA pointer to the beginning of the data in y.
inc_yhow many places to skip to get to next element in y
Returns
the dot product
int psi::C_DGEEV ( int  n,
double **  a,
int  lda,
double *  wr,
double *  wi,
double **  vl,
int  ldvl,
double **  vr,
int  ldvr,
double *  work,
int  lwork,
int  info 
)

C_DGEEV()

This function computes the eigenvalues and the left and right right eigenvectors of a real, nonsymmetric matrix A. For symmetric matrices, refer to C_DSYEV().

Parameters
n= The order of the matrix A. n >= 0.
a= The n-by-n matrix A. As with all other lapack routines, must be allocated by contiguous memory (i.e., block_matrix()).
lda= The leading dimension of matrix A. lda >= max(1,n).
wr= array of length n containing the real parts of computed eigenvalues.
wi= array of length n containing the imaginary parts of computed eigenvalues.
vl= matrix of dimensions ldvl*n. The columns store the left eigenvectors u(j). If the j-th eigenvalues is real, then u(j) = vl(:,j), the j-th column of vl. If the j-th and (j+1)-st eigenvalues form a complex conjugate pair, then u(j) = vl(:,j) + i*vl(:,j+1) and u(j+1) = vl(:,j) - i*vl(:,j+1). Note: this is the Fortran documentation, may need to change cols <-> rows.
ldvl= The leading dimension of matrix vl.
vr= matrix of dimensions ldvr*n. The columns store the right eigenvectors v(j). If the j-th eigenvalues is real, then v(j) = vr(:,j), the j-th column of vr. If the j-th and (j+1)-st eigenvalues form a complex conjugate pair, then v(j) = vr(:,j) + i*vr(:,j+1) and v(j+1) = vr(:,j) - i*vr(:,j+1). Note: this is the Fortran documentation, may need to change cols <-> rows.
ldvr= The leading dimension of matrix vr.
work= Array for scratch computations, of dimension lwork. On successful exit, work[0] returns the optimal value of lwork.
lwork= The dimension of the array work. lwork >= max(1,3*n), and if eigenvectors are required (default for this wrapper at present) then actually lwork >= 4*n. For good performance, lwork must generally be larger. If lwork = -1, then a workspace query is assumed. The routine only calculates the optimal size of the work array, returns this value ans the first entry of the work array, and no error message related to lword is issued by xerbla.
info= On output (returned by C_DGEEV), a status flag. info = 0 for successful exit, if info = -i, the ith argument had an illegal value. If info = i, the QR algorithm failed to compute all the eigenvalues, and no eigenvectors have been computed. Elements i+1:n of wr and wi contain eigenvalues which have converged.
int psi::C_DGESV ( int  n,
int  nrhs,
double *  a,
int  lda,
int *  ipiv,
double *  b,
int  ldb 
)

Purpose

DGESV computes the solution to a real system of linear equations A * X = B, where A is an N-by-N matrix and X and B are N-by-NRHS matrices.

The LU decomposition with partial pivoting and row interchanges is used to factor A as A = P * L * U, where P is a permutation matrix, L is unit lower triangular, and U is upper triangular. The factored form of A is then used to solve the system of equations A * X = B.

Arguments

N (input) INTEGER The number of linear equations, i.e., the order of the matrix A. N >= 0.

NRHS (input) INTEGER The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0.

A (input/output) DOUBLE PRECISION array, dimension (LDA,N) On entry, the N-by-N coefficient matrix A. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored.

LDA (input) INTEGER The leading dimension of the array A. LDA >= max(1,N).

IPIV (output) INTEGER array, dimension (N) The pivot indices that define the permutation matrix P; row i of the matrix was interchanged with row IPIV(i).

B (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS) On entry, the N-by-NRHS matrix of right hand side matrix B. On exit, if INFO = 0, the N-by-NRHS solution matrix X.

LDB (input) INTEGER The leading dimension of the array B. LDB >= max(1,N).

C++ Return value: INFO (output) INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, so the solution could not be computed.


.. External Subroutines ..

C_DGESV()

This function solves a system of linear equations A * X = B, where A is an n x n matrix and X and B are n x nrhs matrices.

Parameters
n= The number of linear equations, i.e., the order of the matrix A. n >= 0.
nrhs= The number of right hand sides, i.e., the number of columns of the matrix B. nrhs >= 0.
A= On entry, the n-by-n coefficient matrix A. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored.

lda = The leading dimension of the array A. lda >= max(1,n).

Parameters
intipiv = An integer array of length n. The pivot indices that define the permutation matrix P; row i of the matrix was interchanged with row ipiv(i).
B= On entry, the n-by-nrhs matrix of right hand side matrix B. On exit, if info = 0, the n-by-nrhs solution matrix X.
ldb= The leading dimension of the array B. ldb >= max(1,n).

Returns: int info. info = 0: successful exit. info < 0: if info = -i, the i-th argument had an illegal value. info > 0: if info = i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, so the solution could not be computed.

int psi::C_DGESVD ( char  jobu,
char  jobvt,
int  m,
int  n,
double *  A,
int  lda,
double *  s,
double *  u,
int  ldu,
double *  vt,
int  ldvt,
double *  work,
int  lwork 
)

C_DGESVD() This function computes the singular value decomposition (SVD) of a real mxn matrix A, ** optionally computing the left and/or right singular vectors. The SVD is written

A = U * S * transpose(V)

where S is an mxn matrix which is zero except for its min(m,n) diagonal elements, U is an M-by-M orthogonal matrix, and V is an N-by-N orthogonal matrix. The diagonal elements of SIGMA are the singular values of A; they are real and non-negative, and are returned in descending order. The first min(m,n) columns of U and V are the left and right singular vectors of A.

Note that the routine returns V^t, not V;

These arguments mimic their Fortran counterparts. See the LAPACK manual for additional information.

Parameters
jobu= 'A' = all m columns of U are returned 'S' = the first min(m,n) columns of U are returned 'O' = the first min(m,n) columns of U are returned in the input matrix A 'N' = no columns of U are returned
jobvt= 'A' = all n rows of VT are returned 'S' = the first min(m,n) rows of VT are returned 'O' = the first min(m,n) rows of VT are returned in the input matrix A 'N' = no rows of VT are returned

Obviously jobu and jobvt cannot both be 'O'.

Parameters
m= The row dimension of the matrix A.
n= The column dimension of the matrix A.
A= On entry, the mxn matrix A with dimensions m by lda. On exit, if jobu='O' the first min(m,n) columns of A are overwritten with the left singular vectors; if jobvt='O' the first min(m,n) rows of A are overwritten with the right singular vectors; otherwise, the contents of A are destroyed.
lda= The number of columns allocated for A.
s= The singular values of A.
u= The right singular vectors of A, returned as column of the matrix u if jobu='A'. If jobu='N' or 'O', u is not referenced.
ldu= The number of columns allocated for u.
vt= The left singular vectors of A, returned as rows of the matrix VT is jobvt='A'. If jobvt='N' or 'O', vt is not referenced.
ldvt= The number of columns allocated for vt.
work= Workspace array of length lwork.
lwork= The length of the workspace array work, which should be at least as large as max(3*min(m,n)+max(m,n),5*min(m,n)). For good performance, lwork should generally be larger. If lwork=-1, a workspace query is assumed, and the value of work[0] upon return is the optimal size of the work array.

Return value: int info. If info=0, successful exit. If info = -i, the i-th argument had an illegal value. If info > 0, Related to failure in the convergence of the upper bidiagonal matrix B. See the LAPACK manual for additiona information.

Interface written by TDC, July 2001, updated April 2004

int psi::C_DGETRF ( int  m,
int  n,
double *  a,
int  lda,
int *  ipiv 
)

Purpose

DGETRF computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges.

The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if m > n), and U is upper triangular (upper trapezoidal if m < n).

This is the right-looking Level 3 BLAS version of the algorithm.

Arguments

M (input) INTEGER The number of rows of the matrix A. M >= 0.

N (input) INTEGER The number of columns of the matrix A. N >= 0.

A (input/output) DOUBLE PRECISION array, dimension (LDA,N) On entry, the M-by-N matrix to be factored. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored.

LDA (input) INTEGER The leading dimension of the array A. LDA >= max(1,M).

IPIV (output) INTEGER array, dimension (min(M,N)) The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i).

C++ Return value: INFO (output) INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, and division by zero will occur if it is used to solve a system of equations.


.. Parameters ..

C_DGETRF(): Compute an LU factorization of a general M-by-N matrix a using partial pivoting with row interchanges.

The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if nrow > ncol), and U is upper triangular (upper trapezoidal if nrow < ncol).

Parameters
nrow= number of rows
ncol= number of colums
a= matrix to factorize
lda= leading dimension of a, lda >= max(1,ncol)
ipiv= output integer array of the pivot indices; for 1 <= i <= min(nrow, ncol), col i of the matrix was interchanged with col ipiv(i).

Returns: int info. info = 0: successful exit. info < 0: if info=-i, the ith-argument had an illegal vale. info>0: if info=i, U(i,i) is exactly zero. The factorization has been completed, but the factor U is exactly singular, and division by zero will occur if it is used to solve a system of equations.

int psi::C_DGETRI ( int  n,
double *  a,
int  lda,
int *  ipiv,
double *  work,
int  lwork 
)

Purpose

DGETRI computes the inverse of a matrix using the LU factorization computed by DGETRF.

This method inverts U and then computes inv(A) by solving the system inv(A)*L = inv(U) for inv(A).

Arguments

N (input) INTEGER The order of the matrix A. N >= 0.

A (input/output) DOUBLE PRECISION array, dimension (LDA,N) On entry, the factors L and U from the factorization A = P*L*U as computed by DGETRF. On exit, if INFO = 0, the inverse of the original matrix A.

LDA (input) INTEGER The leading dimension of the array A. LDA >= max(1,N).

IPIV (input) INTEGER array, dimension (N) The pivot indices from DGETRF; for 1<=i<=N, row i of the matrix was interchanged with row IPIV(i).

WORK (workspace/output) DOUBLE PRECISION array, dimension (MAX(1,LWORK)) On exit, if INFO=0, then WORK(1) returns the optimal LWORK.

LWORK (input) INTEGER The dimension of the array WORK. LWORK >= max(1,N). For optimal performance LWORK >= N*NB, where NB is the optimal blocksize returned by ILAENV.

If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.

C++ Return value: INFO (output) INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, U(i,i) is exactly zero; the matrix is singular and its inverse could not be computed.


.. Parameters ..

C_DGETRI(): computes the inverse of a matrix using the LU factorization computed by DGETRF

This method inverts U and then computes inv(A) by solving the system inv(A)*L = inv(U) for inv(A).

Parameters
n= order of the matrix a. n >= 0.
a= array of dimension (n,lda). On entry, the factors L and U from the factorization A = P*L*U as computed by DGETRF. On exit, if info=0, the inverse of hte original matrix A.
lda= the leading dimension of the array a. lda >= max(1,n).
ipiv= The pivot indices from DGETRF. For 1<=i<=n, row i of the matrix was interchanged with row ipiv(i)
work= workspace of dimension max(1,lwork). On exit, if info=0, then work(0) returns the optimal lwork.
lwork= the dimension of the array work. lwork >= max(1,n). For optimal performance, lwork > n*nb, where nb is the optimal blocksize returned by ilaenv. If lwork=-1, then a workspace query is assumed; the routine only calculates the optimal size of the work array, returns this value as the first entry of the work array, and no error message related to lwork is issued by xerbla.

Returns: int info. info=0: successful exit. info<0: if info=-i, the i-th argument had an illegal value. info>0: if info=i, U(i,i) is exactly zero; the matrix is singular and its inverse could not be computed.

double psi::C_DNRM2 ( size_t  length,
double *  x,
int  inc_x 
)

This function returns the square of the norm of this vector.

Parameters
lengthNumber of elements in x.
xA pointer to the beginning of the data in x. Must be of at least length (1+(N-1)*abs(inc_x).
inc_xhow many places to skip to get to next element in x
Returns
the norm squared product
int psi::C_DPOTRF ( char  uplo,
int  n,
double *  a,
int  lda 
)

Purpose

DPOTRF computes the Cholesky factorization of a real symmetric positive definite matrix A.

The factorization has the form A = U**T * U, if UPLO = 'U', or A = L * L**T, if UPLO = 'L', where U is an upper triangular matrix and L is lower triangular.

This is the block version of the algorithm, calling Level 3 BLAS.

Arguments

UPLO (input) CHARACTER*1 = 'U': Upper triangle of A is stored; = 'L': Lower triangle of A is stored.

N (input) INTEGER The order of the matrix A. N >= 0.

A (input/output) DOUBLE PRECISION array, dimension (LDA,N) On entry, the symmetric matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced.

On exit, if INFO = 0, the factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T.

LDA (input) INTEGER The leading dimension of the array A. LDA >= max(1,N).

C++ Return value: INFO (output) INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, the leading minor of order i is not positive definite, and the factorization could not be completed.


.. Parameters ..

C_DPOTRF(): Compute an Cholesky factorization of a Hermitian positive definite N-by-N matrix A.

The factorization has the form A = U**T* U if uplo = 'U' or A = L*L**T if oplo = 'L'

Parameters
uplo,:uplo = 'U' upper triangle of A is stored. uplo = 'L' lower triangle of A is stored
n= number of rows and columns
A= matrix to factorize
lda= leading dimension of a, lda >= max(1,ncol)

Returns: int info. info = 0: successful exit. info < 0: if info=-i, the ith-argument had an illegal vale. info>0: if info=i, the i-th minor is not positive definite

int psi::C_DPOTRI ( char  uplo,
int  n,
double *  a,
int  lda 
)

Purpose

DPOTRI computes the inverse of a real symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF.

Arguments

UPLO (input) CHARACTER*1 = 'U': Upper triangle of A is stored; = 'L': Lower triangle of A is stored.

N (input) INTEGER The order of the matrix A. N >= 0.

A (input/output) DOUBLE PRECISION array, dimension (LDA,N) On entry, the triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, as computed by DPOTRF. On exit, the upper or lower triangle of the (symmetric) inverse of A, overwriting the input factor U or L.

LDA (input) INTEGER The leading dimension of the array A. LDA >= max(1,N).

C++ Return value: INFO (output) INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, the (i,i) element of the factor U or L is zero, and the inverse could not be computed.


.. External Functions ..

C_DPOTRI(): Compute inverse of a Hermitian positive deifinite matrix A AFTER a Cholesky decomposition is computed by C_DPORF.

Parameters
uplo,:uplo = 'U' upper triangle of A is stored. uplo = 'L' lower triangle of A is stored
n= number of rows and columns
A= matrix to factorize
lda= leading dimension of a, lda >= max(1,ncol)

Returns: int info. info = 0: successful exit. info < 0: if info=-i, the ith-argument had an illegal vale. info>0: if info=i, the i-th minor is no positive definite

int psi::C_DPOTRS ( char  uplo,
int  n,
int  nrhs,
double *  a,
int  lda,
double *  b,
int  ldb 
)

Purpose

DPOTRS solves a system of linear equations A*X = B with a symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF.

Arguments

UPLO (input) CHARACTER*1 = 'U': Upper triangle of A is stored; = 'L': Lower triangle of A is stored.

N (input) INTEGER The order of the matrix A. N >= 0.

NRHS (input) INTEGER The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0.

A (input) DOUBLE PRECISION array, dimension (LDA,N) The triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, as computed by DPOTRF.

LDA (input) INTEGER The leading dimension of the array A. LDA >= max(1,N).

B (input/output) DOUBLE PRECISION array, dimension (LDB,NRHS) On entry, the right hand side matrix B. On exit, the solution matrix X.

LDB (input) INTEGER The leading dimension of the array B. LDB >= max(1,N).

C++ Return value: INFO (output) INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value


.. Parameters ..

C_DPOTRS(): Solve a Hermitian positive definite matrix Ax = B AFTER a Cholesky decomposition is computed for A by C_DPORF.

Parameters
uplo,:uplo = 'U' upper triangle of A is stored. uplo = 'L' lower triangle of A is stored
n= number of rows and columns of A
nrhs= number of vectors to solve in B
A= matrix to solve
lda= leading dimension of a, lda >= max(1,ncol)
B= (input)forcing vector/(output)solution vector
ldb= leading dimension of b, ldb >= max(1,ncol)

Returns: int info. info = 0: successful exit. info < 0: if info=-i, the ith-argument had an illegal value.

void psi::C_DROT ( size_t  length,
double *  x,
int  inc_x,
double *  y,
int  inc_y,
double  costheta,
double  sintheta 
)

Calculates a plane Givens rotation for vectors x, y and angle theta. x = x*cos + y*sin, y = -x*sin + y*cos.

Parameters
xvector x
yvector Y
lengthlength of x,y
inc_xhow many places to skip to get to the next element of x
inc_yhow many places to skip to get to the next element of y
void psi::C_DSCAL ( size_t  length,
double  alpha,
double *  vec,
int  inc 
)

This function scales a vector by a real scalar.

Parameters
lengthlength of array
alphascale factor
vecvector to scale
inchow many places to skip to get to next element in vec
void psi::C_DSWAP ( size_t  length,
double *  x,
int  inc_x,
double *  y,
int  inc_y 
)

Swaps a vector with another vector.

Parameters
lengthSpecifies the number of elements in vectors x and y.
xArray, DIMENSION at least (1 + (n-1)*abs(incx)).
inc_xSpecifies the increment for the elements of x.
yArray, DIMENSION at least (1 + (n-1)*abs(incy)).
inc_ySpecifies the increment for the elements of y.
int psi::C_DSYEV ( char  jobz,
char  uplo,
int  n,
double *  a,
int  lda,
double *  w,
double *  work,
int  lwork 
)

Purpose

DSYEV computes all eigenvalues and, optionally, eigenvectors of a real symmetric matrix A.

Arguments

JOBZ (input) CHARACTER*1 = 'N': Compute eigenvalues only; = 'V': Compute eigenvalues and eigenvectors.

UPLO (input) CHARACTER*1 = 'U': Upper triangle of A is stored; = 'L': Lower triangle of A is stored.

N (input) INTEGER The order of the matrix A. N >= 0.

A (input/output) DOUBLE PRECISION array, dimension (LDA, N) On entry, the symmetric matrix A. If UPLO = 'U', the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A. If UPLO = 'L', the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A. On exit, if JOBZ = 'V', then if INFO = 0, A contains the orthonormal eigenvectors of the matrix A. If JOBZ = 'N', then on exit the lower triangle (if UPLO='L') or the upper triangle (if UPLO='U') of A, including the diagonal, is destroyed.

LDA (input) INTEGER The leading dimension of the array A. LDA >= max(1,N).

W (output) DOUBLE PRECISION array, dimension (N) If INFO = 0, the eigenvalues in ascending order.

WORK (workspace/output) DOUBLE PRECISION array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.

LWORK (input) INTEGER The length of the array WORK. LWORK >= max(1,3*N-1). For optimal efficiency, LWORK >= (NB+2)*N, where NB is the blocksize for DSYTRD returned by ILAENV.

If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.

C++ Return value: INFO (output) INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, the algorithm failed to converge; i off-diagonal elements of an intermediate tridiagonal form did not converge to zero.


.. Parameters ..

C_DSYEV(): Computes all eigenvalues and, optionally, eigenvectors of a real symmetric matrix A.

These arguments mimic their Fortran counterparts.

Parameters
jobz= 'N' or 'n' = compute eigenvalues only; 'V' or 'v' = compute both eigenvalues and eigenvectors.
uplo= 'U' or 'u' = A contains the upper triangular part of the matrix; 'L' or 'l' = A contains the lower triangular part of the matrix.
n= The order of the matrix A.
A= On entry, the two-dimensional array with dimensions n by lda. On exit, if jobz = 'V', the columns of the matrix contain the eigenvectors of A, but if jobz = 'N', the contents of the matrix are destroyed.
lda= The second dimension of A (i.e., the number of columns allocated for A).
w= The computed eigenvalues in ascending order.
work= An array of length lwork. On exit, if the return value is 0, work[0] contains the optimal value of lwork.
lwork= The length of the array work. A useful value of lwork seems to be 3*N.

Returns: 0 = successful exit <0 = the value of the i-th argument to the function was illegal >0 = the algorithm failed to converge.

Interface written by TDC, 10/2002

size_t psi::C_IDAMAX ( size_t  length,
double *  x,
int  inc_x 
)

This function returns the index of the largest absolute value compoment of this vector.

Parameters
lengthNumber of elements in x.
xA pointer to the beginning of the data in x. Must be of at least length (1+(N-1)*abs(inc_x).
inc_xhow many places to skip to get to next element in x
Returns
the index of the largest absolute value
int psi::cc_excited ( const char *  wfn)

cc_excited(): This function takes a WFN string and returns 1 if the WFN is an excited-state method and 0 if the WFN is a ground-state method.

Parameters
*wfn= wavefunction string

Returns: 1 if an excited state method, else 0

int psi::cc_excited ( std::string  wfn)

cc_excited(): This function takes a WFN string and returns 1 if the WFN is an excited-state method and 0 if the WFN is a ground-state method.

Parameters
wfn= wavefunction string

Returns: 1 if an excited state method, else 0

double psi::combinations ( int  n,
int  k 
)

combinations() : Calculates the number of ways to choose k objects from n objects, or "n choose k"

Parameters:

Parameters
n= number of objects in total
k= number of objects taken at a time

Returns: number of combinations of n objects taken k at a time ("n choose k") (returned as a double).

int psi::david ( double **  A,
int  N,
int  M,
double *  eps,
double **  v,
double  cutoff,
int  print 
)

david(): Computes the lowest few eigenvalues and eigenvectors of a symmetric matrix, A, using the Davidson-Liu algorithm.

The matrix must be small enough to fit entirely in core. This algorithm is useful if one is interested in only a few roots of the matrix rather than the whole spectrum.

NB: This implementation will keep up to eight guess vectors for each root desired before collapsing to one vector per root. In addition, if smart_guess=1 (the default), guess vectors are constructed by diagonalization of a sub-matrix of A; otherwise, unit vectors are used.

TDC, July-August 2002

Parameters
A= matrix to diagonalize
N= dimension of A
M= number of roots desired
eps= eigenvalues
v= eigenvectors
cutoff= tolerance for convergence of eigenvalues
print= Boolean for printing additional information

Returns: number of converged roots

void psi::dirprd_block ( double **  A,
double **  B,
int  rows,
int  cols 
)

dirprd_block()

This function takes two block matrices A and B and multiplies each element of B by the corresponding element of A

Parameters
A= block matrix A
B= block matrix B
nrows= number of rows of A and B
ncols= number of columns of A and B

Returns: none

double psi::dot_block ( double **  A,
double **  B,
int  rows,
int  cols,
double  alpha 
)

dot_block(): Find dot product of two block matrices

Parameters
A= block matrix A
B= block matrix B
nrows= number of rows of A and B
ncols= number of columns of A and B
alpha= scale factor by which the dot product is multiplied

Returns: dot product

double psi::factorial ( int  n)

factorial(): Returns n!

Parameters:

Parameters
n= number to take factorial of

Returns: n factorial, as a double word (since n! can get very large).

void psi::fill_sym_matrix ( double **  A,
int  size 
)

fill_sym_matrix(): Fills a symmetric matrix by placing the elements of the lower triangle into the upper triangle.

Parameters
A= matrix to symmetrize
size= number of rows or columns (assume square)

Returns: none

void psi::free_3d_array ( double ***  A,
int  p,
int  q 
)

free_3d_array(): Free a (non-contiguous) 3D array

Parameters
A= triple-star pointer to the 3D array
p= size of first dimension
q= size of second dimension

Returns: none

double *** psi::init_3d_array ( int  p,
int  q,
int  r 
)

3d_array(): Initialize a (non-contiguous) 3D array

Parameters
p= size of first dimension
q= size of second dimension
r= size of third dimension

Returns: triple-star pointer to 3D array

double psi::invert_matrix ( double **  a,
double **  y,
int  N,
std::string  out 
)

INVERT_MATRIX(): The function takes the inverse of a matrix using the C routines in Numerical Recipes in C.

Matt Leininger, Summer 1994

Parameters:

Parameters
a= matrix to take the inverse of (is modified by invert_matrix())
y= the inverse matrix
N= the size of the matrices
outfile= file for error messages

Other variables: col and indx are temporary arrays d is 1 or -1

Returns: double (determinant) Note: The original matrix is modified by invert_matrix()

void psi::mat_print ( double **  matrix,
int  rows,
int  cols,
std::string  out 
)

mat_print(): Prints a matrix to a file in a formatted style

Parameters
matrix= matrix to print
rows= number of rows
cols= number of columns
outfile= output file pointer for printing

Returns: Always returns zero...

void psi::normalize ( double **  A,
int  rows,
int  cols 
)

normalize(): Normalize a set of vectors

Assume we're normalizing the ROWS

Parameters
A= matrix holding vectors to normalize
rows= number of rows in A
cols= number of columns in A

Returns: none

David Sherrill, Feb 1994

void psi::parallel_timer_off ( const std::string &  key,
int  thread_rank 
)

parallel_timer_off(): Turn off the timer with the name given as an argument. Can be turned on and off, time will accumulate while on. Should only be called in OpenMP parallel sections.

Parameters
key= Name of timer
void psi::parallel_timer_on ( const std::string &  key,
int  thread_rank 
)

parallel_timer_on(): Turn on the timer with the name given as an argument. Can be turned on and off, time will accumulate while on. Should only be called in OpenMP parallel sections.

Parameters
key= Name of timer
int psi::pople ( double **  A,
double *  x,
int  dimen,
int  ,
double  tolerance,
std::string  out,
int  print_lvl 
)

POPLE(): Uses Pople's method to iteratively solve linear equations Ax = b

Matt Leininger, April 1998

Parameters
A= matrix
x= initially has vector b, but returns vector x.
dimen= dimension of vector x.
num_vecs= number of vectors x to obtain.
tolerance= cutoff threshold for norm of expansion vector.

Returns: 0 for success, 1 for failure

int psi::ras_set3 ( int  nirreps,
int  nmo,
int *  orbspi,
int *  docc,
int *  socc,
int *  frdocc,
int *  fruocc,
int *  restrdocc,
int *  restruocc,
int **  ras_opi,
int *  core_guess,
int *  order,
int  ras_type,
bool  is_mcscf,
Options &  options 
)

ras_set3()

Updated version of ras_set2() to conform with new DETCI naming of various orbital spaces, particularly to distinguish between restricted orbitals (those that are doubly occupied but can optimized in an MCSCF) and frozen orbitals (those that are doubly occupied but cannot be optimized in an MCSCF). Also remove some deprecated code supporting Mark Hoffmann ordering, etc.

This function sets up the number of orbitals per irrep for each of the RAS subspaces [frozen core, restricted core, RAS I, RAS II, RAS III, RAS IV, restricted virts, frozen virts]. It also obtains the appropriate orbital reordering array. The reordering array takes a basis function in Pitzer ordering (orbitals grouped according to irrep) and gives the corresponding index in the RAS numbering scheme. Orbitals are numbered according to irrep within each of the subspaces.

Guesses for docc and socc should be provided, and they will be left as-is if they are not present in input (maybe auto-guessed by program)

C. David Sherrill

Parameters
nirreps= num of irreps in computational point group
nmo= number of MO's
orbspi= array giving num symmetry orbitals (or MOs) per irrep
docc= array of doubly occupied orbitals per irrep (guess must be provided)
socc= array of singly occupied orbitals per irrep (guess must be provided)
frdocc= array of frozen core per irrep (returned by function, but allocate before call)
fruocc= array of frozen virtuals per irrep (returned by function, but allocate before call)
rstrdocc= array of restricted core per irrep (returned by function, but allocate before call)
rstruocc= array of restricted unoccupied per irrep (returned by function, but allocate before call)
ras_opi= matrix giving num of orbitals per irrep per ras space, addressed as ras_opi[ras_space][irrep] (returned by function, but allocate before call) RAS IV will not be deduced unless RAS III is present in user input (in which case RAS IV needs to make up any remaining orbitals)
core_guess= array of the core orbitals per irrep to drop (must be provided). This is copied into FROZEN_DOCC (if is_mcscf == false) or RESTRICTED_DOCC (if is_mcscf == true), if both of those two arrays are absent from input
order= array nmo big which maps Pitzer to Correlated order (returned by function, but allocate before call)
ras_type= if 1, put docc and socc together in same RAS space (RAS I), as appropriate for DETCI. If 0, put socc in its own RAS space (RAS II), as appropriate for CC.
is_mcscf= true if an MCSCF-type wavefunction (orbital optimization), else false. Used to determine default behavior of FROZEN_DOCC, RESTRICTED_DOCC, FROZEN_UOCC, RESTRICTED_UOCC
options= Options object used to parse user input

Rules for how spaces can be specified:

  1. DOCC and SOCC are not changed by this routine (guesses might have been auto-generated by Psi) unless explicitly given in user input, in which case the guess arrays are overwritten by the user input in this routine
  2. FROZEN_DOCC and RESTRICTED_DOCC will be as given by the user. If the user does not specify either of these keywords, then maybe they specified FREEZE_CORE = TRUE instead, and we want to support that. If neither FROZEN_DOCC nor RESTRICTED_DOCC is specified, then we will use the array core_guess[] to get them. This array is obtained from whatever routine figures out how many orbitals per irrep to freeze when FREEZE_CORE = true (passed down as an argument into this function). If the computation is an MCSCF, then we'll copy core_guess into RESTRICTED_DOCC. If it's not an MCSCF, then we'll copy it into FROZEN_DOCC.
  3. The user can give either ACTIVE (for a CAS type computation, where internally this is basically like RAS 2), or RAS keywords. Not both.
  4. By default, unused virtual orbitals go into FROZEN_UOCC unless we're told this is an MCSCF, in which case they should go by default into RESTRICTED_UOCC
  5. RAS 4 is never used unless explicitly given in input; it's a special feature not normally desired. (Note: this is a change from previous default behavior in which unused orbitals were stuffed into RAS 4)

Returns: 1 for success, 0 otherwise

void psi::reorder_qt ( int *  docc_in,
int *  socc_in,
int *  frozen_docc_in,
int *  frozen_uocc_in,
int *  order,
int *  orbs_per_irrep,
int  nirreps 
)

reorder_qt()

This function constructs a reordering array according to the "Quantum Trio" standard ordering, in which the orbitals are divided into the following sets: frozen core, then doubly occupied, then singly occupied, then virtuals, then deleted (frozen) virtuals. The reordering array takes a basis function in Pitzer ordering (orbitals grouped according to irrep) and gives the corresponding index in the Quantum Trio numbering scheme.

Should give the same reordering array as in the old libread30 routines.

C. David Sherrill Center for Computational Quantum Chemistry University of Georgia, 1995

Parameters
docc_in= doubly occupied orbitals per irrep
socc_in= singly occupied orbitals per irrep
frozen_docc_in= frozen occupied orbitals per irrep
frozen_uocc_in= frozen unoccupied orbitals per irrep
order= reordering array (Pitzer->QT order)
nirreps= number of irreducible representations
void psi::reorder_qt_uhf ( int *  docc,
int *  socc,
int *  frozen_docc,
int *  frozen_uocc,
int *  order_alpha,
int *  order_beta,
int *  orbspi,
int  nirreps 
)

reorder_qt_uhf()

Generalization of reorder_qt() for UHF case

Parameters
docc= doubly occupied orbitals per irrep
socc= singly occupied orbitals per irrep
frozen_docc= frozen occupied orbitals per irrep
frozen_uocc= frozen unoccupied orbitals per irrep
order_alpha= reordering array for alpha (Pitzer->QT order)
order_beta= reordering array for beta (Pitzer->QT order)
nirreps= number of irreducible representations
void psi::schmidt ( double **  A,
int  rows,
int  cols,
std::string   
)

SCHMIDT(): Gram-Schmidt orthogonalize a set of vectors

Assume we're orthogonalizing the ROWS, since in C a vector is usually a row more often than a column.

David Sherrill, Feb 1994

Parameters
A= matrix to orthogonalize (matrix of doubles)
rows= rows of A
cols= columns of A

Returns: none

int psi::schmidt_add ( double **  A,
int  rows,
int  cols,
double *  v 
)

SCHMIDT_ADD(): Assume A is a orthogonal matrix. This function Gram-Schmidt orthogonalizes a new vector v and adds it to matrix A. A must contain a free row pointer for a new row. Don't add orthogonalized v' if norm(v') < NORM_TOL.

David Sherrill, Feb 1994

Parameters
A= matrix to add new vector to
rows= current number of rows in A (A must have ptr for 'rows+1' row.)
cols= columns in A v = vector to add to A after it has been made orthogonal to rest of A

Returns: 1 if a vector is added to A, 0 otherwise

void psi::solve_2x2_pep ( double **  H,
double  S,
double *  evals,
double **  evecs 
)

solve_2x2_pep(): Solve a 2x2 pseudo-eigenvalue problem of the form [ H11 - E H12 - E*S ] [c1] [ H12 - E*S H22 - E ] [c2] = 0

Parameters
H= matrix to get eigenvalues of
S= overlap between states 1 & 2
evals= pointer to array to hold 2 eigenvalues
evecs= matrix to hold 2 eigenvectors

Returns: none

void psi::timer_done ( void  )

timer_done(): Close down all timers and write results to timer.dat

void psi::timer_init ( void  )

timer_init(): Initialize the linked list of timers

void psi::timer_off ( const std::string &  key)

timer_off(): Turn off the timer with the name given as an argument. Can be turned on and off, time will accumulate while on. Should only be called out of OpenMP parallel section.

Parameters
key= Name of timer
void psi::timer_on ( const std::string &  key)

timer_on(): Turn on the timer with the name given as an argument. Can be turned on and off, time will accumulate while on. Should only be called out of OpenMP parallel section.

Parameters
key= Name of timer