Implementation of ConjugateGradients method using C and Nvidia CUDA

p-sto, updated 🕥 2022-06-21 21:10:35

ConjugateGradients

Implementation of Conjugate Gradient method for solving systems of linear equation using Python, C and Nvidia CUDA. Currently only Python implementation is available - it includes Conjugate Gradient Method and Preconditioned Conjugate Gradient with Jacobi pre-conditioner (hopefully others will be added as well).

Road-map:

::

  • [X] Python implementation - [X] Create test matrices generator - [X] Implement pure CG - [ ] Implement PCG - [X] Jacobi preconditioner - [ ] SSOR preconditioner - [ ] Incomplete Cholesky factorization preconditioner

  • [ ] C implementation - [ ] Implement pure CG - [X] Implementation using BLAS library (Intel MKL) - [ ] Custom BLAS implementation using OpenMP - [ ] Implement PCG - [ ] Jacobi preconditioner - [X] Implementation using BLAS library (Intel MKL) - [ ] Custom BLAS implementation using OpenMP

  • [ ] CUDA implementation - [ ] Implement pure CG - [ ] Reference implementation using CUBLAS - [ ] Custom kernels implementation - [ ] Implement PCG - [ ] Jacobi preconditioner - [ ] Reference implementation using CUBLAS - [ ] Custom kernels implementation

  • MKL can be obtained free of charge: https://software.intel.com/en-us/mkl

Getting Started

::

$ git clone https://github.com/stovorov/ConjugateGradients
$ cd ConjugateGradients

Python implementation

Prepare environment ~~~~~~~~~~~~~~~~~~~

::

$ source run_me.sh (sets PYTHONPATH)
$ cd scripts
$ make venv
$ source venv/bin/activate
$ make test

Usage ~~~~~

.. code:: python

from random import uniform
from scripts.ConjugateGradients.test_matrices import TestMatrices
from scripts.ConjugateGradients.utils import get_solver

import numpy as np

matrix_size = 100
# patterns are: quadratic, rectangular, arrow, noise, curve
# pattern='qrana' means that testing matrix will be composition of all mentioned patterns
a_matrix = TestMatrices.get_random_test_matrix(matrix_size)
x_vec = np.vstack([1 for x in range(matrix_size)])
b_vec = np.vstack([uniform(0, 1) for x in range(matrix_size)])
CGSolver = get_solver('CG')             # pylint: disable=invalid-name; get_solver returns Class
PCGJacobiSolver = get_solver('PCG')     # pylint: disable=invalid-name; get_solver returns Class
cg_solver = CGSolver(a_matrix, b_vec, x_vec)
cg_solver.solve()
cg_solver.show_convergence_profile()

pcg_solver = PCGJacobiSolver(a_matrix, b_vec, x_vec)
pcg_solver.solve()

CGSolver.compare_convergence_profiles(cg_solver, pcg_solver)

You can view convergence profile using solver's show_convergence_profile method:

.. image:: doc/cg_conv_visual.png
    :height: 200 px
    :width: 200 px
    :scale: 50 %

You can compare convergence profiles of difference solvers using compare_convergence_profiles method:

.. image:: doc/comparison.png
    :height: 200 px
    :width: 200 px
    :scale: 50 %

Different testing matrices can be generated using TestMatrix class, for more information please refer methods descriptions. Matrices can be viewed using view_matrix function, which can be found in utils.py. Below matrices are symmetric and positively defined.

.. image:: doc/arn_matrix.png
    :height: 200 px
    :width: 200 px
    :scale: 50 %

.. image:: doc/crn_matrix.png
    :height: 200 px
    :width: 200 px
    :scale: 50 %

.. image:: doc/rnqa_matrix.png
    :height: 200 px
    :width: 200 px
    :scale: 50 %

Examples can be found in scripts/ConjugateGradients/demo.py Required Python 3.5+

CPU/GPU implementation

Libraries and compilation ~~~~~~~~~~~~~~~~~~~~~~~~~

Before compiling code, make sure you have installed:

::

1. Intel MKL library
2. Nvidia CUDA with NVCC compiler

Intel MKL library is used for BLAS operations. Implementation was tested on version 2017 though older should work as well. By default MKL will be installed in directory /opt/intel/mkl/. Before compiling make sure prepare_env.sh has proper paths to MKL and CUDA libraries.

In Makefile set accordingly:

::

1. MKLROOT
2. NVCC
3. CUDALIBPATH

