We present a numerically-stable parallel-in-time linear Kalman smoother. The smoother uses a novel highly-parallel QR factorization for a class of structured sparse matrices for state estimation, and an adaptation of the SelInv selective-inversion algorithm to evaluate the covariance matrices of estimated states. Our implementation of the new algorithm, using the Threading Building Blocks (TBB) library, scales well on both Intel and ARM multi-core servers, achieving speedups of up to$47 x$on 64 cores. The algorithm performs more arithmetic than sequential smoothers; consequently it is 1.8x to 2.5x slower on a single core. The new algorithm is faster and scales better than the parallel Kalman smoother proposed by Särkkä and García-Fernández in 2021.

Alternative basis matrix multiplication algorithms are the fastest matrix multiplication algorithms in practice to date. However, are they numerically stable?We obtain the first numerical error bound for alternative basis matrix multiplication algorithms, demonstrating that their error bounds are asymptotically identical to the standard fast matrix multiplication algorithms, such as Strassen’s. We further show that arithmetic costs and error bounds of alternative basis algorithms can be simultaneously and independently optimized. Particularly, we obtain the first fast matrix multiplication algorithm with a 2-by-2 base case that simultaneously attains the optimal leading coefficient for arithmetic costs and optimal asymptotic error bound, effectively beating the Bini and Lotti (1980) speed-stability trade-off for fast matrix multiplication. We provide high-performance parallel implementations of our algorithms with benchmarks that show our algorithm is on par with the best in class for speed and with the best in class for stability. Finally, we show that diagonal scaling stability improvement techniques for fast matrix multiplication are as effective for alternative basis algorithms, both theoretically and empirically. These findings promote the use of alternative basis matrix multiplication algorithms in practical applications.

UltimateKalman is a flexible linear Kalman filter and smoother implemented in three popular programming languages: MATLAB, C, and Java. UltimateKalman is a slight simplification and slight generalization of an elegant Kalman filter and smoother that was proposed in 1977 by Paige and Saunders. Their algorithm appears to be numerically superior and more flexible than other Kalman filters and smoothers, but curiously has never been implemented or used before. UltimateKalman is flexible: it can easily handle time-dependent problems, problems with state vectors whose dimensions vary from step to step, problems with varying numbers of observations in different steps (or no observations at all in some steps), and problems in which the expectation of the initial state is unknown. The programming interface of UltimateKalman is broken into simple building blocks that can be used to construct filters, single or multi-step predictors, multi-step or whole-track smoothers, and combinations. The article describes the algorithm and its implementation as well as a test suite of examples and tests.

We describe the design and implementation of Vildehaye, a family of versatile, widely-applicable, and field-proven tags for wildlife sensing and radio tracking. The family includes 6 distinct hard-ware designs for tags, 3 add-on boards, a programming adapter, and base stations; modular firmware for tags and base stations (both standalone low-power embedded base stations and base stations tethered to a computer running Linux or Windows); and desk-top software for programming and configuring tags, monitoring tags, and downloading and processing sensor data. The tags are versatile: they support multiple packet formats, data rates, and frequency bands; they can be configured for minimum mass (down to less than 1 g), making them applicable to a wide range of flying and terrestrial animals, or for inclusion of important sensors and large memories; they can transmit packets compatible with time-of-arrival transmitter-localization systems, tag identification and state packets, and they can reliably upload sensor data through their radio link. The system has been designed, upgraded, and main-tained as an academic research project, but it has been extensively used by 5 different groups of ecologists in 4 countries over a period of 5 years. More than 7100 tags have been produced and most of these have been deployed. Production used 41 manufacturing runs. The tags have been used in studies that so far resulted in 9 scientific publications in ecology (including in Science). The paper describes innovative design aspects of Vildehaye, field-use experiences, and lessons from the design, implementation, and maintenance of the system. Both the hardware and software of the system are open.

