Here's a compendium of C/C++ libraries that I've put together
to carry out various computational tasks that seem to arise
frequently in my research. Some of these are C/C++ wrappers around
legacy Fortran codes; one of them is a hacked version of another opensource
C/C++ code; and the rest are my own creations. (Which is which is indicated
in the
table
below.)
In this table,

“wrapper” denotes a Fortran code that
I found online and repackaged as a C/C++ library,
sometimes with slight modifications or extensions

“hack” denotes a C/C++ code that I found
online and modified or extended

“original” denotes an original
creation, i.e. a homemade code that I created from
scratch.
Library 
Purpose 
Provenance 
Libraries for fitting and interpolation

libPolyFit

1D polynomial fitting of data

Original 
libMDInterp

Interpolation of functions in 1, 2, 3, or 4 dimensions 
Original 
Libraries for special functions

libAmosBessel

Bessel functions of complex arguments 
Wrapper 
Libraries for numerical integration

libvvqag

Adaptive quadrature of onedimensional functions 
Hack 
libTSQ

Tanhsinh quadrature of onedimensional functions

Original 
libDCUHRE

Adaptive cubature of multidimensional functions
over hypercubic regions

Wrapper 
libDCUTRI

Adaptive cubature of functions over triangles

Wrapper 
libTriInt

Fixedorder cubature of functions over triangles

Original 
Libraries for Matrix and Vector Manipulations

libhmat

Linear algebra operations on matrices and vectors

Original 
Libraries for Electromagnetic Scattering

libRWG

A boundaryelement field solver for 3D scattering problems

Original 
libTDRT

A boundaryelement field solver for 2D scattering problems

Original 
libIncField

Calculation of electromagnetic fields
due to various types of known sources

Original 
libMatProp

Efficient database maintenance for frequencydependent
material properties

Original 
libSpherical

Convenience routines for computations in spherical coordinates

Original 
libSphericalScattering

Electromagnetic scattering from spheres and spherical shells

Original 
Libraries for Quantum Chemistry

libHQC

An infrastructure library for quantum chemistry
calculations with contracted
Gaussian orbitals

Original 
libHHF

A HartreeFock solver for quantum chemistry

Original 
Miscellaneous Utility Libraries

libhrutil

Convenient utility functions for various purposes

Original 
In general, you should be able to just download and unpack the
source packages and do a simple make to compile
each library. (There is no make install ; each
library here has just one .h and one .a
file, and it is up to you whether you want to copy them
somewhere else or leave them where they are.)
In addition to the code snippets below, all source packages
come with a test/demonstration program named
tlibName (where libName is
the library in question; thus, tlibRWG ,
tlibPolyFit , etc.) You can build this program by
doing a make test, and then (if this succeeds)
review its commandline options by doing a
tlibName help.

Many of the libraries here are simple enough that you can understand
what they offer just by looking at the code snippets, the example
programs, and the header files. For the largerscale libraries, I have
provided API documentation.

Some of the libraries on this page depend on other libraries on
this page (“internal” dependencies) or on opensource
libraries available elsewhere (“external” dependencies).
These dependencies are specified in the summaries below.

By default, the Makefile attempts to resolve
internal dependencies by assuming that you have unpacked and
built any other required library packages in the same directory
as the package you are trying to build. Thus, if you are trying
to compile libSpherical (which depends on
libAmosBessel ) and you have unpacked the
libSpherical source package into the directory
~HRLibs/libSpherical , then the Makefile
will attempt to reference the file
~HRlibs/libAmosBessel/libAmosBessel.h
(and also ~HRlibs/libAmosBessel/libAmosBessel.a
if you attempt to build the test program tlibSpherical .)

Send questions, bug reports, job offers, and (particularly) reports of
scenarios in which you found any of these useful to:
.
libPolyFit: Polynomial Fitting of 1D Data
The Setting:  

You have a collection of data points (x_{i},y_{i})
that you want to fit to a onedimensional polynomial y=f(x).
You choose the order of the polynomial, which may be greater than,
equal to, or less than the number of data points.

