NISTIR 7671 

Mathematical and Computational 

Sciences Division 

Summary of Activities for Fiscal Year 2009 

NISTIR 7671 

Mathematical and Computational 

Sciences Division 

Summary of Activities for Fiscal Year 2009 

Ronald F. Boisvert, Editor 

Mathematical and Computational Sciences Division 

Information Technology Laboratory 

February 25, 2010 

U.S. Department of Commerce 

Gary Locke, Secretary 

National Institute of Standards and Technology 

Patrick D. Gallagher, Director 
2 Summary of Activities for Fiscal Year 2009 

 
NIST Mathematical and Computational Sciences Division 3

Abstract

This report summarizes the technical work of the Mathematical and Computational Sciences Division 

(MCSD) of NIST’s Information Technology Laboratory. Part I (Overview) provides a high-level overview of the Division’s activities, including highlights of technical accomplishments during the previous 

year. Part II (Features) provides further details on ten projects of particular note this year. This is followed in Part III (Project Summaries) by brief synopses of all technical projects active during the past 

year. Part IV (Activity Data) provides listings of publications, technical talks, and other professional activities in which Division staff members have participated. The reporting period covered by this 

document is October 2008 through December 2009. 

For further information, contact Ronald F. Boisvert, Mail Stop 8910, NIST, Gaithersburg, MD 20899-

8910, phone 301-975-3812, email boisvert@nist.gov, or see the Division’s web site at 

http://www.nist.gov/itl/math/index.cfm. 

Cover Visualization: A “word cloud” summary of this report. See http://www.wordle.net/. 

Acknowledgement: We are grateful to Robin Bickel for collecting the information and organizing the 

first draft of this report. 

Disclaimer: Certain commercial entities, equipment, or materials may be identified in this document in 

order to describe an experimental procedure or concept adequately. Such identification is not intended to 

imply recommendation or endorsement by the National Institute of Standards and Technology, nor is it 

intended to imply that the entities, materials, or equipment are necessarily the best available for the purpose. 
4 Summary of Activities for Fiscal Year 2009 
NIST Mathematical and Computational Sciences Division 5

Contents

PART I: OVERVIEW......................................................................................................................... 9

Introduction............................................................................................................................... 11 

Highlights.................................................................................................................................. 12 

Technical Accomplishments .................................................................................................. 12

Technology Transfer and Professional Activities................................................................. 15

Staff News.............................................................................................................................. 16

Recognition ........................................................................................................................... 17

Restrictions on Transversal Encoded Quantum Gate Sets........................................................ 23 

High Precision Calculations of Fundamental Properties of Few-electron Systems.................. 25 

Characterization of Electromagnetic Metafilms ....................................................................... 27 

Recovering from Physically Impossible Reconstructions in Blind Deconvolution.................. 29 

A Quality Pre-Processor for Biological Cell Images................................................................ 31 

Computational Modeling of the Flow of Dense Suspensions................................................... 33 

Visualization of Radio Frequency Propagation for Body Area Networks................................ 35 

Mathematical Knowledge Management ................................................................................... 37 

PART III: PROJECT SUMMARIES.................................................................................................. 39

 Information Technology Laboratory Programs ...................................................................... 41 

Enabling Scientific Discovery ............................................................................................... 42

Recovering from Physically Impossible Reconstructions in Variational Blind 

Deconvolution................................................................................................................... 42 

Computational Modeling of the Flow of Dense Suspensions........................................... 42 

Sparse Representations in High Dimensional Geometry.................................................. 42 

Analysis of Diffusion and Stress in Binary Alloys........................................................... 43 

Modeling of Self-Similar Grain Size Distributions .......................................................... 44 

Modeling Fluid Flow in Materials Processing.................................................................. 44 

Materials Data and Metrology for Applications to Machining Processes, Frangible 

Ammunition, and Body Armor ......................................................................................... 45 

Boundary Slip Effects on the Stability of Spiral Poiseuille Flow..................................... 46 

Nanostructures, Nanooptics, and How to Squeeze the Light out of Quantum Dots......... 46 

Computational Modeling and Visualization of Cement Paste Hydration and 

Microstructure Development ............................................................................................ 48 

Micromagnetic Modeling.................................................................................................. 49 

OOF: Finite Element Analysis of Material Microstructures............................................. 50 

Fast Solvers for Models of Fluid Flow ............................................................................. 51 

Parallel Adaptive Refinement and Multigrid Finite Element Methods ............................ 51 
6 Summary of Activities for Fiscal Year 2009 

Information Discovery, Use and Sharing ..............................................................................53

Mathematical Knowledge Management ............................................................................53 

A Quality Pre-Processor for Biological Cell Images.........................................................53 

Estimating Volumes of Simulated Lung Cancer Nodules.................................................53 

Standardization of Lung Cancer Growth Measurements...................................................54 

Complex Systems....................................................................................................................56

Methods for Characterizing Complex Networks...............................................................56 

Standard Reference Data for Complex Network Research ...............................................56 

Community Detection in Network Graphs........................................................................58 

Visualization of Large Complex Networks .......................................................................58 

Fluid-level Modeling of TCP.............................................................................................59 

Analysis of a Distributed Protocol for Network Control...................................................59 

Visualization of Possible Route Probability Distributions................................................60 

Grid Network Simulation Using Piece-Wise Homogeneous Markov Chains...................61 

Quantum Information ............................................................................................................62

Restrictions on Transversal Encoded Quantum Gate Sets ................................................62 

Quantum Computing Theory.............................................................................................62 

Fault-tolerant Quantum Computation................................................................................63 

Quantum Computing Benchmarks ....................................................................................64 

Quantum State Estimation Theory and Applications ........................................................64 

Virtual Measurement Systems................................................................................................66

High Precision Calculations of Fundamental Properties of Few-electron Systems ..........66 

Micromagnetic Modeling ..................................................................................................66 

Uncertainty in Virtual Measurements from Quantum Chemistry Models ........................67 

Pervasive Information Technology........................................................................................68

Visualization of Radio Frequency Propagation for Body Area Networks........................68 

Self-organization of Mobile Sensor Networks..................................................................68 

Trustworthy Information Systems ..........................................................................................70

Automated Combinatorial Testing for Software-systems (ACTS) ...................................70 

Division Programs .....................................................................................................................73 

Digital Library of Mathematical Functions...........................................................................74

Digital Library of Mathematical Functions .......................................................................74 

Visualization of Complex Function Data ..........................................................................75 

Mathematics of Metrology.....................................................................................................77

Multihazard Structural Optimization.................................................................................77 

A Mathematical Model of Atmospheric Retention of Man-made CO2 Emissions ...........77 

Simulation of Bioregulatory Networks Involved in Cell Cycle Control ...........................79 

Modeling the Fluid Dynamics of Tear Films on the Cornea.............................................79 

Modeling Equilibrium, Stability and Transport during Magnetic Fusion .........................80 

Modeling Type Ia Supernova Light Curves ......................................................................81 

High Performance Computing and Visualization..................................................................83

Visualizing Fire and the Resulting Deformation of Structures .........................................83 

Fundamental Math Software Development and Testing........................................................85

SciMark, a Web-based Benchmark for Numerical Computing in Java.............................85 

Division Software and Services.........................................................................................85 
NIST Mathematical and Computational Sciences Division 7

PART IV: ACTIVITY DATA ........................................................................................................... 87

Publications............................................................................................................................... 89 

Appeared ............................................................................................................................... 89

Refereed Journals.............................................................................................................. 89 

Journal of Research of NIST............................................................................................. 90 

Books ................................................................................................................................ 90 

Book Chapters................................................................................................................... 90 

Conference Proceedings.................................................................................................... 90 

Other Formal Publications ................................................................................................ 92 

Technical Reports ............................................................................................................. 92 

Visualizations Published ................................................................................................... 92 

Accepted ................................................................................................................................ 92

In Review............................................................................................................................... 93

Presentations ............................................................................................................................. 94 

Invited Talks .......................................................................................................................... 94

Conference Presentations ..................................................................................................... 95

Patents ....................................................................................................................................... 97 

Software Released..................................................................................................................... 97 

Conferences, Minisymposia, Lecture Series, Shortcourses...................................................... 98 

MCSD Seminar Series........................................................................................................... 98

Local Events Organized ........................................................................................................ 98

External Events Organization............................................................................................... 98

Other Professional Activities .................................................................................................... 99 

Internal.................................................................................................................................. 99

External................................................................................................................................. 99

External Contacts .................................................................................................................... 100 

PART V: APPENDICES ................................................................................................................. 103

Staff......................................................................................................................................... 105 

Glossary of Acronyms ............................................................................................................ 108 
8 Summary of Activities for Fiscal Year 2009 
NIST Mathematical and Computational Sciences Division 9

Part I: Overview

 Part I 

 Overview 
10 Summary of Activities for Fiscal Year 2009 
NIST Mathematical and Computational Sciences Division 11

Introduction

Founded in 1901, NIST is a non-regulatory federal agency within the U.S. Department of Commerce. 

NIST's mission is to promote U.S. innovation and industrial competitiveness by advancing measurement 

science, standards, and technology in ways that enhance economic security and improve our quality of 

life. The NIST Laboratory program is broad-ranging, with research efforts in physics, electronics and 

electrical engineering, chemical science and technology, materials science, building and fire research, 

manufacturing engineering, and information technology. 

The Mathematical and Computational Sciences Division (MCSD) provides leadership within NIST in 

the application of applied and computational mathematics to the solution to science and engineering problems arising in measurement science. In that role Division staff members 

• perform research and development in applied mathematics and computer science and engineering, including analytical methods, numerical algorithms, advanced computing and 

communications architectures and applications, and high performance visualization; 

• engage in peer-to-peer collaborations in the application of mathematical and computational 

technologies to NIST problems; 

• develop and disseminate mathematical reference data, software, and related tools; and 

• work with internal and external groups to develop standards, test procedures, reference implementations, and other measurement technologies for advanced scientific computation on 

current and future architectures. 

Division staff is organized into four groups: 

Mathematical Modeling Group (Geoffrey McFadden, Leader) 

Performs research and maintains expertise in applied mathematics, mathematical modeling, numerical analysis, and scientific computing for application to multidisciplinary problems of interest 

to the NIST measurement science program. 

Mathematical Software Group (Daniel Lozier, Leader) 

Performs research and maintains expertise in the methodology and application of mathematical algorithms and software in support of computational science within NIST as well as in industry and 

academia. 

Discrete Mathematical Analysis Group (Ronald Boisvert, Acting Leader) 

Performs research and maintains expertise in the theory, methodology, and application of discrete 

mathematics for application to the NIST measurement science program, including complex systems and quantum information science. 

Scientific Applications and Visualization Group (Judith Terrill, Leader) 

Develops novel algorithms and implementations for parallel and distributed computation and visual 

data analysis in a variety of hardware environments. 

MCSD is one of six technical Divisions within the NIST Information Technology Laboratory 

(ITL). The work of ITL is organized into 11 cross-cutting programmatic areas, seven of which have significant MCSD participation. See “ITL Programs in Brief” sidebar. MCSD Division-centric work is 

organized into four additional themes. See the “MCSD Programs in Brief” sidebar. Division programs are 

characterized as mature Division-centric efforts in support of NIST goals that do not fit in any of ITL’s 

Program areas, but are nevertheless critical for NIST metrology. Other Division efforts fulfill a need to 

continuously develop and maintain technical competence and facilities in core areas of mathematics and 

computational science which are broadly applicable. 

In Part III of this document we provide descriptions of most ongoing technical projects of the Division. These are organized under the Program areas listed in the sidebars. 
12 Summary of Activities for Fiscal Year 2009 

Highlights

In this section we identify some of the major accomplishments of the Division during the past year. We 

also provide news related to MCSD staff. 

Technical Accomplishments 

MCSD has made significant technical progress on many fronts during the past year. Here we highlight a 

few examples. Further details are provided in Part II (Features) and Part III (Project Summaries) of this 

document. 

ITL Programs In Brief 

* Complex Systems. Complex Systems are composed of large 

interrelated, interacting entities which taken together, exhibit a 

macroscopic behavior which is not predictable by examining the 

components. This Program seeks to understand the fundamental 

science of these systems and develop rigorous descriptions that 

enable prediction and control of their behavior. 

Cyber and Network Security. Cybersecurity is focused on ensuring 

three IT system objectives: confidentiality, integrity, and availability. This Program addresses long-term scientific issues in some of 

the building blocks of IT security ― cryptography, security testing 

and evaluation, security metrics, and security properties ― providing a more scientific foundation for computer and network 

security, while maintaining a focus on near-term issues.

* Enabling Scientific Discovery. Modern scientific research has 

become more and more dependent on mathematical, statistical, and 

computational tools. This Program promotes the use of these tools 

to dramatically advance our ability to predict the behavior of a 

broad range of complex scientific and engineering systems and 

processes. In particular, this Program focuses on inter-disciplinary 

projects that involve novel computational statistics and simulation 

methods and software. 

Identity Management Systems. Identity management systems are 

responsible for the creation, use, and termination of electronic 

identities used to access logical and physical resources. This Program is pursuing the development of common models and metrics 

for identity management, critical standards, and interoperability of 

electronic identities. These efforts will improve the quality, usability, and consistency of identity management systems while 

protecting privacy. 

* Information Discovery, Use, and Sharing. This Program fosters 

innovation throughout the information life cycle by developing the 

measurement infrastructure to enhance knowledge discovery, 

information exchange, and information usability. The Program 

enables novel computational approaches to data collection and 

analysis to be combined with improved interoperability techniques 

to effectively extract needed information from the wealth of available data. 

* Pervasive Information Technologies. Pervasive information 

technology is the trend towards increasingly ubiquitous connected 

computing sensors, devices, and networks that monitor and respond transparently to human needs. This Program facilitates the 

creation of standards for sensor communication, networking interoperability, and sensor information security. 

* Trustworthy Information Systems. Trustworthy information 

systems perform as intended for a specific purpose, with operational resiliency and without unwanted side effects or exploitable 

vulnerabilities. This Program will improve the ability to model, 

produce, measure, and assess trustworthiness through innovative 

technologies, models, measurements, and software tools. 

* Virtual Measurement Systems. A virtual measurement is a 

quantitative result and its uncertainty, obtained by a computer 

simulation or computer-assisted measurements. This Program 

introduces metrology constructs, standard references, uncertainty 

characterization, and traceability into scientific computation and 

computer-assisted measurements to enable predictive computing 

with quantified reliability. 

* Quantum Information. An emerging new discipline at the intersection of quantum mechanics and computer science, quantum 

information science (QIS) is likely to revolutionize science and 

technology in the same way that lasers, electronics, and computing 

did in the 20th century, including phenomenal increases in information storage and processing speed and communication channels 

with extremely high levels of security. However, the needed manipulation and control of quantum states remains a grand 

challenge. We are developing theoretical underpinning to support 

the multi-laboratory experimental program in QIS at NIST. 

Health Information Technology. The demand for online medical 

information and simplified, standardized methods to access healthcare information and services is crucial in making healthcare safe 

and available to all. NIST is collaborating with industry, standards 

organizations, consortia, and government agencies to build tools 

and prototypes to advance the adoption of IT within healthcare 

systems. 

Voting. The 2002 Help America Vote Act gives NIST a key role 

in realizing nationwide improvements in voting systems. To assist 

the Election Assistance Commission with the development of 

voluntary voting system guidelines, HAVA established the Technical Guidelines Development Committee (TGDC) and directs 

NIST to chair the TGDC. NIST research activities include: security 

of computers, computer networks, and computer data storage used 

in voting systems; methods to detect and prevent fraud; protection 

of voter privacy; and the role of human factors and assistive technologies in the design and application of voting systems. 

* indicates MCSD participation 
NIST Mathematical and Computational Sciences Division 13

 Image analysis continues to be a strong emphasis of Division activity. Imagery requiring enhancement and analysis in NIST applications arise in nanotechnology, medicine, biology, and materials 

science. For example, we are working with researchers in the NIST Chemical Sciences and Technology 

Laboratory to develop effective automated analysis of cell imagery. In collaboration with Karen Kafadar 

of Indiana University and Alden Dima of ITL, Adele Peskin has developed a quality metric for cell imagery. The metric is a useful predictor of which of the many available segmentation algorithms will be 

most effective in isolating an individual cell within an image. [Page 31] Images also play an important 

role in the analysis of complex materials. MCSD’s OOF system starts with a micrograph of a complex 

material from which it develops a computational grid to enable a finite element simulation from which 

macroscopic properties can be derived. This year substantial progress has been made in extending OOF’s 

capabilities to three-dimensional image-based analysis. [Page 50] We also take the lead in the development of new image analysis methods for use in metrology. This year Alfred Carasso completed an 

insightful mathematical analysis of a recently proposed blind deconvolution algorithm, demonstrating that 

the visually compelling reconstructions that it produced were actually inadmissible, corresponding to 

physically impossible blurs. He identified the source of the problem and developed a correction procedure. The result has proven to be highly effective. [Page 29] 

Division research in quantum information continues to contribute significantly to the NIST program 

in this area. Two such projects are externally funded efforts: Quantum Optical Metrology with N-Photons 

(joint with EEEL and PL) and Analog Quantum Integrated Circuits (joint with PL). One important contribution we are making to these efforts is the development of novel quantum state estimation procedures 

and associated software. [Page 64] Our efforts to improve our understanding of the potential for the development of quantum computing devices also continue. This year we settled a long-standing conjecture 

regarding the design of quantum gates sets. Bryan Eastin and Manny Knill showed that transversal encoded gate sets, which are ideal for error management, cannot be universal, and hence cannot solely be 

relied upon for the practical implementation of quantum computers. [Page 21] 

We continue our efforts, begun in 2007, to develop foundations of measurement science for information systems. A significant focus of this work is the development and analysis of metrics for 

MCSD Programs In Brief 

Digital Library of Mathematical Functions. The special functions 

of applied mathematics are fundamental tools enabling modeling 

and analysis in all areas of science and engineering. To make effective use of such functions, practitioners must have ready access 

to a reliable source of information on their properties. The goal of 

this work is the development and dissemination of definitive reference data on the special functions of applied mathematics. The 

DLMF will be a freely available interactive and richly linked 

online resource.

Mathematics of Metrology. Mathematics plays an important role 

in the science of metrology. Mathematical models are needed to 

understand how to design effective measurement systems, and to 

analyze the results they produce. Mathematical techniques are used 

to develop and analyze idealized models of physical phenomena to 

be measured, and mathematical algorithms are necessary to find 

optimal system parameters. Finally, mathematical and statistical 

techniques are needed to transform the resulting data into useful 

information. The goal of this work is to develop fundamental 

mathematical methods and analytical tools necessary for NIST to 

continue as a world-class metrology institute, and to apply them 

critical measurement science applications. 

Fundamental Math Software Development, Testing. Modern 

science and engineering in general, and modern measurement 

science in particular require a wide variety of software tools for 

scientific discovery, exploration, and analysis. As scientific inquiry 

becomes deeper and more specialized, so must the supporting 

software tools. The goal of this work is to develop critical software 

tools that support measurement science at NIST, as well as computational science and engineering at large.

High Performance Computing and Visualization. Computational 

capability is advancing rapidly. This means that modeling and 

simulation can be done with greatly increased fidelity (e.g. higher 

resolution, more complex physics). However, developing largescale parallel applications remains highly challenging, requiring 

expertise that application scientists rarely have.In addition, the 

hardware landscape is changing rapidly so new algorithmic techniques must constantly be developed. We are developing and 

applying facilities and expertise of this type for application to 

NIST problems. Large scale computations and laboratory experiments invariably produce large volumes of scientific data, which 

cannot be readily comprehended without some form of visual 

analysis. We are developing the infrastructure necessary for advanced visualization of scientific data, including the use of 3D 

immersive environments and applying this to NIST problems. One 

of our goals is to develop the 3D immersive environment into a 

true interactive measurement laboratory. 

 
14 Summary of Activities for Fiscal Year 2009 

characterizing complex networks. One such example is the reliability polynomial, which gives the probability that a network will remain connected given that each edge can fail independently with a 

probability p. The reliability polynomial of a connected graph G with m edges is an m

th degree polynomial 

in p whose j