Abstract We present robust high‐performance implementations of signal‐processing tasks performed by a high‐throughput wildlife tracking system called ATLAS. The system tracks radio transmitters attached to wild animals by estimating the time of arrival of radio packets to multiple receivers (base stations). Time‐of‐arrival estimation of wideband radio signals is computationally expensive, especially in acquisition mode (when the time of transmission is not known, not even approximately). These computations are a bottleneck that limits the throughput of the system. We developed a sequential high‐performance CPU implementation of the computations a few years back, and more recently a GPU implementation. Both strive to balance performance with simplicity, maintainability, and development effort, as most real‐world codes do. The article reports on the two implementations and carefully evaluates their performance. The evaluations indicates that the GPU implementation dramatically improves performance and power‐performance relative to the sequential CPU implementation running on a desktop CPU typical of the computers in current base stations. Performance improves by more than 50X on a high‐end GPU and more than 4X with a GPU platform that consumes almost 5 times less power than the CPU platform. Performance‐per‐Watt ratios also improve (by more than 16X), and so do the price‐performance ratios.

We develop and evaluate parallel algorithms for a fundamental problem in numerical computing, namely the evaluation of a polynomial of a matrix. The algorithm consists of many building blocks that can be assembled in several ways. We investigate parallelism in individual building blocks, develop parallel implemenations, and assemble them into an overall parallel algorithm. We analyze the effects of both the dimension of the matrix and the degree of the polynomial on both arithmetic complexity and on parallelism, and we consequently propose which variants use in different cases. Our theoretical results indicate that one variant of the algorithm, based on applying the Paterson-Stockmeyer method to the entire matrix, parallelizes very effectively on virtually any matrix dimension and polynomial degree. However, it is not the most efficient from the arithmetic complexity viewpoint. Another algorithm, based on the Davies-Higham block recurrence is much more efficient from the arithmetic complexity viewpoint, but one of its building blocks is serial. Experimental results on a dual-socket 28-core server show that the first algorithm can effectively use all the cores, but that on high-degree polynomials the second algorithm is often faster, in spite of the sequential phase. This indicates that our parallel algorithms for the other phases are indeed effective.

Rendezvous is a fundamental building block in distributed cognitive radio networks (CRNs), where users must find a jointly available channel. Research on the rendezvous problem has focused so far on minimizing the time to rendezvous (to find a suitable channel) or on maximizing the degree (number of channels on which rendezvous can take place). In this paper, we model the rendezvous problem in a more realistic way that acknowledges the fact available channels may suffer from interference, and interference may vary among users in different locations over time. In this setting, CRNs benefit from rendezvous methods that find a quiet channel, which supports high symbol rates and does not suffer much from dropped packets. We propose algorithms that achieve this goal for both initial rendezvous problem (users share no prior information) and continuous rendezvous problem (users who have already established a link must vacate the channel and seek another). We propose both deterministic and randomized methods based on mapping the channel set to a larger set in a way that gives preference to quiet channels. This technique allows us to add interference-awareness to existing rendezvous algorithms. We analyze the new algorithms and substantiate our analyses through extensive simulations.

Most species of birds and bats must be tracked with tracking tags weighing less than 10g and many require tags weighing less than 1g. Tags based on commodity internet-of-things system-on-chips (SoC) can be be mass produced at low-cost, hence allowing many individuals to be tracked. We report on the design and performance of two communication protocols that enable long-range communication with such tags. One is a unidirectional protocol, in which tags transmit unique codes that can be reliably detected from 15km away and that can be used for time-of-arrival and angle-of-arrival localization (tracking). The other is a bidirectional protocol that allows tags to transmit short data packets to low-power low-cost basestations and to receive commands from them. Data packets in this protocol can be reliably received from tags that are 8km away and sometimes from up to 15km, and commands packets can be received by tags from up to 4km away. These protocols have been implemented in low-cost tags that can weigh less than 1g (depending on the choice of battery) and using only about 60uJ per transmission. Our results have been gathered by tagging wild bats. The same tags have been used for time-of-arrival localization of wild bats and birds by several different research groups in 3 countries.

We discuss receiver clock correction and associated performance bounds for passive emitter localization using TDOA measurements from asynchronous sensor networks. In the considered system, passive receiving sensors are augmented with beacons at known locations that perform approximately periodical transmissions, used for calibration of the system synchronization. The precise transmission times as well as the transmission interval of the beacon messages are unknown. Similarly the transmissions of the target are irregular and do not, in general, occur simultaneous with the beacon transmissions. The clocks of each sensor are described by a linear error model. Based on that, we compare different snapshot based approaches for clock correction and derive theoretical limits for the localization of the target based on the modified Cramer-Rao lower bound (MCRLB). Simulation results are presented to illustrate the theoretical findings and show that our proposed estimators perform close to the MCRLB. Additional experimental results verify the analysis and the approach in a realistic large scale scenario.