Download:  
libPolyFit.tar.gz

Dependencies:  
libhrutil

Provenance:  
Original.

Code snippet:  
int NumPts=6;
double X[NumPts], Y[NumPts];
PolyFit *PF4, *PF8;
/* fill in data arrays */
GetExperimentalData(X, Y, NumPts);
/* fit the data to an 4th order polynomial */
PF4=new PolyFit(X,Y,NumPts,4);
/* fit the data to an 8th order polynomial */
PF8=new PolyFit(X,Y,NumPts,8);
/* plot the fits with the data points */
PF4>PlotFit(X, Y, NumPts);
PF8>PlotFit(X, Y, NumPts);
/* use the fits to extrapolate to X=0 */
printf("Extrapolated to X=0 (4th order): %e\n",PF4>f(0.0));
printf("Extrapolated to X=0 (8th order): %e\n",PF8>f(0.0));

libMDInterp: Interpolation of functions in 1, 2, 3, or 4 dimensions
The Setting:  

You have a collection of functions of 1, 2, 3, or 4 variables
that are expensive to compute. (A common example is a definite
integral, which depends on a number of external parameters,
and which cannot be evaluated in closed form but which can
be evaluated numerically for arbitrary values of the parameters.)
You need to evaluate your function(s) at a large number of
points within some 1, 2, 3, or 4 dimensional hypercubic region
(which may be infinite), and you don't want to pay the cost
of recomputing the function at each separate point.
Instead, you want to precompute function values at
the points of a 1, 2, 3, or 4dimensional grid spanning your
hypercube, store these data in a big table, and then obtain
approximate values of your function at arbitrary points
inside the hypercube by interpolating among the known function
values at the grid points.
Often, you will want to do the precomputation step ahead of time
and save the grid data in a big data file, which you can then
subsequently reread into memory whenever you have a code that needs
to compute values of your function(s).

Download:  
libMDInterp.tar.gz

Dependencies:  
libhmat ,
libhrutil

Provenance:  
Original.

Code snippet:  
/***************************************************************/
/* 3D function to be interpolated ******************************/
/***************************************************************/
void My3DFunction(double X1, double X2, double X3, double *PhiVD)
{
PhiVD[0] = exp(X1*X1X2*X2X3*X3); /* value */
PhiVD[1] = 2.0*X1*exp(X1*X1X2*X2X3*X3); /* d/dX */
...
PhiVD[7] = 8.0*X1*X2*X3*exp(X1*X1X2*X2X3*X3); /* d^3/dXdYdZ*/
}
/***************************************************************/
/***************************************************************/
/***************************************************************/
int main()
{
int N1=4;
double X1Points[4]={0.0, 0.1, 0.2, 0.3};
int N1=5;
double X2Points[5]={0.2, 0.4, 0.6, 0.8, 1.0};
int N1=7;
double X3Points[7]={3.0, 2.0, 1.0, 0.0, 1.0, 2.0, 3.0};
/* initialize the interpolator */
Interp3D *I3D = new Interp3D(X1Points, N1,
X2Points, N2,
X3Points, N3,
1, nThread, My3DFunction);
/* get interpolated function values */
double F;
I3D>Evaluate(0.123, 0.234, 0.345, &F);
...
}

libAmosBessel: Bessel functions of general complex arguments
The Setting:  

Download:  
libAmosBessel.tar.bz2

Dependencies:  
None

Provenance:  
C++ wrapper around a
Fortran code
by D. E. Amos. This is the code used by Matlab for
evaluating Bessel and Airy functions.