th coefficient is the number of distinct connected subgraphs of G with exactly j edges. Unfortunately, determining these coefficients is computationally intractable. (In particular, the problem is #Pcomplete.) Thus, the best we can expect to do in general is to approximate the polynomial. Brian 

Cloteaux and Isabel Beichl have developed a new method to do this. Their approach has two advantages. 

First, while there are several algorithms that provide pointwise approximations of the polynomial, their 

approach approximates the polynomial coefficients themselves. Secondly, their method has shown much 

faster rate of convergence than previously proposed methods. A paper has been accepted for publication. 

[Page 56] 

This year also witnessed a significant milestone in our efforts to develop virtual measurement systems 

for fundamental atomic properties. Our recent computations of the nonrelativistic energies for excited 

states of the lithium atom represent the highest level of accuracy ever reached (less than 10-9 hartree) in 

atomic quantum computations with more than two electrons. The achievement was reported in an article 

in Physical Review A by James Sims of MCSD and Stanley Hagstrom of Indiana University. The judicious use of both parallel computing and multiple precision arithmetic (30+ digits) was critical to obtain 

such high accuracy results. In addition, Sims and Hagstrom developed an improved parallel computer 

code for a key computational bottleneck, high-precision solution of the large-scale generalized matrix 

eigenvalue problem. For a 17,190 term wave function they achieved a factor of 30 speedup on 32 processors of a local Linux cluster. The resulting computations obtained the values of four excited states of the 

lithium atom to an accuracy two orders of magnitude greater than has been done before. While the lithium 

calculation is valuable in itself, the deeper import of refining their method is that it should enable the calculation of excitation energies for beryllium, which has four electrons. This achievement would enable 

theorists to predict with greater accuracy values for all of the remaining elements in the second row of the 

periodic table. [Page 25] 

Our efforts with the NIST Building and Fire Research Laboratory to develop a virtual cement and 

concrete testing laboratory reached a new level of maturity this year. Our simulations have enabled the 

exploration of the effects of particle size distribution, particle shape distribution, volume fraction of aggregates, and other relevant system parameters, on critical properties such as the yield stress and viscosity 

of suspensions. In order to perform simulations on systems of sufficient size, density, and simulated time, 

we have developed highly parallel implementation of our models. This implementation, which is based 

upon dissipative particle dynamics, has run simulations efficiently using more than 2,000 processors. In 

these cases we often simulate systems consisting of over 1 million mesoscopic fluid particles and up to 

10,000 aggregates. These computations were enabled by DOE INCITE awards of 750,000 hours of CPU 

time on the IBM Blue Gene/P machine at the Argonne National Laboratory Leadership Computing Facility in each of 2008 and 2009. For 2010 we have been awarded an additional 2 million CPU-hours of 

compute time. In 2009 William George and Judith Terrill of MCSD were part of a team awarded the 

DOC Silver Medal in recognition of this work. [Page 17] 

Visual analysis plays an important role in our modeling of suspensions. Such analysis is supported by 

our immersive visualization facilities in both Gaithersburg and Boulder. Division staff members work 

closely with scientists from a variety of NIST Labs to develop novel visualizations that enable insights 

that would be difficult to obtain otherwise. Examples of projects that have benefited from this facility include ITL efforts to model radio frequency propagation in body area networks [Page 35], work with PL to 

model the optical properties of nanostructures [Page 46], ITL efforts to assess community detection algorithms for large-scale graphs [Page 56], and work with BFRL on modeling the deformation of structures 

resulting from fire [Page 83]. Our visualization team is continuously updating our hardware and software 

capabilities. In particular, we have developed an extensive suite of visualization primitives and tools that 

allow us to be both fast and flexible in our work. These tools run without change on a range of platforms 

from a laptop to our immersive environment. A recent focus has been the use of specialized graphics 

processing unit (GPU) processing to provide novel capabilities for our most challenging visualization 

problems. 
NIST Mathematical and Computational Sciences Division 15

 Finally, we are happy to report that the release of the long awaited Digital Library of Mathematical 

Functions (DLMF) is in sight. A modern successor to the classic Handbook of Mathematical Functions

(M. Abramowitz and I. Stegun, Eds., 1964), the DLMF is a carefully constructed reference on the properties of the special functions of applied mathematics which has been compiled with the cooperation of 

some 50 researchers worldwide. Release is now scheduled for the spring of 2010. The information will be 

made available as both an innovative interactive web resource as well as a traditional printed book produced in cooperation with Cambridge University Press. [Page 74] 

Technology Transfer and Professional Activities 

The volume of technical output of MCSD remains high. During the last 18 months, Division staff members were (co-)authors of 24 articles appearing in peer-reviewed journals, two books, two book chapters, 

26 papers in conference proceedings, and 7 other formal publications. Twenty-five additional papers have 

been accepted for publication, while 19 others are undergoing review. Division staff members gave 27 

invited technical talks and presented 40 others in conferences and workshops. 

MCSD continues to maintain an active Web site with a variety of information and services, including the Guide to Available Mathematical Software, the Matrix Market, and the SciMark Java benchmark. 

During calendar year 2009, the division web server satisfied more than 7 million requests for pages from 

more than 445,000 distinct hosts. On a daily basis, more than 19,000 pages and 1.6 Gbytes of data were 

shipped. In total, there have been more than 217 million “hits” on MCSD Web servers since they went 

online as NIST’s first web servers in 1994. The individual software packages that we produce and distribute also continue to see very high usage. Table 1 lists download counts for some of our most popular 

ones. Another indication of the successful transfer of our technology is references to our software in refereed journal articles. For example, our OOMMF software for nano-magnetic modeling was cited in 120 

such papers which were published in FY 2009 alone (836 such papers have been identified since 2001). 

Members of the Division are also active in professional circles. Isabel Beichl began a term as Editor-in-Chief of Computing in Science & Engineering in January 2009. Staff members hold a total of six 

additional associate editorships of peer-reviewed journals. They are also active in conference organization, serving on five organizing/steering/program committees. Staff members organized five 

minisymposia at technical conferences. 

Service within professional societies is also prevalent among our staff. For example, Ronald Boisvert serves as Co-Chair of the Publications Board of the Association for Computing Machinery (ACM) 

and is a member of the ACM Council, the association’s board of directors. Daniel Lozier serves as Secretary of the Society for Industrial and Applied Mathematics (SIAM) Activity Group on Orthogonal 

Polynomials and Special Functions. Staff members are also active in a variety of working groups. For 

example, Ronald Boisvert serves as Chair of the International Federation for Information Processing 

(IFIP) Working Group 2.5 on Numerical Software, Donald Porter is a member of the Tcl Core Team, and 

Bruce Miller is a member of W3C’s Math Working Group. Judith Terrill represents NIST on the High 

Table 1. Downloads of Selected Division Software Packages During CY 2009

Package Description Downloads 

JAMA Linear algebra in Java 17,067 

TNT Linear algebra using C++ templates 15,312 

LAPACK++ Linear algebra in C++ 6,307 

OOMMF Modeling of nano-magnetic phenomena 5,200 

SparseLib++ Sparse linear algebra in C++ 4,903 

f90gl Fortran 90 interface to OpenGL graphics 3,530 

REGEN Modeling of cryocoolers 1,932 

spblas Sparse basic linear algebra components 1,569 

LaTeXML LaTeX to MathML translator 1,520 

OOF Modeling materials with complex microstructure 1,204 

PHAML Parallel adaptive solution of partial differential equations 767 
16 Summary of Activities for Fiscal Year 2009 

End Computing Interagency Working Group of the Federal Networking and Information Technology Research and Development (NITRD) Program. 

For further details, see Part IV (Activity Data) of this document. 

Staff News 

During FY 2009 MCSD welcomed a new National Research Council (NRC) Postdoctoral Associate, 

James Benson of the University of Missouri – Columbia. A recent PhD in Mathematics, Benson’s thesis 

was in the mathematics of cryobiology, the study of biological systems at low temperatures. In particular, 

Table 2. Student Interns in MCSD

Name Institution Level Program Mentor Project Title 

Matthew 

Du 

Thomas S. 

Wooton 

High School 

High 

School 

Student 

Volunteer 

J. Terrill Modeling Window Breakage due to 

Fire using the Open Dynamics Engine 

Julie 

Eisenberg 

Damascus 

High School 

High 

School 

Student 

Volunteer 

J. Terrill Merging Data Mining with Visualization 

Daniel 

Leifer 

Stanford 

Univ. 

Undergraduate SURF A. Peskin Development of Standardized Lung 

Tumor Data Sets 

Sei 

Masouka 

Winston 

Churchill 

High School 

High 

School 

Student 

Volunteer 

J. Terrill Influence of Particle Size Distribution on Setting of Tricalcium Silicate 

Russell 

McFadden 

St. Olaf College Undergraduate Student 

Volunteer 

W. Mitchell Mathematical Modeling and Software 

Adam 

Meier 

University of 

Colorado 

Graduate PREP E. Knill Quantum Computer Simulation 

Menal 

Modha 

UMBC Undergraduate SURF R. Kacker/ 

R. Kuhn 

Isolating Failure-causing Combinations in Software Discovered via 

Combinatorial Testing

Katherine 

Morrison 

Univ. of 

Nebraska 

Graduate Guest Researcher F. Hunt Modeling the Operation of a Grid 

Computing System 

Michael 

Mullan 

University of 

Colorado 

Graduate PREP E. Knill Quantum Computing Theory 

Sanjeev 

Muralidharan 

Montgomery 

Blair High 

School 

High 

School 

Student 

Volunteer 

J. Terrill Visualization and Measurement in 

Dental Images 

Michael 

Sharpnack 

NYU Undergraduate SURF G. McFadden Sensitivity Analysis for a Biological 

Model of Cellular Response to DNS 

Damage

Nadezhda 

Serova 

UMBC Undergraduate SURF B. Cloteaux Approximation Schemes for the 

Maximal Cover Sum Problem 

Joseph 

Owen 

Montgomery 

Blair High 

School 

High 

School 

Student 

Volunteer 

J. Terrill Language for Generating Geometry 

for the Case of Room and Contents 

in the Project “Impact of Fire on 

Building Structures” 

Jennifer 

Schiller 

Reservoir 

High School 

High 

School 

Student 

Volunteer 

J. Terrill A Mathematical Metric for Stress 

Groupings in Flow of Suspensions 

Yanbao 

Zhang 

University of 

Colorado 

Graduate PREP E. Knill Tests of Local Realism 

SURF: NIST Student Undergraduate Student Fellowship Program. 

 PREP: NIST Professional Research Experience Program
NIST Mathematical and Computational Sciences Division 17

he developed methods for optimizing cryobiological protocols, and in the process defined a new model 

for common cryobiological procedures. At NIST Benson will be working with Anthony Kearsley. Benson’s arrival brings our current number of NRC Postdoctoral Associates to four. 

During FY 2009 MCSD was able to support the work of 15 student interns, including four in the 

Student Undergraduate Research Fellowship (SURF) program, three graduate students in Boulder’s Professional Research Experience Program (PREP), one graduate student guest researcher, and seven high 

school volunteers. See Table 2 for a complete list. 

Recognition 

Division staff garnered a significant number of professional recognitions during the past year. These are 

described below. 

External Awards 

o American Physical Society Outstanding Referee. Manny Knill of MCSD has 

been cited as an Outstanding Referee by the American Physical Society (APS). 

This lifetime APS award was initiated in 2008 “to recognize scientists who have 

been exceptionally helpful in assessing manuscripts for publication in the APS 

journals.” The program is designed to annually recognize about 130 of APS’s approximately 42,000 referees, though larger numbers have been identified during 

2008 and 2009, the program’s first two “catch-up years.” The basis for selection was the quality, 

number and timeliness of their reports, without regard for membership in the APS, country of origin, 

or field of research. 

o SIAM Fellows. Dianne O’Leary and G. W. (Pete) Stewart, Faculty 

Appointees in MCSD, have been named Fellows of the Society for Industrial and Applied Mathematics (SIAM). O’Leary, who is Professor 

of Computer Science at the University of Maryland College Park 

(UMCP), was cited for contributions to linear algebra, regularization, and applications. Stewart, Distinguished University Professor Emeritus of Computer Science at UMCP, was cited for contributions 

to numerical linear algebra. The SIAM Fellows program has just been established; O’Leary and 

Stewart are members of the inaugural class of inductees.

o Excellence in Technology Transfer. Raghu Kacker of MCSD, Rick Kuhn of the ITL Computer Security Division, Jeff Yu Lei of the University of Texas at Arlington (and a guest researcher in MCSD), 

and James Lawrence of George Mason University (and a faculty appointee in MCSD) were selected 

as recipients of the 2009 Federal Laboratory Consortium-Mid Atlantic Regional “Excellence in Technology Transfer” Award for their work in “Advanced 

Combinatorial Test Suites (ACTS) for Testing Software.” This prestigious award recognizes employees 

who have accomplished outstanding work in the process of transferring a technology developed by a Federal 

laboratory to the commercial marketplace. 

Internal Awards 

o Silver Medal. William George and Judith Terrill of 

MCSD are co-winners of a 2009 Silver Medal from the 

Department of Commerce for their work on the modeling and simulation of complex suspensions. The Silver 

Medal recognizes “exceptional performance characterized by noteworthy or superlative contributions which 

have a direct and lasting impact within the DepartWilliam George and Judith Terrill received the Department of Commerce’s Silver Medal 
18 Summary of Activities for Fiscal Year 2009 

ment.” Dale Bentz, Jeffrey Bullard, Edward Garboczi, and Nicos 

Martys of BFRL were co-winners. The award citation credits the 

group wih the creation of the virtual cement and concrete testing 

laboratory (VCCTL), “which uses intensive computer simulation to 

mimic a physical testing laboratory and has proven to be capable of 

revealing behavior beyond the limits of traditional testing instruments.” 

o Allen V. Astin Measurement Science Award. Andrew Dientsfrey of 

MCSD is a co-winner of the 2009 Allen V. Astin Measurement Science Award for work in waveform metrology. The award was 

instituted to recognize “outstanding achievement in the advancement 

of measurement science or in the delivery of measurement services 

developed within NIST, in cooperation with other Government 

agencies or private groups, or in cooperation with international metrology organizations.” Paul Hale and Dylan Williams of EEEL and 

Jack Wang of ITL’s Statistical Engineering Division were named as 

co-winners. The team was cited for “providing a world-first method 

for simultaneously calibrating high-speed electrical test equipment 

in both time and frequency domains with point-by-point uncertainty 

analysis, and for the first time making waveform calibrations traceable to fundamental 

physics.” 

o ITL Outstanding Journal Paper. Bert Rust and Dianne 

O’Leary were honored with the 2008 ITL Outstanding 

Journal Paper Award for their article “Residual periodograms for choosing regularization parameters for illposed problems.” The paper appeared in the journal 

Inverse Problems, Volume 24, Article Number 034055 

(May 2008). The award citation reads "For notable scientific research and exposition in the field of inverse 

problems as applied to measurement science.”

Best Contribution Awards 

o CIRP Conference on Modeling of Machining Operations. A paper presented by Timothy Burns of MCSD 

at the 12th CIRP (The International Academy for Production Engineering) Conference 

on Modeling of Machining Operations, in San Sebastian, Spain, 

May 7-8, 2009 was one of a few 

submissions selected for publication in the international journal Machining Science and Technology. The paper, 

“Dynamic properties for modeling and simulation of machining: An update of results from the NIST pulse-heated Kolsky bar,” discussed some significant implications for finiteelement simulations of new NIST material response measurements on stress vs. strain response under 

conditions of extremely rapid heating of carbon steel in high-speed machining. Co-authors on the paper were Steven Mates of NIST/MSEL, Richard Rhorer and Eric Whitenton of NIST/MEL, and 

Debasis Basak, of Orbital Sciences Corporation. 

Andrew Dienstfrey was a co-winner of 

NIST’s Allen V. Astin Measurement 

Science Award. 

Bert Rust and Dianne O’Leary received ITL’s Outstanding Journal Paper Award for 2009. 
NIST Mathematical and Computational Sciences Division 19

o IEEE Personal, Indoor and Mobile Radio Communications Symposium. A paper co-authored by John 

Hagedorn and Judith Terrill of MCSD received one of the three Best Paper awards at the 20th 

Personal, Indoor and Mobile Radio Communications Symposium (PIMRC) held in Tokyo, Japan on 

September 13-16, 2009. PMIRC is one of the major IEEE Communications Society conferences; 

some 650 papers were presented there this year. The cited paper, 

“A Statistical Path Loss Model for Medical Implant Communication Channels,” describes a new model to characterize RF path 

loss from wireless medical implants. The paper was presented by 

co-author Kamran Sayrafian, ITL Pervasive Information Technology Program Manager. Wen-Bin Yang of the ITL Advanced Networking Technologies Division 

and Kamya Yekeh Yazdandoost of the National Institute of Information and Communications Technology in Yokosuka, Japan were also co-authors. 

o Scientific Discovery through Advanced Computing (SIDAC) Conference. A short movie entitled 

“Simulation of Non-Newtonian Suspensions: Shear Thinning Case” developed by MCSD’s Scientific 

Applications and Visualization Group in collaboration with NIST/BFRL scientists was chosen to receive an OASCR Award at the 2009 Scientific Discovery through Advanced Computing (SciDAC) 

Conference held in San Diego on June 14-18, 

2009. The SciDAC Conference brings together 

about 350 scientists participating in DOE’s SciDAC Program along with other prominent 

researchers from the computational science community. Named for DOE's Office of Advanced Scientific Computing Research, OASCRs were awarded to 10 entries at the conference's Electronic Poster 

and Visualization Night session. The NIST entry illustrates how suspensions such as concrete or paint 

react as strain is applied. The movie was produced by William George, Steven Satterfield, John 

Hagedorn, Marc Olano and Judith Terrill of MCSD, along with Nicos Martys of BFRL's Materials 

and Construction Research Division. This is the second OASCR award garnered by the NIST team. 

Other Significant Recognition 

o Beichl Named Editor-in-Chief of Scientific Computing Magazine. Isabel Beichl of MCSD was selected to serve as Editor-in-Chief of Computing in Science and Engineering (CiSE), a joint 

publication of the IEEE Computer Society and the American Institute of Physics. CiSE is a bimonthly publication that combines peer-reviewed articles, 

columns, tutorials, and book reviews in a magazine format. 

According to its statement of purpose, CiSE “aims to support 

and promote the emerging discipline of computational science and engineering and to foster the use of computers and 

computational techniques in scientific research and education.” Beichl’s appointment as Editor-inChief, which was made by the IEEE Computer Society Board of Governors, is for a two-year term 

which began in January 2009. 
20 Summary of Activities for Fiscal Year 2009 
NIST Mathematical and Computational Sciences Division 21

 Part II 

 Features 
22 Summary of Activities for Fiscal Year 2009 
NIST Mathematical and Computational Sciences Division 23

Restrictions on Transversal Encoded Quantum Gate Sets 

During the last few decades, the theory of computation 

has gradually been expanded to encompass the possibility of quantum-mechanical effects in information 

processing. Quantum computers, as machines capable 

of exploiting such effects are called, have sparked considerable interest because it appears that they should 

be able to efficiently solve certain classically intractable problems. A working quantum computer has yet to 

be built, however, due to the extreme difficulty of precisely controlling complex quantum systems. 

The accuracy required to implement many interesting quantum algorithms is daunting, but this 

accuracy need not be achieved directly in the hardware. As the theory of quantum computation has 

matured, a variety of software approaches have been 

developed to mitigate the effects of hardware imperfections. The foremost of these is quantum error 

correction, which, as in the classical case, relies on 

encoding the data of interest using a suitable code. 

Since the most practical quantum codes cannot reliably 

detect large numbers of errors, it is also important that 

encoded gates (operators used to transform encoded 

data) be designed to limit the spread of errors. One of 

the primary techniques for achieving this is the use of 

transversal operators. 

Transversal operators provide an exceptionally 

simple and robust way of controlling the spread of errors. Thus, ideally, all gates necessary for quantum 

computation would be implemented transversally. 

However, the set of encoded gates that can be implemented transversally is determined by the quantum 

code employed, and despite substantial effort, no nontrivial quantum code has ever been found for which all 

encoded gates could be implemented transversally. 

This observation has led to a long-standing conjecture 

that transversal encoded gate sets cannot be universal. 

Recently, we have shown this conjecture to be true. A 

universal, transversal encoded gate set is incompatible 

with the most basic ability to detect errors. 

Bryan Eastin 

Like early classical devices, the first quantum computers will almost certainly be constructed from rather 

unreliable components. Consequently, strategies for 

managing errors are likely to be central to the operation 

of these machines. While classical error management 

is well understood, however, the requirements of quantum computing introduce new complexities that have 

forced researchers to reconsider the problem. It was, 

for example, not initially clear that the error correction 

of quantum data was even possible. Theoretical work 

on error management thus has an important role to play 

in the construction of a quantum computer. Indeed, in 

the last few years the gap between the error rate that is 

experimentally achievable and that which is proven to 

be sufficiently small as to allow for quantum computation has been narrowed as much by advances in theory 

as by experimental progress. The broader goal of this 

project is to contribute to the effort to build a quantum 

computer by improving our understanding of various 

strategies for managing errors, particularly quantum 

error correction and fault tolerance.

As in the classical case, quantum errors are rendered detectable by encoding the system of interest into 

a subspace of a larger, typically composite, system. A 

quantum code simply specifies which states of a quantum system correspond to which logical (encoded) 

information states. Errors that move states outside of 

the logical subspace can be detected by measuring the 

projector onto that subspace. Of course, not all errors 

can be detected; for any nontrivial code there are operators that act in a nontrivial way within the logical 

subspace. Most commonly, quantum codes are designed to permit the detection of independent, local 

errors and, as a consequence, are incapable of detecting 

some errors that affect many subsystems. 

For quantum computation it is necessary not only 

to detect errors but also to apply operators (gates) that 

transform the logical state of the code. Even when error 

processes are local and independent, however, the operations entailed in computing can generate correlated 

errors from uncorrelated ones. Thus, for error detection 

to be effective, it is important that the logical operators 

employed during a quantum computation be designed 

to limit the spread of errors. It is particularly important 

that operators do not spread errors within code blocks, 

where a block of a quantum code is a collection of subsystems that undergoes error detection independently 

of all other subsystems. Managing the spread of errors 

is the subject of the theory of fault-tolerant quantum 

computing. One of the primary techniques of this theory is the use of transversal encoded gates. 

Roughly, an operator is transversal if it couples no 

subsystem of a code block to any but the corresponding 

subsystem in another code block. (See Figure 1 for an 

illustration and a more detailed description.) Transversal operators are inherently fault tolerant. They can 

spread errors between code blocks, thereby increasing 

the number of locations at which a code block's error 

might have originated, but, since errors on different 

code blocks are treated independently, the total number 

of errors necessary to cause a failure is unchanged. 

This is in contrast to non-transversal operators, where, 

for example, an encoded gate coupling every subsystem in a code block might conceivably convert an error 
24 NIST Mathematical and Computational Sciences Division 

on a single subsystem into an error on every subsystem 

in the block. 

In view of the above, it would be highly desirable 

to carry out quantum computations exclusively using 

transversal encoded gates. To allow for arbitrary computation, it is necessary that the set of gates employed 

be universal, that is, that it be capable of implementing 

(through composition) any encoded operator on the 

logical state space to arbitrarily high accuracy. In spite 

of substantial effort, however, no gate set for a nontrivial quantum code has yet been found that is both 

universal and transversal. Consequently, a longstanding question in quantum information theory is 

whether there exist nontrivial quantum codes for which 

all logical gates can be implemented transversally. 

For an important subclass of quantum codes 

known as stabilizer codes, this question was recently 

answered in the negative. Zeng et al. [1] showed that 

transversal unitary operators are not universal for stabilizer codes on two-level subsystems (qubits); the 

companion result for the case of d-level subsystems 

(qudits) was proven by Chen et al. [2]. In reference 

[3], Emanuel Knill and I used a different approach to 

prove a stronger version, showing that the ability of a 

quantum code to detect arbitrary errors on component 

subsystems is incompatible with the existence of a universal, transversal, and unitary encoded gate set. 

The proof of this result is too involved to be included here, but it can be motivated as follows. 

Consider a quantum code such that 

=⊗⊗⊗= ∏

j

Z Z Z Z j Z θ θ θ θ

θ

i ii i ee...ee

is a valid logical gate for any θ. For sufficiently small 

θ, 

∏ +≈+≈ ∑

j

j

j

θ θ j i)i( θ ZIZIZ

where each Z j acts on a single subsystem. Thus, the 

code cannot detect arbitrary single-subsystem errors 

since otherwise it would sometimes mistake the encoded gate Zθ for an error. There are several flaws in 

this argument, but using the theory of Lie groups it can 

be generalized and made rigorous, as we demonstrate 

in [3]. 

Our result is fairly broadly applicable. The proof 

makes no assumptions about the dimensions of the 

quantum subsystems beyond requiring that they be 

finite. The quantum system encoded is assumed only to 

be nontrivial, that is, to have dimension greater than 

one. The precise structure of the quantum code and its 

initialization state are left unspecified. Nonetheless, 

like any impossibility proof, ours can be circumvented 

by violating the assumptions. 

In [3] we discuss several possible circumventions, 

some of which might make interesting projects for future research. While further study could be fruitful, 

however, it is important to note that the success of 

quantum computing does not depend on the existence 

of universal, transversal encoded gate sets. There are 

well-known schemes for achieving fault-tolerant universal quantum computation using ancillary states, 

where the preparation, testing, and coupling of the ancillary states makes extensive use of measurements and 

classical feed-forward. Thus, we have ruled out a desirable path to quantum computing but not an 

indispensable one. 

References 

[1] B. Zeng, A. Cross, and I. L. Chuang, Transversality 

versus universality for additive quantum codes, 

arXiv:0706.1382 (2007). 

[2] X. Chen, H. Chung, A. W. Cross, B. Zeng, and L. 

Chuang, Subsystem stabilizer codes cannot have a 

universal set of transversal gates for even one encoded 

qudit, arXiv:0801.2360, Physical Review A 78 (2008), 

012353. 

[3] B. Eastin and E. Knill, Restrictions on Transversal 

Encoded Quantum Gate Sets, arXiv:0811.4262, Physical 

Review Letters 102 (2009), 110502. 

Participants 

Bryan Eastin and Emanuel Knill 

Figure 1. Illustration of a transversal partitioning of a quantum 

code. The blue spheres represent quantum subsystems, which are 

arranged into four code blocks, each identified by a gray rectangle. 

The ellipses divide the code into transversal parts, that is, parts 

containing one subsystem from each code block. Given such a 

transversal partitioning, a transversal quantum gate is one that 

only couples subsystems in the same transversal part. 
NIST Mathematical and Computational Sciences Division 25

High Precision Calculations of Fundamental Properties of Fewelectron Systems

Computational scientists at MCSD and Indiana University have achieved record levels of accuracy in the 

development of computational methods for determining 

fundamental properties of three electron atoms. Their 

recent result for the nonrelativistic energies for four 

excited states of the lithium (Li) atom represents the 

highest level of accuracy ever reached (10-9 hartree) in 

atomic quantum computations with more than two electrons. Familiar chemical electron affinities and 

ionization energies involve these nonrelativistic energies, so today we may be witnessing again a situation 

in which quantum mechanical calculations yield more 

accurate determinations of these fundamental properties than can be measured experimentally. 

James Sims 

NIST has long been involved in supplying critically 

evaluated data on atomic and molecular properties such 

as the atomic properties of the elements contained in 

the Periodic Table and the vibrational and electronic 

energy level data for neutral and ionic molecules contained in the NIST Chemistry WebBook. Fundamental 

to this endeavor is the ability to predict, theoretically, a 

property more accurately than even the most accurate 

experiments. Our goal is to accomplish this for fewelectron atomic and molecular systems. 

Impressive advances have been made in the study 

of atomic and molecular structure, at both the experimental and theoretical levels. For atomic hydrogen and 

other equivalent two-body systems, exact analytical 

solutions to the nonrelativistic Schrödinger equation 

are known. It is now possible to calculate essentially 

exact nonrelativistic energies for atomic helium (He) 

and other three-body systems as well. Even for properties other than the nonrelativistic energy, the precision 

of the calculation has been referred to as “essentially 

exact for all practical purposes” [1], i.e., well beyond 

what can be achieved experimentally. Notwithstanding 

this theoretical advance, the data (experimental and 

theoretical) on atomic energy levels remains scarce, 

especially for highly ionized atoms. On the theoretical 

side, the availability of high precision results tails off 

as the state of ionization increases, not to mention 

higher angular momentum states. In addition, atomic 

anions have more diffuse electronic distributions, and 

therefore represent more challenging computational 

targets. On the diatomic molecular side, the presence of 

two fixed nuclei rather than one makes the integrals 

harder, hence exact analytical solutions to the molecular nonrelativistic Schrödinger equation are known (in 

the Born-Oppenheimer approximation) only for the 

one-electron H2

+ ion and other equivalent one electron 

systems. The challenge for computational scientists is 

to extend the phenomenal accomplishments on atomic 

helium to three, four, and more electron atomic states 

and to molecular systems. 

Dihydrogen, H2, is the first molecule whose dissociation energy was correctly predicted by quantum 

mechanical calculation (1968) before this same quantity was measured reliably by experiment. Today we 

may be witnessing again a situation in which quantum 

mechanical calculations yield more accurate determinations of this fundamental property than can be 

measured experimentally. Very high precision approximations are now available for molecular hydrogen 

(a two electron system) as a result of our calculated 

Born-Oppenheimer (BO) energies of 1

Σg

+

 states of H2. 

With an error of 1 in the 13th digit, our results [2] were 

the most accurate energy values ever obtained for a 

molecule of that size, 100 times better than the best 

previous value, calculated or experimental. 

The BO energies are used to provide the X1

Σg

+

 potential energy curve which is critical to determining, 

theoretically, the dissociation energy D0 of H2, as 

shown in Figure 2. Recently [3], a new theoretical bond 

dissociation energy was obtained for the hydrogen 

molecule by fitting the BO potential needed to our 

computed energies over 47 inter-nuclear distances from 

0.4-6.0 bohr. The resulting dissociation energy is accurate to 9 significant digits, two orders of magnitude 

better than the best previous experimental result. 

While these results are important, the methods 

used are perhaps equally so. The calculation requires 

solving an approximation of the Schrödinger equation, 

a central equation in quantum mechanics. Its solution 

can be approximated as an infinite sum, each additional 

term contributing a bit more to the accuracy of the result. For all but the simplest systems or a relative 

handful of terms, however, the calculation rapidly becomes impossibly complex. While very precise 

calculations have been done for systems of just three 

components such as helium (nucleus and two electrons), we were the first to attain such precision for H2, 

a system with two nuclei and two electrons. 

To make the problem computationally practical, 

we merged two earlier algorithms – one which has advantages in ease of calculation, and one which more 

rapidly achieves accurate results – into a hybrid, the 

Hylleraas-Configuration Interaction (Hy-CI) method, 

with some of the advantages of both. Our most recent 

work has shifted from the two-electron molecular case 

to the three-electron atomic case, using this hybrid 
26 NIST Mathematical and Computational Sciences Division 

method to extend the high precision accomplishments 

on atomic helium to atomic lithium, computing four 

excited states of the lithium atom to two orders of 

magnitude greater than has been done before. In calculations published in Physical Review A [4], we showed 

how Hy-CI calculations can achieve sub-nanohartree 

accuracy for the three-electron lithium atom. 

There are two reasons why sub-nanohartree nonrelativistic energies are desirable. First, the physical 

energies of interest, e.g., a transition energy or familiar 

chemical ionization potentials and electron affinities, 

involve these nonrelativistic energies. So the nonrelativistic energies must be calculated very accurately to 

guarantee the accuracy of downstream results. Second, 

when a nonrelativistic energy is obtained, a wave function is also obtained which can be used to compute 

other physical properties like oscillator strengths, electric polarizabilities, Fermi contact terms, or nuclear 

mass polarization corrections. Quantum theory tells us 

that a first order error in the wave function manifests 

itself as a second order error in the energy. This can be 

understood from the fact that the energy is a sum of 

kinetic and potential energy contributions, so there can 

be cancellation of errors between the two energy contributions leading to a precision in the energy 

determination that does not carry over to other properties. Hence the drive for high accuracy calculations. 

All results reported in our paper were obtained using quadruple precision (~30+ digits) floating point 

subroutines written in Fortran 90, and required removing several bottlenecks to highly accurate Hy-CI 

calculations: the three electron triangle integrals [5] 

and treating the remaining three-electron kinetic energy 

and nuclear attraction integrals [6]. We solved the 

secular equation, a generalized eigenvalue problem, 

using our own portable parallel inverse iteration solver. 

The assembly of the necessary matrices was also parallelized since that allows the total memory needed to be 

spread across the processors and eliminates a need to 

communicate matrix elements between processors. For 

a 17,190 term wave function we achieved a 30-fold 

speedup on 32 processors for the order N3 step running 

on the NIST cluster of 411 Intel Xeon, Intel Pentium 4, 

and AMD Opteron processors running Red Hat Linux. 

In going from He (two electrons) to Li (lithium, 

three electrons) to Be (beryllium, four electrons), the 

high precision calculation situation degrades to the 

point that already at four electrons (Be) there are no 

calculations of the ground or excited states with an 

error of less than 10-6 a.u. We are currently building the 

infrastructure for a calculation on the ground state energy level of the Be atom which has involved our 

dealing with four-electron integrals [7]. The beryllium 

calculation will determine whether really accurate calculations on atoms with N¥4 become a reality. 

References 

[1] G. W. F. Drake and W. C. Martin, Ionization Energies 

and Quantum Electrodynamic Effects in the Lower 1sns 

and 1snp Levels of Neutral Helium (4He I), Canadian 

Journal of Physics 76 (1998), pp. 679–698. 

[2] J. S. Sims and S. A. Hagstrom, High Precision Variational Calculations for the Born-Oppenheimer Energies 

of the Ground State of the Hydrogen Molecule, Journal 

of Chemical Physics 124 (2006), 094101. 

[3] K. Piszczatowski, G. Lach, B. Jeziorski, M. Przybytek, 

J. Komasa, and K. Pachucki, Theoretical Determination 

of the Dissociation Energy of Molecular Hydrogen, J. 

Chemical Theory and Composition, in review. 

[4] J. S. Sims and S. A. Hagstrom, Hy-CI Study of the 2 2

S 

Ground State of Neutral Lithium and the First Five Excited 2

S States, Physical Review A 80 (2009), 052507. 

[5] J. S. Sims and S. A. Hagstrom, Mathematical and Computational Science Issues in High Precision Hy-CI 

Variational Calculations I. Three-electron integrals, 

Journal of Physics B 37 (2004), pp. 1519–1540. 

[6] J. S. Sims and S. A. Hagstrom, Mathematical and Computational Science Issues in High Precision Hy-CI 

Variational Calculations II. Kinetic Energy and Nuclear 

Attraction Integrals, Journal of Physics B 40 (2007), pp. 

1575–1587. 

[7] J. S. Sims and S. A. Hagstrom, Mathematical and Computational Science Issues in High Precision Hy-CI 

Variational Calculations III. Four-electron integrals, in 

process. 

Participants 

James Sims. Stanley Hagstrom (Indiana University). 

Figure 2. Potential energy curves of H2 showing transitions 

from the ground state to the B′

1

Σu

+

 state in photodissociation. 
NIST Mathematical and Computational Sciences Division 27

Characterization of Electromagnetic Metafilms 

Increasing reliability and control of nanometer-scale 

fabrication has created new possibilities for engineering materials with novel optical scattering properties. 

These devices will have tremendous impact in a variety 

of commercial, defense-related, and scientific applications such as magneto-optic storage, microantenna 

arrays, and integrated optical circuits. Metafilms are 

novel two-dimensional structures in which microscatterers are arranged as a planar lattice. The result 

is a new class of flexible, subwavelength coatings that 

exhibit strong EM scattering signatures. Novel applications proposed for metafilms include controllable 

“smart surfaces,” microscale cavity resonators, novel 

waveguiding structures, and biomedical devices. Specialized computational tools capable of accurately 

guiding design parameters are essential for successful 

development of such devices. This past year NIST scientists made significant contributions to metafilm 

research by comparing a variety of theoretical predictions with a new high-accuracy simulation of welldefined vector EM scattering problems. 

Andrew Dienstfrey 

Metamaterials comprise one class of nano-scale electromagnetic (EM) scattering structures that have been 

the subject of ongoing study at NIST. A metamaterial

consists of a three-dimensional lattice of tailor engineered scatterers embedded in a uniform substrate, the 

later referred to as “free-space.” The scatterers exhibit 

geometry and composition such that their size is a fraction of a wavelength of light in free-space, nevertheless 

their internal structures admit resonant scattering. The 

resulting tension between length scales suggests the 

possibility that EM characteristics of the more extreme 

resonant variety may be observed in a bulk homogenized sense. Over the past decade many theories have 

emerged to predict the homogenized propagation characteristics of metamaterials. Many of these theories 

predict effective medium properties with exotic “double-negative” or “left-handed” EM propagation 

character, hence their considerable interest for both 

applications and theory. Metafilms are 2-D analogs of 

metamaterials in which micro-scatterers are arranged 

as a planar lattice. In the following the metafilm is 

composed of magneto-dielectric spherical particles in a 

planar, square-lattice configuration as in Figure 3. 

Generalized Sheet Transition Conditions. There is 

disagreement in the theoretical community on the “correct” structural form for modeling EM scattering by a 

metafilm. One approach models the metafilm scattering 

problem as a bulk metamaterial with a small effective

thickness, d. With fixed d reflection and transmission 

coefficients (obtained by either computation or measurement) may be inverted for effective material

properties using classical thin-film inversion techniques. Alternatively, the NIST approach models a 

metafilm as an infinitesimal sheet of electric and magnetic polarizabilities. Electromagnetic scattering across 

such a sheet is described by generalized sheet transition 

conditions (GSTC)

[ ]

[ ]

[ ] ( )

[ ] ( ) )4(

)3(

ˆ ˆ )2(

ˆ ˆ )1(

z tH

z tE

z

zz

tE Et

z

zz

tE Ht

B

E

i E

i H

H

E

Ez zH

Hz zE

⋅⋅∇=

⋅⋅−∇=

∇×−⋅=×

∇×+⋅=×

μ μα

ε εα

ωμα α

ωεα α

Here the subscript t refers to transverse components, 

the polarizability tensors αE and αH are assumed diagonal, and the brackets <···> and [···] represent the 

average and difference across the boundary z = 0. Details can be found in [3]. In forward or predictive 

modeling, expressions for the tensors αE and αH are 

substituted into (1)–(4) and the reflection and transmission coefficients for the metafilm are computed. 

Alternatively, assuming the GSTC’s, reflection and 

transmission coefficients may be inverted for the underlying tensors and comparisons can between theory 

and “experiment” can be performed. 

This year we performed a detailed analysis to 

compare these competing theories. We observed that 

the artificial thickness parameter required of the metamaterial thin-film approach results in non-unique 

recovery of the effective medium properties, ε = ε(d) 

and similarly for μ(d). An example of this non-trivial 

dependence is shown in Figure 4. As there is no a priori way to fix d, this non-uniqueness poses a 

fundamental problem for this approach. Although similar problems exist in classical thin film analysis, the 

Figure 3. Metafilm composed of magneto-dielectric spherical particles with radius a, permittivity ε2, permeability μ2, and lattice 

spacing p. (a = 10mm, p = 25.59mm, ε2 = 2ε0 and μ2 = 900μ0). 
28 NIST Mathematical and Computational Sciences Division 

quantitative extent in the analysis of metafilms was 

previously unexplored. By contrast, the NIST approach 

retrieves unambiguous dielectric tensors that compare 

qualitatively and quantitatively with theoretical predictions. These results are documented in [2]. 

All Angle Total Reflection. The above analysis was 

performed assuming planewave incident light normal 

to the plane of scatterers, i.e. = kzk ˆ . In an attempt to 

further validate the GSTC theory, scattering under offaxis illumination was computed and compared with the 

full angle-dependent GSTC theory. In the course of this 

work a surprising observation was made. Fullyresolved EM computations predicted that the lossless 

metafilm structure Figure 3 would exhibit angleindependent total reflection at specific resonant frequencies, for example see Figure 5. This strong, angleindependent scattering by a sub-wavelength structure 

was entirely unexpected. Subsequent analysis of the 

model polarizability tensors coupled with the GSTC’s 

revealed that the NIST model also predicted this effect. 

The results of this analysis appear in [1]. 

Next Steps. We plan to continue investigating the 

GSTC model. In our study of angular-dependent scattering, detailed analysis of the GSTC’s revealed 

surprising symmetries required for solvability. Briefly, 

for general planewave incidence all E and H fields may 

be expressed with respect to complex scattering coefficients in the canonical orthogonal frames defined by 

the scattering plane, i.e., the Rte, Rtm, Tte, and Ttm reflection and transmission coefficients. In doing this, the 

GSTC equations yield a linear system that may be 

solved either for scattering coefficients (R’s and T’s) in 

terms of polarizability tensors (αE, αH) or conversely. 

In attempting to develop formulas for this inversion we 

discovered that the GSTC system is degenerate. Careful factorization of the linear system revealed four 

independent equations for the six unknown diagonal 

components of the polarizability tensors. What is more, 

we derived a solvability condition that physically corresponds to rotational symmetry of the scattering 

structure about the z-axis. For the spherical scattering 

structure depicted in Figure 3 this symmetry is reasonable. In fact, since our analysis assumed it from the 

outset this symmetry (solvability) constraint went unrecognized. That the GSTC's require such symmetry in 

general is unexpected, and suggests that fundamental 

terms may be missing from the equation. We also plan 

to explore applications for metafilms to problems in 

optical waveguiding. Finally, recent results have stimulated reanalysis of the computation of the Casimir force 

in the presence of metafilm scattering structures. We 

hope to consider this problem as well. 

References 

[1] J. A. Gordon, C. L. Holloway, and A. Dienstfrey, “A 

physical explanation of angle-independent reflection and 

transmission properties of metafilms/metasurfaces,” 

IEEE Antennas and Wireless Propagation Letters 8

(2009), pp. 1127–1130. 

[2] C. L. Holloway, A. Dienstfrey, E. F. Kuester, J. F. 

O’Hara, A. K. Azad, and A. J. Taylor, “A discussion on 

the interpretation and characterization of metafilms/metasurfaces: The two-dimensional equivalent of 

metamaterials,” Metamaterials 3 (2009) pp. 100–112. 

[3] E. F. Kuester, M.A. Mohamed, M. Piket-May, and C. L. 

Holloway, “Averaged transition conditions for electromagnetic fields at a metafilm,” IEEE Transactions on 

Antennas and Propagation 51(10) (October 2003) pp. 

2641–2651. 

Participants 

Andrew Dienstfrey. C. Holloway and J. Gordon (NIST 

EEEL). 

Figure 5. TE reflection coefficient shows all-angle, total reflection 

at resonance near f ≈ 3.5GHz. 

Figure 4. The retrieved effective permittivity of a metafilm modeled

as a thin film exhibits a strong dependence on assumed layer thickness d. Same geometry and composition as above. 
NIST Mathematical and Computational Sciences Division 29

Recovering from Physically Impossible Reconstructions in 

Blind Deconvolution

A pervasive problem in imaging science is that most 

recorded imagery is blurred and degraded by noise. In 

many cases, the cause of the blur is only poorly known, 

and a reliable mathematical formula for the point 

spread function (psf) describing the blur is not available. Blind deconvolution seeks to reconstruct both the 

unknown psf k(x,y), and the unknown sharp image 

f(x,y), given only the recorded blurred image g(x,y). 

This is a daunting problem fraught with serious 

mathematical difficulties, yet it is one that is of vital 

interest in numerous scientific contexts. Considerable 

effort continues to be directed at finding the right 

mathematical and computational formulation that 

would result in reliable blind deconvolution of large 

size scientific imagery. 

Alfred Carasso 

An important theoretical approach to blind deconvolution seeks to formulate that problem in terms of a 

variational principle. Using a-priori mathematical constraints on the unknown blur k(x,y), and on the 

unknown sharp image f(x,y), that approach seeks to 

obtain the pair {k, f} as the unique minimizer of an 

appropriate cost functional. In recent years, several 

such linear and nonlinear functionals have been proposed. However, it may be anticipated that formulating 

a variational principle that would result in full deblurring of images with unknown blurs, remains a 

challenging problem in mathematical analysis. Indeed, 

the great diversity of possible image blurs, together 

with the difficulty of encapsulating the subtleties of the 

unknown sharp image into a mathematical statement 

that can serve as an a-priori constraint, make it unlikely 

that a useful all-encompassing principle can be found. 

Recently, a highly efficient nonlinear variational 

approach was developed at the University of Bremen, 

[1]. This method is based on seeking that solution 

which differs from an educated first guess {k †

(x,y), 

f †

(x,y)} by the least amount in the L2 norm. The 

blurred image g(x,y) is used for f †

(x,y). The strength of 

this so-called Minimum Norm Solution (MNS) technique is that it can be formulated and calculated in 

Fourier space, and efficiently implemented using FFT 

algorithms. Indeed, almost instantaneous blind deconvolution ensues, even for 1024ä1024 imagery. 

The MNS method is critically examined in a recent 

paper by A. Carasso in the SIAM Journal on Applied 

Mathematics [2]. Fourier space behavior is used to 

identify and document three basic flaws in the MNS 

procedure. First, classical theorems of Bochner, 

Schoenberg, and Pólya are used to show that MNSdetected optical transfer functions (otf) are generally 

spurious, and are not characteristic functions. Such 

detected otfs must therefore correspond to non-physical 

psfs with negative values. Next, it is shown that only 

limited partial deblurring is theoretically possible in the 

MNS method, full deblurring being precluded by the 

MNS variational formulation. Lastly, it is shown that 

the MNS approach induces an erroneous relationship at 

high frequencies between the detected optical transfer 

function and the Fourier space reconstructed image. 

Each of these two objects becomes delusive at high 

frequencies. Unexpected spurious relationships between detected blurs and reconstructed images also 

occur in other variational deconvolution procedures. 

Figure 6. Failure of MNS approach in synthetically blurred 

1024ä1024 USS Kittyhawk image. Top: MNS detected otf (solid 

line) displays non-monotone behavior, at variance with classical theorems, and is not a characteristic function. This is 

confirmed in the bottom image showing that the defective otf 

above corresponds to a physically impossible, highly oscillatory psf k † (x,y), exhibiting sustained negative values. While k †

integrates to 1, its negative part integrates to -395. 
30 NIST Mathematical and Computational Sciences Division 

In the latter part of Carasso's SIAM paper [2], it is 

shown how to extract useful information from defective MNS reconstructions in the important special case 

of class L blurs, which includes Gaussians, Lorentzians, and other heavy-tailed densities, together with 

their convolutions. This intervention requires analytical 

considerations extraneous to the variational framework 

in [1], and involves several steps. First, the high frequency portion of the non-physical MNS-detected blur 

is discarded and a class L blur is associated with the 

low frequency part using least squares fitting. Next, 

this rectified blur is interpreted as the p

th convolution 

root of the true blur for some postulated real number p 

¥ 2. Infinite divisibility of this candidate true blur allows deconvolution to be performed in slow motion, by 

marching backwards in time in an associated parabolic 

equation, with the blurred image as data at t=1. Behavior of the image evolution as t∞0 can be monitored and 

used to readjust the value of p. Previously developed 

APEX/SECB methodologies [3] provide the necessary 

computational tools to make such an approach feasible. 

Figure 6 illustrates one fundamental defect in the 

MNS method. In the top graph, the dashed curve depicts an educated guess k †

 for the psf in the blurred 

USS Kitty Hawk image. The solid curve is the MNS 

detected otf. The oscillatory non-monotone behavior in 

that curve contradicts classical theorems of Schoenberg 

and Pólya on characteristic functions. This is confirmed in the bottom graph, which shows that the 

corresponding psf has sustained negative values extending over the entire range, and hence cannot

represent a physical blur. Thus, the corresponding 

MNS image reconstruction is of questionable validity. 

The next two figures illustrate the remarkable effectiveness of Carasso's intervention in blind 

restoration of a blurred image of Uranus and its moons. 

In Figure 7, the non-physical oscillatory MNS detected 

otf (solid line), restricted to the frequency range |ξ| < 

50, is approximated by a class L otf (dashed curve), 

with the latter now assumed valid over entire frequency 

range. In Figure 8, the partially deblurred MNS reconstruction is shown to be equivalent to marching 

backwards in time from t=1 to t=0.6 in the governing 

associated parabolic equation. Such MNS deblurring 

reveals 5 Uranus moons. However, using the new class 

L otf from Figure 7, Carasso's SECB intervention can 

continue the backwards solution much further, to 

t=0.05, and can recover 11 Uranus moons. (Uranus has 

27 known moons). 

Coupled with Carasso's recovery methodology, 

MNS blind deconvolution is a valuable tool, one that 

can be used to verify the results of other blind deconvolution procedures. Future work will apply this 

technique to new state of the art scanning electron microscope imagery. 

References 

[1] L. Justen and R. Ramlau, A non-iterative regularization 

approach to blind deconvolution, Inverse Problems 22

(2006), pp. 771-800. 

[2] A. S. Carasso, False characteristic functions and other 

pathologies in variational blind deconvolution. A 

method of recovery, SIAM Journal on Applied Mathematics 70 (2009), pp. 1097-1119. 

[3] A. S. Carasso, The APEX method in image sharpening 

and the use of low exponent Lévy stable laws, SIAM 

Journal on Applied Mathematics 63 (2002), pp. 593-

618. 

Figure 8. Blind deconvolution of Uranus image by marching 

backwards in time, from t=1.0 to t=0, using the class L otf detected 

in Figure 7. Partially deblurred MNS image corresponds to solution 

at t=0.6, and recovers 5 moons. However, Carasso’s SECB methodology can continue solution to t=0.05, and recover 11 moons. 

Figure 7. Blurred Uranus image experiment. Least squares fit of 

false MNS otf (solid line), with physically valid class L otf (dashed 

line) on restricted frequency interval |ξ| § 50. Dashed curve is 

postulated to be pth convolution root of true otf, with p=2.5. This 

tentative value can be readjusted by marching backwards in time.
NIST Mathematical and Computational Sciences Division 31

A Quality Pre-Processor for Biological Cell Images 

High content screening (HCS) has become a critical 

method for rapid large scale analyses in cell biology 

and drug discovery. HCS uses automated microscopes 

and sample preparation and includes the acquisition 

and analysis of cellular images without human intervention. Typical HCS-based experiments can involve 

the analysis of more than a million cells [1, 2, 3]. The 

process of automatically identifying objects and their 

boundaries in an image, known as image segmentation, 

is the most important part of analyzing biological image data. Consequently, many segmentation methods 

have been published, including histogram-based, edgedetection-based, watershed, morphological, and stochastic techniques. However, segmenting the same 

image with different methods can lead to different 

masks for a given cell, and hence different estimates of 

cell characteristics (e.g., area, perimeter, etc.). Accurate segmentation is one of the critical challenges 

facing high content screening. Therefore it is important 

to obtain a predictive, objective and efficient measure 

of the segmentation quality resulting from a particular 

method, on a cell-by-cell basis. Not all cells in the 

same image will lead to segmentations with the same 

accuracy. We have developed a simple, yet reasonably 

faithful, measure of the quality of an individual cell in 

an image to facilitate the automated choice of segmentation method and focus setting. 

Adele Peskin 

The goal of this project is to develop a quality index for 

the representation of a cell in a biological image. This 

index is determined on a cell-by-cell basis by assessing 

the characteristics of the edge of each cell. As such, it 

is different from other quality indices proposed previously which focus on overall image noise and 

distortion [3, 4]. 

Data Description. To define a quality index for each 

cell, we first examine a series of images and their 

ranges in quality, prepared by our collaborators in the 

NIST Chemical Sciences and Technology Laboratory. 

These reveal the types of features associated with cell 

images acquired in practice. These images are of A10 

rat smooth vascular muscle cells and 3T3 mouse fibroblasts stained with Texas Red cell body stain. For each 

set, we compared images acquired under five sets of 

imaging conditions, combinations of low, medium, and 

high illumination, as well as with optimal and nonoptimal filters. Figure 9 provides a visual rendition of 

the five conditions for two example cells.

Quality Calculation. For each cell in an image, we 

look at pixel intensities within an isolated region containing the cell. We fit this distribution using a 3-

component Gaussian mixture via the EM (ExpectationMaximization) algorithm (background, edge region, 

cell) [6]. For each pixel intensity value between the 

mean of the background distribution and the mean of 

the edge distribution, we also calculate an average 

value for the magnitude of the gradient at that intensity. 

We consider the gridded data within a bounding box on 

the image that contains a single cell. At each pixel location, pij, with pixel intensity, Iij, we find a local 

derivative of that intensity and its magnitude, using a 

Sobel mask, and then find an average gradient magnitude for each intensity value.

The average gradient magnitude curve, G(I) vs. I, 

reflects two cell characteristics. The initial slope of the 

curve indicates the sharpness of the edge: the steeper 

the slope, the sharper the edge. The shape of this curve 

also indicates a feature of the cell edges. For high illumination settings with very sharp edges, the average 

gradient magnitude increases monotonically across the 

range from the mean of the background to the mean of 

the edge distribution. For the low and medium illumination images, the gradient increases, reaches a 

maximum, and then falls, as seen in Figure 10. 

We find the 3-component pixel intensity distribution, find and smooth the gradient in the region of the 

cell edge, and then find the expected fraction of the 

edge that lies in the 2-pixel-length region outward from 

the edge point of maximum intensity gradient. As an 

Figure 9. Five images acquired under different conditions, along with a ground truth mask (right). Individual pixel intensities are color-coded

over the range in each image, to show differences in edge sharpness. 
32 NIST Mathematical and Computational Sciences Division 

example, look at a cell with a very thin connection between the main part of the cell and pseudopodia (see 

Figure 11). These five cells have quality indices of 

1.015, 1.000, 0.852, 0.467, and 0.429 respectively. 

Physically, that represents cell edges that cover from 2 

to 4 pixel to pixel length distances. 

Data Creation. We then synthetically altered cell image examples to examine the effects of different 

imaging parameters on the quality index. Figure 12 

shows the effects of altering the spread of the background pixel intensities and of shifting the edge pixel 

intensities closer to the background pixel intensities.

Conclusions. The goal of our quality index is to provide an indication of the segmentation method that 

should be used to yield the highest accuracy in derived 

cell quantities (e.g., area, perimeter, etc.). Ideally, we 

would like to be able to distinguish between those 

segmentation methods that over- or under-estimate cell 

area from those that provide unbiased estimates of cell 

area, based on the quality of the edge of the cell. Because our quality index as defined here measures only 

edge clarity, we need to modify it to better capture cell 

image features. Our future work will includes a more 

thorough investigation of the properties of this quality 

index, particularly its relation to accuracy of cell area. 

We hope to use this quality index (modified to include 

perimeter to area ratios) to automatically select segmentation methods during rapid analysis of large 

numbers of biological images. 

References 

[1] V. C. Abraham, D. L. Taylor, and J. R. Haskins, High 

content screening applied to large-scale cell biology, 

Trends in Biotechnology 22 (2004), pp. 15–22 . 

[2] K. A. Giuliano, J. R. Haskins, and D. L. Taylor, Advances in High Content Screening for Drug Discovery, 

Assay and Drug Development Technologies 1(4) (2003) 

pp. 565–577. 

[3] X. Zhou and S. T. C. Wong, High content cellular imaging for drug development, IEEE Signal Processing 23(2) 

(2006), pp. 170–174. 

[4] Z. Wang, A. C. Bovik, H. R. Sheikh, and E. P. 

Simoncelli, Image quality assessment: from error 

visibility to structural similarity, IEEE Transactions on 

Image Processing 13(4) (2004), pp. 600–612. 

[5] Z. Wang and A. C. Bovik, A universal image quality 

index, IEEE Signal Processing Letters 9(3) (2002), pp 

81-84. 

[6] A. P. Peskin, K. Kafadar, A. M. Santos, and G. G. 

Haemer, Robust Volume Calculations of Tumors of 

Various Sizes, in 2009 International Conference on 

Image Processing, Computer Vision and Pattern 

Recognition (2009). 

Participants 

Adele Peskin. Karen Kafadar (Indiana University). 

Alden Dima (NIST ITL). Anne Plant, John Elliot, and 

Michael Halter (NIST/CSTL). 

Figure 10. Left: Averaged magnitudes of the pixel gradient for

low (red), medium (blue), and high (purple) illumination settings

with the optimal filter. Right: For low, using optimal (red), vs. nonoptimal (blue) filters. 

Figure 11. A cell mask and enlarged pictures of the lower right hand marked section, in order of decreasing quality index.

Figure 12. From left to right: (a) Original cell imaged with medium illumination, optimal filter. (b) Synthesized cell of the same type as (a), 

made as a test case. (c) and (d) Two cells identical except for the spread of the background pixels; standard deviations are 10 and 12 compared

to original 5. (e) and (f) Two cells identical except that the edge distribution has been shifted from a mean of 225 to a mean of 215 and 205. 
Yearly Report, Fiscal Year 2009 33

Computational Modeling of the Flow of Dense Suspensions

Understanding the flow properties of complex fluids 

such as suspensions is of technological importance, but 

presents significant theoretical challenges. Computational models of such systems require the tracking of 

boundaries between different fluid/fluid and fluid/solid 

phases, which can be extremely difficult. In our models 

we use methods from Dissipative Particle Dynamics 

and Smoothed Particle Hydrodynamics to simulate the 

flow of cement and concrete. By using these methods, 

we avoid the traditional difficulties in tracking boundaries by instead tracking individual mesoscopic scale 

particles. We have developed efficient parallel implementations of these methods, allowing us to simulate 

larger systems, with realistically shaped particles (instead of spherical particles), and for longer simulation 

times, than was previously possible. The results from 

our simulations have been validated against both theoretical results and results from physical experiments. 

Analysis of recent simulations has led to an advance in 

the conceptual understanding of the mechanisms that 

control the onset of flow in suspensions. 

William George 

While concrete may sound like a low tech material, 

modeling and predicting its rheological properties, that 

is, its physical flow properties, remains a great challenge. Concrete is actually a very complex fluid 

involving phenomena taking place over many length 

and time scales. For example, concrete is typically 

composed of cement (micrometer scale), sand (millimeter scale) and aggregates (centimeter scale). 

Clearly, no computer simulation can simultaneously 

include all these length scales, so we have developed a 

multi-scale approach linking phenomena at different 

length scales. For example the rheological properties of 

cement paste can be used as input to simulations of a 

mortar (cement paste and sand) and mortar properties 

as input into simulations of concrete (mortar and coarse 

aggregate). However, even at each representative 

length scale there can be considerable variation and 

other factors to consider such as the addition of chemical admixtures in cement paste, variation by over a 

factor of 100 in the size of aggregates, and shape variation of aggregates (rounded or crushed). 

There is a wide body of literature on the rheology 

of suspensions. Unfortunately, most analytical results 

can only describe the properties of simple suspensions, 

most typically very dilute suspensions. Understanding 

the flow of more complex suspensions (e.g., dense suspensions, suspensions composed of particles with 

different size and shapes, etc.) represents a great challenge. Current prediction of the rheological properties 

of dense suspensions is usually based on semiempirical models such as the Kreiger-Doughtry equation which is used to predict the relative viscosity as a 

function of maximum packing fraction and knowledge 

of the intrinsic viscosity of a representative particle of 

the suspension. Rigid bodies having shapes like 

crushed aggregates are poorly accounted for by such 

approaches and there is no accepted theoretical understanding. Direct measurement of the rheological 

properties of such systems is itself a significant challenge, and there is no standard practice or database for 

validation of measurements. 

As a result, we have turned to computational modeling of suspensions to better understand the flow of 

such systems and to give interpretation of measurements a better scientific basis. Modeling complex 

suspensions using computational fluid dynamics, while 

promising, remains a challenge because it is difficult to 

track boundaries between different fluid and solid 

phases. Thus, its practical applicability to modeling 

dense suspension remains in question. 

Recently, a new computational method called dissipative particle dynamics (DPD) has been developed 

which holds promise for modeling such complex fluids. DPD has the advantage that it can naturally 

accommodate many boundary conditions while not 

requiring meshing (or remeshing) of the computational 

domain. On the surface, DPD looks very much like a 

molecular dynamics algorithm where particles, subject 

to inter-atomic forces, move according to Newton's 

laws. However, the particles in DPD are not atomistic, 

but instead are a mesoscopic representation of the fluid. 

The interactions between DPD particles are typically 

produced by three classes of forces: conservative, dissipative, and random. The conservative force is simply 

a central force, derivable from some potential. The 

dissipative force is proportional to the difference in 

velocity between particles and acts to slow down their 

relative motion. The dissipative force can be shown to 

produce a viscous effect. The random force (usually 

Figure 13. A snapshot of a simulation. Only the aggregates with 

the highest stresses are shown in this display, with the other 

aggregates shown as outlines. Actual stress values are painted 

on the surface of the aggregates.
34 NIST Mathematical and Computational Sciences Division 

Gaussian noise) helps maintain the temperature of the 

system while producing some viscous effects. In order 

to maintain a well defined temperature consistent with 

a fluctuation-dissipation theorem, coefficients describing the strength of the dissipative and random forces 

must be coupled. We have adopted this approach as it 

has the potential of modeling rigid bodies with a wide 

variety of shapes. 

In addition to this DPD-based model, which simulates suspensions with Newtonian fluid matrices, we 

have also implemented a Smoothed Particle Hydrodynamics (SPH) model of dense suspensions. SPH, which 

also incorporates a mesoscopic particle model similar 

to DPD, simulates suspensions with non-Newtonian 

fluid matrices. This is important because while cement 

paste has a Newtonian fluid matrix (water), mortar 

(cement + sand) and concrete (mortar + aggregates) 

have non-Newtonian fluid matrices. In non-Newtonian 

fluids, viscosity varies with the local shear rate. 

In order to perform simulations on systems of sufficient size, density, and simulation time, we have 

developed highly parallel implementation of our models. This implementation, called QDPD, has run 

simulations efficiently using more than 2,000 processors. In these cases we often simulate systems 

consisting of over 1 million mesoscopic fluid particles 

and up to 10,000 aggregates. The output from these 

simulations includes the overall system stress over 

time, from which the viscosity can be computed, as 

well as periodic snapshots of the system that captures 

the location, orientation, and velocity of each of the 

system's aggregates. Post-processing of these snapshots 

generates an animation of the simulation which can be 

displayed in our immersive visualization system. One 

such simulation snapshot is shown in Figure 13, where 

only the high-stress aggregates are shown fully. 

Through awards granted by the DOE INCITE program, we have been able to run simulations on Intrepid, 

the IBM Blue Gene/P machine at the Argonne National 

Laboratory Leadership Computing Facility. Intrepid, as 

of November 2009, is 8th in the Top 500 supercomputing list1

. In 2008 and 2009 we were granted (and used) 

750,000 CPU-hours on Intrepid. (1 CPU-hour = 1 hour 

of processing time on 1 Intrepid CPU). On Intrepid our 

simulations use 512, 1000, or 4096 CPUs, depending 

on the size of the system being simulated. For 2010 we 

have been awarded 2 million CPU-hours of compute 

time on Intrepid, which will allow us to complete more 

large simulations than in previous years. 

We have made several improvements to QDPD in 

the past year. One is a complete re-write of the file 

input/output routines to use parallel I/O to and from 

single files, rather than to per-process files. This 

greatly improves post-processing of QDPD output, 

especially when running on more than 1000 processors, 

and makes simulation restarts from checkpoints much 

simpler. We have also added the ability to change the 

simulation time step dynamically to avoid instabilities 

in the simulations without the need to stop and restart. 

This allows the simulator to find the largest time step 

that it can use safely, thus reducing the total time 

needed to complete a simulation. 

In general, our simulations have been designed to 

explore the effects of particle size and shape distribution, the volume fraction of aggregates in the 

suspension, and other relevant system parameters, on 

the yield stress and viscosity of suspensions. These 

have led to a variety of very useful insights. For example, we discovered that, for dense suspensions, a 

previously unidentified physical mechanism plays a 

role in yield stress, which is a measure of the internal 

resistance to the startup of flow. This discovery was 

made by analyzing the pattern of contacts between the 

aggregates in snapshots from the simulation [1]. 

Figure 14 shows the location of high stress regions in 

these simulations before the aggregates begin to move 

(left image) and after (right image). 

References 

[1] N. Martys, D. Lootens, W. George, P. Hebraud, Contact 

and stress anisotropies in start-up flow of colloidal suspensions, Physical Review E 80 (2009), 031401. 

Participants 

W. George, S. Satterfield, J. Hagedorn, M. Olano, and 

J. Terrill. N. Martys (NIST BFRL). Virtual Cement and 

Concrete Testing Laboratory Consortium. 

http://math.nist.gov/mcsd/savg/parallel/dpd

http://math.nist.gov/mcsd/savg/vis/concrete

 

1 http://www.top500.org/

 

Figure 14. Left: This image illustrates the correlation of shear 

stress in a suspension. The green represents high stress regions and 

the red represents lower stress regions. The white arrow highlights 

the stress that is highly correlated in the “compression quadrant”, 

that is, with aggregates being pushed together as we begin to strain 

the system. The black arrows indicate the direction of the applied 

shearing force. Right: Image is from the same simulation, showing 

the stress in the system after is has been sheared. Note the emergence of an ordered six-fold symmetry.
Yearly Report, Fiscal Year 2009 35 

Visualization of Radio Frequency Propagation for Body Area 

Networks 

A Body Area Network (BAN) is a network intended to 

be used in or around the human body. While this is an 

emerging field, networks of medical sensors are anticipated to be a primary application. Such sensors would 

either be attached to or implanted in the human body 

and would communicate wirelessly both within the 

body and to devices outside. For this technology to 

develop, greater understanding of radio frequency 

propagation through the human body is needed. In this 

project, we seek to provide such insight. This work is 

assisting IEEE standards committees to develop communications standards and will help RF engineers to 

design effective communication devices for BANs. 

John Hagedorn 

Body Area Networks (BANs) consisting of radioenabled wearable and implantable sensory nodes are 

poised to be a promising technology with novel uses in 

pervasive health information technology, such as smart 

pills for precision drug delivery, glucose monitors, 

blood pressure sensing systems, 

and eye pressure sensors for 

glaucoma patients. However, 

numerous challenges remain, 

including size, cost, energy 

source, sensing/actuator technology, and transceiver design. 

Some recent advances in microelectronics indicate that the 

technology to achieve ultrasmall and ultra low power devices is within reach. 

For this technology to develop, greater understanding of 

radio frequency (RF) propagation through the human body is 

needed. In this project, we seek 

to provide such insight. In particular, we are providing 

technical support to IEEE standards committees to develop 

communications standards and 

to RF engineers who seek to 

design effective communication 

devices for BANs. 

Because experimentation 

on human subjects is not feasible, RF propagation through the 

human body is being modeled 

using a 3D full-wave electromagnetic field simulator. The software we are using is 

HFSS, a product of the Ansoft Corporation. Simulations are based on a 3D human body model which 

includes frequency dependent dielectric properties of 

300+ parts in a male body. 

The output data produced by this simulation software is brought into our custom 3D immersive 

visualization system, which enables researchers to 

study the modeled RF propagation through direct interactions with the data. Figure 15 shows a visualization 

of a simulation of a transmitter implanted in the right 

forearm with color coding on the surface of the body 

based on signal loss. 

Because communication with an implanted or 

wearable device can be done from any direction inside 

or outside the body, a true 3D environment is needed to 

capture, visualize and understand RF propagation to 

and from implants. The 3D immersive visualization 

platform that we have developed is shown in Figure 16. 

It includes several key components: three orthogonal 

screens that provide the visual display, motion tracked 

stereoscope glasses, and a hand-held motion tracked 

input device. The screens are 

large projection video displays 

that are placed edge-to-edge in a 

corner configuration. These three 

screens are used to display a single 3D stereo scene. The scene is 

updated based on the position of 

the user as determined by the 

motion tracker. This allows the 

system to present to the user a 

3D virtual world within which 

the user can move and interact 

with the virtual objects. The 

main interaction device is a 

hand-held three button motiontracked wand with a joystick. 

This virtual environment allows for natural interaction 

between experts with different 

backgrounds such as engineering 

and medical sciences. The user 

can look at data representations 

at any scale and position, move 

through data, change orientation, 

and control the elements of the 

virtual world using a variety of 

interaction techniques including 

measurement and analysis. 

All of these capabilities are 

extremely useful when studying 

Figure 15. Simulation of a transmitter implanted in 

the right forearm. Color denotes signal loss. 
36 NIST Mathematical and Computational Sciences Division 

RF propagation to and from medical implants. For example, we have implemented interactive tools for 

probing the 3D data fields. One tool enables the researcher to move the motion-tracked wand through the 

virtual scene, yielding a continuously updated display 

of the value of the data field at the position of the 3D 

wand. Another tool enables the user to interactively 

stretch a line segment through virtual body, and to generate graphs of the 3D data fields along that path (see 

Figure 17). We have applied these tools to the results 

of RF simulations with antennas of various characteristics at various locations throughout the body. 

Using the insight gained from interacting with the 

results of these simulations in the immersive environment, together with additional analysis, we have 

constructed a statistical model of signal strength path 

loss as a function of distance. These results are described in reference [1], which received the Best Paper 

Award at the IEEE 2009 Personal, Indoor, and Mobile 

Radio Communications Symposium. Further results 

were presented in an invited paper [2] at the 2009 International Symposium on Applied Sciences in 

Biomedical and Communication Technologies. Additionally, implant channel models that have been 

obtained as a result of our work have been adopted by 

the IEEE 802.15.6 Task Group on BAN that is developing an international communication standard for 

implantable and wearable devices [3,4,5]. 

We plan to continue this work in several areas. 

We will be using our simulation and visualization system to compare antennas of various placements and 

designs. We also plan to compare our simulation results with experimental measurements and to use our 

immersive visualization system to investigate deviations of simulations from experiments. In addition, the 

3D visualization platform is an ideal environment to 

observe and study creeping wave propagation on the 

body surface. We believe that the immersive visualization environment provides a flexible platform to 

visualize and understand RF propagation for BANs. 

References 

[1] K. Sayrafian-Pour, W. Yang, J. Hagedorn, J. Terrill, K. 

Y. Yazdandoost, A Statistical Path Loss Model for 

Medical Implant Communication Channels, Proceedings of the IEEE PIMRC 2009. 

[2] W. Yang, K. Sayrafian-Pour, J. Hagedorn, J. Terrill, 

Simulation Study of Body Surface RF Propagation for 

UWB Wearable Medical Sensors, IEEE International 

Symposium on Applied Sciences in Biomedical and 

Communication Technologies, November 2009. 

[3] J. Hagedorn, J. Terrill, W. Yang, K. Sayrafian-Pour, K. 

Y. Yazdandoost and R. Kohno, A Statistical Path Loss 

Model for MICS , September 2008, Report to the IEEE 

P802.15. ID: IEEE 802.15-08-0519-01-0006. 

[4] J. Hagedorn, J. Terrill, W. Yang, K. Sayrafian-Pour, K. 

Y. Yazdandoost and R. Kohno, MICS Channel Characteristic: Preliminary Results, May 2008, Report to the 

IEEE P802.15. ID: IEEE 802.15-08-0351-00-0006. 

[5] K. Y. Yazdandoost and K. Sayrafian-Pour, Channel 

Model for Body Area Network (BAN), April 2009, Report to the IEEE P802.15. ID: IEEE 802.15-08-0780-02-

0006. 

Participants 

John Hagedorn, Judith Terrill, and Adele Peskin. 

Kamran Sayrafian-Pour and Wen-Bin Yang (NIST 

ITL). 

http://math.nist.gov/mcsd/savg/vis/BodyAreaNetwork

Figure 16. Our 3D immersive visualization system. Several important components of the system are labeled and the user is 

shown interacting with 3D RF data produced by our simulations. 

Figure 17. Simulation of a transmitter implanted in the right forearm. The user has interactively drawn a line segment in the 3D 

virtual scene. A graph showing signal strength along that line is 

displayed. 
Yearly Report, Fiscal Year 2009 37 

Mathematical Knowledge Management 

The field of knowledge management deals with the 

practices used to create, represent, distribute and identify “knowledge.” The ubiquity of word processing, 

the web and the fact that the name of a successful web 

search engine is commonly treated as a verb reveal the 

broad impact of these techniques on society. 

But these successes are primarily achieved for 

knowledge that naturally takes a textual form. The 

characteristics of mathematical knowledge and its publishing present a unique set of challenges. The special 

symbols and notations require special means to write 

and represent. The ambiguities created by the disconnect between the familiar notations (mathematical 

presentation) and the underlying meaning (mathematical semantics) inhibit reuse and accessibility. Both 

effects create difficulties for finding mathematical material in the creation of search indices as well as in 

formulating search queries. 

As part of the Digital Library of Mathematical 

Functions project – a web oriented successor to the 

classic NBS Handbook of Mathematical Functions [1] 

– we began addressing these issues. In that context we 

developed a tool for the conversion from TeX, the de 

facto standard for writing in many scientific fields, to 

an abstract XML format with MathML – a W3C standard for representing mathematics with growing 

support in web browsers and applications. 

Beyond the flexibility provided for use on the web 

and accessibility, this format enables further developments for knowledge management. We explore the 

inference of mathematical semantics from the markup, 

enhancing the reuse of the formula. Techniques for 

search within document collections specialized to 

matching fragments of mathematics have been developed. Additionally, the half-million article corpus of 

scientific documents at the arXiv is being converted, 

providing a large set of documents for further use and 

experimentation. 

Bruce Miller 

In the future, interactive computing environments may 

combine computer algebra, theorem proving, numerical 

computation, graphics and word processing to create 

scientific documents containing complete, verified and 

semantically rich mathematical content. Presumably, 

such a format would easily support the finding and 

reuse of such mathematical material for further proof 

or computation. Meanwhile, most mathematics is at 

least published, if not created, in a less rigorous, printoriented, fashion. Still, one would like to maximize the 

accessibility, reuse and finding of such materials. 

Moreover, there exists a rich legacy of printed scientific literature that we can hardly afford to abandon. 

We will focus here on TeX, a powerful, but 

quirky, document formatting language, and the de facto

standard for authoring in most mathematically oriented 

scientific fields, due to its support of complex mathematics and high-quality typesetting. It was the clear 

choice for authoring in Digital Library of Mathematical 

Functions2

 project and is the preferred format for the 

preprint archive at arXiv3

. 

While conversion of TeX to simple HTML was 

possible, conversion to a more structured XML, particularly supporting the MathML standard was not. We 

thus developed LaTeXML4

, a tool for converting LaTeX to XML. The XML format used here preserves the 

document structure, but is easily transformed to HTML 

for use on the web. The mathematical content within 

the document receives special treatment; it is parsed 

with the goal of at least recognizing the mathematical 

structure of the expressions. With author supplied declarations more of the semantic intent of the formula 

may also be inferred. Further work in this area will 

allow more semantics to be synthesized with less author input, although given the inherent ambiguities of 

mathematical notation, perfect enhancement cannot be 

expected. In any case, the better the semantics obtained, the better the reuse in other contexts such as 

copying to computer algebra systems. Enhanced semantics also enhances search capabilities. 

The method used to enable searching for mathematics leverages the advanced state of text search 

engines. The idea is to convert the mathematics to a 

textual form that parallels the MathML markup but 

contains both the mathematical symbols and keywords 

reflecting the mathematical structure [2]. Search queries are similarly transformed into this textual form. 

Search accuracy and expressiveness is improved by colocating alternative forms for symbols including aliases, semantic terms and any presentational variations. 

Since mathematics expressions often contain complicated sub-expressions that are nevertheless not critical 

to a given search, a technique called query relaxation 

can be used to allow the match to skip over subexpressions in the quest for a reasonable number of 

matching documents. By textualizing in parallel to the 

MathML markup, the matching terms within the search 

result can be marked for highlighting providing for 

useful feedback. Figure 18 shows these several search 

features. 

 

2 http://dlmf.nist.gov/

3 http://arXiv.org/

4 http://dlmf.nist.gov/LaTeXML/
38 NIST Mathematical and Computational Sciences Division 

Meanwhile, a project at Jacobs University in Bremen has been carrying out its own research in 

mathematically aware search, as well as enhancement 

of semantics. A large body of mathematically oriented 

documents was needed for experimentation. Thus, we 

set about to convert the papers in the preprint archive at 

arXiv from TeX to XML [3]. A cluster of workstations 

carries out document conversions in parallel. 

Currently, 400,000 documents have been processed with the conversion of 60% being of good 

quality. This corpus is being used for several experiments in document understanding which attempts to 

integrate the grammar based parsing of the mathematics with natural language processing of the textual 

content; clearly both must be understood to fully understand a document. 

This exercising of LaTeXML has driven improved 

coverage and robustness, leading to a more generally 

useful tool. The rich variety of mathematical notations 

encountered is key to improving the mathematical 

grammar for semantic enhancement. 

References 

[1] M. Abramowitz and I. Stegun, Eds., Handbook of 

Mathematical Functions, National Bureau of Standards, 

1964. 

[2] B. Miller and A. Youssef, Technical Aspects of the 

Digital Library of Mathematical Functions, Annals of 

Mathematics and Artificial Intelligence 38 (2003), pp. 

121-136. 

[3] H. Stammerjohanns, M. Kohlhase, D. Ginev, C. David, 

and B. Miller, Transforming Large Collections of Scientific Publications to XML, Mathematics in Computer 

Science, Special Issue on Authoring, Digitalization and 

Management of Mathematical Knowledge, to appear. 

Participants 

Bruce Miller. Abdou Youssef (MCSD and George 

Washington University). Michael Kohlhase andHeinrich Stamerjohanns (Jacobs University, Bremen, 

Germany). 

Figure 18. A search with a query with both semantic (“Hankel”) and presentational (“Ai”) terms; the “~” stands for asymptotically equal. The 

result shows how the query was relaxed to skip over function arguments and leading factors to find the best matching asymptotic equation. The 

term-by-term highlighting gives feedback on how the match was made. 
Yearly Report, Fiscal Year 2009 39

Part III: Project Summaries 

Part III 

 Project Summaries
40 NIST Mathematical and Computational Sciences Division 
Yearly Report, Fiscal Year 2009 41

 Information Technology Laboratory Programs Information Technology Laboratory Programs

 

 

 

Enabling 

Scientific 

Discovery 

Information 

Discovery, Use 

and Sharing 

Complex 

Systems 

Quantum 

Information 

Virtual 

Measurement 

Systems 

Pervasive 

Information 

Technology 

Trustworthy 

Information 

Systems 
42 NIST Mathematical and Computational Sciences Division 

Enabling Scientific Discovery

Modern scientific research has become more and more dependent on mathematical, statistical, and computational tools for enabling discovery. The Enabling Scientific Discovery Program promotes the use of these 

tools to dramatically advance our ability to predict the behavior of a broad range of complex scientific and 

engineering systems and enhance our ability to explore fundamental scientific processes. This Program focuses on interdisciplinary scientific projects that involve novel computational statistics and the development 

of simulation methods and software. These efforts will have a foundational impact on scientific discovery 

throughout U.S. industry, government, and academia. – Anthony Kearsley, Program Manager

Recovering from Physically Impossible Reconstructions in Variational 

Blind Deconvolution 

Alfred Carasso 

See feature article, page 29. 

Computational Modeling of the Flow 

of Dense Suspensions 

William George 

Steve Satterfield 

John Hagedorn 

Marc Olano 

Judith Terrill 

Nicos Martys (NIST BFRL) 

Virtual Cement and Concrete Testing Lab Consortium 

See feature article, page 31. 

Sparse Representations in High Dimensional Geometry 

Bradley Alpert 

Yu Chen (New York University) 

Laboratory and clinical diagnostic systems often produce images that lack quantitative information or 

comparability across imaging systems, which can lead 

to difficulties in diagnosis and long-term management 

of disease that requires patient monitoring. Elaborate 

measurements are taken, digitally processed, then reported in essentially qualitative form. For MRI, 

functional MRI, and some other processes, formation 

of the images relies on filtering and Fourier transformation. 

 Although the computational recovery of functions 

from their Fourier representations is considered routine, aside from possible efficiency concerns arising 

from non-grid data, standard approaches assume that 

the functions are smooth and can be recovered with the 

discrete Fourier transform (computed via the FFT). 

When this assumption fails, as is evident from Fourier 

data that are not small at the highest frequencies measured, some sort of attenuation, or belling, scheme is 

typically used. This procedure tends to blur sharp features and, being more art than science, can lead to 

different functions (images) in different systems. Although this problem, and the Gibbs phenomenon, have 

received considerable attention from mathematicians in 

recent years, recently developed methods tend to work 

well only in an asymptotic sense; they do not exploit 

available data efficiently. 

 Despite inherent limits on resolution that can be 

obtained from truncated Fourier data, a change of assumptions from smooth to piecewise smooth can lead 

to significantly improved recoveries. Procedures implementing this idea are not completely established 

even for functions (signals) in one dimension; they are 

yet more challenging in two and three dimensions, 

where discontinuities may be expected to occur along 

mostly smooth curves or surfaces. This constraint, 

which magnifies the advantage of the piecewise 

smoothness assumption, must be appropriately reflected in the methods used. Bradley Alpert and Yu 

Chen are conducting extensive numerical experiments 

to understand this environment and to develop reliable 

procedures for these problems. 

 In recent months, a technique was developed for 

exact recovery of piecewise polynomial functions, in 

which the discontinuities between pieces are located by 

examining the empirical translation matrix of the Fourier data, which has eigenvalues of multiplicity 

corresponding to the order of the discontinuities. This 

structure, which is not invariant in the presence of 

noise or removal of the piecewise polynomial assumption, nevertheless is suggestive of the structure when 

these stipulations are removed. The eigenvalues then 

split into clusters, in which the mean of each cluster 

approaches the corresponding discontinuity with high 

order convergence. Current numerical experiments and 

implementations are aimed at turning these observations into a robust method for locating discontinuities. 
Yearly Report, Fiscal Year 2009 43

Once the discontinuity locations are known, recovery 

of the discontinuous function is a generally wellconditioned procedure consisting of solving a linear 

system of equations. 

 A related, yet more general, challenge is parsimonious, or sparse, representation and recovery of functions 

under assumptions appropriate to an application. Although there is considerable current interest in the 

mathematical community in these problems (L1

-norm 

minimization, compressive sensing, sparse representation in high-dimensional geometry), and much recent 

progress (initiated in pivotal work by Candes, Romberg, and Tao) most of the methods being explored are 

limited to linear spaces. While linearity is a natural 

starting point, strong evidence suggests that image recovery cannot be done this way yet it nevertheless is 

within reach. It is the goal of this project to develop 

procedures for robust recovery of piecewise smooth, or 

otherwise constrained, functions from Fourier data. 

Analysis of Diffusion and Stress in 

Binary Alloys 

Geoffrey B. McFadden 

William Boettinger (NIST MSEL) 

The Kirkendall effect refers to the motion of lattice 

planes in a crystalline binary alloy that is caused by a 

difference in the rates of diffusion of the two species. 

The difference in diffusion rates leads to a stress-free 

stain in the sample, and the process is an example of 

the coupling between diffusion and elasticity during 

phase transformations. Current interest in stress generation in thin films due to diffusion processes has 

stimulated new interest in this classical problem. The 

effect can lead to measurable shape changes perpendicular to the direction of diffusion [1]. It is often 

detected by a shift of inert markers parallel to the diffusion direction; that is, deformation with pure axial 

displacement [2,3]. 

In current work we are analyzing the timedependent bending of a bimetallic strip due to interdiffusion in a binary alloy. In this case the deformation 

is due to diffusion across the short dimension of the 

strip, so that the deformation and diffusion are colinear. We have considered two cases: first, diffusion in 

a single-phase material, in which the results, which 

employ a Fourier method, couple simple beam theory 

and diffusion. The second case considers a two-phase 

material, in which the interface separating the two 

phases moves in response to the concentration gradients induced in the sample. In that case we employ a 

similartity solution for short times, when end effects of 

the sample are negligible, followed by a numerical 

solution using finite differences for the long time solution representing the approach to steady-state behavior 

(see Figure 19). The results are in good agreement with 

previously reported experiments [4]. 

[1] W. J. Boettinger, G. B. McFadden, S. R. Coriell, R. F. 

Sekerka, and J. A. Warren, A Model for the Lateral Deformation of Diffusion Couples, Acta Materialia 53

(2005), pp. 1995-2008. 

[2] W. J. Boettinger, J. E. Guyer, C. E. Campbell, and G. B. 

McFadden, “Computation of the Kirkendall Velocity 

and Displacement Fields in a One-dimensional Binary 

Diffusion Couple with a Moving Interface,” Proceedings of the Royal Society A 463 (2007), pp. 3347–3373. 

[3] J. A. Dantzig, W. J. Boettinger, J. A. Warren, G. B. 

McFadden, S. R. Coriell, and R. F. Sekerka, “Numerical 

Modeling of Diffusion-induced Deformation,” Metallurgical and Materials Transactions 37A (2006), pp. 

2701-2714. 

[4] W. J. Boettinger and G. B McFadden, ”Bending of a 

Bimetallic Beam Due to the Kirkendall Effect,” Journal 

of Phase Equilibrium and Diffusion, in press (2010). 

Figure 19. The solid curves show the dimensionless curvature δκ of 

the beam vs. the dimensionless time tDβδ2

 for various initial positions of the interphase boundary z/δ= 0.3, 0.5, 0.649, 0.7, where δ, 

where is the beam thickness and Dβ is the diffusion coefficient. The 

dashed lines are the initial approximation and the long-time limit. 
44 NIST Mathematical and Computational Sciences Division 

Modeling of Self-Similar Grain Size 

Distributions 

Geoffrey. B. McFadden 

C.S. Pande (Naval Research Laboratory) 

Polycrystalline solids typically consist of a large number of grains that differ in the orientation of their 

crystallographic axes. These grains have a distribution 

of sizes, and are separated by curved grain boundaries 

that in turn meet at grain boundary junctions. Driven 

by combinations of elasticity, diffusion, or capillary 

forces the grain boundaries migrate in time, and the 

distribution of grain sizes consequently evolves dynamically. A statistical description of the process is 

given by a grain size distribution function, which is 

assumed to depend on a mean grain radius and time. 

Grain growth under usual circumstances is known to 

approach a quasi-stationary distribution of grain sizes 

after a transient period. The quasi-stationary state in a 

wide variety of materials exhibits a scaling property, 

such that the grain size distribution has an invariant 

form when expressed in terms of the grain size scaled 

by its mean value. An accurate description of the 

spatial and temporal evolution of a polycrystal from an 

initial stage, through the transient period, and finally to 

the quasi-stationary state is still only poorly understood. In this work, a size-based continuum stochastic 

formulation is presented based on topological considerations. This analysis leads to a Fokker-Planck 

equation for the size distribution which admits a similarity solution that describes the long time behavior of 

the grain growth. The solution for the case of twodimensional grain growth recently appeared [1]. The 

three-dimensional version of this problem is has been 

successfully addressed this year [2]. The resulting grain 

size distributions are in good agreement with those 

obtained from computer simulations (see Figure 20), 

indicating the validity of the stochastic approach. 

[1] C. S. Pande, K. P. Cooper, and G. B. McFadden, “Grain 

Size Distribution in Two Dimensions in the Long Time 

Limit,” Acta Materialia 56 (2008), pp. 5304–5311. 

[2] C. S. Pande and G. B. McFadden, “Self-Similar Grain 

Size Distribution in Three Dimensions,” Acta Materialia, in press.

Modeling Fluid Flow in Materials 

Processing 

Geoffrey B. McFadden 

P. Aaron Lott 

Sam R. Coriell (NIST MSEL) 

Daniel Anderson (George Mason University) 

Bruce Murray (SUNY Binghampton) 

The study of the stability of a fluid-fluid interface is 

important in a number of scientific and technological 

applications. In this project we consider two fluid layers separated by a horizontal planar interface subject to 

a vertical temperature gradient. The effects of various 

driving forces on the stability of the system can be 

taken into account, including buoyancy (known as 

Rayleigh-Benard convection), the effects of bulk density differences (known as Rayleigh-Taylor 

instabilities), and the effects of surface tension gradients along the interface (known as Marangoni 

instabilities). If the two layers represent different 

phases of the same material the stability results for a 

two-phase bi-layer system are quantitatively and even 

qualitatively different than for those of an immiscible 

system. To estimate the relative importance of these 

types of instabilities for a two-phase system, we have 

previously considered a bi-layer geometry in which a 

horizontal fluid-fluid interface separates two semiinfinite layers of a single component fluid [1,2]. Recently we studied [3] a novel version of the two-phase 

system in which the onset of the instability is oscillatory in time. This instability was first noted by 

researchers at the University of Minnesota, although 

the physical origins of the instability are poorly understood. We were able obtain both numerical and 

analytical approximations which help explain the 

physical mechanism of the instability. This year we 

have worked on an extension of the two-phase, bi-layer 

problem to the case of a binary alloy system [4]. This 

work is relevant to experimental studies of “monotectic” systems (in our case the aluminum-indium system) 

in which a liquid parent phase undergoes a phase transition leading to two daughter phases, one a liquid and 

one a solid. The geometry of the two-phase system is 

often composed of rod-like inclusions of one phase 

within the other, forming a composite that has useful 

technological applications. An open problem is to pre0.0 0.5 1.0 1.5 2.0 0.0

0.2

0.4

0.6

0.8

1.0

 Normalized Grain Size

 Normalized Grain Size Distribution

 Wakai et al.

α = 7

α = 8

 Hillert

Figure 20. Comparison of numerical simulation results (data 

points) with three self-similar grain size distributions from Ref. 

[2]. 
Yearly Report, Fiscal Year 2009 45

dict the inter-rod spacing in such systems, which is 

found to be sensitive to the flow conditions at the liquid-liquid phase boundary. In another related study, the 

effect of buoyant convection on the stability of a 

mushy zone in a ternary alloy has been considered [5]. 

The mushy zone is an idealized model of a two-phase 

region where liquid and solid co-mingle. References 

[4] and [5] are invited papers for a special edition of 

the Journal of Fluid Mechanics celebrating the seventieth birthday of its editor, S. H. Davis, who has served 

as the Ph.D. advisor for several MCSD postdocs. G. 

McFadden is acting as a guest co-editor for this volume. 

[1] G. B. McFadden, S. R. Coriell, K. F. Gurski and D. L. 

Cotrell, “Onset of Convection in Two Liquid Layers 

with Phase Change,” Physics of Fluids 19 (2007), 

104109. 

[2] G. B. McFadden, S. R. Coriell, K. F. Gurski and D. L. 

Cotrell, “Convective Instabilities in Two Liquid Layers,” Journal of Research of NIST 112 (2007), pp. 271-

281. 

[3] G. B. McFadden and S. R. Coriell, “Onset of Oscillatory 

Convection in Two Liquid Layers with Phase Change,” 

Physics of Fluids 21 (2009), 034101. 

[4] G.B. McFadden, S.R. Coriell, and P.A. Lott, “Onset of 

Convection in Two Layers of a Binary Liquid,” Journal 

of Fluid Mechanics, in press (2010). 

[5] D. M. Anderson, G. B. McFadden, S. R. Coriell, and B. 

T. Murray, “Convective Instabilities During the Solidification of an Ideal Ternary Alloy in a Mushy Layer,” 

Journal of Fluid Mechanics, in press (2010). 

Materials Data and Metrology for 

Applications to Machining Processes, 

Frangible Ammunition, and Body 

Armor 

Timothy Burns 

Steven Mates (NIST MSEL) 

Richard Rhorer (NIST MEL) 

Eric Whitenton (NIST MEL) 

Debasis Basak (Orbital Sciences Corporation) 

High-speed metal machining processes can cause extremely rapid heating of the workpiece material. If this 

material is a carbon steel, a small region of the workpiece is deformed plastically at a strain rate on the 

order of 10,000 s-1, on a time interval on the order of 

one millisecond. On this time interval, the material 

undergoes a change in temperature from ambient (~25 

oC) to a temperature on the order of magnitude of 1000 

oC, as a result of dissipation due to rapid shear in the 

primary and secondary shear zones. This means that, in 

a fairly routine cutting operation on a modern machining center, a heating rate on the order of one million 

degrees centigrade per second is not uncommon for 

iron-carbon alloys of interest in manufacturing; see 

Figure 21 [1]. This poses a challenge for modeling the 

constitutive response, i.e., the material flow stress vs. 

the strain, strain rate, temperature, etc., of such materials, for use in computer simulations of machining 

operations. 

As part of a research program on the measurement and modeling of machining operations, we have 

developed the NIST Pulse-Heated Kolsky Bar Laboratory for the study of material response at high 

temperatures and rates of loading [2]. The Laboratory 

is a unique testing facility in which the strength of metals can be measured under conditions of rapid DC 

resistive heating, followed by rapid compression. The 

Figure 21. Thermal image of high-speed orthogonal machining of

AISI 1045 steel; temperature is a maximum along the chip-tool interface (A-B) due to friction. 

Figure 22. Pulse-heated Kolsky bar data on shear flow stress vs. 

temperature, at a true compressive strain of 0.1, and a true strain 

rate of 3500 s-1 [3]; sample microstructures were a uniform fine 

pearlite prior to testing; error bars denote 2σ [2]; upper left data 

point corresponds to room-temperature test in which sample was not

pre-heated. 
46 NIST Mathematical and Computational Sciences Division 

flow stress can be measured in metals that have been 

pre-heated to uniform temperatures of up to 1200 oC, at 

heating rates up to 6000 oC s

-1, and then dynamically 

loaded at strain rates up to the order of 10,000 s-1. The 

facility is also being used for more conventional roomtemperature testing of materials, such as frangible ammunition made from sintered powders, in support of 

the NIST Office of Law Enforcement Standards 

(OLES). 

In recent work [3], new pulse-heated compression 

test results on AISI 1075 steel were presented, in which 

a rapid phase transformation from pearlite to austenite 

was observed to take place near the austenization 

temperature (723 oC) of the material. We are currently 

studying the implications of this behavior for finiteelement simulations of high-speed machining of this 

material. 

[1] M. A. Davies, H. Yoon, T. L. Schmitz, T. J. Burns, and 

M. D. Kennedy, “Calibrated Thermal Microscopy of 

the Tool-Chip Interface in Machining,” Journal of Machining Science and Technology 7(2) (2003), pp. 166-

190. 

[2] S. P. Mates, R. L. Rhorer, E. P. Whitenton, T. J. Burns, 

and D. Basak, “A Pulse-heated Kolsky Bar Technique 

for Measuring Flow Stress of Metals Subjected to High 

Loading and Heating Rates,” Experimental Mechanics 

48 (2008), pp. 799-807. 

[3] T. J. Burns, S. P. Mates, R. L. Rhorer, E. P. Whitenton, 

and D. Basak, “Dynamic Properties for Modeling and 

Simulation of Machining: an Update of Results from 

the NIST Pulse-heated Kolsky Bar,” in 12th CIRP Conference on Modeling of Machining Operations, 2009. 

San Sebastian, Spain. 

Boundary Slip Effects on the Stability of Spiral Poiseuille Flow 

P. Aaron Lott 

Geoffrey B. McFadden 

David L. Cotrell (Lawrence Livermore National Lab) 

Advances in micro-fabrication techniques have 

provided significant understanding into the behavior of 

fluids at small scales. Small scale experiments, with 

flow driven by pressure gradients, drainage, shear, or 

an electric field have reported a breakdown of the noslip boundary condition, with observed slip lengths as 

large as microns. In addition, numerical studies of 

fluids are often performed with no-slip boundary 

conditions, and the effects of introducing slip as a 

model of surface roughness is of current interest. 

 We consider the effect boundary slip has on the 

stability of the spiral Poiseuille flow, which is an annular flow with an applied pressure gradient [1]. We have 

found that under certain conditions, slip has a moderate 

effect on the linear stability of the flow. Numerical 

results indicate that linear stability is sensitive to the 

ratio between the inner and outer cylinder radii, as well 

as the slip length. Changes in these parameters cause 

the flow to be stabilized or destabilized in reference to 

a flow where no-slip conditions are applied; this 

change in stability is indicated by change in critical 

Taylor number. The dependence of stability on radius 

ratio becomes more complicated as the Reynolds number is increased; in particular, the flow becomes 

destabilized by small ratios, stabilized by mid-range 

ratios and de-stabilized again for large ratios. We have 

computed estimates of the critical radius ratios where 

slip has no effect on the stability of the flow. We are 

currently examining whether changes in the angular 

velocity of the inner and outer cylinders exhibits similar variations in stability of the flow. 

[1] D.L. Cotrell, P.A. Lott and G.B. McFadden, “Boundary 

Slip Effects on the Stability of Spiral Poiseuille Flow,” 

in preparation. 

Nanostructures, Nanooptics, and 

How to Squeeze the Light out of 

Quantum Dots 

James S. Sims 

John Kelso 

Garnett W. Bryant (NIST PL) 

Research and development of nanotechnology, with 

applications ranging from solar cells and light emitting 

diodes to smart materials to quantum computation to 

biolabs on a chip, is of high national priority. Semiconductor nanoparticles, also known as nanocrystals 

and quantum dots, are one of the most intensely studied 

nanotechnology paradigms. Nanoparticles are typically 

1 nm to 10 nm in size, with 103

 to 106 atoms. Precise 

control of particle size, shape, and composition allows 

one to tailor charge distributions and control quantum 

effects to tailor properties completely different from 

the bulk or from small clusters. As a result of enhanced 

quantum confinement effects, nanoparticles act as artificial man-made atoms with discrete electronic spectra 

that can be exploited as light absorbers for solar cells, 

light sources for novel enhanced lasers and light emitting diodes, discrete components in nanoelectronics, 

qubits for quantum information processing, and enhanced ultrastable fluorescent labels for biosensors to 

detect, for example, cancers, malaria or other pathogens, and to do cell biology. 
Yearly Report, Fiscal Year 2009 47

We are working 

with the NIST Physics 

Lab to develop computationally efficient large 

scale simulations of 

such nanostructures. We 

are also working to develop real-time 

(immersive and desktop) 

visualization techniques 

and tools to enable 

analysis of highly complex computational 

results of this type. In 

particular, we are studying the electrical and 

optical properties of 

semiconductor 

nanocrystals and quantum dots. In the most 

complex structures this 

entails modeling structures with on the order 

of a million atoms. 

Highly parallel computational and visualization platforms are critical for obtaining the computational speeds required for a 

systematic, comprehensive study. 

While it often is easy to define the simple subsystems that make up a complex, heterogeneous 

nanosystem [1], it may be difficult to explicitly define 

an entire large-scale structure. A novel feature of our 

code is the ability to link together heterogeneous nanosubsystems. Since we can do multiprocessor runs 

routinely we have the basic building blocks for making 

larger structures by “stitching” together disparate subsystems into composite structures, each separate 

component subsystem being a smaller multiprocessor 

run. We are currently focusing on modulating and controlling the optical properties of self-assembled 

quantum dots using external strain [2-5]. 

 To understand the coupling between externally imposed strains and the electronic states of the selfassembled quantum dots, we have developed an atomistic tight-binding theory of the confined states in 

quantum dots that incorporates the local strain due to 

atomistic lattice mismatch and the externally imposed 

strain due to surface acoustic waves to applied stressors 

or to the bend in a vibrating nanomechanical oscillator. 

Both local strain and the externally imposed strain are 

important, so we include them on an equal footing via 

an atomistic valence force field approach, with the externally imposed strain modeled as a distortion of the 

lattice on the boundary of structure. A full tightbinding model including an sp3s* orbital model and 

spin-orbit effects is used [2-5]. To understand how 

applied stress can be used to actively control dot optical properties, dots buried at points with 

different strain in a 

nanobridge oscillator or 

cantilever are considered. The dependence 

of the quantum dot electronic states and optical 

transitions on the coupling to bending modes 

of the oscillator is studied. Both isolated and 

coupled dots are considered. Applied strain can 

be used to transfer electrons and holes between 

vertically or laterally 

coupled dots. Level 

crossings for electron 

and hole states are 

found. State symmetries 

can be significantly distorted by the applied 

stress. Transitions can 

be strongly polarized 

parallel or perpendicular to the external strain or suppressed, depending on how the dots are stressed. The 

dependence on applied stress can be seen as a competition between the internal and external stress that can 

either enhance or suppress the local strain at the dot, 

depending on how the external stress is applied. 

 We have found that the applied strain can make 

large changes in the exchange splitting between bound 

electron-hole pairs (excitons) optically bright states to 

modify the exciton fine structure induced by QD 

asymmetry or atomistic effects, can induce crossing 

between the bright states for certain strains, and can be 

used to rotate the polarized response of the bright 

states. These capabilities should be useful for applications such as entangled photon sources, where the 

elimination of asymmetric exchange splitting is essential, or cavity-coupled QDs, where control of the 

exchange-induced dark and bright states is needed. 

Analysis shows that the applied strain tunes the size 

and phase of the exchange coupling which, in turn controls the phase of the exciton state and the resulting 

polarization. Such phase control is critical for applications in quantum science. 

 As the computational model has been extended to 

handle more complex and larger systems by including 

not only the nanostructures but also the substrate and 

environment around them, parallel processing has become even more necessary. We use both shared 

memory OpenMP parallelization and distributed memory message passing. 

 To study the increasingly larger systems we are 

working with, we have also expanded the visualization 

Figure 23. The nano visualization interface.
48 NIST Mathematical and Computational Sciences Division 

capability with multiple selection tools that can sum up 

charge density in different types of regions. This year 

we made improvements in: speed – much faster building, loading and higher graphic frame rate for a given 

number of orbitals; capacity – a million orbitals can be 

displayed with an interactive frame rate; and dynamic 

updates – all displays are dynamically updated as orbitals and selection criteria are modified. Dynamic 

updates give immediate feedback to the researcher. 

 New capabilities in the visualization (see Figure 23 

for a snapshot of the interface) include: 

• Selection of which orbitals to display. Orbitals are 

currently represented using “fatpoints,” i.e., 

spheres drawn by a GPU shader program. A P orbital is represented by two spheres and an S orbital 

by one sphere. The spheres are colored and sized 

as specified by the orbital data. 

• Dynamic visual updating of which orbitals lie inside or outside specified geometric objects. 

• Tools to dynamically exclude orbitals that are 

greater or less than a specified size. 

• Assistance of a gnomon object to display orbital 

orientation, the option to display all orbitals at a 

fixed size, and to magnify all orbitals. 

• A dynamically displayed status window showing 

the number of visible orbitals, and the exact position and orientation of the selection geometry. 

• Dynamic summing of visible orbitals, with optional orbital coefficients. 

[1] G. W. Bryant and W. Jaskolski, “Electronic Properties 

of Quantum-dot Molecules,” Physica E 17 (2003), pp. 

40-41. 

[2] G. W. Bryant, N. Malkova, J. Sims, M. Zielinski, W. 

Jaskolski and J. Aizpurua, “Re-engineering the Optics of 

Quantum Dots with Mechanical Strain,” in Optics of 

Excitons in Confined Systems (OECS11), Madrid, Spain, 

September 7-11, 2009. 

[3] G. W. Bryant, N. Malkova, J. Sims, M. Zielinski, W. 

Jaskolski and J. Aizpurua, “Re-engineering the Optics of 

Quantum Dots with Mechanical Strain,” in Proceedings 

of the Electronic Materials Conference, University Park, 

PA, June 24-26, 2009. 

[4] G. W. Bryant and N. Malkova, “Engineering the Optics 

of Quantum Dots Using Strain,” in Proceedings of the 

2009 Nanoelectronic Devices for Defense & Security 

(NANO-DDS) Conference, Fort Lauderdale, FL, September 21-25, 2009. 

[5] N. Malkova, G. W. Bryant, and J. Sims, “Straining 

Quantum Dot Excitons,” in Proceedings of the Materials 

Research Society Fall Meeting, Boston, MA, November 

30-December 4, 2009. 

[6] J. Sims, G. W. Bryant, W. L. George, T. G. Griffin, J. G. 

Hagedorn, J. T. Kelso, T. M. Olano, A. P. Peskin, S. G. 

Satterfield, and J. D. Terrill, “Computation and Visualization of Nano-structures and Nano-optics,” Computing 

and Visualization in Science. 

Computational Modeling and Visualization of Cement Paste Hydration 

and Microstructure Development 

William George 

Steven Satterfield 

Judith Terrill 

Jeffrey Bullard (NIST BFRL) 

Edith Enjolras (NIST BRFL) 

Clarissa Ferraris (NIST BRFL) 

Edward Garboczi (NIST BFRL) 

Nicos S. Martys (NIST BFRL) 

Paul E. Stutzman (NIST BRFL) 

When cement powder is mixed with water, the hydration process that transforms the paste from a fluid 

suspension into a hardened solid involves complex 

chemical and microstructure changes. Understanding 

and predicting the rates of these changes is a longstanding but extremely challenging technological goal. 

Fundamental computational modeling of the hydration 

of cement is difficult because it involves a large number of coupled nonlinear rate equations that must be 

solved in an irregular three-dimensional spatial domain. 

To address these challenges, we are applying a 

new computational model called HydratiCA, which has 

several advantages over previous attempts to model 

cement paste hydration. HydratiCA is based on the 

principles of mass action and detailed balances and 

uses stochastic cellular automaton algorithms to simultaneously model 3-D reaction and transport 

phenomena. This allows us to track the detailed kinetics and equilibria that occur in a diverse range of 

cementitious systems. At the length scales required to 

finely resolve the reaction mechanisms and microstructure changes in cement paste, HydratiCA must take 

small time steps (approximately 10-5 per second) to 

remain numerically stable. In other words, tens of millions of time steps are required to simulate the behavior 

of cement paste for just one hour. Therefore, parallelization of the model is important so that we can 

simulate systems that are large enough to be realistic, 

avoiding finite size effects, and still be able to complete the simulations in a reasonable amount of time. 

Visualization of the output of this model is important 

both for validation and to understand the distribution of 

phases in 3D. The output is a 3D volume of data with 

percentage values for each of multiple material phases 

at each lattice site. Over the course of the simulation 

time, a series of data volumes is produced at the time 

values of interest. 

This year we modified the HydratiCA simulator 

to use XML encoded parameter files. These changes 

have simplified the automatic validation of input and 

has allowed for the easy development of a GUI inter-
Yearly Report, Fiscal Year 2009 49

face for creating and editing such files. We have improved the checkpoint/restart handling of HydratiCA 

such that no post-processing or pre-processing of files 

is needed when restarting a simulation. This greatly 

simplifies the job of running long simulations within 

batch systems that enforce relatively short runtimes 

compared to the runtime of a complete simulation. 

 We were awarded 540,000 CPU-hours on the 

NASA Advanced Supercomputing machine Columbia 

for this project. This time will be used to run HydratiCA for the more compute-intensive simulations 

required in this project. A sample visualization of our 

results is shown in Figure 24. We have completed a 

large number of validation runs and written a paper on 

the result, which has been accepted for publication. We 

will be using HydratiCA to help resolve two prominent 

hypotheses for reduced alite dissolution rates during 

the induction period. 

[1] J. W. Bullard, E. Enjolras, W. L. George, S. G. 

Satterfield, and J. E. Terrill, “Validating a Scalable Parallel Computer Model of Cement Hydration and 

Microstructure Development,” Modeling and Simulation 

in Materials Science and Engineering, accepted. 

[2] J. W. Bullard, C. F. Ferraris, E. J. Garboczi, N. S. Martys, P. E. Stutzman, and J. E. Terrill, “Virtual Cement 

and Concrete,” Innovations in Portland Cement Manufacturing, in press. 

Micromagnetic Modeling 

Michael Donahue 

Donald Porter 

Robert McMichael (NIST CNST) 

June Lau (NIST MSEL) 

http://math.nist.gov/oommf/

Advances in magnetic devices such as recording heads, 

field sensors, magnetic nonvolatile memory (MRAM), 

and magnetic logic devices are dependent on an understanding of magnetization processes in magnetic 

materials at the nanometer level. Micromagnetics, a 

mathematical model used to simulate magnetic behavior, is needed to interpret measurements at this scale. 

MCSD is working with industrial and academic partners, as well as with colleagues in NIST CNST, MSEL, 

and EEEL, to improve the state-of-the-art in micromagnetic modeling. 

 Michael Donahue and Donald Porter have developed a widely used public domain computer code for 

doing computational micromagnetics, the ObjectOriented Micromagnetic Modeling Framework 

(OOMMF). OOMMF serves as an open, welldocumented environment in which algorithms can be 

evaluated on benchmark problems. OOMMF has a 

modular structure that allows independent developers 

to contribute extensions that add to the basic functionality of OOMMF. OOMMF also provides a fully 

functional micromagnetic modeling system, handling 

both two and three-dimensional problems, with sophisticated extensible input and output mechanisms. In 

fiscal year 2009 alone, the software was downloaded 

more than 5,200 times, and use of OOMMF was acknowledged in more than 120 peer-reviewed journal 

articles. OOMMF has become an invaluable tool in the 

magnetics research community. 

 Key developments in OOMMF last year include: 

1) Addition of multi-threaded computation in the 

