Similarly, if you have a frequently recurring code fragment within asingle file you can define a static or static inline function for it.This is typesafe and saves writing out everything in full.@node Comments, Minimal structs, Object-orientation, Design@section CommentsFollow the GNU Coding Standards.  A relevant quote is,``Please write complete sentences and capitalize the first word.  If alower-case identifier comes at the beginning of a sentence, don'tcapitalize it!  Changing the spelling makes it a different identifier.If you don't like starting a sentence with a lower case letter, writethe sentence differently (e.g., "The identifier lower-case is ...").''@node Minimal structs, Algorithm decomposition, Comments, Design@section Minimal structsWe prefer to make structs which are @dfn{minimal}.  For example, if acertain type of problem can be solved by several classes of algorithm(e.g. with and without derivative information) it is better to makeseparate types of struct to handle those cases.  i.e. run time typeidentification is not desirable.@node Algorithm decomposition, Memory allocation and ownership, Minimal structs, Design@section Algorithm decompositionIterative algorithms should be decomposed into an INITIALIZE, ITERATE,TEST form, so that the user can control the progress of the iterationand print out intermediate results.  This is better than usingcall-backs or using flags to control whether the function prints outintermediate results.  In fact, call-backs should not be used -- if theyseem necessary then it's a sign that the algorithm should be broken downfurther into individual components so that the user has complete controlover them.  For example, when solving a differential equation the user may need tobe able to advance the solution by individual steps, while tracking arealtime process.  This is only possible if the algorithm is broken downinto step-level components.  Higher level decompositions would not givesufficient flexibility.@node Memory allocation and ownership, Memory layout, Algorithm decomposition, Design@section Memory allocation and ownershipFunctions which allocate memory on the heap should end in _alloc(e.g. gsl_foo_alloc) and be deallocated by a corresponding _free function(gsl_foo_free).Be sure to free any memory allocated by your function if you have toreturn an error in a partially initialized object.Don't allocate memory 'temporarily' inside a function and then free itbefore the function returns.  This prevents the user from controllingmemory allocation.  All memory should be allocated and freed throughseparate functions and passed around as a "workspace" argument.  Thisallows memory allocation to be factored out of tight loops.@node Memory layout, Linear Algebra Levels, Memory allocation and ownership, Design@section Memory layoutWe use flat blocks of memory to store matrices and vectors, not C-stylepointer-to-pointer arrays.  The matrices are stored in row-major order-- i.e. the column index (second index) moves continuously through memory. @node Linear Algebra Levels, Exceptions and Error handling, Memory layout, Design@section Linear Algebra LevelsFunctions using linear algebra are divided into two levels:For purely "1d" functions we use the C-style arguments (double *,stride, size) so that it is simpler to use the functions in a normal Cprogram, without needing to invoke all the gsl_vector machinery.The philosophy here is to minimize the learning curve. If someone onlyneeds to use one function, like an fft, they can do so without havingto learn about gsl_vector.  This leads to the question of why we don't do the same for matrices.In that case the argument list gets too long and confusing, with(size1, size2, tda) for each matrix and potential ambiguities over rowvs column ordering. In this case, it makes sense to use gsl_vector andgsl_matrix, which take care of this for the user.So really the library has two levels -- a lower level based on C typesfor 1d operations, and a higher level based on gsl_matrix andgsl_vector for general linear algebra.Of course, it would be possible to define a vector version of thelower level functions too. So far we have not done that because it wasnot essential -- it could be done but it is easy enough to get byusing the C arguments, by typing v->data, v->stride, v->size instead.A gsl_vector version of low-level functions would mainly be aconvenience.Please use BLAS routines internally within the library whenever possiblefor efficiency.@node Exceptions and Error handling, Persistence, Linear Algebra Levels, Design@section Exceptions and Error handlingThe basic error handling procedure is the return code (see gsl_errno.hfor a list of allowed values).  Use the GSL_ERROR macro to mark anerror.  The current definition of this macro is not ideal but it can bechanged at compile time. You should always use the GSL_ERROR macro to indicate an error, ratherthan just returning an error code. The macro allows the user to traperrors using the debugger (by setting a breakpoint on the functiongsl_error).  The only circumstances where GSL_ERROR should not be used are where thereturn value is "indicative" rather than an error -- for example, theiterative routines use the return code to indicate the success orfailure of an iteration. By the nature of an iterative algorithm"failure" (a return code of GSL_CONTINUE) is a normal occurrence andthere is no need to use GSL_ERROR there.Be sure to free any memory allocated by your function if you return anerror (in particular for errors in partially initialized objects).@node Persistence, Using Return Values, Exceptions and Error handling, Design@section PersistenceIf you make an object foo which uses blocks of memory (e.g. vector,matrix, histogram) you can provide functions for reading and writingthose blocks,@smallexampleint gsl_foo_fread (FILE * stream, gsl_foo * v);int gsl_foo_fwrite (FILE * stream, const gsl_foo * v);int gsl_foo_fscanf (FILE * stream, gsl_foo * v);int gsl_foo_fprintf (FILE * stream, const gsl_foo * v, const char *format);@end smallexample@noindentOnly dump out the blocks of memory, not any associated parameters suchas lengths.  The idea is for the user to build higher level input/outputfacilities using the functions the library provides.  The fprintf/fscanfversions should be portable between architectures, while the binaryversions should be the "raw" version of the data. Use the functions@smallexampleint gsl_block_fread (FILE * stream, gsl_block * b);int gsl_block_fwrite (FILE * stream, const gsl_block * b);int gsl_block_fscanf (FILE * stream, gsl_block * b);int gsl_block_fprintf (FILE * stream, const gsl_block * b, const char *format);@end smallexample@noindentor@smallexampleint gsl_block_raw_fread (FILE * stream, double * b, size_t n, size_t stride);int gsl_block_raw_fwrite (FILE * stream, const double * b, size_t n, size_t stride);int gsl_block_raw_fscanf (FILE * stream, double * b, size_t n, size_t stride);int gsl_block_raw_fprintf (FILE * stream, const double * b, size_t n, size_t stride, const char *format);@end smallexample@noindentto do the actual reading and writing.@node Using Return Values, Variable Names, Persistence, Design@section Using Return ValuesAlways assign a return value to a variable before using it.  This allowseasier debugging of the function, and inspection and modification of thereturn value.  If the variable is only needed temporarily then encloseit in a suitable scope.For example, instead of writing,@examplea = f(g(h(x,y)))@end example@noindentuse temporary variables to store the intermediate values,@example@{  double u = h(x,y);  double v = g(u);  a = f(v);@}@end example@noindentThese can then be inspected more easily in the debugger, and breakpointscan be placed more precisely.  The compiler will eliminate the temporaryvariables automatically when the program is compiled with optimization.@node Variable Names, Datatype widths, Using Return Values, Design@section Variable NamesTry to follow existing conventions for variable names,@table @code@item dimnumber of dimensions@item wpointer to workspace @item statepointer to state variable (use @code{s} if you need to save characters)@item resultpointer to result (output variable)@item abserrabsolute error @item relerrrelative error@item epsabsabsolute tolerance @item epsrelrelative tolerance@item sizethe size of an array or vector e.g. double array[size] @item stridethe stride of a vector@item size1the number of rows in a matrix@item size2the number of columns in a matrix@item ngeneral integer number, e.g. number of elements of array, in fft, etc@item rrandom number generator (gsl_rng)@end table@node Datatype widths, size_t, Variable Names, Design@section Datatype widthsBe aware that in ANSI C the type @code{int} is only guaranteed toprovide 16-bits. It may provide more, but is not guaranteed to.Therefore if you require 32 bits you must use @code{long int}, whichwill have 32 bits or more.  Of course, on many platforms the type@code{int} does have 32 bits instead of 16 bits but we have to code tothe ANSI standard rather than a specific platform.@node size_t, Arrays vs Pointers, Datatype widths, Design@section size_tAll objects (blocks of memory, etc) should be measured in terms of a@code{size_t} type.  Therefore any iterations (e.g. @code{for(i=0; i<N;i++)}) should also use an index of type @code{size_t}.  Don't mix @code{int} and @code{size_t}.  They are @emph{not}interchangeable.If you need to write a descending loop you have to be careful because@code{size_t} is unsigned, so instead of@examplefor (i = N - 1; i >= 0; i--) @{ ... @} /* DOESN'T WORK */@end example@noindentuse something like@examplefor (i = N; i > 0 && i--;) @{ ... @}@end example@noindentto avoid problems with wrap-around at @code{i=0}.  If you really want to avoid confusion use a separate variable to invertthe loop order,@examplefor (i = 0; i < N; i++) @{ j = N - i; ... @}@end example@node Arrays vs Pointers, Pointers, size_t, Design@section Arrays vs PointersA function can be declared with either pointer arguments or arrayarguments.  The C standard considers these to be equivalent. However, itis useful to distinguish between the case of a pointer, representing asingle object which is being modified, and an array which represents aset of objects with unit stride (that are modified or not depending onthe presence of @code{const}).  For vectors, where the stride is notrequired to be unity, the pointer form is preferred.@smallexample/* real value, set on output */int foo (double * x);                           /* real vector, modified */int foo (double * x, size_t stride, size_t n);  /* constant real vector */int foo (const double * x, size_t stride, size_t n);  /* real array, modified */int bar (double x[], size_t n);                 /* real array, not modified */int baz (const double x[], size_t n);           @end smallexample@node  Pointers, Constness, Arrays vs Pointers, Design@section PointersAvoid dereferencing pointers on the right-hand side of an expression wherepossible.  It's better to introduce a temporary variable.  This iseasier for the compiler to optimise and also more readable since itavoids confusion between the use of @code{*} for multiplication anddereferencing.@examplewhile (fabs (f) < 0.5)@{  *e = *e - 1;  f  *= 2;@}@end example@noindentis better written as,@example@{   int p = *e;  while (fabs(f) < 0.5)    @{     p--;     f *= 2;    @}  *e = p;@}@end example@node Constness, Pseudo-templates, Pointers, Design@section ConstnessUse @code{const} in function prototypes wherever an object pointed to bya pointer is constant (obviously).  For variables which are meaningfullyconstant within a function/scope use @code{const} also.  This preventsyou from accidentally modifying a variable which should be constant(e.g. length of an array, etc).  It can also help the compiler do