July 20, 1996 This distribution consists of the following files A. manpak.f B. daepak.shar C. drivers.f D. manaux.f The packages were written by W. C. Rheinboldt Dept. of Mathematics, University of Pittsburgh Pittsburgh, PA 15260 wcrhein@vms.cis.pitt.edu Brief descriptions follow: A. manpak.f ======== is a package of utility programs for computations with submanifolds of R^NVAR that are implicitly defined by a system of nonlinear equations F(x) = 0. Here F is a mapping from R^NVAR to R^NALG, NALG < NVAR, which is sufficiently smooth on an open set E in R^NVAR and satisfies rank DF(x) = NALG for all x in M = { u in E; F(u) = 0 } Then M is a submanifold of R^NVAR of dimension MDIM = NVAR-NALG. The routines in the package establish several types of local parametrizations (coordinate systems) and -- once available -- compute points on the manifold with given local coordinates. In addition, there are routines for computing the first and second derivatives of these local parametrizations and some other quantities, such as sensitivity measures and the second fundamental tensor. For a discussion of the algorithms see W. C. Rheinboldt, MANPAK: A Set of Algorithms for Computations on Implicitly Defined Manifolds Inst. for Comp. Math. and Appl., Univ. of Pittsburgh, Tech. Report. TR-ICMA-96-198, July 1996 J. Comp. and Math. Applic. submitted Routines in the package ----------------------- AUGM Routine for generating a specific augmented matrix COBAS General-purpose routine for computing an orthonormal basis matrix of the nullspace of a given NALG x NVAR matrix with rank = NALG CURVT Computes approximations of the curvature and principal normal of a path on the manifold passing through three consecutive points. This also provides an approximatation of the diagonal terms of the second fundamental tensor. DGPHI Computes the first derivative of the local parametrization of a general local coordinate system. D2GPHI Computes the second derivative of the local parametrization of a general local coordinate system. DTPHI Differentiating the local parametrization of a tangential local coordinate system GCBAS Computes an orthonormal basis matrix of a local coordinate space consisting of MDIM canonical basis vectors. GLOB Globalizes a vector of local coordinates GNBAS Computes an orthonormal basis matrix of a local coordinate space that contains the NVAR-th canonical basis vector GPHI Computes a point on M with specified local coordinates in a general local coordinate system. MOVFR Uses a moving frame algorithm for ordering a local coordinate basis such that its orientation agrees with that of another basis ORIENT Uses a combinatorial algorithm for ordering a local coordinate basis such that its orientation agrees with that of another basis PROJ Computes the orthogonal projection of a point onto a given local coordinate space SENMAP Computes the sensitivity map at a given point with respect to specified natural parameters SENSNR Computes the Euclidean norm of the sensitivity map at a given point with respect to specified natural parameters TPHI Computes a point on M with specified local coordinates in a tangential local coordinate system TSFT Computes a component of the second fundamental tensor in a tangential local coordinate system. B. daepak.shar =========== is a package of subroutines based on MANPAK and MANAUX for solving the following six types of algebraically explicit differential algebraic equations (DAEs); that is, of DAEs in which either the algebraic equations and/or the algebraic variables are explicitly specified. DAEN1: Solver for Nonlinear, index 1 problems F(u,u',w,t) = 0, u' = du/dt with explicitly given algebraic variable, subject to the consistent initial conditions u(t0) = U0, u'(t0) = UP0, w(t0) = W0, F(U0,UP0,W0,T0) = 0 DAEN2: Solver for nonlinear, index-2 problems G(u,u',w,t) = 0 F(u,t) = 0 subject to the consistent initial conditions u(t0) = u0, u'(t0) = p0, w(t0) = w0, G(u0,p0,t0,w0) = 0, F(u0,t0) = 0 DAEQ2: Solver for Quasi-Linear, index-2 problems A(u,t)u' + B(u,t)w = G(u,t) F(u,t) = 0 subject to the partial initial condition u(t0) = u0, F(u0,t0) = 0 The routine determines the remaining initial conditions up(t0) = up0, w(t0) = w0 such that A(u0,t0)up0 + B(u0,t0)w0 = G(u0,t0) DAEQ3: Solver for quasi-linear, index-3 problems of second order A(u,u',t)u" + B(u,u',t)w = G(u,u',t) F(u,t) = 0 subject to the partial initial condition u(t0) = u0, u'(t0) = up0, such that F(u0,t0) = 0 The routine determines the remaining initial conditions u"(t0) = upp0, w(t0) = w0, such that A(u0,up0,t0)upp0 + B(u0,up0,t0)w0 = G(u0,up0,t0) DAESQ1: Solver for Quasi-linear, autonomous, index-1 problems A(u)u' = G(u), F(u) = 0, u(0) = u0, F(u0) = 0 with optional facility for handling Singular (impasse) points. DAEUL3: Solver for the EUler Lagrange problem of index 3 M(u,t)u" + D_uF(u,t)^T w = G(u,u',t) F(u,t) = 0 subject to the consistent initial conditions u(t0) = u0, u'(t0) = du0, F(u0,t0) = 0 For existence results on the above classes of DAEs and algorithms for the six solvers see W. C. Rheinboldt, Solving Algebraically Explicit DAEs with the MANPAK - Manifold - Algorithms Inst. for Comp. Math. and Appl., Univ. of Pittsburgh, Tech. Reportt. TR-ICMA-96-199, July 1996 J. Comp. and Math. Applic. submitted For an analysis of impasse points and the algorithms used in DAESQ1 see P. Rabier and W. C. Rheinboldt, On Impasse Points of Quasilinear Differential Algebraic Equations J. Math. Anal. and Appl. 181, 1994, 429-454 P. Rabier and W. C. Rheinboldt, On the Computation of Impasse Points of Quasilinear Differential Algebraic Equations Math. of Comput. 62, 1994, 133-154 For the algorithm used in DAEUL3 see P. Rabier and W. C. Rheinboldt, Numerical Solution of Euler-Lagrange Equations for Constrained Mechanical Systems SIAM J. Num. Anal. 32, 1995, 318-329 W. C. Rheinboldt and B. Simeon, Performance Analysis for Solving Euler-Lagrange Equations Appl. Math. Letters 8, 1995, 77-82 SUBROUTINES IN THE PACKAGE DAEN1 Driver for solving the coresponding DAE liested above DRVN1 Rk-step driver for DAEN1 DYN1 Evaluates the local ODE for DAEN1 INTN1 Output interpolation for DAEN1 WSTN1 Writes statistics of a run with DAEN1 DAEN2 Driver for solving the coresponding DAE liested above DRVN2 Rk-step driver for DAEN2 DYN2 Evaluates the local ODE for DAEN2 INTN2 Output interpolation for DAEN2 MATN2 Evaluates the iteration matrix for the iterative solver used in DAEN2 EVXN2 Evaluates the next global point in DAEN2 NWTN2 Chord Newton process for DAEN2 WSTN2 Writes statistics of a run with DAEN2 DAEQ2 Driver for solving the coresponding DAE liested above DRVQ2 Rk-step driver for DAEQ2 DYQ2 Evaluates the local ODE for DAEQ2 INTQ2 Output interpolation for DAEQ2 WSTQ2 Writes statistics of a run with DAEQ2 DAEQ3 Driver for solving the coresponding DAE liested above DRVQ3 Rk-step driver for DAEQ3 DYQ3 Evaluates the local ODE for DAEQ3 INTQ3 Output interpolation for DAEQ3 EVXQ3 Evaluates the next global point in DAEQ3 SLVQ3 Solver for the reduced linear system in DAEQ3 WSTQ3 Writes statistics of a run with DAEQ3 DAESQ1 Driver for solving the coresponding DAE liested above DRVSQ1 Rk-step driver for DAEQ2 DYSQ1 Evaluates the local ODE for DAESQ1 at points not near a singularity DYSYN Evaluates the local ODE for DAESQ1 at points near a singularity SNAUG Evaluates an augmented matrix used in the scaling procedure near singularities. WSTSQ1 Writes statistics of a run with DAESQ1 DAEUL3 Driver for solving the coresponding DAE liested above DRVEUL Rk-step driver for DAEUL3 DYEUL Evaluates the local ODE for DAEUL3 EVXEUL Evaluates the next global point in DAEUL3 WSTEUL Writes statistics of a run with DAEUL3 C. drivers.shar ============ A set of sample drivers for the DAE solvers in daepak.shar For DAEN1: P1N1 'Figure Eight' problem with an exact solution P2N1 Modified orgenator problem P3N1 A batch reactor model For DAEN2: P1N2 A problem by R. Maerz and C. Tischendorf with exact solution P2N2 A problem by U. Ascher and L. Petzold with exact solution P3N2 K. Brenan's Trajectory-Prescribed-Path Control Problem P4N2 K. Brenan's index-two test problem with exact solution For DAEQ2: P1Q2 A problem by U. Ascher and L. Petzold with exact solution P2Q2 K. Brenan's index-two test problem with exact solution P3Q2 A version of the pendulum problem by C. W.Gear, B. Leimkuhler, G.K. Gupta For DAEQ3: P1Q3 A linear test problem with exact solution P2Q3 A nonlinear test problem with exact solution For DAESQ1: P1SQ1 A simple test problem with exact solution P2SQ1 A problem by F. Takens with impasse points For DAEUL3: P1EUL An autonomous pendulum version P2EUL A nonautonomous double-pendulum with rotating pivot-point D. manaux.f ======== is a package of subroutines for uniform support of MANPAK and DAEPAK. The package contains the following subroutines: ODE Solvers DOPSTA Dormand-Prince order 5 RK-step routine for autonomous ODEs using reverse communication for function calls DOPSTN Dormand-Prince order 5 RK-step routine for non-autonomous ODEs using reverse communication for function calls Linear algebra routines BIDIA Reduces a general M by N matrix A, M >= N to an M by N upper bidiagonal matrix DDIST2 Funcion for computing either the Euclidean distance between two vectors X and Y or the Euclidean norm of one such vector X HOUSG Generates a Householder reflector HOUSL Householder reflector times matrix HOUSR Matrix times Householder reflector LQF LQ factorization with column pivoting LQAS For an M x N matrix A with rank A = M <= N and given M-vector y, computes the solution x of A x = y orthogonal to ker A LUF LU factorization with row pivoting LUS1 Solves a linear system for an LU factored matrix A and a given vector X LUSK Solves a linear system for an LU factored matrix A and a given matrix B on the right side QRF QR factorization with column pivoting QORG Generate a matrix Q with orthonormal columns as the product of given Householder reflectors QRS Least squares solver ROTG Generates a plane rotation ROTML Plane rotation times matrix ROTMR Matrix times plane rotation SVD Singular value decomposition of upper bidiagonal matrix SVD2D singular value decomposition of a 2-by-2 triangular matrix Other support routines MSGPRT Uniform message-printing routine DMACH Returns relative machine precision and safe minimum