modules of OOMMF for use on multi-processor, 

multi-core computers, enabling faster (observed 

speed-up of four to seven times) and larger simulations, and more efficient utility of common 

computing platforms [1]. 

2) Creation of additional modules tailored for the 

simulation of polycrystalline materials, common 

in practical device designs [4]. 

 In addition to the continuing development of 

OOMMF, the project also does collaborative research 

using OOMMF [2,3]. The project plays an instrumental role in the 2007 Innovation in Measurement Science 

project “Magnetic nanostructures for post-CMOS electronics,” with CNST, MSEL and EEEL. We also 

provide modeling and simulation consultation for work 

in EEEL on spin momentum transfer. In collaboration 

Figure 24. Snapshot from a visualization of formation of slaked 

lime. The number in the figure is generated with the virtual cement 

analysis probe and indicates the concentration of the Ca++ ion at 

the point marked by the 3D cross glyph. 
50 NIST Mathematical and Computational Sciences Division 

with the University of Konstanz we have used 

OOMMF to study magnetic structures in nanowires. 

The MCSD project produced four journal papers and 

four conference presentations this past year. 

 Work in progress includes a continuing effort to 

develop techniques to minimize or control for errors 

introduced when the spatial resolution of the computation grid does not exactly represent the boundaries of 

the simulated material. Further progress on multithreaded computation is focused on designing computations to more effectively use limited memory 

