%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%% Methods for Large-Scale Correlated Electronic Structure Calculations
%%%%% on Parallel Computers
%%%%% for NSF/MPS/CHEM/Theoretical  
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\chapter{Project Description}
\section{Introduction}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Computational chemistry has, in recent years, been transformed into a
virtually indispensable tool used by a large cross-section of the
discipline.  The ability to model ``real world'' systems with the
necessary sophistication to obtain chemically meaningful results has
helped give rise to a remarkable level of synergy between experimental
and theoretical treatments of chemical problems.  This, in turn, has
fueled interest in further expanding the role of computational
chemistry to even larger, more sophisticated, and more demanding
simulations.  In order to meet these challenges, we must insure our
ability to effectively utilize cutting edge high-performance computer
systems, and simultaneously look for new approaches which reduce the
computational cost of these simulations.

As the combination of hardware and software improvements enable ever
larger simulations, it has become clear that we must also examine the
traditional electronic structure methods used in these calculations
with an eye to improving their performance. Somewhat ironically, while
the cost of electronic structure calculations increases (often quite
rapidly) as the chemical systems increase in size, the opportunities
for methodological improvement also increase -- as the number of
interactions (represented by integrals, excitation amplitudes, etc.)
increase, the importance of any individual contribution decreases, and
of course the magnitude of these interactions fall off with distance
as well.

In density functional theory (DFT) and Hartree-Fock self-consistent
field (SCF) methodology, there has been a great deal of excitement
recently over so-called ``linear scaling'' techniques such as those
based on the Fast Multipole
Method\cite{greengard:1988:rapid,saunders:1992:electrostatic,panas:1992:fragment,Petersen:1994:PRSA-448-389,Petersen:1994:PRSA-448-401,Petersen:1994:JCP-101-8870,White:1994:JCP-101-6593,White:1994:CPL-230-8,Kutteh:1995:CPL-238-173,Burant:1996:CPL-248-43,Burant:1996:CPL-258-45,Challacombe:1996:CCR,Challacombe:1996:JCP-104-4685,Strain:1996:Science-271-51,White:1996:JCP-105-5061,White:1996:CPL-253-268,White:1996:CPL-257-647},
which promise calculations at a cost which depends linearly on the
size of the basis set ($N$) rather than the formal $O(N^{3})$ or $O(N^{4})$
for these methods.  While DFT and SCF methods are quite popular,
higher level (correlated) methods are crucial in providing reliable
results in many chemical applications.  Such methods, including
many-body perturbation theory, coupled cluster, and configuration
interaction, scale formally as $O(N^{5})$, $O(N^{6})$, and higher as the
sophistication and complexity of the method increase.  So far,
correlated methods have received relatively little attention aimed at
improving their performance, no doubt due, at least in part, to their
complexity relative to SCF and DFT.  Although multipole methods and
other approaches now being developed for SCF and DFT will undoubtedly
be applied to correlated methods, their utility in the more complex
methods is unclear at present.  Moreover, correlated methods offer
both challenges and opportunities which the linear scaling techniques
do not address.  It is therefore important to examine other promising
techniques which can provide alternatives to ideas imported from SCF
and DFT usage, or better yet work with them in a complementary fashion
to provide additional performance improvements.  This proposal
considers two approaches to correlated electronic structure
calculations which are well-suited for large-scale calculations, in a
sense, ``taking advantage'' of certain features of large chemical
systems to improve performance: the ``resolution of the identity''
integral approximation (RI), and the orbital-invariant formulation
(OI).

When considering methods appropriate to large-scale correlated
calculations, it is also appropriate to consider how such calculations
will be performed.  Pushing the boundaries of problem size usually
entails the use of high-performance computers.  Despite the seemingly
constant introductions of faster and faster CPUs, the performance of
``real'' applications not keeping pace -- performance is in practice
often limited by the speed of memory, bus, or I/O rather than by the
CPU.  This is especially true of electronic structure calculations,
which make intensive use of \textit{all} resources, not just CPU.
Consequently, virtually all modern high-performance computers employ
multiple processors to provide improved performance, and this trend
will continue.  In fact, parallelism is not limited to high-end
machines.  In the product lines of most computer vendors, one can find
multi-processor systems suitable for research groups, or even as
individual desktop machines.  There is enough accumulated experience
to show that automatic parallelization by compilers, and related
approaches are not sufficient to produce computational chemistry
software capable of effectively exploiting these architectures, but at
the same time, with appropriate algorithms, it is possible to drive
them quite efficiently.  Explicit consideration of parallel computing
is, therefore, necessary in order to achieve the largest calculations
and beneficial even for problems to be tackled with the local
computational resources of a research group or individual.

The resolution of the identity approximation and the orbital-invariant
formulation addressed in this proposal can be applied to a wide range
of correlated methods, and some preliminary work has been done with
both approaches in methods such as second- and fourth-order
perturbation theory, coupled cluster, and configuration interaction.
However neither approach is currently considered a generally
applicable alternative to the traditional formulation of these methods
for a variety of reasons. A key goal of my research program is to to
advance the state of the art in large-scale correlated calculations,
through the development of novel approaches such as RI and OI, as well
as by furthering the practical use of parallel computing.  As a
vigorous first step towards that goal, this proposal considers a
deeper investigation of the resolution of the identity approximation
and the orbital invariant formulation in the context of second-order
many-body perturbation theory (abbreviated MBPT(2) or MP2).  This
builds upon and extends my previous work with RI-MP2 and OI-MP2, and
will provide valuable insights which will greatly facilitate extending
the proposed work to higher-level correlated methods.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Second-Order Many-Body Perturbation Theory}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

MP2 is the simplest and computationally least expensive {\it ab
initio} electronic structure method which includes the effects of
dynamic electron correlation.  In many cases the MP2 method provides
acceptable results for properties for which correlation is important,
including reaction energetics, transition states, and electronic and
vibrational spectra.  Where MP2 itself is not sufficient, it provides
an inexpensive starting point which can be improved through the use of
higher levels of perturbation theory, coupled cluster, or
configuration interaction approaches.  As a consequence of its
computational and scientific efficiency, MP2 is the most commonly used
{\it ab initio} method after Hartree-Fock.  In the context of this
proposal, MP2 provides the most readily accessible framework in which
to refine the RI and OI approaches, while sharing many characteristics
with higher-level methods so that the proposed work has ramifications
well beyond just MP2.

