DSJM: A Software Toolkit for Direct Determination of Sparse Jacobian Matrices ⋆ Mahmudul Hasan 1 , Shahadat Hossain 1 , Trond Steihaug 2 1 University of Lethbridge, Department of Mathematics and Computer Science Lethbridge, AB T1K 3M4, Canada shahadat.hossain@uleth.ca, hasan@uleth.ca 2 University of Bergen, Department of Informatics Bergen, N-5020 Norway Trond.Steihaug@ii.uib.no Abstract. We describe DSJM, a software toolkit written in portable C++ that enables the direct determination of sparse Jacobian matrices whose sparsity pattern is a priori known. Using the seed matrix S ∈ℜ n×p defined by the partitioning information the nonzero un- knowns of A ∈ℜ m×n can be obtained in AS = B where B ∈ℜ m×p has been obtained via finite difference approximation or forward automatic differentiation. Our preliminary imple- mentation includes well-known as well as new column ordering heuristics. Early numerical testing is highly promising both in terms of running time and the number of matrix-vector products needed to determine A. Keywords. Structural Orthogonality, Direct Determination, Software Implementation. 1 Efficient Software Implementation for Direct Determination of Sparse Jacobian Matrices The determination of a sparse Jacobian matrix with a priori known sparsity pattern of at least once continuously differentiable mapping F : ℜ n →ℜ m can be viewed as a computation of p matrix-vector products AS ≡ B where A is an approximation of the Jacobian matrix F ′ (x) but has the same sparsity pattern and S ∈ℜ n×p is a seed matrix of p directions. The nonzero entries of A are obtained using finite difference approximation ∂F (x + ts) ∂t     t=0 = F ′ (x)s ≈ As = 1 ε [F (x + εs) − F (x)] ≡ b, (1) or via the forward mode of automatic differentiation at a computational effort equal to a small mul- tiple of the computational effort for evaluating F at x. If the ρ i rows of S defined by the ρ i nonzero unknowns of A(i, :) are linearly independent then A(i, :) is uniquely determined from B(i, :). We have direct determination if the nonzero unknowns can be recovered from B without any extra arithmetic operations. The Curtis, Powell, and Reid method, henceforth the CPR [2] method, ex- ploits A’s sparsity to define S(:,k)= ∑ j∈C k e j ,k =1,...,p where C = {C 1 ,..., C p } is a partition of column indices such that for each pair of indices j, l ∈C k , the product a ij a il ,i =1,...,m is identically 0, thereby yielding direct determination of the unknowns. The DSJM currently im- plements the well-known column ordering algorithms “largest first order (LFO)”, “smallest last order (SLO)”, and “incidence degree order (IDO)”, as well as a new column partitioning algorithm M(atrix)RLF based on “recursive largest first (RLF)” heuristic [8] for direct determination of sparse Jacobian matrices. Inspired by the highly acclaimed software DSM [1] our implementation uses general purpose data structures compressed row storage (CRS) and compressed column storage (CCS) for input matrices emphasizing efficiency and simplicity of use. Results from preliminary numerical experiments are given in Tables 1 and 2. Table 1 displays timing results for selected large problems [4](m rows, n columns, nnz nonzero entries) in seconds obtained on an IBM PC with 2.8 GHz Intel Pentium CPU, 1 GB RAM, and 512 KB L2 cache running Linux. The timing shown is that of IDO column ordering (ot) and greedy partitioning (pt) [5]. For both ColPack and DSJM the column ordering step is computationally more expensive than the greedy partitioning for almost all the test problems. This is expected since ⋆ This research was supported in part by the Natural Sciences and Engineering Research Council of Canada (NSERC) and the Research Council of Norway (NFR).