We present a novel algorithm for self-localization in sensor networks without any prior knowledge on the locations of the sensors. We assume that all sensors in the network can receive and transmit, thus we obtain time difference of arrival measurements for all combinations of sensors. Using the full set of these differences in arrival times in the network we are able to obtain the relative location of the sensors nodes, the shape of the network. This leaves us with the problem of anchoring the network to its absolute location, which we solve using additional transmitting beacons at known locations. Experimental results from numerical simulation demonstrate the performance of our approach under various conditions.

We propose a novel high-level abstraction for real-time control, called Proper Timed I/O (PTIO). The abstraction allows user-space programs running on a stock operating system (without real-time extensions) to perform high-resolution real-time digital I/O (setting pins high or low, responding to input transitions, etc.). PTIO programs express their real-time I/O behavior in terms of a timed automaton that can communicate with the user-space program. Simple behaviors are encoded in the timed automaton; complex behaviors are implemented by the user-space program. We present two implementations of the PTIO abstraction, both for Linux. One utilizes a deterministic co-processor that is available on some ARM-based system-on-a-chip processors. This implementation can achieve timing accuracy of 100ns or better and can perform millions of finite-state transitions per second. The other implementation uses hardware timers that are available on every system-on-a-chip; it achieves a timing accuracy of 6µs or better, but it is limited to about 2000 state transitions per second. Both implementations guarantee that the PTIO never fails silently: if the mechanism missed a deadline, the user space program is always notified. In many cases, PTIOs eliminate the need for bare-metal programming or for specialized real-time operating systems.

We present a fast algorithm for approximate Canonical Correlation Analysis (CCA). Given a pair of tall-and-thin matrices, the proposed algorithm first employs a randomized dimensionality reduction transform to reduce the size of the input matrices, and then applies any standard CCA algorithm to the new pair of matrices. The algorithm computes an approximate CCA to the original pair of matrices with provable guarantees, while requiring asymptotically less operations than the state-of-the-art exact algorithms.

Factorization of a dense symmetric indefinite matrix is a key computational kernel in many scientific and engineering simulations. However, there is no scalable factorization algorithm that takes advantage of the symmetry and guarantees numerical stability through pivoting at the same time. This is because such an algorithm exhibits many of the fundamental challenges in parallel programming like irregular data accesses and irregular task dependencies. In this paper, we address these challenges in a tiled implementation of a blocked Aasen's algorithm using a dynamic scheduler. To fully exploit the limited parallelism in this left-looking algorithm, we study several performance enhancing techniques; e.g., parallel reduction to update a panel, tall-skinny LU factorization algorithms to factorize the panel, and a parallel implementation of symmetric pivoting. Our performance results on up to 48 AMD Opteron processors demonstrate that our implementation obtains speedups of up to 2.8 over MKL, while losing only one or two digits in the computed residual norms.

Parallel algorithms for sparse matrix-matrix multiplication typically spend most of their time on inter-processor communication rather than on computation, and hardware trends predict the relative cost of communication will only increase. Thus, sparse matrix multiplication algorithms must minimize communication costs in order to scale to large processor counts.

High performance for numerical linear algebra often comes at the expense of stability. Computing the LU decomposition of a matrix via Gaussian Elimination can be organized so that the computation involves regular and efficient data access. However, maintaining numerical stability via partial pivoting involves row interchanges that lead to inefficient data access patterns. To optimize communication efficiency throughout the memory hierarchy we confront two seemingly contradictory requirements: partial pivoting is efficient with column-major layout, whereas a block-recursive layout is optimal for the rest of the computation. We resolve this by introducing a shape morphing procedure that dynamically matches the layout to the computation throughout the algorithm, and show that Gaussian Elimination with partial pivoting can be performed in a communication efficient and cache-oblivious way. Our technique extends to QR decomposition, where computing Householder vectors prefers a different data layout than the rest of the computation.

This paper considers the problem of scheduling streaming applications on uniprocessors in order to minimize the number of cache-misses. Streaming applications are represented as a directed graph (or multigraph), where nodes are computation modules and edges are channels. When a module fires, it consumes some data-items from its input channels and produces some items on its output channels. In addition, each module may have some state (either code or data) which represents the memory locations that must be loaded into cache in order to execute the module. We consider synchronous dataflow graphs where the input and output rates of modules are known in advance and do not change during execution. We also assume that the state size of modules is known in advance.