The traditional form of the MP2 energy expression (in spin-orbital
form, following the usual notation\footnote{$i, j, \ldots m$ represent
occupied molecular orbitals, $a, b, \ldots e$ virtuals, and the
$\epsilon_p$ are the canonical SCF eigenvalues.}) can be given as,
\begin{equation} \label{eq:mo-e2-compact}
E^{(2)} = \frac{1}{2} \sum_{i,j,a,b} \frac{(ia|jb) [ (ia|jb) - (ib|ja) ]}
{\epsilon_i + \epsilon_j - \epsilon_a - \epsilon_b} .
\end{equation}
The energy evaluation itself has a formal $O(N^{4})$ scaling, but the
four-index transformation required to produce the molecular orbital
(MO) integrals from the atomic orbital (AO) integrals scales as
$O(N^{5})$.

For large systems which are spatially rather extended and use modest
atomic basis sets, one intuitively expects that distant parts of the
molecule do not interact strongly, leading to sparsity in the two-electron
integrals.  Empirical observation
\cite{clementi:1990:kgnmol,almlof:private-communication} and more
rigorous theoretical analysis \cite{dyczmons:1973:no} suggest that the
number of non-negligible two-electron integrals over Gaussian-like
atomic orbitals (AOs) has an asymptotic dependence of $O(N^2 (\ln
N)^2)$ compared to the formal $O(N^{4})$.  For such systems, an
approach which can take advantage of the locality of interactions,
such as ``orbital-invariant'' formulations (OI-MP2 when applied to the
MP2 method) has the potential to reduce the scaling significantly.  At
the other end of the spectrum are spatially compact systems with
extensive atomic basis sets.  Both factors mitigate against any useful
sparsity in such cases.  But an approach such as the ``resolution of
the identity'' integral approximation which takes advantage of the
near-redundancy in the basis set to simplify integral evaluation,
offers significant performance improvements.

It has been estimated that using a large, modern parallel computer
(512 node IBM SP) and an efficient parallel code, it will be possible
to perform exact MP2 gradient calculations with up to 1200 basis
functions ($C_{1}$ symmetry) before multi-passing through the integrals
is necessary \cite{Harrison:private-communication} (slightly larger for
MP2 energies only).  By targeting both ends of the spectrum of
chemical problems, the novel methods discussed in this proposal will
enable calculations on a broad variety of systems which would be
considered beyond the reach of traditional methods.  In the spatially
compact/large basis regime are ``benchmark'' calculations, where large
basis sets are applied to obtain high-quality results with a given
method, including extrapolations to the complete basis set limit.
Computational thermodynamics and other applications which require
accurate reaction energies would fall into this category. On the
spatially extended/modest basis end of the spectrum are, for example,
a wide range of systems of interest from the materials standpoint--for
molecular electronics, optical properties, nanomaterials, and numerous
other applications.  Precisely where along this spectrum the RI and OI
approaches are applicable is not understood at present, and this is
one of the questions the proposed work addresses.  However there is
reason to believe that the range of applicability of both approaches is
large, and perhaps partially overlapping.

In earlier research, Robert Harrison and I demonstrated parallel
implementations of both of these methods for energy evaluation
\cite{Bernholdt:1995:JCP-102-9582,Bernholdt:1996:CPL-250-447}.  This
proposal builds on that earlier work to develop RI-MP2 and OI-MP2 to
the point at which they can be put to routine use as serious
alternatives to traditional MP2 calculations for large-scale
calculations, especially those requiring high-performance (parallel)
computers.  After introducing the methods in more detail, I will
describe the necessary extensions to that previous work which
constitute the heart of this proposal.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{The Resolution of the Identity Approximation to MP2}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

The so-called ``resolution of the identity'' approach
\cite{Feyereisen:1993:CPL-208-359,Vahtras:1993:CPL-213-514,Kendall:1997:TCA-97-158}
has been applied to SCF
\cite{Vahtras:1993:CPL-213-514,Fruchtl:1997:IJQC-64-63}, MCSCF
\cite{Ten-No:1995:CPL-240-578,Ten-No:1996:IJQC-60-1319,Ten-No:1996:JCP-105-3604},
DFT
\cite{Eichkorn:1995:CPL-240-283,Eichkorn:1995:CPL-242-652,Eichkorn:1997:TCA-97-119},
MP2
\cite{Feyereisen:1993:CPL-208-359,Bernholdt:1996:CPL-250-447,Feller:1996:JCP-105-1940,Weigend:1997:TCA-97-331},
and coupled cluster methods \cite{Rendell:1994:JCP-101-400}.  It is
very closely related to the fitting of the density commonly used in
density functional theory
\cite{Dunlap:1979:JCP-71-3396,Mintmire:1982:PRA-25-88} and to the
multiplicative integral approximation of Van Alsenoy
\cite{vanalsenoy:1988:abinitio}.  The essence of the RI approach is the
replacement of pairwise products of the atomic orbital (AO) basis
functions, $\{\chi_p\chi_q\}$, by linear combinations of a ``fitting''
basis, $\{\alpha_t\}$. If $\{\alpha_t\}$ is the Gaussian Product
Theorem (GPT) result of the AO basis $\{\chi_p\}$ with itself, the
substitution is exact, while smaller fitting bases result in
approximations.  The success of this approach depends on the
observation that the AO product space is usually nearly linearly dependent
and can be spanned (or nearly so) by a smaller set of functions
\cite{beebe:1977:simplification,oneal:1989:application}.  Such an
approximation should work particularly well for calculations involving
large atomic basis sets, which provide more redundancy in the AO
product space and at the same time would provide little sparsity that
might be exploited in other approaches to reduce the computational
effort.

When applied to the familiar the four-center two-electron repulsion
integrals (ERIs),
\begin{equation}
	(pq|rs) = \int \chi_p(r_1) \chi_q(r_1) \frac{1}{|r_1 - r_2|}
		\chi_r(r_2) \chi_s(r_2) d^3 r_1 d^3 r_2 ,
\end{equation}
the RI approximation results in combinations of two- and three-center
one- and two-electron integrals.  Several criteria can be chosen for
optimizing the fitting basis representation of the AO product space,
resulting in different expansions of the four-center ERIs
\cite{Vahtras:1993:CPL-213-514}.  The ``V approximation'' 
\begin{equation}\label{eq:v-approx}
   (pq|rs) = \sum_{t,u} (pq|t) V^{-1}_{tu} (u|rs) ,