By default MKL will be compiled as a static library. CUDA is linked dynamically. LDFLAGS are used to set dependencies for MKL, please refer to MKL link line advisor to be sure to have it set properly:

https://software.intel.com/en-us/articles/intel-mkl-link-line-advisor

CUDAFLAGS are responsible for setting CUDA libraries.

GCC is used for compiling .c files, NVCC is used for .cu files. Whole project is linked by GCC.

To compile:

::

$ source prepare_env.sh
$ make

Use make clean command to delete compiled build.

Running ConjugateGradient ~~~~~~~~~~~~~~~~~~~~~~~~~

Running single core CPU MKL implementation:

./ConjugateGradient -i input_matrix.txt

Running multiple core CPU MKL implementation:

./ConjugateGradient -i input_matrix.txt -mt 4

Running GPU implementation (single device only available):

./ConjugateGradient -i input_matrix.txt --gpu

::

If there are no CUDA devices, CPU implementation will be launched.

input_matrix.txt is expected to be CSR formatted matrix, various examples can be generated by Python scripts.

Conjugate Gradients description

A bit about Conjugate Gradients and when it actually works (collection of information found over internet):

CG will work when is applied on symmetrical and positively defined matrix.

CG is equivalent to applying the Lanczos algorithm on the given matrix with the starting vector given by the (normalized) residual of the initial approximation. source: https://math.stackexchange.com/questions/882713/application-of-conjugate-gradient-method-to-non-symmetric-matrices

Resources: ~~~~~~~~~~

General overview and derivation is described on Wiki: https://en.wikipedia.org/wiki/Derivation_of_the_conjugate_gradient_method

Though this description has a lot of shortcuts and will probably leave you with a more questions then before reading it...

A good description can be found in Painless Conjugate Gradient:

https://www.cs.cmu.edu/~quake-papers/painless-conjugate-gradient.pdf A bit complex work but worth reading (but requires a lot of focus...at least from me...).

A lot about preconditioning could be found here: http://netlib.org/linalg/html_templates/node51.html haven't read everything but may explain a lot (still, will probably leave you with a lot of questions...).

Issues

Bump numpy from 1.13.0 to 1.22.0 in /scripts

opened on 2022-06-21 21:10:30 by dependabot[bot]

Bumps numpy from 1.13.0 to 1.22.0.

Release notes

Sourced from numpy's releases.

v1.22.0

NumPy 1.22.0 Release Notes

NumPy 1.22.0 is a big release featuring the work of 153 contributors spread over 609 pull requests. There have been many improvements, highlights are:

  • Annotations of the main namespace are essentially complete. Upstream is a moving target, so there will likely be further improvements, but the major work is done. This is probably the most user visible enhancement in this release.
  • A preliminary version of the proposed Array-API is provided. This is a step in creating a standard collection of functions that can be used across application such as CuPy and JAX.
  • NumPy now has a DLPack backend. DLPack provides a common interchange format for array (tensor) data.
  • New methods for quantile, percentile, and related functions. The new methods provide a complete set of the methods commonly found in the literature.
  • A new configurable allocator for use by downstream projects.

These are in addition to the ongoing work to provide SIMD support for commonly used functions, improvements to F2PY, and better documentation.

The Python versions supported in this release are 3.8-3.10, Python 3.7 has been dropped. Note that 32 bit wheels are only provided for Python 3.8 and 3.9 on Windows, all other wheels are 64 bits on account of Ubuntu, Fedora, and other Linux distributions dropping 32 bit support. All 64 bit wheels are also linked with 64 bit integer OpenBLAS, which should fix the occasional problems encountered by folks using truly huge arrays.

Expired deprecations

Deprecated numeric style dtype strings have been removed

Using the strings "Bytes0", "Datetime64", "Str0", "Uint32", and "Uint64" as a dtype will now raise a TypeError.

(gh-19539)

Expired deprecations for loads, ndfromtxt, and mafromtxt in npyio

numpy.loads was deprecated in v1.15, with the recommendation that users use pickle.loads instead. ndfromtxt and mafromtxt were both deprecated in v1.17 - users should use numpy.genfromtxt instead with the appropriate value for the usemask parameter.

(gh-19615)

... (truncated)

Commits


Dependabot compatibility score

Dependabot will resolve any conflicts with this PR as long as you don't alter it yourself. You can also trigger a rebase manually by commenting @dependabot rebase.