We present a design for a high-performance low-cost solid-state disk (SSD). Ignoring garbage-collection costs, our SSD performs only 1 + ε physical accesses to NAND flash pages for every request of a page-size block by the host, for some small ε. This is true for all access patterns, including random writes, which are usually slow on low-cost SSDs. Garbage collection in all SSDs is determined primarily by how full the SSD is, and its cost is similar in most SSDs. The unique feature in our design is that it achieves high performance even with when the SSD contains only a small amount of RAM. In most SSD designs, this would imply low performance; in ours, it does not. A small RAM lowers the cost of an SSD with a given flash array. Our design achieves high performance with a small RAM using two innovative ideas. One is the use of a clever mapping data structure. The second is a host-assisted hinting mechanism that uses RAM on the host to compensate for the small amount of RAM within the SSD. This mechanism is implemented as an enhanced SCSI driver (kernel module). Our prototyping methodology is also a significant contribution. We simulate the SSD in software, using files to represent the flash array, but the resulting prototype is a working SCSI device that file systems can be mounted on.

We analyze the convergence of randomized trace estimators. Starting at 1989, several algorithms have been proposed for estimating the trace of a matrix by 1/MΣ i =1 M z i T Az i , where the z i are random vectors; different estimators use different distributions for the z i s, all of which lead to E (1/MΣ i =1 M z i T Az i ) = trace( A ). These algorithms are useful in applications in which there is no explicit representation of A but rather an efficient method compute z T Az given z . Existing results only analyze the variance of the different estimators. In contrast, we analyze the number of samples M required to guarantee that with probability at least 1-δ, the relative error in the estimate is at most ϵ. We argue that such bounds are much more useful in applications than the variance. We found that these bounds rank the estimators differently than the variance; this suggests that minimum-variance estimators may not be the best. We also make two additional contributions to this area. The first is a specialized bound for projection matrices, whose trace (rank) needs to be computed in electronic structure calculations. The second is a new estimator that uses less randomness than all the existing estimators.

Flash memories are widely used in computer systems ranging from embedded systems to workstations and servers to digital cameras and mobile phones. The memory cells of flash devices can only endure a limited number of write cycles, usually between 10,000 and 1,000,000. Furthermore, cells containing data must be erased before they can store new data, and erasure operations erase large blocks of memory, not individual cells. To maximize the endurance of the device (the amount of useful data that can be written to it before one of its cells wears out), flash-based systems move data around in an attempt to reduce the total number of erasures and to level the wear of the different erase blocks. This data movement introduces an interesting online problem called the wear-leveling problem . Wear-leveling algorithms have been used at least since 1993, but they have never been mathematically analyzed. In this article we analyze the two main wear-leveling problems. We show that a simple randomized algorithm for one of them is essentially optimal both in the competitive sense and in the absolute sense (our competitive result relies on an analysis of a nearly-optimal offline algorithm). We show that deterministic algorithms cannot achieve comparable endurance. We also analyze a more difficult problem and show that offline algorithms for it can improve upon naive approaches, but that online algorithms essentially cannot.

We present a partitioned algorithm for reducing a symmetric matrix to a tridiagonal form, with partial pivoting. That is, the algorithm computes a factorization PAPT = LTLT , where, P is a permutation matrix, L is lower triangular with a unit diagonal and entries’ magnitudes bounded by 1, and T is symmetric and tridiagonal. The algorithm is based on the basic (nonpartitioned) methods of Parlett and Reid and of Aasen. We show that our factorization algorithm is componentwise backward stable (provided that the growth factor is not too large), with a similar behavior to that of Aasen’s basic algorithm. Our implementation also computes the QR factorization of T and solves linear systems of equations using the computed factorization. The partitioning allows our algorithm to exploit modern computer architectures (in particular, cache memories and high-performance blas libraries). Experimental results demonstrate that our algorithms achieve approximately the same level of performance as the partitioned Bunch-Kaufman factor and solve routines in lapack .

For a 256GB Silicon Blue SSD, a table that maps each 4KB host sector to an arbitrary flash address requires about 256MB of RAM; alas, the SSD only has 64MB of RAM.