\end{equation} 
where $V_{tu}$ is the ERI $(t|u)$, is representative and currently
seems to be the fit of choice based on results of Vahtras and
coworkers, and similarity to the variational density fitting procedure
widely used in DFT.  This choice can be thought of as a least squares
fit of the electrostatic potential of the product space density as a
linear combination of the electrostatic potential arising from the
fitting density.

When the RI approximation is applied to MP2, the working equations are
obtained by the straightforward insertion of
Eq. \ref{eq:v-approx} into the MP2 energy expression
(Eq. \ref{eq:mo-e2-compact}), resulting in the equation
\begin{equation} \label{eq:ri-e2}
E^{(2)} = \frac{1}{2} \sum_{i,j,a,b,t^{\prime}} 
\frac{(ia|t^{\prime})(t^{\prime}|jb) [ (ia|t^{\prime})(t^{\prime}|jb) 
- (ib|t^{\prime})(t^{\prime}|ja) ]}
{\epsilon_i + \epsilon_j - \epsilon_a - \epsilon_b} ,
\end{equation}
where for convenience, the $V^{-1}$ term has been subsumed
symmetrically into the three-center integrals, following Rendell and
Lee \cite{Rendell:1994:JCP-101-400}.  

At the heart of RI's advantage is the is the reduced volume of
integrals to be calculated and manipulated--down from $O(N^{4})$ to
$O(N^{3})$, taking $N_{fit}$ to be $O(N)$.  The transformation of
integrals from the AO to MO basis is thus $O(N^{4})$, while a new
$O(N^{5})$ step is introduced, forming approximate four-center
integrals from products of three-center integrals using
Eq. \ref{eq:v-approx}, which is a large matrix multiplication.
When the RI approximation is applied to more complex methods,
intermediates are typically produced with fewer indices than in the
exact equations (i.e. three instead of four indices), and later
contracted or combined to produce approximations of exact quantities,
in analogy with the RI-MP2 case, but with the additional advantage
that smaller intermediates are carried through substantial parts of
the method.

Even with relatively small chemical systems (4 atoms, 140 basis
functions, 396 fitting functions) RI-MP2 calculations have shown as
much as a 12.4$\times$ performance improvement over exact MP2 on a
single processor workstation \cite{Bernholdt:submitted:FBS}.  When the
parallel implementation was applied to the K$^{+}$:12-crown-4 complex
(461 AO and 1227 fitting basis functions), the RI-MP2 calculation was
approximately 17 times faster than the exact MP2
\cite{Bernholdt:1996:CPL-250-447}, while accurately reproducing the
exact MP2 binding energy.  Although a sizable system for an exact MP2
calculation, the parallel speedup curves shown in
Fig. \ref{fig:rimp2-speedup} suggest that as an RI-MP2 calculation, it
is in fact somewhat small.  The distinctive features of the RI-MP2
method are well suited to the strengths of modern parallel processors.
Distributed matrix multiplication
\cite{choi:1993:pumma,husslederman:1994:CPE-6-571,vandeGeijn:1995:SSU}
performs much more efficiently than most parallel four-index
transformation algorithms which have been developed
\cite{clementi:1986:lecture,whiteside:1987:parallel,hurley:1988:optimized,bauschlicher:1989:efficient,watts:1989:vector,covick:1990:fourindex,wiest:1991:program,harrison:1992:abinitio,limaye:1994:general,dupuis:1994:modern,windus:1994:parallel,Wong:1996:TCA-93-317,Schutz:1997:TCA-95-13},
and the reduced number of integrals means disk-based implementations
can be used much more readily than for the exact MP2.  

\begin{figure}[tbh]
\begin{minipage}[b]{3.5in}
%\psfig{file=k12c4-speedup.eps,angle=00.0,width=3in}
\end{minipage}
\begin{minipage}[b]{3in}
\caption{Parallel speedup for the K$^{+}$:12-crown-4 RI-MP2 energy on
three different MPP platforms.  Error bars resulting from the
uncertainty associated with the load on the remainder of the machine
were obtained from analysis of repeated runs of the same job.  Error
bars for the Cray-T3D data are too small to be distinguished at this
scale.}
\label{fig:rimp2-speedup}
\end{minipage}
\end{figure}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{The Orbital-Invariant Formulation of MP2}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

While the resolution of the identity approximation is well-suited to
problems in which sizable AO basis sets are applied to molecules which
are spatially compact, problems at the other end of the spectrum are
equally interesting: those involving modest AO basis sets on large,
spatially extended systems.  As noted above, in such cases one expects
substantial sparsity in the interactions between distant parts of the
molecule.  It is well known that screening AO integral blocks prior to
evaluation will provide some of the benefits of this sparsity, but the
basis in which the interactions are represented can play an important
role in determining how much of this reduction will be present in
practical calculations.  Canonical Hartree-Fock MOs, normally used to
formulate correlated methods, tend to be delocalized.  This means that
even if there is significant sparsity in the AO basis representation
of the integrals, it will be ``filled in'' by the delocalization of
the MOs, and the number of non-negligible integrals will not approach
the asymptotic limit of $O(N^2 (\ln N)^2)$.

An orbital-invariant formulation is more general because it allows the
use of localized MOs (from
Edmiston-Ruedenberg\cite{edmiston:1963:localized},
Foster-Boys\cite{boys:1966:localized}, or Pipek-Mezey
\cite{pipek:1989:fast} procedures, for example) as occupied
orbitals (or even AOs projected into the occupied space). Following
the framework of the Method of Self-Consistent Electron Pairs
(SCEP)\cite{meyer:1976:theory,pulay:1984:efficient}, the virtual
orbitals may in fact be represented in any (non-orthogonal) basis that
spans the virtual space.  Thus instead of HF molecular orbitals, which
do not localize very well, it is possible to use the AO basis itself,
greatly improving the spatial locality of the individual functions
spanning the virtual space and thereby improving sparsity.  Orbital
invariant formulations have been reported for a number of perturbation
theory, coupled cluster, and configuration interaction-based
energies\cite{meyer:1976:theory,pulay:1984:efficient,Saebo:1987:JCP-86-914,Hampel:1996:JCP-104-6286},
but the approach cannot not be described as being widely used.