Code snippet:  
cdouble z;
double MinOrder = 0.25;
int NumOrders=3;
cdouble JValues[NumOrders], iValues[NumOrders], AValue;
z = 1.2 + 3.4fi; // arbitrary complex number
AmosBessel('J', z, MinOrder, NumOrders, 0, JValues);
// at this point i have
// JValues[0] = regular cylindrical bessel func J_{0.25}(z)
// JValues[1] = regular cylindrical bessel func J_{1.25}(z)
// JValues[2] = regular cylindrical bessel func J_{2.25}(z)
AmosBessel('i', z, MinOrder, NumOrders, 0, iValues);
// at this point i have
// iValues[0] = regular modified spherical bessel func i_{0.25}(z)
// iValues[1] = regular modified spherical bessel func i_{1.25}(z)
// iValues[2] = regular modified spherical bessel func i_{2.25}(z)
AmosAiry('A', z, 0, &AValue);
// at this point i have AValue = airy function Ai(z)

libhmat: Manipulation of Matrices and Vectors
The Setting:  

Download:  
libhmat.tar.bz2

Dependencies:  
libC2ML

Provenance:  
Original.

Code snippet:  
// create a matrix and fill in its entries by hand
int m, n;
HMatrix *M1=new HMatrix(DIM, DIM, LHM_SYMMETRIC);
for(m=0; m<DIM; m++)
for(n=m; n<DIM; n++)
M1>SetEntry(m,n, MyMatrixEntryFunction(m,n) );
// read in a matrix from a matlabstyle ascii data file
HMatrix *M2=new HMatrix("MyMatrixDataFile.mat");
// create and initialize a vector (with random entries in this case)
HVector *V=new HVector(DIM);
for(m=0; m<DIM; m++)
V>SetEntry(m, lrand48() );
// now factorize the matrix and solve the linear system
Success=M>CholFactorize();
if ( Success )
{
// if M is positivedefinite, then i can solve the linear
// system using cholesky factorization
M>CholSolve(V);
}
else
{ // otherwise, i fall back on solving the system
// by LU factorization
M>LUFactorize();
M>LUSolve(V);
};

libMatProp: Database Maintenance for FrequencyDependent
Material Properties
The Setting:  

You want to do computational electromagnetism (possibly
with one or more of my software packages) involving
materials with frequencydependent material properties
(permittivity ε and permeability μ).
For each of the various materials you will be working with,
you have either (a) tabulated data for ε and
μ at a range of frequencies, or (b)
functional forms (such as a sum of lorentzians) for
ε(ω) and μ(ω), which you would like
to store somewhere in a universal database file which can be
used over and over again by various software tools.

Download:  
libMatProp.tar.bz2

Dependencies:  
libhrutil,
muParser

Provenance:  
Original.

Code snippet:  
//
// material properties specified by functional forms
// tabulated in a global database file
//
MatProp *MP1 = new MatProp("Gold");
MatProp *MP2 = new MatProp("Bromobenzene");
// permittivity of gold at omega = 1e14 radians / second
cdouble Eps = MP1>GetEps(1.0e14, REAL_FREQ);
// permeability of bromobenzene at \omega = 1e14 on imaginary axis
double Mu = MP2>GetMu(1.0e14, IMAG_FREQ);
//
// material properties obtained by interpolating among
// values in a data file
//
MatProp *MP3 = new MatProp("FILE_MyInterpolatedDataFile.dat");
// often it's more convenient to work in terms of
// units in which c==1 and the default length unit is 1 um
MatProp::SetFreqUnit( SPEED_OF_LIGHT / 1.0e6 );
MP3>GetEpsMu(0.12, &Eps, &Mu);

libSpherical: Utilities for Computations in 3D Spherical Polar Coordinates
The Setting:  

You need to convert coordinates, vectors,
finitedifference expressions for div, grad, and curl,
and other entities back and forth between cartesian and
spherical polar coordinates.
You need to do computations involving solutions of
the scalar or vector Helmholtz equation, with either
real or imaginary values of the k parameter,
i.e. spherical Bessel functions and spherical harmonics.
You want to do electromagnetic scattering from
spherical objects using my
libSphericalScattering code.

Download:  
libSpherical.tar.bz2

Dependencies:  
libC2ML

Provenance:  
Original.