NANDFS is a flash file system that exposes a memory-performance tradeoff to system integrators. The file system can be configured to use a large amount of RAM, in which case it delivers excellent performance. In particular, when NANDFS is configured with the same amount of RAM that YAFFS2 uses, the performance of the two file systems is comparable (YAFFS2 is a file system that is widely used in embedded Linux and other embedded environments). But YAFFS2 and other state-of-the-art flash file systems allocate RAM dynamically and do not provide the system builder with a way to limit the amount ofmemory that they allocate. NANDFS, on the other hand, allows the system builder to configure it to use a specific amount of RAM. The performance of NANDFS degrades when the amount of RAM it uses shrinks, but the degradation is graceful, not catastrophic. NANDFS is able to provide this flexibility thanks to a novel data structure that combines a coarsegrained logical-to-physical mapping with a log-structured file system.

Traffic delays and congestion are a major source of inefficiency, wasted fuel, and commuter frustration. Measuring and localizing these delays, and routing users around them, is an important step towards reducing the time people spend stuck in traffic. As others have noted, the proliferation of commodity smartphones that can provide location estimates using a variety of sensors---GPS, WiFi, and/or cellular triangulation---opens up the attractive possibility of using position samples from drivers' phones to monitor traffic delays at a fine spatiotemporal granularity. This paper presents VTrack, a system for travel time estimation using this sensor data that addresses two key challenges: energy consumption and sensor unreliability. While GPS provides highly accurate location estimates, it has several limitations: some phones don't have GPS at all, the GPS sensor doesn't work in "urban canyons" (tall buildings and tunnels) or when the phone is inside a pocket, and the GPS on many phones is power-hungry and drains the battery quickly. In these cases, VTrack can use alternative, less energy-hungry but noisier sensors like WiFi to estimate both a user's trajectory and travel time along the route. VTrack uses a hidden Markov model (HMM)-based map matching scheme and travel time estimation method that interpolates sparse data to identify the most probable road segments driven by the user and to attribute travel times to those segments. We present experimental results from real drive data and WiFi access point sightings gathered from a deployment on several cars. We show that VTrack can tolerate significant noise and outages in these location estimates, and still successfully identify delay-prone segments, and provide accurate enough delays for delay-aware routing algorithms. We also study the best sampling strategies for WiFi and GPS sensors for different energy cost regimes.

We present a new parallel sparse LU factorization algorithm and code. The algorithm uses a column-preordering partial-pivoting unsymmetric-pattern multifrontal approach. Our baseline sequential algorithm is based on UMFPACK 4, but is somewhat simpler and is often somewhat faster than UMFPACK version 4.0. Our parallel algorithm is designed for shared-memory machines with a small or moderate number of processors (we tested it on up to 32 processors). We experimentally compare our algorithm with SuperLU_MT, an existing shared-memory sparse LU factorization with partial pivoting. SuperLU_MT scales better than our new algorithm, but our algorithm is more reliable and is usually faster. More specifically, on matrices that are costly to factor, our algorithm is usually faster on up to 4 processors, and is usually faster on 8 and 16. We were not able to run SuperLU_MT on 32. The main contribution of this article is showing that the column-preordering partial-pivoting unsymmetric-pattern multifrontal approach, developed as a sequential algorithm by Davis in several recent versions of UMFPACK, can be effectively parallelized.

Abstract We present ATSL, an automatically‐tuned sorting library. ATSL generates a sorting routine optimized to the target machine for a specific data type. ATSL finds a high‐performance sorting routine by searching an algorithmic space that we have defined. The search space includes basic sorting algorithms and automatically‐generated compositions of sorting algorithms. Performance measurements are used both for ranking candidate algorithms and for characterizing the behavior of candidates in specific settings (e.g. ranges of input sizes). These characterizations allow ATSL to generate hybrid algorithms that intelligently exploit the strengths of particular algorithms, such as high speed at specific input‐size ranges. Many sorting algorithms can be tuned using numeric parameters and ATSL searches these parameter spaces to find values that yield high performance on the target machine. The building blocks from which ATSL synthesizes sorting algorithms include adaptations of many of the most effective hand‐tuned sorting routines, including several that are tuned for cache efficiency. An extensive experimental evaluation shows that ATSL generates high‐performance codes that are well tuned for the target machine and data type. The experiments were conducted on six different machines, of several architectures, and with three different compilers. The algorithms that are generated are fast; in particular, they beat the hand‐tuned building blocks and the compiler's C++ built‐in sorting routine. The algorithms that ATSL generates on different machines and using different compilers are different from each other. Copyright © 2007 John Wiley & Sons, Ltd.