Rather than repeat the detailed derivation of the orbital invariant
MP2 (OI-MP2) equations, which have already been published
\cite{pulay:1986:orbitalinvariant,Bernholdt:1995:JCP-102-9582}, I will
briefly sketch the relationship between OI and traditional
formulations.  The MP2 energy of Eq. \ref{eq:mo-e2-compact} can also
be written with the excitation amplitudes, $t_{ij}^{ab}$ made
explicit, as
%\begin{subeqnarray}\label{eq:mo-e2}
%	E^{(2)} & = & 
%            \frac{1}{4} \sum_{i,j,a,b} ((ia|jb) - (ib|ja)) t_{ij}^{ab} \\
%	0 & = & (ia|jb) - (ib|ja)
%	+ (t_{ij}^{ab} f_{bb} + f_{aa} t_{ij}^{ab})
%	- (t_{ij}^{ab} f_{jj} + f_{ii} t_{ij}^{ab}),
%	\forall i, j, a, b.
%\end{subeqnarray}
Note that the SCF eigenvalues $\epsilon_{p}$ of
Eq. \ref{eq:mo-e2-compact} have been replaced with the equivalent
notation $f_{pp}$, denoting the diagonal elements of the MO basis Fock
matrix.

To accommodate localized orbitals, which do not diagonalize the
occupied-occupied (o-o) and virtual-virtual (v-v) blocks of the Fock
matrix, the amplitude equation above must include all elements of the
o-o and v-v blocks of the Fock matrix instead of just the diagonals
(the energy equation is unaffected),
\begin{equation}\label{eq:mo-t2}
	0 = r_{ij}^{ab} = (ia|jb) - (ib|ja)
	+ \sum_e (t_{ij}^{ae} f_{be} + f_{ae} t_{ij}^{eb})
	- \sum_m (t_{im}^{ab} f_{mj} + f_{im} t_{mj}^{ab}),
	\forall i, j, a, b.
\end{equation}
The residual $r_{ij}^{ab}$ is introduced at this stage as a reminder
that in this more general form, the equation must be solved iteratively
for $t$, with the residual measuring ``how far'' the current $t$ is
from the final solution.

The equations are now recast in the form of coupled matrix equations.
The amplitudes are collected into a set of matrices, $\{{\bf t}_{ij}\}$
which are labeled by virtual orbitals $a$ and $b$ in their row and
column dimensions.  Similarly, the integrals are collected as a set of
exchange operators, ${\bf K}_{ij}^{ab} = (ia|jb)$.  The MP2 energy
is now represented as 
\begin{equation}\label{eq:mo-matrix-e2}
	E^{(2)} = \frac{1}{4} \sum_{i,j} ( {\bf t}_{ij} \cdot 
	({\bf K}_{ij} - {\bf K}_{ij}^{\dagger}) ),
\end{equation}
and the amplitude equation becomes
\begin{equation}\label{eq:mo-matrix-t2}
	0 = {\bf r}_{ij} = {\bf K}_{ij} - {\bf K}_{ij}^{\dagger}
	+ {\bf F} {\bf t}_{ij} + {\bf t}_{ij} {\bf F}
	- \sum_m f_{im} {\bf t}_{mj} - \sum_m {\bf t}_{im} f_{mj}, \forall i, j.
\end{equation}

Any transformation ${\bf C}$ of the MOs which satisfies the
orthonormality constraint ${\bf C}^\dagger{\bf S}{\bf C} = {\bf 1}$ for
some metric ${\bf S}$ may be used to transform the excitation space.
Orthogonality of the occupied and excitation spaces can be guaranteed by
appropriate use of projection operators during solution of the amplitude
equations \cite{meyer:1976:theory}. This results in the most general
form of the MBPT(2) equations,
%\begin{subeqnarray}
%	E^{(2)} & = & \frac{1}{4} \sum_{i,j} ( {\bf t}_{ij} \cdot
%		( {\bf K}_{ij} - {\bf K}_{ij}^{\dagger} ) ) , \\
%	0 = {\bf S} {\bf r}_{ij} {\bf S}
%	& = &      {\bf K}_{ij} - {\bf K}_{ij}^{\dagger} 
%	+       {\bf F} {\bf t}_{ij} {\bf S} 
%	+ {\bf S} {\bf t}_{ij} {\bf F}
%	-       \sum_m f_{im} {\bf S} {\bf t}_{mj} {\bf S} 
%	- \sum_m {\bf S} {\bf t}_{im} {\bf S} f_{mj} .
%\end{subeqnarray}
In this form, the occupied orbitals can be localized (noncanonical) and
the virtual space can be any appropriate set of orbitals.  The AO basis
itself provides a very localized representation which can be used
rather than localized virtual orbitals, which are often hard to obtain.

The orbital invariant formulation makes it possible to retain a great
deal more sparsity than is possible using the traditional formulation.
This lays the groundwork for approaches that can take advantage of the
localization and sparsity.  One such idea is the ``local correlation''
approximation, originated by Pulay and Saeb{\o}, which extends the
localization concept to the excitation amplitudes by neglecting
excitations which go beyond some local neighborhood of the originating
orbital.  Although relatively little work has been done based on the
orbital-invariant approach, it should yield an MP2 procedure, for
example, with an $O(N^{3})$ scaling compared to the traditional
approach's $O(N^{5})$\/.  Saeb{\o} and coworkers have reported good
performance for the several applications of this approach
\cite{sellers:1986:ring,sabo:1987:effect,saebo:1993:local} using their
sequential implementation.  Very recently, Maslen and
Head-Gordon\cite{Maslen:inpress:NLS} have given another derivation of the
orbital-invariant MP2 which leads naturally to a local correlation
ansatz which is slightly different from that used by Pulay and
Saeb{\o}, though the concept is very much the same.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Relationship to Other Methods}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Methods specifically targeted at large-scale correlated calculations
have not yet received the kind of attention that has been given to
reducing the scaling of SCF and DFT calculations in recent years.  The
orbital invariant and resolution of the identity approaches represent
rather different ideas from those used in linear scaling methods, and
I anticipate that as the ``linear scaling revolution'' is extended to
correlated methods, these various approaches will be combined in
synergistic fashion to provide even greater benefit.

When considering the spatially compact/large basis set regime, where
sparsity is negligible, there do not seem to be any alternatives to
the RI-MP2 approach at this time.  When sparsity does come into play,
there are several approaches have been put forward.