Dependabot commands and options
You can trigger Dependabot actions by commenting on this PR: - `@dependabot rebase` will rebase this PR - `@dependabot recreate` will recreate this PR, overwriting any edits that have been made to it - `@dependabot merge` will merge this PR after your CI passes on it - `@dependabot squash and merge` will squash and merge this PR after your CI passes on it - `@dependabot cancel merge` will cancel a previously requested merge and block automerging - `@dependabot reopen` will reopen this PR if it is closed - `@dependabot close` will close this PR and stop Dependabot recreating it. You can achieve the same result by closing it manually - `@dependabot ignore this major version` will close this PR and stop Dependabot creating any more for this major version (unless you reopen the PR or upgrade to it yourself) - `@dependabot ignore this minor version` will close this PR and stop Dependabot creating any more for this minor version (unless you reopen the PR or upgrade to it yourself) - `@dependabot ignore this dependency` will close this PR and stop Dependabot creating any more for this dependency (unless you reopen the PR or upgrade to it yourself) - `@dependabot use these labels` will set the current labels as the default for future PRs for this repo and language - `@dependabot use these reviewers` will set the current reviewers as the default for future PRs for this repo and language - `@dependabot use these assignees` will set the current assignees as the default for future PRs for this repo and language - `@dependabot use this milestone` will set the current milestone as the default for future PRs for this repo and language You can disable automated security fix PRs for this repo from the [Security Alerts page](https://github.com/p-sto/ConjugateGradients/network/alerts).

Bump ipython from 6.1.0 to 7.16.3 in /scripts

opened on 2022-01-21 18:57:19 by dependabot[bot]

Bumps ipython from 6.1.0 to 7.16.3.

Commits


Dependabot compatibility score

Dependabot will resolve any conflicts with this PR as long as you don't alter it yourself. You can also trigger a rebase manually by commenting @dependabot rebase.

Dependabot commands and options
You can trigger Dependabot actions by commenting on this PR: - `@dependabot rebase` will rebase this PR - `@dependabot recreate` will recreate this PR, overwriting any edits that have been made to it - `@dependabot merge` will merge this PR after your CI passes on it - `@dependabot squash and merge` will squash and merge this PR after your CI passes on it - `@dependabot cancel merge` will cancel a previously requested merge and block automerging - `@dependabot reopen` will reopen this PR if it is closed - `@dependabot close` will close this PR and stop Dependabot recreating it. You can achieve the same result by closing it manually - `@dependabot ignore this major version` will close this PR and stop Dependabot creating any more for this major version (unless you reopen the PR or upgrade to it yourself) - `@dependabot ignore this minor version` will close this PR and stop Dependabot creating any more for this minor version (unless you reopen the PR or upgrade to it yourself) - `@dependabot ignore this dependency` will close this PR and stop Dependabot creating any more for this dependency (unless you reopen the PR or upgrade to it yourself) - `@dependabot use these labels` will set the current labels as the default for future PRs for this repo and language - `@dependabot use these reviewers` will set the current reviewers as the default for future PRs for this repo and language - `@dependabot use these assignees` will set the current assignees as the default for future PRs for this repo and language - `@dependabot use this milestone` will set the current milestone as the default for future PRs for this repo and language You can disable automated security fix PRs for this repo from the [Security Alerts page](https://github.com/stovorov/ConjugateGradients/network/alerts).

Bump py from 1.4.34 to 1.10.0 in /scripts

opened on 2021-04-20 17:07:04 by dependabot[bot]

Bumps py from 1.4.34 to 1.10.0.

Changelog

Sourced from py's changelog.

1.10.0 (2020-12-12)

  • Fix a regular expression DoS vulnerability in the py.path.svnwc SVN blame functionality (CVE-2020-29651)
  • Update vendored apipkg: 1.4 => 1.5
  • Update vendored iniconfig: 1.0.0 => 1.1.1

1.9.0 (2020-06-24)

  • Add type annotation stubs for the following modules:

    • py.error
    • py.iniconfig
    • py.path (not including SVN paths)
    • py.io
    • py.xml

    There are no plans to type other modules at this time.

    The type annotations are provided in external .pyi files, not inline in the code, and may therefore contain small errors or omissions. If you use py in conjunction with a type checker, and encounter any type errors you believe should be accepted, please report it in an issue.

1.8.2 (2020-06-15)

  • On Windows, py.path.locals which differ only in case now have the same Python hash value. Previously, such paths were considered equal but had different hashes, which is not allowed and breaks the assumptions made by dicts, sets and other users of hashes.

1.8.1 (2019-12-27)

  • Handle FileNotFoundError when trying to import pathlib in path.common on Python 3.4 (#207).

  • py.path.local.samefile now works correctly in Python 3 on Windows when dealing with symlinks.

1.8.0 (2019-02-21)

  • add "importlib" pyimport mode for python3.5+, allowing unimportable test suites to contain identically named modules.

  • fix LocalPath.as_cwd() not calling os.chdir() with None, when being invoked from a non-existing directory.

... (truncated)

Commits
  • e5ff378 Update CHANGELOG for 1.10.0
  • 94cf44f Update vendored libs
  • 5e8ded5 testing: comment out an assert which fails on Python 3.9 for now
  • afdffcc Rename HOWTORELEASE.rst to RELEASING.rst
  • 2de53a6 Merge pull request #266 from nicoddemus/gh-actions
  • fa1b32e Merge pull request #264 from hugovk/patch-2
  • 887d6b8 Skip test_samefile_symlink on pypy3 on Windows
  • e94e670 Fix test_comments() in test_source
  • fef9a32 Adapt test
  • 4a694b0 Add GitHub Actions badge to README
  • Additional commits viewable in compare view


Dependabot compatibility score

Dependabot will resolve any conflicts with this PR as long as you don't alter it yourself. You can also trigger a rebase manually by commenting @dependabot rebase.

Dependabot commands and options
You can trigger Dependabot actions by commenting on this PR: - `@dependabot rebase` will rebase this PR - `@dependabot recreate` will recreate this PR, overwriting any edits that have been made to it - `@dependabot merge` will merge this PR after your CI passes on it - `@dependabot squash and merge` will squash and merge this PR after your CI passes on it - `@dependabot cancel merge` will cancel a previously requested merge and block automerging - `@dependabot reopen` will reopen this PR if it is closed - `@dependabot close` will close this PR and stop Dependabot recreating it. You can achieve the same result by closing it manually - `@dependabot ignore this major version` will close this PR and stop Dependabot creating any more for this major version (unless you reopen the PR or upgrade to it yourself) - `@dependabot ignore this minor version` will close this PR and stop Dependabot creating any more for this minor version (unless you reopen the PR or upgrade to it yourself) - `@dependabot ignore this dependency` will close this PR and stop Dependabot creating any more for this dependency (unless you reopen the PR or upgrade to it yourself) - `@dependabot use these labels` will set the current labels as the default for future PRs for this repo and language - `@dependabot use these reviewers` will set the current reviewers as the default for future PRs for this repo and language - `@dependabot use these assignees` will set the current assignees as the default for future PRs for this repo and language - `@dependabot use this milestone` will set the current milestone as the default for future PRs for this repo and language You can disable automated security fix PRs for this repo from the [Security Alerts page](https://github.com/stovorov/ConjugateGradients/network/alerts).

Bump pygments from 2.2.0 to 2.7.4 in /scripts

opened on 2021-03-29 17:01:53 by dependabot[bot]

Bumps pygments from 2.2.0 to 2.7.4.

Release notes

Sourced from pygments's releases.

2.7.4

  • Updated lexers:

    • Apache configurations: Improve handling of malformed tags (#1656)

    • CSS: Add support for variables (#1633, #1666)

    • Crystal (#1650, #1670)

    • Coq (#1648)

    • Fortran: Add missing keywords (#1635, #1665)

    • Ini (#1624)

    • JavaScript and variants (#1647 -- missing regex flags, #1651)

    • Markdown (#1623, #1617)

    • Shell

      • Lex trailing whitespace as part of the prompt (#1645)
      • Add missing in keyword (#1652)

    • SQL - Fix keywords (#1668)

    • Typescript: Fix incorrect punctuation handling (#1510, #1511)

  • Fix infinite loop in SML lexer (#1625)

  • Fix backtracking string regexes in JavaScript/TypeScript, Modula2 and many other lexers (#1637)

  • Limit recursion with nesting Ruby heredocs (#1638)

  • Fix a few inefficient regexes for guessing lexers

  • Fix the raw token lexer handling of Unicode (#1616)

  • Revert a private API change in the HTML formatter (#1655) -- please note that private APIs remain subject to change!

  • Fix several exponential/cubic-complexity regexes found by Ben Caller/Doyensec (#1675)

  • Fix incorrect MATLAB example (#1582)

Thanks to Google's OSS-Fuzz project for finding many of these bugs.

2.7.3

... (truncated)

Changelog

Sourced from pygments's changelog.

Version 2.7.4

(released January 12, 2021)

  • Updated lexers:

    • Apache configurations: Improve handling of malformed tags (#1656)

    • CSS: Add support for variables (#1633, #1666)

    • Crystal (#1650, #1670)

    • Coq (#1648)

    • Fortran: Add missing keywords (#1635, #1665)

    • Ini (#1624)

    • JavaScript and variants (#1647 -- missing regex flags, #1651)

    • Markdown (#1623, #1617)

    • Shell

      • Lex trailing whitespace as part of the prompt (#1645)
      • Add missing in keyword (#1652)

    • SQL - Fix keywords (#1668)

    • Typescript: Fix incorrect punctuation handling (#1510, #1511)

  • Fix infinite loop in SML lexer (#1625)

  • Fix backtracking string regexes in JavaScript/TypeScript, Modula2 and many other lexers (#1637)

  • Limit recursion with nesting Ruby heredocs (#1638)

  • Fix a few inefficient regexes for guessing lexers

  • Fix the raw token lexer handling of Unicode (#1616)

  • Revert a private API change in the HTML formatter (#1655) -- please note that private APIs remain subject to change!

  • Fix several exponential/cubic-complexity regexes found by Ben Caller/Doyensec (#1675)

  • Fix incorrect MATLAB example (#1582)

Thanks to Google's OSS-Fuzz project for finding many of these bugs.

Version 2.7.3

(released December 6, 2020)

... (truncated)

Commits
  • 4d555d0 Bump version to 2.7.4.
  • fc3b05d Update CHANGES.
  • ad21935 Revert "Added dracula theme style (#1636)"
  • e411506 Prepare for 2.7.4 release.
  • 275e34d doc: remove Perl 6 ref
  • 2e7e8c4 Fix several exponential/cubic complexity regexes found by Ben Caller/Doyensec
  • eb39c43 xquery: fix pop from empty stack
  • 2738778 fix coding style in test_analyzer_lexer
  • 02e0f09 Added 'ERROR STOP' to fortran.py keywords. (#1665)
  • c83fe48 support added for css variables (#1633)
  • Additional commits viewable in compare view


Dependabot compatibility score

Dependabot will resolve any conflicts with this PR as long as you don't alter it yourself. You can also trigger a rebase manually by commenting @dependabot rebase.

Dependabot commands and options
You can trigger Dependabot actions by commenting on this PR: - `@dependabot rebase` will rebase this PR - `@dependabot recreate` will recreate this PR, overwriting any edits that have been made to it - `@dependabot merge` will merge this PR after your CI passes on it - `@dependabot squash and merge` will squash and merge this PR after your CI passes on it - `@dependabot cancel merge` will cancel a previously requested merge and block automerging - `@dependabot reopen` will reopen this PR if it is closed - `@dependabot close` will close this PR and stop Dependabot recreating it. You can achieve the same result by closing it manually - `@dependabot ignore this major version` will close this PR and stop Dependabot creating any more for this major version (unless you reopen the PR or upgrade to it yourself) - `@dependabot ignore this minor version` will close this PR and stop Dependabot creating any more for this minor version (unless you reopen the PR or upgrade to it yourself) - `@dependabot ignore this dependency` will close this PR and stop Dependabot creating any more for this dependency (unless you reopen the PR or upgrade to it yourself) - `@dependabot use these labels` will set the current labels as the default for future PRs for this repo and language - `@dependabot use these reviewers` will set the current reviewers as the default for future PRs for this repo and language - `@dependabot use these assignees` will set the current assignees as the default for future PRs for this repo and language - `@dependabot use this milestone` will set the current milestone as the default for future PRs for this repo and language You can disable automated security fix PRs for this repo from the [Security Alerts page](https://github.com/stovorov/ConjugateGradients/network/alerts).

src/utils.h:3:10: fatal error: mkl.h: No such file or directory

opened on 2020-10-19 21:15:47 by kamadforge

When runinng make in the home directory I get the error:

src/utils.h:3:10: fatal error: mkl.h: No such file or directory

Reworking project

opened on 2018-12-22 11:21:09 by p-sto None
Paweł Stoworowicz
GitHub Repository

conjugate-gradient nvidia-cuda numerical-methods gpgpu cuda-kernels numpy linear-equations mkl-pardiso c