bandwidth on non-uniform memory access architectures. Another focus area is the use of OOMMF to 

explore spintronic designs for magnetic memory and 

logic devices, including a new collaboration with University of Alabama. 

[1] M. J. Donahue, “Parallelizing a Micromagnetic Program 

for Use on Multi-processor Shared Memory Computers,” IEEE Transactions on Magnetics 45 (2009), pp. 

3923-3925. 

[2] W. C. Uhlig, M. J. Donahue, D. T. Pierce and J. Unguris, “Direct Imaging of Current-driven Domain Walls 

in Ferromagnetic Nanostripes,” Journal of Applied 

Physics 105 (2009), 103902. 

[3] M. D. Stiles, W. M. Saslow, M. J. Donahue, and A. 

Zangwill, ‘Reply to “Comment on ‘Adiabatic domain 

wall motion and Landau-Lifshitz Damping,” Physical 

Review B 78 (Dec. 2008), 216402. 

[4] J. W. Lau, R. D. McMichael, and M. J. Donahue, “Implementation of Two-dimensional Polycrystalline Grains 

in Object Oriented Micromagnetic Framework,” Journal 

of Research of the National Institute of Standards and 

Technology 114(Jan 2009), pp. 57-67. 

OOF: Finite Element Analysis of 

Material Microstructures 

Stephen A. Langer 

Andrew C.E. Reid 

Valerie R. Coffman 

Gunay Dogan (Cornell University) 

http://www.ctcms.nist.gov/oof/

The OOF Project, a collaboration between MCSD and 

the NIST Material Science and Engineering Laboratory, is developing software tools for analyzing real 

material microstructure. The microstructure of a material is the (usually) complex ensemble of 

polycrystalline grains, second phases, cracks, pores, 

and other features occurring on length scales large 

compared to atomic sizes. The goal of OOF is to use 

data from a micrograph of a real material to predict the 

macroscopic behavior of the material via finite element 

analysis. OOF is intended to be a general-purpose tool, 

applicable to a wide variety of microstructures in a 

wide variety of physical situations. 

 During FY 2009 the OOF team concentrated on 

improving the non-linear and time dependent solvers in 

the 2D program, OOF2, and on developing the 3D version, OOF3D. Numerous bugs and glitches delayed the 

release of OOF2 version 2.1, which will be the first 

version to solve time dependent problems. There were, 

however, three bug-fix releases of 2.0.5 this past year. 

Including all versions, OOF2 was downloaded 1,204 

times from the OOF website this past year (10% more 

than last year, during which seven minor versions were 

released). OOF2's predecessor, OOF1, is still available 

and was downloaded 265 times. This is significantly 

less than last year, which is good. 

As usual, the OOF customer base is broad and difficult to define, since users aren't required to register 

before using the software. 42 users wrote to the developers for help on installing and running the program 

this year. OOF2 is still available to be run on-line at 

Purdue University's nanoHUB site, where it had 259 

users in 2,420 interactive sessions this year, up about 

50% from last year. 

 Externally visible progress on the time-dependent 

and three dimensional versions has been slow this year. 

However, almost everything that was not working at 

the end of last year has now been fixed. 

 OOF3D can now solve simple 3D diffusion problems 

from start to finish, beginning with a 3D image and 

ending with a finite element solution. It can load orientation data from 3D electron backscatter diffraction 

(EBSD) experiments. More 3D mesh refinement methods have been added. Most of the 3D development 

effort this year went towards optimization and corner 

cases in the element homogeneity calculation, which 

computes the fraction of each mesh element that overlaps each voxel type. This surprisingly subtle 

Figure 25. Preliminary data shows that OOF2 can achieve up to 

a 100 fold reduction in the number of elements required to mesh a 

3D micrograph of trabecular bone, compared to the simple approach of creating one element per voxel. 
Yearly Report, Fiscal Year 2009 51

calculation is fundamental to many OOF computations, 

so it needs to be fast and robust. 

 OOF2 gained a module for solving RambergOsgood plasticity, which is a form of nonlinear elasticity that can be used to model the plastic deformation of 

materials. The nonlinear solvers were restructured to be 

easier to use and to give users more control over the 

solution process. The Newton and Picard solvers now 

share infrastructure. 

 It is now possible to animate the solutions of timedependent problems. This and other aspects of the 

time-dependence code required major under-the-hood 

changes in the way in which OOF2 handles progress 

bars and manages execution threads. 

 The OOF developers hope that version 2.1, with 

time-dependence, and version 3.0, with 3D solutions, 

will both be made available in FY 2010. 

[1] A. C. E. Reid, R. C. Lua, R. E. Garcia, V. R. Coffman, 

and S. A. Langer, “Modeling Microstructures with 

OOF2,” International Journal of Materials and Product 

Technology 35 (2009), pp. 361-373. 

Fast Solvers for Models of Fluid 

Flow 

P. Aaron Lott 

Geoffrey B. McFadden 

Sam Coriell (NIST MSEL) 

Howard Elman (University of Maryland) 

Katherine Evans (Oak Ridge National Laboratory) 

Numerical simulation of fluids enables improved prediction and design of a wide range of natural and 

engineered systems. Fluid systems often involve dynamics that occur on disparate length and time scales 

due to variations in inertial and viscous forces. In order 

to numerically simulate these dynamics, advanced 

mathematical modeling and scalable computational 

methods are required. The principal goal of this project 

is to develop fast numerical methods for conducting 

high resolution fluid flow simulations. 

 In numerical fluid simulations, the majority of time 

is spent by linear system solvers that reside within the 

inner loop of the simulation. Reducing the cost of these 

linear system solves is a challenging problem in computational mathematics. We are developing 

preconditioners to accelerate the convergence of iterative linear system solvers for linear systems that arise 

from spectral element discretizations. Preconditioners 

are chosen as approximations of the linear system that 

mimic system properties, but provide an efficient 

means of applying the inverse operation. We take advantage of underlying structure, multiscale features and 

parallelism present in the discrete system and use domain decomposition coupled with a fast 

diagonalization procedure to efficiently apply the preconditioner. These techniques have resulted in 

significant speedup in simulating steady flow models 

including convection-diffusion systems [1] and incompressible fluids [2]. We are working to develop 

improved preconditioning strategies that will enable 

more challenging flow problems to be solved at higher 

resolution and on modern supercomputers. 

 We are also extending the use of our preconditioning strategies to unsteady flow models including 

thermosolutal convection and the shallow-water equations. Many numerical methods for fluid simulation are 

based on explicit and semi-implicit time integration 

schemes. However, in order for these methods to be 

numerically stable, the maximum time step is required 

to be bounded by the minimum grid space divided by 

the maximum flow velocity. This causes time steps to 

become prohibitively small when simulating fast moving flows at high resolution. As modern computing 

facilities are increasing the ability to perform higher 

resolution studies, numerical methods must be developed to provide a scalable means of performing these 

simulations. Fully implicit methods would allow for 

significantly larger time steps that are independent of 

the grid resolution. Similar to steady state fluid simulations, the cost of fully implicit time integration is 

dependent on the expense of subsidiary linear system 

solves. By developing efficient preconditioning strategies, we are working to lower the cost of these linear 

system solves to enable scalable high-resolution flow 

studies. 

[1] P. A. Lott and H. C. Elman, “Fast Iterative Solvers for 

Convection-diffusion Systems with Spectral Elements,” 

Numerical Methods in Partial Differential Equations, 

accepted. 

[2] P. A. Lott and H. C. Elman, “Fast Iterative Solvers for 

Incompressible Navier-Stokes Equations with Spectral 

Elements,” in preparation. 

Parallel Adaptive Refinement and 

Multigrid Finite Element Methods 

William F. Mitchell 

Marjorie A. McClain 

Eite Tiesinga (NIST PL) 

Thomas Hanna (NIST PL) 

http://math.nist.gov/phaml

Finite element methods using adaptive refinement and 

multigrid techniques have been shown to be very efficient for solving partial differential equations. 

Adaptive refinement reduces the number of grid points 

by concentrating the grid in the areas where the action 

is, and multigrid methods solve the resulting linear 
52 NIST Mathematical and Computational Sciences Division 

systems in an optimal number of operations. Recent 

research has focused on hp-adaptive methods where 

adaptivity is in both the grid size and the polynomial 

order of approximation, resulting in exponential rates 

of convergence. W. Mitchell has been developing a 

code, PHAML, to apply these methods on parallel 

computers. The expertise and software developed in 

this project are useful for many NIST laboratory programs, including material design, semiconductor 

device simulation, and the quantum physics of matter. 

 This year the effort was focused on two areas. The 

first is preparation for an experiment to study the effectiveness of several competing strategies for hp-adaptive 

refinement, which involves implementing these methods in PHAML and defining a set of test problems. 

The second is the application of PHAML to solve the 

time-dependent Schrödinger Equation for eigenvalues 

and eigenstates relevant to optical traps for neutral atoms, in collaboration with the Quantum Processes 

group of NIST's Atomic Physics Division. Understanding the interactions of adjacent atoms corresponding to 

qubits of a quantum gate involves computing multiple 

eigenvalues in the middle of the spectrum, with eigenstates that have sharp gradients, which is a very 

challenging computation. Major accomplishments for 

FY 2009 include the following. 

• Implementation of the final three hp-adaptive 

strategies in PHAML. 

• Completion of preliminary experiments to determine parameters to use in the hp-adaptive 

strategies comparative experiment. 

• Completion of the definition and implementation 

of a set of test problems. 

• Improvement of the MPI communication layer of 

PHAML, to improve parallel performance. 

• Extension of the class of problems solved by 

PHAML to include first order and cross derivative 

terms, and many other minor improvements. Two 

releases of the code occurred as the code evolved. 

• Design and implementation of a method to solve 

the time-dependent Schrödinger equation with 

PHAML. 

 Future work will continue to enhance PHAML with 

additional capabilities and robustness, complete the 

experiment to compare the effectiveness of hp-adaptive 

strategies, add OpenMP directives to PHAML for a 

hybrid MPI/OpenMP parallel implementation for clusters of multi-core computers, apply PHAML to a timedependent model of two dipolar molecules held in a 

cylindrically symmetric trap, and perform further physics experiments using the code. 

[1] W. F. Mitchell, “The hp-Multigrid Method Applied to 

hp-Adaptive Refinement of Triangular Grids,” Numerical Linear Algebra with Applications, to appear. 

[2] W. F. Mitchell and M. A. McClain, “A Survey of hpAdaptive Strategies for Elliptic Partial Differential 

Equations,” Annals of the European Academy of Sciences, invited paper, submitted. 

[3] W.F. Mitchell, “hp-Adaptive Finite Elements for the 

Schrödinger Equation,” University of Nevada at Reno, 

April 2009, invited talk. 

[4] W. F. Mitchell, “hp-Adaptive Finite Elements in a Parallel PDE Solver,” SIAM Conference on Computational 

Science and Engineering, March 2009, invited. 

[5] W. F. Mitchell, “Application of a Parallel Adaptive Finite Element Code to Confined Interacting Atoms,” 

SIAM Annual Meeting, July 2008. 

. 
Yearly Report, Fiscal Year 2009 53

Information Discovery, Use and Sharing

Society is awash in data – our ability to amass data has outpaced our ability to use it. Extracting knowledge, 

information, and relationships from this data is one of the greatest challenges faced by the scientists in the 

twenty-first century. The data can be as diverse as biological research data, medical images, automated 

newswire, speech, or video. The Information Discovery, Use, and Sharing Program fosters innovation 

throughout the information life cycle by developing the measurement infrastructure to enhance knowledge 

discovery, information exchange, and information usability. The Program enables novel computational approaches to data collection and analysis to be combined with improved interoperability techniques to 

effectively extract needed information from the wealth of available data. – Mary Brady, Program Manager. 

Mathematical Knowledge Management 

Bruce Miller 

Abdou Youssef 

Michael Kohlhase (Jacobs University, Bremen) 

Heinrich Stamerjohanns (Jacobs University, Bremen) 

See feature article, page 37. 

A Quality Pre-Processor for Biological Cell Images 

Adele Peskin 

Karen Kafadar (Indiana University) 

Alden Dima (NIST ITL) 

Anne Plant (NIST CSTL) 

John Elliot (NIST CSTL) 

Michael Halter (NIST CSTL) 

See feature article, page 31. 

Estimating Volumes of Simulated 

Lung Cancer Nodules 

David Gilsinn 

Amelia Tebbe 

Bruce Borchardt (NIST MEL) 

Lung cancer is a disease of uncontrolled cell growth in 

tissues of the lung. It is the most common cause of 

cancer related deaths in men and second most in 

women and is responsible for an estimated 560,000 

deaths in the US in 2007. Lung cancer lesions may be 

seen in Computerized Tomography (CT) images. CT is 

a medical imaging method that employs digital image 

processing techniques to generate three dimensional 

images from a large sequence of two-dimensional Xrays made around a single axis. 

 The Food and Drug Administration (FDA) is conducting research on the development of reference 

cancer lesions, called phantoms, to test CTs and their 

underlying proprietary software. Two sample phantoms 

were loaned to NIST to estimate their volumes. The 

material of which the phantoms were composed simulates lung cancer material. The FDA can insert these 

phantoms into a simulated body torso for testing CT 

scanning and image reconstruction. The two phantoms 

are shown in Figure 26. Although they seem spherical 

they are slightly non-spherical. The phantom on the 

right in Figure 26 is referred to as Green and the one on 

the left is referred to as Pink. 

The approach chosen to determine the volumes of 

the phantoms was based on a fundamental theorem in 

calculus, the Divergence Theorem. In this theorem a 

volume integral is shown to be equal to a particular 

surface integral. In order to develop a surface model of 

the phantoms, surface data was measured by a coordinate measuring machine (CMM) in the Manufacturing 

Engineering Laboratory (MEL) at NIST; see Figure 27. 

This machine produced a set of precisely measured (x, 

y, z) points on the surface of each phantom. As a reference, surface data was obtained on a well calibrated 

metal sphere in order to test the surface modeling process. In all cases, two sets of data points were measured, 

which we refer to as coarse (122 points) and dense 

(181 points). 

Figure 26. Physical phantoms used to simulate presence of lung 

cancers. Left: Pink phantom. Right: Green phantom. 
54 NIST Mathematical and Computational Sciences Division 

 The data were then transformed to spherical coordinates and modeled using a set of basis functions, called 

B-splines. The initial best fit, in the sense of least 

squares, of linear combinations of bivariate forms of 

the B-splines exhibited undesirable oscillations when 

the models were used to extrapolate values to points on 

the phantom surfaces. The least squares objective function was modified to include smoothing terms. This 

process, called regularization, produced very satisfactory extrapolation values at triangle vertices on the 

phantom surfaces. The surface triangulation was then 

used to compute the necessary surface integrals and 

finally the volumes using the Divergence Theorem. 

 In order to test the B-spline modeling process two 

approaches were taken to model the surface data from 

the calibrated metallic sphere. A spherical model was 

fit to the calibrated sphere by least squares parameter 

fitting and the B-spline model was fit to the sphere data 

using the regularization form of the objective function. 

The predicted volumes differed by approximately 

0.01%. This confirmed that the B-spline modeling 

process, using regularization, was a valid process and 

volume estimates of the phantoms would be acceptable. 

 As a check on the potential non-spherical shape of 

the phantoms a spherical model was fit to both of the 

phantoms. These produced volume estimates for the 

coarse and dense data from the Green phantom of 

4331.3 mm3

 and 4296.5 mm3 respectively. The expanded uncertainties for the center point and radius 

were estimated as order 10-1 and 10-2 respectively. For 

the Pink phantom the volumes were estimated as 

3892.1 mm3

 for the coarse data and 3912.5 mm3 for the 

dense data. The expanded uncertainties for the center 

point and radius were of order 1 and 10-1 respectively. 

This suggested a potentially spherical shape for the 

Green phantom and a slight non-spherical shape for the 

Pink phantom. When these results were estimated using the B-spline model, the volumes for the Green 

phantom were computed as 4331.8 mm3

 and 4296.8 

mm3

 for the coarse and dense data respectively. Notice 

the closeness to the spherical model fit results. The 

expanded volume uncertainties were estimated as approximately 12 mm3

 for the coarse data and 8 mm3 for 

the dense data. For the Pink phantom the volume results were 3854.7 mm3

 and 3906.5 mm3 with expanded 

uncertainties of approximately 28 mm3

 and 18 mm3. 

These differed slightly from the spherical model fit, 

again indicating a potential non-spherical form for the 

Pink phantom. The final report of this project has been 

submitted for publication in the Journal of Research of 

NIST. 

[1] D. E. Gilsinn, B. R. Borchardt, and A. Tebbe, “Estimating Volumes of Simulated Lung Cancer Nodules,” 

NISTIR 7571, July 2009. 

[2] D. E. Gilsinn, B. R. Borchardt, and A. Tebbe, “Volume 

Estimation of Molded Artifacts by B-Splines,” Proceedings of the SIAM Conference on Mathematics for 

Industry, 2009. 

Standardization of Lung Cancer 

Growth Measurements 

Adele Peskin 

Javier Bernal 

David Gilsinn 

Alden Dima (NIST ITL) 

Karen Kafadar (Indiana University) 

Lung tumor malignancy is assessed by measuring the 

increase in tumor size over time in computed tomography (CT) scans. Current methods for estimating tumor 

size from CT scans are inaccurate, with results varying 

greatly from lab to lab. Standardized and accurate 

methods for estimating tumor volumes are needed. In 

order to assess the quality of volume measurements, we 

are preparing synthetic lung tumors of known volume 

for insertion into clinical CT data. In this regard, it is 

important that we reconstruct the complex geometries 

of real tumors that make them challenging to measure. 

Figure 27. NIST Coordinate measureing machine. 
Yearly Report, Fiscal Year 2009 55

 Lung tumors are difficult to measure for a number 

of reasons. They are often small compared to the discretization of the CT data, so that large portions of the 

tumors only fill a fraction of the each voxel of the 3D 

data set, In addition, there are no standard methods to 

measure the partial volumes. Finally, the tumors are 

hard to separate from blood vessels that often surround 

them and from the pleural lining of the lung. By combining different geometric shapes, such as spheres, 

ellipsoids, cylinders, and cones, of different sizes and 

at different angles with respect to one another, we can 

recreate common shapes of lung tumor nodules inside 

clinical data at positions in the data that give rise to 

challenging volumetric issues. Several different tumor 

geometries are shown in Figure 28 that shows 3D isosurfaces of a region of lung data that contains both a 

clinical tumor (left) and several different synthetic tumors (right). 

[1] A. Peskin, K. Kafadar, A. Dima, J. Bernal, and D. 

Gilsinn, “Synthetic Lung Tumor Data Sets for Comparison of Volumetric Algorithms,” in Proceedings of the 

2009 International Conference on Image Processing, 

Las Vegas, NV, July 13-16, 2009. 

 

 

Figure 28. Isosurfaces of lung CT data taken at 400 Hounsfield units. From left to right: (a) The original lung CT data. (b) Synthetic tumor added 

at right. (c) Same synthetic tumor attached to the pleural lining of the lung using a cylindrical attachment. 
56 NIST Mathematical and Computational Sciences Division 

Complex Systems

Complex Systems are composed of large interrelated, interacting entities which taken together, exhibit a 

macroscopic behavior which is not predictable by examination of the individual entities. The Complex Systems Program seeks to understand the fundamental science of these systems and develop rigorous 

descriptions (analytic, statistical, or semantic) that enable prediction and control of their behavior. Initially 

focused on the Internet and Grid Computing, this Program will facilitate predictability and reliability in 

these areas and other complex systems such as biotechnology, nanotechnology, semiconductors, and complex 

engineering. – Sandy Ressler, Program Manager

Methods for Characterizing Complex Networks 

Isabel Beichl 

Brian Cloteaux 

David Harris (NSA) 

Francis Sullivan (IDA Center for Computing Sci.) 

Over the past decade, it has been realized that many of 

the networks inherent to natural and engineered systems share similar characteristics. These complex 

networks have been shown in systems such as the 

power grid, the Internet, communications and transportation networks, and software components. Since these 

networks are at the heart of our commercial infrastructure, it is important to develop a measurement science 

for such objects based on sound mathematical methods. 

Our project looks at both the modeling and measurement of complex networks of these types. 

 The development of predictive models of network 

characteristics is an important prerequisite. There have 

been a number of recent papers proposing possible 

mechanisms for explaining how complex networks 

could evolve. Unfortunately, the resulting models have 

not been able to realistically capture many of the structural characteristics that are observed in real life 

networks. Thus, we are currently unable to accurately 

model the time progression of complex networks. As 

an alternative, we have proposed a method for building 

random models through a completely structural approach [2]. 

 This structural approach does not try to mimic how 

the network evolves, but rather it generates random 

networks with a series of characteristics that we choose 

as being important to our model. In order for us to produce dynamic models of these complex networks we 

have the additional constraint that our models must 

contain a core that reflects the network from an earlier 

time period. In other words, our method produces 

models that simulate the evolution of a network over 

time while maintaining certain critical characteristics, 

but without having to specify the exact generating 

mechanism. This approach has shown considerable 

promise in tests to predict the evolution of the autonomous system topology of the Internet. 

 A second research direction has been in the measuring of characteristics of these networks. One result has 

been in the area of computing the reliability polynomial of a network [3]. The reliability polynomial 

computes the probability that a network stays connected given that each edge in the network has an 

identical and independent chance of failure. Since it 

has been shown that this problem is computationally 

difficult, we have produced a new algorithm for approximating this polynomial. When compared to the 

best known published algorithm, our method has 

shown a significant speedup in its computation time. 

 We have also proposed an enhanced algorithm for 

computing reliability using a new importance function. 

A new technique for reducing variance is also being 

investigated. This technique involves adding a term to 

the estimator whose expected value is zero. 

[1] I. Beichl and B. Cloteaux, “Measuring The Effectiveness of the s-Metric to Produce Better Network 

Models,” in Proceedings of the Winter Simulation Conference (2008), Miami, FL, December 7-10, 2008, pp. 

1020-1028. 

[2] I. Beichl and B. Cloteaux, “A Structural Approach to the 

Temporal Modeling of Networks,” in Proceedings of the 

2009 Winter Simulation Conference. 

[3] I. Beichl and B. Cloteaux, “An Approximation Algorithm for the Coefficients of the Reliability Polynomial,” 

Congressus Numerantium, accepted. 

Standard Reference Data for Complex Network Research 

Roldan Pozo 

http://math.nist.gov/~RPozo/complex_datasets.html

The study of complex or self-organizing networks is a 

relatively young, yet active research area. Much of the 

focus is to determine how graphs from real phenomena 

(such as protein interactions, social networks, structure 

of the World Wide Web, etc.) differ from synthetic or 
Yearly Report, Fiscal Year 2009 57

engineered graphs that are constructed by design. 

Clearly they are different, but how? Fundamental metrics, such as the clustering coefficient, graph diameter, 

and degree distributions, are commonly used in the 

literature, and these have served as encouraging first 

steps, but further research is needed to yield more 

meaningful distinctions and truly understand the nature 

of complex networks. 

 One important question to ask is: what are the precise real graphs that are being used to drive these 

theories? Clearly, the conclusions being reached can be 

only as meaningful as the quality and availability of the 

original data. Most papers provide little such information, making comparisons or verification of approaches 

difficult. For example, a commonly studied collaboration network in the literature is the Movie Actor graph, 

where two actors are connected if they have worked 

together on the same movie. While, this information is 

usually culled from the Internet Movie Database, there 

is no single Movie Actor graph because this database is 

a changing snapshot of the entertainment industry at 

any given time. Furthermore, what one chooses to include in the search (e.g., foreign films, documentaries, 

TV movies, specific epochs or genres) will generate 

different graphs. In fact, there are several versions of 

the Movie Actor graph cited in the literature, so a citation is not meaningful unless it contains a time-stamp 

and a record of the specific queries used. A similar 

argument can be made for other association networks, 

such as co-author or citation graphs. 

 Another example is networks based on the structure 

of the World Wide Web. In this case, each web page is 

a node and two pages are connected if there a hypertext 

link from one to the other. These are typically generated by web crawlers and can vary greatly depending 

on where they crawl and what decisions they make 

about what constitutes a link. For example: Do different anchors to the same page count as separate links? 

How are non-static pages handled? What file types are 

processed? How robust it is with broken links and 

sloppy HTML pages? This list can be quite long, and 

each decision path will yield a different graph. Aside 

from this, the Web is, of course, not a static entity. 

Many of the Web graphs cited provide little of these 