The Laplace Transformation (LT) approach, proposed by Alml{\"o}f and
H{\"a}ser
\cite{almlof:1991:elimination,haser:1992:laplace,haser:1993:mollerplesset,Rauhut:1996:CPL-248-223}
allows the use of localized orbitals in a fashion similar to the
OI-MP2 equations.  Although I believe the approach to be very
promising for higher-level correlated methods (where it allows
high-complexity terms to be separated into less expensive components,
thus reducing scaling), it suffers from several deficiencies at the
MP2 level when compared with an orbital invariant formulation.
\begin{itemize}
\item The LT approach is approximate in the sense that it relies on a
numerical quadrature, while the OI-MP2 is exact. (The quadrature can be
evaluated to matching accuracy, but this can be quite expensive.)
\item As published, LT-MP2 lacks the flexibility in representing the
virtual space present in OI-MP2.  The limitation to localized HF
virtual orbitals would significantly reduce the sparsity which can be
obtained.
\item The LT-MP2's numerical quadrature requires multiple integral
evaluation and transformation steps, usually the least efficient steps
in a parallel code suitable for such large-scale calculations.  By
contrast, the OI-MP2 only requires a partial transformation, which can
be much more efficient, and can utilize memory and disk to reduce or
eliminate the need for integral recomputation.
\end{itemize}
As a result, I believe that the OI-MP2 approach is better suited for
application to spatially extended/modest basis calculations.

Pseudospectral approximations have been applied to a number of methods
ranging from SCF to full-CI including MP2 and MP3
\cite{ringnalda:1990:pseudospectral,langlois:1990:pseudospectral,martinez:1992:pseudospectral,martinez:1994:pseudospectral},
and might be considered as yet another approach to large-scale
calculations.  Pseudospectral approximations combine the usual atomic
orbital basis with a numerical grid, expanding the Coulomb potential
due to AO basis function products onto the grid.  Use of a physical
grid simplifies taking advantage of the locality of electronic
interactions.  However it must be pointed out that pseudospectral MP2
relies in essentially the same formulation as presented here to
provide localized orbitals which gives the pseudospectral approach its
performance advantage.  Essentially, pseudospectral MP2 is an integral
approximation applied to OI-MP2, and the effort described in this
proposal will apply equally with or without the additional
approximation.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{The Role of High-Performance Computing}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

These methods are designed push back the size limitations of
correlated calculations by ``taking advantage'' of certain features of
large chemical systems. This naturally suggests that high-performance
computers should be taken into account during their development. As
noted earlier, the trend in computing is such that
``high-performance'' must include parallel computers, and
increasingly, parallel computing extends all the way down to the
desktop.

The NWChem parallel computational chemistry package \cite{nwchem3.1},
developed by the Environmental Molecular Sciences Laboratory (EMSL) at
the Pacific Northwest National Laboratory provides an excellent
framework in which to perform the work proposed here.  It is one of
the few packages in existence which was designed with the specific
intention of providing fully-distributed highly-scalable computational
chemistry tools. The Global Array programming model
\cite{ga-tools:anon-ftp,nieplocha:1994:global,Nieplocha:1996:JSC-10-169},
which underlies NWChem, provides much of the ease-of-programming
common to shared-memory-only programming models, while retaining the
generality of the message-passing paradigm, thus greatly simplifying
parallel software development.  In addition, NWChem embodies a great
deal of experience in developing parallel MP2 codes in particular
using this model.  It presently includes independent
fully-direct\cite{Wong:1996:TCA-93-317} and
semi-direct\cite{Deegan:in-preparation} implementations of exact MP2,
the latter including gradients, as well as the parallel RI-MP2
energy\cite{Bernholdt:1996:CPL-250-447}; the first parallel
implementation of the OI-MP2 method\cite{Bernholdt:1995:JCP-102-9582}
served as a prototype for the design of NWChem and is being
reimplemented within the package.

Consequently, the proposed work can be carried out incorporating
parallelism with almost no additional effort or risk.  The benefits of
such an approach, however, extend well beyond this proposal, in that
they will lead to further improvements to the Global Array toolkit,
and provide additional examples and experience to help the
computational chemistry community assimilate parallel computing, now
available from desktop to high-end machines, into the standard
``tool box''.

As a continuing member of the NWChem development team, I have the full
support and cooperation of the EMSL for the proposed work, including
the use of NWChem (see supporting letter in Section
\ref{chap:supplementary-docs}).  Distribution of NWChem is managed by
the EMSL, and the intent is to make the package widely available for
use on parallel supercomputers.  In addition, certain components of
NWChem, such as the Global Array Toolkit
\cite{ga-tools:anon-ftp,nieplocha:1994:global,Nieplocha:1996:JSC-10-169}
and PeIGS parallel eigensolver \cite{peigs2.1} have been placed in the
public domain and are freely distributed as stand-alone packages to
facilitate their use in other applications.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Proposed Work}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Both OI- and RI-MP2 energies have been implemented already and their
promise has been demonstrated, both in terms of parallel scalability and
improved performance compared to traditional MP2 calculations.  The goal
of this project is to develop both methods to the point at which they
are sufficiently well understood, and provide the capabilities
necessary, to offer a useful (more efficient) alternative to traditional
MP2 calculations on large systems. In order to accomplish this goal,
a number of different issues must be addressed.
\begin{itemize}
\item For the RI-MP2 method, better understanding of the fitting basis
set is required, so that reliable results can be obtained without
a large number of expensive ``calibration'' runs.
\item For the OI-MP2 method, effective use of sparsity is crucial, but
little is known about how sparsity grows in for real chemical systems
or how to design parallel algorithms to use it efficiently.
\item For any routinely used method, analytic gradients and other
first-order properties are a practical necessity.  They allow efficient
optimization of structures and can be differentiated numerically to
produce vibrational frequencies.
\end{itemize}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{RI-MP2 Fitting Basis Sets}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

In addition to some work towards fitting basis sets for use in RI-DFT
\cite{Eichkorn:1995:CPL-240-283,Eichkorn:1995:CPL-242-652,Eichkorn:1997:TCA-97-119}
and RI-SCF and -MCSCF
\cite{Ten-No:1995:CPL-240-578,Ten-No:1996:JCP-105-3604,Ten-No:1996:IJQC-60-1319}
calculations, there has been some preliminary work on RI-MP2 fitting
sets.  Weigend and H\"aser \cite{Weigend:1997:TCA-97-331} have
employed fitting basis sets created without any optimization in
conjunction with the SVP AO basis \cite{Schafer:1992:JCP-97-2571} for
a series of calculations of reaction energies.  A more detailed study,
partly carried out in collaboration with Robert Harrison, has lead to
RI-MP2 fitting basis sets
\cite{Bernholdt:submitted:FBS,Bernholdt:to-be-published} for use in
conjunction with the correlation consistent AO basis sets
\cite{Dunning:1989:JCP-90-1007,Woon:1994:JCP-100-2975,Woon:li-be-cc-pvxz,Kendall:1992:JCP-96-6796}
cc-pVDZ, cc-pVTZ, and aug-cc-pVTZ.  The ``-fit2-1'' fitting sets
(i.e. cc-pVDZ-fit2-1 and aug-cc-pVTZ-fit2-1), are designed to
reproduce a set of 32 test cases involving molecules comprised of
zeroth- and first-row atoms to better than 1\% on average and a
maximum of 2\% error in any individual test case.  These tests were
chosen to represent a range of valance electron properties, including
reaction energies involving single and multiple bonds, hydrogen
bonding, van der Waals interactions, and electrostatic properties.
Table
\ref{tbl:fit2-1-summary} briefly describes the fitting sets and their
numerical performance.

\begin{table}
\caption{Performance RI-MP2 calculations with the ``-fit2-1'' fitting
basis sets.  ``Reactions'' refer to the error (compared to exact MP2)
results for the 32 test cases covering a range of interactions an
properties.  ``Total Energies'' refer to the error (compared to exact
MP2) of the 55 separate calculations which contribute to the test
cases.}
\label{tbl:fit2-1-summary}
\begin{tabular}{llccc}
\hline\hline
 & & \multicolumn{3}{c}{AO and Fitting Basis Set} \\
 & & cc-pVDZ & cc-pVTZ & aug-cc-pVTZ \\
 & & cc-pVDZ-fit2-1 & cc-pVTZ-fit2-1 & aug-cc-pVTZ-fit2-1 \\
\hline
AO Descr. & H-He (size) 
  & 2$s$1$p$ (5) & 3$s$2$p$1$d$ (15) & 4$s$3$p$2$d$ (25) \\
  & Li-Ne (size) 
  & 3$s$2$p$1$d$ (15) & 4$s$3$p$2$d$1$f$ (35) & 5$s$4$p$3$d$2$f$ (55) \\
\\
Fit Descr. & H-He (size) 
  & 9$s$6$p$4$d$ (51) & 7$s$6$p$5$d$4$f$ (95) & 7$s$6$p$5$d$4$f$ (95) \\
  & Li-Ne (size) 
  & 10$s$5$p$4$d$1$f$ (59) & 11$s$6$p$4$d$3$f$1$g$ (98) & 11$s$6$p$6$d$3$f$2$g$ (125) \\
\\
Reactions & Avg. Error (\%) & 0.13 & 0.06 & 0.11 \\
 & Max. Error (\%) & 1.3 & 0.4 & 0.5 \\
\\
Total Energies & Min. Error ($\mu E_{h}$) & 0.3 & $<$ 0.1 & $<$ 0.1 \\
 & Avg. Error (m$E_{h}$) & 0.2 & 0.3 & 0.3 \\
 & Max. Error (m$E_{h}$) & 0.6 & 1.5 & 1.4 \\
\hline\hline
\end{tabular}
\end{table}

The practical experience of this work agrees with the theoretical
expectation that as the AO basis increases in size, the AO product
space grows increasingly redundant, and therefore the fitting basis
required to fit it adequately is smaller in relative terms.  This
helps to alleviate concerns that calculations involving large fitting
basis sets will suffer from linear dependencies in the fitting basis
itself.  Small eigenvalues were observed in the aug-cc-pVTZ case, but
were easily handled without effecting the quality of the calculation.
While our work thus far has not addressed the use of this method on
potential energy surfaces, preliminary calculations by Weigend and
H\"aser \cite{Weigend:1997:TCA-97-331} are favorable.

The next logical step in the development of RI-MP2 fitting basis sets
is to develop a more detailed and systematic understanding of the role
of the fitting basis--both where errors arise and how to improve the
approximation in a cost-effective manner.  This will involve direct
comparison of exact and approximate integrals, categorization
according to the number of unique atomic centers involved in the
integrals, and numerical optimization procedures as well as
theoretical examination of the approximations targeted at bounding
their errors.

This work will consider not only energies and finite field property
calculations, but also gradients, potential energy surfaces, and a
broader range of other properties in conjunction with the availability
of relaxed density matrices (as a by-product of the gradient
implementation).

In addition, we will investigate the effectiveness of ideas for
constructing fitting basis sets beyond the usual atom-centered
approach.  On the one hand, it may be possible to draw on the work of
Beebe and Linderberg \cite{beebe:1977:simplification} and O'Neal and
Simons \cite{oneal:1989:application} to develop numerically-oriented
approaches to the creation of fitting sets. On the other hand we will
consider the use of off-atom fitting functions, and more
physically-based approaches -- for example, fitting each product using
only a local neighborhood of functions, which would effectively lead
to a fitting basis of constant size, thus eliminating a factor of
$N_{fit}$ from the cost of the method.

The work described here on understanding and developing fitting basis
sets will result both in improvements to the simple recipes, developed
in the preliminary work, as well as entirely new approaches to the
construction of fitting basis sets which I expect to be more accurate
and more efficient.  The detailed numerical analysis of the fitting
procedure will also help define the error bars for this method in
comparison with the exact MP2.  Moreover, much of the proposed work is
largely independent of the particular electronic structure method to
which the RI approximation is being applied.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{Locality and Sparsity in OI-MP2}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Orbital invariant formulations of correlated methods facilitate taking
advantage of the sparsity or locality of interactions in large quantum
mechanical systems.  Theoretically, we understand that the number of
two-electron integrals and double excitation amplitudes must
asymptotically approach $O(N^{2})$; practically, information about how
real chemical systems approach this is quite limited and spotty.
Because of the importance of sparsity and localized representations in
quantum chemical methods (not just at the MP2 level), this area needs
to be examined systematically and in some detail.

Thus far, studies of locality have focused largely on SCF-level
quantities, as in the recent work by Maslen and coworkers
\cite{Maslen:submitted:LSA} (though this particular work does briefly
consider sparsity of the MP2 density matrix as well), or have relied
solely on physical arguments to enforce a locality \textit{a priori},
as in the local correlation work of S{\ae}b{\o} and Pulay
\cite{pulay:1983:localizability,saebo:1985:local,saebo:1993:local},
and Maslen and Head-Gordon \cite{Maslen:inpress:NLS}.  In this work,
we propose to study sparsity and locality in correlated methods,
starting with the MP2 level, in order to connect the physically-based
approximations with a numerical understanding of the sparsity in
integrals, excitation amplitudes and related quantities.  The goal of
this work is error bounds for the physically-based local correlation
approximations, and approximations, controllable by numerical
thresholds and tolerances, which can bridge between exact MP2 and the
physically-based local correlation approximations.  The development of
such new approximations will involve controlled neglect of integrals,
amplitudes, and orbital pairs, and will involve a combination of
theoretical work and numerical experiments.  This work will be
extended to include gradients, once that capability is implemented for
OI-MP2 (see Section \ref{sect:gradients}).

Collateral to this study of the physical and numerical aspects of
sparsity are considerations of how it can be incorporated into
high-performance (parallel) algorithms efficiently, so that storage
requirements are minimized while performance is maximized.  Sparsity
is present in a wide variety of problems in computational science,
though typically with somewhat different characteristics than that in
chemical problems -- logically, a kind of ``block sparsity'' fits the
structure of many electronic structure applications, but this
matrix structure is not commonly used in other fields.

Work will initially focus on efficient storage and manipulation of
block sparse matrices, and later expand to emphasize improving
algorithms for computations on such matrices (i.e. matrix
multiplication, diagonalization, etc.) and providing interfaces to
other commonly used sparse data structures in order to take advantage
of existing packages.  This work will be done within the context of
the Global Array (GA) Toolkit
\cite{ga-tools:anon-ftp,nieplocha:1994:global,Nieplocha:1996:JSC-10-169},
the parallel programming model used in NWChem, and in close
collaboration with researchers at EMSL.  This allows our work to take
advantage of the extensive capabilities of the GAs for dense linear
algebra and facilitates rapid adoption of this work in a broad range
of other applications, both within computational chemistry and in
other areas, which already use GAs as their programming model.  It
will also be possible to use the lessons learned in this work to
further the development of a general, efficient parallel programming
language because of the conceptual similarity between the GAs and High
Performance Fortran \cite{HPF1.1,HPF2.0}, which has been accepted by
the computer industry as the {\em de facto} standard version of
Fortran for parallel computers.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{Analytic Gradients in OI- and RI-MP2}
\label{sect:gradients}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

These days, for a method to be practical for routine application to
chemical problems, it must have analytic gradients and other first-order
properties available.  This capability allows the efficient
determination of structures as well as offering the possibility of
computing vibrational spectra by finite difference of gradients.  Theory
and implementation of traditional, exact MP2 gradients are well
established and widely used.  A number of parallel implementations have
been developed as well
\cite{Nielsen:1996:CPL-255-210,Fletcher:1997:MP-91-431,Deegan:in-preparation}.
This proposal includes the corresponding development and implementation
which is required for both OI- and RI-MP2. 

The theory of OI-MP2 gradients was outlined by Pulay and Saeb{\o}
\cite{pulay:1986:orbitalinvariant}, and has been
implemented in the Jaguar software package from Schr\"odinger, Inc.~
in the context of the pseudospectral approximation, but neither
detailed equations nor a description of the implementation have been
published \cite{Murphy:private-communication}.  I have previously
derived detailed equations for the orbital-invariant form of the
coupled-perturbed Hartree-Fock (CPHF) equations
\cite{Bernholdt:unpublished}, which are required for MP2 gradients,
and am currently in the process of producing detailed equations for
the OI-MP2 gradient itself.  There is reason to believe that OI-MP2
gradients may in fact be more efficient than the traditional
formulation even without taking advantage of sparsity.  This stems
from the fact that in the traditional formulation, a four-index
back-transformation of the two-particle density is usually performed
in order to contract with the AO derivative integrals, while in the
OI-MP2 framework, we can use the AO space representation of the
virtual space directly which should simplify this step.  Our
implementation of OI-MP2 gradients will accommodate the use of UHF
reference functions for calculations on open-shell systems, and will
run efficiently on parallel machines.

As with the RI-MP2 energy, the RI-MP2 gradient equations can be
obtained straightforwardly by inserting the RI approximation into the
familiar exact gradient equations, though some care must be taken to
distinguish quantities arising from the RI-MP2, which involve
approximate integrals, from those arising from the SCF reference,
which require exact integrals.  This has been done recently by Weigend
and H\"aser\cite{Weigend:1997:TCA-97-331}, who have also reported a
preliminary implementation, as well as by
myself\cite{Bernholdt:to-be-published}.  Our implementation will be
based on NWChem's parallel semidirect MP2 energy and gradient
module\cite{Deegan:in-preparation}, and will be distinguished from the
work of Weigend and H\"aser by the ability to handle open-shell
systems via a UHF reference, and extremely large systems through a
parallel implementation.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{RI-MP2 Applications}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Because of the performance improvements they offer over the
traditional MP2 formulation, the RI- and OI-MP2 methods make feasible
calculations and entire studies which would otherwise be quite
impractical.  For instance, the RI-MP2 method has already demonstrated
an order of magnitude performance improvement in calculations of the
binding energy of K$^{+}$ to 12-crown-4 ether
\cite{Bernholdt:1996:CPL-250-447}.

In order to showcase the full capabilities of the RI-MP2 method on a
problem of environmental interest, I plan to extend the earlier
K$^{+}$:12-crown-4 calculations
\cite{Bernholdt:1996:CPL-250-447,Feller:1996:JCP-105-1940}.  Crown
ether systems are interesting because of their selective binding to
metal cations, and have been proposed for use in separating high-level
nuclear wastes such as $^{90}$Sr and $^{137}$Cs
\cite{Horwitz:1991:SEIE-9-1} among other applications.  In order to
design appropriately selective crown ether systems for specific
situations, it is necessary to have a thorough understanding of the
nature and origin of the selectivity.  Recently it was shown that
theory and experimental results for the gas-phase binding energy of
the K$^{+}$ and Cs$^{+}$:18-crown-6 complexes were in marked
disagreement
\cite{Katritzky:1990:JACS-112-2471,Katritzky:1992:RCMS-6-25,Glendenning:1994:JACS-116-10657}.
This kind of disagreement must be resolved if we are to have any hope
of achieving a detailed understanding of the nature of the
preferential binding and ultimately of designing new chemical systems
tailored to specific applications.  The preferred theoretical
treatment of these systems involves the use of progressively larger
basis sets, allowing extrapolation to the complete basis set (CBS)
limit for the level of theory applied.  Because of the size of these
systems, only relatively small basis set MP2 calculations can be
performed, making credible extrapolation difficult.  A comparison of
DFT, MP2 and RI-MP2 treatments of the analogous 12-crown-4 system
\cite{Feller:1996:JCP-105-1940} showed the MP2 approaches were
preferable to DFT, and that the RI-MP2 calculation was capable of
reproducing the binding energy in a small fraction of the time
required for exact MP2.

\begin{figure}[tbh]
\begin{minipage}[b]{3.5in}
%\vspace{2.5in}
%\psfig{file=k12c4_and_k18c6.eps,angle=00.0,width=3in}
\end{minipage}
\begin{minipage}[b]{3in}
\caption{Representative structures of bare 12-crown-4 and 18-crown-6
ethers and their complexes with K$^{+}$.}
\label{fig:crown-ethers}
\end{minipage}
\end{figure}

With a better understanding of the fitting basis and analytic
gradient, as proposed here, RI-MP2 will be in a position to replace
the traditional MP2, allowing investigations of these crown ether
systems to be expanded to larger basis sets and larger systems.  In
particular, I propose to use RI-MP2 calculations to complete the CBS
extrapolation of the structures and energetics of K$^{+}$:18-crown-6
(see Fig. \ref{fig:crown-ethers}).  Calculations on the Cs$^{+}$
system, utilizing relativistic effective core potentials, will be the
second phase, and extension to other cations is planned for the
future, as well as more complex crowns. This work will be coordinated
with related activities at EMSL to insure that this work extends
theirs.  The key advantage of using the RI-MP2 is that it will allow
calculations using larger basis sets, and on larger crowns, than are
possible with the traditional methods being used by EMSL researchers.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\subsection{OI-MP2 Applications}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Photoluminescent polymers are of great interest for use in organic
light emitting diodes and related applications.  In these materials,
the role of interchain interactions is of great importance.  A recent
paper by Blatchford and coworkers\cite{Blatchford:1996:PRB-54-9180}
reports for three pyridine-based polymers, the photoluminescence (PL)
spectrum of samples cast as films is strongly red-shifted ($\sim$0.5
eV) relative to the solution and powder form spectra and the films
show a much lower quantum yield.  These effects are attributed to
interchain interactions, which would be much stronger in the film than
in either the solution or powder sample, and for practical
applications of such polymers, which typically involve a solid-state
form, need to be controlled or eliminated.  Blatchford and coworkers
propose a mechanism involving the formation of \textit{aggregate
states} which tend to delocalize the ground and excited states across
several chains.  The energy of the ground state is lowered, and the
excited state is split into several distinct levels.  The OI-MP2
method is well-suited to a theoretical investigation of this effect,
as the systems involved are large, and the effect in question is large
enough that extremely high-level (i.e. large basis set) calculations
should not be required to model it.  We will consider various
oligomers of poly(\textit{p}-pyridine) and poly(\textit{p}-pyridyl
vinylene), and the effect of bringing two or more oligomers into
proximity on the ground and excited electronic states of the system.
An important question to be answered is how many chains must interact
in order to produce this effect, as the answer may contribute to the
development of ways to avoid the formation of aggregate states through
design and processing of the polymers.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Educational Aspects}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