The reconstruction of a complete watertight model from scan data is still a difficult process. In particular, since scanned data is often incomplete, the reconstruction of the expected shape is an ill-posed problem. Techniques that reconstruct poorly-sampled areas without any user intervention fail in many cases to faithfully reconstruct the topology of the model. The method that we introduce in this paper is topology-aware: it uses minimal user input to make correct decisions at regions where the topology of the model cannot be automatically induced with a reasonable degree of confidence. We first construct a continuous function over a three-dimensional domain. This function is constructed by minimizing a penalty function combining the data points, user constraints, and a regularization term. The optimization problem is formulated in a mesh-independent manner, and mapped onto a specific mesh using the finite-element method. The zero level-set of this function is a first approximation of the reconstructed surface. At complex under-sampled regions, the constraints might be insufficient. Hence, we analyze the local topological stability of the zero level-set to detect weak regions of the surface. These regions are suggested to the user for adding local inside/outside constraints by merely scribbling over a 2D tablet. Each new user constraint modifies the minimization problem, which is solved incrementally. The process is repeated, converging to a topology-stable reconstruction. Reconstructions of models acquired by a structured-light scanner with a small number of scribbles demonstrate the effectiveness of the method.

We present the design and implementation of TINYSTORE, a per-sistent, transactional, garbage-collected memory-management sys-tem, designed to be called from the Java virtual machine of a Java Card. The system is designed for flash-based implementations of Java Card, a variant of the Java platform for smart cards. In the Java Card platform, objects are persistent by default. The platform supports transactions: a sequence of accesses to objects can be ex-plicitly declared to constitute a transaction. TINYSTORE supports explicit transactions and atomically executes individual accesses that are not part of transactions; it also supports garbage collection, even on systems with a small constant amount of RAM. TINYS-TORE uses a novel approach and specialized data structures to effi-ciently manage flash memory. We demonstrate its effectiveness by comparing it to a traditional EEPROM-based memory management system for Java Cards.

We present a new out-of-core sparse symmetric-indefinite factorization algorithm. The most significant innovation of the new algorithm is a dynamic partitioning method for the sparse factor. This partitioning method results in very low I/O traffic and allows the algorithm to run at high computational rates, even though the factor is stored on a slow disk. Our implementation of the new code compares well with both high-performance in-core sparse symmetric-indefinite codes and a high-performance out-of-core sparse Cholesky code.

This article presents an algebraic analysis of a mesh-compression technique called high-pass quantization [Sorkine et al. 2003]. In high-pass quantization, a rectangular matrix based on the mesh topological Laplacian is applied to the vectors of the Cartesian coordinates of a polygonal mesh. The resulting vectors, called δ-coordinates, are then quantized. The applied matrix is a function of the topology of the mesh and the indices of a small set of mesh vertices (anchors) but not of the location of the vertices. An approximation of the geometry can be reconstructed from the quantized δ-coordinates and the spatial locations of the anchors. In this article, we show how to algebraically bound the reconstruction error that this method generates. We show that the small singular value of the transformation matrix can be used to bound both the quantization error and the rounding error which is due to the use of floating-point arithmetic. Furthermore, we prove a bound on this singular value. The bound is a function of the topology of the mesh and of the selected anchors. We also propose a new anchor-selection algorithm, inspired by this bound. We show experimentally that the method is effective and that the computed upper bound on the error is not too pessimistic.

We describe a new out-of-core sparse Cholesky factorization method. The new method uses the elimination tree to partition the matrix, an advanced subtree-scheduling algorithm, and both right-looking and left-looking updates. The implementation of the new method is efficient and robust. On a 2 GHz personal computer with 768 MB of main memory, the code can easily factor matrices with factors of up to 48 GB, usually at rates above 1 Gflop/s. For example, the code can factor audikw, currenly the largest matrix in any matrix collection (factor size over 10 GB), in a little over an hour, and can factor a matrix whose graph is a 140-by-140-by-140 mesh in about 12 hours (factor size around 27 GB).