details and usually present only the topological structure (i.e., URLs removed) so it is difficult to confirm or 

validate them independently. 

To address these issues, we have developed and 

deployed our own set of web crawlers that report not 

only topological structure, but also record the node 

URLs for verification and validation, as well as generate a unique checksum to ensure that pages which 

contain identical material (e.g., aliased, or mirrored 

pages) are uniquely identified. Furthermore, these web 

crawls can be fined-tuned to record and characterize 

data along a precise decision path, thus allowing us to 

see what differences are generated when we alter the 

definition of what is a link and what is a node. 

 Similar concerns exist for graphs taken from other 

application domains. Often the data is incomplete, 

poorly documented, or have missing metadata crucial 

to its characterization. There are ongoing debates about 

whether particular networks are really scale-free or not, 

or whether the conclusions of some papers are actually 

valid. Usually these concerns lie at the heart of how 

data was collected and analyzed. 

 To help improve this situation we are developing 

standard reference data sets for complex networks. 

This collection of public and well-documented “realworld” graphs will serve as a test bed of reasonably 

validated data which researchers can utilize to 

1) verify and compare the algorithms and analysis of 

various approaches in the literature, 

2) ensure that everyone is using the same reference 

data in their studies, 

3) have one convenient location where one can 

browse and search for network graphs with various 

characteristics, 

4) provide a focal point for the research community 

to contribute and exchange network graphs from 

various application domains, and 

5) provide a testbed for development of software 

analysis tools (e.g., viewers, partitioners, clustering algorithms) that will aid in the further 

development and research of complex networks. 

 In further aiding the understanding of such network 

systems, we have developed a set of software tools for 

working with mathematical graphs (Ngraph++) and 

processing and transforming data collections into various formats (e.g., Matrix Market, GraphViz, Matlab) 

for output and analysis. Ngraph++ computes clustering 

coefficients and degree distributions for directed/undirected graphs. Although a small and simple 

library, it provides intuitive graph operations (subgraphs, unions, differences) and complex-network 

features not commonly found in more complex packages. Together with various Python and scripting tools, 

and a reference dataset of network graphs, this forms a 

kernel package for analyzing complex network data in 

a portable and universal fashion. 
58 NIST Mathematical and Computational Sciences Division 

Community Detection in Network 

Graphs 

Roldan Pozo 

Clement Rey 

Detecting community structure in large graphs is one 

of the fundamental challenges in network science. This 

line of research probes at the underlying structure of 

complex systems, and is necessary in understanding 

how such systems behave over time. A community is 

loosely defined as a subset of nodes which are highly 

connected to each other, but less connected to the remaining graph. Identifying such node groups is 

important in understanding graph evolution, synchronization, and the dynamics of networks. In various 

application contexts, community detection can classify 

customers with similar purchasing habits, determine 

graph layouts for visualization, model the spreading of 

infectious diseases, and help in identifying crime cells. 

 Extracting such structure reliably from large-scale 

networks, however, is a difficult problem. First, there 

must be some mathematical metric for what constitutes 

a graph community; secondly we need an algorithm for 

finding such communities (related solutions are NPcomplete); lastly, we have to have some way of verifying that the computed groups are indeed meaningful. 

Various metrics have been proposed and several optimization algorithms introduced in the literature. 

 The fundamental question remains: How well do 

these approaches work in practice? So far, there has 

been little work in verifying and validating the results 

on a wide array of application domains. To aid this 

study, we are currently investigating eigenvalue-based 

methods and fast greedy algorithm approaches which 

run in O(n log2 n) for sparse graphs, and studying their 

behavior on specific data sets for which we have 

ground truth community information available. Initial 

results suggest that current state-of-the-art algorithms 

provide only mixed results and further research is 

needed to provide truly reliable solutions. 

Visualization of Large Complex 

Networks 

Roldan Pozo 

Clement Rey 

Terence Griffin 

One of the fundamental challenges in complex systems 

is to understand the nature and structure of large-scale 

networks. These can be information networks (such as 

the World Wide Web) biological networks (like protein-protein interactions, food webs) social networks 

(e.g., online social sites, personal interactions, criminal 

rings, terrorist cells) and engineered networks. Such 

networks typically exhibit thousands to millions of 

nodes and about 10 times as many edges. Given the 

rapid growth of electronic data, networks with billions 

of nodes soon will be readily available for analysis, 

providing greater insight into the behavior of largescale systems. 

 To the aid the understanding of such systems, we 

have developed visualization techniques that utilize the 

latest immersive and 3D environments, together with 

state-of-the-art graph layout algorithms, so that researchers can walk around and examine these 

structures interactively. The visualizations can also be 

run on independent desktops. 

 For example, Figure 29 is a visualization of a web 

graph, were nodes are individual web pages of George 

Washington University. In this instance, we have highlighted the community structure by using a fatpoint 

shader (with an interactive size control) to render the 

nodes of the graph. The largest groups (7 in this case) 

are colored individually. The set of small groups (179 

in this case) is colored a transparent grey. An individual node can be selected using a cross-hair glyph and 

its URL displayed in a box tethered to the node with 

lines. In addition, the web page that the node represents 

can be brought up, so researchers can investigate its 

content. 

Figure 29. 3D visualization of community structure of the World 

Wide Web graph of the George Washington University Computer 

Science Department (www.cs.gwu.edu). Each node is a separate 

web page and colors indicate different community groups. The 

transparent grey color indicates the subset of very small groups.
Yearly Report, Fiscal Year 2009 59

Fluid-level Modeling of TCP 

Vladimir Marbukh 

Daniel Genin (NIST ITL) 

While modeling of Transmission Control Protocol 

(TCP) at the packet level yields accurate results, its 

computational complexity makes packet-level modeling infeasible even for moderate size networks. Flowlevel models allow for significant reduction in computational and simulation complexity at the cost of loss in 

accuracy. Conventional flow-level models [1], which 

assume that at equilibrium each flow has fixed rate, do 

not take into account TCP flow burstiness. This simplifying assumption results in underestimation of the 

packet loss probability and thus overestimation of the 

TCP throughput which is inversely related to the 

packet losses. Quantification and mitigation of this 

inaccuracy is an open issue of significant practical and 

theoretical importance due to the central role played by 

TCP in the Internet. 

In [2] we assess that conventional fluid models [1] 

result in unacceptably high underestimation of the 

numbers of flows in progress and file transfer times in 

TCP-controlled network under arriving/departing 

flows. To alleviate this deficiency [2] proposes modeling TCP at flow-level by an on-off fluid process [3]. 

Model [3], often named the “AMS model” after authors 

D. Anick, D. Mitra and M. M. Sondhi, accounts for the 

flow burstiness while retaining the computational tractability. Our extensive ns2 simulations [2] suggest that 

the AMS flow-level model of TCP may have substantially better accuracy than the conventional fluid model 

with fixed rate [1], potentially leading to faster yet accurate Internet simulation. More efforts are needed to 

improve accuracy of the proposed AMS model under 

high load by accounting for the TCP slow start phase. 

[1] R. Srikant, Mathematics of Internet Congestion Control, 

Willey, 2001. 

[2] D. Genin and V. Marbukh, “Bursty Fluid 

Approximation of TCP for Modeling Internet 

Congestion at the Flow Level,” Forty Seventh Annual 

Allerton Conference on Communication, Control, and 

Computing, Monticello, Illinois, USA, 2009. 

[3] M.M. Sondhi D. Anick, D. Mitra. Stochastic theory of a 

data handling system with multiple sources. ICC'80; 

International Conference on Communications, Seattle, 

USA, 1980. 

Analysis of a Distributed Protocol 

for Network Control 

Fern Y. Hunt 

Vladimir Marbukh 

This year we made progress in studying the properties 

of a discrete dynamical system that models TCP congestion protocol with multipath routing. Here, traffic is 

split over several routes according to a probability distribution that is constrained by an entropy parameter. 

This parameter measures the robustness or path diversity of the splitting scheme. F. Kelly and others have 

observed that many protocols currently in use can be 

viewed as distributed algorithms for solving an optimization problem. In two sample network topologies we 

show that the corresponding discrete system does tend 

to a stable equilibrium that is the solution of a network 

utility optimization problem when the entropy of the 

route distribution is greater than a critical value. This 

value is a function of the network link capacities and 

network topology. For our model we show that there is 

no stable equilibrium for entropy values less than the 

critical value. A NIST Internal Report describing this 

research has been released [1]. 

 The existence of an explicit parameter measuring 

path diversity allows us to evaluate the trade-off between path robustness and network utility. We looked 

again at the two sample networks which in both cases 

have a single source and destination. We computed the 

the long term average utility as a function of entropy. 

Starting with the critical value, as the entropy increases, the network utility (which we show depends 

on the minimum average or mean cost) decreases. Thus 

increased robustness leads to decreased utility. This 

work is also described in the internal report [1] and in a 

later paper that focuses on the tradeoff between utility 

and stability [2]. 

[1] F. Hunt and V. Marbukh, “A Model of Joint Congestion 

Control and Routing Through Random Assignment of 

Paths,” NIST Interagency Report 7590. 

[2] F. Hunt and V. Marbukh, “Measuring the Utility/Path 

Diversity Trade-Off in Multipath Protocols,” Proceedings of the Fourth International ICST Conference on 

Performance Evaluation Methodologies and Tools, ValueTools09, Pisa, Italy. 
60 NIST Mathematical and Computational Sciences Division 

Visualization of Possible Route 

Probability Distributions 

Fern Y. Hunt 

Terence Griffin 

In order to control congestion and maximize network 

utility, a variety of routing protocols can be used to 

allocate network traffic. When multiple paths linking a 

source and destination are available, but only a single 

path, e.g., a shortest path or path of minimum cost, is 

selected, we have the Open Shortest Path First (OSPF) 

protocol. On the other hand, traffic could be allocated 

uniformly when the paths have the same cost. A spectrum of path specification strategies lies in between 

these extremes. We visualize these possibilities in the 

case of a single source destination pair with 3 possible 

connecting paths. If β1, β2, and β3 are the fractions 

(probabilities) of traffic assigned to paths 1, 2, and 3 

respectively, then each allocation can be associated 

with a point in or on the equilateral triangle with an 

altitude height of 1. Here, β1, β2, and β3 are the lengths 

of perpendiculars from points on edges of the triangle. 

The entropy associated with an allocation is a function 

of β1, β2, and β3, and provides an important quantitative 

description of its degree of randomness. A point in the 

triangle is colored according to the value of the entropy 

of the allocation represented by that point. The surface 

above the triangle was created by plotting a point 

above the data point in the triangle where the height is 

equal to the Entropy. OSPF is associated with H = 0, 

and traffic is allocated according to one of the vertices 

of the triangle. This is the least robust allocation. Equal 

cost multipath allocation, assigned to the center of the 

triangle has the maximum value of H = log(3) is the 

most robust but it has reduced network utility. (See the 

magenta dot in the center of the triangle). Bands of 

points with the same value of H (up to accuracy .001) 

are given the same color. As H decreases, an oval 

emerges from the center and grows until it touches the 

sides of the triangle at the maximum entropy for a network with 2 paths (depicted by the magenta oval in the 

triangle and band in the surface above it). Theoretical 

studies predict a change in network dynamics for H 

smaller than this critical value. 

Our visualization runs immersively and on the 

desktop. Entropy is displayed as the third dimension 

and uses an RGB color ramp on entropy. The visualization allows selection of entropy values, and groups of 

related points. The structure can also be examined with 

clipping planes of various types, and also by selection 

of two dimensional or three dimensional representations. The data can be displayed in three ways: points 

with fatpoint shaders with adjustable point size; polygonal surface; and wireframe surface. 

Figure 30. This figure shows the surface over the triangle defined by the equation: H(a1

,a2,a3) = −a1*log(a1)-a2log(a2)-

a3

log(a3) , where the height of a point on the surface above the 

triangle is the value of H for the point in the triangle representing 

(a1

,a2,a3). H is the entropy function. The top of the surface is the 

maximum value of the entropy, i.e. log 3 which corresponds to the 

point (1/3,1/3,1/3) whose image in the triangle is the pink dot in 

the center. The surface illustrates the changes in the level curves 

of H (i.e. points in the triangle with a constant entropy value) as 

H decreases from log 3 to 0. Colors indicate different values of 

entropy. On the left is the fatpoint shader representation. On the 

right (top) are two dimensional and polygonal (bottom) representations. 
Yearly Report, Fiscal Year 2009 61

Grid Network Simulation Using 

Piece-Wise Homogeneous Markov 

Chains 

Fern Y. Hunt 

Christopher A. Dabrowski (NIST ITL) 

The growing use of services provided by cloud and 

grid computing networks underlines the need for analytical and computational tools for understanding and 

predicting the behavior of such complex systems. Our 

work centers on the analysis of grid networks. Here the 

challenge is that under decentralized control, the interactions of many service providers and customers could 

result in emergent behaviors resulting in unpredictable 

and undesirable outcomes. A straightforward approach 

is to perform a discrete event simulation to model customers and service providers entering into service level 

agreements for a variety of computer services. Simulations that seek to accurately reproduce system structure 

and component behavior are computationally expensive however. Even simulating an entire day takes 

many hours of CPU time (see [4]). As analysts seek to 

identify system execution paths and scenarios that lead 

to system breakdown and other anomalous behavior, 

there is a critical need for models that can provide 

rapid scalable simulation of large scale grid systems. 

One should also be able to quickly perturb these models in a search for such scenarios. 

 In the previous year, [2,3], we developed a piecewise homogeneous Markov chain model of a grid network system that compared very favorably with 

discrete event simulations (i.e., large-scale simulation) 

including several failure scenarios arising from the 

perturbation of transfer rates in the grid system. Moreover, the Markov chain simulation produced these 

results in hours while the large-scale simulation required days of processing. 

 This year, we took another step in our goal, introducing an analytical approximation of the Markov 

chain that eliminated the need for simulation when 

computing the task completion rate. These computations are much faster than Markov chain simulation. 

Secondly we introduced a method for finding failure 

scenarios in the Markov chain model based on identifying cut point sets in the graph of the Markov chain. 

Identifying the most critical cut sets is an alternative to 

the brute force perturbation method used in [2,3], and 

in a forthcoming paper, we demonstrate that a heuristic 

node contraction algorithm can very quickly identify 

many of these sets in some sample large graphs. The 

results of this algorithm compare quite favorably with 

the algorithm developed by Tsukiyama et al. that performs a complete enumeration of all cutsets. The speed 

of the node contraction algorithm relative to the algorithm of Tsukiyama et al. depends on the connectivity 

of the graph. In the examples tested, node contraction 

is faster on graphs with moderate connectivity but it is 

dramatically slower when graph connectivity is high. 

We are currently developing a report on our work this 

year [1]. 

[1] C. Dabrowski, K. Morrison, and F. Hunt, “A Fast and 

Scalable Representation of Grid Network Failure,” in 

preparation. 

[2] C. Dabrowski and F. Hunt, “Markov Chain Analysis for 

Large-Scale Grid Systems,” NIST Interagency Report 

7566. 

[3] C. Dabrowski and F. Hunt, “Using Markov Chain 

Analysis to Study Dynamic Behavior in Large-Scale 

Grid Systems,” in Proceedings Seventh Australasian 

Symposium on Grid Computing and e-Research (Ausgrid2009) , Wellington, New Zealand, CRPIT,99, P. 

Roe and W. Kelly, Eds. ACS, pp.29-40. 

[4] K. Mills and C. Dabrowski, “Investigating Global Behavior in Computing Grids,” in Self Organizing Systems, 

Lecture Notes in Computer Science, volume 4124,

Springer-Verlag, pp. 120-136. 
62 NIST Mathematical and Computational Sciences Division 

Quantum Information

An emerging new discipline at the intersection of quantum mechanics and computer science, quantum information science (QIS) is likely to revolutionize science and technology in much the same way that lasers, 

transistors, electronics, and computing did in the 20th century. At the very least, QIS has the potential to provide phenomenal increases in information storage and processing speed and communication channels with 

extremely high levels of security guaranteed by the laws of physics. However, the needed manipulation and 

control of quantum states remains a grand challenge. We are developing theoretical underpinning to support 

the multi-laboratory experimental program in quantum information at NIST. This work is aimed at the development of a measurement science for quantum information processing systems. — Ronald Boisvert, Program 

Coordinator

Restrictions on Transversal Encoded 

Quantum Gate Sets 

Bryan Eastin 

Emanuel Knill 

See feature article, page 21. 

Quantum Computing Theory 

Bryan Eastin 

Scott Glancy 

Emanuel Knill 

Adam Meier 

Mike Mullan 

Yanbao Zhang 

R. Somma (Perimeter Institute and LANL) 

S. Boixo (California Institute of Technology) 

G. Ortiz (University of Indiana) 

Quantum computing involves the exploitation of quantum mechanical effects to improve our ability to solve 

computational problems such as factoring and physics 

simulation. Our contributions to quantum computing 

theory include work on better understanding the limitations as well as the power of quantum computers. 

 To improve our ability to exploit the advantages of 

quantum computing, we are working on better understanding and applying adiabatic processes in 

algorithms. An adiabatic process normally requires 

slowly varying the parameters of the Hamiltonian determining the evolution of a quantum system. An initial 

ground state is thus induced to follow the ground states 

of the changing Hamiltonian. In this way one can prepare complex states from simple initial states. The 

conditions for this following are quite subtle, but we 

found that by randomizing the evolution and making 

energy measurements, the desired state preparation can 

be performed more efficiently, with the conditions having a simpler dependence on the parameters. 

 The advantages of quantum computing are easiest to 

see in the context of so-called oracle problems. In these 

problems, we are given an oracle or “black box” that 

we can apply to a quantum register. The task is to infer 

a property of the oracle, such as which computational 

basis state's phase it changed (this is the “search” problem). For the purpose of proving advantages, one must 

determine upper bounds on the probabilities of success 

of classical algorithms and compare them to probabilities of success of quantum algorithms. Perhaps the 

most powerful general technique for this purpose is 

called the “adversary method.” An exact version of 

this method can be based on semidefinite programming, and we have been investigating its application to 

determine exact optimal success probabilities of both 

classical and quantum algorithms for small instances of 

oracle problems. An important generalization that we 

are developing involves oracles that interpolate between quantum and classical by containing varying 

amounts of decoherence. We intend to apply this to 

characterize the power of experimentally available 

quantum registers. 

 One of the often-mentioned features of spatially 

separated quantum systems is that they have states that 

cannot be described with so-called local realistic models in terms of independent states of each system, as 

one might expect to be able to do classically. Violation 

of local realism can be detected from the statistics of 

independent measurements of the systems. By making 

such measurements, one can obtain very strong security guarantees for protocols such as secret key 

exchange. Experiments that demonstrate such violations have been performed, but have one or more 

loopholes, the best known of which are the “detection 

loophole” (introduced by the inefficiencies of typical 

quantum measurements) and the “locality loophole” 

(associated with the possibility of hidden communication between the quantum systems consistent with 

relativistic causality). To help close the detection 

loophole, we previously initiated a study of an experimental configuration involving photon sources to 

determine tradeoffs between detector inefficiency, 

source parameters, and the rate at which confidence in 
Yearly Report, Fiscal Year 2009 63

violation is gained. We have now completed this work. 

A possible loophole not often considered is introduced 

by possible dependencies of each data point on what 

happened before in an experiment. This is called the 

memory loophole. Previously, several researchers 

have pointed out that it can be closed by a more sophisticated statistical analysis of the data based on 

martingale inequalities from probability theory and 

statistics. We found that a better and asymptotically 

optimal approach can be derived from prediction-based 

likelihood ratios. This is insensitive to experimental 

stability, possibly simplifying experimental control. 

 Violations of local realism, in the sense just described, appear not to be necessary for quantum 

speedups over classical computation. We are investigating whether a different kind of non-classicality can 

be shown to be necessary for the efficient solution of 

any classically intractable problem. We have demonstrated that this is indeed the case for quantum 

computations composed of two-qubit gates, though a 

similar result may not hold for more general gate sets. 

[1] S. Boixo, E. Knill, and R. D. Somma, “Quantum State 

Preparation by Phase Randomization,” arXiv: 

0903.1652, Journal of Quantum Information and Computation 9 (2009), pp. 833-855. 

[2] E. Knill, “Quantum Computing Q&A,” in review. 

[3] R. Somma, S. Boixo, and E. Knill, “Measurement-Based 

Quantum Computation,” in preparation. 

[4] Y. Zhang, E. Knill, and S. Glancy, “The Statistical 

Strength of Experiments to Reject Local Realism with 

Photon Pairs and Inefficient Detectors,” in review. 

Fault-tolerant Quantum Computation 

Bryan Eastin 

Emanuel Knill 

Adam Meier

Advances in the theory of fault-tolerant quantum computation are important for the development of quantum 

computers as the hardware employed in the first generation of such devices is likely to be highly prone to 

error. The theory of fault tolerance addresses the problem of implementing encoded gates and error 

correction in a way that precludes a catastrophic failure 

from resulting from the malfunction of a few isolated 

components. This can be achieved by preventing errors from spreading within individual blocks, i.e., 

within independent error-correcting units, of a quantum 

code. 

 Previously, we settled a long-standing question in 

quantum information theory, namely, whether there 

exist nontrivial quantum codes for which all logical 

gates can be implemented transversally. Transversal 

encoded gates never cause subsystems within a code 

block to be coupled and are, consequently, fault tolerant. Their simplicity and robustness makes transversal 

encoded gates valuable tools for quantum computation, 

but, for all known codes, such gates are inadequate for 

achieving full computational universality. We showed 

that no nontrivial transversal encoded gate set could be 

universal. 

 We have also begun investigating new procedures 

for magic-state distillation, a technique for improving 

the quality of prepared magic states. At present, the 

most commonly considered method of achieving universal quantum computation is through (magic) state 

injection and distillation. Magic states are those that, 

given the ability to perform Clifford gates (a class of 

operations that are frequently amenable to being made 

fault tolerant), can be both distilled and used to perform universal quantum computation. Though it is a 

widely known technique, many interesting questions 

regarding magic state distillation remain. We hope to 

begin to address some of these by developing a general 

theory of magic state distillation. In addition, we plan 

to calculate the output error probability for a novel 

distillation procedure that we have already developed 

and believe to have some superior qualities with regard 

to resource usage. 

 Finally, we have made some progress on extending 

the techniques of DiVincenzo and Aliferis for using 

unverified ancillae in quantum error correction and 

measurement. They showed that a fault-tolerant 

scheme based on concatenating the 7-qubit Steane code 

could be implemented, with little additional overhead, 

using unverified ancillae and extremely slow measurements. We have extended their result to faulttolerant schemes based on a variety of other codes and 

hope to further expand the range of applicability of this 

intriguing technique. 

 In addition to continuing our work on state preparation, we plan, in the coming year, to investigate 

techniques for passive quantum error correction. Specifically, we plan to study the efficacy of proposed 

passive quantum error correcting architectures and to 

seek new architectures supporting passive correction. 

The success of any of these research projects could 

contribute to the goal of constructing a quantum computer by reducing the overhead and/or precision 

required. 

[1] B. Eastin and E. Knill, “Restrictions on Transversal 

Encoded Quantum Gate Sets,” arXiv: 0811.4262, Physical Review Letters 102 (2009), 110502. 
64 NIST Mathematical and Computational Sciences Division 

Quantum Computing Benchmarks 

Emanuel Knill 

Bryan Eastin 

Adam Meier 

Michael Mullan 

Dietrich Leibfried (NIST PL) 

David Wineland (NIST PL) 

Quantum computers promise to greatly speed up the 

solutions to some problems in mathematics and physics. However, there are substantial obstacles to our 

building quantum computers. Today, we are able to 

see quantum effects with up to 8 quantum bits (qubits), 

but barely able to properly compute with two. We are 

working on ways to consistently determine the amount 

and quality of experimentally available quantum resources by developing standard, device-independent 

benchmarks and participating in experiments to realize 

them, particularly in ion-trap quantum information 

processing systems. 

 In previous years we developed and implemented 

randomized benchmarking strategies for one qubit, a 

strategy that has now been tried with nuclear magnetic 

resonance and superconducting qubits. We are currently working on extending this strategy to more 

qubits. The basic idea is to perform experiments with 

many random sequences of so-called Clifford gates and 

compare the measurement outcomes with the ones predicted from simulations. The results are processed to 

infer an effective probability of errors for each step of a 

typical computation. Clifford gates allow efficient 

prediction of measurement outcomes while being sufficiently powerful to enable universal quantum 

computing given the additional capability to prepare a 

special but robust state. We have the necessary simulation capabilities and developed principles for 

randomizing multi-qubit gates and inferring not only 

overall error probabilities, but given enough data and a 

simplifying assumption, the error behavior of each type 

of gate. Verification of the principles are in progress 

and possible ion-trap experiments are in the planning 

stages. 

 Part of this project is to develop measurement 

strategies for quantifying the amount of quantum resources contained in a given experimental quantum 

information processing device. One approach is to 

obtain standard criteria for when one can claim to have 

a one, two or more qubit quantum register. One of the 

criteria should involve demonstration of a quantum 

advantage, namely that the device can do better than 

any quantum register of smaller size. Randomized oracle problems, such as the search problem, can be used 

for this purpose, but to do so requires a rigorous way of 

comparing not only classical and quantum algorithms, 

but algorithms that use restricted quantum resources. 

We anticipate using the semidefinite programming 

techniques developed in the course of the Quantum 

Computing Theory project for this purpose. 

Quantum State Estimation Theory 

and Applications 

Scott Glancy 

Emanuel Knill 

Bryce Calkins (NIST EEEL) 

Tracy Clement (NIST EEEL) 

Thomas Gerrits (NIST EEEL) 

Adriana Lita (NIST EEEL) 

Alan Migdall (NIST PL) 

Richard Mirin (NIST EEEL) 

Sae Woo Nam (NIST EEEL) 

Gene Hilton (NIST EEEL) 

Kent Irwin (NIST EEEL) 

Konrad Lehnert (NIST PL) 

Leila Vale (NIST EEEL) 

Many emerging technologies will exploit quantum mechanical effects to enhance metrology, computation, 

and communication. Developing these technologies 

requires improved methods to measure the state of 

quantum systems. Quantum state estimation is a statistical problem of estimating an underlying quantum 

state using a collection of measurements made on independently prepared copies of the state. Accurate 

quantum state estimation allows experimentalists to 

answer the questions “What quantum state does my 

device prepare?” and “How confident am I that it has 

prepared that state?” In the last year we have made 

several improvements to the algorithms used to answer 

these questions. These include a criterion for stopping 

iterations of the estimation algorithm and a new iterative algorithm with faster convergence. We have also 

implemented a parametric bootstrap method to assign 

uncertainties in the quantum state estimate and discovered that conventional algorithms for quantum state 

estimation are biased toward states of lower purity. 

 We have been applying quantum state estimation in 

two experiments. The first of these is an experiment to 

prepare Schrodinger cat states of light, which is part of 

the IMS project “Quantum Optical Metrology with NPhotons.” This experiment uses high efficiency transition edge sensors to subtract three photons from a pulse 

of squeezed vacuum light. The photon subtraction prepares the Schrodinger cat state, which is measured by 

homodyne detection. Our quantum state estimation 

software is used to reconstruct the state from the measurements. Our state estimation algorithm corrects for 

the photon loss and electronic noise of the homodyne 

detector to infer the originally created state. The state 
Yearly Report, Fiscal Year 2009 65

created in the laboratory has a fidelity of 04.0 59.0 14.0

+

− with 

an ideal Schrodinger cat state with a mean of 3.08 photons. 

 The second experiment where we have been applying quantum state estimation is the DARPA funded 

collaboration with JILA titled “Analog Quantum Integrated Circuits.” The goal of this project is to create 

entangled states in superconducting microwave circuits 

and to use that entanglement for quantum teleportation. 

We have been developing methods to calibrate the amplifiers used to measure the quantum states' of 

microwaves. We have also used our quantum state 

estimation software for preliminary analysis of 

squeezed states of microwaves. 

[1] S. Glancy, E. Knill, T. Gerrits, T. Clement, B. Calkins, 

A. Lita, A. Miller, A. Migdall, S. W. Nam, and R. Mirin, 

“Innovations in Maximum Likelihood Quantum State 

Tomography,” Perimeter Institute for Theoretical Physics Quantum Information Seminar, Perimeter Institute, 

Waterloo, Canada, September 9, 2009. 

[2] T. Gerrits, S. Glancy, T. Clement. B. Calkins, A. Lita, 

A. Miller, A. Migdall, S. W. Nam, and R. Mirin, “Generation of Optical Schrodinger Cat States by Numberresolved Photon Subtraction from Squeezed Vacuum,” 

in preparation. 
66 NIST Mathematical and Computational Sciences Division 

Virtual Measurement Systems

A virtual measurement is a quantitative result and its uncertainty, obtained primarily by a nontrivial computer simulation or computer-assisted measurements. Examples of virtual measurements include 

computational models of physical systems and visualizations of the results. The Virtual Measurement Systems 

Program introduces metrology constructs ― standard references, uncertainty characterization and propagation, and traceability ― into scientific computation and computer-assisted measurement technologies. As 

with physical measurement systems, development of a virtual metrology infrastructure will result in predictive computing with quantified reliability. In turn, this will enable improved decision making contingent on 

virtual measurement systems. – Andrew Dienstfrey, Program Manager

High Precision Calculations of 

Fundamental Properties of Fewelectron Systems 

James Sims 

Stanley Hagstrom (Indiana University) 

See feature article, page 25. 

Micromagnetic Modeling 

Michael Donahue 

Donald Porter 

Robert McMichael (NIST CNST) 

June Lau (NIST MSEL) 

http://math.nist.gov/oommf/

Advances in magnetic devices such as recording 

heads, field sensors, magnetic nonvolatile memory 

(MRAM), and magnetic logic devices are dependent 

on an understanding of magnetization processes in 

magnetic materials at the nanometer level. Micromagnetics, a mathematical model used to simulate 

magnetic behavior, is needed to interpret measurements at this scale. MCSD is working with industrial 

and academic partners, as well as with colleagues in 

NIST CNST, MSEL, and EEEL, to improve the stateof-the-art in micromagnetic modeling. 

 Michael Donahue and Donald Porter have developed a widely used public domain computer code for 

doing computational micromagnetics, the ObjectOriented Micromagnetic Modeling Framework 

(OOMMF). OOMMF serves as an open, welldocumented environment in which algorithms can be 

evaluated on benchmark problems. OOMMF has a 

modular structure that allows independent developers 

to contribute extensions that add to the basic functionality of OOMMF. OOMMF also provides a fully 

functional micromagnetic modeling system, handling 

both two and three-dimensional problems, with sophisticated extensible input and output mechanisms. 

In fiscal year 2009 alone, the software was 

downloaded more than 5,200 times, and use of 

OOMMF was acknowledged in more than 120 peerreviewed journal articles. OOMMF has become an 

invaluable tool in the magnetics research community. 

 Key developments in OOMMF last year include: 

1) Addition of multi-threaded computation in the 

modules of OOMMF for use on multiprocessor, multi-core computers, enabling 

faster (observed speed-up of four to seven 

times) and larger simulations, and more efficient utility of common computing platforms 

[2]. 

2) Creation of additional modules tailored for the 

simulation of polycrystalline materials, common in practical device designs [3]. 

 OOMMF is part of a larger activity, the Micromagnetic Modeling Activity Group (μMAG), formed 

to address fundamental issues in micromagnetic 

modeling through two activities: the development of 

public domain reference software, and the definition 

and dissemination of standard problems for testing 

modeling software. MCSD staff members are involved in development of the standard problem suite 

as well. There are currently four standard problems in 

the suite, testing both static and dynamic magnetization properties. A new standard problem proposed 

jointly by several authors [1] aims at validation of 

simulations of spin momentum transfer in the "current parallel to plane" geometry. We are in the 

process of publishing this proposed problem as the 

fifth μMAG standard problem. Work continues to 

develop additional standard problems to address the 

spin momentum transfer phenomenon in other geometries. 

[1] M. Najafi, et al., “Proposal for a Standard Problem for 

Micromagnetic Simulations Including Spin-transfer 

Torque,”' Journal of Applied Physics 105 (2009), pp. 

113914-113921. 

[2] M. J. Donahue, “Parallelizing a Micromagnetic Program for Use on Multi-processor Shared Memory 
Yearly Report, Fiscal Year 2009 67

Computers,” IEEE Transactions on Magnetics 45

(2009), pp. 3923-3925. 

[3] J. W. Lau, R. D. McMichael, and M. J. Donahue, 

“Implementation of Two-dimensional Polycrystalline 

Grains in Object Oriented Micromagnetic Framework,” Journal of Research of the National Institute of 

Standards and Technology 114(Jan 2009), pp. 57-67. 

Uncertainty in Virtual Measurements from Quantum Chemistry 

Models 

Raghu Kacker 

Rüdiger Kessel 

Karl Irikura (NIST CSTL) 

Russell Johnson (NIST CSTL) 

By a virtual measurement we mean a prediction 

along with its associated uncertainty for the value of 

a measurand determined from a computational model 

as an alternative to a physical measurement. An important application is quantum chemistry, where 

carefully determined uncertainties have not been reported. As the technology improves, the need and 

importance of reliable uncertainties in virtual measurements is being recognized. This project is focused 

on developing and applying methods for quantifying 

the uncertainty associated with a virtual measurement 

in quantum chemistry. The benefits accrue to research and development of chemical processes, 

materials development, and drug discovery. 

 Predictions from computational quantum chemistry models seldom agree with the corresponding 

high-quality physical measurements. The differences 

are not random but systematic. Therefore, a common 

practice is to apply an empirical scaling factor to 

computational predictions to bring them closer to the 

true quantity values. The empirical scaling factor 

carries uncertainty. We have developed a methodology to quantify the uncertainty associated with a scaling factor. This approach is based on the Guide to the 

Expression of Uncertainty in Measurement, which is 

an international standard. The uncertainties for scaling factors lead to the corresponding uncertainties for 

virtual predictions. In 2005 we reported uncertainties 

in the scaling factors for ab initio vibrational frequencies from 40 models (combinations of theory 

and basis sets) intended for predicting fundamental 

frequencies from computed harmonic frequencies. 

The uncertainties were larger than generally acknowledged. Subsequently in 2006-2007, we 

determined a new reference data set of 60 diatomic 

and 15 polyatomic experimental vibrational zeropoint energies (ZPEs) that includes estimated uncertainties. In 2008 we refined the methodology for 

calculating uncertainties and reported uncertainties 

for 32 models intended for predicting vibrational 

ZPEs form computed harmonic frequencies. The 

ZPE paper appeared this year in the Journal of 

Chemical Physics [1]. An erratum was published to 

correct an insignificant error in one measurement. 

 In FY 2009 we determined uncertainties for virtual measurements from computed an-harmonic 

vibrational frequencies. A paper on this topic has 

been submitted for publication [2]. Subsequently, 

some small changes were implemented in the latest 

version of a research-software (Gaussian) used in our 

calculations. So we are recalculating to verify our 

results. 

[1] K. Irikura, R. Johnson, R. Kacker, and R. Kessel, 

“Uncertainties in Scaling Factors for ab-initio Vibrational Zero-Point Energies,” Journal of Chemical 

Physics 130 (11)(2009), pp. 114102-114111. 

[2] R. Johnson III, K. Irikura, R. Kacker and R. Kessel, 

“Scaling factors and uncertainties for ab-initio anharmonic vibrational frequencies,” in review. 
68 NIST Mathematical and Computational Sciences Division 

Pervasive Information Technology

Pervasive information technology is the trend towards increasingly ubiquitous connected computing sensors, 

devices, and networks that monitor and respond transparently to human needs. The Pervasive Information 

Technologies Program facilitates the creation of standards for sensor communication, networking interoperability, and sensor information security. The Program enables the use of pervasive information 

technologies to enhance personal and professional productivity and quality of life. – Kamran Sayrafian-Pour, 

Program Manager

Visualization of Radio Frequency 

Propagation for Body Area Networks 

John Hagedorn 

Judith Terrill 

Adele Peskin 

Kamran Sayrafian-Pour (NIST/ITL) 

Wen-Bin Yang (NIST/ITL) 

See feature article, page 35. 

Self-organization of Mobile Sensor 

Networks 

Vladimir Marbukh 

Kamran Sayrafian-Pour (NIST ITL) 

Mobile sensor networks are envisioned for detecting 

and tracking potential targets and events for civilian as 

well as military purposes. Sensor locations in a mobile 

network affect both the sensor’s ability to acquire information on the intended target(s) and event(s) as well 

as the sensor’s ability to communicate this information 

to the intended recipient(s). Information acquisition 

needs, which require proximity to the target(s), often 

compete with communication needs, which require 

proximity to the recipient(s) of the sensor information. 

Communication ability can be improved if sensors are 

capable of forming and optimizing a multi-hop Mobile 

Ad-hoc NETwork (MANET) with some sensors relaying the other sensors’ information possibly in addition 

to transmitting their own information. The Network 

Utility Maximization (NUM) framework for crosslayer network optimization has been successfully applied in various networking domains from wire-line 

network to MANET [1]. NUM optimizes trade-offs 

between demands of different users in a networked 

system by a distributed, price-based algorithm. 

In [2] we propose NUM extension intended to balance competing mobile sensor network information 

acquisition and communication goals. Since achieving 

both these goals depends on the sensor locations relatively to the target(s), to the access point(s), and to 

each other, the system utility to be maximized is a 

complex function of sensor rates and locations. The 

appeal of this unified view of mobile sensor network 

optimization is that while system utility maximization 

over sensor information flow rates in “fast” time scale 

yields the optimal cross-layer design of the sensor 

MANET, system utility maximization over sensor locations in “slow” time scale guides sensor repositioning. 

Realization of this vision, however, requires ability to 

quantify the sensor rate and location dependent system 

utility. Due to possible redundancy of information acquired by different sensors, the aggregate utility of 

sensor information streams may not be the sum of utilities of information streams from different sensors. 

It is natural to assume that the sensor information 

recipient, who fuses all sensor information into a coherent picture, is in the best position to estimate the 

marginal utility of each sensor information stream. On 

the other hand, each mobile sensor has direct knowledge of its remaining battery power supply and 

surrounding terrain which affect sensor ability to communicate and relocate. In [3] we propose to resolve 

this information asymmetry by assuming that mobile 

sensor network self-organization is guided by the sensor information fusion center. The proposed framework 

for self-organization of mobile sensor networks includes cross-layer network optimization as well as 

controlled sensor mobility. Given sensor locations, 

cross-layer network optimization allocates resources 

and configures protocols to ensure the delivery of the 

highest utility of the sensor information to the access 

point taking into account battery energy conservation 

needs. This framework overcomes deficiencies of numerous existing phenomenological approaches to 

sensor mobility control based on phenomenologically 

defined “potential” fields and corresponding virtual 

forces, which do not account for the specific sensor 

information acquisition and communication needs, e.g., 