The proposed work touches on a number of very significant areas of
correlated electronic structure methods and high-performance
computing.  It is rare to find candidates for postdoctoral research
positions with the particular combination of experience required for
this project; this is not a problem, but rather an excellent
opportunity for the researcher to extend their training into new
areas.  It is perhaps most likely that the researcher will have little
or no background in parallel computing.  Our experience in the High
Performance Computational Chemistry group at EMSL indicates that the
Global Array programming model is fairly easy to learn and use --
postdoctoral researchers in that group were quickly able to make
substantive contributions to the NWChem package.

In addition, my own unique position of being affiliated with both the
Northeast Parallel Architectures Center and the Chemistry Department
of Syracuse University (which are housed in the same building), will
provide the researcher with a flourishing interdisciplinary
environment including both chemists and computational scientists in
other fields.  In the course of my own work at NPAC, I have already
gotten a number of computer science graduate students and research
staff interested in problems on the computational side of
computational chemistry.  I also have active or nascent collaborations
with a number of the Chemistry faculty, coupling theoretical and
experimental work.

The nature of the proposed work and the environment here at Syracuse
provide an excellent opportunity for the postdoctoral researcher to
broaden their background and learn new skills, such as parallel
computing, which can be applied well outside the confines of the
proposed work.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Summary}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

This proposal describes a plan to develop the OI-MP2 and RI-MP2
methods, placing them on a solid footing as viable alternatives to
traditional MP2 calculations on large chemical systems.  The impact of
the work, however, will extend far beyond these two methods.  The
resolution of the identity approximation has already been applied to a
variety of methods, and is especially promising for higher-level
correlated methods.  The ideas of sparsity and locality associated
with orbital invariant formulations are also much more widely
applicable than just the MP2 method.  In addition, the work of this
proposal will be implemented for parallel computers, and will result
in tools (i.e. NWChem and the Global Array Toolkit), algorithms, and
practical experience that will benefit the entire field of
computational chemistry in the course of adapting to the latest
directions in high performance computing.  At the same time, this
proposal builds on my previous work in the areas of correlated
electronic structure methods and high-performance (parallel)
computing, and would provide a clear step forward in my overall
research goal of making large-scale correlated electronic structure
calculations more efficient and more accessible to the chemistry
community.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Bibliography is a separate chapter in the NSF style
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\renewcommand{\bibname}{Bibliography}
%\bibliographystyle{cv}
\bibliographystyle{plain}
\bibliography{\jobname}
%\bibliography{\jobname,new-research,vle2,bernholdt,nmr}