This article describes the design, implementation, and evaluation of a parallel algorithm for the Cholesky factorization of symmetric banded matrices. The algorithm is part of IBM's parallel engineering and scientific subroutine library version 1.2 and is compatible with ScaLAPACK's banded solver. Analysis, as well as experiments on an IBM SP2 distributed-memory parallel computer, shows that the algorithm efficiently factors banded matrices with wide bandwidth. For example, a 31-mode SP2 factors a large matrix more than 16 times faster than a single node would factor it using the best sequential algorithm, and more than 20 times faster than a single node would using LAPACK's DPBTRF. The algorithm uses novel ideas in the area of distributed dense-matrix computations that include the use of a dynamic schedule for a blocked systolic-like algorithm and the separation of the input and output layouts from the layout the algorithm uses internally. The algorithm alson uses known techniques such as blocking to improve its communication-to-computation ratio and its data-cache behavior.

The communication complexity of a linear transformation is the number of scalars which must be transferred between two processors, each of which holds half the input vector in order to perform the transformation such that each processor holds half the output vector. This letter shows that the communication complexity of the the discrete Fourier transform of size n is at most n/2. That is, given a suitable partition of the input and output vectors, only n/4 complex numbers must be sent in each direction. In contrast, fast Fourier transform (FFT) algorithms transfer at least n/2 complex numbers in each direction.

When a numerical computation fails to fit in the primary memory of a serial or parallel computer, a so-called "out-of-core" algorithm must be used which moves data between primary and secondary memories. In this paper, we study out-of-core algorithms for sparse linear relaxation problems in which each iteration of the algorithm updates the state of every vertex in a graph with a linear combination of the states of its neighbors. We give a general method that can save substantially on the I/O traffic for many problems. For example, our technique allows a computer with M words of primary memory to perform T=/spl Omega/(M/sup 1/5/) cycles of a multigrid algorithm for a two-dimensional elliptic solver over an n-point domain using only /spl Theta/(nT/M/sup 1/5/) I/O transfers, as compared with the naive algorithm which requires /spl Omega/(nT) I/O's.>

Consider a convex set P in R/sup d/ and a piece wise polynomial concave function F: P to R. Let A be an algorithm that given a point x in IR/sup d/ computes F(x) if x in P, or returns a concave polynomial p such that p(x)or= 0. The author assumes that d is fixed and that all comparisons in A depend on the sign of polynomial functions of the input point. He shows that under these conditions, one can find max/sub P/ F in time which is polynomial in the number of arithmetic operations of A. Using this method he gives the first strongly polynomial algorithms for many nonlinear parametric problems in fixed dimension, such as the parametric max flow problem, the parametric minimum s-t distance, the parametric spanning tree problem and other problems. In addition he shows that in one dimension, the same result holds even if one only knows how to approximate the value of F. Specifically, if one can obtain an alpha -approximation for F(x) then one can alpha -approximate the value of maxF. He thus obtains the first polynomial approximation algorithms for many NP-hard problems such as the parametric Euclidean traveling salesman problem.>

In this paper we study fault tolerance in area-universal networks and provide both positive and negative results. We present a new network layout, the mesh of ladders, which is area-universal and fault-tolerant. This network, laid out in an n n area, can simulate any other network layout laid out in the same area with O(log 4 n) slowdown, even if cn/2 log n blocks of size 2 log n 2 log n become faulty (for some constant c < 1). Furthermore, it can simulate any other network layout which has at most f(n) bends in each wire with slowdown O(f(n) log 2 n) even after any number of such blocks become faulty in both networks. Our results are tight, in the sense that if one of our assumptions is removed, it is no longer possible to obtain similar simulation results. We show for example that the width of a layout for a network with at least n nodes and diameter at most n 1-# for any fixed 1/2 < # < 1 must be at least # n). Therefore, if faults can happen in any pattern the remaining non...

Given a convex polygonal object P and an environment consisting of polygonal obstacles, we seek a placement for the largest copy of P that does not intersect any of the obstacles, allowing translation, rotation and scaling. We employ the parametric search technique of Megiddo [Me], and the fixed size polygon placement algorithms developed by Leven and Sharir [LS, LS1], to obtain an algorithm that runs in time O(k 2 n# 4 (kn) log 3 (kn) log log(kn)). We also present several other e#cient algorithms for restricted variants of the extremal polygon containment problem, using the same ideas. These variants include: placement of the largest homothetic copies of one or two convex polygons in another convex polygon and placement of the largest similar copy of a triangle in a convex polygon. 1 Introduction Let P be a convex polygon having k vertices and edges, and let Q be a closed two dimensional space bounded by a collection of polygonal obstacles (the "environment") having altogether n...