different “power costs” of mobile sensors with different remaining battery energy levels. 

[1] X. Lin, N. B. Shroff, R. Srikant, “A tutorial on crosslayer optimization in wireless networks,” IEEE Journal 

of Selected Areas of Communications, Vol. 24, Issue 8, 

pp. 1452-1463, 2006. 
Yearly Report, Fiscal Year 2009 69

[2] V. Marbukh and K. Sayrafian-Pour, “A Framework for 

Joint Cross-Layer and Node Location Optimization in 

Mobile Sensor Networks,” in Proceedings of AdhocNow’09, Nice, France, 2008. 

[3] V. Marbukh and K. Sayrafian-Pour, “Mobile Sensor 

Network Self-organization for System Utility 

Maximization: Work in Progress,” in Proceedings of the 

International Conference on Wireless and Mobile 

Communications 2009 (ICWMC’09), Cannes, France, 

2009. 
70 NIST Mathematical and Computational Sciences Division 

Trustworthy Information Systems

A trustworthy system is one that performs as intended for a specific purpose, when needed, with operational 

resiliency and without unwanted side effects, behaviors, or exploitable vulnerabilities. The Trustworthy Information Systems (TIS) Program conducts research, development, and testing to improve the ability to 

model, build, test, measure, and assess information system trustworthiness through the development and application of new and innovative technologies, models, measurement methods, and tools. The aim of the TIS 

program is to reduce the risk and uncertainty associated with information systems that must be trusted, and 

to improve the ability to build in and evaluate trustworthiness in information applications and systems. – 

Tom Rhodes, Program Manager

Automated Combinatorial Testing 

for Software-systems (ACTS) 

Raghu Kacker 

Rüdiger Kessel 

James Lawrence 

Stephen Langer 

Andrew Reid 

Jeff Yu Lei (U. Texas, Arlington) 

Rick Kuhn (NIST ITL) 

Vincent Hu (NIST ITL) 

Thomas Rhodes (NIST ITL) 

Renee Bryce (Utah State University) 

Sreedevi Sampath (U. Maryland Baltimore County) 

Tao Xie (North Carolina State University) 

A combinatorial interaction fault is a defect that may 

cause a complex software system to fail when particular values of certain parameters combine with 

particular values of other parameters. Combinatorial 

interaction faults often remain dormant until the unfortunate combination of values is encountered in 

practice. Many critical system failures in recent times 

have been traced to interaction faults in the underlying 

software. Combinatorial testing is an approach to detect such interaction faults. It is a critical technology 

for ensuring trustworthy software in business, industrial, medical, scientific, and transport systems. 

 It is often difficult or impossible to characterize all 

possible combinatorial interactions among software 

input parameters, let alone test them. Combinatorial 

testing is based on the observation that most interaction 

faults involve only a few parameters. Several studies of 

actual faults show that most involved a single parameter, a smaller proportion resulted from interactions 

between values of two parameters, progressively fewer 

interaction faults involved 3, 4, 5, and 6 parameters. So 

far a fault involving more than six parameters has not 

been seen. Most combinatorial testing has so far been 

limited to pair-wise (2-way) testing in which all interactions involving two parameters are detected because 

development of efficient test methods to detect higher 

order interaction faults have not been available. The 

goal of this project is to advance the technology from 

pair-wise to higher order (strength) testing and to demonstrate successful applications. 

 In 2008 we issued version 1 of a publicly available 

research tool, FireEye, to generate test suites for high 

strength testing. The name subsequently was changed 

to ACTS (Advanced Combinatorial Test Suites) because FireEye was trademarked. In September 2009 we 

released version 9. Major improvements include the 

following. 

1) Base-Choice Testing allows the user to designate 

a base value for each parameter (which is typically the default or most frequently used value), 

and then it combines other values with these base 

values in a well-defined manner. 

2) Test Output Specification allows the user to 

specify the output parameters (in addition to input parameters) and their expected values for 

each test case. Expected outcome is needed for 

test evaluation. 

3) New Interface for Constraint Support allows the 

user to free-type constraints, which is much more 

user-friendly than button-clicks. The new release 

also adds constraint support to the command line 

interface. Significant validation support has also 

been added in the backend to ensure the legality 

of constraints. 

4) Miscellaneous Improvements have been continually made in response to user feedback. 

Examples include better error handling, additional export formats, better mixed-strength 

support, and various GUI enhancements such as 

in-/out-synchronization (between system configurations and test sets) indicators, new 

relation/constraint display format in the tree 

view, and buttons that are made more contextsensitive, and so on. 

We have distributed over 265 copies of ACTS tool 

to US government agencies and global corporations 

including Microsoft, Lockheed Martin, Cisco, SAP, 

AT&T, Intel, and SAIC. The Mid-Atlantic section of 

the Federal Laboratory Consortium recognized ACTS 
Yearly Report, Fiscal Year 2009 71

tool with its Regional Excellence in Technology Transfer Award on September 17, 2009. 

 We have published eight papers in refereed publications and three magazine articles. We have ongoing 

collaborations with the Food and Drug Administration 

and the Applied Physics Lab of Johns Hopkins University. We are investigating applications in computer 

security, health care, GUI based applications, webapplications, and for optimization and testing of simulations of complex systems. 

[1] R. Kuhn, R. Kacker, Y. Lei, and J. Hunter, “Combinatorial Software Testing,” Computer, August 2009, pp. 94-

96. 

[2] W. Wang, S. Sampath, Y. Lei, and R. Kacker, “An Interaction-based Test Sequence Generation Approach for 

Testing Web-applications,” in Proceedings of 11-th 

IEEE High Assurance Systems Engineering Symposium, 

Nanjing, China, December 3-5, 2008. 

[3] B. Johnasson, R. Kacker, R. Kessel, C. McLean, and R. 

Sriram, “Utilizing Combinatorial Testing to Detect Interactions and Optimize a Discrete Event Simulation 

Model for Sustainable Manufacturing,” in Proceedings 

of the 14th Design for Manufacturing and the Life Cycle 

Conference, San Diego, CA, August 30-September 2, 

2009. 

[4] R. Bryce, Y. Lei, D. R. Kuhn, and R. Kacker, “Combinatorial Testing” in Handbook of Research on Software 

Engineering and Productivity Technologies: Implications of Globalization, to appear. 

[5] W. Wang, Y. Lei, S. Sampath, R. Kacker, R. Kuhn and 

J. Lawrence, “A Combinatorial Approach to Building 

Navigation Graphs for Dynamic Web Applications,” 

Proceedings of 25th IEEE International Conference on 

Software Maintenance, to appear. 

[6] J. Lawrence, R. Kacker, Y. Lei, D. R. Kuhn, and M. 

Forbes. “Binary Covering Arrays,” in review. 
72 NIST Mathematical and Computational Sciences Division 
Yearly Report, Fiscal Year 2009 73

 

Division Programs

 

 

 

Division Programs

 

 

 

 

Digital 

Library of 

Mathematical 

Functions 

Mathematics of 

Metrology 

High 

Performance 

Computing and 

Visualization 

Fundamental 

Mathematical 

Software 

Development 
74 NIST Mathematical and Computational Sciences Division 

Digital Library of Mathematical Functions

The special functions of applied mathematics are fundamental tools enabling modeling and analysis in all 

areas of science and engineering. To make effective use of such functions, practitioners must have ready access to a reliable source of information on their properties. The goal of this work is the development and 

dissemination of definitive reference data on the special functions of applied mathematics. The DLMF will be 

a freely available interactive and richly linked online resource. – Daniel Lozier, Program Coordinator

Digital Library of Mathematical 

Functions 

Daniel Lozier 

Ronald Boisvert 

Peter Ketcham 

Marjorie McClain 

Bruce Miller 

Frank Olver 

Bonita Saunders 

Abdou Youssef 

Brian Antonishek (NIST ITL) 

Charles Clark (NIST PL) 

Eight associate editors, 28 authors, and 23 validators 

at universities and research institutions 

http://dlmf.nist.gov/ 

Mathematics and science go hand in hand. Purely 

mathematical developments find application in the description of scientific phenomena, seemingly by 

serendipity; in the other direction, cutting-edge science 

is a major motivator and driver of mathematical research. Often the mathematical objects at the 

intersection of math and science are mathematical 

functions. Thus there has always been an interest in, 

and a practical need for, accurate and comprehensive 

compilations of the many intricate properties of 

mathematical functions. 

Traditionally, massive published handbooks have 

been the accepted method of supplying compilations of 

scientific reference data. An outstanding example for 

mathematical reference data is the 1,000-page NBS 

Handbook of Mathematical Functions, created in 1964 

under the editorship of Milton Abramowitz and Irene 

Stegun. The impact of this handbook was tremendous 

and attracted immediate attention, especially among 

physical scientists and engineers. It remains a standard 

reference today, even though it was never updated. 

The first goal of the current project is to produce a 

new handbook that serves the same readership and 

takes into account the many advances of the last halfcentury in the theory and computation of mathematical 

functions. The second goal is to develop a novel Web 

service that takes advantage of, and builds upon, the 

powerful capabilities of the Web to disseminate essentially the same information as is contained in the new 

handbook. The new handbook is called the NIST 

Handbook of Mathematical Functions (HMF), and the 

Web service is called the NIST Digital Library of 

Mathematical Functions (DLMF). 

After a decade of effort by the project leaders, assisted by NIST personnel and many external experts, 

the full HMF and DLMF will be publicly available 

early in 2010. External authors under contract to NIST 

produced the 36 chapters of the HMF and transferred 

their copyrights for retention by NIST. These chapters 

were independently checked by validators, and they 

were also checked and edited for consistent style and 

format by the NIST editors (Olver, Lozier, Boisvert, 

and Clark). Highly qualified NIST experts (Miller, 

Saunders, and Youssef) led the development of the 

DLMF, taking advantage of the latest Web technology 

and making significant research contributions to 

mathematical representation and semantics, mathematical graphics, and mathematical search. 

The HMF is licensed for production, marketing, 

and distribution by Cambridge University Press, which 

has already begun its sales campaign; see Figure 31. 

The DLMF will be a Web service5

 operated by NIST. 

The HMF contains 

o Three methodology chapters covering parts of algebra, real and complex analysis, asymptotics, and 

numerical analysis which are important for work 

with mathematical functions. 

o Thirty-three chapters covering the elementary 

functions and selected higher, or special, functions. 

o Definitions, notations, and properties useful in current or potential applications of mathematical 

functions in scientific and related fields. 

o Graphics, including many line graphs and surface 

plots for real and complex arguments and parameters. 

o Methods of computation and approximation. 

o Examples of special functions applications in 

physics, other scientific and engineering fields, 

and applied mathematics. 

 

5 http://dlmf.nist.gov/
Yearly Report, Fiscal Year 2009 75

o Large bibliography and index. 

The DLMF provides for 

o Delivery of the full content of the HMF in a carefully designed Web-oriented format, that contains 

links to external reviews and full texts of bibliographic items, links to external sources of 

mathematical software, and links to the HMF definitions of symbols and HMF cross-references. 

o Delivery of text and math in XML and MathML to 

Web browsers that support these technologies, or 

of text in HTML and math in GIF images to less 

capable Web browsers, and detection of browser 

capabilities automatically at the Web server. 

o Zoomable and rotatable color visualizations of surfaces of real or complex functions of two real 

variables in browsers that support VRML or X3D. 

o Search based on mathematical queries as well as 

text queries (for example, the query “sin^2 + 

cos^2” locates many identities that contain the 

sum of squares of sines and cosines). 

Several future developments related to the DLMF 

project are planned or under consideration. One is to 

further the semantic capabilities of the XML tools and 

techniques developed for the DLMF. Another is to 

develop tools and techniques for generating validated 

numerical values and graphs on demand. A third is to 

extend the power and applicability of the DLMF search 

tools and techniques. Also, resumed development of 

chapters on Probability Functions and Computer Algebra is planned. 

Visualization of Complex Function 

Data 

Bonita Saunders 

Qiming Wang 

Brian Antonishek (NIST ITL) 

Sandy Ressler (NIST ITL) 

Daniel Lozier 

Frank Olver 

http://dlmf.nist.gov

Understanding of high level mathematical functions 

arising from the solution of complex mathematical or 

physical problems can be greatly enhanced by manipulating accurate visualizations of those functions. 

Consequently, we have put a substantial effort into 

finding ways to create precise and informative visualizations of complex function surfaces that can be 

manipulated inside a web environment. While we 

originally thought this to be a somewhat minor task, we 

quickly discovered that our options for creating such 

visualizations were limited both in number and quality. 

Using techniques from numerical grid generation and 

web-based visualization technologies such as the Virtual Reality Modeling Language (VRML) and X3D, 

we have successfully created accurate visualizations 

that capture key function features such as zeros, poles, 

and branch cuts. We have produced over 200 interactive 3D visualizations for the NIST Digital Library of 

Mathematical Functions (DLMF). 

Figure 31. Announcement of the NIST Handbook in the Cambridge University Press online catalog. 
76 NIST Mathematical and Computational Sciences Division 

Having completed all the visualizations for the 

first release of the DLMF, our primary focus this year 

was on improving the presentation quality of the visualizations for the website, validating plot data, and 

making the visualizations accessible on multiple platforms. The validation of the data has been completed 

and substantial progress has been made in refining the 

underlying computational grids in order to improve 

surface color maps. 

Ensuring that our visualizations are accessible on 

multiple platforms has been one of the more problematic issues. While the DLMF VRML files work well on 

Windows systems, where most VRML browsers are 

provided as browser plugins, problems appear in Linux 

and MacOS where VRML browsers are usually standalone applications. Such applications can be made to 

operate like plugins by connecting them to the browser 

as helper applications, but they may not be able to detect the location of support routines located outside the 

main VRML file. Fortunately, we found a way to 

imbed all the routines and images needed by a particular figure into its main VRML file. This key 

breakthrough made it possible to view the DLMF 

VRML files on Mac systems and helped us make substantial progress on Linux systems. We have been able 

to view a stripped down VRML file, containing only 

surface geometry information without our customized 

cutting plane, color map, and scaling features, on a 

Linux system. We plan to provide this type of file for 

users whose systems cannot handle the full featured 

DLMF VRML files. 

 Also, in record time we successfully converted the 

DLMF VRML files into X3D (Extensible 3D) using a 

VRML to X3D translator, with minor post-processing, 

designed by Qiming Wang several years ago. X3D, 

which uses an XML type format, is becoming more 

prevalent for the manipulation and viewing of 3D 

graphics on the web. Access to these files will also be 

put on the website shortly. 

 As our visualization work for the current phase of 

the NIST DLMF winds down, we will concentrate 

more on the on-demand generation of mathematical 

graphs and tables. This may be a stand-alone capability 

or ultimately be incorporated as part of the DLMF. We 

are researching websites that offer some of these capabilities and looking at the possibility of extending work 

on the live generation of tables started by Marjorie 

McClain, Joyce Conlon and former NIST postdoc and 

guest worker Bruce Fabijonas. 

[1] B. Saunders and Q. Wang, “Tensor Product Grid Generation for Complex Surface Visualizations in a Digital 

Library of Mathematical Functions,” in Proceedings of 

the 11th International Society for Grid Generation 

(ISGG) Conference, Montreal, Canada, May 25-28, 

2009, p. SV.3. 

[2] B. Saunders and Q. Wang, “Tensor Product B-Spline 

Mesh Generation for Accurate Surface Visualizations in 

the NIST Digital Library of Mathematical Functions,” in 

Proceedings of the Seventh International Conference on 

Mathematical Methods for Curves and Surfaces, Tonsberg, Norway, June 26-July 1, 2008, to appear. 

[3] Q. Wang, B. Saunders, and S. Ressler, “Dissemination 

of 3D Visualizations of Complex Function Data for the 

NIST Digital Library of Mathematical Functions,” CODATA Data Science Journal 6 (2007) (Supplement on 

the Proceedings of the 20th International CODATA Conference, Beijing, China, October 2006), pp. S146- S154. 

[4] B. Saunders and Q. Wang, “From B-Spline Mesh Generation to Effective Visualizations for the NIST Digital 

Library of Mathematical Functions,” in Curve and Surface Design: Avignon 2006, Sixth International 

Conference on Curves and Surfaces, Avignon, France, 

Neshkoro Press, 2007, pp. 235-243. 

Figure 32. Plot of Jacobian elliptic function sn(x,k) where k = 

1-exp(-n). 

Figure 33. Phase density plot of Dedekind’s eta function η(x+iy) . 
Yearly Report, Fiscal Year 2009 77

Mathematics of Metrology

Mathematics plays an important role in the science of metrology. Mathematical models are needed to understand how to design effective measurement systems, and to analyze the results they produce. Mathematical 

techniques are used to develop and analyze idealized models of physical phenomena to be measured, and 

mathematical algorithms are necessary to find optimal system parameters. Finally, mathematical and statistical techniques are needed to transform the resulting data into useful information. Our goal is to develop 

fundamental mathematical methods and analytical tools necessary for NIST to continue as a world-class metrology institute, and to apply them critical measurement science applications. 

Multihazard Structural Optimization 

Florian Potra 

Emil Simiu (NIST BFRL) 

The objective of multihazard structural engineering is 

to develop methodologies for achieving designs that 

are safe and cost-effective under multiple hazards. Optimization is a natural tool for achieving such designs. 

In general, its aim is to determine a vector of design 

variables subjected to a given set of constraints, such 

that an objective function of those variables is minimized. In the particular case of structural design, the 

design variables may be member sizes; the constraints 

pertain to structural strength and serviceability (e.g., 

keeping the load-induced stresses and deflections below specified thresholds); and the objective function is 

the structure’s cost or weight. In a multihazard context, the design variables are subjected to the 

constraints imposed by all the hazards to which the 

structure is exposed. In particular, the design variables 

should be determined in such a way that the loadinduced stresses and deflections are kept below specified thresholds at all points of the resulting structure. 

Since there are infinitely many such points, the optimization problem becomes a semi-infinite programming 

problem. 

In this project, we formulated the multihazard 

structural design problem in nonlinear programming 

terms, and presented some simple illustrative examples 

involving two hazards: earthquake and strong winds. 

We have also developed efficient bootstrap strategies 

for solving the semi-infinite programming problem. 

Our numerical experiments show that interior-point 

methods are significantly more efficient than classical 

optimization methods in solving the nonlinear programming problem associated with our illustrative 

example. 

 The results concerning the semidefinite programming aspect were presented at the Eleventh 

International Conference of Optimum Design of Structures and Materials in Engineering (OPTI 2009), 

Algarwe, Portugal, June 8-10, 2009, and published in 

the proceedings in that conference [1]. Other results 

were published in [2] and [3]. 

[1] F. A. Potra and E. Simiu, “Semi-infinite Programming 

in Multi-hazard Structural Design,” Computer Aided 

Optimum Design in Engineering XI, in Proceedings of 

the Eleventh International Conference of Optimum Design of Structures and Materials in Engineering (OPTI 

2009), Editors: S. Hernandez and C.A. Brebbia, WIT 

Press, 2009, pp. 13-19. 

[2] F. A. Potra and E. Simiu, “Optimization and Multihazard Structural Design,” Journal of Engineering Mechanics 135 (1)2 (2009), pp. 1472-1475. 

[3] F. A. Potra and E. Simiu, “Multihazard Design: Structural Optimization Approach,” Journal of Optimization 

Theory and Applications, 2009, Online First, DOI 

10.1007/s10957-009-9586-4. 

A Mathematical Model of Atmospheric Retention of Man-made CO2

Emissions 

Bert Rust 

Climate change skeptics have now mostly conceded 

that the Earth is indeed warming, but they still strongly 

deny that the warming is caused by man's activities. In 

2008, Rust and Thijsse showed that from 1856 to the 

present, global annual average temperatures have varied linearly with atmospheric CO2 concentrations. That 

relationship is somewhat obscured by a 71 year cycle 

in the temperature data which, though quite visible (see 

the lower plot in Figure 36), accounts for only about 

8% of the total variance in the measured record. By 

contrast, the accelerating baseline obtained from the 

CO2 relationship accounts for about 77% of that variance. 

The present study shows how the increases in atmospheric CO2 concentration can be derived from 

increases in man-made CO2 emissions, thus completing 

the connection between man's activities and global 

warming. This relationship is established without using a complicated global climate model, but rather by 

using simple mathematical modeling techniques to 
78 NIST Mathematical and Computational Sciences Division 

connect and synthesize six different measured data sets 

containing 5 different kinds of measured quantities. 

Figure 34 shows that man-made CO2 emissions to 

the atmosphere have increased quasi-exponentially 

since 1850. The estimated fraction remaining in the 

atmosphere was calculated by solving a linear inverse 

problem. Estimates obtained by two different approaches are shown in Figure 35. The forward problem 

assumed a harmonic expansion for the unknown function, using a basic period of 150 years. The number of 

harmonics was chosen by optimizing a fit of the emissions remaining in the atmosphere to the record of 

atmospheric CO2 concentrations shown in the top plot 

in Figure 36. The inverse problem approach was more 

difficult because the regression model was highly underdetermined, with an 83ä153 matrix. Estimation was 

possible only because the solution is constrained to be 

between 0 and 1. There is good agreement between the 

two estimates in spite of the high uncertainties in the 

first half of the record. The similarity of the two 

dashed curves in Figure 36 strongly suggests that manmade CO2 emissions are driving the observed global 

warming. 

[1] B. W. Rust, “Atmospheric Retention of Man-made CO2

Emissions,” in MAMERN'09: Proceedings of the 3rd International Conference on Approximation Methods and 

Numerical Modeling in Environment and Natural Resources, B. Amaziane, D. Barrera, M.A. Fortes, et. al., 

(eds.), Universidad de Granada (June 2009), pp. 813-

818. 

Figure 34. The uppermost plot gives yearly total CO2 emissions to 

the atmosphere. The smooth curve, scaled on the right, is an estimate of the fraction of those emissions that remain in the 

atmosphere. The lower plot, the product of the other two, estimates 

the emissions remaining in the atmosphere. 

Figure 35. Two estimates of the fraction of CO2 emissions remaining in the atmosphere. The smooth curve is the forward problem 

estimate, and the dashed piecewise constant curve is an inverse 

problem estimate for 6 nonoverlapping 25-year averages. The solid 

piecewise constant curves are 95% confidence bounds for those 

averages. 

Figure 36. The relationship between atmospheric carbon dioxide 

and global temperatures. The dashed curve in the upper plot is the 

estimate of atmospheric concentrations of CO2 given by the retention 

fraction plotted in Figure 34. The dashed curve in the lower plot is 

the temperature response to those concentrations. 
Yearly Report, Fiscal Year 2009 79

Simulation of Bioregulatory Networks Involved in Cell Cycle Control 

Geoffrey McFadden 

S. Kim (National Institutes of Health) 

Augustin Luna (NIH and Boston University) 

Michael Sharpnack (New York University) 

Mirit Aladjem (National Institutes of Health) 

Kurt Kohn (National Institutes of Health) 

G. McFadden of MCSD is a co-advisor to S. Kim, a 

postdoctoral fellow, and A. Luna, a graduate student 

from Boston University, who are both doing research at 

the Laboratory of Molecular Pharmacology in the National Cancer Institute at NIH; their co-advisors are 

NIH researchers M. Aladjem and K. Kohn. The team is 

developing models of bioregulatory networks that are 

involved in cell cycle control. The models consist of 

systems of nonlinear ordinary differential equations or 

delay differential equations that typically exhibit 

switching behavior, limit cycles, and other types of 

bifurcations. 

Proper cell growth depends on a network of interacting molecules that monitors cellular metabolism and 

environmental signals. This network ensures that cells 

halt their growth in response to unfavorable conditions 

such as the absence of sufficient nutrients or the presence of potentially damaging agents. When cells 

escape these controls, the results are developmental 

abnormalities, genomic instability, and cancer. Much 

current work focuses on a protein known as p53, a tumor suppressor that causes cell cycle arrest or 

programmed cell death in response to stress signals 

such as DNA damage. Regulating the appropriate levels of p53 is essential for cell survival. Two associated 

proteins, Mdm2 and Mdmx, are known regulators of 

p53: Mdm2 can facilitate degradation of p53 whereas 

the mechanism of the regulatory interaction of Mdmx 

with p53 is not clear. It is also not obvious how those 

three molecules will operate together under various 

conditions. 

To address those questions a mathematical model 

has been developed [1] to investigate the interactions 

of these three partner molecules by numerical simulations. An interesting feature of this system is the 

experimental observation of time-periodic behavior of 

the measured amounts of p53 and Mdm2 in the system 

under certain conditions. These results show the stabilizing effect that Mdmx has on the system at long 

times: for large enough amounts of Mdmx in the system the entire branch of steady state solutions is found 

to be linearly stable. This modeling work is intended to 

guide experimental investigations of the role of Mdmx 

in the cell cycle that are being performed by Dr. Kim at 

NIH. 

A. Luna has recently joined the team to conduct 

his Ph.D. research though a combination of stochastic 

and deterministic modeling of cell cycle regulation. M. 

Sharpnack, an undergraduate at New York University, 

was a SURF student in MCSD during the summer of 

2009, and worked on a project entitled “A mathematical analysis of the switch-like and oscillatory behaviors 

of the tumor suppressor p53 feedback loop.” An example of one of the simplified networks considered by 

Mr. Sharpnack in shown in Figure 37. 

[1] S. Kim, M. I. Aladjem, G. B. McFadden, and K. W. 

Kohn, “Predicted Functions of MdmX in Fine Tuning the 

Response of p53 to DNA Damage,” PLoS Computational Biology, in press.

Modeling the Fluid Dynamics of 

Tear Films on the Cornea 

Geoffrey B. McFadden 

Richard J. Braun (University of Delaware) 

R. Usha (Indian Institute of Technology) 

P. E. King-Smith (Ohio State University) 

Dry eye disorders refer to eye problems associated with 

insufficient amounts of fluids coating the surface of the 

cornea. They can be caused by reduced production of 

aqueous tears, and/or increased evaporation due to a 

deficient lipid layer in the fluid. In both of these conditions, tear film breakup is more rapid than normal. 

Breakup can be observed either as dark spots in fluorescein-stained tear films, or, noninvasively, by 

Figure 37. A molecular interaction map showing the reactions 

between p53 and mdm2 in a simple model. The red dots indicate 

formation of a complex p53:mdm2, which can degrade to form 

products indicated by the red null set symbols; the reactants ('in') 

forming p53 and mdm2 are also indicated. Each arrow denotes a 

reaction taking place at a specified rate. 
80 NIST Mathematical and Computational Sciences Division 

distortions in the reflected image of a grid. In a recent 

study [1] we evaluated the contributions of evaporation 

and tangential flow to thinning of the pre-corneal tear 

film between blinks. It was argued that tangential flow 

of the aqueous tears is largely determined by the 

movement of the lipid layer. The results indicated that 

the observed movement of the lipid layer was too slow 

to explain the observed thinning rate of the tear film, 

suggesting that most of the observed tear thinning between blinks is due to evaporation. 

 In order to better understand the dynamics of the 

spreading of tear films on the human eye we have gone 

on to develop in more detail two hydrodynamic models 

based on lubrication theory [2]. The eye is modeled by 

an ellipsoid of revolution, using a prolate spheroidal 

coordinate system. One model is a self-consistent leading-order hyperbolic PDE for relatively large substrate 

curvatures; the other retains another order resulting in a 

fourth order parabolic PDE for the film dynamics. We 

consider both Newtonian viscous flow models and nonNewtonian models that take into account shearthinning of the tear fluid. We explore a wide parameter 

range of shear thinning and find a significant effect on 

finite-time singularities present in the model, representing discontinuous solutions allowed by the underlying 

hyperbolic PDE. Approximate results illustrating the 

nature of the flows are obtained both numerically and 

analytically. We are able to use this approach to evaluate the influence of the cornea's shape in providing a 

driving force for fluid motion through the LaplaceYoung boundary condition that relates pressure and 

interface curvature, an effect that has been neglected in 

previous studies. 

[1] P. E. King-Smith, B. A. Fink, J. J. Nichols, K. K. Nichols, R. J. Braun, and G. B. McFadden, “The 

Contribution of Lipid Layer Movement to Tear Film 

Thinning and Breakup,” Investigative Ophthalmology & 

Visual Science 50 (2009), pp. 2747–2756. 

[2] R. J. Braun, R. Usha, G. B. McFadden, T. A. Driscoll, 

L. P. Cook, and P. E. King-Smith, “Thin Film Dynamics 

on a Prolate Ellipsoid with Application to the Cornea,” 

submitted to AIMS Journal on Discrete and Continuous 

Dynamical Systems A (2009).

Modeling Equilibrium, Stability and 

Transport during Magnetic Fusion 

Geoffrey B. McFadden 

Paul R. Garabedian (Courant Institute, NYU) 

A hot plasma of hydrogen isotopes can be confined in a 

strong magnetic field with toroidal geometry so the 

ions fuse to form helium and release energetic neutrons. Models of such magnetic fusion devices 

developed in computational science have been implemented previously by P. Garabedian and his students 

and colleagues at New York University in codes that 

comprise an effective numerical simulation of the most 

essential features of modern tokamak and stellarator 

experiments. This has led to the discovery of advanced 

concepts that make fusion reactors a future prospect for 

a commercial source of energy. In recent work [1] we 

have applied computer codes [2,3] to issues of equilibrium, stability and transport that provide an example of 

the theory. 

The equilibrium code employs a divergence-free 

representation of the magnetic field that minimizes a 

discretization of the energy using a combination of 

spectral techniques in the toroidal and poloidal directions, and finite differences in the radial direction 