Code snippet:  
/*
* given the cartesian coordinates of a point in space,
* compute the polar coordinates of the same point
* note: C2S = 'cartesian to spherical'
* S2C = 'spherical to cartesian'
*/
double xC[3], r, Theta, Phi;
xC[0] = 0.1; // x coord
xC[1] = 0.2; // y coord
xC[2] = 0.3; // z coord
CoordinateC2S(xC, &r, &Theta, &Phi);
/*
* given the x, y, z components of a vector at some
* point in space, compute the r, theta, phi components
* of the same vector
*
* note: on entry to VectorC2S i have
* vC[0,1,2] = x, y, z components of vector
* on return i have
* vS[0,1,2] = r, theta, phi components of vector
*/
double vC[3], vS[3];
vC[0] = 0.1; // x component
vC[1] = 0.2; // y component
vC[2] = 0.3; // z componet
VectorC2S(Theta, Phi, vC, vS);
/*
* get the values of the 'M' and 'N' vectorvalued
* Helmholtz solutions for a range of l,m values
*/
int lMax=3;
int NAlpha=(lMax+1)*(lMax+1); // total # of functions with l<= lMax
cdouble MInt[3*NAlpha]; // interior 'M' functions
cdouble NInt[3*NAlpha]; // interior 'N' functions
cdouble MExt[3*NAlpha]; // exterior 'M' functions
cdouble NExt[3*NAlpha]; // exterior 'N' functions
GetMNlmArray(lMax, k, 1, r, Theta, Phi, INTERIOR, MInt, NInt);
GetMNlmArray(lMax, k, 1, r, Theta, Phi, EXTERIOR, MExt, NExt);

libSphericalScattering:
Electromagnetic Scattering from Spherical Geometries
The Setting:  
You need to solve problems involving the scattering of
electromagnetic radiation from a spherical object
(either a solid sphere or a spherical shell).
You set up the problem by specifying

the incident field (which may be a plane wave,
the field of a point source, or any other arbitrary
incident field you like)

the frequency (which may be real or imaginary)

the radius (for a sphere) or the inner
and outer radii (for a spherical shell)

the material properties (permittivity and permeability)
of all regions: the exterior region, the scatterer,
and the innermost region in the case of a spherical
shell. Lossy materials are allowed;
permittivities may have negative real parts and
nonzero imaginary parts.
Then, libSphericalScattering solves the problem
for you, after which point you can call a GetFields()
routine to compute the scattered E and H fields at
arbitrary points in space (both inside and outside the scatterer).

Download:  
libSphericalScattering.tar.bz2

API Documentation:  
libSphericalScattering.pdf

Dependencies:  
libSpherical,
libAmosBessel,
libhmat,
libC2ML,
libIncField,
libSGJV,
libhrutil

Provenance:  
Original.

Code snippet:  
// create an instance of the SSGeometry class
SSGeometry *SSG = new SSGeometry();
// set frequency, geometry, material properties, etc.
SSG>SetParameters("freq 0.1 R 0.2 rEpsI 13.0 MuI 0.4");
// prepare to solve scattering problems using the given parameters
SSG>Setup();
// solve a scattering problem involving a particular incident field
SSG>Solve(MyIncidentFieldRoutine, MyUserData, EXTERIOR);
// now compute scattered E and H fields at arbitrary points
SSG>GetFields(X1, EH1);
SSG>GetFields(X2, EH2);
...

libhrutil: Miscellaneous Utility Functions
The Setting:  

Download:  
libhrutil.tar.bz2

Dependencies:  
None.

Provenance:  
Original.

Code snippet:  
#include
int main(int argc, char *argv[])
{
// Vararg versions of common functions
char *FB=GetFileBase(argv[0]);
FILE *f=vfopen("%s.out","w",FB);
int err=vsystem("grep iq %s.out 'MySymbol'",FB);
// Matlabstyle timing
Tic();
DoSomeExpensiveCalculation();
double SecondsElapsed=Toc();
// ... and many more similar utility functions
}

