=================================== Solving Systems of Linear Equations [TOC] =================================== Compare this example with __flens/examples/lapack-gesv__. Like in __flens/examples/lapack-gesv__ we solve a system of linear equations $Ax = b.$ But instead of using the driver function __lapack::sv__ we first compute the $LU$ factorization of $A$ with __lapack::trf__ and then use the triangular solver __blas::sm__ to finish the job. So basicaly we do manually what __lapack::sv__ does internally. :links: __lapack::sv__ -> file:flens/lapack/gesv/sv.h __lapack::trf__ -> file:flens/lapack/gesv/trf.h __blas::sm__ -> file:flens/blas/level3/sm.h __(.+)__ -> doc:$1 Example Code ============ :import: flens/examples/lapack-getrf.cc [stripped, downloadable] Comments on Example Code ======================== :import: flens/examples/lapack-getrf.cc Compile ======= *--[SHELL]----------------------------------------------------------------* | | | cd flens/examples | | clang++ -std=c++0x -Wall -I../.. -o lapack-getrf lapack-getrf.cc | | | *-------------------------------------------------------------------------* Run === *--[SHELL]----------------------------------------------------------------* | | | cd flens/examples | | ./lapack-getrf | | | *-------------------------------------------------------------------------* Using QD (Double-Double and Quad-Double Package) ================================================ The following description is taken from the __QD Library__ page: *--[BOX]------------------------------------------------------------------* | | | The __QD Library__ supports both a double-double datatype (approx. 32 | | decimal digits) and a quad-double datatype (approx. 64 decimal digits). | | The computational library is written in C++. Both C++ and Fortran-90 | | high-level language interfaces are provided to permit one to use convert| | an existing C++ or Fortran-90 program to use the library with only minor| | changes to the source code. In most cases only a few type statements and| | (for Fortran-90 programs) read/write statements need to be changed. PSLQ| | and numerical quadrature programs are included. | | | *-------------------------------------------------------------------------* :links: __QD Library__ -> http://crd-legacy.lbl.gov/~dhbailey/mpdist/ Modified Example Code --------------------- :import: flens/examples/qd-lapack-getrf.cc [stripped, downloadable] Comments on Modifications ------------------------- Only a few modifications were made which will be pointed out in the following: :import: flens/examples/qd-lapack-getrf.cc [brief] Compile and Run --------------- *--[SHELL]----------------------------------------------------------------------* | | | cd flens/examples | | clang++ -Wall -std=c++0x -I../.. -I/opt/local/include/ -c qd-lapack-getrf.cc | | clang++ -o qd-lapack-getrf qd-lapack-getrf.o /opt/local/lib/libqd.a | | ./qd-lapack-getrf | | | *-------------------------------------------------------------------------------* Using mpfr::real (C++ Interface for mpfr) ========================================= The following description is taken from the __mpfr::real__ page: *--[BOX]------------------------------------------------------------------* | | | The class __mpfr::real__ is a high-level C++ wrapper for the | | __GNU MPFR library__, _a C library for multiple-precision floating-point| | computations with correct rounding_. | | | | The objects of `mpfr::real` can (almost) be used like doubles or other | | fundamental floating-point types, thus avoiding the use of C functions | | to manipulate `MPFR`'s low-level data type directly. In addition to | | all arithmetic and comparison operators available for fundamental | | floating-point types, `mpfr::real` also supports the set of | | mathematical functions for `double` from `math.h/cmath`. Finally, | | `std::istream operator >>` and `std::ostream operator <<` are | | implemented for `mpfr::reals`. This allows to substitute `double` with | | `mpfr::real` with no further modifications of the code in many cases. | | | *-------------------------------------------------------------------------* :links: __mpfr::real__ -> http://chschneider.eu/programming/mpfr_real/ __GNU MPFR library__ -> http://www.mpfr.org/ Modified Example Code --------------------- :import: flens/examples/mpfr-real-lapack-getrf.cc [stripped, downloadable] Comments on Modifications ------------------------- :import: flens/examples/mpfr-real-lapack-getrf.cc [brief] Compile and Run --------------- *--[SHELL]----------------------------------------------------------------* | | | cd flens/examples | | clang++ -Wall -std=c++0x -I../.. -I/opt/local/include +++ | | -o mpfr-real-lapack-getrf mpfr-real-lapack-getrf.cc +++ | | -L /opt/local/lib -lmpfr -lgmp | | ./mpfr-real-lapack-getrf | | | *-------------------------------------------------------------------------* For Performance: Compile with ATLAS =================================== For high performance you have to link with an optimized BLAS implementation like __GotoBLAS__ or __ATLAS__. For ATLAS pass `-DWITH_ATLAS` to the compiler: *--[SHELL]----------------------------------------------------------------* | | | cd flens/examples | | clang++ -DWITH_ATLAS +++| | -Wall -std=c++0x -I../.. +++| | -I/opt/local/include -o atlas-lapack-getrf lapack-getrf.cc +++| | -latlas -lcblas | | ./atlas-lapack-getrf | | | *-------------------------------------------------------------------------* If you want to check which BLAS function gets called also pass `-DDEBUG_CXXBLAS` to the compiler and grep and sort for `CXXBLAS` (note that the problem size in this example is so small that no LEVEL 3 BLAS function gets used): *--[SHELL]----------------------------------------------------------------* | | | cd flens/examples | | clang++ -DWITH_ATLAS -DCXXBLAS_DEBUG +++| | -Wall -std=c++0x -I../.. +++| | -I/opt/local/include -o atlas-lapack-getrf lapack-getrf.cc +++| | -latlas -lcblas | | ./atlas-lapack-getrf 2>&1 | grep CXXBLAS | sort -u | | | *-------------------------------------------------------------------------* For optimal performance also pass `-DNDEBUG` and `-O3`. :links: __ATLAS__ -> http://math-atlas.sourceforge.net __GotoBLAS__ -> http://www.tacc.utexas.edu/tacc-projects/gotoblas2