(outward from the plasma's magnetic axis). The stability of the nonlinear equilibrium states is assessed by 

applying appropriate perturbations to the system and 

determining whether they lead to multiple steady states 

that indicate the occurrence of bifurcated solutions. 

Transport properties of the charged species are estimated by following the orbits of test particles that track 

the magnetic field lines between particle collisions. 

The codes have been used to evaluate possible stellarator designs as candidates for follow-up (“DEMO”) 

experiments that may eventually be performed after the 

completion of planned international experiments 

(“ITER”) that are based on tokamak configurations. 

 Tokamaks have axisymmetric symmetry, which is 

an advantage for purposes of particle confinement, but 

they rely on strong toroidal currents in the plasma to 

generate the magnetic fields that can lead to stability 

issues. The magnetic fields for stellarators are generated by external coils which can alleviate many of the 

Figure 38. Magnetic fusion of hydrogen forms helium and emits 

energetic neutrons from a torus of plasma shaped to optimize confinement. Modular coils generate a field keeping the ions from 

hitting walls. 
Yearly Report, Fiscal Year 2009 81

stability problems of tokamaks, but generally lack the 

high symmetry of tokamaks associated with enhanced 

confinement. The stellarator designs reported recently 

[1,4] offer the advantages of a type of quasiaxisymmetry (see Figure 38) that combines the best 

features of each, with the possibility of both good stability and transport properties. 

[1] P. R. Garabedian and G. B. McFadden, “Design of the 

DEMO Fusion Reactor following ITER,” Journal of Research of the NIST 114 (2009), pp. 229-236. 

[2] M. Taylor, “A High Performance Spectral Code for 

Nonlinear MHD Stability,” Journal of Computational 

Physics 110 (1994), pp. 407-418. 

[3] N. Kuhl, “Monte Carlo Simulation of Transport,” Journal of Computational Physics 129 (1996), pp. 170-180. 

[4] G. B. McFadden, “QAS Design of the DEMO Reactor,” 

17th International Stellarator/Heliotron Workshop, 

Princeton Plasma Physics Laboratory, 12-16 October 

2009; Communications in Plasma Physics, submitted. 

Modeling Type Ia Supernova Light 

Curves 

Bert Rust 

Dianne O’Leary 

Katharine Mullen (NIST MSEL)

The regularity and uniformity of their light curves 

make Type Ia supernovae the best available candidates 

for the role of “standard candle” in estimating extragalactic distances. An empirically established correlation 

between their peak absolute magnitudes and the rates 

of decline of their light curves provides a means of 

estimating the peak absolute magnitude of a distant 

supernova from measurements of its rate of decline in 

brightness. The estimate thus obtained can then be 

combined with a measurement of its peak apparent 

magnitude to give an estimate of its distance. We have 

developed a mathematical model for the luminosity 

which, using only six adjustable parameters, gives superior fits to measured light curves and should prove 

useful in calibrating the correlation described above. 

An example is given in Figure 39 where the model 

explains 99.79 % of the variance in the record. 

 Type Ia supernova luminosity is thought to be 

powered by the radioactive decay chain Ni → Co →

Fe, even though observed luminosity decline rates do 

not match the terrestrial decay rates of Ni and Co. In 

1976 Rust, Leventhal and McCall [1] presented evidence supporting the hypothesis that those decay rates 

are accelerated in a supernova because the decays occur in the interior of a carbon-oxygen white dwarf star 

at extremely high temperatures and densities. Our 

model [2] revives this idea and uses a three-parameter 

Weibull probability density function to model the 

 

Figure 39. The B-magnitude light curve for the supernova 

SN2003du. The time variable is measured in Julian Day units. The 

discrete data points are the observed B-magnitudes and the smooth 

curve is the model fit. 

Figure 40. The Weibull pulse powering the Ni deposition in 

SN1999dq. Here it is normalized to subtend unit area, but it would

have to be multiplied by a scaling constant to model the amount of Ni 

that is actually deposited. The time variable is days since the beginning of the fusion reactions that produced the nickel. 

Figure 41. The relative abundances of Ni, Co, and Fe during the 

evolution of the light curve for SN1999dq. The time scale here is the 

same as in Figure 40. 
82 NIST Mathematical and Computational Sciences Division 

deposition of the Ni which is the product of fusion reactions between the carbon, oxygen and lighter 

elements. An example of this Weibull pulse is plotted 

in Figure 40. 

Plots of the corresponding relative abundances of 

Ni, Co, and Fe are given in Figure 41. These abundances are modeled by a system of three coupled 

ordinary differential equations involving four adjustable parameters. Three of those parameters specify the 

Weibull pulse, and the fourth is the common acceleration that must be applied to the terrestrial decay rates of 

Ni and Co. Two more adjustable fitting parameters are 

required to convert the relative abundances to observed 

luminosities. The final six-parameter model gives 

much better fits to observed light curves than any previous theoretical model. The previous models have all 

posited a complete disruption of the star by the supernova explosion, but this model requires the star to 

remain intact to maintain the high densities needed to 

produce the stable accelerated decay rate observed in 

the long tail of the light curve. None of the observed 

luminosity is produced by the fusion reactions that deposit the Ni. 

[1] B. W. Rust, M. Leventhal and S. L. McCall. “Evidence 

for a radioactive decay hypothesis for supernova luminosity,” Nature 262 (1976), pp. 118-120 

[2] B. W. Rust, D. P. O’Leary, and K. M. Mullen, “Modeling Type Ia Supernova Light Curves,” Exponential Data 

Fitting and its Applications, Victor Pereyra and Godela 

Scherer, Eds., to appear. 
Yearly Report, Fiscal Year 2009 83

High Performance Computing and Visualization 

Computational capability is advancing rapidly. This means that modeling and simulation can be done with 

greatly increased fidelity (e.g. higher resolution, more complex physics). However, developing large-scale 

parallel applications remains highly challenging, requiring expertise that application scientists rarely have. 

In addition, the hardware landscape is changing rapidly so new algorithmic techniques must constantly be 

developed. We are developing and applying facilities and expertise of this type for application to NIST problems. Large scale computations and laboratory experiments invariably produce large volumes of scientific 

data, which cannot be readily comprehended without some form of visual analysis. We are developing the infrastructure necessary for advanced visualization of scientific data, including the use of 3D immersive 

environments and applying this to NIST problems. One of our goals is to develop the 3D immersive environment into a true interactive measurement laboratory. 

Visualizing Fire and the Resulting 

Deformation of Structures 

William Hess 

Marc Olano 

Judith Terrill 

Dilip Banerjee (NIST BFRL) 

Researchers in the NIST Building and Fire Research 

Laboratory (BFRL) are currently studying the effects 

of fire on building structures to improve structural engineering and building safety. We are using a 

sequential process in which the BFRL Fire Dynamics 

Simulator (FDS) is used to simulate the start and 

propagation of fire in a room, after which a second 

code calculates how the heat computed by FDS diffuses into the structural members. A third program 

takes the output of the second code and computes how 

the structure deforms over time due to combined effects of thermal and mechanical loads. We use a test 

case with a single beam in a room. We link the fire, 

thermal, and structural data with a separate real-time 

visualization capability that can be run on the desktop 

and in the NIST 3D immersive environment. This enables interactive analysis of the resulting data. 

 This year we applied volume visualization to the 

beam in the fire, and created measurement tools for it. 

We handle all preprocessing of data with a Python 

script simplifying data and program management. Another Python program reads FDS input files and 

automatically creates the geometry for the visualization 

environment. We implemented a Python script for interpolation, simplifying the transformation of output 

between programs. This script is very general, with the 

ability to handle any number of datasets relating to 

beam. We modified the open source physics engine 

ODE, to get windows to blow out of their frames. We 

created a multi-channel volume visualization capability 

and applied it to the fire and smoke. This visualization 

Figure 42. Multichannel volume visualization of fire and smoke.
84 NIST Mathematical and Computational Sciences Division 

technique draws the bounding box of the volume and 

any clipping planes using a ray-casting shader. The 

shader integrates the ray through the volume until it 

exits the volume through the box, a clipping plane, or 

by hitting one of the solid objects. The colors are computed using a shader function. The color of the flame is 

selected from a pseudo-black body spectrum based on 

the fire temperature. The opacity of the flame comes 

from the heat release rate per unit volume (a measure 

of where the burning is actually happening). The 

smoke is black with an opacity based on the soot density. All three data values are available at each voxel in 

each time step, and combined to create the final color. 

See Figure 42. Individual points can be selected on the 

beam and the graph of temperature vs. time displayed. 

See Figure 43. We developed a specialty GUI to interact with the data that allows all of the above including 

cutting planes. 

[1] D. Banerjee, W. Hess, M. Olano, J. Terrill, and J. Gross, 

“Visualization of Structural Behavior under Fire,” 

NISTIR 7619, August 2009. 

Figure 43. Graph of Temperature Vs Frame number for a selected

point on beam 
Yearly Report, Fiscal Year 2009 85

Fundamental Math Software Development and Testing 

Modern science and engineering in general, and modern measurement science in particular require a wide 

variety of software tools for scientific discovery, exploration, and analysis. As scientific inquiry becomes 

deeper and more specialized, so must the supporting software tools. The goal of this work is to develop critical software tools that support measurement science at NIST, as well as computational science and 

engineering at large.

SciMark, a Web-based Benchmark 

for Numerical Computing in Java 

Roldan Pozo 

Bruce Miller 

http://math.nist.gov/scimark

The NIST SciMark benchmark continues to be one of 

the most widely used Java scientific benchmarks and 

recently became an official component of the industry 

standard SPECjvm2008. In December 2009, a Google 

search on “java benchmark” reveals SciMark as the #1 

web site out of 2.6 million results. The SciMark suite 

consists of computational kernels for FFTs, SOR, 

Monte Carlo integration, sparse matrix multiply, and 

dense LU factorization; these comprise a representative 

set of computational styles commonly found in numerical applications. 

 SciMark can be run interactively from Web browsers, or can be downloaded and compiled for standalone Java platforms. Full source code is provided in 

Java and C programming languages for comparison 

under different compilers and execution environments. 

The SciMark result is recorded as a series of megaflop 

rates for the numerical kernels, as well as an aggregate 

score for the complete benchmark. The current results 

database lists results of more than 3,000 submissions 

from computational platforms ranging from palm devices to high-end servers, representing nearly every 

operating system and virtual machine environment 

currently in use, including Solaris, FreeBSD, MacOS, 

Sun OS, IRIX, OSF1, Linux, OS/2, and Windows 95, 

98, 2K, ME, NT, XP, Vista, and Windows 7. 

 SciMark and its kernel components have become a 

de facto standard in industry and academia. They have 

been adopted by the Java Grande Benchmark Forum, 

and Sun Microsystems has used SciMark 2.0 to demonstrate continued floating-point improvements to their 

Java Virtual Machine (JVM). SciMark is one of 10 

benchmarks included in SPECjvm2008 to represent the 

wide range of existing Java applications. Other 

benchmarks included in the suite assess performance in 

areas such as cryptographic operations, audio decoding, graphics processing, database logic, and XML 

processing. SPECjvm2008 was released on May 7, 

2008 by the System Performance Evaluation Corporation13 (SPEC), a non-profit open industry cooperative 

which develops benchmarks that ensure that the marketplace has a fair and useful set of metrics to 

differentiate system performance. For example, results 

of the well-know SPEC CPU benchmarks (such as 

SPECint, SPECfp, and SPECrate) are typically cited on 

specification sheets used in the marketing of computer 

processors. Members of SPEC include AMD, Dell, 

HP, IBM, Intel, Microsoft, and Sun.

Division Software and Services 

Staff 

MCSD maintains a number of significant services 

in support of the computational science community at 

large. This includes the following: 

o The Guide to Available Mathematical Software6

(GAMS): a cross-index and virtual repository of 

mathematical software components. 

o The Matrix Market7

: a database of matrices for use 

in testing software for linear algebra. 

o The Digital Library of Mathematical Functions8

: 

an online handbook of the properties of the special 

functions of applied mathematics 

o SciMark9

: an benchmark for assessing the performance of Java Virtual Machines for core 

numerical computing tasks 

o The JavaNumerics Web Page10: a directory of resources for numerical computing in Java 

Each such resource requires attention of Division staff 

to add information, make corrections, and respond to 

user queries. 

In addition, MCSD has developed a wide variety 

of software which is made available for public 

download. Among the most popular of these are the 

following (numbers in parenthesis indicate number of 

source/executable downloads for calendar 2009): 

 

6 http://gams.nist.gov/

7 http://math.nist.gov/MatrixMarket/

8 http://dlmf.nist.gov/

9 http://math.nist.gov/scimark/

10 http://math.nist.gov/javanumerics/
86 NIST Mathematical and Computational Sciences Division 

o JAMA: elementary numerical linear algebra in 

Java (17,067). 

o TNT: Template Numerical Toolkit, linear algebra 

components based on C++ templates (15,312). 

o OOMMF: micromagnetic modeling system 

(5,200). Under current development. 

o LAPACK++: C++ version of the popular package 

LAPACK for dense linear algebra (6,307). 

o SparseLib++: sparse linear algebra in C++ 

(4,903). 

o F90gl: standard Fortran interface to OpenGL 

graphics (3,530). 

o REGEN: cryocooler modeling (1,932). 

o Spblas: sparse basic linear algebra subprograms 

(1,569). 

o LaTeXML: LaTeX to XML translator (1,520). 

o OOF: modeling of materials with complex microstructures (1,204). Under current development. 

o PHAML: Parallel Hierarchical Adaptive Multilevel, adaptive multigrid solution of partial differential equations (767). Under current 

development. 

Each such software package requires attention of Division staff to add information, make corrections, and 

respond to user queries, even when the software is no 

longer under current development. 

MCSD maintains a Web server to enable distribution of its products and services. The server continues 

to see high usage. During calendar year 2009, the 

server math.nist.gov satisfied more than 7M requests 

for pages, or some 19,000 per day. More than 1.6 

Gbytes of data was shipped each day. More than 

445,000 distinct hosts were served. There have been 

more than 217M “hits” on MCSD Web servers since 

they went online as NIST first web servers in 1994. 
Yearly Report, Fiscal Year 2009 87

Part IV: Activity Data 

 Part IV 

 Activity Data 
88 NIST Mathematical and Computational Sciences Division 
Yearly Report, Fiscal Year 2009 89

Publications

Appeared 

Refereed Journals 

1. S. Boixo, E. Knill, and R. D. Somma, “Quantum 

State Preparation by Phase Randomization,” Journal of Quantum Information and Computation 9

(2009), pp. 833-855. 

2. S. S. Bullock and D. P. O’Leary, “Locality Bounds 

on Hamiltonians for Stabilizer Codes,” Quantum 

Information and Computation 9 (2009), pp. 0487-

0499. 

3. A. Carasso, “False Characteristic Functions and 

Other Pathologies in Variational Blind Deconvolution. A Method of Recovery,” SIAM Journal of 

Applied Math 70 (4) (2009), pp. 1097–1119. 

4. A. Dienstfrey, T. L. Oreskovic, L. T. Hudson, and 

H. S. Bennett, “Analysis of ISCD-NIST Survey for 

Bone Health,” Journal of Clinical Densitometry 12

(2) (2009), pp. 162-169. 

5. M. J. Donahue, “Parallelizing a Micromagnetic 

Program for Use on Multi-processor Shared Memory Computers,” IEEE Transactions on Magnetics

45 (2009), pp. 3923-3925. 

6. B. Eastin and E. Knill, “Restrictions on Transversal Encoded Quantum Gate Sets,” Physical Review 

Letters 102 (11) (2009), article no. 110502 (4 

pages). 

7. D. E. Gilsinn, “On Algorithms for Estimating 

Computable Error Bounds for Approximate Periodic Solutions of an Autonomous Delay 

Differential Equation,” Communications in 

Nonlinear Science and Numerical Simulation 14

(2009), pp. 1526-1550. 

8. P. D. Hale, A. Dienstfrey, J. C. M. Wang, D. F. 

Williams, A. Lewandowski, D. A. Keenan, and T. 

Clement, “Traceable Waveform Calibration with a 

Covariance-Based Uncertainty Analysis,” IEEE 

Transactions on Instrumentation and Measurement, 58 (10) (2009), pp. 3554-3568. 

9. K. Irikura, R. Johnson, R. Kacker, and R. Kessel, 

“Uncertainties in Scaling Factors for ab-initio Vibrational Zero-Point Energies,” Journal of 

Chemical Physics 130 (11)(2009), pp. 114102-

114111. 

10. R. Kacker, A. Forbes, R. Kessel, and K.-D. Sommer, “Bayesian Posterior Predictive p-value of 

Statistical Consistency in Interlaboratory Evaluations,” Metrologia 45 (2008), pp. 512-523. 

11. R. Kacker and J. Lawrence, “Rectangular Distribution Whose Width is Not Exactly Known: 

Isocurvilinear Trapezoidal Distribution,” Metrologia 46 (2009), pp. 254-260. 

12. R. Kessel and R. Kacker, “A Test of Linearity Using Covering Arrays for Evaluating Uncertainty in 

Measurement,” Advanced Mathematical and 

Computational Tools in Metrology and Testing 

VIII, Series on Advances in Mathematics for Applied Sciences; 78 (2009), Editors F. Pavese, M. 

Bdr, A. B. Forbes, J.-M. Linares, C. Perruchet, N. 

F. Zhang World Scientific, Singapore, pp. 194-

202. 

13. P. E. King-Smith, B. A. Fink, J. J. Nichols, K. K. 

Nichols, R. J. Braun, and G. B. McFadden, “The 

Contribution of Lipid Layer Movement to Tear 

Film Thinning and Breakup,” Investigations in 

Ophthalmology and Visual Science 50 (2009), pp. 

2747-2756. 

14. R. Kuhn, R. Kacker, Y. Lei, and J. Hunter, “Combinatorial Software Testing,” Computer, August 

2009, pp. 94-96. 

15. Z. H. Levine, M. Li, A. P. Reeves, D. F. Yankelevitz, J. J. Chen, E. L. Siegel, A. Peskin, and D. N. 

Zeiger, “A Low-Cost Density Reference Phantom 

for Computed Tomography,” Medical Physics 36

(2009), pp. 286-288. 

16. N. S. Martys, D. Lootens, W. L. George, and P. 

Hebraud, “Contact and Stress Anisotropies in 

Start-up Flow of Colloidal Suspensions,” Physical 

Review E 80 (2009), article no. 031401. 

17. G. B. McFadden and S. R. Coriell, “Onset of Oscillatory Convection in Two Liquid Layers with 

Phase Change,” Physics of Fluids 21 (2009), article no. 034101. 

18. Y. Mishin, W. J. Boettinger, J. A. Warren, and 

G.B. McFadden, “Thermodynamics of Grain 

Boundary Premelting in Alloys. I. Phase-field 

Modeling,” Acta Materialia 57 (2009), pp. 3771-

3785. 

19. M. J. O’Hara and D. P. O’Leary, “The Adiabatic 

Theorem in the Presence of Noise,” Physical Review A 77 (2008), article no. 042319, 20 pages. 

20. C. S. Pande, K. P. Cooper, and G. B. McFadden, 

“Grain Size Distribution in Two Dimensions in the 

Long Time Limit,” Acta Materialia 56 (2008), pp. 

5304-5311. 

21. A. C. E. Reid, R. C. Lua, R. E. Garcia, V. R. 

Coffman, and S. A. Langer, “Modeling Microstructures with OOF2,” International Journal of 

Materials and Product Technology 35 (2009), pp. 

361-373. 
90 NIST Mathematical and Computational Sciences Division 

22. B. W. Rust and D. P. O’Leary, “Residual Periodograms for Choosing Regularization Parameters for 

Ill-Posed Problems,” Inverse Problems 24 (2008), 

article no. 034005 (30 pages). 

23. J. S. Sims and S. A. Hagstrom, “HylleraasConfiguration-Interaction Study of the 22

S Ground 

State of Neutral Lithium and the First Five Excited 

2

S States,” Physical Review A 80 (2009), article 

no. 052507 (8 pages). 

24. W. C. Uhlig, M .J. Donahue, D. T. Pierce and J. 

Unguris, “Direct Imaging of Current-driven Domain Walls in Ferromagnetic Nanostripes,” 

Journal of Applied Physics 105 (2009), article no. 

103902. 

Journal of Research of NIST 

1. P. R. Garabedian and G. B. McFadden, “Design of 

the DEMO Fusion Reactor Following ITER,” 

Journal of Research of the National Institute of 

Standards and Technology 114 (2009), pp. 229-

236. 

2. J. W. Lau, R. D. McMichael, M. J. Donahue, “Implementation of Two-Dimensional Polycrystalline 

Grains in Object Oriented Micromagnetic Framework,” Journal of Research of the National 

Institute of Standards and Technology 114 (2009), 

pp. 57-67. 

Books 

1. B. Balachandran, T. Kalmar-Nagy, D. E. Gilsinn 

(Ed.), Delay Differential Equations: Recent Advances and New Directions, Springer, New York, 

2009. 

2. D. P. O’Leary, “Scientific Computing with Case 

Studies,” SIAM Press, Philadelphia, 2009. 

Book Chapters 

1. D. E. Gilsinn, “Bifurcations, Center Manifolds, 

and Periodic Solutions,” in Delay Differential 

Equations: Recent Advances and New Directions, 

 Balakumar Balachandran, Tamas Kalmar-Nagy 

and David E. Gilsinn (Ed.), Springer, New York, 

2009, pp. 155-202. 

2. J. Terrill, W. George, T. Griffin, J. Hagedorn, J. 

Kelso, M. Olano, A. Peskin, S. Satterfield, J. 

Sims, J. Bullard, J. Dunkers, N. Martys, A. 

O’Gallagher, and G. Haemer “Extending Measurement Science to Interactive Visualization 

Environments,” in Trends in Interactive Visualization: A-State-of-the-Art Survey, pp. 287-302 

(Elena Zudilova-Seinstra, Tony Adriaansen and 

Robert van Liere, Eds.), Springer, UK. 

Conference Proceedings 

1. I. Beichl and B. Cloteaux, “Measuring The Effectiveness Of The s-Metric To Produce Better 

Network Models,” in Proceedings of the Winter 

Simulation Conference (2008), Miami, FL, December 7-10, 2008, pp. 1020-1028. 

2. G. W. Bryant, N. Malkova, J. Sims, M. Zielinski, 

W. Jaskolski, and J. Aizpurua, “Re-engineering the 

Optics of Quantum Dots with Mechanical Strain,” 

in Proceedings of Optics of Excitons in Confined 

Systems (OECS11), Madrid, Spain, Sept. 7-11, 

2009.

3. T. J. Burns, S. P. Mates, R. L. Rhorer, E. P. 

Whitenton, D. Basak, “Dynamic Properties for 

Modeling and Simulation of Machining: An Update of Results from the NIST Pulse-Heated 

Kolsky Bar,” in Proceedings of the 12th CIRP 

Conference on Modeling of Machining Operations, San Sebastian, Spain, May 7-8, 2009, pp. 

255-262. 

4. B. Cloteaux and D. Ranjan, “Lower Bounds for 

Accessing Information on Pure Pointer Machines,” 

in Proceedings of the International Conference on 

the Foundations of Computer Science, Las Vegas, 

NV, July 13-16, 2009, pp. 103-107.

5. C. Dabrowski and F. Hunt, “Using Markov Chain 

Analysis to Study Dynamic Behavior in LargeScale Grid Systems,” in Proceedings of the Seventh Australasian Symposium on Grid Computing 

and e-Research (AusGrid 2009), CRPIT, 99, Wellington, New Zealand, Jan. 20-23, 2009, pp. 29-40. 

6. J. T. Fong, R. deWit, P. V. Marcal, J. J. Filliben, 

and N. A. Heckert, “A Design-of-Experiments 

Plug-In for Estimating Uncertainties in Finite 

Element Simulations,” in Proceedings of 2009 International SIMULIA Conference, London, U.K., 

May 18-21, 2009, pp. 828-842. 

7. J. T. Fong, and P. V. Marcal, “Artificial Intelligence Tools for Failure Event Data Management 

and Probabilistic Risk Analysis for Failure Prevention,” in Proceedings of the Materials Science & 

Technology 2009 Conference, Pittsburgh, PA, October 25-29, 2009, pp. 529-540. 

8. F. Hunt and V. Marbukh, “Measuring the Utility/Path Diversity Trade-Off in Multipath 

Protocols,” in Proceedings of the 4th International 

Conference on Performance Evaluation, Methods 

and Tools (ValueTools09), Pisa, Italy, October 20-

22, 2009 . 
Yearly Report, Fiscal Year 2009 91

9. V. Marbukh and D. Genin, “Towards Understanding Metastability in Large-scale Loss Network 

with Mobile Users: Emergence and Implications 

for Performance,” in Proceedings of the IEEE 

Globecom 2008 Conference, New Orleans, Nov. 

30 – Dec. 4 , 2008. 

10. V. Marbukh, “From Network Microeconomics to 

Network Infrastructure Emergence,” in Proceedings of the First IEEE International Workshop on 

Network Science for Communication Networks 

(NetSciCom’09), Rio de Janeiro, Brazil, April 24, 

2009. 

11. V. Marbukh and K. Sayrafian-Pour, “A Framework for Joint Cross-Layer and Node Location 

Optimization in Mobile Sensor Networks,” in Proceedings of Adhoc-Now’09, Nice, France, 2008. 

12. V. Marbukh and K. Sayrafian-Pour, “Mobile Sensor Network Self-organization for System Utility 

Maximization: Work in Progress,” in Proceedings 

of the International Conference on Wireless and 

Mobile Communications, Cannes, France, 2009. 

13. P. V. Marcal, J. T. Fong, and N. Yamagata, “Artificial Intelligence (AI) Tools for Data Acquisition 

and Probability Risk Analysis of Nuclear Piping 

Failure Databases,” in Proceedings of the International Symposium on PVP Inservice Inspection 

and Failure Event Database," Prague, Czech Republic July 26-30, 2009. 

14. N. S. Martys, W. L. George, S. G. Satterfield, D. 

Lootens, P. Hebraud, and J. Hagedorn, “A 

Smoothed-particle-hydrodynamics-based Fluid 

Model with a Local Shear-rate Dependent Viscosity: Application to Flow of a Suspension with a 

Non-Newtonian Fluid Matrix,” in Proceedings of 

the 2009 Society of Rheology Meeting, Madison, 

WI, October 18-22, 2009. 

15. N. S. Martys, D. Lootens, W. L. George, and P. 

Hebraud, “Contact and Stress Anisotropies in the 

Start-up Flow of Colloidal Suspensions,” in Proceedings of the 2009 Society of Rheology Meeting, 

Madison, WI, October 18-22, 2009. 

16. W. F. Mitchell, “Strategies for hp-Adaptive Refinement,” in Proceedings of the 2008 

International Conference on Numerical Analysis 

and Applied Mathematics, AIP Conference Proceedings 1048, Kos, Greece, Sept. 16-20, 2008, 

pp. 23-25. 

17. A. Peskin, K. Kafadar, A. M. Santos, and G. 

Haemer, “Robust Volume Calculations of Tumors 

of Various Sizes,” in Proceedings of the 2009 International Conference on Image Processing, 

Computer Vision, and Pattern Recognition, Las 

Vegas, NV, July 13-16, 2009. 

18. A. Peskin, K. Kafadar, A. Dima, J. Bernal, and D. 

Gilsinn, “Synthetic Lung Tumor Data Sets for 

Comparison of Volumetric Algorithms,” in Proceedings of the 2009 International Conference on 

Image Processing, Las Vegas, NV, July 13-16, 

2009. 

19. A. Peskin, K. Kafadar, and A. Dima, “A Quality 

Pre-Processor for Biological Cell Images,” Computational Bioimaging 2009, 5th International 

Symposium on Visual Computing, Las Vegas, 

NV, November 30-December 2, 2009. 

20. B. W. Rust, “Atmospheric Retention of Man-made 

CO2 Emissions,” in MAMERN'09: Proceedings of 

the 3rd International Conference on Approximation Methods and Numerical Modeling in 

Environment and Natural Resources, B. 

Amaziane, D. Barrera, M.A. Fortes, et. al., (eds.), 

Universidad de Granada (June 2009), pp. 813-818. 

21. B. Saunders and Q. Wang, “Tensor Product Grid 

Generation for Complex Surface Visualizations in 

a Digital Library of Mathematical Functions,” in 

Proceedings of the 11th International Society for 

Grid Generation (ISGG) Conference, Montreal, 

Canada, May 25-28, 2009. 

22. K. Sayrafian, W. Yang, J. Hagedorn, J. Terrill, and 

K. Y. Yazdandoost, “A Statistical Path Loss 

Model for Implant Communication Channels,” in 

Proceedings of the IEEE 2009 Personal, Indoor, 

and Mobile Radio Communications Symposium, 

Tokyo, Japan, Sept. 13-16, 2009. 

23. X. Tang, L. Ma, A. Mink, T. Chang, H. Xu, O. 

Slattery, A. Nakassis, B. Hershman, D. Su, and R. 

F. Boisvert, “High-Speed Quantum Key Distribution System for Optical Fiber Networks in Campus 

and Metro Areas,” in Quantum Communications 

and Quantum Imaging VI (R. E. Meyers, Y. Shih, 

and K. S. Deacon, eds.), Proceedings of the SPIE

7092, pp. 70920I-1 – 70920I-15 (2008). 

24. H. M. Vasconcelos, L. Sanz, and S. Glancy, “Generation of GKP States with Optical States,” in 

Proceedings of the XXXII Brazilian Meeting on 

Condensed Matter Physics, Aguas de Lindoia, 

Brazil, May 11 - 15, 2009. 

25. W. Wang, S. Sampath, Y. Lei, and R. Kacker, “An 

Interaction-based Test Sequence Generation Approach for Testing Web-applications,” in 

Proceedings of 11-th IEEE High Assurance Systems Engineering Symposium, Nanjing, China, 

December 3-5, 2008. 

26. W. Yang, K. Sayrafian-Pour, J. Hagedorn, J. 

Terrill, “Simulation Study of Body Surface RF 

Propagation for UWB Wearable Medical Sensors,” 

in Proceedings of the IEEE International Sympo-
92 NIST Mathematical and Computational Sciences Division 

sium on Applied Sciences in Biomedical and 

Communication Technologies, November 2009. 

Other Formal Publications 

1. I. Beichl, “Begin at the Beginning,” Computing in 

Science & Engineering 11 (1) (2009), p. 3. 

2. I. Beichl, “We’ll Take Care of All Your Problems,” Computing in Science & Engineering 11 (2) 

(2009), pp. 3-4. 

3. I. Beichl, “Where the Rubber Meets the Road,” 

Computing in Science and Engineering 11 (5) 

(2009), p 3. 

4. I. Beichl and F. Sullivan, “Cut It Out,” Computing 

in Science and Engineering 11 (3) (2009), pp. 74-

79. 

5. W. George, “What Insights Are Revealed When 

Modeling Concrete’s Rheological Properties?,” 

Argonne Leadership Computing Facility Brochure

(2008), p. 32. 

6. W. George, “Modeling the Rheological Properties 

of Concrete,” in Transforming the World through 

Breakthrough Science, Department of Energy, Argonne Leadership Computing Facility Annual 

Report (2008), p. 21. 

7. W. George, “Modeling the Rheological Properties 

of Concrete,” in High End Computing at NASA - 

Highlights of Science and Engineering Achievements on NASA’s Supercomputers 2007-2008, pp. 

112-113. 

Technical Reports 

1. D. Banergee, W. Hess, M. Olano, and J. Terrill, 

“Visualization of Structural Behavior under Fire,” 

NISTIR 7619, August 2009. 

2. C. Dabrowski and F. Hunt, “Markov Chain Analysis for Large Scale Grid Systems” NISTIR 7566, 

April 2009. 

3. D. E. Gilsinn, B. R. Borchardt, and A. Tebbe, “Estimating Volumes of Simulated Lung Cancer 

Nodules,” NISTIR 7571, July 2009. 

4. C. Witzgall and M. Franaszek, “Sensitivities and 

Variances for Fitted Parameters of Spheres,” 

NISTIR 7616, July 2009. 

Visualizations Published 

1. W. George, N. Martys, S. Satterfield, J. Hagedorn, 

M. Olano, and J. Terrill, “Simulation of NonNewtonian Suspensions: Shear Thinning Case,” 

movie presented at the 2009 SciDAC Conference 

Visualization Night., San Diego, CA, June 14-18, 

2009. 

2. W. George, N. Martys, S. Satterfield, J. Hagedorn, 

M. Olano, and J. Terrill, “Simulation of NonNewtonian Suspensions: Shear Thinning Case11,” 

movie, in Hadley Leggett, “Best Science Visualizations Videos of 200912,” Wired Science, on-line 

magazine, August 19, 2009. 

Accepted 

1. I. Beichl and B. Cloteaux, “A Structural Approach 

to the Temporal Modeling of Networks,” in Proceedings of the 2009 Winter Simulation 

Conference. 

2. I. Beichl and B. Cloteaux, “An Approximation 

Algorithm for the Coefficients of the Reliability 

Polynomial,” Congressus Numerantium. 

3. W. J. Boettinger and G. B. McFadden, “Bending 

of a Bimetallic Beam Due to the Kirkendall effect,” Journal of Phase Equilibria and Diffusion. 

4. R. Bryce, Y. Lei, D. R. Kuhn, and R. Kacker, 

“Combinatorial Testing” in Handbook of Research 

on Software Engineering and Productivity Technologies: Implications of Globalization. 

5. J. W. Bullard, E. Enjolras, W. L. George, S. G. 

Satterfield, and J.E. Terrill, “Validating a Scalable 

Parallel Computer Model of Cement Hydration 

and Microstructure Development,” Modeling and 

Simulation in Materials Science and Engineering. 

6. J. W. Bullard, C. F. Ferraris, E. J. Garboczi, N. S. 

Martys, P. E. Stutzman, and J. E. Terrill, “Virtual 

Cement and Concrete,” in Innovations in Portland 

Cement Manufacturing. 

7. R. Datla, R. Kessel, A. Smith, R. Kacker, and D. 

Pollock “Uncertainty Analysis of Remote Sensing 

Optical Sensor Data: Guiding Principles to 

Achieve Metrological Consistency,” International 

Journal of Remote Sensing. 

8. J. T. Fong, R. deWit, P. V. Marcal, J. J. Filliben, 

N. A. Heckert, and S. R. Gosselin, “Design of a 

Python-based Plug-In for Benchmarking Probabilistic Fracture Mechanics Computer Codes with 

Failure Event Data,” Proceedings of the ASME 

Pressure Vessels & Piping Conference. 

9. D. E. Gilsinn, B. Borchardt, and A. Tebbe, “Volume Estimation of Molded Artifacts by B-

 

11 http://www.wired.com/video/simulation-of-nonnewtoniansuspensions/32972797001

12 http://www.wired.com/wiredscience/2009/08/visualizations/6/
Yearly Report, Fiscal Year 2009 93

Splines,” SIAM Conference on Mathematics for 

Industry. 

10. J. A. Gordon, C. L. Holloway, and A. Dienstfrey, 

“A Physical Explanation of Angle-Independent 

Reflection and Transmission Properties of Metafilms/Metasurfaces,” IEEE Transactions on 

Antennas and Propagation. 

11. C. L. Holloway, A. Dienstfrey, E. F. Kuester, J. F. 

O’Hara, A. E. Azad, and A. J. Taylor, “A Discussion on the Interpretation and Characterization of 

Metafilms/Metasurfaces: The Two-dimensional 

Equivalent of Metamaterials,” Metamaterials. 

12. B. Johnasson, R. Kacker, R. Kessel, C. McLean, 

and R. Sriram, “Utilizing Combinatorial Testing to 

Detect Interactions and Optimize a Discrete Event 

Simulation Model for Sustainable Manufacturing,” 

in Proceedings of the 14th Design for Manufacturing and the Life Cycle Conference, San Diego, 

CA, August 30-September 2, 2009. 

13. R. N. Kacker, R. Kessel, K.-D. Sommer, and X. 

Bian, “Comparison of Statistical Consistency and 

Metrological Consistency,” Proceedings of the 19th

World Congress on Fundamental and Applied Metrology. 

14. R. Kessel and R. N. Kacker, “Correlation in Uncertainty of Measurement - A Discussion of State 

of the Art Techniques,” in Proceedings of the 19th

World Congress on Fundamental and Applied Metrology. 

15. S. Kim, M. I. Aladjem, G. B. McFadden, and K. 

W. Kohn, “Predicted Functions of MdmX in Fine 

Tuning the Response of p53 to DNA Damage,” 

PLoS Computational Biology. 

16. P. A. Lott and H. Elman, “Fast Iterative Solver for 

Convection-Diffusion Systems with Spectral Elements,” Numerical Methods in Partial Differential 

Equations. 

17. V. Marbukh and K. Sayrafian, “Mobile Sensor 

Networks Self-organization for System Utility 

Maximization: Work in Progress,” The Fifth International Conference on Wireless and Mobile 

Communications (ICWMC’09). 

18. V. Marbukh and D. Genin, “Bursty Fluid Approximation of TCP for Modeling Internet 

Congestion at the Flow Level,” in Forty-Seventh 

(2009) Annual Allerton Conference on Communication, Control, and Computing. 

19. W. F. Mitchell, “The hp-Multigrid Method Applied to hp-Adaptive Refinement of Triangular 

Grids,” Numerical Linear Algebra with Applications. 

20. J. Obrzut, A. Anopchenko, J. F. Douglas and B. 

W. Rust, “Relaxation and Antiplasticization 

Measurements in Trehalose Glycerol Mixtures -- 

A Model for Protein Preservation,” Journal of 

Non-Crystalline Solids. 

21. R. Radebaugh, Y. Huang, A. O’Gallagher, and J. 

Gary, “Optimization Calculations for a 30 HZ, 4 K 

Regenerator with Helium-3 Working Fluid,” in 

Proceedings of the Cryogenic Engineering Conference, Tucson, AZ, June 28 - July 2, 2009. 

22. B. W. Rust, D. P. O’Leary, and K. M. Mullen, 

“Modeling Type Ia Supernova Light Curves,” Exponential Data Fitting and its Applications, Victor 

Pereyra and Godela Scherer, Eds. 

23. B. Saunders and Q. Wang, “Tensor Product BSpline Mesh Generation for Accurate Surface 

Visualization in the NIST Digital Library of 

Mathematical Functions,” Proceedings of the Seventh International Conference on Mathematical 

Methods for Curves and Surfaces, Tonsberg, Norway, June 26-July 1, 2008. 

24. J. S. Sims, G. W. Bryant, W. L. George, T. G. 

Griffin, J. G. Hagedorn, J. T. Kelso, T. M. Olano, 

A. P. Peskin, S. G. Satterfield, and J. D. Terrill, 

“Computation and Visualization of Nanostructures and Nano-optics,” Computing and Visualization in Science. 

25. W. Wang, Y. Lei, S. Sampath, R. Kacker, R. Kuhn 

and J. Lawrence, “A Combinatorial Approach to 

Building Navigation Graphs for Dynamic Web 

Applications,” Proceedings of 25th IEEE International Conference on Software Maintenance. 

In Review 

1. D. M. Anderson, G. B. McFadden, S. R. Coriell, 

and B.T. Murray, “The Solidification of an Ideal 

Ternary Alloy in a Mushy Layer.” 

2. R. Braun, R. Usha, G. McFadden, T. Driscoll, L. 

P. Cook and P. King-Smith, “Thin Film Dynamics 

on a Prolate Ellipsoid with Application to the Cornea.” 

3. T. Burns, S. Mates, R. Rhorer, E. Whitenton, and 

D. Basak, Dynamic Properties for Modeling and 

Simulation of Machining: Effect of Pearlite to 

Austenite Phase Transition on Flow Stress in AISI 

1075 Steel.” 

4. A. Carasso, “Bochner Subordination, Logarithmic 

Diffusion Equations, and Blind Deconvolution of 

Hubble Space Telescope Imagery and Other Scientific Data.” 
94 NIST Mathematical and Computational Sciences Division 

5. B. Cloteaux, “Approximating the Number of Bases 

for Almost All Matroids.” 

6. D. E. Gilsinn, B. R. Borchardt, and A. Tebbe, 

“Volume Estimation of Molded Artifacts by BSplines.” 

7. R. Kessel, R. Kacker, and K.-D. Sommer, “Metrological Consistency and Consensus Result of 

Multiple Evaluations.” 

8. E. Knill, “Quantum Computing Q&A.” 

9. J. Lawrence, R. Kacker, Y. Lei, D. R. Kuhn, and 

M. Forbes. “Binary Covering Arrays.” 

10. H. Mahboubi, K. Moezzi, A. G. Aghdam, K. Sayrafian-Pour, and V. Marbukh, “Distributed 

Deployment Algorithms for Field Coverage in 

Mobile Sensor Networks.” 

11. H. Mahboubi, K. Moezzi, A. G. Aghdam, K. Sayrafian-Pour, and V. Marbukh, “Self-Deployment 

Algorithms in Coverage Problem for a Network of 

Mobile Sensors with Unidentical Sensing Range.” 

12. H. Mahboubi, A. Momeni, A. G. Aghdam, K. Sayrafian-Pour, and V. Marbukh, “Optimal Target 

Tracking Strategy with Controlled Node Mobility 

in Mobile Sensor Networks.” 

13. N. S. Martys, W. L. George, S. G. Satterfield, and 

M. Olano, “Spatial-temporal Correlations in Startup Flows of Colloidal Suspensions.” 

14. N. S. Martys, W. L. George, B.-W. Chun, and D. 

Lootens, “A Smoothed Particle Hydrodynamics 

Based Fluid Model with a Spatially Dependent 

Viscosity: Application to Flow of a Suspension 

with a Non-Newtonian Fluid Matrix.” 

15. G. B. McFadden, S. R. Coriell, and P. A. Lott, 

“Onset of Convection in Binary Liquid Layers.” 

16. G. B. McFadden and P. R. Garabedian, “Quasiaxisymmetric Design of the DEMO Reactor.” 

17. W. F. Mitchell and M. A. McClain, “A Survey of 

hp-Adaptive Strategies for Elliptic Partial Differential Equations.” 

18. C. S. Pande and G. B. McFadden, “Self-similar 

Grain Size Distribution in Three Dimensions.” 

19. Y. Zhang, E. Knill, and S. Glancy, “The Statistical 

Strength of Experiments to Reject Local Realism 

with Photon Pairs and Inefficient Detectors.” 

Presentations

Invited Talks 

1. R. Boisvert, “Succeeding in a Federal Lab,” Future 

Practitioner Colloquium, INFORMS Annual Meeting, Washington, DC, October 11, 2008. 

2. R. Boisvert, “A Handbook of Mathematical Functions for the Digital Age,” Department of 

Computer Engineering, University of Palermo, Italy, July 17, 2009. 

3. R. Boisvert, “A Handbook of Mathematical Functions for the Digital Age,” International Federation 

of Information Processing Working Group 2.5 

Meeting, Raleigh, NC, August 29, 2009. 

4. G. W. Bryant, N. Malkova, J. Sims, M. Zielinski, 

W. Jaskolski, and J. Aizpurua, “Reengineering the 

Optics of Quantum Dots,” Condensed Matter and 

Surface Science Seminar, Ohio University, Athens, OH, April 2009. 

5. G. W. Bryant, N. Malkova, J. Sims, M. Zielinski, 

W. Jaskolski, and J. Aizpurua, “Reengineering the 

Optics of Quantum Dots,” Atomic Physics Division Seminar, NIST, May 2009. 

6. G. W. Bryant, N. Malkova, J. Sims, M. Zielinski, 

W. Jaskolski, and J. Aizpurua, “Reengineering the 

Optics of Quantum Dots with Mechanical Strain,” 

Electronic Materials Conference, Penn State University, University Park, PA, June 2009. 

7. B. Cloteaux, “Measurements on Networks,” Computer Science Colloquium, New Mexico State 

University, Las Cruces, NM, April 15, 2009. 

8. A. Dienstfrey, “NIST and Virtual Measurement 

Systems,” International Federation of Information 

Processing Working Group 2.5, Raleigh, NC, August 31, 2009. 

9. J. Fong, “Engineers in the Global Economy: Innovative and Symphonic,” Mechanical Engineering 

Department Graduate Course Seminar on Nonlinear Mechanics, Stanford University, Stanford, CA, 

May 15, 2009. 

10. J. Fong, “Engineers in the Global Economy: Innovative and Symphonic,” Mechanical Engineering 

Department Undergraduate Course Seminar on 

Statics, Stanford University, Stanford, CA, May 

15, 2009. 

11. J. Fong, “An Uncertainty Estimation Plug-in for 

Applications in Computational Engineering and 

Science Modeling,” Mathematics and Scientific 

Computing Group Seminar, National Physical 
Yearly Report, Fiscal Year 2009 95

Laboratory, Teddington, England, UK, May 20, 

2009. 

12. J. Fong, “Engineers in the Global Economy: Innovative and Symphonic,” Mechanical Engineering 

Department Seminar, University of South Carolina, Columbia, SC, May 29, 2009. 

13. J. Fong, “Engineers in the Global Economy: Innovative and Symphonic,” Mechanical Engineering 

and Mechanics Department Seminar, Drexel University, Philadelphia, PA, June 5, 2009. 

14. J. Fong, “Engineers in the Global Economy: Innovative and Symphonic,” Civil and Environmental 

Engineering Department Seminar, University of 

British Columbia, Vancouver, BC, Canada, June 

16, 2009. 

15. J. Fong, “Aging Plant Risk-Informed In-Service 

Inspection (RI-ISI) Methods and Artificial Intelligence Tools for Probabilistic Risk Assessment,” 

Applied Physics Division Seminar, Pacific Northwest National Laboratory, Richland, WA, Sept. 

18, 2009. 

16. S. Glancy, “Use and Creation of Superpositions of 

Optical Coherent States,” JILA Atomic Molecular 

Optical Physics Seminar, University of Colorado, 

Boulder, CO, November 12, 2008. 

17. S. Glancy, “Innovations in Maximum Likelihood 

Quantum State Tomography,” Perimeter Institute 

Quantum Information Seminar, Perimeter Institute, 

Waterloo, Canada, Sept. 9, 2009. 

18. F. Hunt, “The Dynamics of a Model of a Computer Protocol,” Mathematical Problems in 

Industry Workshop 2009, University of Delaware, 

Newark, DE, June 15-19, 2009. 

19. J. Kelso, “DIVERSE and OpenSceneGraph,” Virginia Tech, Blacksburg, VA, October 28, 2008. 

20. E. Knill, “How Are Your Qubits?” Time and Frequency Division Seminar, NIST Boulder, Feb. 26, 

2009. 

21. S. Langer, “OOF: Analyzing Material Microstructure,” SIAM Conference on Computational 

Science and Engineering, Costa Mesa, CA, Feb 

19-23, 2007. 

22. W. F. Mitchell, “hp-Adaptive Finite Elements for 

the Schrödinger Equation,” University of Nevada, 

Reno, NV, April 24, 2009. 

23. D. P. O’Leary, “A Noisy Adiabatic Theorem and 

Implications for Quantum Computing,” AMS Joint 

Math Meetings, Washington, DC, January 2009. 

24. D. P. O’Leary “Validating Approximate Solutions 

to Ill-Posed Problems,” George Mason University, 

January 30, 2009. 

25. D. P. O’Leary “Optimization Methods for Image 

Deblurring,” NASA Goddard Space Flight Center 

IS&T Colloquium, Greenbelt, MD, April, 2009. 

26. D. P. O’Leary “Validating Approximate Solutions 

to Ill-Posed Problems,” Washington-Baltimore 

Section of the Society for Industrial and Applied 

Mathematics Meeting, College Park, Maryland, 

May 2009. 

27. R. Radebaugh, Y. Huang, A. O’Gallagher, and J. 

Gary, “Optimization Calculations for a 30 HZ, 4 K 

Regenerator with Helium-3 Working Fluid,” The 

Cryogenic Engineering Conference, Tucson, AZ, 

June 28 - July 2, 2009. 

Conference Presentations 

1. B. Alpert, “Compressive Sampling Reconstruction 

by Subspace Refinement,” Institute for Mathematics and Its Applications, University of Minnesota, 

Workshop on Multi-Manifold Data Modeling and 

Applications, Minneapolis, MN, October 27-30, 

2008. 

2. G. W. Bryant, N. Malkova, J. Sims, M. Zielinski, 

W. Jaskolski and J. Aizpurua, “Re-engineering the 

Optics of Quantum Dots with Mechanical Strain,” 

Electronic Materials Conference, University Park, 

PA, June 24-26, 2009. 

3. G. W. Bryant, N. Malkova, J. Sims, M. Zielinski, 

W. Jaskolski, and J. Aizpurua, “Re-engineering the 

Optics of Quantum Dots with Mechanical Strain,” 

Optics of Excitons in Confined Systems 

(OECS11), Madrid, Spain, Sept. 7-11, 2009. 

4. T. J. Burns, “Dynamic Properties for Modeling 

and Simulation of Machining: An Update of Results from the NIST Pulse-Heated Kolsky Bar,” 

12th CIRP Conference on Modeling of Machining 

Operations, San Sebastian, Spain, May 7-8, 2009. 

5. B. Cloteaux, “Measuring The Effectiveness of the 

s-Metric to Produce Better Network Models,” 

Winter Simulation Conference, Miami, FL, December 9, 2008. 

6. B. Cloteaux, “An Approximation Algorithm For 

The Coefficients Of The Reliability Polynomial,” Southeastern International Conference on 

Combinatorics, Graph Theory, and Computing, 

Boca Raton, FL, March 5, 2009. 

7. V. R. Coffman, “Validating Virtual Polycrystals 

with OOF,” 10th U.S. National Congress on Com-
96 NIST Mathematical and Computational Sciences Division 

putational Mechanics, Columbus, OH, July 17, 

2009. 

8. M. J. Donahue, “Parallelizing a Micromagnetic 

Program for Use on Multi-processor Shared Memory Computers,” Intermag 2009, Sacramento, CA, 

May 5, 2009. 

9. M. J. Donahue and K. M. Lebecki, “Quasi-stable 

Vortex Magnetization Structures in Nanowires 

with Perpendicular Anisotropy,” 7th International 

Symposium on Hysteresis Modeling and Micromagnetics (HMM 2009), Gaithersburg, MD, May 

13, 2009. 

10. B. Eastin, “Restrictions on Transversal Encoded 

Quantum Gate Sets,” Quantum Information Processing (QIP) Conference, Santa Fe, NM, 

January 14, 2009. 

11. B. Eastin, “Restrictions on Transversal Encoded 

Quantum Gate Sets,” Southwest Quantum Information and Technology (SQuInT) 

Conference, Seattle, WA, February 20, 2009. 

12. J. Fong, “A Design-of-Experiments Plug-In for 

Estimating Uncertainties in Finite Element Simulations,” 2009 Simulia Customers Conference, 

London, UK, May 18, 2009. 

13. J. Fong, “A Dataplot-Python-Anlap (DPA) Plug-in 

for High Temperature Mechanical Property Databases to Facilitate Stochastic Modeling of FireStructure Interactions,” ASME 2009 Pressure Vessels and Piping Division Conference, Prague, The 

Czech Republic, July 27, 2009. 

14. J. Fong, “A Web-based Uncertainty Plug-in 

(WUPI) for Fatigue Life Prediction based on NDE 

Data and Fracture Mechanics Analysis,” ASME 

2009 Pressure Vessels and Piping Division Conference, Prague, The Czech Republic, July 27, 

2009 

15. J. Fong, “Artificial Intelligence (AI) Tools for 

Data Acquisition and Probability Risk Analysis of 

Nuclear Piping Failure Databases,” ASME 2009 

Pressure Vessels and Piping Division Conference, 

Prague, The Czech Republic, July 30, 2009 

16. J. Fong, “Design of a Python-based Plug-in for 

Benchmarking Probabilistic Fracture Mechanics 

Computer Codes with Failure Event Data,” ASME 

2009 Pressure Vessels and Piping Division Conference, Prague, The Czech Republic, July 30, 

2009 

17. W. George, S. Satterfield, M. Olano, and J. Terrill, 

“Flow of Suspensions,” movie at the Argonne National Laboratory booth, Supercomputing 2008, 

Austin, TX, November 15-21, 2008. 

18. J. Hagedorn and S. Satterfield, “An Immersive 

Visualization Environment to Study RF Propagation from Medical Implants,” DIVERSE BOF, 

SIGGRAPH 2009, New Orleans, LA, August 3-7, 

2009. 

19. J. Kelso, “Introduction to DIVERSE,” DIVERSE 

Birds of a Feather, SIGGRAPH 2009, New Orleans, LA, Aug. 3-7, 2009. 

20. P. A. Lott, “Computing Steady Flow States using 

Fast Iterative Solvers,” 5th Annual Symposium of 

the Burgers Program for Fluid Dynamics, College 

Park, MD, November 14, 2008. 

21. P. A. Lott, “Computing Steady Flow States using 

Fast Iterative Solvers,” 16th Annual Sigma Xi 

poster session, NIST, Gaithersburg, MD, Feb. 11, 

2009. 

22. P. A. Lott, “Fast Solver for Models of Fluid 

Flow,” MD-DC-VA MAA Section Meeting, Fredericksburg, VA, April 18, 2009. 

23. V. Marbukh, “From Network Microeconomics to 

Network Infrastructure Emergence,” First IEEE 

International Workshop on Network Science for 

Communication Networks (NetSciCom’09), Rio 

de Janeiro, Brazil, April 24, 2009. 

24. V. Marbukh and K. Sayrafian, “Mobile Sensor 

Networks Self-organization for System Utility 

Maximization: Work in Progress,” Fifth International Conference on Wireless and Mobile 

Communications. (ICWMC 2009), Cannes/La 

Bocca, France, August 23-29, 2009. 

25. N. S. Martys, W. George, S. Satterfield, D. 

Lootens, P. Hebraud, and J. Hagedorn, “A 

Smoothed-particle-hydrodynamics-based Fluid 

Model with a Local Shear-rate Dependent Viscosity: Application to Flow of a Suspension with a 

Non-Newtonian Fluid Matrix,” Society for Rheology 81st Annual Meeting, Madison, WI, October 

18-22, 2009. 

26. N.S. Martys, D. Lootens, W.L. George, and P. 

Hebraud, “Contact and Stress Anisotropies in the 

Start-up Flow of Colloidal Suspensions,” 2009 Society of Rheology Meeting, Madison, WI, October 

18-22, 2009. 

27. N. S. Martys, W.L. George, M. Olano, J. Kelso, J. 

Hagedorn, S.G. Satterfield, and J. Terrill, “VCCTL 

Rheology Update,” VCCTL Bi-Annual Meeting, 

Gaithersburg, MD, December 3-4, 2008. . 

28. W. F. Mitchell, “hp-Adaptive Finite Elements in a 

Parallel PDE Solver,” SIAM Conference on Computational Science and Engineering, Miami, FL, 

March 2-6, 2009. 
Yearly Report, Fiscal Year 2009 97

29. W F. Mitchell, “The hp-Multigrid Method Applied 

to hp-Adaptive Finite Elements,” Fourteenth Copper Mountain Conference on Multigrid Methods, 

Copper Mountain, CO, March 22-27, 2009. 

30. M. Olano, “Multi-Channel Volume Visualization 

of Fire and Smoke,” DIVERSE BOF, SIGGRAPH 

2009, New Orleans, LA, Aug. 3-7, 2009. 

31. A. Peskin, K. Kafadar, A. M. Santos, and G. 

Haemer, “Robust Volume Calculations of Tumors 

of Various Sizes,” The 2009 International Conference on Image Processing, Computer Vision, and 

Pattern Recognition, Las Vegas, NV, July 13-16, 

2009. 

32. A. Peskin, K. Kafadar, A. Dima, J. Bernal, and D. 

Gilsinn, “Synthetic Lung Tumor Data Sets for 

Comparison of Volumetric Algorithms,” The 2009 

International Conference on Image Processing, Las 

Vegas, NV, July 13-16, 2009. 

33. D. Porter and M. Donahue, “Exploiting Effective 

Field Time Derivative Information to Improve Accuracy of a Norm-preserving Landau-Lifshitz 

Solver,” MMM2008, Austin, TX, November 14, 

2008. 

34. D. Porter (panelist), “Open Discussion with the 

Tcl/Tk Core Team,” Fifteenth Annual Tcl/Tk Conference, Manassas, VA, October 23, 2008. 

35. D. Porter and M. Donahue, “Semi-analytical LLG 

Solver Using Effective Field Time Derivative,” 7th

International Symposium on Hysteresis Modeling 

and Micromagnetics (HMM 2009), Gaithersburg, 

Maryland, May 12, 2009. 

36. B. Rust, “Atmospheric Retention of Man-made 

CO2 Emissions,” 3rd International Conference on 

Approximation Methods and Numerical Modeling 

in Environment and Natural Resources, Pau, 

France, June 8-11, 2009. 

37. B. Rust, "Using Confidence Ellipsoids to Choose 

the Regularization Parameter for Ill-Posed Problems," SIAM Conference on Applied Linear 

Algebra, Monterey Bay - Seaside, California, October 26-29, 2009. 

38. B. Saunders, “Tensor Product Grid Generation for 

Complex Surface Visualizations in a Digital Library of Mathematical Functions,” 11th 

International Society for Grid Generation (ISGG) 

Conference, Montreal, Canada, May 25-28, 2009. 

39. K. Sayrafian, W. Yang, J. Hagedorn, J. Terrill, K. 

and Y. Yazdandoost, “A Statistical Path Loss 

Model for Implant Communication Channels,” 

IEEE 2009 Personal, Indoor, and Mobile Radio 

Communications Symposium, Tokyo, Japan, Sept 

13-16, 2009. 

40. J. Terrill, “Opportunities at NIST,” Supercomputing 2008 Birds-of-a-Feather: Opportunities in 

Federal NITRD High End Computing Interagency 

Working Group (HEC-IWG) Agencies, Austin, 

TX, November 15-21, 2008. 

Patents

1. Alfred S. Carasso, Singular Integral Image Deblurring Method, Patent Number 7,437,012, October 

14, 2008. 

Software Released

1. DiVisa13: A multi-dimensional visualization tool 

developed for researchers to understand the behavior of their data. Version 8.10.30 – C. Houard and 

J. Hagedorn. 

2. f90gl14: Fortran 90 interface for the OpenGL 

graphics library. Version 1.2.14 (June 9, 2009). – 

W. Mitchell. 

3. NGraph++15: a simplified C++ graph library for 

analysis of complex networks. Version 1.0.1. – R. 

Pozo. 

4. OOF16: Object Oriented Finite Elements, modeling of materials with complex microstructures. 

Versions 2.0.2 (Feb. 16, 2007), 2.0.5a7 (Dec. 1, 

2008), 2.0.5a8 (Dec. 16, 2008), and 2.0.5a9 (Jan. 

26, 2009). – S. Langer. 

5. Pcrawler17: a Python-based web crawler which 

generates a web graph of hyperlinks, together with 

time stamps, SHA checksums, and annotated links 

so that validation of graph can be independently 

performed. Version 0.9 – R. Pozo. 

6. PHAML18: Parallel Hierarchical Adaptive Multi 

Level, general purpose software for 2D elliptic 

partial differential equations. Versions 1.5.0 (October 9, 2008) and 1.6.0 (August 2009). – W. 

Mitchell

7. Rendering Error Metrics19: Methods for assessing 

geometric errors introduced by the rendering of 

three 3D geometric primitive forms: points, line 

 

13 http://math.nist.gov/mcsd/savg/software/divisa/ 14 http://math.nist.gov/f90gl 15 http://math.nist.gov/pozo/ngraph 16 http://www.ctcms.nist.gov/oof/ 17 http://math.nist.gov/~RPozo/complex_datasets.html 18 http://math.nist.gov/phaml/ 19 http://math.nist.gov/mcsd/savg/vis/metrology/index.html
98 NIST Mathematical and Computational Sciences Division 

segments, and triangles. – C. Houard, J. Hagedorn, J. Terrill, A. Peskin, and J. Filliben. 

8. Sparse BLAS20: with added testing functionality 

for upper and lower triangular sparse matrices. 

Version v. 1.0.2 – R. Pozo. 

9. SparseLib++21: a C++ library for sparse matrix 

computations, which included updating various 

modules to the latest ANSI/ISO C++ coding standards, and improved compatibility with multiplatform C++ compilers. Version 1.7 ( ). – R. 

Pozo. 

10. Tcl/Tk22: Tool Command Language and Toolkit. 

Versions 8.6a3 (October 10, 2008), 8.5.5 (October 

15, 2008), 8.6b1 (December 23, 2008), 8.5.6 (December 23, 2008), and 8.5.7 (April 15, 2009). – D. 

Porter.

Conferences, Minisymposia, 

Lecture Series, Shortcourses

MCSD Seminar Series 

1. O. Marques (Lawrence Berkeley National Laboratory), “Eigensolvers for Large Electronic Structure 

Calculations,” December 16, 2008. 

2. J. M. Burke (Boehringer-Ingelheim Pharmaceuticals), “Mathematical and Experimental Single-cell 

Analysis of Caspase Amplification in the Death 

Receptor Network,” December 23, 2008. 

3. A. Patra (University of Buffalo), “Parallelism and 

Adaptivity – A Marriage Not Made in Heaven,” 

February 17, 2009. 

4. G. Dogan, “Fast Methods to Find Optimal Shapes 

in Images,” April 21, 2009. 

5. S. McMaster (University of Maryland), “Advances 

in Coverage-Based Test Suite Reduction,” April 

24, 2009. 

6. R. Barnard (Louisiana State University), “Modeling Background Noise for Denoising in Chemical 

Spectroscopy,” June 30, 2009. 

7. J. Harlim (New York University), “Mathematical 

Strategies for Filtering Turbulent Systems: Sparse 

Observations, Model Errors, and Stochastic Parameter Estimation,” July 1, 2009. 

8. R. C. Peralta (NIST/ITL), “Small Circuits for Linear Functions,” August 24, 2009. 

 

20 http://math.nist.gov/spblas 21 http://math.nist.gov/sparselib 22 http://tcl.sourceforge.net/

9. S. G. Nash (George Mason University), “Multilevel Optimization Methods for Engineering 

Design and PDE-Constrained Optimization,” September 8, 2009. 

10. J. Chung (University of Maryland), “Numerical 

Methods for Large-Scale Ill-Posed Inverse Problems,” September 15, 2009. 

Local Events Organized 

1. M. Donahue, Member, Local Organizing Committee, Hysteresis and Micromagnetic Modeling 

Symposium, Gaithersburg, MD, May 11-14, 2009. 

2. M. Donahue, Member, Program Committee, Hysteresis and Micromagnetic Modeling Symposium, 

Gaithersburg, MD, May 11-14, 2009. 

3. D. Porter, Member, Local Organizing Committee, 

Hysteresis and Micromagnetic Modeling Symposium, Gaithersburg, MD, May 11-14, 2009. 

External Events Organization 

1. R. Boisvert, Member, Program Committee, PSE 

Workshop08: Innovative and Collaborative Problem Solving Environments in Distributed 

Resources, 4th IEEE International Conference on 

e-Science, Indianapolis, December 7-12, 2008. 

2. R. Boisvert, Co-Organizer, IFIP Working Group 

2.5 Symposium on Advances in Numerical Computing and Problem Solving, Raleigh, NC, August 

31 – September 1, 2009. 

3. M. Donahue, Member, Organizing Committee for 

MS10, the 6th SIAM meeting on Mathematical 

Aspects of Materials Science. 

4. J. Fong, Minisymposium Organizer, “Materials 

Property Databases,” ASME Pressure Vessel and 

Piping Conference, Prague, Czech Republic, July 

26-30, 2009. 

5. J. Fong, Minisymposium Organizer, “PVP Inservice Inspection (ISI) and Failure Event Databases. (NDE),” ASME Pressure Vessel and Piping 

Conference, Prague, Czech Republic, July 26-30, 

2009. 

6. J. Fong, Minisymposium Organizer, “Field & 

Laboratory Non-destructive Evaluation (NDE) of 

Residual Stresses,” ASME Pressure Vessel and 

Piping Conference, Prague, Czech Republic, July 

26-30, 2009. 

7. J. Fong, Minisymposium Organizer, “Reliability of 

Continuous NDE Monitoring and Uncertainty Estimation of Defect Detection, Location, and 
Yearly Report, Fiscal Year 2009 99

Sizing,” ASME Pressure Vessel and Piping Conference, Prague, Czech Republic, July 26-30, 

2009. 

8. F. Hunt, Member, Program Planning Committee, 

Mathematics Association of America. 

9. J. Kelso, Organizer, DIVERSE Bird’s of a Feather 

(BOF), SIGGRAPH 2009, New Orleans, LA. 

10. V. Marbukh, Member, Technical Program Committee, Advances in Learning for Networks, ACM 

SIGMETRICS, Seattle, WA, June 15, 2009. 

11. Y. Parker and J. Terrill, NIST Exhibit and Demonstration Booth, Supercomputing 2008, Austin, TX, 

November 15-21, 2008. 

12. Y. Parker, J. Terrill and J. Hagedorn, NIST Exhibit 

and Demonstration Booth, Supercomputing 2009, 

Portland, OR, November 16-20, 2009. 

Other Professional Activities 

Internal 

1. I. Beichl, Member, NIST Research Advisory 

Committee. 

2. R. Boisvert, Member, NIST Scientific Computing 

Steering Group. 

3. R. Boisvert, Member, ITL Diversity Committee. 

4. R. Kacker, Member, NIST Library Advisory 

Committee. 

5. W. Mitchell, Member, ITL Awards Committee. 

6. A. O’Gallagher, Member, ITL Diversity Committee. 

7. A. Peskin, Member, Boulder Research Advisory 

Council. 

External 

Editorial 

1. I. Beichl, Editor-in-Chief, IEEE Computing in 

Science & Engineering. 

2. R. Boisvert, Associate Editor, ACM Transactions 

on Mathematical Software. 

3. R. Boisvert, Editor, Numerical Analysis, Mathematical Software, and Computational Engineering, 

Finance, and Science areas of the Computing Research Repository (CoRR) preprint service hosted 

by www.arXiv.org. 

4. R. Boisvert, Member, Editorial Committee, Special Issue on Innovative and Collaborative 

Problem Solving Environments in Distributed Resources, Journal of Convergence Information 

Technology. 

5. D. Lozier, Member, Editorial Board, Journal of 

Numerical Analysis, Industrial and Applied 

Mathematics. 

6. G. McFadden, Associate Editor, Journal of Crystal 

Growth. 

7. G. McFadden, Associate Editor, Interfaces and 

Free Boundaries. 

8. W. Mitchell, Associate Editor, International Journal of Applied Mathematics and Computational 

Science. 

9. W. Mitchell, Associate Editor, Journal of Numerical Analysis, Industrial and Applied Mathematics. 

10. D. Porter, Editor, November 2009 issue of IEEE 

Transactions on Magnetics (Selected papers from 

the 2009 Hysteresis and Micromagnetic Modeling 

Symposium). 

Boards and Committees 

1. I. Beichl, Member, Electronic Publications and 

Services Committee, IEEE Computer Society. 

2. I. Beichl, Member, IEEE Magazine Operations 

Committee. 

3. R. Boisvert, Chair, International Federation for 

Information Processing’s (IFIP) Working Group 

2.5 on Numerical Software. 

4. R. Boisvert, Co-chair, Publication Board, Association for Computing Machinery (ACM). 

5. R. Boisvert, Ex-Officio Member, ACM Council. 

6. R. Boisvert, Member, ACM Awards Committee. 

7. R. Boisvert, Member, Program Review Committee, Center for Computing Sciences, Institute for 

Defense Analysis. 

8. R. Boisvert, Member, External Review Committee, Computer Science Department, George 

Washington University. 

9. R. Boisvert, Member, Ad-hoc External Review 

Panel, Faculty of Engineering, University of Palermo (Italy). 

10. F. Hunt, Member, MAA Carriage House Program 

Committee. 

11. F. Hunt, Member, Advisory Board, Georgetown 

University’s Masters Program. 
100 NIST Mathematical and Computational Sciences Division 

12. F. Hunt, Emeritus Member, Board of Trustees, 

Bryn Mawr College. 

13. D. Lozier, Secretary, Activity Group on Orthogonal Polynomials and Special Functions, Society for 

Industrial and Applied Mathematics (SIAM). 

14. B. Miller, Member, Math Working Group, World 

Wide Web Consortium (W3C). 

15. D. Porter, Member, Tcl Core Team. 

16. B. Saunders, Member, Selection Committee, Etta 

Z. Falconer Lecture, Mathfest 2009. 

17. J. Terrill, Member, Federal High End Computing 

Implementation Task Force. 

18. J. Terrill, High End Computing Research and Development, and Infrastructure Interagency 

Working Groups, NITRD. 

Reviewing 

1. Division staff members referee manuscripts for a 

wide variety of journals and conferences including 

ACM Transactions on Mathematical Software;

Communications in Nonlinear Science and Numerical Simulations, Computing in Science and 

Engineering, European Physics Letters, IEEE 

Symposium on Information Theory, IEEE Transactions on Magnetics, International Symposium on 

Information Theory, Inverse Problems, Journal of 

Applied Physics, Journal of Computing and Information Science in Engineering, Journal of 

Materials in Civil Engineering, Journal of Physics 

A, Mathematical Reviews, Nature, Measurement 

Science and Technology, New Journal of Physics, 

Optical Engineering, Optics Communications, 

Physics Letters A), Physical Review A, Physical 

