From c77eb59bc575876fb3dfd475f96e8716cb6b5762 Mon Sep 17 00:00:00 2001 From: Adam Furman Date: Fri, 29 Sep 2023 14:51:24 -0700 Subject: [PATCH] Added paper draft compiler for JOSS --- .github/workflows/draft-pdf.yml | 23 +++++ paper/paper.md | 98 ++++++++++++--------- paper/pysci.bib | 146 ++++++++++++++++++++++++++++++++ 3 files changed, 225 insertions(+), 42 deletions(-) create mode 100644 .github/workflows/draft-pdf.yml create mode 100644 paper/pysci.bib diff --git a/.github/workflows/draft-pdf.yml b/.github/workflows/draft-pdf.yml new file mode 100644 index 0000000..6563593 --- /dev/null +++ b/.github/workflows/draft-pdf.yml @@ -0,0 +1,23 @@ +on: [push] + +jobs: + paper: + runs-on: ubuntu-latest + name: Paper Draft + steps: + - name: Checkout + uses: actions/checkout@v3 + - name: Build draft PDF + uses: openjournals/openjournals-draft-action@master + with: + journal: joss + # This should be the path to the paper within your repo. + paper-path: paper/paper.md + - name: Upload + uses: actions/upload-artifact@v1 + with: + name: paper + # This is the output path where Pandoc will write the compiled + # PDF. Note, this should be the same directory as the input + # paper.md + path: paper/paper.pdf \ No newline at end of file diff --git a/paper/paper.md b/paper/paper.md index 0bfc3c9..640d29b 100644 --- a/paper/paper.md +++ b/paper/paper.md @@ -1,10 +1,20 @@ --- author: -- Adam Furman +- name: Adam Furman + orcid: 0000-0002-6675-1064 + equal-contrib: true + affiliation: 1 +affiliations: + - name: University of Oxford Mathematical Institute, Oxford, United Kingdom + index: 1 bibliography: - pysci.bib -date: September 29, 2023 -title: "[MathMat]{.smallcaps}: A Mathematically Aware Matrix Library" +date: 29 September 2023 +title: "MathMat: A Mathematically Aware Matrix Library" +tags: + - Python + - Linear Algebra + - Matrix --- # Introduction @@ -39,7 +49,7 @@ of any mathematical properties associated with that data. It is left to the programmer to keep track of such properties in a specific application. -The [MathMat]{.smallcaps} library provides a mathematically aware +The MathMat library provides a mathematically aware toolbox for working with numerical linear algebra problems. The core idea is to represent a matrix in two parts: the data, comprising the entries of the matrix, and a collection of known properties about the @@ -49,11 +59,11 @@ addition, the methods for computing properties are stateful, in that they respond to what other properties have already been computed. For example, if a matrix is triangular, the eigenvalues are the entries on the diagonal. Thus, when computing eigenvalues of a matrix known to be -triangular, [MathMat]{.smallcaps} returns the diagonal instead of +triangular, MathMat returns the diagonal instead of performing an additional expensive computation. In addition to mathematical awareness of matrix properties, -[MathMat]{.smallcaps} contains implementations of some matrix +MathMat contains implementations of some matrix algorithms, specifically Krylov subspace methods. The package has visualization capabilities as well, focusing on generating plots that are aesthetically pleasing and clearly visible within academic papers, @@ -61,7 +71,7 @@ with focus again towards those useful in numerical linear algebra. ## The Library Structure -[MathMat]{.smallcaps} is organized into one library with five +MathMat is organized into one library with five submodules, including the main module. Each submodule covers a specific range of functionality, outlined here. @@ -191,13 +201,13 @@ their methods have the prefix `is_*`, those of type `V` are *values*. | `triangular_U` | `C` | `True` if $M$ is upper (right) triangular. | | `unitary` | `C` | `True` if $M$ is unitary, $M^H = M^{-1}$. | -Some [MathMat]{.smallcaps} functions may store additional properties +Some MathMat functions may store additional properties using the `set` method of the `Matrix`, but these are not accessible using methods. ## Efficiency through Computable Properties -The philosophy of [MathMat]{.smallcaps} is to increase efficiency for +The philosophy of MathMat is to increase efficiency for multiple computations on the same matrix. Efficiency of computations is preferred over memory consumption: storing properties will never be more efficient in terms of storage requirements. Under this philosophy, it is @@ -205,7 +215,7 @@ acceptable for an individual initial computation to be more computationally complex, as long as that increased complexity produces knowledge about the matrix that can be reused later. -An illustrative example of the [MathMat]{.smallcaps} philosophy is the +An illustrative example of the MathMat philosophy is the `eigenvalues` method: ```python @@ -247,22 +257,21 @@ SciPy to find eigenvalues, one of which is optimized for Hermitian matrices. Thus, another check is performed, once again yielding valuable information for future tasks. -A call to [MathMat]{.smallcaps}'s `eigenvalues` may take substantially +A call to MathMat's `eigenvalues` may take substantially longer than directly invoking `scipy.linalg.eigvals`, but it is also much more informative. For example, if the matrix was found to be upper triangular then computing the determinant, checking for invertability, and solving the linear equation $Mx = b$ will all be automatically optimized, without the check for triangularity being performed again. -Another consistent advantage of [MathMat]{.smallcaps} is the equivalence +Another consistent advantage of MathMat is the equivalence of sparse and dense matrix representations. While NumPy arrays and SciPy sparse objects are often compatible, many NumPy specific functions -cannot handle a mix of the two. [MathMat]{.smallcaps} abstracts the +cannot handle a mix of the two. MathMat abstracts the distinction away, automatically adapting when a sparse representation is used, and employing the appropriate sparse algorithms to take full advantage of the sparse format. An example is the `is_diagonal` method -(Appendix [8.2](#code:is_diagonal){reference-type="ref" -reference="code:is_diagonal"}), which uses a check based on the +(Appendix \ref{code:is_diagonal}), which uses a check based on the `count_nonzero()` method of sparse objects, or one of two different checks for dense matrices. If the matrix is found to be diagonal, it is automatically stored as a `DiagonalMatrix` in sparse format. @@ -397,7 +406,7 @@ solve the system. The implemented solvers are: triangular $L$ or $U$. The `automatic` method is available in Appendix -[8.3](#code:automatic){reference-type="ref" reference="code:automatic"}. +\ref{code:automatic}. It prefers using an inverse if one exists[^3], and will then check triangularity and the existence of a QR or LU decomposition before falling back on using GMRES. Note that using `automatic` on a `Matrix` @@ -418,7 +427,9 @@ The Generalized Mininimum Residual Method is an example of a Krylov subspace method for solving the (possibly overdetermined) linear system $Ax = b$. It is an iterative algorithm where at iteration $k$, $x_k$ satisfies + $$x_k = \underset{x_k \in \mathcal{K}_k}{\text{argmin}} \lVert \underbrace{Ax_k - b}_{\text{residual } r_k} \rVert_2$$ + where $r_k = Ax_k - b$ is called the residual and is minimized over all $x_k$ within the Krylov subspace $\mathcal{K}_k$ of dimension $k$ [@Saad1986]. As the dimension of the subspace searched over grows, $x_k$ @@ -428,14 +439,16 @@ The effectiveness of GMRES depends on the properties of the Krylov subspace, named after applied mathematician Alexey Krylov [@Krylov1931]. A Krylov subspace is the span of multiplications of increasing powers of a matrix $A$ and a vector $b$: -$$\mathcal{K}_k = \text{span} \left\{ b, Ab, A^2b, \hdots, A^{k-1}b \right\}$$ + +$$\mathcal{K}_k = \text{span} \left\{ b, Ab, A^2b, ..., A^{k-1}b \right\}$$ + which forms a linearly independent set up to some critical dimension $k = k_\text{crit.}$, at which iteration the solution achieved via GMRES is exact [@Saad1986]. Searching for solutions $x_k$ in $\mathcal{K}_k$ yields increasingly accurate approximations for most matrices. While the vectors -$\left\\{ b, Ab, A^2b, \hdots, A^{k_\text{crit.}-1}b \right\\}$ are +$\left\\{ b, Ab, A^2b, ... , A^{k_\text{crit.}-1}b \right\\}$ are mathematically linearly independent, in floating-point arithmetic they quickly lose independence. The limitation can be bypassed by instead constructing an orthonormal basis for $\mathcal{K}_k$. The algorithm for @@ -447,7 +460,9 @@ $q_{k} = A q_{k-1}$ and orthogonalizing it with respect to all previous $q_j$s, recording the orthogonalization dot products $q_j^H q_k$ in an upper Hessenberg matrix $H_k$ [@Saad1986]. The minimization problem is then rewritten as: -$$y_k = \underset{y_k \in \mathbb{R}^N}{\text{argmin}} \Big\Vert H_k y - \lVert b \rVert_2 \mathbf{e}_1 \Big\Vert_2 && \text{and} && x_k = Q_k y_k$$ + +$$y_k = \underset{y_k \in \mathbb{R}^N}{\text{argmin}} \Big\Vert H_k y - \lVert b \rVert_2 \mathbf{e}_1 \Big\Vert_2 \text{ and } x_k = Q_k y_k$$ + where $\mathbf{e}_1$ is the first unit vector of size $N$. The upper Hessenberg structure yields a much more efficiently solvable least-squares problem for $y_k$ [@Saad1986]. @@ -471,8 +486,7 @@ Gram-Schmidt (MGS) still fails, repeated MGS, achieved by re-orthogonalizing the resulting vector, is employed [@Daniel1976]. Both MGS and re-orthogonalizing MGS are implemented, as `mathmat.solve.Krylov.arnoldi_mgs` (available in Appendix -[8.4](#code:arnoldi_mgs){reference-type="ref" -reference="code:arnoldi_mgs"}) and `arnoldi_mgs_2x` respectively. +\ref{code:arnoldi_mgs}) and `arnoldi_mgs_2x` respectively. A novel GMRES-based algorithm, developed by Nakatsukasa and Tropp [@Nakatsukasa2022], is Sketched GMRES. sGMRES differs from GMRES in two @@ -489,7 +503,7 @@ GMRES[^4], and is currently being further studied. # Visualizations -The final major submodule of [MathMat]{.smallcaps} is `mathmat.vis`, the +The final major submodule of MathMat is `mathmat.vis`, the visualizations toolkit. The goal of this submodule is to provide aesthetically pleasing plots that fit well within mathematical publications with minimum adjustment required. `matplotlib` @@ -503,7 +517,7 @@ and math is changed to Computer Modern, the default LaTeX typeface. Finally, the thickness of lines is slightly increased to yield better visibility when embedded into papers. All figures included in this report have been exported as `.svg` files using the -[MathMat]{.smallcaps} configuration. +MathMat configuration. Each plot is an instance of the `Plot` class. A generic `Plot` object displays no data, but is capable of setting the title, $x$ and $y$ axis @@ -555,24 +569,22 @@ the visualization submodule without displaying any of the plots. To see the plots, `demo_visualizations.py` can be run to see plots representative of the library's functionality, mirroring tasks that might be done in the course of normal use. A selection of these plots is -reproduced in Appendix [8.5](#visualizations){reference-type="ref" -reference="visualizations"}. +reproduced in Appendix \ref{visualizations}. -To run the demonstrations and tests, first install [MathMat]{.smallcaps} +To run the demonstrations and tests, first install MathMat using the instructions in Appendix -[8.1](#install_instructions){reference-type="ref" -reference="install_instructions"}. Make sure to additionally install +\ref{install_instructions}. Make sure to additionally install `pytest`. Then, ensuring that the `tests` directory and `demo_visualizations.py` file exists, run `pytest`. Note that it is expected that the visualization tests will raise warnings, as the figures are not shown. -Interactive documentation for [MathMat]{.smallcaps} is also available in +Interactive documentation for MathMat is also available in the `docs` folder. Simply open the `index.html` file in any web browser. # Summary and Conclusion -The [MathMat]{.smallcaps} package is a purpose-specific Python library +The MathMat package is a purpose-specific Python library for numerical linear algebra. It aims to provide a mathematically-aware suite of tools for working with matrices, utilizing the results of past computations on the same `Matrix` object to identify which algorithms @@ -590,9 +602,7 @@ work, several key methods are reproduced in the Appendix. Effort has been taken to render the code legible and thoroughly commented, and to adhere to standard Python guidelines. -After following the installation instructions (Appendix -[8.1](#install_instructions){reference-type="ref" -reference="install_instructions"}), try the following to see the library +After following the installation instructions (Appendix \ref{install_instructions}), try the following to see the library at work: in a new Python shell, run the commands ```python @@ -615,9 +625,10 @@ width="0.9\\columnwidth"} # Appendix -## Installation Instructions {#install_instructions} +## Installation Instructions +\label{install_instructions} -[MathMat]{.smallcaps} can be installed in several ways. In all cases, it +MathMat can be installed in several ways. In all cases, it is recommended to create a new virtual environment using the command `py -m venv venv` and then activating the virtual environment using the appropriate command for the operating system. @@ -628,13 +639,14 @@ can be installed by running code of the library can be unpacked from `mathmat-1.0.0.tar.gz`, giving easy access to the tests and demo visualizations. -[MathMat]{.smallcaps} depends on `numpy`, `scipy`, `matplotlib`, and +MathMat depends on `numpy`, `scipy`, `matplotlib`, and `lazy-import`, which must be installed for the library to work. `pip` will do so automatically if installing from the `.whl` file. The full development environment can be recreated by installing the packages listed in `requirements_dev.txt`. -## The `is_diagonal` Method {#code:is_diagonal} +## The `is_diagonal` Method +\label{code:is_diagonal} ```python @_computable_property("diagonal") @@ -664,7 +676,8 @@ listed in `requirements_dev.txt`. return check ``` -## The `automatic` Method {#code:automatic} +## The `automatic` Method +\label{code:automatic} ```python def automatic(A: Matrix, b: Vector, get_method=False): @@ -701,7 +714,8 @@ def automatic(A: Matrix, b: Vector, get_method=False): return A.lin_solver(b) ``` -## The `arnoldi_mgs` Method {#code:arnoldi_mgs} +## The `arnoldi_mgs` Method +\label{code:arnoldi_mgs} ```python def arnoldi_mgs(A: Matrix, b: Vector, k: int, tolerance=1e-4, **_): @@ -751,7 +765,8 @@ def automatic(A: Matrix, b: Vector, get_method=False): return K_k, H ``` -## Selected Visualizations {#visualizations} +## Selected Visualizations +\label{visualizations} **Eigenvalues of the `west132` Matrix** @@ -786,5 +801,4 @@ convention. which $M^{-1}$ can be computed in a numerically stable way. [^4]: Such as when applied to the `FS 680` matrix used as an example in -Appendix [8.5](#visualizations){reference-type="ref" -reference="visualizations"}. +Appendix \ref{visualizations}. diff --git a/paper/pysci.bib b/paper/pysci.bib new file mode 100644 index 0000000..f07f19e --- /dev/null +++ b/paper/pysci.bib @@ -0,0 +1,146 @@ +@Electronic{Carbonnelle2023, + author = {Pierre Carbonnelle}, + month = jun, + note = {https://pypl.github.io/PYPL.html}, + title = {{PYPL PopularitY of Programming Language Index}}, + url = {https://pypl.github.io/PYPL.html}, + year = {2023}, + creationdate = {2023-06-27T16:55:00}, +} + +@Electronic{Wang2022, + author = {Peter Wang}, + month = feb, + note = {https://www.anaconda.com/blog/why-python}, + organization = {Anaconda}, + title = {{Why Python: The Factors Leading to the Language’s Ascendancy}}, + url = {https://www.anaconda.com/blog/why-python}, + year = {2022}, + creationdate = {2023-06-27T17:06:35}, +} + +@Electronic{NumPy, + note = {https://numpy.org/}, + title = {{NumPy}}, + url = {https://numpy.org/}, + creationdate = {2023-06-27T17:08:53}, +} + +@Electronic{SciPy, + note = {https://scipy.org/}, + title = {{SciPy}}, + url = {https://scipy.org/}, + creationdate = {2023-06-27T17:18:54}, +} + + + + + + +@Article{Saad1986, + author = {Saad, Youcef and Schultz, Martin H.}, + journal = {SIAM Journal on Scientific and Statistical Computing}, + title = {{GMRES: A Generalized Minimal Residual Algorithm for Solving Nonsymmetric Linear Systems}}, + year = {1986}, + number = {3}, + pages = {856-869}, + volume = {7}, + abstract = {We present an iterative method for solving linear systems, which has the property of minimizing at every step the norm of the residual vector over a Krylov subspace. The algorithm is derived from the Arnoldi process for constructing an \$l\_2 \$-orthogonal basis of Krylov subspaces. It can be considered as a generalization of Paige and Saunders’ MINRES algorithm and is theoretically equivalent to the Generalized Conjugate Residual (GCR) method and to ORTHODIR. The new algorithm presents several advantages over GCR and ORTHODIR.}, + creationdate = {2023-06-28T17:41:01}, + doi = {10.1137/0907058}, + eprint = {https://doi.org/10.1137/0907058}, + url = {https://doi.org/10.1137/0907058}, +} + +@Article{Nakatsukasa2022, + author = {Yuji Nakatsukasa and Joel A. Tropp}, + title = {{Fast and Accurate Randomized Algorithms for Linear Systems and Eigenvalue Problems}}, + year = {2022}, + archiveprefix = {arXiv}, + creationdate = {2023-06-28T17:45:40}, + eprint = {2111.00113}, + primaryclass = {math.NA}, +} + +@Article{Krylov1931, + author = {Aleksey Krylov}, + journal = {Bulletin de l'Acad\'emie des Sciences de l'URSS. Classe des sciences math\'ematiques et na}, + title = {De la r\'esolution num\'erique de l'\'equation servant \`a d\'eterminer dans des questions de m\'ecanique appliqu\'ee les fr\'equences de petites oscillations des syst\`emes mat\'eriels}, + year = {1931}, + number = {4}, + pages = {491--539}, + creationdate = {2023-06-29T12:03:21}, + date = {1931}, +} + + +@Article{Arnoldi1951, + author = {Arnoldi, W. E.}, + journal = {Quarterly of Applied Mathematics}, + title = {{The Principle of Minimized Iterations in the Solution of the Matrix Eigenvalue Problem}}, + year = {1951}, + issn = {0033-569X, 1552-4485}, + number = {1}, + pages = {17--29}, + volume = {9}, + abstract = {Advancing research. Creating connections.}, + creationdate = {2023-06-29T14:30:47}, + doi = {10.1090/qam/42792}, + file = {Full Text PDF:https\://www.ams.org/qam/1951-09-01/S0033-569X-1951-42792-9/S0033-569X-1951-42792-9.pdf:application/pdf}, + language = {en}, + url = {https://www.ams.org/qam/1951-09-01/S0033-569X-1951-42792-9/}, + urldate = {2023-06-29}, +} + +@Article{Matinfar2012, + author = {M. Matinfar and H. Zareamoghaddam and M. Eslami and M. Saeidy}, + journal = {Computers and Mathematics with Applications}, + title = {{GMRES Implementations and Residual Smoothing Techniques for Solving Ill-Posed Linear Systems}}, + year = {2012}, + issn = {0898-1221}, + number = {1}, + pages = {1-13}, + volume = {63}, + abstract = {There are verities of useful Krylov subspace methods to solve nonsymmetric linear system of equations. GMRES is one of the best Krylov solvers with several different variants to solve large sparse linear systems. Any GMRES implementation has some advantages. As the solution of ill-posed problems are important. In this paper, some GMRES variants are discussed and applied to solve these kinds of problems. Residual smoothing techniques are efficient ways to accelerate the convergence speed of some iterative methods like CG variants. At the end of this paper, some residual smoothing techniques are applied for different GMRES methods to test the influence of these techniques on GMRES implementations.}, + creationdate = {2023-06-29T15:46:59}, + doi = {https://doi.org/10.1016/j.camwa.2011.09.022}, + keywords = {Least square, GMRES, Givens rotations, Householder}, + url = {https://www.sciencedirect.com/science/article/pii/S0898122111007905}, +} + +@Article{Daniel1976, + author = {James W. Daniel and William B. Gragg and Linda Kaufman and G. W. Stewart}, + journal = {Mathematics of Computation}, + title = {{Reorthogonalization and Stable Algorithms for Updating the Gram-Schmidt QR Factorization}}, + year = {1976}, + pages = {772-795}, + volume = {30}, + creationdate = {2023-06-29T15:53:11}, +} + +@Electronic{Matplotlib, + note = {https://matplotlib.org/}, + title = {MatPlotLib}, + url = {https://matplotlib.org/}, + creationdate = {2023-06-29T16:56:13}, +} + +@Electronic{PyTest, + note = {https://docs.pytest.org/en/7.3.x/}, + title = {{PyTest}}, + url = {https://docs.pytest.org/en/7.3.x/}, + creationdate = {2023-06-30T14:29:34}, +} + +@Article{Nakatsukasa2020, + author = {Yuji Nakatsukasa}, + title = {{Fast and Stable Randomized Low-Rank Matrix Approximation}}, + year = {2020}, + archiveprefix = {arXiv}, + creationdate = {2023-06-30T15:33:54}, + eprint = {2009.11392}, + primaryclass = {math.NA}, +} + +@Comment{jabref-meta: databaseType:bibtex;}