00001 /* 00002 ------------------------------------------------------------------- 00003 00004 Copyright (C) 2006, 2007, 2008, Andrew W. Steiner 00005 00006 This file is part of O2scl. 00007 00008 O2scl is free software; you can redistribute it and/or modify 00009 it under the terms of the GNU General Public License as published by 00010 the Free Software Foundation; either version 3 of the License, or 00011 (at your option) any later version. 00012 00013 O2scl is distributed in the hope that it will be useful, 00014 but WITHOUT ANY WARRANTY; without even the implied warranty of 00015 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00016 GNU General Public License for more details. 00017 00018 You should have received a copy of the GNU General Public License 00019 along with O2scl. If not, see <http://www.gnu.org/licenses/>. 00020 00021 ------------------------------------------------------------------- 00022 */ 00023 /** 00024 \mainpage O2scl User's Guide 00025 00026 \htmlonly 00027 <div align="right"><a href="dl_page.html"><b>Download 00028 O<span style='position: relative; top: 0.3em; 00029 font-size: 0.8em'>2</span>scl</b></a></div> 00030 \endhtmlonly 00031 00032 \o2 is a C++ class library for object-oriented numerical 00033 programming. It includes 00034 - Classes based on numerical routines from GSL and CERNLIB 00035 - Vector and matrix classes which are fully compatible with 00036 \c gsl_vector and \c gsl_matrix, yet offer indexing with 00037 \c operator[] and other object-oriented features 00038 - The CERNLIB-based classes are rewritten in C++ and are often 00039 faster than their GSL counterparts 00040 - Classes which require function inputs are designed to accept 00041 (public or private) member functions, even if they are virtual. 00042 - Classes use templated vector types, which allow the use of 00043 object-oriented vectors or C-style arrays. 00044 - Highly compatible - Recent versions have been tested on 00045 Linux (32- and 64-bit systems, with Intel and AMD chips), 00046 Windows XP with Cygwin, and MacOSX. 00047 - Free! \o2 is provided under Version 3 of the GNU 00048 Public License 00049 - Two mini-libraries 00050 <ul> 00051 <li>Thermodynamics of ideal and nearly-ideal particles with 00052 quantum statistics</li> 00053 <li>Equations of state for finite density relevant for 00054 neutron stars</li> 00055 </ul> 00056 00057 This is a beta version. The library should install and 00058 test successfully, and most of the classes are ready 00059 for production use. Some of the interfaces may change 00060 slightly in future versions. There are a few classes 00061 which are more experimental, and this is clearly stated 00062 at the top of the documentation for these classes. 00063 00064 See licensing information at \ref license_section. 00065 00066 \htmlonly 00067 <hr /> 00068 <a href="dl_page.html">Download 00069 O<span style='position: relative; top: 0.3em; font-size: 00070 0.8em'>2</span>scl</a> 00071 <a href="../latex/refman.pdf">PDF documentation</a> 00072 \endhtmlonly 00073 00074 \hline 00075 \section toc_section Quick Reference to User's Guide 00076 00077 - \ref install_section 00078 - \ref usage_section 00079 - \ref examples_section 00080 - \ref related_section 00081 - \ref cplxnum_section 00082 - \ref vecmat_section 00083 - \ref permute_section 00084 - \ref linalg_section 00085 - \ref intp_section 00086 - \ref const_section 00087 - \ref funct_section 00088 - \ref table_section 00089 - \ref string_section 00090 - \ref diff_section 00091 - \ref inte_section 00092 - \ref poly_section 00093 - \ref solve_section 00094 - \ref min_section 00095 - \ref conmin_section 00096 - \ref mcarlo_section 00097 - \ref anneal_section 00098 - \ref fit_section 00099 - \ref ode_section 00100 - \ref rng_section 00101 - \ref tintp_section 00102 - \ref gslcheb_section 00103 - \ref unitconv_section 00104 - \ref other_section 00105 - \ref lset_section 00106 - \ref io_section 00107 - \ref source_section 00108 - \ref design_section 00109 - \ref license_section 00110 - \ref ack_section 00111 - \ref ref_section 00112 - \ref todo 00113 - \ref bug 00114 00115 \htmlonly 00116 <UL><LI>Other mini-libraries 00117 <UL> 00118 <LI><B><a href="../part/html/index.html">Particles</a></b> 00119 <LI><B><a href="../eos/html/index.html">Equations of State</a></b> 00120 </UL> 00121 </LI> 00122 <br /> 00123 <li>Add-on library: There is also a related library, 00124 <a href='../../o2scl_ext/html/index.html'>O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl_ext</a>, 00125 which has a separate source distribution. 00126 </li> 00127 </UL> 00128 \endhtmlonly 00129 00130 \hline 00131 \section install_section Installation 00132 00133 The rules for installation are generally the same as that for 00134 other GNU libraries. The file \c INSTALL has some details on this 00135 procedure. Generally, you should be able to run 00136 <tt>./configure</tt> and then type \c make and \c make \c 00137 install. More information on the \c configure command can also be 00138 obtained from <tt>./configure --help</tt>. \o2 requires the GSL 00139 libraries. If the <tt>configure</tt> script cannot find them, you 00140 may have to specify their location in the <tt>CPPFLAGS</tt> and 00141 <tt>LDFLAGS</tt> environment variables (<tt>./configure 00142 --help</tt> shows some information on this). The documentation is 00143 included in the distribution and automatically installed by 00144 \c make \c install. 00145 00146 After \c make \c install, you may test the library with \c make \c 00147 o2scl-test. At the end, the phrase <tt>All O2scl tests passed.</tt> 00148 indicates that all the testing was successful. 00149 00150 This library requires GSL and is designed to work with GSL 00151 versions 1.12 or greater. Some classes may work with older 00152 versions of GSL, but this cannot be guaranteed. 00153 00154 Range-checking for vectors and matrices is performed similar to 00155 the GSL approach, and is turned on by default. You can 00156 disable range-checking by defining -DO2SCL_NO_RANGE_CHECK 00157 \code 00158 CPPFLAGS="-DO2SCL_NO_RANGE_CHECK" ./configure 00159 \endcode 00160 00161 The separate libraries \o2e and \o2p are installed by default. To 00162 disable the installation of these libraries and their associated 00163 documentation, run <tt>./configure</tt> with the flags \c 00164 --disable-eoslib or \c --disable-partlib. Note that \o2e depends 00165 on \o2p so using \c --disable-partlib without \c --disable-eoslib 00166 will not work. 00167 00168 There are several warning flags that are useful when configuring 00169 and compiling with \o2. See the GSL documentation for an 00170 excellent discussion, and also see the generic installation 00171 documentation in the file <tt>INSTALL</tt> in the \o2 top-level 00172 directory. For running \c configure, for example, if you do 00173 not have privileges to write to <tt>/usr/local</tt>, 00174 \code 00175 CPPFLAGS="-O3" -I/home/asteiner/install/include" \ 00176 LDFLAGS="-L/home/asteiner/install/lib" ./configure -C \ 00177 --prefix=/home/asteiner/install 00178 \endcode 00179 In this example, specifying 00180 <tt>-I/home/asteiner/install/include</tt> and 00181 <tt>-L/home/asteiner/install/lib</tt> above ensures that the GSL 00182 libraries can be found (this is where they are installed on my 00183 machine). The <tt>--prefix=/home/asteiner/install</tt> argument to 00184 <tt>./configure</tt> ensures that \o2 is installed there as well. 00185 00186 The documentation is generated with \doxygen. In principle, the 00187 documentation can be regenerated by the end-user, but this is 00188 not supported and requires several external applications 00189 not included in the distribution. 00190 00191 \b Un-installation: While there is no explicit "uninstall" 00192 procedure, there are only a couple places to check. Installation 00193 creates directories named <tt>o2scl</tt> in the include, doc and shared 00194 files directory (which default to \c /usr/local/include, \c 00195 /usr/local/doc, and \c /usr/local/share) which can be 00196 removed. Finally, all of the libraries are named with the prefix 00197 \c libo2scl and are created by default in /usr/local/lib. As 00198 configured with the settings above, the files are in \c 00199 /home/asteiner/install/include/o2scl, \c /home/asteiner/install/lib, 00200 \c /home/asteiner/install/share/o2scl, and \c 00201 /home/asteiner/install/doc/o2scl. 00202 00203 \hline 00204 \section usage_section General Usage 00205 00206 \subsection namespace_subsect Namespaces 00207 00208 Most of the classes reside in the namespace \c o2scl (this 00209 namespace has been removed from the documentation for clarity). 00210 Numerical constants (many of them based on the GSL constants) are 00211 placed in separate namespaces (gsl_cgs, gsl_cgsm, gsl_mks, 00212 gsl_mksa, gsl_num, o2scl_const, and o2scl_fm). There are also two 00213 namespaces which hold integration coefficients, 00214 o2scl_inte_qag_coeffs and o2scl_inte_qng_coeffs. There are also 00215 some namespaces for the linear algebra functions, see \ref 00216 linalg_section for more information. 00217 00218 \subsection docconv_subsect Documentation conventions 00219 00220 In the following documentation, function parameters are denoted by 00221 \c parameter, except when used in mathematical formulas as in \f$ 00222 \mathrm{variable} \f$ . 00223 00224 \subsection nomen_subsect Nomenclature 00225 00226 Classes directly derived from the GNU Scientific Library are 00227 preceeded by the prefix \c gsl_ and classes derived from CERNLIB 00228 are preceeded by the prefix \c cern_. Some of those classes 00229 derived from GSL and CERNLIB operate slightly differently 00230 from the original versions. The differences are detailed 00231 in the corresponding class documentation. 00232 00233 \subsection errorhand_subsect Basic error handling 00234 00235 Error handling is similar to GSL. When an error occurs, functions 00236 and/or classes call a GSL-like error handler and (when 00237 appropriate) return a non-zero value. When functions succeed they 00238 return 0 (\ref gsl_success). The error handler, \ref err_hnd is a 00239 global pointer to an object of type \ref err_base. There is a 00240 global default error handler of type \ref err_class. The list of 00241 error codes is given in the documentation for the file \ref 00242 err_hnd.h. The default error handler can be replaced by simply 00243 assigning the address of a descendant of \ref err_base to \ref 00244 err_hnd. Most of the time, if a function returns an integer value, 00245 and the returned value is nonzero, then the error handler was 00246 called (at least once). Functions can return success even when the 00247 error handler was called during execution, and this happens 00248 particularly in template classes where templated methods are 00249 beyond the direct control of the function using them. Also, \o2 00250 functions never reset the error handler. 00251 00252 The default behavior for all errors is to store the error 00253 information, print the information to <tt>cout</tt> and call 00254 <tt>exit(error_number)</tt>. You can modify the behavior of the 00255 default error handler with the err_class::set_mode() function, 00256 Mode '2' is the default, mode '1' just prints the information, 00257 and mode '0' ignores the error. 00258 00259 Although \o2 does not yet have this functionality by default, it 00260 is striaghtforward to define an error handler which throws a C++ 00261 exception. Object destructors do not generally call the error 00262 handler. Internally, \o2 does not use <tt>try</tt> blocks, but 00263 these can easily be effectively employed by an \o2 user. 00264 00265 Errors can be set by the user through the macros \ref O2SCL_ERR, 00266 which sets an error, and \ref O2SCL_ERR_RET, which sets an error 00267 and returns the error number. 00268 00269 Functionality similar to assert() is provided with the macro 00270 \ref O2SCL_ASSERT, which exits if its argument is non-zero, and 00271 \ref O2SCL_BOOL_ASSERT which exits if its argument is false. 00272 00273 Note that the GSL function <tt>gsl_set_error_handler()</tt> is 00274 called in the abstract base class err_base to set the GSL error 00275 handler any time a new instance of \ref err_base is created. If 00276 the user creates a new handler, that new handler will be used by 00277 GSL functions as well. 00278 00279 \subsection errorhand2_subsect What is an error? 00280 00281 \o2 assumes that errors are events which should happen 00282 infrequently. Error handling strategies are often time-consuming 00283 and they are not a replacement for normal code flow. However, even 00284 with this in mind, one can still distinguish a large spectrum of 00285 posibillities from "fatal" errors, those likely to corrupt the 00286 stack and/or cause a dreaded "segmentation fault" and "non-fatal" 00287 errors, those errors which might cause incorrect results, but 00288 might be somehow recoverable. One of the purposes of error 00289 handling is to decide if and how these different types of errors 00290 should be handled differently. Most errors in \o2 result in a call 00291 to the error handler by default. Some of the classes which attempt 00292 to reach numerical convergence have an option (e.g. 00293 mroot::err_nonconv) to turn this default behavior off for 00294 convergence errors (which are not necessarily "fatal" errors). To 00295 set these "convergence" errors, the macros \ref O2SCL_CONV and 00296 \ref O2SCL_CONV_RET can be used. 00297 00298 Another related issue is that \o2 often calls functions which are 00299 supplied by the user, these user-designed functions may create 00300 errors, and the library needs to decide how to deal with them, 00301 even though it knows little about what is actually happening 00302 inside these user-defined functions. Most of the time, \o2 assumes 00303 that if a function returns a nonzero value, then an error has 00304 occured and the present calculation should abort. 00305 00306 \subsection objscope_subsect Objects and scope 00307 00308 \o2 objects frequently take inputs which are of the form of a 00309 reference to a smaller object. This is particularly convenient 00310 because it allows a lot of flexibility, while providing a certain 00311 degree of safety. However, the user retains the responsibility of 00312 ensuring that input objects do not go out of scope before they are 00313 utilized by objects which require them. (This is actually no 00314 different than the requirements on the user imposed by GSL, for 00315 example.) 00316 00317 A simple example of this is provided by the simulated annealing 00318 class \ref gsl_anneal. Simulated annealing works by providing a 00319 temperature annealing schedule, which slowly decreases the 00320 temperature. However, if the temperature schedule object goes out 00321 of scope before the annealing is done, a segmentation fault will 00322 result. For example, the following code will fail: 00323 \include scope_doc.cpp 00324 00325 This will fail because the temperature schedule object goes out of 00326 scope (it's a local variable in the function and its destructor is 00327 called before the <tt>set_schedule()</tt> function exits) before 00328 the <tt>mmin()</tt> function is called. 00329 00330 \hline 00331 \section examples_section Compiling examples 00332 00333 The example programs are in the \c examples directory. After 00334 installation, they can be compiled and executed by running \c make 00335 \c o2scl-examples in that directory. This will also test the 00336 output of the examples to make sure it is correct and summarize 00337 these tests in \c examples-summary.txt. The output for each 00338 example is placed in the corresponding file with a <tt>.scr</tt> 00339 extension. 00340 00341 Alternatively, you can make the executable for each example in the 00342 \c examples directory individually using, e.g. <tt>make 00343 ex_mroot</tt>. 00344 00345 See \ref source_section for the documented source code 00346 of the individual examples 00347 00348 Also, the testing code for each class is occasionally useful for 00349 providing examples of their usage. The testing source code for 00350 each source file is named with an \c _ts.cpp prefix in the same 00351 directory as the class source. 00352 00353 \hline 00354 \section related_section Related projects 00355 00356 Several noteworthy related projects: 00357 00358 - GSL - The GNU Scientific Library \n The first truly free, 00359 ANSI-compliant, fully tested, and well-documented scientific 00360 computing library. The GSL is located at 00361 http://www.gnu.org/software/gsl . Manual is available at 00362 http://www.gnu.org/software/gsl/manual/html_node/ . Many GSL 00363 routines are included and reworked in \o2 (the corresponding 00364 classes begin with the \c gsl_ prefix) and \o2 was specifically 00365 designed to be used with GSL. GSL is a must-have for most serious 00366 numerical work in C or C++ and is required for installation of 00367 \o2. 00368 00369 - CERNLIB - The gold standard in FORTRAN computing libraries \n 00370 Several CERNLIB routines are rewritten completely in C++ and 00371 included in \o2 (they begin with the \c cern_ prefix). CERNLIB is 00372 located at http://cernlib.web.cern.ch/cernlib/mathlib.html 00373 00374 - LAPACK and ATLAS - The gold standard for linear algebra \n 00375 Available at http://www.netlib.org/lapack and 00376 http://www.netlib.org/atlas . See also 00377 http://www.netlib.org/clapack . 00378 00379 - QUADPACK - FORTRAN adaptive integration library \n This is the 00380 library on which the GSL integration routines are based (prefix \c 00381 gsl_inte_). It is available at http://www.netlib.org/quadpack . 00382 00383 - TNT - Template numerical toolkit (http://math.nist.gov) \n TNT 00384 provides vector and matrix types and basic arithmetic operations. 00385 Most of the classes in \o2 which use vector and matrix types are 00386 templated to allow compatbility with TNT. (Though there are a few 00387 small differences.) This software is in the public domain. 00388 00389 - Blitz++ - http://www.oonumerics.org/blitz \n Another linear 00390 algebra library designed in C++ which includes vector and matrix 00391 types and basic arithmetic. As the name implies, Blitz++ has the 00392 capability to be particularly fast. Distributed with a Perl-like 00393 artistic license "or the GPL". 00394 00395 - Boost - http://www.boost.org \n Free (license is compatible with 00396 GPL) peer-reviewed portable C++ source libraries that work well 00397 with the C++ Standard Library. Boost also contains uBlas, 00398 for linear algebra computing. For examples of how to use 00399 \o2 with uBlas vector and matrices, see \ref ex_ublas_sect. 00400 00401 - FLENS - http://flens.sourceforge.net \n A C++ library with 00402 object-oriented linear algebra types and an interface to BLAS 00403 and LAPACK. 00404 00405 - MESA - Modules for Experiments in Stellar Astrophysics 00406 (http://mesa.sourceforge.net) \n An excellent FORTRAN library with 00407 accurate low-density equations of state, interpolation, opacities 00408 and other routines useful in stellar physics. Work is currently 00409 under way to rewrite some of the MESA routines in C/C++ for \o2. 00410 Licensed with LGPL (not GPL). 00411 00412 - OOL - Open Optimization Library (http://ool.sourceforge.net) \n 00413 Constrained minimization library designed with GSL in mind. The 00414 \o2 constrained minimization classes are derived from this 00415 library. See \ref conmin_section . 00416 00417 - Root - CERN's new C++ analysis package (http://root.cern.ch) \n 00418 A gargantuan library for data analysis, focused mostly on 00419 high-energy physics. Their histograms, graphics, file I/O and 00420 support for large data sets is particuarly good. 00421 00422 \hline 00423 \section cplxnum_section Complex Numbers 00424 00425 Several arithmetic operations for \c gsl_complex are 00426 defined in cx_arith.h, but no constructor has been written. The 00427 object \c gsl_complex is still a \c struct, not a \c class. For 00428 example, 00429 \include cx_arith_doc.cpp 00430 00431 In case the user needs to convert between \c gsl_complex and \c 00432 std::complex<double>, two conversion functions gsl_to_complex() 00433 and complex_to_gsl() are also provided in \ref cx_arith.h. 00434 00435 \hline 00436 \section vecmat_section Arrays, Vectors, Matrices and Tensors 00437 00438 \subsection vmintro_subsect Introduction 00439 00440 The \o2 library uses a standard nomenclature to distinguish 00441 a couple different concepts. The word "array" is always used to 00442 refer to C-style arrays, i.e. \c double[]. If there are two 00443 dimensions in the array, it is a "two-dimensional array", i.e. \c 00444 double[][] . The word "vector" is reserved generic objects with 00445 array-like semantics, i.e. any type of object or class which can 00446 be treated similar to an array in that it has an function \c 00447 operator[] which gives array-like indexing. Thus arrays are 00448 vectors, but not all vectors are arrays. There are a couple 00449 specific vector types defined in \o2 like \ref ovector and \ref 00450 uvector . The STL vector <tt>std::vector<double></tt> is also a 00451 vector type in this language. The word "matrix" is reserved for 00452 the a generic object which has matrix-like semantics and can be 00453 accessed using either \c operator() or two successive applications 00454 of \c operator[] (or sometimes both). The \o2 objects \ref omatrix 00455 and \ref umatrix are matrix objects, as is a C-style 00456 two-dimensional array, \c double[][] . The header files are named 00457 in this spirit also: functions and classes which operate on 00458 generic vector types are in \ref vector.h and functions and 00459 classes which work only with C-style arrays are in \ref array.h . 00460 The word "tensor" is used for a generic object which has 00461 rank \c n and then has \c n associated indices. A vector is just 00462 a \tensor of rank 1 and a matrix is just a \tensor of rank 2. 00463 00464 Most of the classes in \o2 which use vectors and/or matrices are 00465 designed so that they can be used with any appropriately-defined 00466 vector or matrix types. This is a major part of the design goals 00467 for \o2 and most of the classes are compatible with matrix and 00468 vector objects from GSL, TNT, MV++, uBlas, and Blitz++. 00469 00470 The first index of a matrix type is defined always to be the index 00471 associated with "rows" and the second is associated with 00472 "columns". The \o2 matrix output functions respect this notation 00473 as well, so that all of the elements of the first row is sent to 00474 the screen, then all of the elements in the second row, and so on. 00475 With this in mind, one can make the distinction between 00476 "row-major" and "column-major" matrix storage. C-style 00477 two-dimensional arrays are "row-major" in that the elements 00478 of the first row occupy the first locations in memory as 00479 opposed "column-major" where the first column occupies the 00480 first locations in memory. It is important to note that the 00481 majority of the classes in \o2 do not care about the details 00482 of the underlying memory structure, so long as two 00483 successive applications of <tt>operator[]</tt> (or 00484 in some cases <tt>operator(,)</tt> ) works properly. The storage 00485 format used by \ref omatrix and \ref umatrix is row-major, and there 00486 are no column-major matrix classes in \o2 (yet). 00487 00488 \subsection index_subsect Matrix indexing with [][] vs. (,) 00489 00490 While vector indexing with \c operator[] is very compatible with 00491 almost any vector type, matrix indexing is a bit more difficult. 00492 There are two options: assume matrix objects provide an \c 00493 operator[] method which can be applied twice, i.e. 00494 <tt>m[i][j]</tt>, or assume that matrix elements should be 00495 referred to with <tt>m(i,j)</tt>. Most of the \o2 classes use the 00496 former approach so that they are also compatible with 00497 two-dimensional arrays. However, there are sometimes good reasons 00498 to want to use <tt>operator()</tt> for matrix-intensive operations 00499 from linear algebra. For this reason, some of the functions given 00500 in the <tt>linalg</tt> directory are specified in two forms: the 00501 first default form which assumes <tt>[][]</tt>, and the second 00502 form with the same name, but in a namespace which has a suffix 00503 <tt>_paren</tt> and assumes matrix types use <tt>(,)</tt>. 00504 00505 \subsection rowcol_subsect Rows and columns vs. x and y 00506 00507 Sometimes its useful to think about the rows and columns in a 00508 matrix as referring to a elements of a grid, and the matrix 00509 indices refer to points in a grid in \f$ (x,y) \f$. It might seem 00510 intuitive to think of a matrix as \c A[ix][iy] where \c ix and \c 00511 iy are the \f$ x \f$ and \f$ y \f$ indices because the ordering of 00512 the indices is alphabetical. However, it is useful to note that 00513 because functions like \ref matrix_out print the first "row" first 00514 rather than the first column, a matrix constructed as \c A[ix][iy] 00515 will be printed out with x on the "vertical axis" and y on the 00516 "horizontal axis", so it is sometimes useful to store data 00517 in the form \c A[iy][ix] (for example, in the two dimensional 00518 interpolation class, twod_intp). In any case, all classes 00519 which take matrix data as input will document how the matrix 00520 ought to be arranged. 00521 00522 \subsection genvec_subsection Generic vector functions 00523 00524 There are a couple functions which operate on generic vectors of 00525 any type in \ref vector.h . They perform sorting, summing, 00526 rotating, copying, and computations of minima and maxima. For 00527 more statistically-oriented operations, see also \ref vec_stats.h . 00528 00529 \subsection nativetypes_subsect Native matrix and vector types 00530 00531 The rest of this section in the User's guide is dedicated to 00532 describing the \o2 implementations of vector, matrix, 00533 and \tensor types. 00534 00535 Vectors and matrices are designed using the templates \ref 00536 ovector_tlate and \ref omatrix_tlate, which are compatible with \c 00537 gsl_vector and \c gsl_matrix. Vectors and matrices with unit 00538 stride are provided in \ref uvector_tlate and \ref umatrix_tlate. 00539 The most commonly used double-precision versions of these template 00540 classes are \ref ovector, \ref omatrix, \ref uvector and \ref 00541 umatrix, and their associated "views" (analogous to GSL vector and 00542 matrix views) which are named with a \c _view suffix. The classes 00543 \ref ovector_tlate and \ref omatrix_tlate offer the syntactic 00544 simplicity of array-like indexing, which is easy to apply to 00545 vectors and matrices created with GSL. 00546 00547 The following sections primarily discuss the operation objects of 00548 type \ref ovector. The generalizations to objects of type \ref 00549 uvector, \ref omatrix, and the complex vector and matrix objects 00550 \ref ovector_cx, \ref omatrix_cx, \ref uvector_cx, and \ref 00551 umatrix_cx are straightforward. 00552 00553 \subsection view_subsect Vector and matrix views 00554 00555 Vector and matrix views are provided as parents of the vector and 00556 matrix classes which do not have methods for memory allocation. As 00557 in GSL, it is simple to "view" normal C-style arrays or pointer 00558 arrays (see \ref ovector_array), parts of vectors (\ref 00559 ovector_subvector), and rows and columns of matrices (\ref 00560 omatrix_row and \ref omatrix_col). Several operations are defined, 00561 including addition, subtraction, and dot products. 00562 00563 \subsection typedef_subsect Vector and matrix typedefs 00564 00565 Several typedefs are used to give smaller names to often used 00566 templates. Vectors and matrices of double-precision numbers all 00567 begin with the prefixes \ref ovector and \ref omatrix. Integer 00568 versions begin with the prefixes \ref ovector_int and \ref 00569 omatrix_int. Complex versions have an additional \c "_cx", 00570 e.g. \ref ovector_cx. See \ref ovector_tlate.h, \ref 00571 ovector_cx_tlate.h, \ref omatrix_tlate.h, and \ref 00572 omatrix_cx_tlate.h. 00573 00574 \subsection unit_subsect Unit-stride vectors 00575 00576 The \ref uvector_tlate objects are naturally somewhat faster 00577 albeit less flexible than their finite-stride counterparts. 00578 Conversion to GSL vectors and matrices is not trivial for \ref 00579 uvector_tlate objects, but demands copying the entire 00580 vector. Vector views, operators, and the nomenclature is similar 00581 to that of \ref ovector. 00582 00583 \subsection memalloc_subsect Memory allocation 00584 00585 Memory for vectors can be allocated using ovector::allocate() and 00586 ovector::free(). Allocation can also be performed by the 00587 constructor, and the destructor automatically calls <tt>free()</tt> if 00588 necessary. In contrast to \c gsl_vector_alloc(), 00589 ovector::allocate() will call ovector::free(), if necessary, to 00590 free previously allocated space. Allocating memory does not clear 00591 the recently allocated memory to zero. You can use \ref 00592 ovector::set_all() with a zero argument to clear a vector 00593 (and similarly for a matrix). 00594 00595 If the memory allocation fails, either in the constructor or in 00596 allocate(), then the error handler will be called, partially 00597 allocated memory will be freed, and the size will be reset to 00598 zero. 00599 00600 Although memory allocation in \o2 vectors is very similar to 00601 that in GSL, the user must not mix allocation and deallocation 00602 between GSL and \o2. 00603 00604 \subsection getset_subsect Get and set 00605 00606 Vectors and matrices can be modified using ovector::get() and 00607 ovector::set() methods analogous to \c gsl_vector_get() and \c 00608 gsl_vector_set(), or they can be modified through 00609 ovector::operator[] (or ovector::operator() ), e.g. 00610 \comment 00611 For some reason the references to ovector::operator have been 00612 failing in the Doxygen compilation. 00613 \endcomment 00614 \code 00615 ovector a(4); 00616 a.set(0,2.0); 00617 a.set(1,3.0); 00618 a[2]=4.0; 00619 a[3]=2.0*a[1]; 00620 \endcode 00621 00622 If you want to set all of the values in an \ref ovector or an \ref 00623 omatrix at the same time, then use ovector::set_all() or 00624 omatrix::set_all(). 00625 00626 \subsection rangecheck_subsect Range checking 00627 00628 Range checking is performed depending on whether or not \c 00629 O2SCL_NO_RANGE_CHECK is defined. It can be defined in the 00630 arguments to \c ./configure upon installation to turn off range 00631 checking. Note that this is completely separate from the GSL range 00632 checking mechanism, so range checking may be on in \o2 even if 00633 it has been turned off in GSL. Range checking is used primarily 00634 in the vector, matrix, and \tensor \c get() and \c set() methods. 00635 00636 To see if range checking was turned on during installation 00637 (without calling the error handler), use 00638 lib_settings_class::range_check(). 00639 00640 Note that range checking in \o2 code is most often present in 00641 header files, rather than in source code. This means that range 00642 checking can be turned on or off in user-defined functions 00643 separately from whether or not it was used in the library classes 00644 and functions. 00645 00646 \subsection copy_subsect Shallow and deep copy 00647 00648 Copying \o2 vectors using constructors or the \c = operator 00649 is performed according to what kind of object is 00650 on the left-hand side (LHS) of the equals sign. If the LHS is a 00651 view, then a shallow copy is performed, and if the LHS is a \ref 00652 ovector, then a deep copy is performed. If an attempt is made to 00653 perform a deep copy onto a vector that has already been allocated, 00654 then that previously allocated memory is automatically freed. 00655 The user must be careful to ensure that information is not lost 00656 this way, even though no memory leak will occur. 00657 00658 For generic deep vector and matrix copying, you can use the 00659 template functions \ref vector_copy(), \ref matrix_copy(). These 00660 would allow you, for example, to copy an \ref ovector to a \c 00661 std::vector<double> object. These functions do not do any 00662 memory allocation so that must be handled beforehand by the user. 00663 00664 \subsection vecarith_subsect Vector and matrix arithmetic 00665 00666 Several operators are available as member functions of the 00667 corresponding template: 00668 00669 Vector_view unary operators: 00670 - vector_view += vector_view 00671 - vector_view -= vector_view 00672 - vector_view += scalar 00673 - vector_view -= scalar 00674 - vector_view *= scalar 00675 - scalar = norm(vector_view) 00676 00677 Matrix_view unary operators: 00678 - matrix_view += matrix_view 00679 - matrix_view -= matrix_view 00680 - matrix_view += scalar 00681 - matrix_view -= scalar 00682 - matrix_view *= scalar 00683 00684 Binary operators like addition, subtraction, and matrix 00685 multiplication are also defined for \ref ovector, \ref uvector, 00686 and related objects. The generic template for a binary operator, 00687 e.g. 00688 \code 00689 template<class vec_t> vec_t &operator+(vec_t &v1, vec_t &v2); 00690 \endcode 00691 is difficult because the compiler has no way of distinguishing vector 00692 and non-vector classes. At the moment, this is solved by 00693 creating a define macro for the binary operators. In addition 00694 to the predefined operators for native classes, the user may 00695 also define binary operators for other classes using the same 00696 macros. For example, 00697 \code 00698 O2SCL_OP_VEC_VEC_ADD(o2scl::ovector,std::vector<double>, 00699 std::vector<double>) 00700 \endcode 00701 would provide an addition operator for \ref ovector and vectors 00702 from the Standard Template Library. A full list of the 00703 operators which are defined by default and the associated macros 00704 are detailed in the documentation for \ref vec_arith.h. 00705 00706 The GSL BLAS routines can also be used directly with \ref ovector 00707 and \ref omatrix objects. 00708 00709 Note that some of these arithmetic operations succeed even with 00710 non-matching vector and matrix sizes. For example, adding a 3x3 00711 matrix to a 4x4 matrix will result in a 3x3 matrix and the 7 outer 00712 elements of the 4x4 matrix are ignored. 00713 00714 \subsection gslconv_subsect Converting to and from GSL forms 00715 00716 Because of the way \ref ovector is constructed, you may use 00717 type conversion to convert to and from objects of type \c gsl_vector. 00718 \code 00719 ovector a(2); 00720 a[0]=1.0; 00721 a[1]=2.0; 00722 gsl_vector *g=(gsl_vector *)(&a); 00723 cout << gsl_vector_get(g,0) << " " << gsl_vector_get(g,1) << endl; 00724 \endcode 00725 Or, 00726 \code 00727 gsl_vector *g=gsl_vector_alloc(2); 00728 gsl_vector_set(0,1.0); 00729 gsl_vector_set(1,2.0); 00730 ovector &a=(ovector &)(*g); 00731 cout << a[0] << " " << a[1] << endl; 00732 \endcode 00733 This sort of type-casting is discouraged among unrelated classes, 00734 but is permissible here because ovector_tlate is a descendant of 00735 gsl_vector. In particular, this will not generate "type-punning" 00736 warnings in later gcc versions. If this bothers your 00737 sensibilities, however, then you can use the following approach: 00738 \code 00739 ovector a(2); 00740 gsl_vector *g=a.get_gsl_vector(); 00741 \endcode 00742 The ease of converting between these two kind of objects makes it 00743 easy to use gsl functions on objects of type \ref ovector, i.e. 00744 \code 00745 ovector a(2); 00746 a[0]=2.0; 00747 a[1]=1.0; 00748 gsl_vector_sort((gsl_vector *)(&a)); 00749 cout << a[0] << " " << a[1] << endl; 00750 \endcode 00751 00752 \subsection stlconv_subsect Converting from STL form 00753 00754 To "view" a \c std::vector<double>, you can use \ref ovector_array 00755 \code 00756 std::vector<double> d; 00757 d.push_back(1.0); 00758 d.push_back(3.0); 00759 ovector_array aa(d.size,&(d[0])); 00760 cout << aa[0] << " " << aa[1] << endl; 00761 \endcode 00762 00763 However, you should note that if the memory for the \c std::vector 00764 is reallocated (for example because of a call to \c push_back()), 00765 then a previously created \ref ovector_view will be incorrect. 00766 00767 \subsection pushpop_subsect push_back() and pop() methods 00768 00769 These two functions give a behavior similar to the corresponding 00770 methods for \c std::vector<>. This will work in \o2 classes, 00771 but may not be compatible with all of the GSL functions. This 00772 will break if the address of a \ref ovector_tlate is given 00773 to a GSL function which accesses the \c block->size parameter 00774 instead of the \c size parameter of a \c gsl_vector. Please 00775 contact the author of \o2 if you find a GSL function with this 00776 behavior. 00777 00778 \subsection view_subsect Views 00779 00780 Views are slightly different than in GSL in that they are now 00781 implemented as parent classes. Effectively, this means that vector 00782 views are just the same as normal vectors, except that they have 00783 no memory allocation functions (because the memory is presumably 00784 owned by a different object or scope). The code 00785 \code 00786 double x[2]={1.0,2.0}; 00787 gsl_vector_view_array v(2,x); 00788 gsl_vector *g=&(v.vector); 00789 gsl_vector_set(g,0,3.0); 00790 cout << gsl_vector_get(g,0) << " " << gsl_vector_get(g,1) << endl; 00791 \endcode 00792 can be replaced by 00793 \code 00794 double x[2]={1.0,2.0}; 00795 ovector_array a(2,x); 00796 a[0]=3.0; 00797 cout << a << endl; 00798 \endcode 00799 00800 \subsection vecparm_subsect Passing ovector parameters 00801 00802 It is often best to pass an \ref ovector as a const reference 00803 to an \ref ovector_view, i.e. 00804 \code 00805 void function(const ovector_view &a); 00806 \endcode 00807 If the function may change the values in the \ref ovector, 00808 then just leave out \c const 00809 \code 00810 void function(ovector_view &a); 00811 \endcode 00812 This way, you ensure that the function is not allowed to modify 00813 the memory for the vector argument. 00814 00815 If you intend for a function (rather than the user) to handle 00816 the memory allocation, then some care is necessary. The 00817 following code 00818 \include ovector_doc1.cpp 00819 is confusing because the user may have already allocated memory 00820 for \c a. To avoid this, you may want to ensure that the user 00821 sends an empty vector. For example, 00822 \include ovector_doc2.cpp 00823 In lieu of this, it is often preferable to use a local variable for 00824 the storage and offer a \c get() function, 00825 \include ovector_doc3.cpp 00826 The \o2 classes run into this situation quite frequently, but 00827 the vector type is specified through a template 00828 \include ovector_doc4.cpp 00829 00830 \subsection eqoper_subsect Vectors and operator=() 00831 00832 An "operator=(value)" method for setting all vector elements to 00833 the same value is not included because it leads to confusion 00834 between, \c ovector_tlate::operator=(const data_t &val) and 00835 ovector_tlate::ovector_tlate(size_t val) For example, after 00836 implementing \c operator=() and executing the following 00837 \code 00838 ovector_int o1=2; 00839 ovector_int o2; 00840 o2=2; 00841 \endcode 00842 \c o1 will be a vector of size two, and \c o2 will be an 00843 empty vector! 00844 00845 To set all of the vector elements to the same value, use 00846 ovector_tlate::set_all(). Because of the existence of 00847 constructors like ovector_tlate::ovector_tlate(size_t val), the 00848 following code 00849 \code 00850 ovector_int o1=2; 00851 \endcode 00852 still compiles, and is equivalent to 00853 \code 00854 ovector_int o1(2); 00855 \endcode 00856 while the code 00857 \code 00858 ovector_int o1; 00859 o1=2; 00860 \endcode 00861 will not compile. As a matter of style, \c ovector_int \c o1(2); 00862 is preferable to \c ovector_int \c o1=2; to avoid confusion. 00863 00864 \subsection matstruct_subsect Matrix structure 00865 00866 The matrices from \ref omatrix_tlate 00867 are structured in exactly the same way as in GSL. 00868 For a matrix with 2 rows, 4 colums, and a "tda" or "trailing 00869 dimension" of 7, the memory for the matrix is structured in the 00870 following way: 00871 00872 \verbatim 00873 00 01 02 03 XX XX XX 00874 10 11 12 13 XX XX XX 00875 \endverbatim 00876 00877 where \c XX indicates portions of memory that are unreferenced. The tda 00878 can be accessed through, for example, the method 00879 omatrix_view_tlate::tda(). The \c get(size_t,size_t) methods 00880 always take the row index as the first argument and the column 00881 index as the second argument. The matrices from \ref umatrix_tlate 00882 have a trailing dimension which is always equal to the number of 00883 columns. 00884 00885 \subsection rev_subsect Reversing the order of vectors 00886 00887 You can get a reversed vector view from \ref 00888 ovector_reverse_tlate, or uvector_reverse_tlate. For these 00889 classes, \c operator[] and related methods are redefined to 00890 perform the reversal. If you want to make many calls to these 00891 indexing methods for a reversed vector, then simply copying the 00892 vector to a reversed version may be faster. 00893 00894 \subsection constvec_subsect Const-correctness with vectors 00895 00896 Vector objects, like \ref ovector, can be const in the usual way. 00897 However because vector views, like \ref ovector_view are like 00898 pointers to vectors, there are two ways in which they can be 00899 const. The view can point to a const vector, or the view can be a 00900 'const view' which is fixed to point only to one vector (whether 00901 or not the vector pointed to happens to be const). This is exactly 00902 analogous to the distinction between const double *p, double * 00903 const p, and const double * const p. The first is a non-const 00904 pointer to const data, the second is a const pointer to non-const 00905 data, and the third is a const pointer to const data. The 00906 equivalent expressions in \o2 are 00907 \code 00908 ovector_const_view v1; 00909 const ovector_view v2; 00910 const ovector_const_view v3; 00911 \endcode 00912 Explicitly, 00913 - <tt>ovector_const_view</tt> is a view of a const vector, the 00914 view may change, but the vector may not. 00915 - <tt>const ovector_const_view</tt> is a const view of a const 00916 vector, the view may not point to a different vector and the 00917 vector may not change. 00918 - <tt>const ovector_view</tt> is a const view of a normal vector, 00919 the view may not change, but the vector can. 00920 This same distinction is also present, for example, in ublas vectors 00921 views within boost. 00922 00923 A reference of type \ref ovector_base is often used as a function 00924 argument, and can hold either a \ref ovector or a \ref 00925 ovector_view. The important rule to remember with \ref ovector_base 00926 is that, a const reference, i.e. 00927 \code 00928 const ovector_base &v; 00929 \endcode 00930 is always a const reference to data which cannot be changed. A 00931 normal \ref ovector_base reference does not refer to const data. 00932 00933 \comment 00934 00935 Vector objects, like \ref ovector, can be const in the usual 00936 way. However, in order to properly prevent the user from 00937 making a non-const view of a const vector, A separate 00938 view class \ref ovector_const_view was created. To demonstrate, 00939 the following code is allowed 00940 \code 00941 ovector x(2); 00942 const ovector_view y=x; 00943 \endcode 00944 However, if the original vector is <tt>const</tt>, then the 00945 \ref ovector_const_view class should be used 00946 \code 00947 const ovector x(2); 00948 const ovector_const_view y=x; 00949 \endcode 00950 00951 \endcomment 00952 00953 \subsection vmout_subsect Vector and matrix output 00954 00955 Both the \ref ovector_view_tlate and \ref omatrix_view_tlate 00956 classes have an << operator defined for each class. For 00957 writing generic vectors to a stream, you can use \ref vector_out() 00958 which is defined in \ref vector.h . Pretty matrix output is 00959 performed by global template functions \ref matrix_out(), \ref 00960 matrix_cx_out_paren(), and \ref matrix_out_paren() which are 00961 defined in \ref columnify.h . 00962 00963 \subsection tensor_subsect Tensors 00964 00965 Some preliminary support is provided for tensors of arbitrary rank 00966 and size in the class \ref tensor. Classes \ref tensor1, \ref 00967 tensor2, \ref tensor3, and \ref tensor4 are rank-specific versions 00968 for 1-, 2-, 3- and 4-rank tensors. For n-dimsional data defined on 00969 a grid, \ref tensor_grid provides a space to define a hyper-cubic 00970 grid in addition to the the tensor data. This class \ref 00971 tensor_grid also provides n-dimensional interpolation of the data 00972 defined on the specified grid. 00973 00974 \hline 00975 \section permute_section Permutations 00976 00977 Permutations are implemented through the \ref permutation class. 00978 This class is fully compatible with gsl_permutation objects since 00979 it is inherited from gsl_permutation_struct. The class also 00980 contains no new data members, so upcasting and downcasting can 00981 always be performed. It is perfectly permissible to call 00982 GSL permutation functions from \ref permutation objects by 00983 simply passing the address of the permutation, i.e. 00984 \code 00985 permutation p(4); 00986 p.init(); 00987 gsl_permutation_swap(&p,2,3); 00988 \endcode 00989 00990 The functions which apply a permutation to a 00991 user-specified vector are member template functions in the \ref 00992 permutation class (see \ref permutation::apply() ). 00993 00994 Memory allocation/deallocation between the class and the 00995 gsl_struct is compatible in many cases, but mixing these forms is 00996 strongly discouraged, i.e. avoid using 00997 <tt>gsl_permutation_alloc()</tt> on a \ref permutation object, but 00998 rather use \ref permutation::allocate() instead. The use of \ref 00999 permutation::free() is encouraged, but any remaining memory will 01000 be deallocated in the object destructor. 01001 01002 \hline 01003 \section linalg_section Linear algebra 01004 01005 There is a small set of linear algebra routines. These are not 01006 intended to be a replacement for higher performance linear algebra 01007 libraries, but are designed for use by \o2 routines so that they 01008 work for generic matrix and vector types. For vector and matrix 01009 types using \c operator[], the BLAS and linear algebra routines 01010 are inside the \ref o2scl_cblas and \ref o2scl_linalg namespaces. 01011 For vector and matrix types using \c operator(), the BLAS and 01012 linear algebra routines routines are inside the \ref 01013 o2scl_cblas_paren and \ref o2scl_linalg_paren namespaces. 01014 01015 The linear algebra classes and functions include: 01016 - Householder transformations (\ref householder.h) 01017 - Householder solver (\ref hh.h) 01018 - LU decomposition and solver (\ref lu.h) 01019 - QR decomposition and solver (\ref qr.h) 01020 - Solve tridiagonal systems (\ref tridiag.h) 01021 - Lanczos diagonalization is inside class \ref lanczos, which 01022 also can compute the eigenvalues of a tridiagonal matrix. 01023 - Givens rotations (\ref givens.h) 01024 01025 There is also a set of linear solvers for generic matrix and 01026 vector types which descend from \ref o2scl_linalg::linear_solver. 01027 These classes provide GSL-like solvers, but are generalized so 01028 that they are compatible with vector and matrix types which allow 01029 access through <tt>operator[]</tt>. 01030 01031 For users who require high-performance linear algebra, the \ref 01032 ovector and \ref omatrix objects can be used to call LAPACK 01033 routines directly, just as can be done with GSL. For an example of 01034 how to do this, see 01035 http://sourceware.org/ml/gsl-discuss/2001/msg00326.html . For 01036 uBlas users, there is an example of how to use uBlas vectors with 01037 \o2 in \ref ex_ublas_sect. Finally, there are also a couple of 01038 examples, <tt>gesvd.cpp</tt> and <tt>zheev.cpp</tt> in the 01039 <tt>src/internal</tt> directory which show how to call LAPACK with 01040 \o2 objects which may be adaptable for your platform and 01041 configuration. 01042 01043 \hline 01044 \section intp_section Interpolation 01045 01046 The classes o2scl_interp and o2scl_interp_vec allow basic 01047 interpolation, lookup, differentiation, and integration of data 01048 given in two ovectors or ovector views. In contrast to the GSL 01049 routines, data which is presented with a decreasing independent 01050 variable is handled automatically. For interpolation with 01051 arrays rather than ovectors, use array_interp or array_interp_vec. 01052 01053 For fast interpolation of arrays where the independent variable 01054 is strictly increasing and without error-checking, you can directly 01055 use the children of \ref base_interp. 01056 01057 \subsection twoint_subsect The two interpolation interfaces 01058 01059 The difference between the two classes, o2scl_interp and 01060 o2scl_interp_vec, analogous to the difference between using \c 01061 gsl_interp_eval() and \c gsl_spline_eval() in \c GSL. You can 01062 create a o2scl_interp object and use it to interpolate among any 01063 pair of chosen vectors, i.e. 01064 \code 01065 ovector x(20), y(20); 01066 // fill x and y with data 01067 o2scl_interp oi; 01068 double y_half=oi.interp(0.5,20,x,y); 01069 \endcode 01070 Alternatively, you can create a o2scl_interp_vec object which can be 01071 optimized for a pair of vectors that you specify in advance 01072 \code 01073 ovector x(20), y(20); 01074 // fill x and y with data 01075 o2scl_interp_vec oi(20,x,y); 01076 double y_half=oi.interp(0.5); 01077 \endcode 01078 01079 \subsection search_subsect Lookup and binary search 01080 01081 The class search_vec contains a searching functions for objects of 01082 type \ref ovector which are monotonic. Note that if you want to 01083 find the index of an \ref ovector where a particular value is 01084 located without any assumptions with regard to the ordering, you 01085 can use ovector::lookup() which performs an exhaustive search. 01086 01087 \subsection smintp_subsect "Smart" interpolation 01088 01089 The classes smart_interp and smart_interp_vec allow interpolation, 01090 lookup, differentiation, and integration of data which is 01091 non-monotonic or multiply-valued outside the region of interest. 01092 As with o2scl_interp above, the corresponding array versions are 01093 given in sma_interp and sma_interp_vec. 01094 01095 \subsection lowinterp_subsect Lower-level interpolation objects 01096 01097 There are five lower-level interpolation objects which 01098 specify the basic \o2 interpolation types. They are 01099 - \ref linear_interp 01100 - \ref cspline_interp 01101 - \ref cspline_peri_interp 01102 - \ref akima_interp 01103 - \ref akima_peri_interp 01104 01105 \subsection interp_mgr_subsect Interpolation manager objects 01106 01107 Many \o2 classes require the ability to create and destroy 01108 interpolation objects at will, given one of the lower-level 01109 interpolation types described above. For this purpose, 01110 interpolation managers have been created, which allow the user to 01111 provide a class with an interpolation manager associated with a 01112 low-level interpolation type (and thus free the user from having 01113 to worry about the associated memory management). The abstract 01114 base interpolation manager type is \ref base_interp_mgr, and the 01115 default interpolation manager object is \ref def_interp_mgr . The 01116 \ref def_interp_mgr class has two template arguments: the first is 01117 the vector type to be interpolated, and the second is the 01118 lower-lever interpolation type. 01119 01120 For an example usage of the default interpolation manager, see 01121 the \ref ex_contour_sect , which specifies an interpolation manager 01122 for the \ref contour class. The \ref table class also works with 01123 interpolation manager objects (see \ref table::set_interp() ). 01124 01125 \subsection highdintp_subsect Two and higher dimensional interpolation 01126 01127 Support for two-dimensional interpolation (by successive 01128 one-dimensional interpolations) is given in \ref twod_intp (see 01129 \ref tintp_section, and a simple version of n-dimensional 01130 interpolation in \ref tensor_grid. 01131 01132 \hline 01133 \section const_section Physical constants 01134 01135 The constants from GSL are reworked with the type \c const \c 01136 double and placed in namespaces called gsl_cgs, gsl_cgsm, gsl_mks, 01137 gsl_mksa, and gsl_num. Some additional constants are given in the 01138 namespace o2scl_const. Some of the numerical values have been 01139 updated from recently released data from NIST. 01140 01141 \hline 01142 \section funct_section Function Objects 01143 01144 Functions are passed to numerical routines using template-based 01145 function classes, sometimes called "functors". There are several 01146 basic kinds of function 01147 objects: 01148 - \ref funct : One function of one variable 01149 - \ref multi_funct : One function of several variables 01150 - \ref mm_funct : \c n functions of \c n variables 01151 - \ref fit_funct : One function of one variable with \c n 01152 fitting parameters 01153 - \ref ode_funct : \c n derivatives as a function of \c n function 01154 values and the value of the independent variable 01155 - \ref jac_funct : Jacobian function for solver 01156 - \ref grad_funct : Gradient function for minimizers 01157 - \ref ool_hfunct : Hessian product for constrained minimization 01158 01159 For each of these classes (except \ref funct), there is a version 01160 named \c _vfunct instead of \c _funct which is designed to be used 01161 with C-style arrays instead of vector classes. 01162 01163 The class name suffixes denote children of a generic function type 01164 which are created using different kinds of inputs: 01165 - _fptr: function pointer for a static or global function 01166 - _gsl: GSL-like function pointer 01167 - _mfptr: function pointer template for a class member function 01168 - _cmfptr: function pointer template for a class member function which 01169 is const 01170 - _strings: functions specified using strings, e.g. "x^2-2", 01171 which are in the separate \o2x package. 01172 - _noerr: (for \ref funct and multi_funct) a function which 01173 directly returns the function value rather than returning 01174 an integer error value 01175 - _nopar: (for \ref funct) a function which has no parameters and 01176 directly returns the function value. 01177 01178 See the \ref ex_fptr_sect and the \ref ex_mroot_sect which 01179 provide detailed examples of how functions can be specified 01180 to classes through these function objects. 01181 01182 There is a small overhead associated with the indirection: a "user 01183 class" accesses the function class which then calls function which 01184 was specified in the constructor of the function class. In many 01185 problems, the overhead associated with the indirection is small. 01186 Some of this overhead can always be avoided by inheriting directly 01187 from the function class and thus the user class will make a direct 01188 virtual function call. To eliminate the overhead entirely, one can 01189 specify a new type for the template parameter in the user class. 01190 01191 Note that virtual functions can be specified through 01192 this mechanism as well. For example, if \ref cern_mroot 01193 is used to solve a set of equations specified as 01194 \include fptr_doc.cpp 01195 Then the solver will solve the member function in the derived 01196 type, not the parent type. 01197 01198 Note also that providing a user access to a function object 01199 instantiatied with a protected or private member function is 01200 (basically) the same as providing them access to that function. 01201 01202 \hline 01203 \section table_section Data tables 01204 01205 The class \ref table is a container to hold and perform operations 01206 on related columns of data. It supports column operations, 01207 interpolation, column reference by either name or index, binary 01208 searching (in the case of ordered columns), sorting, and fitting 01209 two columns to a user-specified function. 01210 01211 \hline 01212 \section string_section String manipulation 01213 01214 There are a couple classes and functions to help manipulate 01215 strings of text. Conversion routines for \c std::string 01216 objects are given in \ref string_conv.h and include 01217 - \ref ptos() - pointer to string 01218 - \ref itos() - integer to string 01219 - \ref dtos() - double to string 01220 - \ref stoi() - string to integer 01221 - \ref stod() - string to double 01222 01223 See also \ref size_of_exponent(), \ref format_float, 01224 and \ref double_to_ieee_string(). 01225 01226 A class called \ref columnify converts a set of 01227 strings into nicely formatted columns by padding with the 01228 necessary amount of spaces. This class operates on string objects 01229 of type \c std::string, and also works will for formatting columns 01230 of floating-point numbers. This class is used to provide output 01231 for matrices in the functions \ref matrix_out(), \ref 01232 matrix_out_paren(), and \ref matrix_cx_out_paren(). For output 01233 of vectors, see \ref vector_out() in \ref array.h. 01234 01235 A related function, \ref screenify(), reformats a column of 01236 strings into many columns stored row-by-row in a new string 01237 array. It operates very similar to the way the classic Unix 01238 command \c ls organizes files and directories in multiple columns 01239 in order to save screen space. 01240 01241 The function \ref count_words() counts the number of "words" 01242 in a string, which are delimited by whitespace. 01243 01244 \hline 01245 \section diff_section Differentiation 01246 01247 Differentiation is performed by descendants of \ref deriv and the 01248 classes are provided. These allow one to calculate either first, 01249 second, and third derivatives. The GSL approach is used in \ref 01250 gsl_deriv, and the CERNLIB routine is used in \ref 01251 cern_deriv. Both of these compute derivatives for a function 01252 specified using a descendant of \ref funct. For functions which 01253 are tabulated over equally-spaced abscissas, the class \ref 01254 eqi_deriv is provided which applies the formulas from Abramowitz 01255 and Stegun at a specified order. 01256 01257 \b Warning: For gsl_deriv and cern_deriv, the second and third 01258 derivatives are calculated by naive repeated application of the 01259 code for the first derivative and can be particularly troublesome 01260 if the function is not sufficiently smooth. Error estimation is 01261 also incorrect for second and third derivatives. 01262 01263 \hline 01264 \section inte_section Integration 01265 01266 Integration is performed by descendants of \ref inte and is 01267 provided in the library \c o2scl_inte. 01268 01269 There are several routines for one-dimensional integration. 01270 - General integration over a finite interval: cern_adapt, 01271 cern_gauss, cern_gauss56, gsl_inte_qag, and gsl_inte_qng. 01272 - General integration from 0 to \f$ \infty \f$ : gsl_inte_qagiu 01273 - General integration from \f$ -\infty \f$ to 0: gsl_inte_qagil 01274 - General integration from \f$ -\infty \f$ to \f$ \infty \f$ : 01275 gsl_inte_qagi 01276 - General integration over a finite interval for a function with 01277 singularities: gsl_inte_qags and gsl_inte_qagp 01278 - Cauchy principal value integration over a finite interval: 01279 cern_cauchy and gsl_inte_qawc 01280 - Integration over a function weighted by \c cos(x) or \c sin(x): 01281 \ref gsl_inte_qawo_cos and \ref gsl_inte_qawo_sin 01282 - Fourier integrals: \ref gsl_inte_qawf_cos and \ref 01283 gsl_inte_qawf_sin 01284 - Integration over a weight function 01285 \f[ 01286 W(x)=(x-a)^{\alpha}(b-x)^{\beta}\log^{\mu}(x-a)\log^{\nu}(b-x) 01287 \f] 01288 is performed by \ref gsl_inte_qaws. 01289 01290 For the GSL-based integration routines, the variables inte::tolx 01291 and inte::tolf have the same role as the quantities usually 01292 denoted in the GSL integration routines by \c epsabs and \c 01293 epsrel. In particular, the integration classes attempt to ensure 01294 that 01295 \f[ 01296 |\mathrm{result}-I| \leq \mathrm{Max}(\mathrm{tolx}, 01297 \mathrm{tolf}|I|) 01298 \f] 01299 and returns an error to attempt to ensure that 01300 \f[ 01301 |\mathrm{result}-I| \leq \mathrm{abserr} \leq 01302 \mathrm{Max}(\mathrm{tolx},\mathrm{tolf}|I|) 01303 \f] 01304 where \c I is the integral to be evaluated. Even when the 01305 corresponding descendant of inte::integ() returns success, these 01306 inequalities may fail for sufficiently difficult functions. All 01307 of the GSL integration routines except for gsl_inte_qng use 01308 a workspace given in gsl_inte_table which holds the results 01309 of the various subdivisions of the original interval. The 01310 size of this workspace (and thus then number of subdivisions) can 01311 be controlled with gsl_inte_table::set_wkspace(). 01312 01313 The GSL routines were originally based on QUADPACK, which is 01314 available at http://www.netlib.org/quadpack . 01315 01316 \note The GSL integration routines can sometimes lose precision 01317 if the integrand is everywhere much smaller than unity. Some 01318 rescaling is required in these cases. 01319 01320 Multi-dimensional hypercubic integration is performed by \ref 01321 composite_inte, the sole descendant of multi_inte. \ref 01322 composite_inte allows you to specify a set of one-dimensional 01323 integration routines (objects of type \ref inte) and apply them to 01324 a multi-dimensional problem. 01325 01326 General multi-dimensional integration is performed by \ref 01327 comp_gen_inte, the sole descendant of \ref gen_inte. The user is 01328 allowed to specify a upper and lower limits which are functions of 01329 the variables for integrations which have not yet been performed, 01330 i.e. the n-dimensional integral 01331 \f[ 01332 \int_{x_0=a_0}^{x_0=b_0} f(x_0) \int_{x_1=a_1(x_0)}^{x_1=b_1(x_0)} 01333 f(x_0, x_1) ... 01334 \int_{x_{\mathrm{n}-1}=a_{\mathrm{n}-1}(x_0,x_1,..,x_{\mathrm{n}-2})}^ 01335 {x_{\mathrm{n}-1}=b_{\mathrm{n}-1}(x_0,x_1,..,x_{\mathrm{n}-2})} 01336 f(x_0,x_1,...,x_{\mathrm{n-1}})~d x_{\mathrm{n}-1}~...~d x_1~d x_0 01337 \f] 01338 Again, one specifies a set of \ref inte objects to apply to 01339 each variable to be integrated over. 01340 01341 Monte Carlo integration is also provided (see \ref mcarlo_section). 01342 01343 \hline 01344 \section poly_section Roots of Polynomials 01345 01346 Classes are provided for solving quadratic, cubic, and quartic 01347 equations as well as general polynomials. There is a standard 01348 nomenclature: classes which handle polynomials with real 01349 coefficients and real roots end with the suffix <tt>_real</tt> 01350 (quadratic_real, cubic_real and quartic_real), classes which 01351 handle real coefficients and complex roots end with the suffix 01352 <tt>_real_coeff</tt> (quadratic_real_coeff, cubic_real_coeff, 01353 quartic_real_coeff, and poly_real_coeff), and classes which handle 01354 complex polynomials with complex coefficients end with the suffix 01355 <tt>_complex</tt> (quadratic_complex, cubic_complex, 01356 quartic_complex, and poly_complex). As a reminder, complex roots 01357 may not occur in conjugate pairs if the coefficients are not 01358 real. Most of these routines will return an error if the leading 01359 coefficient is zero. 01360 01361 In the public interfaces to the polynomial solvers, the 01362 complex type <tt>std::complex<double></tt> is used. These can 01363 be converted to and from the GSL complex type using the 01364 \ref complex_to_gsl() and \ref gsl_to_complex() functions. 01365 01366 At present, the polynomial routines work with complex numbers as 01367 objects of type \c std::complex<double> and are located in library 01368 \c o2scl_other. 01369 01370 For quadratics, \ref gsl_quadratic_real_coeff is the best if the 01371 coefficients are real, while if the coefficients are complex, use 01372 \ref quadratic_std_complex. For cubics with real coefficients, 01373 \ref cern_cubic_real_coeff is the best, while if the coefficients 01374 are complex, use \ref cubic_std_complex. 01375 01376 For a quartic polynomial with real coefficients, \ref 01377 cern_quartic_real_coeff is the best, unless the coefficients of 01378 odd powers happen to be small, in which case, \ref 01379 gsl_quartic_real2 tends to work better. For quartics, generic 01380 polynomial solvers such as \ref gsl_poly_real_coeff can provide 01381 more accurate (but slower) results. If the coefficients are 01382 complex, then you can use \ref simple_quartic_complex. 01383 01384 \hline 01385 \section solve_section Equation Solving 01386 01387 \subsection onedsolve_subsect One-dimensional solvers 01388 01389 Solution of one equation in one variable is accomplished by 01390 children of the class \ref root. 01391 01392 For one-dimensional solving, use cern_root or gsl_root_brent if 01393 you have the root bracketed, or gsl_root_stef if you have the 01394 derivative available. If you have neither a bracket or a 01395 derivative, you can use cern_mroot_root. 01396 01397 The \ref root base class provides the structure for three 01398 different solving methods: 01399 - \ref root::solve() which solves a function given an initial 01400 guess \c x 01401 - \ref root::solve_bkt() which solves a function given a solution 01402 bracketed between \c x1 and \c x2. The values of the function at 01403 \c x1 and \c x2 should have different signs. 01404 - \ref root::solve_de() which solves a function given an initial 01405 guess \c x and the derivative of the function \c df. 01406 01407 If not all of these three functions are overloaded, then the 01408 source code in the \ref root base class is designed to try to 01409 automatically provide the solution using the remaining 01410 functions. Most of the one-dimensional solving routines, in their 01411 original form, are written in the second or third form above. For 01412 example, gsl_root_brent is originally a bracketing routine of the 01413 form root::solve_bkt(), but calls to either root::solve() or 01414 root::solve_de() will attempt to automatically bracket the 01415 function given the initial guess that is provided. Of 01416 course, it is frequently most efficient to use the solver in the 01417 way it was intended. 01418 01419 \subsection multisolve_subsect Multi-dimensional solvers 01420 01421 Solution of more than one equation is accomplished by descendants 01422 of the class \ref mroot . The higher-level interface is provided 01423 by the function \ref mroot::msolve() . 01424 01425 For multi-dimensional solving, you can use either cern_mroot or 01426 gsl_mroot_hybrids. While cern_mroot cannot utilize user-supplied 01427 derivatives, gsl_mroot_hybrids can use user-supplied derivative 01428 information (as in the GSL hybridsj method) using the function 01429 \ref gsl_mroot_hybrids::msolve_de() . 01430 01431 \hline 01432 \section min_section Minimization 01433 01434 One-dimensional minimization is performed by descendants of \ref 01435 minimize . There are two one-dimensional minimization algorithms, 01436 \ref cern_minimize and \ref gsl_min_brent, and they are both 01437 bracketing algorithms type where an interval and an initial guess 01438 must be provided. If only an initial guess and no bracket is 01439 given, these two classes will attempt to find a suitable bracket 01440 from the initial guess. While the \ref minimize base class is 01441 designed to allow future descendants to optionally use derivative 01442 information, this is not yet supported for any one-dimensional 01443 minimizers. 01444 01445 Multi-dimensional minimization is performed by descendants of 01446 multi_min: \ref gsl_mmin_simp2, \ref gsl_mmin_conp, \ref 01447 gsl_mmin_conf, and \ref gsl_mmin_bfgs2. (The class \ref 01448 gsl_mmin_simp2 has been updated with the new "simplex2" method 01449 from GSL-1.12. The older "simplex" method is also available in 01450 \ref gsl_mmin_simp .) The classes \ref gsl_mmin_simp and \ref 01451 gsl_mmin_simp2 do not require any derivative information. The 01452 remaining minimization classes are intended for use when the 01453 gradient of the function is available, but they can also 01454 automaticallly compute the \gradient numerically. The standard way 01455 to provide the \gradient is to use a child of \ref grad_funct (or 01456 \ref grad_vfunct). Finally, the user may specify the automatic 01457 gradient object of type \ref gradient (or \ref gradient_array) 01458 which is used by the minimizer to compute the gradient numerically 01459 when a function is not specified. 01460 01461 See an example for the usage of the multi-dimensional minimizers 01462 in \ref ex_mmin_sect . 01463 01464 Simulated annealing methods are also provided for multi-dimensional 01465 minimization (see \ref anneal_section). 01466 01467 It is important to note that not all of the minimization routines 01468 test the second derivative to ensure that it doesn't vanish to 01469 ensure that we have indeed found a true minimum. 01470 01471 The class \ref multi_min_fix provides a convenient way of 01472 fixing some of the parameters and minimizing over others, 01473 without requiring a the function interface to be rewritten. An 01474 example is given in \ref ex_mmin_fix_sect. 01475 01476 A naive way of implementing constraints is to add a function to 01477 the original which increases the value outside of the allowed 01478 region. This can be done with the functions \ref constraint() and 01479 \ref lower_bound(). There are two analogous functions, \ref 01480 cont_constraint() and \ref cont_lower_bound(), which continuous 01481 and differentiable versions. Where possible, it is better to use 01482 the constrained minimization routines described below. 01483 01484 \hline 01485 \section conmin_section Constrained Minimization 01486 01487 \o2 reimplements the Open Optimization Library (OOL) available at 01488 http://ool.sourceforge.net. The associated classes allow 01489 constrained minimization when the constraint can be expressed as a 01490 hyper-cubic constraint on all of the independent variables. The 01491 routines have been rewritten and reformatted for C++ in order to 01492 facilitate the use of member functions and user-defined vector 01493 types as arguments. The base class is ool_constr_mmin and there 01494 are two different constrained minimzation algorithms implemented 01495 in ool_mmin_pgrad, ool_mmin_spg. (The ool_mmin_gencan minimizer is 01496 not yet finished). The \o2 implementation should be essentially 01497 identical to the most recently released version of OOL. 01498 01499 The constrained minimization classes operate in a similar way 01500 to the other multi-dimensional minimization classes (which are 01501 derived from multi_min). The constraints are specified with 01502 the function 01503 \code 01504 ool_constr_mmin::set_constraints(size_t nc, vec_t &lower, 01505 vec_t &upper); 01506 \endcode 01507 and the minimization can be performed by calling either 01508 multi_min::mmin() or multi_min::mmin_de() (if the \gradient is 01509 provided by the user). The method in \ref ool_mmin_gencan requires 01510 a Hessian vector product and the user can specify this product for 01511 the minimization by using \ref ool_constr_mmin::mmin_hess(). The 01512 Hessian product function can be specified as an object of type 01513 \ref ool_hfunct or \ref ool_hvfunct in a similar way to the other 01514 function objects in \o2. 01515 01516 There are five error codes defined in \ref ool_constr_mmin 01517 which are specific to the OOL classes. 01518 01519 \hline 01520 \section mcarlo_section Monte Carlo Integration 01521 01522 Monte Carlo integration is performed by descendants of mcarlo_inte 01523 (\ref gsl_monte, \ref gsl_miser, and \ref gsl_vegas). These 01524 routines are generally superior to the direct methods for 01525 integrals over regions with large numbers of spatial dimensions. 01526 01527 \hline 01528 \section anneal_section Simulated Annealing 01529 01530 Minimization by simulated annealing is performed by descendants of 01531 sim_anneal (see \ref gsl_anneal). Because simulated annealing is 01532 particularly well-suited to parallelization, a multi-threaded 01533 minimizer analogous to \ref gsl_anneal is given in \ref anneal_mt. 01534 This header-only class uses the Boost libraries and requires it 01535 for use. 01536 01537 \hline 01538 \section fit_section Non-linear Least-Squares Fitting 01539 01540 Fitting is performed by descendants of \ref fit_base and fitting 01541 functions can be specifed using \ref fit_funct. The GSL fitting 01542 routines (scaled and unscaled) are implemented in gsl_fit. A 01543 generic fitting routine using a minimizer object specified as a 01544 child of multi_min is implemented in min_fit. When the multi_min 01545 object is (for example) a sim_anneal object, min_fit can avoid 01546 local minima which can occur when fitting noisy data. 01547 01548 \hline 01549 \section ode_section Solution of Ordinary Differential Equations 01550 01551 Classes for non-adaptive integration are provided as descendents 01552 of \ref odestep and classes for adaptive integration are 01553 descendants of \ref adapt_step. To specify a set of functions 01554 to these classes, use a child of ode_funct for a generic 01555 vector type or a child of ode_vfunct when using arrays. 01556 01557 Solution of simple initial value problems is performed by \ref 01558 ode_iv_solve. This class contains a couple different methods, 01559 depending on whether the user needs only the final value or the 01560 solution on a fixed grid. 01561 01562 The solution of boundary-value problems is based on the abstract 01563 base class \ref ode_bv_solve. At the moment, a simple shooting method 01564 is the only implementation of this base class and is given 01565 in \ref ode_bv_shoot . An experimental multishooting class is 01566 given in \ref ode_bv_multishoot . 01567 01568 An experimental application of linear solvers to solve the 01569 finite-difference equations for a boundary value problem 01570 is given in \ref ode_it_solve . 01571 01572 \hline 01573 \section rng_section Random Number Generation 01574 01575 Random number generators are descendants of \ref rnga . While the 01576 base object \ref rnga is created to allow user-defined random 01577 number generators, the only random number generator presently 01578 included are from GSL. The GSL random number generator code is 01579 reimplemented in the class \ref gsl_rnga to avoid an additional 01580 performance penalty. This may not be a truly "object-oriented" 01581 interface in that it does not use virtual functions, but it avoids 01582 any possible performance penalty. Random number generators are 01583 implemented as templates in \ref sim_anneal and \ref mcarlo_inte. 01584 In these classes, the random number generator is a template type, 01585 rather than a member data pointer, in order to ensure fast 01586 execution. 01587 01588 \hline 01589 \section tintp_section Two-dimensional Interpolation 01590 01591 Successive use of \ref smart_interp is implemented in \ref 01592 twod_intp, which is useful for interpolating data presented 01593 on a two-dimensional grid (though the spacings between grid points 01594 need not be equal). Also, one can compute \contour lines from 01595 data represented on a grid using the \ref contour class. 01596 01597 If data is arranged without a grid, then \ref planar_intp can 01598 be used. At present, the only way to compute \contour lines on 01599 data which is not defined on a grid is to use \ref planar_intp to 01600 recast the data on a grid and then use \ref contour afterwards. 01601 01602 Higher-dimensional interpolation is possible with \ref tensor_grid. 01603 01604 \hline 01605 \section gslcheb_section Chebyshev Approximation 01606 01607 A class implementing the Chebyshev approximations based on GSL is 01608 given in \ref gsl_chebapp. This class has its own copy 01609 constructor, so that Chebyshev approximations can be copied and 01610 passed as arguments to functions. 01611 01612 \hline 01613 \section unitconv_section Unit Conversions 01614 01615 There is a class which performs conversion between units specified 01616 with strings in \ref convert_units. The \ref convert_units 01617 class also uses a system call to the GNU units command (if it's 01618 installed in the path) to obtain the proper conversion factors. 01619 01620 \hline 01621 \section other_section Other classes and functions 01622 01623 The \o2 library contains several classes which are still under 01624 development and are to be considered experimental. While it is 01625 expected that the testing code associated with these classes will 01626 succeed on most any system, the interface, structure, and basic 01627 approaches used by these classes may change in future versions of 01628 \o2. 01629 01630 \b Three-dimensional data tables - \ref table3d \n 01631 \b Two-dimensional Gaussian distribution - \ref gaussian_2d \n 01632 \b Series \b acceleration - \ref gsl_series \n 01633 \b Command-line interface - \ref cli \n 01634 \b Automatic bin sizing - \ref bin_size \n 01635 \b Fourier \b transforms - \ref gsl_fft \n 01636 \b Polylogarithms - \ref polylog \n 01637 01638 \hline 01639 \section lset_section Library settings 01640 01641 There are a couple library settings which are handled by 01642 a global object \ref lib_settings of type \ref lib_settings_class. 01643 01644 There are several data files that are used by various classes in 01645 the library. The installation procedure should ensure that these 01646 files are automatically found. However, if these data files are 01647 moved after installation, then a call to 01648 lib_settings_class::set_data_dir() can adjust the library to use 01649 the new directory. It is assumed that the directory structure 01650 within the data directory has not changed. 01651 01652 \hline 01653 \section io_section Object I/O 01654 01655 The I/O portion of the library is still experimental. 01656 01657 Collections of objects can be stored in a \ref collection class, 01658 and these collections can be written to or read from text or 01659 binary files. User-defined classes may be added to the collections 01660 and may be read and written to files as long as a descendant of 01661 \ref io_base is provided. 01662 01663 Every type has an associated I/O type which is a descendant of 01664 \ref io_base. In order to perform any sort of input/output on any 01665 type, an object of the corresponding I/O type must be instantiated 01666 by the user. This is not done automatically by the library. (Since 01667 it doesn't know which objects are going to be used ahead of time, 01668 the library would have to instantiate <em>all</em> of the I/O 01669 objects, which is needlessly slow.) This makes the I/O slightly 01670 less user-friendly, but much more efficient. For convenience, each 01671 subsection of the library has a class (named with an \c _ioc 01672 suffix) which will automatically allocate all I/O types for that 01673 subsection. 01674 01675 <b>Level 1 functions:</b> Functions that input/output data from 01676 library-defined objects and internal types from files and combine 01677 these objects in collections. These are primarily member functions 01678 of the class \ref collection. 01679 01680 <b>Level 2 functions:</b> Functions which are designed to allow 01681 the user to input or output data for user-generated objects. These 01682 are primarily member functions of classes \ref cinput and \ref 01683 coutput. 01684 01685 <b>Level 3 functions:</b> Functions which allow low-level 01686 modifications on how input and output is performed. Usage of 01687 level 3 functions is not immediately recommended for the casual 01688 user. 01689 01690 <b>Level 1 usage</b>: 01691 01692 For adding an object to a collection when you have a pointer 01693 to the I/O object for the associated type: 01694 \code 01695 int collection::add(std::string name, io_base *tio, void *vec, 01696 int sz=0, int sz2=0, bool overwrt=true, bool owner=false); 01697 \endcode 01698 For adding an object to a collection otherwise: 01699 \code 01700 int collection::add(std::string name, std::string stype, 01701 void *vec, int sz=0, int sz2=0, 01702 bool overwrt=true, bool owner=false); 01703 \endcode 01704 01705 To retrieve an object as a <pre>void *</pre> from a collection 01706 use one of: 01707 \code 01708 int get(std::string tname, void *&vec); 01709 int get(std::string tname, void *&vec, int &sz); 01710 int get(std::string tname, void *&vec, int &sz, int &sz2); 01711 int get(std::string tname, std::string &stype, void *&vec); 01712 int get(std::string tname, std::string &stype, void *&vec, int &sz); 01713 int get(std::string tname, std::string &stype, void *&vec, int &sz, 01714 int &sz2); 01715 \endcode 01716 When retrieving a scalar object without error- and type-checking you 01717 can use the shorthand version: 01718 \code 01719 void *get(std::string name); 01720 \endcode 01721 01722 To output one object to a file: 01723 \code 01724 int collection::out_one(out_file_format *outs, std::string stype, 01725 std::string name, void *vp, int sz=0, int sz2=0); 01726 \endcode 01727 01728 To input one object from a file with a given type and name: 01729 \code 01730 int collection::in_one_name(in_file_format *ins, std::string stype, 01731 std::string name, void *&vp, int &sz, int &sz2); 01732 \endcode 01733 To input the first object of a given type from a file: 01734 \code 01735 int collection::in_one(in_file_format *ins, std::string stype, 01736 std::string &name, void *&vp, int &sz, int &sz2); 01737 \endcode 01738 01739 <b>Level 2 usage (string-based)</b>: 01740 01741 If you don't have a pointer to the io_base child object corresponding 01742 to the type of subobject that you are manipulating, then you can use 01743 the following functions, which take the type name as a string. 01744 01745 To input a sub-object in an io_base template for which memory 01746 has already been allocated use one of: 01747 \include io_doc1.cpp 01748 01749 To automatically allocate memory and input a sub-object of a 01750 io_base template use one of: 01751 \include io_doc2.cpp 01752 01753 To output a subobject in an io_base template use: 01754 \code 01755 int collection::object_out(std::string type, out_file_format *outs, 01756 void *op, int sz=0, int sz2=0, std::string name=""); 01757 \endcode 01758 01759 <b>Level 2 usage (with io_base pointer)</b>: 01760 01761 To input a sub-object in an io_base template for which memory 01762 has already been allocated use one of: 01763 \include io_doc4.cpp 01764 01765 To automatically allocate memory and input a sub-object of a 01766 io_base template use one of: 01767 \include io_doc5.cpp 01768 01769 To output a subobject in an io_base template use: 01770 \include io_doc6.cpp 01771 01772 To automatically allocate/deallocate memory for an object, use: 01773 \code 01774 virtual int mem_alloc(object *&op); 01775 virtual int mem_alloc_arr(object *&op, int sz); 01776 virtual int mem_alloc_2darr(object **&op, int sz, int sz2); 01777 virtual int mem_free(object *op); 01778 virtual int mem_free_arr(object *op); 01779 virtual int mem_free_2darr(object **op, int sz); 01780 \endcode 01781 01782 \subsection iotlate_subsect Usage of io_tlate 01783 01784 The functions io_tlate::input() and io_tlate::output() need to be 01785 implemented for every class has information for I/O. For 01786 subobjects of the class, cinput::object_in() and 01787 cinput::object_out() can be called to input or output the 01788 information associated with the subobject. For input, 01789 cinput::object_in_name(), cinput::object_in_mem(), and 01790 cinput::object_in_mem_name() allow the freedom to input an object 01791 with a name or with memory allocation. The function 01792 coutput::object_out_name() allows one to output an object with a 01793 name. If the class contains a pointer to the subobject, then 01794 io_base::pointer_in() or io_base::pointer_out() can be used. 01795 01796 \hline 01797 \section source_section Example source code 01798 01799 \subsection exlist_subsect Example list 01800 01801 - \ref ex_fptr_sect shows how member functions and external 01802 parameters are supplied to \o2 numerical routines. In this case, a 01803 member function representing an equation with one unknown 01804 is solved using a one-dimesional solver. 01805 - \ref ex_mroot_sect demonstrates the multidimensional function 01806 solver and several different methods of specifying the function to 01807 be solved. 01808 - \ref ex_mmin_sect 01809 - \ref ex_mmin_fix_sect 01810 - \ref ex_deriv_sect 01811 - \ref ex_inte_sect 01812 - \ref ex_ode_sect 01813 - \ref ex_anneal_sect demonstrates multidimensional 01814 minimization by simulated annealing 01815 - \ref ex_minte_sect demonstrates several ways to compute a 01816 multidimensional integral including Monte Carlo and direct methods 01817 - \ref ex_contour_sect 01818 - \ref ex_twod_intp_sect 01819 - \ref ex_ublas_sect demonstrates the use of interpolation 01820 and solver classes with ublas vectors 01821 01822 \comment 01823 \c ex_diff.cpp computes the derivative of a tabulated function 01824 01825 \c ex_table.cpp demonstrates the use of the \ref table class 01826 01827 \c ex_file.cpp demonstrates the file I/O 01828 \endcomment 01829 01830 \subsection ex_fptr_sect Function and solver example 01831 01832 \dontinclude ex_fptr.cpp 01833 \skip Example: 01834 \until End of example 01835 01836 The image below shows how the solver progresses to the 01837 solution of the example function. 01838 \image html ex_fptr.png "ex_fptr.png" 01839 \image latex ex_fptr.eps "ex_fptr.eps" width=9cm 01840 01841 \subsection ex_mroot_sect Multidimensional solver example 01842 01843 \dontinclude ex_mroot.cpp 01844 \skip Example: 01845 \until End of example 01846 01847 \subsection ex_mmin_sect Multidimensional minimizer example 01848 01849 \dontinclude ex_mmin.cpp 01850 \skip Example: 01851 \until End of example 01852 01853 \subsection ex_mmin_fix_sect Minimizer fixing variables example 01854 01855 \dontinclude ex_mmin_fix.cpp 01856 \skip Example: 01857 \until End of example 01858 01859 \subsection ex_deriv_sect Numerical differentiation example 01860 01861 \dontinclude ex_deriv.cpp 01862 \skip Example: 01863 \until End of example 01864 01865 \subsection ex_inte_sect Numerical integration example 01866 01867 \dontinclude ex_inte.cpp 01868 \skip Example: 01869 \until End of example 01870 01871 \subsection ex_ode_sect Ordinary differential equations example 01872 01873 \dontinclude ex_ode.cpp 01874 \skip Example: 01875 \until End of example 01876 01877 \subsection ex_anneal_sect Simulated annealing example 01878 01879 \dontinclude ex_anneal.cpp 01880 \skip Example: 01881 \until End of example 01882 01883 \subsection ex_minte_sect Multidimensional integration example 01884 01885 \dontinclude ex_minte.cpp 01886 \skip Example: 01887 \until End of example 01888 01889 \subsection ex_contour_sect Contour lines example 01890 01891 \dontinclude ex_contour.cpp 01892 \skip Example: 01893 \until End of example 01894 01895 \subsection ex_twod_intp_sect Two-dimensional interpolation example 01896 01897 \dontinclude ex_twod_intp.cpp 01898 \skip Example: 01899 \until End of example 01900 01901 \subsection ex_ublas_sect ublas vector example 01902 01903 \dontinclude ex_ublas.cpp 01904 \skip Example: 01905 \until End of example 01906 01907 \hline 01908 \section design_section Design Considerations 01909 01910 The design goal is to create an object-oriented computing library 01911 with classes that perform common numerical tasks. The most 01912 important principle is that the library should add functionality 01913 to the user while at the same time retaining as much freedom for 01914 the user as possible and allowing for ease of use and extensibility. 01915 To that end, 01916 - The classes which utilize user-specified functions 01917 should be able to operate on member functions without requiring 01918 a particular inheritance structure, 01919 - The interfaces ought to be generic so that the user 01920 can create new classes which perform related numerical 01921 tasks through inheritance, 01922 - The classes should not use static variables or functions 01923 - Const-correctness and type-safety should be used wherever possible, and 01924 - The design should be somewhat compatible with GSL. 01925 Also, the library provides higher-level routines 01926 for situations which do not require lower-level access. 01927 01928 <b>Header file dependencies</b> 01929 01930 For reference, it's useful to know how the top-level header files 01931 depend on each other, since it can be difficult to trace 01932 everything down. In the \c base directory, the following are some 01933 of the most "top-level" header files and their associated 01934 dependencies within \o2 (there are other dependencies on GSL and 01935 the C standard library not listed here). 01936 \code 01937 err_hnd.h : (none) 01938 lib_settings.h : (none) 01939 array.h: err_hnd.h 01940 vector.h: err_hnd.h 01941 sring_conv.h : lib_settings.h 01942 misc.h : err_hnd.h lib_settings.h 01943 test_mgr.h : string_conv.h 01944 uvector_tlate.h: err_hnd.h string_conv.h array.h vector.h 01945 ovector_tlate.h: err_hnd.h string_conv.h 01946 uvector_tlate.h array.h vector.h 01947 \endcode 01948 01949 <B>The use of templates</b> 01950 01951 Templates are used extensively, and this makes for longer 01952 compilation times so any code that can be removed conveniently 01953 from the header files should be put into source code files 01954 instead. 01955 01956 \subsection errdesign_subsect Error handling 01957 01958 <b>Thread safety</b> 01959 01960 Two approaches to thread-safe error handling which are worth 01961 comparing: the first is GSL which uses return codes and global 01962 function for an error handler, and the second is the Math/Special 01963 Functions section of Boost, which uses a separate policy type for 01964 each function. One issue is thread safety: the GSL approach is 01965 thread safe only in the sense that one can in principle use the 01966 return codes in different threads to track errors. What one cannot 01967 do in GSL is use different user-defined error handlers for 01968 different threads. The Special Functions library allows one to 01969 choose a different Policy for every special function call, and 01970 thus allows quite a bit more flexibility in designing 01971 multi-threaded error handling. 01972 01973 \subsection vecdesign_subsect Vector design 01974 01975 \o2 vector and matrix types are a hybrid approach: creating objects 01976 compatibile with GSL, while providing syntactic simplicity and 01977 object-oriented features common to C++ vector classes. In terms 01978 of their object-oriented nature, they are not as elegant as 01979 the ublas vector types from ublas, but for many applications they 01980 are also faster (and they are always at least as fast). 01981 01982 However, ensuring const-correctness makes the design a bit thorny, 01983 and this is still in progress. 01984 01985 \subsection typecase_subsect Type-casting in vector and matrix design 01986 01987 \o2 uses a GSL-like approach where viewing const double * 01988 arrays is performed by explicitly casting aways const'ness 01989 internally and then preventing the user from changing the data. 01990 01991 In GSL, the preprocessor output for \c vector/view_source.c is: 01992 \include gsl_vector_doc.cpp 01993 Note the explicit cast from const double * to double *. This 01994 is similar to what is done in \c src/base/ovector_tlate.h. 01995 01996 \subsection define_subsect Define constants and macros 01997 01998 There are a couple define constants and macros that \o2 01999 understands, they are all in upper case and begin with the prefix 02000 <tt>O2SCL_</tt>. 02001 02002 Range-checking for arrays and matrices is turned on by default, 02003 but can be turned off by defining <tt>O2SCL_NO_RANGE_CHECK</tt> 02004 during the initial configuration of the library. To see how the 02005 library was configured at runtime, use the \ref lib_settings 02006 class. 02007 02008 There are several macros for error handling defined in \ref 02009 err_hnd.h, and several for vector/matrix arithmetic in \ref 02010 vec_arith.h . 02011 02012 There is a define constant O2SCL_NO_SYSTEM_FUNC which permanently 02013 disables the shell command <tt>'!'</tt> in \ref cli (when the 02014 constant is defined, the shell command doesn't work even if 02015 \ref cli::shell_cmd_allowed is <tt>true</tt>). 02016 02017 The constant O2SCL_DATA_DIR is defined internally to provide 02018 the directory which contains the \o2 data files. After installation, 02019 this can be accessed in \ref lib_settings. 02020 02021 All of the header files have their own define constant of 02022 the form <tt>O2SCL_HEADER_FILE_NAME</tt> which ensures that 02023 the header file is only included once. 02024 02025 Finally, I sometimes comment out sections of code with 02026 \code 02027 #ifdef O2SCL_NEVER_DEFINED 02028 ... 02029 #endif 02030 \endcode 02031 This constant should not be defined by the user as it will cause 02032 compilation to fail. 02033 02034 \comment 02035 These are makefile constants not source code define constants 02036 02037 The two define constants O2SCL_PARTLIB and O2SCL_EOSLIB are used 02038 internally to control which sublibraries are compiled together 02039 with the main library (see \ref install_section ). The end-user 02040 shouldn't have to worry about these. 02041 \endcomment 02042 02043 \subsection global_subsect Global objects 02044 02045 There are three global objects that are created in 02046 libo2scl: 02047 - \ref def_err_hnd is the default error handler 02048 - \ref err_hnd is the pointer to the error handler (points to 02049 def_err_hnd by default) 02050 - \ref lib_settings to control a few library settings 02051 02052 All other global objects are to be avoided. 02053 02054 \subsection thread_subsect Thread safety 02055 02056 Most of the classes are thread-safe, meaning that two instances of 02057 the same class will not clash if their methods are called 02058 concurrently since static variables are only used for compile-time 02059 constants. However, two threads cannot, in general, safely access 02060 the same instance of a class. In this respect, \o2 is no 02061 different from GSL. 02062 02063 \subsection docdesign_subsect Documentation design 02064 02065 The commands \\comment and \\endcomment delineate comments about 02066 the documentation that are present in the header files but don't 02067 ever show up in the HTML or LaTeX documentation. 02068 02069 \subsection crightfoo_subsect Copyright notices 02070 02071 For files where it is appropriate to do so, I have followed the 02072 prescription suggested in 02073 http://lists.gnu.org/archive/html/help-gsl/2008-11/msg00017.html 02074 retaining the GSL copyright notices and putting the \o2 notices at 02075 the top. CERNLIB has no such standard, but their licensing information 02076 is outlined at 02077 http://cernlib.web.cern.ch/cernlib/conditions.html . 02078 02079 \subsection futurework_subsect Design plans 02080 02081 <b>Exception classes:</b> \n Eventually I'd like to create a sensible 02082 hierarchy of exception classes to be thrown by the handler 02083 if the user would prefer (this is already in progress). 02084 02085 <b>Boost and linear algebra:</b> \n I would like to ensure this 02086 class is compatible with boost, and start integrating things 02087 accordingly. IMHO object-oriented linear algebra is in a rather 02088 sad state at the moment. uBlas and MTL are both promising, 02089 however, and I'd like to start implementing some sort of 02090 compatibility with uBlas vectors and matrices soon. The uBlas 02091 documentation is pretty sparse, but that's the pot calling the 02092 kettle a cheap piece of metal. 02093 02094 <b>Other Improvements:</b> \n I'm particularly interested in 02095 improving the ODE and fitting classes, as well as updating the 02096 BFGS2 minimizer. Of course, more examples and better documentation 02097 are also a must. 02098 02099 <b>Algorithms to include</b> 02100 - Method of lines for PDEs 02101 - Some of the MESA interpolation routines. 02102 - C++ translation of MINUIT (done already by ROOT, but quite difficult). 02103 - Creating closed regions from contour lines (I have no idea how to 02104 do this at the moment, though I'm sure someone has solved this 02105 problem already somewhere.) 02106 02107 <b>Complex numbers</b> \n I'm not sure where to go with complex 02108 numbers. My guess is that <tt>std::complex</tt> is not 02109 significantly slower (or is faster) than <tt>gsl_complex</tt>, but 02110 it would be good to check this. Then there's the C99 standard, 02111 which is altogether different. Unfortunately the interfaces may be 02112 sufficiently different that it's impossible to make templated 02113 classes which operate on generic complex number types. 02114 02115 <b>I/O</b> \n The I/O classes are (admittedly) not structured very 02116 well, and there are a lot of things which could be done there. 02117 There are a lot of alternative solutions out there, but most of 02118 them work on binary file formats (which I find difficult to use). 02119 - On a shorter time scale I'd like to replace the inner workings 02120 of the \ref collection class with a mechanism using 02121 <tt>boost::any</tt>, which will avoid the silly <tt>void *</tt>s 02122 all over the place. This is already in progress. 02123 - On a longer time scale, I'd like to replace the file format I/O 02124 classes with something using Boost pipes. 02125 02126 \hline 02127 \section license_section License Information 02128 02129 \o2 (as well as CERNLIB and the Gnu Scientific Library (GSL)) is 02130 licensed under version 3 of the GPL as provided in the files \c 02131 COPYING and in \c doc/o2scl/extras/gpl_license.txt. After 02132 installation, it is included in the documentation in \c 02133 PREFIX/doc/extras/gpl_license.txt where the default \c PREFIX is 02134 \c /usr/local. 02135 02136 \verbinclude gpl_license.txt 02137 02138 This documentation is provided under the GNU Free Documentation 02139 License, as given below and provided in \c 02140 doc/o2scl/extras/fdl_license.txt. After installation, it is included 02141 in the documentation in \c PREFIX/doc/extras/fdl_license.txt where 02142 the default \c PREFIX is \c /usr/local. 02143 02144 \verbinclude fdl_license.txt 02145 02146 \hline 02147 \section ack_section Acknowledgements 02148 02149 I would like to thank the creators of GSL and Doxygen for their 02150 excellent work! Thanks also to Julien Garaud for contributing the 02151 ODE multishooting class and for the <tt>ex_hydrogen</tt> example. 02152 02153 \hline 02154 \section ref_section Bibliography 02155 02156 Some of the references which contain links should direct you to 02157 the work referred to directly through dx.doi.org. 02158 02159 \anchor Bader83 Bader83: 02160 G. Bader and P. Deuflhard, Numer. Math. \b 41 (1983) 373. 02161 02162 \anchor Bus75 Bus75: 02163 \htmlonly 02164 <a href="http://doi.acm.org/10.1145/355656.355659"> 02165 J.C.P. Bus and T.J. Dekker</a>, 02166 \endhtmlonly 02167 \latexonly 02168 \href{http://doi.acm.org/10.1145/355656.355659}{ 02169 J.~C.~P. Bus and T.~J. Dekker}, 02170 \endlatexonly 02171 ACM Trans. Math. Soft. \b 1 (1975) 330. 02172 02173 \anchor Cash90 Cash90: 02174 \htmlonly 02175 <a href="http://doi.acm.org/10.1145/79505.79507"> 02176 J. R. Cash and A. H. Karp</a>, 02177 \endhtmlonly 02178 \latexonly 02179 \href{http://doi.acm.org/10.1145/79505.79507}{ 02180 J. R. Cash and A. H. Karp}, 02181 \endlatexonly 02182 ACM Trans. Math. Soft. \b 16 (1990) 201. 02183 02184 \anchor Fletcher87 Fletcher87: 02185 R. Fletcher, Practical methods of optimization (John Wiley & Sons, 02186 Chichester 1987) p. 39. 02187 02188 \anchor Hairer00 Hairer00: 02189 Hairer, E., Norsett S.P., Wanner, G. Solving ordinary 02190 differential equations I, Nonstiff Problems, 2nd revised edition, 02191 Springer, 2000. 02192 02193 \anchor Krabs83 Krabs83: 02194 W. Krabs, Einfuhrung in die lieare und nictlineare Optimierung 02195 fur Ingenieure (BSB B.G. Teubner, Leipzig 1983) p. 84. 02196 02197 \anchor Lepage78 Lepage78: 02198 \htmlonly 02199 <a href="http://dx.doi.org/10.1016/0021-9991(78)90004-9"> 02200 G.P. Lepage</a>, 02201 \endhtmlonly 02202 \latexonly 02203 \href{http://dx.doi.org/10.1016/0021-9991(78)90004-9}{ 02204 G.P. Lepage}, 02205 J. Comp. Phys. \b 27 (1978) 192. 02206 02207 \anchor Lewin83 Lewin83: 02208 L. Lewin, Polylogarithms and Associated Functions (North-Holland, 02209 New York, 1983). 02210 02211 \anchor Longman58 Longman58: 02212 I.M. Longman, MTAC (later renamed Math. Comp.) \b 12 (1958) 205. 02213 02214 \anchor More79 More79: 02215 J.J. More' and M.Y. Cosnard, ACM Trans. Math. Software, 02216 \b 5 (1979) 64-85. 02217 02218 \anchor More80 More80: 02219 J.J. More' and M.Y. Cosnard, Algoright 554 BRENTM, 02220 Collected Algorithms from CACM (1980). 02221 02222 \anchor Nelder65 Nelder65: 02223 Nelder, J.A., Mead, R., Computer Journal 7 (1965) pp. 308-313. 02224 02225 \anchor Press90 Press90: 02226 W.H. Press, G.R. Farrar, "Recursive Stratified Sampling for 02227 Multidimensional Monte Carlo Integration", Computers in Physics, 02228 v4 (1990), pp. 190-195. 02229 02230 \anchor Prince81 Prince81: 02231 P.J. Prince and J.R. Dormand J. Comp. Appl. Math., 7 (1981) 67. 02232 02233 \anchor Rutishauser63 Rutishauser63: 02234 H. Rutishauser, Ausdehnung des Rombergschen Prinzips 02235 (Extension of Romberg's Principle), Numer. Math. \b 5 (1963) 48-54. 02236 02237 */ 02238 /* 02239 --------- 02240 This section has been removed from the documentation since the 02241 vector classes have been updated. 02242 --------- 02243 02244 \page cvn_page Notes on const vectors 02245 02246 The code 02247 02248 \code 02249 const ovector x; 02250 //ovector_const_base &cb=x; 02251 \endcode 02252 02253 doesn't complile, which is a bit weird (you have to use a 02254 const ovector_const_base reference instead). 02255 02256 Also, ovector_base has a function like 02257 02258 \code 02259 double &operator[](size_t i) 02260 \endcode 02261 02262 which in principle could be actually const, i.e. 02263 02264 \code 02265 double &operator[](size_t i) const 02266 \endcode 02267 02268 but isn't. Also annoying is that "const ovector_base &x" is not a 02269 perfectly generic base, since it's not a base for 02270 "ovector_const_base &y" or "const ovector_const_base &y2" (last one 02271 I'm not sure about). I had thought about making ovector_const_base a 02272 child of ovector_base instead, but that means that the interpolation 02273 routines which have ovector_const_base as a default template 02274 parameter will have to be changed. (Is this a problem?) Even worse, 02275 ovector_const_base classes could be upcasted easily to ovector_base 02276 removing const correctness. 02277 02278 Finally, the meaning of "const ovector_view" vs. "ovector_view" is 02279 kind of confusing, though maybe we just have to explain to the 02280 user that "const" means const where-ever it appears. 02281 02282 Happily, this new branch removes the need for "ovector_view_assign", 02283 in the trunk which was kind of frustrating. 02284 02285 Probably the best new solution is to rework the hierarchy altogether 02286 and admit that ovector and ovector_view aren't really related 02287 objects. There is no such thing as a generic base for both ovector 02288 and ovector_view. This of course might require a redesign of the 02289 solver classes, for example, as they treat these two classes as 02290 necessarily related (alloc_vec_t is supposed to be related to 02291 vec_t). Maybe not though, I need to think about this more. 02292 02293 */ 02294 /** 02295 \page dl_page Download O2scl 02296 02297 The current version is 0.904. The source distribution can 02298 be obtained from \n 02299 - \c http://sourceforge.net/project/showfiles.php?group_id=206918 02300 02301 You may also download \o2 from version from the Subversion 02302 repository. To obtain the bleeding-edge developer version, use 02303 something like 02304 \verbatim 02305 svn co https://o2scl.svn.sourceforge.net/svnroot/o2scl/trunk o2scl 02306 \endverbatim 02307 The checkout process is a bit more time consuming than directly 02308 downloading the source distribution. 02309 02310 Version 0.904 corresponds with version 40 in the svn repository. 02311 02312 Version 0.903 does not exist (this number was skipped to help 02313 coordinate numbering with \o2x). 02314 02315 Version 0.902 corresponds with version 37 in the svn respository. 02316 02317 Version 0.901 corresponds with version 31 in the svn respository. 02318 02319 Version 0.9 corresponds with version 26 in the svn respository. 02320 02321 Version 0.806 was a quick release to fix some issues with 02322 compilations in Mac OS X and does not correspond exactly 02323 with a version in the repository. 02324 02325 Version 0.805 corresponds with version 16 in the svn respository. 02326 02327 */ 02328 /* 02329 \htmlonly 02330 <table border=0> 02331 <tr><td>Release version number</td><td>SVN repository number</td></tr> 02332 <tr><td>0.901</td><td>31</td></tr> 02333 <tr><td>0.9</td><td>26</td></tr> 02334 <tr><td>0.806</td><td>none</td></tr> 02335 <tr><td>0.805</td><td>16</td></tr> 02336 </table> 02337 \endhtmlonly 02338 02339 */ 02340 /** 02341 \brief An empty class to add some items to the todo and bug lists 02342 02343 \todo 02344 - The o2scl-test and o2scl-examples targets require grep, awk, 02345 tail, cat, and wc. It would be good to reduce this list 02346 to ensure better compatibility. 02347 - More examples and benchmarks 02348 - There are a couple classes which Doxygen doesn't yet parse 02349 properly for the class hierarchy: gsl_root_stef, 02350 and ool_mmin_gencan. Also should double check 02351 smart_interp, etc. 02352 - Make sure abs(-double) isn't allowed by the O2scl and 02353 GSL headers with the correct flags 02354 02355 \future 02356 - Fix the PNG images so that they're smaller and repair the 02357 transparency issue (probably can be done just using different 02358 arguments to 'convert' in the pngfix target) 02359 - Make sure we have a uvector_cx_alloc, 02360 ovector_cx_const_reverse, ovector_cx_const_subvector_reverse, 02361 uvector_reverse, uvector_const_reverse, uvector_subvector_reverse, 02362 uvector_const_subvector_reverse, omatrix_cx_diag, blah_const_diag, 02363 umatrix_diag, and umatrix_cx_diag 02364 - ovector_cx_view::operator=(uvector_cx_view &) is missing 02365 - ovector_cx::operator=(uvector_cx_view &) is missing 02366 - uvector_c_view::operator+=(complex) is missing 02367 - uvector_c_view::operator-=(complex) is missing 02368 - uvector_c_view::operator*=(complex) is missing 02369 02370 \todo At present, the testing code assumes that shared libraries are 02371 installed. Can this be improved? (12/3/08 - I'm not sure what the 02372 status is on this...this should be checked.) 02373 02374 \todo Make sure default template parameters are documented 02375 in each class, and make sure default template parameters exist 02376 where they should. 02377 02378 \future Consider breaking documentation up into sections? 02379 02380 \bug 02381 - BLAS libraries not named \c libblas or \c libgslblas are not 02382 properly detected in ./configure and will have to be added 02383 manually. 02384 - The -lm flag may not be added properly by ./configure 02385 (1/4/09 - I'm not sure I remember why this is a problem, 02386 but it probably has to do with the detection of libraries 02387 which is done in configure.ac). 02388 02389 */ 02390 class other_todos_and_bugs { 02391 } 02392
Documentation generated with Doxygen and provided under the GNU Free Documentation License. See License Information for details.
Project hosting provided by
,
O2scl Sourceforge Project Page