Review E, Physical Review Letters, Quantum Information and Computation, Science, SIAM 

Journal on Matrix Analysis and Applications, 

SIAM Journal on Scientific Computing.

2. MCSD staff reviewed proposals for the following 

agencies/competitions during FY 2009: Binational 

Science Foundation, Defense Threat Reduction 

Agency, Department of Energy, Louisiana Faculty 

Pilot Funding Program, National Science Foundation, NIST ARRA Measurement Science and 

Engineering Grants Program, NIST Technology 

Innovation Program, National Science and Engineering Council of Canada, NIST Construction 

Grant Program, U.S.-Israel Binational Science 

Foundation

External Contacts

MCSD staff members make contact with a wide variety 

of organizations in the course of their work. Examples 

of these follow. 

Industrial Labs 

Avia Training 

Boehringer-Ingelheim Pharmaceuticals 

ISI Foundation (Italy) 

Northrop-Grumman 

Numerical Algorithms Group (NAG), Ltd. (UK) 

NVIDIA 

Orbital Sciences Corporation 

Pfizer Global R&D 

RedPrairie 

SGI 

2Wire, Inc. 

Government/Non-profit Organizations 

Air Force Research Laboratory 

Army Research Laboratory 

Association for Computing Machinery 

Centre National de la Recherche Scientifique (France) 

Centrun voor Wiskunde en Informatica (Holland) 

Department of Energy 

Food and Drug Administration 

IEEE 

Instituto Superior Tecnico (Portugal) 

International Institute of IT (India) 

Institute for Defense Analysis 

International Federation for Information Processing 

Lawrence Berkeley National Laboratory 

Lawrence Livermore National Laboratory 

Los Alamos National Laboratory 

Mathematical Association of America 

NASA Goddard Space Flight Center 

National Academy of Sciences 

National Institutes of Health 

National Physical Laboratory (UK) 

National Security Agency 

National Science Foundation 

Naval Research Laboratory 

Oak Ridge National Laboratory 

Perimeter Institute for Theoretical Physics (Canada) 

Sandia National Laboratory 

Society for Industrial and Applied Mathematics 

SP Technical Research Institute (Sweden) 

World Wide Web Consortium 

Universities 

American University of Sharjah (UAE) 

Beloit College 

Boston University 
Yearly Report, Fiscal Year 2009 101

Bryn Mawr College 

California Institute of Technology 

Centre College (Kentucky) 

City University of Hong Kong 

Collège de Polytechnique (France) 

Cornell University 

Courant Institute of Mathematical Sciences 

Dalhousie University (Canada) 

Delft University of Technology (The Netherlands) 

Drexel University 

George Mason University 

George Washington University 

Georgetown University 

Giessen University (Germany) 

Hampton University 

Hood College 

Indiana University 

Indian Institute of Technology (India) 

Instituto Superior Tecnico (Portugal) 

Iowa State University 

Johannes Kepler University (Austria) 

Johns Hopkins University 

Louisiana State University 

Loyola College (Maryland) 

Macalester College 

Millersville University 

North Carolina State University 

New Mexico State University 

New York University 

North Carolina State University 

The Ohio State University 

Penn State University 

Seoul National University (Korea) 

Sichuan University (China) 

Stanford University 

SUNY Binghampton 

Technische Universität Berlin (Germany) 

Tulane University 

Universidad de Cantabria (Spain) 

Universidad de Valladolid (Spain) 

Universidade Federal do Rio Grande do Sul (Brazil) 

Universität Duisberg-Essen (Germany) 

University of Albertay Dundee (UK) 

University of Amsterdam (The Netherlands) 

University of Bristol (UK) 

University of Buffalo 

University of Cambridge (UK) 

University of California at Berkeley 

University of California at Los Angeles 

University of California at San Diego 

University of Catalonia (Spain) 

University of Clermont-Ferrand II (France) 

University of Colorado 

University of Delaware 

University of the District of Columbia 

University of Edinburgh (UK) 

University of Groningen (Germany) 

University of Illinois at Carbondale 

University of Iowa 

University of Karlsruhe (Germany) 

University of Kent (UK) 

University of Leeds 

University of Linköping (Sweden) 

University of Manitoba 

University of Maryland 

University of Maryland – Baltimore County 

University of Missouri – Columbia 

University of Nebraska 

University of Nevada at Reno 

University of New Mexico 

University of Palermo (Italy) 

University of Puerto Rico, Bayamon 

University of Sheffield (UK) 

University of St. Andrews (UK) 

University of Texas at Arlington 

University of Utah 

University of Waikato (New Zealand) 

University of Washington 

University of Wisconsin 

University of Wisconsin-Milwaukee 

Utah State University 

Vrije Universiteit Amsterdam (The Netherlands) 
102 NIST Mathematical and Computational Sciences Division 
Yearly Report, Fiscal Year 2009 103

Part V: Appendices 

 Part V 

 Appendices 
104 NIST Mathematical and Computational Sciences Division 
Yearly Report, Fiscal Year 2009 105

Staff

MCSD consists of full time permanent staff located at NIST laboratories in Gaithersburg, MD and Boulder, CO. This is supplemented with a variety of faculty appointments, guest researchers, postdoctoral 

appointments, and student appointments. The following list reflects all appointments held during any portion of the reporting period. (*) Denotes staff at the NIST Boulder Laboratory. 

Division Staff 

Ronald Boisvert, Chief, Ph.D. (Computer Science), Purdue University, 1979 

Robin Bickel, Secretary 

Jeffrey Fong, Ph. D. (Applied Mechanics and Mathematics), Stanford University, 1966 

Roldan Pozo, Ph.D. (Computer Science), University of Colorado at Boulder, 1991 

Christopher Schanzle, B.S. (Computer Science), University of Maryland – Baltimore County, 1989 

 Guest Researchers 

Barry Bernstein / Illinois Institute of Technology 

Mathematical Modeling Group 

Geoffrey McFadden, Leader, Ph.D. (Mathematics), New York University, 1979 

*Bradley Alpert, Ph.D. (Computer Science), Yale University, 1990 

Timothy Burns, Ph.D. (Mathematics), University of New Mexico, 1977 

Alfred Carasso, Ph.D. (Mathematics), University of Wisconsin, 1968 

*Andrew Dienstfrey, Ph.D. (Mathematics), New York University, 1998 

Michael Donahue, Ph.D. (Mathematics), Ohio State University, 1991 

Fern Hunt, Ph.D. (Mathematics), New York University, 1978 

Raghu Kacker, Ph.D. (Statistics), Iowa State University, 1979 

Anthony Kearsley, Ph.D. (Computational and Applied Mathematics), Rice University, 1996 

Peter Ketcham. M.S. (Mathematics), University of Minnesota, 1997 

Stephen Langer, Ph.D. (Physics), Cornell University, 1989 

*Agnes O’Gallagher, M.S. (Applied Math), University of Colorado at Boulder, 1991 

Donald Porter, Ph.D. (Electrical Engineering), Washington University, 1996 

 NRC Postdoctoral Associates 

James Benson, Ph.D. (Mathematics), University of Missouri – Columbia, 2009 

Valerie Coffman, Ph.D. (Physics), Cornell University, 2006 

Aaron Lott, Ph.D. (Mathematics), University of Maryland, 2008 

Contractors 

Andrew C.E. Reid, Ph.D. (Physics), Queen’s University, Kingston, Ontario, 1994 

 Faculty Appointees 

Daniel Anderson / George Mason University 

Dianne O’Leary / University of Maryland College Park 

Florian Potra / University of Maryland Baltimore County 

Guest Researchers 

Mirit Aladjem / National Institutes of Health 

Richard Braun / University of Delaware 

David Cotrell / Lawrence Livermore National Laboratory 

Gunay Dogan / Cornell University 

*John Gary / NIST (retired) 
106 NIST Mathematical and Computational Sciences Division 

Katharine Gurski / George Washington University 

Rüdiger Kessel / Germany 

Sohyoung Kim / National Institutes of Health 

Yu (Jeff) Lei / University of Texas at Arlington 

Bruce Murray / SUNY Binghampton

Students 

Katherine Morrison / University of Nebraska 

Michael Sharpnack / NYU 

Mathematical Software Group 

Daniel Lozier, Leader, Ph.D. (Applied Mathematics), University of Maryland, 1979 

Marjorie McClain, M.S. (Mathematics), University of Maryland College Park, 1984 

Bruce Miller, Ph.D. (Physics), University of Texas at Austin, 1983 

William Mitchell, Ph.D. (Computer Science), University of Illinois at Urbana-Champaign, 1988 

Bert Rust, Ph.D. (Astronomy), University of Illinois at Urbana-Champaign, 1974 

Bonita Saunders, PhD (Mathematics), Old Dominion University, 1985 

Faculty Appointees 

Frank Olver / University of Maryland College Park 

G.W. Stewart / University of Maryland College Park 

Abdou Youssef / George Washington University 

Guest Researchers 

Leonard Maximon / George Washington University 

Qiming Wang / NIST (retired) 

Students 

Russell McFadden / St. Olaf College 

Menal Modha / UMBC 

Discrete Mathematical Analysis Group 

Ronald Boisvert, Acting Leader 

 Isabel Beichl, Ph.D. (Mathematics), Cornell University, 1981 

 Javier Bernal, Ph.D. (Mathematics), Catholic University, 1980 

Brian Cloteaux, Ph.D. (Computer Science), New Mexico State University, 2007 

David Gilsinn, Ph.D. (Mathematics), Georgetown University, 1969 

*Scott Glancy, Ph.D. (Physics), University of Notre Dame, 2003 

 *Emanuel (Manny) Knill, Ph.D., (Mathematics), University of Colorado at Boulder, 1991 

 Vladimir Marbukh, Ph.D. (Mathematics) Leningrad Polytechnic University, 1986 

NRC Postdoctoral Associates 

*Bryan Eastin, Ph.D. (Physics), University of New Mexico, 2007 

Faculty Appointees 

Saul Gass / University of Maryland College Park 

James Lawrence / George Mason University 

Guest Researchers 

Stephen Bullock / IDA Center for Computing Sciences 

Sita Ramamurti / Trinity University, DC 

Francis Sullivan / IDA Center for Computing Sciences 

Christoph Witzgall / NIST Scientist Emeritus 
Yearly Report, Fiscal Year 2009 107

Students 

*Michael Mullan / University of Colorado at Boulder 

*Adam Meier / University of Colorado at Boulder 

Nadezhda Serova / UMBC 

*Yanbao Zhang / University of Colorado at Boulder 

Scientific Applications and Visualization Group 

Judith Terrill, Leader, Ph.D. (Information Technology), George Mason University, 1998 

Yolanda Parker, Office Manager

Robert Bohn, Ph.D. (Physical Chemistry), University of Virginia, 1991 

William George, Ph.D. (Computer/Computational Science), Clemson University, 1995 

Terence Griffin, B.S. (Mathematics), St. Mary’s College of Maryland, 1987 

John Hagedorn, M.S. (Mathematics), Rutgers University, 1980 

John Kelso, M.S. (Computer Science), George Washington University, 1984 

*Adele Peskin, Ph.D. (Chemical Engineering), University of Colorado at Boulder, 1985 

Steven Satterfield, M.S. (Computer Science), North Carolina State University, 1975 

James Sims, Ph.D. (Chemical Physics), Indiana University, 1969 

Faculty Appointees 

Marc Olano / University of Maryland Baltimore County 

Guest Researchers 

Aurelei Canale / ISIMA Campus de Clermont-Ferrand, France 

Dong Yeon Cho / Seoul National University, Korea 

Edith Enjolras / Université Blaise Pascal, France 

Pradeep Gaddam / Core Projects and Technologies PVT LTD, India 

Cedric Houard / Université Blaise Pascal Clermont Ferrand II, France 

Julien Lancien / Institut National des Telecommunications, France 

Audrey Lemoussu / ISIMA Campus de Clermont-Ferrand, France 

Alexandre Marie / ISIMA Campus de Clermont-Ferrand, France 

Jean-Loup Traore / Université Blaise Pascal Clermont Ferrand II, France 

Students 

Matthew Du / Thomas S. Wooton High School, Maryland 

Julie Eisenberg / Damascus High School, Maryland 

*Daniel Leifer / Stanford University 

Sei Masouka / Winston Churchill High School, Maryland 

Sanjeev Muralidharan / Montgomery Blair High School, Maryland 

Joseph Owen / Montgomery Blair High School, Maryland 

Jennifer Schiller / Reservoir High School 

* NIST Boulder Laboratory 
108 NIST Mathematical and Computational Sciences Division 

Glossary of Acronyms

ACM Association for Computing Machinery 

ACTS Advanced Combinatorial Test Suites (software) 

AISI American Iron and Steele Institute 

ANSI American National Standards Institute 

APS American Physical Society 

ASME American Society of Mechanical Engineers 

BAN body area network 

BFRL NIST Building and Fire Research Laboratory 

BLAS Basic Linear Algebra Subprograms 

CIRP International Academy for Production Engineering 

CiSE Computing in Science & Engineering 

CMM coordinate measuring machine 

CMOS complementary metal-oxide semiconductor 

CNST NIST Center for Nanoscale Science and Technology 

CODATA Committee on Data for Science and technology 

CPU central processing unit 

CSTL NIST Chemical Science and Technology Laboratory 

CT computed tomography 

DARPA DOD Defense Advanced Research Projects Agency 

DIVERSE Device Independent Virtual Environments ― Reconfigurable, Scalable, Extensible (visualization 

software) 

DLMF Digital Library of Mathematical Functions 

DNA deoxyribonucleic acid 

DOD U.S. Department of Defense 

DOE U.S. Department of Energy 

DPD dissipative particle dynamics 

EEEL NIST Electronics and Electrical Engineering Laboratory 

EM electromagnetic 

EM expectation maximization 

FDA Food and Drug Administration 

FDS Fire Dynamics Simulator 

FFT fast Fourier transform 

FY fiscal year 

F90gl Fortran 90 interface to OpenGL graphics standard 

GAMS Guide to Available Mathematical Software 

GPU graphics processing unit 

GSTC generalized sheet transition conditions 

GUI graphical user interface 

HAVA Help America Vote Act 

HEC high-end computing 

HMF Handbook of Mathematical Functions 

HTML hypertext markup language 

Hy-CI Hylleraas-Configuration Interaction technique 

IDA Institute for Defense Analysis 

IEEE Institute of Electronics and Electrical Engineers 

IFIP International Federation for Information Processing 

IMS Innovations in Measurement Science 

INCITE Innovative and Novel Computational Impact on Theory and Experiment (DOE Program) 

ISO International Standards Organization 

IT information technology 

ITL NIST Information Technology Laboratory 

IWG Interagency Working Group 

JAMA Java Matrix package 

LANL Los Alamos National Laboratory 
Yearly Report, Fiscal Year 2009 109

LAPACK++ software for numerical linear algebra in C++ 

LaTeXML LaTeX to MathML translator 

MANET mobile ad-hoc network 

MCSD ITL Mathematical and Computational Sciences Division 

MEL NIST Manufacturing Engineering Laboratory 

MICS medical implant communication system 

MKM mathematical knowledge management 

MMM magnetism and magnetic materials 

MNS minimum norm solution 

MPI Message Passing Interface 

MRAM magneto-resistive random access memory 

MRI magnetic resonance imaging 

MSEL NIST Materials Science and Engineering Laboratory 

μmag Micromagnetics Activity Group 

NASA National Aeronautics and Space Administration 

NBS National Bureau of Standards 

NIH National Institutes of Health 

NIST National Institute of Standards and Technology 

NISTIR NIST Internal Report 

NITRD Networking and Information Technology Research and Development 

NRC National Research Council 

NSF National Science Foundation 

NUM network utility maximization 

OASCR DOE Office of Advanced Scientific Computing Research 

ODE ordinary differential equation 

OLES NIST Office of Law Enforcement Standards 

OOF Object-Oriented Finite Elements 

OOMMF Object-Oriented Micromagnetic Modeling Framework 

OSPF open shortest path first 

PDE partial differential equation 

PHAML Parallel Hierarchical Adaptive Multi Level (software) 

PIMRC Personal, Indoor and Mobile Radio Communications Symposium 

PL NIST Physics Laboratory 

PREP Professional Research Experience Program 

QDPD quarternion-based dissipative particle dynamics 

QIS quantum information science 

R&D research and development 

PDE partial differential equation 

REGEN software for modeling of cryocoolers 

RF radio frequency 

SAVG MCSD Scientific Applications and Visualization Group 

SciDAC Scientific Discovery Through Advanced Computing (DOE Program) 

SECB slow evolution from the continuation boundary 

SIAM Society for Industrial and Applied Mathematics 

SIGGRAPH ACM Special Interest Group on Graphics 

SIMA NIST Systems Integration for Manufacturing Applications Program 

SparseLib++ software for spare linear algebra in C++ 

Spblas sparse basic linear algebra subprograms 

SPEC System Performance Evaluation Corporation 

SPIE International Society for Optical Engineering 

SOR successive over-relaxation 

SUNY State University of New York 

SURF Student Undergraduate Research Fellowship 

TCP ` transmission control protocol 

TNT Template Numerical Toolkit 

TGDC Technical Guidelines Development Committee 

UCLA University of California at Los Angeles 
110 NIST Mathematical and Computational Sciences Division 

UMBC University of Maryland Baltimore County 

UMCP University of Maryland College Park 

URL universal resource locator 

VCCTL Virtual Cement and Concrete Testing Laboratory 

VRML virtual reality modeling language 

W3C World Wide Web Consortium 

XML Extensible Markup Language 

ZPE zero-point energy