SET(NB_KERNEL_DIRS_TO_IGNORE_IN_DOXYGEN
"${NB_KERNEL_DIRS_TO_IGNORE_IN_DOXYGEN} \\\n ${NB_KERNEL_DIR}")
ENDFOREACH(NB_KERNEL_DIR)
+ CONFIGURE_FILE(${CMAKE_SOURCE_DIR}/Doxyfile-common.cmakein
+ ${CMAKE_BINARY_DIR}/Doxyfile-common)
CONFIGURE_FILE(${CMAKE_SOURCE_DIR}/Doxyfile.cmakein
${CMAKE_BINARY_DIR}/Doxyfile)
CONFIGURE_FILE(${CMAKE_SOURCE_DIR}/Doxyfile-lib.cmakein
--- /dev/null
+PROJECT_NAME = @CMAKE_PROJECT_NAME@
+PROJECT_NUMBER = @PROJECT_VERSION@
+OUTPUT_DIRECTORY = doxygen-doc
+INPUT = @CMAKE_SOURCE_DIR@/src \
+ @CMAKE_SOURCE_DIR@/share/template
+EXAMPLE_PATH = @CMAKE_SOURCE_DIR@
+RECURSIVE = YES
+EXCLUDE = @CMAKE_SOURCE_DIR@/src/contrib \
+ @CMAKE_SOURCE_DIR@/src/external \
+ @CMAKE_SOURCE_DIR@/src/gromacs/legacyheaders/thread_mpi/atomic \
+ @CMAKE_SOURCE_DIR@/src/gromacs/selection/scanner.cpp @NB_KERNEL_DIRS_TO_IGNORE_IN_DOXYGEN@
+EXCLUDE_SYMBOLS = YY* yy* _gmx_sel_yy*
+FULL_PATH_NAMES = YES
+STRIP_FROM_PATH = @CMAKE_SOURCE_DIR@
+STRIP_FROM_INC_PATH = @CMAKE_SOURCE_DIR@/src
+INCLUDE_PATH = @CMAKE_SOURCE_DIR@/src \
+ @CMAKE_SOURCE_DIR@/src/gromacs/legacyheaders
+HAVE_DOT = @DOXYGEN_DOT_FOUND@
+DOT_PATH = @DOXYGEN_DOT_PATH@
+
+# This is for thread_mpi to #ifdef some code out that should not be documented.
+PREDEFINED = DOXYGEN
+# This is for parser.cpp to make it produce code that Doxygen understands
+# and that does not have unnecessary function declarations.
+PREDEFINED += __STDC__ YYMALLOC=malloc YYFREE=free
+
+JAVADOC_AUTOBRIEF = YES
+BUILTIN_STL_SUPPORT = YES
+INLINE_INHERITED_MEMB = YES
+SORT_BY_SCOPE_NAME = YES
+ALPHABETICAL_INDEX = YES
+SHOW_DIRECTORIES = YES
+HTML_DYNAMIC_SECTIONS = YES
+GENERATE_LATEX = NO
-PROJECT_NAME = @CMAKE_PROJECT_NAME@
-PROJECT_NUMBER = @PROJECT_VERSION@
-OUTPUT_DIRECTORY = doxygen-doc
-INPUT = @CMAKE_SOURCE_DIR@/src \
- @CMAKE_SOURCE_DIR@/share/template
-EXAMPLE_PATH = @CMAKE_SOURCE_DIR@
-RECURSIVE = YES
-EXCLUDE = @CMAKE_SOURCE_DIR@/src/contrib \
- @CMAKE_SOURCE_DIR@/src/external @NB_KERNEL_DIRS_TO_IGNORE_IN_DOXYGEN@
-EXCLUDE_SYMBOLS = YY* yy* _gmx_sel_yy*
-FULL_PATH_NAMES = YES
-STRIP_FROM_PATH = @CMAKE_SOURCE_DIR@
-STRIP_FROM_INC_PATH = @CMAKE_SOURCE_DIR@/src
-INCLUDE_PATH = @CMAKE_SOURCE_DIR@/src \
- @CMAKE_SOURCE_DIR@/src/gromacs/legacyheaders
-HAVE_DOT = @DOXYGEN_DOT_FOUND@
-DOT_PATH = @DOXYGEN_DOT_PATH@
+@INCLUDE = Doxyfile-common
ENABLED_SECTIONS = libapi
INTERNAL_DOCS = NO
WARN_LOGFILE = doxygen-doc/doxygen-lib.log
HTML_OUTPUT = html-lib
-JAVADOC_AUTOBRIEF = YES
-BUILTIN_STL_SUPPORT = YES
-INLINE_INHERITED_MEMB = YES
-SORT_BY_SCOPE_NAME = YES
-ALPHABETICAL_INDEX = YES
-SHOW_DIRECTORIES = YES
-HTML_DYNAMIC_SECTIONS = YES
-GENERATE_LATEX = NO
-
ALIASES += inpublicapi="\ingroup group_publicapi"
ALIASES += inlibraryapi="\ingroup group_libraryapi"
ALIASES += addtopublicapi="\addtogroup group_publicapi"
-PROJECT_NAME = @CMAKE_PROJECT_NAME@
-PROJECT_NUMBER = @PROJECT_VERSION@
-OUTPUT_DIRECTORY = doxygen-doc
-INPUT = @CMAKE_SOURCE_DIR@/src \
- @CMAKE_SOURCE_DIR@/share/template
-EXAMPLE_PATH = @CMAKE_SOURCE_DIR@
-RECURSIVE = YES
-EXCLUDE = @CMAKE_SOURCE_DIR@/src/contrib \
- @CMAKE_SOURCE_DIR@/src/external @NB_KERNEL_DIRS_TO_IGNORE_IN_DOXYGEN@
-EXCLUDE_SYMBOLS = YY* yy* _gmx_sel_yy*
-FULL_PATH_NAMES = YES
-STRIP_FROM_PATH = @CMAKE_SOURCE_DIR@
-STRIP_FROM_INC_PATH = @CMAKE_SOURCE_DIR@/src
-INCLUDE_PATH = @CMAKE_SOURCE_DIR@/src \
- @CMAKE_SOURCE_DIR@/src/gromacs/legacyheaders
-HAVE_DOT = @DOXYGEN_DOT_FOUND@
-DOT_PATH = @DOXYGEN_DOT_PATH@
+@INCLUDE = Doxyfile-common
ENABLED_SECTIONS =
INTERNAL_DOCS = NO
WARN_LOGFILE = doxygen-doc/doxygen-user.log
HTML_OUTPUT = html-user
-JAVADOC_AUTOBRIEF = YES
-BUILTIN_STL_SUPPORT = YES
-INLINE_INHERITED_MEMB = YES
-SORT_BY_SCOPE_NAME = YES
-ALPHABETICAL_INDEX = YES
-SHOW_DIRECTORIES = YES
-HTML_DYNAMIC_SECTIONS = YES
-GENERATE_LATEX = NO
-
ALIASES += inpublicapi="\ingroup group_publicapi"
ALIASES += inlibraryapi="\ingroup group_libraryapi"
ALIASES += addtopublicapi="\addtogroup group_publicapi"
-PROJECT_NAME = @CMAKE_PROJECT_NAME@
-PROJECT_NUMBER = @PROJECT_VERSION@
-OUTPUT_DIRECTORY = doxygen-doc
-INPUT = @CMAKE_SOURCE_DIR@/src \
- @CMAKE_SOURCE_DIR@/share/template
-EXAMPLE_PATH = @CMAKE_SOURCE_DIR@
-RECURSIVE = YES
-EXCLUDE = @CMAKE_SOURCE_DIR@/src/contrib \
- @CMAKE_SOURCE_DIR@/src/external @NB_KERNEL_DIRS_TO_IGNORE_IN_DOXYGEN@
-EXCLUDE_SYMBOLS = YY* yy* _gmx_sel_yy*
-FULL_PATH_NAMES = YES
-STRIP_FROM_PATH = @CMAKE_SOURCE_DIR@
-STRIP_FROM_INC_PATH = @CMAKE_SOURCE_DIR@/src
-INCLUDE_PATH = @CMAKE_SOURCE_DIR@/src \
- @CMAKE_SOURCE_DIR@/src/gromacs/legacyheaders
-HAVE_DOT = @DOXYGEN_DOT_FOUND@
-DOT_PATH = @DOXYGEN_DOT_PATH@
+@INCLUDE = Doxyfile-common
MACRO_EXPANSION = YES
EXPAND_ONLY_PREDEF = YES
-PREDEFINED = F77_FUNC(name,NAME)=name
+PREDEFINED += F77_FUNC(name,NAME)=name
ENABLED_SECTIONS = libapi internal
INTERNAL_DOCS = YES
HIDE_UNDOC_CLASSES = YES
WARN_LOGFILE = doxygen-doc/doxygen.log
-JAVADOC_AUTOBRIEF = YES
-BUILTIN_STL_SUPPORT = YES
INLINE_INHERITED_MEMB = NO # Makes it easier to go through all documentation
-SORT_BY_SCOPE_NAME = YES
-ALPHABETICAL_INDEX = YES
-SHOW_DIRECTORIES = YES
-HTML_DYNAMIC_SECTIONS = YES
-GENERATE_LATEX = NO
ALIASES += inpublicapi="\ingroup group_publicapi"
ALIASES += inlibraryapi="\ingroup group_libraryapi"
#include "futil.h"
int continuing(char *s)
-/* strip trailing spaces and if s ends with a CONTINUE remove that too.
- * returns TRUE if s ends with a CONTINUE, FALSE otherwise.
- */
{
int sl;
assert(s);
char *fgets2(char *line, int n, FILE *stream)
-/* This routine reads a string from stream of max length n
- * and zero terminated, without newlines
- * line should be long enough (>= n)
- */
{
char *c;
if (fgets(line,n,stream) == NULL) {
return dest;
}
-/*!
- * \param[in] pattern Pattern to match against.
- * \param[in] str String to match.
- * \returns 0 on match, GMX_NO_WCMATCH if there is no match.
- *
- * Matches \p str against \p pattern, which may contain * and ? wildcards.
- * All other characters are matched literally.
- * Currently, it is not possible to match literal * or ?.
- */
int
gmx_wcmatch(const char *pattern, const char *str)
{
#define __has_feature(x) 0 // Compatibility with non-clang compilers.
#endif
-/** \def GMX_ATTRIBUTE_NORETURN \brief Indicate that a function is not
- * expected to return.
+/*! \def GMX_ATTRIBUTE_NORETURN
+ * \brief
+ * Indicate that a function is not expected to return.
+ *
* WARNING: In general this flag should not be used for compiler
* optimizations, since set_gmx_error_handler can be set to a
* handler which does not quit.
gmx_fft_flag flags);
-
-/*! \brief Setup a 2-dimensional complex-to-complex transform
- *
- * \param fft Pointer to opaque Gromacs FFT datatype
- * \param nx Length of transform in first dimension
- * \param ny Length of transform in second dimension
- * \param flags FFT options
- *
- * \return status - 0 or a standard error message.
- *
- * \note Since some of the libraries (e.g. MKL) store work array data in their
- * handles this datatype should only be used for one thread at a time,
- * i.e. you should create one copy per thread when executing in parallel.
- */
-int
-gmx_fft_init_2d (gmx_fft_t * fft,
- int nx,
- int ny,
- gmx_fft_flag flags);
-
-
/*! \brief Setup a 2-dimensional real-to-complex transform
*
* \param fft Pointer to opaque Gromacs FFT datatype
gmx_fft_flag flags);
-/*! \brief Setup a 3-dimensional complex-to-complex transform
- *
- * \param fft Pointer to opaque Gromacs FFT datatype
- * \param nx Length of transform in first dimension
- * \param ny Length of transform in second dimension
- * \param nz Length of transform in third dimension
- * \param flags FFT options
- *
- * \return status - 0 or a standard error message.
- *
- * \note Since some of the libraries (e.g. MKL) store work array data in their
- * handles this datatype should only be used for one thread at a time,
- * i.e. you should create one copy per thread when executing in parallel.
- */
-int
-gmx_fft_init_3d (gmx_fft_t * fft,
- int nx,
- int ny,
- int nz,
- gmx_fft_flag flags);
-
-
-/*! \brief Setup a 3-dimensional real-to-complex transform
- *
- * \param fft Pointer to opaque Gromacs FFT datatype
- * \param nx Length of transform in first dimension
- * \param ny Length of transform in second dimension
- * \param nz Length of transform in third dimension
- * \param flags FFT options
- *
- * The normal space is assumed to be real, while the values in
- * frequency space are complex.
- *
- * \return status - 0 or a standard error message.
- *
- * \note Since some of the libraries (e.g. MKL) store work array data in their
- * handles this datatype should only be used for one thread at a time,
- * i.e. you should create one copy per thread when executing in parallel.
- */
-int
-gmx_fft_init_3d_real (gmx_fft_t * fft,
- int nx,
- int ny,
- int nz,
- gmx_fft_flag flags);
-
-
-
/*! \brief Perform a 1-dimensional complex-to-complex transform
*
* Performs an instance of a transform previously initiated.
void * in_data,
void * out_data);
-
-/*! \brief Perform a 2-dimensional complex-to-complex transform
- *
- * Performs an instance of a transform previously initiated.
- *
- * \param setup Setup returned from gmx_fft_init_1d()
- * \param dir Forward or Backward
- * \param in_data Input grid data. This should be allocated with gmx_new()
- * to make it 16-byte aligned for better performance.
- * \param out_data Output grid data. This should be allocated with gmx_new()
- * to make it 16-byte aligned for better performance.
- * You can provide the same pointer for in_data and out_data
- * to perform an in-place transform.
- *
- * \return 0 on success, or an error code.
- *
- * \note Data pointers are declared as void, to avoid casting pointers
- * depending on your grid type.
- */
-int
-gmx_fft_2d (gmx_fft_t setup,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data);
-
-
/*! \brief Perform a 2-dimensional real-to-complex transform
*
* Performs an instance of a transform previously initiated.
void * in_data,
void * out_data);
-
-/*! \brief Perform a 3-dimensional complex-to-complex transform
- *
- * Performs an instance of a transform previously initiated.
- *
- * \param setup Setup returned from gmx_fft_init_1d()
- * \param dir Forward or Backward
- * \param in_data Input grid data. This should be allocated with gmx_new()
- * to make it 16-byte aligned for better performance.
- * \param out_data Output grid data. This should be allocated with gmx_new()
- * to make it 16-byte aligned for better performance.
- * You can provide the same pointer for in_data and out_data
- * to perform an in-place transform.
- *
- * \return 0 on success, or an error code.
- *
- * \note Data pointers are declared as void, to avoid casting pointers
- * depending on your grid type.
- */
-int
-gmx_fft_3d (gmx_fft_t setup,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data);
-
-
-/*! \brief Perform a 3-dimensional real-to-complex transform
- *
- * Performs an instance of a transform previously initiated.
- *
- * \param setup Setup returned from gmx_fft_init_1d_real()
- * \param dir Real-to-complex or complex-to-real
- * \param in_data Input grid data. This should be allocated with gmx_new()
- * to make it 16-byte aligned for better performance.
- * \param out_data Output grid data. This should be allocated with gmx_new()
- * to make it 16-byte aligned for better performance.
- * You can provide the same pointer for in_data and out_data
- * to perform an in-place transform.
- *
- * \return 0 on success, or an error code.
- *
- * \note If you are doing an in-place transform, the last dimension of the
- * array MUST be padded up to an even integer length so n/2 complex numbers can
- * fit. Thus, if the real grid e.g. has dimension 7*5*3, you must allocate it as
- * a 7*5*4 array, where the last element in the second dimension is padding.
- * The complex data will be written to the same array, but since that dimension
- * is 7*5*2 it will now fill the entire array. Reverse complex-to-real in-place
- * transformation will produce the same sort of padded array.
- *
- * The padding does NOT apply to out-of-place transformation. In that case the
- * input will simply be 7*5*3 of real, while the output is 7*5*2 of complex.
- *
- * \note Data pointers are declared as void, to avoid casting pointers
- * depending on transform direction.
- */
-int
-gmx_fft_3d_real (gmx_fft_t setup,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data);
-
-
/*! \brief Release an FFT setup structure
*
* Destroy setup and release all allocated memory.
int nx,
int ny);
-
-/*! \brief Transpose 2d multi-element matrix
- *
- * This routine is very similar to gmx_fft_transpose_2d(), but it
- * supports matrices with more than one data value for each position.
- * It is extremely useful when transposing the x/y dimensions of a 3d
- * matrix - in that case you just set nelem to nz, and the routine will do
- * and x/y transpose where it moves entire columns of z data
- *
- * This routines works when the matrix is non-square, i.e. nx!=ny too,
- * without allocating an entire matrix of work memory, which is important
- * for huge FFT grid.
- *
- * For performance reasons you need to provide a \a small workarray
- * with length at least 2*nelem (note that the type is char, not t_complex).
- *
- * \param in_data Input data, to be transposed
- * \param out_data Output, transposed data. If this is identical to
- * in_data, an in-place transpose is performed.
- * \param nx Number of rows before transpose
- * \param ny Number of columns before transpose
- * \param nelem Number of t_complex values in each position. If this
- * is 1 it is faster to use gmx_fft_transpose_2d() directly.
- * \param work Work array of length 2*nelem, type t_complex.
+/*! \brief Cleanup global data of FFT
*
- * \return GMX_SUCCESS, or an error code from gmx_errno.h
- */
-int
-gmx_fft_transpose_2d_nelem(t_complex * in_data,
- t_complex * out_data,
- int nx,
- int ny,
- int nelem,
- t_complex * work);
-
-
+ * Any plans are invalid after this function. Should be called
+ * after all plans have been destroyed.s
+ * */
+void gmx_fft_cleanup();
#ifdef __cplusplus
#include "gmxcomplex.h"
#include "gmx_fft.h"
+#ifdef __cplusplus
+extern "C" {
+#endif
+
typedef struct gmx_parallel_3dfft *
gmx_parallel_3dfft_t;
*
* \param pfft_setup Pointer to parallel 3dfft setup structure, previously
* allocated or with automatic storage.
- * \param ngridx Global number of grid cells in the x direction. Must be
- * divisible by the number of nodes.
- * \param ngridy Global number of grid cells in the y direction. Must be
- * divisible by the number of nodes.
- * \param ngridz Global number of grid cells in the z direction.
- * \param node2slab Node id to slab index array, can be NULL.
- * \param slab2grid_x Slab index to grid_x array (nnodes+1), can be NULL.
- * \param comm MPI communicator, must have been initialized.
+ * \param ndata Number of grid cells in each direction
+ * \param real_data Real data. Input for forward and output for backward.
+ * \param complex_data Complex data.
+ * \param comm MPI communicator for both parallelization axis.
+ * Needs to be either initialized or MPI_NULL for
+ * no parallelization in that axis.
+ * \param slab2index_major Not used
+ * \param slab2index_minor Not used
* \param bReproducible Try to avoid FFT timing optimizations and other stuff
* that could make results differ for two runs with
* identical input (reproducibility for debugging).
int
gmx_parallel_3dfft_destroy(gmx_parallel_3dfft_t pfft_setup);
-
+#ifdef __cplusplus
+}
+#endif
#endif /* _gmx_parallel_3dfft_h_ */
* And Hey:
* Gnomes, ROck Monsters And Chili Sauce
*/
-#ifndef _GMX_SORT_H_
-#define _GMX_SORT_H_
-
-/** @file gmx_sort.h
- *
- * @brief Portable implementation of thread-safe sort routines.
+/*! \internal \file
+ * \brief
+ * Portable implementation of thread-safe sort routines.
*
+ * This module provides a Gromacs version of the qsort() routine defined.
+ * It is not highly optimized, but it is thread safe, i.e. multiple threads
+ * can simultaneously call gmx_qsort() with different data.
*
- * This module provides a Gromacs version of the qsort() routine defined.
- * It is not highly optimized, but it is thread safe, i.e. multiple threads
- * can simultaneously call gmx_qsort with different data.
+ * The rational is that some implementations of qsort() are not threadsafe.
+ * For instance qsort in glibc contains a bug which makes it non-threadsafe:
+ * http://sources.redhat.com/bugzilla/show_bug.cgi?id=11655
+ * On the other hand, system qsort might be faster than our own.
*/
+#ifndef GMX_SORT_H
+#define GMX_SORT_H
#include <stdlib.h>
} /* fixes auto-indentation problems */
#endif
-
-/*
- * @param base Pointer to first element in list to sort
- * @param nmemb Number of elements in list
- * @param size Size in bytes of each element
- * @param compar Comparison function that takes two pointers to elements
- * being compared as arguments. The function should return an
- * integer less than, equal to, or greater than zero if the
- * first argument is considered to be respectively less than,
- * equal to, or greater than the second.
+/*! \brief
+ * Portable thread-safe sort routine.
+ *
+ * \param base Pointer to first element in list to sort
+ * \param nmemb Number of elements in list
+ * \param size Size in bytes of each element
+ * \param compar Comparison function that takes two pointers to elements
+ * being compared as arguments. The function should return an
+ * integer less than, equal to, or greater than zero if the
+ * first argument is considered to be respectively less than,
+ * equal to, or greater than the second.
*/
void
-gmx_qsort(void * base,
- size_t nmemb,
- size_t size,
+gmx_qsort(void * base,
+ size_t nmemb,
+ size_t size,
int (*compar)(const void *, const void *));
-#ifdef GMX_THREAD_MPI
-/* Some implementations of qsort are not threadsafe.
- * For instance qsort in glibc contains a bug which makes it non-threadsafe:
- * http://sources.redhat.com/bugzilla/show_bug.cgi?id=11655
+/*! \def qsort_threadsafe
+ * \brief
+ * Thread-safe qsort.
+ *
+ * Expands to gmx_qsort() if Gromacs is built with threading, or system qsort()
+ * otherwise.
*/
+#ifdef GMX_THREAD_MPI
#define qsort_threadsafe gmx_qsort
#else
-/* System qsort might be faster than our own */
#define qsort_threadsafe qsort
#endif
-
#ifdef __cplusplus
}
#endif
-
-#endif /* _GMX_SORT_H_ */
+#endif /* GMX_SORT_H */
typedef t_complex cvec[DIM];
-static t_complex cnul = { 0.0, 0.0 };
-
static t_complex rcmul(real r,t_complex c)
{
t_complex d;
* of the rotation output files.
* \param cr Pointer to MPI communication data.
* \param x The positions of all MD particles.
+ * \param box Simulation box, needed to make group whole.
* \param mtop Molecular topology.
* \param oenv Needed to open the rotation output xvgr file.
+ * \param bVerbose Whether to print extra status information.
* \param Flags Flags passed over from main, used to determine
* whether or not we are doing a rerun.
*/
#include <errno.h>
#include "../utility/gmx_header_config.h"
-/*#include "typedefs.h"*/
#include "types/simple.h"
/* Suppress Cygwin compiler warnings from using newlib version of
extern "C" {
#endif
+/** Continuation character. */
#define CONTINUE '\\'
+/** Comment sign to use. */
#define COMMENTSIGN ';'
+/*! \brief
+ * Strip trailing spaces and if s ends with a ::CONTINUE remove that too.
+ *
+ * \returns TRUE if s ends with a CONTINUE, FALSE otherwise.
+ */
int continuing(char *s);
+/*! \brief
+ * Reads a line from a stream.
+ *
+ * This routine reads a string from stream of max length n
+ * and zero terminated, without newlines.
+ * \p s should be long enough (>= \p n)
+ */
char *fgets2(char *s, int n, FILE *stream);
-void strip_comment (char *line);
-
-int break_line (char *line,
- char *variable,
- char *value);
+/** Remove portion of a line after a ::COMMENTSIGN. */
+void strip_comment(char *line);
-void upstring (char *str);
+/** Make a string uppercase. */
+void upstring(char *str);
-void ltrim (char *str);
+/** Remove leading whitespace from a string. */
+void ltrim(char *str);
-void rtrim (char *str);
+/** Remove trailing whitespace from a string. */
+void rtrim(char *str);
-void trim (char *str);
+/** Remove leading and trailing whitespace from a string. */
+void trim(char *str);
-void nice_header (FILE *out,const char *fn);
+/** Prints creation time stamp and user information into a file as comments. */
+void nice_header(FILE *out, const char *fn);
+/** Version of gmx_strcasecmp() that also ignores '-' and '_'. */
int gmx_strcasecmp_min(const char *str1, const char *str2);
+/** Version of gmx_strncasecmp() that also ignores '-' and '_'. */
int gmx_strncasecmp_min(const char *str1, const char *str2, int n);
-/* This funny version of strcasecmp, is not only case-insensitive,
- * but also ignores '-' and '_'.
- */
+/** Case-insensitive strcmp(). */
int gmx_strcasecmp(const char *str1, const char *str2);
+/** Case-insensitive strncmp(). */
int gmx_strncasecmp(const char *str1, const char *str2, int n);
+/** Creates a duplicate of \p src. */
char *gmx_strdup(const char *src);
+/** Duplicates first \p n characters of \p src. */
char *gmx_strndup(const char *src, int n);
-
-/** Pattern matcing with wildcards. */
-int gmx_wcmatch(const char *pattern, const char *src);
+
+/*! \brief
+ * Pattern matcing with wildcards.
+ *
+ * \param[in] pattern Pattern to match against.
+ * \param[in] str String to match.
+ * \returns 0 on match, GMX_NO_WCMATCH if there is no match.
+ *
+ * Matches \p str against \p pattern, which may contain * and ? wildcards.
+ * All other characters are matched literally.
+ * Currently, it is not possible to match literal * or ?.
+ */
+int gmx_wcmatch(const char *pattern, const char *str);
/** Return value for gmx_wcmatch() when there is no match. */
#define GMX_NO_WCMATCH 1
-
-/* this is our implementation of strsep, the thread-safe replacement for
- strtok */
+/** Our implementation of strsep, the thread-safe replacement for strtok. */
char *gmx_strsep(char **stringp, const char *delim);
-
-char *wrap_lines(const char *buf,int line_width, int indent,
- gmx_bool bIndentFirst);
-/* wraps lines at 'linewidth', indenting all following
- * lines by 'indent' spaces. A temp buffer is allocated and returned,
- * which can be disposed of if no longer needed.
- * If !bIndentFirst, then the first line will not be indented, only
+/*! \brief
+ * Wraps lines, optionally indenting lines.
+ *
+ * Wraps lines at \p linewidth, indenting all following lines by \p indent
+ * spaces. A temp buffer is allocated and returned, which can be disposed of
+ * if no longer needed.
+ * If \p bIndentFirst is FALSE, then the first line will not be indented, only
* the lines that are created due to wapping.
*/
+char *wrap_lines(const char *buf,int line_width, int indent,
+ gmx_bool bIndentFirst);
-
+/** Implementation of the well-known Perl function split. */
char **split(char sep,const char *str);
-/* Implementation of the well-known Perl function split */
+/*! \brief
+ * Convert a string to gmx_large_int_t.
+ *
+ * This method works as the standard library function strtol(), except that it
+ * does not support different bases.
+ *
+ * \attention
+ * The following differences are present from the standard behavior:
+ * - \p endptr cannot be NULL.
+ * - If an overflow occurs, returns zero and \p *endptr will equal \p str.
+ * errno is still set to ERANGE.
+ */
gmx_large_int_t str_to_large_int_t(const char *str, char **endptr);
#ifdef GMX_NATIVE_WINDOWS
#error No atomic operations implemented for this cpu/compiler combination.
#endif
+/** Indicates that no support for atomic operations is present. */
#define TMPI_NO_ATOMICS
*/
typedef struct tMPI_Atomic
{
- int value;
+ int value; /**< The atomic value. */
}
tMPI_Atomic_t;
*/
typedef struct tMPI_Atomic_ptr
{
- void* value;
+ void* value; /**< The atomic pointer value. */
}
tMPI_Atomic_ptr_t;
tMPI_Comm comm;
} tMPI_Reduce_req;
+/** Allocate data structure for asynchronous reduce. */
tMPI_Reduce_req *tMPI_Reduce_req_alloc(tMPI_Comm comm);
#if 0
-/** Execute fast a asynchronious reduce over comm.
+/** Execute fast a asynchronous reduce over comm.
Reduces array input with supplied funtion. This function may return before
the input array is ready to be written to again; to check for its completion,
Sets the number of events that had occurred during the wait in N.
\param ev The event structure to wait on.
- \ret The number of events that have occurred at function
+ \returns The number of events that have occurred at function
return time. */
int tMPI_Event_wait(tMPI_Event *ev);
Returns the total number of cores and SMT threads that can run.
- \ret The maximum number of threads that can run simulataneously. If this
- number cannot be determined for the current architecture, -1 is
- returned.
+ \returns The maximum number of threads that can run simulataneously.
+ If this number cannot be determined for the current architecture,
+ -1 is returned.
*/
int tMPI_Get_hw_nthreads(void);
equal to the number of hardware threads available, or 1 if the number
can't be determined, or if there are no atomics for this platform.
- \ret The maximum number of threads to run on.
+ \returns The maximum number of threads to run on.
*/
int tMPI_Get_recommended_nthreads(void);
Is a list with push, pop and detach operations */
typedef struct
{
- tMPI_Atomic_ptr_t head;
+ tMPI_Atomic_ptr_t head; /**< Pointer to the top stack element. */
} tMPI_Stack;
/** A single element in stack */
typedef struct tMPI_Stack_element
{
- struct tMPI_Stack_element *next; /*< pointer to the next stack element */
- void *data; /*< pointer to data */
+ struct tMPI_Stack_element *next; /**< Pointer to the next stack element. */
+ void *data; /**< Pointer to data. */
} tMPI_Stack_element;
-
-/** Lock-free double-ended queue (LIFO)
+#if 0
+/** Lock-free double-ended queue (FIFO)
Is a list with enqueue and dequeue operations */
typedef struct
/** A single element in a queue */
typedef struct tMPI_Queue_element
{
- struct tMPI_Queue_element *next,*prev; /*< pointer to the next, prev queue
- element */
- void *data; /*< pointer to data */
+ struct tMPI_Queue_element *next; /**< Pointer to the next queue element. */
+ struct tMPI_Queue_element *prev; /**< Pointer to the prev queue element. */
+ void *data; /**< Pointer to data. */
} tMPI_Queue_element;
/** Initialize a queue */
struct tMPI_List_element *next, *prev;
void *data;
} tMPI_List_element;
-
+
+/** Initialize a list */
void tMPI_List_init(tMPI_List *l);
+/** Deallocates a list */
void tMPI_List_destroy(tMPI_List *l);
tMPI_List_element* tMPI_List_first(tMPI_List *l);
void tMPI_List_insert(tMPI_List *l, tMPI_List_element *after,
tMPI_List_element *le);
void tMPI_List_remove(tMPI_List *l, tMPI_List_element *le);
+#endif
#ifdef __cplusplus
#ifndef TMPI_NUMA_MALLOC_H_
#define TMPI_NUMA_MALLOC_H_
-/*! \file numa_alloc.h
+/*! \file
\brief NUMA aware memory allocators.
*/
typedef struct
{
- tMPI_Atomic_t initialized;
- struct tMPI_Mutex* mutex;
+ tMPI_Atomic_t initialized; /*!< Whether \a mutex has been initialized. */
+ struct tMPI_Mutex* mutex; /*!< Actual mutex data structure. */
} tMPI_Thread_mutex_t;
/*! \brief Static initializer for tMPI_Thread_mutex_t
*
*/
typedef struct
{
- tMPI_Atomic_t initialized;
- struct tMPI_Thread_key *key;
+ tMPI_Atomic_t initialized; /*!< Whether \a key has been initialized. */
+ struct tMPI_Thread_key *key; /*!< Actual key data structure. */
} tMPI_Thread_key_t;
*/
typedef struct
{
- tMPI_Atomic_t once;
+ tMPI_Atomic_t once; /*!< Whether the operation has been performed. */
} tMPI_Thread_once_t;
/*! \brief Static initializer for tMPI_Thread_once_t
*
*/
typedef struct
{
- tMPI_Atomic_t initialized;
- struct tMPI_Thread_cond* condp;
+ tMPI_Atomic_t initialized; /*!< Whether \a condp has been initialized. */
+ struct tMPI_Thread_cond* condp; /*!< Actual condition variable data structure. */
} tMPI_Thread_cond_t;
/*! \brief Static initializer for tMPI_Thread_cond_t
*
*/
typedef struct
{
- tMPI_Atomic_t initialized;
- struct tMPI_Thread_barrier* barrierp;
+ tMPI_Atomic_t initialized; /*!< Whether \a barrierp has been initialized. */
+ struct tMPI_Thread_barrier* barrierp; /*!< Actual barrier data structure. */
volatile int threshold; /*!< Total number of members in barrier */
volatile int count; /*!< Remaining count before completion */
volatile int cycle; /*!< Alternating 0/1 to indicate round */
\param message format string for error message.
*/
void tMPI_Fatal_error(const char *file, int line, const char *message, ...);
+/** Convenience macro for the first two arguments to tMPI_Fatal_error(). */
#define TMPI_FARGS __FILE__,__LINE__
Returns the total number of cores and SMT threads that can run.
- \ret The maximum number of threads that can run simulataneously. If this
- number cannot be determined for the current architecture, 0 is
- returned. */
+ \returns The maximum number of threads that can run simulataneously.
+ If this number cannot be determined for the current architecture,
+ 0 is returned.
+ */
int tMPI_Thread_get_hw_number(void);
/** Wait until some of several messages are transferred. Waits until at least
one message is transferred.
- \param[in] count The number of requests
+ \param[in] incount The number of requests
\param[in,out] array_of_requests List of count requests obtained with
tMPI_Isend()/tMPI_Irecv().
\param[out] outcount Number of completed requests
/** Test whether some of several messages are transferred.
- \param[in] count The number of requests
+ \param[in] incount The number of requests
\param[in,out] array_of_requests List of count requests obtained with
tMPI_Isend()/tMPI_Irecv().
\param[out] outcount Number of completed requests
-file(GLOB MDLIB_SOURCES *.c)
+file(GLOB MDLIB_SOURCES *.c *.cpp)
if(GMX_FFT_FFTPACK)
list(APPEND MDLIB_SOURCES ${CMAKE_SOURCE_DIR}/src/external/fftpack/fftpack.c)
endif()
#include <config.h>
#endif
+#include <algorithm>
#include <stdio.h>
#include <stdlib.h>
#ifdef GMX_FFT_FFTW3
-#ifdef GMX_THREAD_MPI
+#include "thread_mpi/mutex.h"
+#include "gromacs/utility/exceptions.h"
/* none of the fftw3 calls, except execute(), are thread-safe, so
we need to serialize them with this mutex. */
-static tMPI_Thread_mutex_t big_fftw_mutex=TMPI_THREAD_MUTEX_INITIALIZER;
-
-#define FFTW_LOCK tMPI_Thread_mutex_lock(&big_fftw_mutex)
-#define FFTW_UNLOCK tMPI_Thread_mutex_unlock(&big_fftw_mutex)
-#else /* GMX_THREAD_MPI */
-#define FFTW_LOCK
-#define FFTW_UNLOCK
-#endif /* GMX_THREAD_MPI */
+static tMPI::mutex big_fftw_mutex;
+#define FFTW_LOCK try { big_fftw_mutex.lock(); } GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
+#define FFTW_UNLOCK try { big_fftw_mutex.unlock(); } GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
#endif /* GMX_FFT_FFTW3 */
-static double fft5d_fmax(double a, double b){
- return (a>b)?a:b;
-}
-
/* largest factor smaller than sqrt */
static int lfactor(int z) {
int i;
- for (i=sqrt(z);;i--)
+ for (i=static_cast<int>(sqrt(static_cast<double>(z)));;i--)
if (z%i==0) return i;
return 1;
}
static int lpfactor(int z) {
int f = l2factor(z);
if (f==1) return z;
- return fft5d_fmax(lpfactor(f),lpfactor(z/f));
+ return std::max(lpfactor(f),lpfactor(z/f));
}
#ifndef GMX_MPI
*/
/* int lsize = fmax(N[0]*M[0]*K[0]*nP[0],N[1]*M[1]*K[1]*nP[1]); */
- lsize = fft5d_fmax(N[0]*M[0]*K[0]*nP[0],fft5d_fmax(N[1]*M[1]*K[1]*nP[1],C[2]*M[2]*K[2]));
+ lsize = std::max(N[0]*M[0]*K[0]*nP[0],std::max(N[1]*M[1]*K[1]*nP[1],C[2]*M[2]*K[2]));
/* int lsize = fmax(C[0]*M[0]*K[0],fmax(C[1]*M[1]*K[1],C[2]*M[2]*K[2])); */
if (!(flags&FFT5D_NOMALLOC)) {
snew_aligned(lin, lsize, 32);
}
static void print_localdata(const t_complex* lin, const char* txt, int s, fft5d_plan plan) {
- int x,y,z,l,t;
+ int x,y,z,l;
int *coor = plan->coor;
int xs[3],xl[3],xc[3],NG[3];
int ll=(plan->flags&FFT5D_REALCOMPLEX)?1:2;
#ifdef GMX_MPI
MPI_Comm *cart=plan->cart;
#endif
-
+#ifdef NOGMX
double time_fft=0,time_local=0,time_mpi[2]={0},time=0;
+#endif
int *N=plan->N,*M=plan->M,*K=plan->K,*pN=plan->pN,*pM=plan->pM,*pK=plan->pK,
*C=plan->C,*P=plan->P,**iNin=plan->iNin,**oNin=plan->oNin,**iNout=plan->iNout,**oNout=plan->oNout;
int s=0,tstart,tend,bParallelDim;
{
FFTW(destroy_plan)(plan->mpip[s]);
}
+#endif /* FFT5D_MPI_TRANSPOS */
if (plan->p3d)
{
FFTW(destroy_plan)(plan->p3d);
}
-#endif /* FFT5D_MPI_TRANSPOS */
+ FFTW_UNLOCK;
#endif /* GMX_FFT_FFTW3 */
if (!(plan->flags&FFT5D_NOMALLOC))
#include <config.h>
#endif
+#ifdef __cplusplus
+extern "C" {
+#endif
+
#ifdef NOGMX
/*#define GMX_MPI*/
/*#define GMX_FFT_FFTW3*/
void fft5d_destroy(fft5d_plan plan);
fft5d_plan fft5d_plan_3d_cart(int N, int M, int K, MPI_Comm comm, int P0, int flags, t_complex** lin, t_complex** lin2, t_complex** lout2, t_complex** lout3, int nthreads);
void fft5d_compare_data(const t_complex* lin, const t_complex* in, fft5d_plan plan, int bothLocal, int normarlize);
+
+#ifdef __cplusplus
+}
+#endif
#endif /*FFTLIB_H_*/
}
}
-#endif
+#endif //not GMX_FFT_FFTW3
int gmx_fft_transpose_2d(t_complex * in_data,
t_complex * out_data,
-/* Same as the above, but assume each (x,y) position holds
- * nelem complex numbers.
- * This is used when transposing the x/y dimensions of a
- * 3D matrix; set nelem to nz in this case.
- */
-int
-gmx_fft_transpose_2d_nelem(t_complex * in_data,
- t_complex * out_data,
- int nx,
- int ny,
- int nelem,
- t_complex * work)
-{
- int i,j,k,im,n,ncount,done1,done2;
- int i1,i1c,i2,i2c,kmi,max,ncpy;
-
- t_complex *tmp1,*tmp2,*tmp3;
- t_complex *data;
-
- /* Use 500 bytes on stack to indicate moves.
- * This is just for optimization, it does not limit any dimensions.
- */
- char move[500];
- int nmove = 500;
-
- ncpy = nelem*sizeof(t_complex);
-
- if(nx<2 || ny<2)
- {
- if(in_data != out_data)
- {
- memcpy(out_data,in_data,ncpy*nx*ny);
- }
- return 0;
- }
-
- /* Out-of-place transposes are easy */
- if(in_data != out_data)
- {
- for(i=0;i<nx;i++)
- {
- for(j=0;j<ny;j++)
- {
- memcpy(out_data + (j*nx+i)*nelem,
- in_data + (i*ny+j)*nelem,
- ncpy);
- }
- }
- return 0;
- }
-
-
- /* In-place transform. in_data=out_data=data */
- data = in_data;
-
-
- /* Check the work array */
- if(work == NULL)
- {
- gmx_fatal(FARGS,
- "No work array provided to gmx_fft_transpose_2d_nelem().");
- return EINVAL;
- }
-
- tmp1 = work;
- tmp2 = work + nelem;
-
- if(nx==ny)
- {
- /* trivial case, just swap elements */
- for(i=0;i<nx;i++)
- {
- for(j=i+1;j<nx;j++)
- {
- memcpy(tmp1,data+(i*nx+j)*nelem,ncpy);
- memcpy(data+(i*nx+j)*nelem,data+(j*nx+i)*nelem,ncpy);
- memcpy(data+(j*nx+i)*nelem,tmp1,ncpy);
- }
- }
- return 0;
- }
-
- for(i=0;i<nmove;i++)
- {
- move[i]=0;
- }
-
- ncount = 2;
-
- if(nx>2 && ny>2)
- {
- i = nx-1;
- j = ny-1;
- do
- {
- k = i % j;
- i = j;
- j = k;
- }
- while(k);
- ncount += i-1;
- }
-
- n = nx*ny;
- k = n - 1;
- i = 1;
- im = ny;
-
- done1=0;
- do
- {
- i1 = i;
- kmi = k-i;
- i1c=kmi;
- memcpy(tmp1,data+i1*nelem,ncpy);
- memcpy(tmp2,data+i1c*nelem,ncpy);
-
- done2=0;
- do
- {
- i2 = ny*i1-k*(i1/nx);
- i2c = k-i2;
- if(i1<nmove)
- {
- move[i1]=1;
- }
- if(i1c<nmove)
- {
- move[i1c]=1;
- }
- ncount += 2;
- if(i2==i)
- {
- done2 = 1;
- }
- else if(i2 == kmi)
- {
- /* Swapping pointers instead of copying data */
- tmp3=tmp1;
- tmp1=tmp2;
- tmp2=tmp3;
- done2=1;
- }
- else
- {
- memcpy(data+i1*nelem,data+i2*nelem,ncpy);
- memcpy(data+i1c*nelem,data+i2c*nelem,ncpy);
- i1=i2;
- i1c = i2c;
- }
- }
- while(!done2);
-
- memcpy(data+i1*nelem,tmp1,ncpy);
- memcpy(data+i1c*nelem,tmp2,ncpy);
-
- if(ncount>=n)
- {
- done1=1;
- }
- else
- {
- done2=0;
- do
- {
- max=k-i;
- i++;
- im+=ny;
- if(im>k)
- {
- im-=k;
- }
- i2=im;
- if(i!=i2)
- {
- if(i>=nmove)
- {
- while(i2>i && i2<max)
- {
- i1=i2;
- i2=ny*i1-k*(i1/nx);
- }
- if(i2==i)
- {
- done2=1;
- }
- }
- else if(!move[i])
- {
- done2=1;
- }
- }
- }
- while(!done2);
- }
- }
- while(!done1);
-
- return 0;
-}
-
-
-
+++ /dev/null
-/* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
- *
- *
- * Gromacs 4.0 Copyright (c) 1991-2003
- * David van der Spoel, Erik Lindahl, University of Groningen.
- *
- * This program is free software; you can redistribute it and/or
- * modify it under the terms of the GNU General Public License
- * as published by the Free Software Foundation; either version 2
- * of the License, or (at your option) any later version.
- *
- * To help us fund GROMACS development, we humbly ask that you cite
- * the research papers on the package. Check out http://www.gromacs.org
- *
- * And Hey:
- * Grim Ragnarok Overthrows Midgard Amidst Cursing Silence
- */
-#ifdef HAVE_CONFIG_H
-#include <config.h>
-#endif
-
-#ifdef GMX_FFT_ACML
-
-#include <errno.h>
-#include <stdlib.h>
-#include <memory.h>
-#include <math.h>
-#include <assert.h>
-
-#include "acml.h"
-
-
-#include "gmx_fft.h"
-#include "gmx_fatal.h"
-
-/* Since c has no support of templates, we get around the type restrictions
- with macros. Our custom macro names obscure the float vs double interfaces
- */
-#ifdef GMX_DOUBLE
-#define ACML_FFT1DX zfft1dx
-#define ACML_FFT2DX zfft2dx
-#define ACML_FFT3DX zfft3dy
-#define ACML_FFT1DMX zfft1mx
-
-#define ACML_RCFFT1D dzfft
-#define ACML_CRFFT1D zdfft
-#define ACML_RCFFT1DM dzfftm
-#define ACML_CRFFT1DM zdfftm
-
-#define acmlComplex doublecomplex
-#else
-
-#define ACML_FFT1DX cfft1dx
-#define ACML_FFT2DX cfft2dx
-#define ACML_FFT3DX cfft3dy
-#define ACML_FFT1DMX cfft1mx
-
-#define ACML_RCFFT1D scfft
-#define ACML_CRFFT1D csfft
-#define ACML_RCFFT1DM scfftm
-#define ACML_CRFFT1DM csfftm
-
-#define acmlComplex complex
-#endif
-
-void cPrintArray( t_complex* arr, int outer, int mid, int inner )
-{
- int i, j, k;
-
- printf("\n");
- for( i = 0; i < outer; ++i )
- {
- for( j = 0; j < mid; ++j )
- {
- for( k = 0; k < inner; ++k )
- {
- printf("%f,%f ", arr[i*inner*mid+j*inner+k].re, arr[i*inner*mid+j*inner+k].im );
- }
- printf("\n");
- }
- printf("\n");
- }
-
- return;
-}
-
-void rPrintArray( real* arr, int outer, int mid, int inner )
-{
- int i, j, k;
-
- printf("\n");
- for( i = 0; i < outer; ++i )
- {
- for( j = 0; j < mid; ++j )
- {
- for( k = 0; k < inner; ++k )
- {
- printf("%f ", arr[i*inner*mid+j*inner+k] );
- }
- printf("\n");
- }
- printf("\n");
- }
-
- return;
-}
-
-/* Contents of the AMD ACML FFT fft datatype.
- *
- * Note that this is one of several possible implementations of gmx_fft_t.
- *
- * The ACML _API_ supports 1D,2D, and 3D transforms, including real-to-complex.
- * Unfortunately the actual library implementation does not support 2D&3D real
- * transforms as of version 4.2 In addition, ACML outputs their hermitian data
- * in such a form that it can't be fed straight back into their complex API's.
- *
- * To work around this we roll our own 2D and 3D real-to-complex transforms,
- * using separate X/Y/Z handles defined to perform (ny*nz), (nx*nz), and
- * (nx*ny) transforms at once when necessary. To perform strided multiple
- * transforms out-of-place (i.e., without padding in the last dimension)
- * on the fly we also need to separate the forward and backward
- * handles for real-to-complex/complex-to-real data permutation.
- *
- * So, the handles are enumerated as follows:
- *
- * 1D FFT (real too): Index 0 is the handle for the entire FFT
- * 2D complex FFT: Index 0 is the handle for the entire FFT
- * 3D complex FFT: Index 0 is the handle for the entire FFT
- * 2D, real FFT: 0=FFTx, 1=FFTrc, 2=FFTcr handle
- * 3D, real FFT: 0=FFTx, 1=FFTy, 2=FFTrc, 3=FFTcr handle
- *
- * AMD people reading this: Learn from FFTW what a good interface looks like :-)
- *
- */
-struct gmx_fft
-{
- /* Work arrays;
- * c2c = index0
- * r2c = index0, c2r = index1
- */
- void* comm[4];
- void* realScratch;
-
- int ndim; /**< Number of dimensions in FFT */
- int nx; /**< Length of X transform */
- int ny; /**< Length of Y transform */
- int nz; /**< Length of Z transform */
- int real_fft; /**< 1 if real FFT, otherwise 0 */
-};
-
-/* ACML's real FFT support leaves the data in a swizzled 'hermitian' format. The
- problem with this is that you can't just feed this data back into ACML :( A
- pre-processing step is required to transform the hermitian swizzled complex
- values into unswizzled complex values
- -This function assumes that the complex conjugates are not expanded; the
- calling function should not assume that they exist.
- -This function pads all real numbers with 0's for imaginary components.
- */
-void hermitianUnPacking( float scale, int row, int col, void* src, int srcPitch, void* dst )
-{
- int mid = col/2;
- int loopCount = (col-1)/2;
- int hermLength = mid + 1;
-
- int nFFT;
-
- /* Currently this function is not written to handle aliased input/output pointers */
- assert( src != dst );
-
- /* These two pointers could be aliased on top of each other; be careful */
- real* realData = (real*)src;
- t_complex* hermData = (t_complex*)dst;
-
- /* We have to expand the data from real to complex, which will clobber data if you start
- * from the beginning. For this reason, we start at the last array and work backwards,
- * however, we skip the very first row (nFFT == 0) because we still have data collision
- * on that row. We copy the first row into a temporary buffer.
- */
- for( nFFT = row-1; nFFT >= 0; --nFFT )
- {
- realData = (real*)src + (nFFT*srcPitch);
- hermData = (t_complex*)dst + (nFFT*hermLength);
-
- /* The first complex number is real valued */
- t_complex tmp = { realData[0]*scale, 0 };
- hermData[0] = tmp;
-
- int i;
- for(i = 1; i <= loopCount; ++i)
- {
- t_complex tmp = { realData[i]*scale, realData[col-i]*scale };
- hermData[i] = tmp;
- }
-
- /* A little cleanup for even lengths */
- if( loopCount != mid )
- {
- /* The n/2 number is real valued */
- t_complex tmp = { realData[mid]*scale, 0 };
- hermData[mid] = tmp;
- }
- }
-
- return;
-}
-
-/* ACML's real FFT support requires the data in a swizzled 'hermitian' format.
- This is a non-standard format, and requires us to manually swizzle the data :(
- A pre-processing step is required to transform the complex values into the
- swizzled hermitian format
- */
-void hermitianPacking( float scale, int row, int col, void* src, void* dst, int dstPitch )
-{
- int mid = col/2;
- int loopCount = (col-1)/2;
- int hermLength = mid + 1;
-
- int nFFT;
-
- real* realData;
- t_complex* hermData;
-
- /* Currently this function is not written to handle aliased input/output pointers */
- assert( src != dst );
-
- for( nFFT = 0; nFFT < row; ++nFFT )
- {
- realData = (real*)dst + (nFFT*dstPitch);
- hermData = (t_complex*)src + (nFFT*hermLength);
-
- /* The first complex number is real valued */
- realData[0] = hermData[0].re * scale;
-
- /* ACML's complex to real function is documented as only handling 'forward'
- * transforms, which mean we have to manually conjugate the data before doing
- * the backwards transform. */
- int i;
- for(i = 1; i <= loopCount; ++i)
- {
- realData[i] = hermData[i].re * scale;
- realData[col-i] = -hermData[i].im * scale;
- }
-
- /* A little cleanup for even lengths */
- if( loopCount != mid )
- {
- /* The n/2 number is real valued */
- realData[mid] = hermData[mid].re * scale;
- }
-
- }
-
- return;
-}
-
-/* Gromacs expects results that are returned to be in a packed form; No
- complex conjugate's. This routine takes a 2D array of complex data,
- and compresses it to fit into 'real' form by removing all the complex
- conjugates and any known imaginary 0's. Data is expected in row-major
- form, with ny describing the length of a row.
- */
-void compressConj2D( int nX, int nY, void* src, void* dst )
-{
- /* These two pointers are aliased on top of each other; be careful */
- real* realData = (real*)dst;
- t_complex* complexData;
-
- int halfRows = nX/2 + 1;
- int halfCols = nY/2 + 1;
- int rowOdd = nX & 0x1;
- int colOdd = nY & 0x1;
-
-
- int nRow;
- for( nRow = 0; nRow < nX; ++nRow )
- {
- complexData = (t_complex*)src + (nRow*halfCols);
-
- int nCol;
- for( nCol = 0; nCol < halfCols; ++nCol )
- {
- if(( nRow == 0) && (nCol == 0 ))
- {
- /* The complex number is real valued */
- realData[0] = complexData[nCol].re;
- realData += 1;
- }
- /* Last column is real if even column length */
- else if( (nRow == 0) && ( !colOdd && ( nCol == (halfCols-1) ) ) )
- {
- /* The complex number is real valued */
- realData[0] = complexData[nCol].re;
- realData += 1;
- }
- /* The middle value of the first column is real if we have an even number of rows and
- * the last column of middle row is real if we have an even number of rows and columns */
- else if( !rowOdd && ( ( nRow == (halfRows-1) ) && ( ( nCol == 0 ) || ( !colOdd && ( nCol == (halfCols-1) ) ) ) ) )
- {
- /* The complex number is real valued */
- realData[0] = complexData[nCol].re;
- realData += 1;
- }
- else if( (nCol == 0) || ( !colOdd && ( nCol == (halfCols-1) ) ) )
- {
- /* Last half of real columns are conjugates */
- if( nRow < halfRows )
- {
- /* Copy a complex number, which is two reals */
- realData[0] = complexData[nCol].re;
- realData[1] = complexData[nCol].im;
- realData += 2;
- }
- }
- else
- {
- /* Copy a complex number, which is two reals */
- realData[0] = complexData[nCol].re;
- realData[1] = complexData[nCol].im;
- realData += 2;
- }
- }
- }
-
- return;
-}
-
-/* Gromacs expects results that are returned to be in a packed form; No
- complex conjugate's. This routine takes a 2D array of real data,
- and uncompresses it to fit into 'complex' form by creating all the complex
- conjugates and any known imaginary 0's. Data is expected in row-major
- form, with ny describing the length of a row.
- */
-void unCompressConj2D( int nX, int nY, void* src, void* dst )
-{
- /* These two pointers are aliased on top of each other; be careful */
- real* realData = (real*)src;
- t_complex* complexData = (t_complex*)dst;
-
- int halfRows = nX/2 + 1;
- int halfCols = nY/2 + 1;
- int rowOdd = nX & 0x1;
- int colOdd = nY & 0x1;
-
- int nRow;
- for( nRow = 0; nRow < nX; ++nRow )
- {
- int nCol;
- for( nCol = 0; nCol < halfCols; ++nCol )
- {
- int iIndex = nRow*halfCols + nCol;
-
- if(( nRow == 0) && (nCol == 0 ))
- {
- /* The complex number is real valued */
- complexData[iIndex].re = realData[0];
- complexData[iIndex].im = 0;
- realData += 1;
- }
- /* Last column is real if even column length */
- else if( (nRow == 0) && ( !colOdd && ( nCol == (halfCols-1) ) ) )
- {
- /* The complex number is real valued */
- complexData[iIndex].re = realData[0];
- complexData[iIndex].im = 0;
- realData += 1;
- }
- /* The middle value of the first column is real if we have an even number of rows and
- * the last column of middle row is real if we have an even number of rows and columns */
- else if( !rowOdd && ( ( nRow == (halfRows-1) ) && ( ( nCol == 0 ) || ( !colOdd && ( nCol == (halfCols-1) ) ) ) ) )
- {
- /* The complex number is real valued */
- complexData[iIndex].re = realData[0];
- complexData[iIndex].im = 0;
- realData += 1;
- }
- else if( (nCol == 0) || (!colOdd && ( nCol == (halfCols-1) ) ) )
- {
- /* Last half of real columns are conjugates */
- if( nRow < halfRows )
- {
- /* Copy a complex number, which is two reals */
- complexData[iIndex].re = realData[0];
- complexData[iIndex].im = realData[1];
- realData += 2;
- }
- else
- {
- int oIndex = (nX-nRow)*halfCols + nCol;
- complexData[iIndex].re = complexData[oIndex].re;
- complexData[iIndex].im = -complexData[oIndex].im;
- }
- }
- else
- {
- /* Copy a complex number, which is two reals */
- complexData[iIndex].re = realData[0];
- complexData[iIndex].im = realData[1];
- realData += 2;
- }
- }
- }
-
- return;
-}
-
-/* Support routine for the 1D array tranforms. ACML does not support a MODE
- * flag on the real/complex interface. It assumes a 'forward' transform both
- * directions, so that requires a manual conjugation of the imaginary comps. */
-void negateConj( int len, void* src )
-{
- real* imag = (real*)src;
- int half = len/2 + 1;
- int i;
-
- for( i = half; i < len; ++i)
- {
- imag[i] = -imag[i];
- }
-
- return;
-}
-
-int
-gmx_fft_init_1d(gmx_fft_t * pfft,
- int nx,
- enum gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- int info = 0;
- acmlComplex* comm = NULL;
- int commSize = 0;
-
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- if( (fft = malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- // Single precision requires 5*nx+100
- // Double precision requires 3*nx+100
- if( sizeof( acmlComplex ) == 16 )
- commSize = (3*nx+100)*sizeof( acmlComplex );
- else
- commSize = (5*nx+100)*sizeof( acmlComplex );
-
- // Allocate communication work array
- if( (comm = (acmlComplex*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- // Initialize communication work array
- ACML_FFT1DX( 100, 1.0f, TRUE, nx, NULL, 1, NULL, 1, comm, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- fft->ndim = 1;
- fft->nx = nx;
- fft->real_fft = 0;
- fft->comm[0] = comm;
-
- *pfft = fft;
- return 0;
-}
-
-int
-gmx_fft_init_1d_real(gmx_fft_t * pfft,
- int nx,
- enum gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- int info = 0;
- real* commRC = NULL;
- real* commCR = NULL;
- int commSize = 0;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- if( (fft = malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
-
- commSize = (3*nx+100)*sizeof( float );
-
- // Allocate communication work array, r2c
- if( (commRC = (real*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- // Allocate communication work array, c2r
- if( (commCR = (real*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- // Initialize communication work array
- ACML_RCFFT1D( 100, nx, NULL, (real*)commRC, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- // Initialize communication work array
- ACML_CRFFT1D( 100, nx, NULL, (real*)commCR, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- /* Allocate scratch work array that ACML uses to splat from hermitian complex format to
- * full complex format */
- if( (fft->realScratch = (acmlComplex*)malloc( nx*sizeof( acmlComplex ) ) ) == NULL )
- {
- return ENOMEM;
- }
-
- fft->ndim = 1;
- fft->nx = nx;
- fft->real_fft = 1;
- fft->comm[0] = commRC;
- fft->comm[1] = commCR;
-
- *pfft = fft;
- return 0;
-}
-
-int
-gmx_fft_init_2d(gmx_fft_t * pfft,
- int nx,
- int ny,
- enum gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- int info = 0;
- acmlComplex* comm = NULL;
- int commSize = 0;
-
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- if( (fft = malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- // Single precision requires nx*ny+5*(nx+ny)
- // Double precision requires nx*ny+3*(nx+ny)
- if( sizeof( acmlComplex ) == 16 )
- commSize = (nx*ny+3*(nx+ny)+200)*sizeof( acmlComplex );
- else
- commSize = (nx*ny+5*(nx+ny)+200)*sizeof( acmlComplex );
-
- // Allocate communication work array
- if( (comm = (acmlComplex*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- // Initialize communication work array
- ACML_FFT2DX( 100, 1.0f, TRUE, TRUE, nx, ny, NULL, 1, nx, NULL, 1, nx, (acmlComplex*)comm, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- fft->ndim = 2;
- fft->nx = nx;
- fft->ny = ny;
- fft->real_fft = 0;
- fft->comm[0] = comm;
-
- *pfft = fft;
- return 0;
-}
-
-int
-gmx_fft_init_2d_real(gmx_fft_t * pfft,
- int nx,
- int ny,
- enum gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- int info = 0;
- acmlComplex* comm = NULL;
- real* commRC = NULL;
- real* commCR = NULL;
- int commSize = 0;
- int nyc = 0;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- if( (fft = malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- nyc = (ny/2 + 1);
-
- /* Roll our own 2D real transform using multiple transforms in ACML,
- * since the current ACML versions does not support our storage format,
- * and all but the most recent don't even have 2D real FFTs.
- */
-
- // Single precision requires 5*nx+100
- // Double precision requires 3*nx+100
- if( sizeof( acmlComplex ) == 16 )
- commSize = (3*nx+100)*sizeof( acmlComplex );
- else
- commSize = (5*nx+100)*sizeof( acmlComplex );
-
- // Allocate communication work array
- if( (comm = (acmlComplex*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- // Initialize communication work array
- ACML_FFT1DMX( 100, 1.0f, FALSE, nyc, nx, NULL, nyc, 1, NULL, nyc, 1, comm, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- commSize = (3*ny+100)*sizeof( real );
- // Allocate communication work array
- if( (commRC = (real*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- // TODO: Is there no MODE or PLAN for multiple hermetian sequences?
- // Initialize communication work array
- ACML_RCFFT1D( 100, ny, NULL, commRC, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- commSize = (3*ny+100)*sizeof( real );
- // Allocate communication work array
- if( (commCR = (real*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- // TODO: Is there no MODE or PLAN for multiple hermetian sequences?
- // Initialize communication work array
- ACML_CRFFT1D( 100, ny, NULL, commCR, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- /* Allocate scratch work array that ACML uses to splat from hermitian complex format to
- * full complex format */
- if( (fft->realScratch = (acmlComplex*)malloc( (nx*ny)*sizeof( acmlComplex ) ) ) == NULL )
- {
- return ENOMEM;
- }
-
- fft->ndim = 2;
- fft->nx = nx;
- fft->ny = ny;
- fft->real_fft = 1;
- fft->comm[0] = comm;
- fft->comm[1] = commRC;
- fft->comm[2] = commCR;
-
- *pfft = fft;
- return 0;
-}
-
-int
-gmx_fft_init_3d(gmx_fft_t * pfft,
- int nx,
- int ny,
- int nz,
- enum gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- int info = 0;
- acmlComplex* comm = NULL;
- int commSize = 0;
-
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- if( (fft = malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- commSize = (nx*ny*nz+4*(nx+ny+nz)+300)*sizeof( acmlComplex );
-
- // Allocate communication work array
- if( (comm = (acmlComplex*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- ACML_FFT3DX( 100, 1.0f, TRUE, nx, ny, nz, NULL, 1, nx, nx*ny, NULL, 1, nx, nx*ny, comm, commSize, &info );
-
- fft->ndim = 3;
- fft->nx = nx;
- fft->ny = ny;
- fft->nz = nz;
- fft->real_fft = 0;
- fft->comm[0] = comm;
-
- *pfft = fft;
- return 0;
-}
-
-int
-gmx_fft_init_3d_real(gmx_fft_t * pfft,
- int nx,
- int ny,
- int nz,
- enum gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- int info = 0;
- acmlComplex* commX = NULL;
- acmlComplex* commY = NULL;
- real* commRC = NULL;
- real* commCR = NULL;
- int commSize = 0;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- /* nzc = (nz/2 + 1); */
-
- if( (fft = malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- /* Roll our own 3D real transform using multiple transforms in ACML,
- * since the current ACML versions does not support 2D
- * or 3D real transforms.
- */
-
- /* In-place X FFT.
- * ny*nz complex-to-complex transforms, length nx
- * transform distance: 1
- * element strides: ny*nz
- */
-
- /* Single precision requires 5*nx+100
- Double precision requires 3*nx+100
- */
- if( sizeof( acmlComplex ) == 16 )
- commSize = (3*nx+100)*sizeof( acmlComplex );
- else
- commSize = (5*nx+100)*sizeof( acmlComplex );
-
- /* Allocate communication work array */
- if( (commX = (acmlComplex*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- /* Initialize communication work array */
- ACML_FFT1DMX( 100, 1.0f, TRUE, ny*nz, nx, NULL, ny*nz, 1, NULL, ny*nz, 1, commX, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- /* In-place Y FFT.
- * We cannot do all NX*NZ transforms at once, so define a handle to do
- * NZ transforms, and then execute it NX times.
- * nz complex-to-complex transforms, length ny
- * transform distance: 1
- * element strides: nz
- */
- /* Single precision requires 5*nx+100
- Double precision requires 3*nx+100
- */
- if( sizeof( acmlComplex ) == 16 )
- commSize = (3*ny+100)*sizeof( acmlComplex );
- else
- commSize = (5*ny+100)*sizeof( acmlComplex );
-
- /* Allocate communication work array */
- if( (commY = (acmlComplex*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- /* Initialize communication work array */
- /* We want to do multiple 1D FFT's in z-y plane, so we have to loop over x
- * dimension recalculating z-y plane for each slice.
- */
- ACML_FFT1DMX( 100, 1.0f, TRUE, nz, ny, NULL, nz, 1, NULL, nz, 1, commY, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- /* In-place Z FFT:
- * nx*ny real-to-complex transforms, length nz
- * transform distance: nzc*2 -> nzc*2
- * element strides: 1
- */
-
- commSize = (3*nz+100)*sizeof( real );
- /* Allocate communication work array */
- if( (commRC = (real*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- /* TODO: Is there no MODE or PLAN for multiple hermetian sequences? */
- // Initialize communication work array
- ACML_RCFFT1D( 100, nz, NULL, commRC, &info );
-
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- /* Out-of-place complex-to-real (affects distance) Z FFT:
- * nx*ny real-to-complex transforms, length nz
- * transform distance: nzc*2 -> nz
- * element STRIDES: 1
- */
- commSize = (3*nz+100)*sizeof( real );
- /* Allocate communication work array */
- if( (commCR = (real*)malloc( commSize ) ) == NULL )
- {
- return ENOMEM;
- }
-
- // Initialize communication work array
- ACML_CRFFT1D( 100, nz, NULL, commCR, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error initializing ACML FFT; status=%d", info);
- gmx_fft_destroy( fft );
- return info;
- }
-
- /* Allocate scratch work array that ACML uses to splat from hermitian complex format to
- * full complex format */
- if( (fft->realScratch = (acmlComplex*)malloc( (nx*ny*nz)*sizeof( acmlComplex ) ) ) == NULL )
- {
- return ENOMEM;
- }
-
- fft->ndim = 3;
- fft->nx = nx;
- fft->ny = ny;
- fft->nz = nz;
- fft->real_fft = 1;
- fft->comm[0] = commX;
- fft->comm[1] = commY;
- fft->comm[2] = commRC;
- fft->comm[3] = commCR;
-
- *pfft = fft;
- return 0;
-}
-
-int
-gmx_fft_1d(gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int inpl = (in_data == out_data);
- int info = 0;
- int mode = (dir == GMX_FFT_FORWARD) ? -1 : 1;
-
- if( (fft->real_fft == 1) || (fft->ndim != 1) ||
- ((dir != GMX_FFT_FORWARD) && (dir != GMX_FFT_BACKWARD)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- ACML_FFT1DX( mode, 1.0f, inpl, fft->nx, in_data, 1, out_data, 1, fft->comm[0], &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error executing AMD ACML FFT.");
- info = -1;
- }
-
- return info;
-}
-
-int
-gmx_fft_1d_real(gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * inPtr,
- void * outPtr)
-{
- int info = 0;
- int nx = fft->nx;
-
- if( (fft->real_fft != 1) || (fft->ndim != 1) ||
- ((dir != GMX_FFT_REAL_TO_COMPLEX) && (dir != GMX_FFT_COMPLEX_TO_REAL)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- /* I apply a correction scale to the automatic scaling that was done in the real-complex step */
- float recipCorrection = sqrtf( nx );
-
- /* ACML needs to do complex2complex intermediate transforms, which will not fit in the amount
- * of memory allocated by the gromacs program, which assumes real. The real2complex transforms
- * are also only in-place, so we manually do a memcpy() first */
- if(dir==GMX_FFT_REAL_TO_COMPLEX)
- {
- memcpy( fft->realScratch, inPtr, nx*sizeof( real ) );
- ACML_RCFFT1D( 1, nx, fft->realScratch, fft->comm[0], &info );
-
- hermitianUnPacking( recipCorrection, 1, nx, fft->realScratch, 0, outPtr );
- }
- else
- {
- memcpy( fft->realScratch, inPtr, nx*sizeof( t_complex ) );
- hermitianPacking( recipCorrection, 1, nx, fft->realScratch, outPtr, 0 );
-
- ACML_CRFFT1D( 1, nx, outPtr, fft->comm[1], &info );
- }
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error executing AMD ACML FFT.");
- info = -1;
- }
-
- return info;
-}
-
-int
-gmx_fft_2d(gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int inpl = (in_data == out_data);
- int info = 0;
- int mode = (dir == GMX_FFT_FORWARD) ? -1 : 1;
-
- if( (fft->real_fft == 1) || (fft->ndim != 2) ||
- ((dir != GMX_FFT_FORWARD) && (dir != GMX_FFT_BACKWARD)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- ACML_FFT2DX( mode, 1.0f, TRUE, inpl, fft->nx, fft->ny,
- in_data, 1, fft->nx, out_data, 1, fft->nx, fft->comm[0], &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error executing AMD ACML FFT.");
- info = -1;
- }
-
- return info;
-}
-
-
-int
-gmx_fft_2d_real(gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * inPtr,
- void * outPtr )
-{
- int inpl = (inPtr == outPtr);
- int info = 0;
- int i = 0;
-
- int nx = fft->nx;
- int ny = fft->ny;
- int nyc = (fft->ny/2 + 1);
- int nyLen = 0;
-
- /* Depending on whether we are calculating the FFT in place or not, the rows of the
- * 2D FFT will be packed or not.
- */
- if( inpl )
- {
- /* Not packed; 1 or 2 extra reals padding at end of each row */
- nyLen = nyc*2;
- }
- else
- {
- nyLen = ny;
- }
-
- if( (fft->real_fft != 1) || (fft->ndim != 2) ||
- ((dir != GMX_FFT_REAL_TO_COMPLEX) && (dir != GMX_FFT_COMPLEX_TO_REAL)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- /* I apply a correction scale to the automatic scaling that was done in the real-complex step */
- float recipCorrection = sqrtf( ny );
-
- if(dir==GMX_FFT_REAL_TO_COMPLEX)
- {
- /* ACML needs to do complex2complex intermediate transforms, which will not fit in the amount
- * of memory allocated by the gromacs program, which assumes real. The real2complex transforms
- * are also only in-place, so we manually do a memcpy() first */
- if( !inpl )
- {
- memcpy( outPtr, inPtr, nx*ny*sizeof( real ) );
- }
-
- /* real-to-complex in Y dimension, in-place
- * SCFFTM is not valid to call here. SCFFTM does not take any stride information, and assumes that
- * the rows are tightly packed. GROMACS pads rows with either 1 or 2 extra reals depending
- * on even or odd lengths.
- */
- for( i = 0; i < nx; ++i )
- {
- if ( info == 0 )
- ACML_RCFFT1D( 1, ny, (real*)outPtr+i*nyLen, fft->comm[1], &info );
- }
-
- hermitianUnPacking( 1.0f, nx, ny, outPtr, nyLen, fft->realScratch );
-
- /* complex-to-complex in X dimension */
- if ( info == 0 )
- ACML_FFT1DMX( -1, recipCorrection, FALSE, nyc, nx, fft->realScratch, nyc, 1,
- outPtr, nyc, 1, fft->comm[0], &info );
- }
- else
- {
- /* complex-to-complex in X dimension, in-place */
- ACML_FFT1DMX( 1, recipCorrection, FALSE, nyc, nx, inPtr, nyc, 1,
- fft->realScratch, nyc, 1, fft->comm[0], &info );
-
- hermitianPacking( 1.0f, nx, ny, fft->realScratch, outPtr, nyLen );
-
- /* complex-to-real in Y dimension
- * CSFFTM is not valid to call here. SCFFTM does not take any stride information, and assumes that
- * the rows are tightly packed. GROMACS pads rows with either 1 or 2 extra reals depending
- * on even or odd lengths.
- */
- if ( info == 0 )
- {
- for( i = 0; i < nx; ++i )
- {
- if ( info == 0 )
- ACML_CRFFT1D( 1, ny, (real*)outPtr+i*nyLen, fft->comm[2], &info );
- }
- }
- }
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error executing AMD ACML FFT.");
- info = -1;
- }
-
- return info;
-}
-
-
-int
-gmx_fft_3d(gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int mode = (dir == GMX_FFT_FORWARD) ? -1 : 1;
- int inpl = ( in_data == out_data );
-
- int commSize = 0;
- int nx = fft->nx;
- int ny = fft->ny;
- int nz = fft->nz;
- int info = 0;
-
- if( (fft->real_fft == 1) || (fft->ndim != 3) ||
- ((dir != GMX_FFT_FORWARD) && (dir != GMX_FFT_BACKWARD)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- commSize = (nx*ny*nz+4*(nx+ny+nz)+300)*sizeof( acmlComplex );
-
- ACML_FFT3DX( mode, 1.0f, inpl, nx, ny, nz, in_data, 1, nx, nx*ny,
- out_data, 1, nx, nx*ny, fft->comm[0], commSize, &info );
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error executing AMD ACML FFT.");
- info = -1;
- }
-
- return info;
-}
-
-int
-gmx_fft_3d_real(gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * inPtr,
- void * outPtr)
-{
- int inpl = (inPtr == outPtr);
- int info = 0;
- int i;
- int nx,ny,nz, nzLen = 0;
-
- nx = fft->nx;
- ny = fft->ny;
- nz = fft->nz;
- int nzc = (nz/2 + 1);
-
- /* Depending on whether we are calculating the FFT in place or not, the rows of the
- * 3D FFT will be packed or not.
- */
- if( inpl )
- {
- /* Not packed; 1 or 2 extra reals padding at end of each row */
- nzLen = nzc*2;
- }
- else
- {
- nzLen = nz;
- }
-
- if( (fft->real_fft != 1) || (fft->ndim != 3) ||
- ((dir != GMX_FFT_REAL_TO_COMPLEX) && (dir != GMX_FFT_COMPLEX_TO_REAL)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- /* I apply a correction scale to the automatic scaling that was done in the real-complex step */
- float recipCorrection = sqrtf( nz );
-
- if(dir==GMX_FFT_REAL_TO_COMPLEX)
- {
- /* ACML needs to do complex2complex intermediate transforms, which will not fit in the amount
- * of memory allocated by the gromacs program, which assumes real. The real2complex transforms
- * are also only in-place, so we manually do a memcpy() first */
- if( !inpl )
- {
- memcpy( outPtr, inPtr, nx*ny*nz*sizeof( real ) );
- }
-
- /* real-to-complex in Z dimension, in-place
- * SCFFTM is not valid to call here. SCFFTM does not take any stride information, and assumes that
- * the rows are tightly packed. GROMACS pads rows with either 1 or 2 extra reals depending
- * on even or odd lengths.
- */
- for( i = 0; i < nx*ny; ++i )
- {
- if ( info == 0 )
- ACML_RCFFT1D( 1, nz, (real*)outPtr+i*nzLen, fft->comm[2], &info );
- }
-
- hermitianUnPacking( 1.0f, nx*ny, nz, outPtr, nzLen, fft->realScratch );
-
- /* complex-to-complex in Y dimension, in-place */
- for(i=0;i<nx;i++)
- {
- if ( info == 0 )
- ACML_FFT1DMX( -1, 1.0f, TRUE, nzc, ny, (acmlComplex*)fft->realScratch+i*ny*nzc, nzc, 1,
- (acmlComplex*)fft->realScratch+i*ny*nzc, nzc, 1, fft->comm[1], &info );
- }
-
- /* complex-to-complex in X dimension, in-place */
- if ( info == 0 )
- ACML_FFT1DMX( -1, recipCorrection, FALSE, ny*nzc, nx, fft->realScratch, ny*nzc, 1, outPtr, ny*nzc, 1, fft->comm[0], &info );
- }
- else
- {
- /* complex-to-complex in X dimension, from inPtr to work */
- ACML_FFT1DMX( 1, recipCorrection, FALSE, ny*nzc, nx, inPtr, ny*nzc, 1, fft->realScratch, ny*nzc, 1, fft->comm[0], &info );
-
- /* complex-to-complex in Y dimension, in-place */
- for( i=0; i < nx; i++ )
- {
- if ( info == 0 )
- ACML_FFT1DMX( 1, 1.0f, TRUE, nzc, ny, (acmlComplex*)fft->realScratch+i*ny*nzc, nzc, 1,
- (acmlComplex*)fft->realScratch+i*ny*nzc, nzc, 1, fft->comm[1], &info );
- }
-
- hermitianPacking( 1.0f, nx*ny, nz, fft->realScratch, outPtr, nzLen );
-
- /* complex-to-real in Z dimension, in-place
- * CSFFTM is not valid to call here. CSFFTM does not take any stride information, and assumes that
- * the rows are tightly packed. GROMACS pads rows with either 1 or 2 extra reals depending
- * on even or odd lengths.
- */
- for( i = 0; i < nx*ny; ++i )
- {
- if ( info == 0 )
- ACML_CRFFT1D( 1, nz, (real*)outPtr+i*nzLen, fft->comm[3], &info );
- }
-
- }
-
- if( info != 0 )
- {
- gmx_fatal(FARGS,"Error executing AMD ACML FFT.");
- info = -1;
- }
-
- return info;
-}
-
-void
-gmx_fft_destroy(gmx_fft_t fft)
-{
- int d;
-
- if(fft != NULL)
- {
- for(d=0;d<4;d++)
- {
- if(fft->comm[d] != NULL)
- {
- free(fft->comm[d]);
- }
- }
-
- if( fft->realScratch != NULL)
- free( fft->realScratch );
-
- free( fft );
- }
-}
-
-#else
-int
-gmx_fft_acml_empty;
-#endif /* GMX_FFT_ACML */
return 0;
}
-
-
-int
-gmx_fft_init_2d(gmx_fft_t * pfft,
- int nx,
- int ny,
- int flags)
-{
- gmx_fft_t fft;
- int rc;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid FFT opaque type pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- /* Create the X transform */
- if( (rc = gmx_fft_init_1d(&fft,nx,flags)) != 0)
- {
- return rc;
- }
-
- /* Create Y transform as a link from X */
- if( (rc=gmx_fft_init_1d(&(fft->next),ny,flags)) != 0)
- {
- free(fft);
- return rc;
- }
-
- *pfft = fft;
- return 0;
-};
-
-
int
gmx_fft_init_2d_real(gmx_fft_t * pfft,
int nx,
int
-gmx_fft_init_3d(gmx_fft_t * pfft,
- int nx,
- int ny,
- int nz,
- int flags)
-{
- gmx_fft_t fft;
- int rc;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid FFT opaque type pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- /* Create the X transform */
-
- if( (fft = (struct gmx_fft *)malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- fft->n = nx;
-
- /* Need 4*nx storage for 1D complex FFT, and another
- * 2*nz elements for gmx_fft_transpose_2d_nelem() storage.
- */
- if( (fft->work = (real *)malloc(sizeof(real)*(4*nx+2*nz))) == NULL)
- {
- free(fft);
- return ENOMEM;
- }
-
- fftpack_cffti1(nx,fft->work,fft->ifac);
-
-
- /* Create 2D Y/Z transforms as a link from X */
- if( (rc=gmx_fft_init_2d(&(fft->next),ny,nz,flags)) != 0)
- {
- free(fft);
- return rc;
- }
-
- *pfft = fft;
- return 0;
-};
-
-
-int
-gmx_fft_init_3d_real(gmx_fft_t * pfft,
- int nx,
- int ny,
- int nz,
- int flags)
-{
- gmx_fft_t fft;
- int nzc = (nz/2 + 1);
- int rc;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid FFT opaque type pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- /* Create the X transform */
- if( (fft = (struct gmx_fft *)malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- fft->n = nx;
-
- /* Need 4*nx storage for 1D complex FFT, another
- * 2*nx*ny*nzc elements to copy the entire 3D matrix when
- * doing out-of-place complex-to-real FFTs, and finally
- * 2*nzc elements for transpose work space.
- */
- if( (fft->work = (real *)malloc(sizeof(real)*(4*nx+2*nx*ny*nzc+2*nzc))) == NULL)
- {
- free(fft);
- return ENOMEM;
- }
- fftpack_cffti1(nx,fft->work,fft->ifac);
-
- /* Create 2D real Y/Z transform as a link from X */
- if( (rc=gmx_fft_init_2d_real(&(fft->next),ny,nz,flags)) != 0)
- {
- free(fft);
- return rc;
- }
-
- *pfft = fft;
- return 0;
-}
-
-
-int
gmx_fft_1d (gmx_fft_t fft,
enum gmx_fft_direction dir,
void * in_data,
}
-int
-gmx_fft_2d (gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int i,nx,ny;
- t_complex * data;
-
- nx = fft->n;
- ny = fft->next->n;
-
- /* FFTPACK only does in-place transforms, so emulate out-of-place
- * by copying data to the output array first.
- * For 2D there is likely enough data to benefit from memcpy().
- */
- if( in_data != out_data )
- {
- memcpy(out_data,in_data,sizeof(t_complex)*nx*ny);
- }
-
- /* Much easier to do pointer arithmetic when base has the correct type */
- data = (t_complex *)out_data;
-
- /* y transforms */
- for(i=0;i<nx;i++)
- {
- gmx_fft_1d(fft->next,dir,data+i*ny,data+i*ny);
- }
-
- /* Transpose in-place to get data in place for x transform now */
- gmx_fft_transpose_2d(data,data,nx,ny);
-
- /* x transforms */
- for(i=0;i<ny;i++)
- {
- gmx_fft_1d(fft,dir,data+i*nx,data+i*nx);
- }
-
- /* Transpose in-place to get data back in original order */
- gmx_fft_transpose_2d(data,data,ny,nx);
-
- return 0;
-}
-
-
-
int
gmx_fft_2d_real (gmx_fft_t fft,
enum gmx_fft_direction dir,
return 0;
}
-
-
-int
-gmx_fft_3d (gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int i,nx,ny,nz,rc;
- t_complex * data;
- t_complex * work;
- nx=fft->n;
- ny=fft->next->n;
- nz=fft->next->next->n;
-
- /* First 4*nx positions are FFTPACK workspace, then ours starts */
- work = (t_complex *)(fft->work+4*nx);
-
- /* FFTPACK only does in-place transforms, so emulate out-of-place
- * by copying data to the output array first.
- * For 3D there is likely enough data to benefit from memcpy().
- */
- if( in_data != out_data )
- {
- memcpy(out_data,in_data,sizeof(t_complex)*nx*ny*nz);
- }
-
- /* Much easier to do pointer arithmetic when base has the correct type */
- data = (t_complex *)out_data;
-
- /* Perform z transforms */
- for(i=0;i<nx*ny;i++)
- gmx_fft_1d(fft->next->next,dir,data+i*nz,data+i*nz);
-
- /* For each X slice, transpose the y & z dimensions inside the slice */
- for(i=0;i<nx;i++)
- {
- gmx_fft_transpose_2d(data+i*ny*nz,data+i*ny*nz,ny,nz);
- }
-
- /* Array is now (nx,nz,ny) - perform y transforms */
- for(i=0;i<nx*nz;i++)
- {
- gmx_fft_1d(fft->next,dir,data+i*ny,data+i*ny);
- }
-
- /* Transpose back to (nx,ny,nz) */
- for(i=0;i<nx;i++)
- {
- gmx_fft_transpose_2d(data+i*ny*nz,data+i*ny*nz,nz,ny);
- }
-
- /* Transpose entire x & y slices to go from
- * (nx,ny,nz) to (ny,nx,nz).
- * Use work data elements 4*n .. 4*n+2*nz-1.
- */
- rc=gmx_fft_transpose_2d_nelem(data,data,nx,ny,nz,work);
- if( rc != 0)
- {
- gmx_fatal(FARGS,"Cannot transpose X & Y/Z in gmx_fft_3d().");
- return rc;
- }
-
- /* Then go from (ny,nx,nz) to (ny,nz,nx) */
- for(i=0;i<ny;i++)
- {
- gmx_fft_transpose_2d(data+i*nx*nz,data+i*nx*nz,nx,nz);
- }
-
- /* Perform x transforms */
- for(i=0;i<ny*nz;i++)
- {
- gmx_fft_1d(fft,dir,data+i*nx,data+i*nx);
- }
-
- /* Transpose back from (ny,nz,nx) to (ny,nx,nz) */
- for(i=0;i<ny;i++)
- {
- gmx_fft_transpose_2d(data+i*nz*nx,data+i*nz*nx,nz,nx);
- }
-
- /* Transpose from (ny,nx,nz) to (nx,ny,nz)
- * Use work data elements 4*n .. 4*n+2*nz-1.
- */
- rc = gmx_fft_transpose_2d_nelem(data,data,ny,nx,nz,work);
- if( rc != 0)
- {
- gmx_fatal(FARGS,"Cannot transpose Y/Z & X in gmx_fft_3d().");
- return rc;
- }
-
- return 0;
-}
-
-
-int
-gmx_fft_3d_real (gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int i,j,k;
- int nx,ny,nz,nzc,rc;
- t_complex * data;
- t_complex * work_transp;
- t_complex * work_c2r;
- real * p1;
- real * p2;
-
- nx=fft->n;
- ny=fft->next->n;
- nz=fft->next->next->n;
- nzc=(nz/2+1);
-
-
- /* First 4*nx positions are FFTPACK workspace, then ours starts.
- * We have 2*nx*ny*nzc elements for temp complex-to-real storage when
- * doing out-of-place transforms, and another 2*nzc for transpose data.
- */
- work_c2r = (t_complex *)(fft->work+4*nx);
- work_transp = (t_complex *)(fft->work+4*nx+2*nx*ny*nzc);
-
- /* Much easier to do pointer arithmetic when base has the correct type */
- data = (t_complex *)out_data;
-
- if(dir==GMX_FFT_REAL_TO_COMPLEX)
- {
- /* FFTPACK only does in-place transforms, so emulate out-of-place
- * by copying data to the output array first. This is guaranteed to
- * work for real-to-complex since complex data is larger than the real.
- * For 3D there is likely enough data to benefit from memcpy().
- */
- if( in_data != out_data )
- {
- p1 = (real *)in_data;
- p2 = (real *)out_data;
-
- for(i=0;i<nx;i++)
- {
- for(j=0;j<ny;j++)
- {
- for(k=0;k<nz;k++)
- {
- p2[(i*ny+j)*2*nzc+k] = p1[(i*ny+j)*nz+k];
- }
- }
- }
- }
- data = (t_complex *)out_data;
-
- /* Transform the Y/Z slices real-to-complex */
- for(i=0;i<nx;i++)
- {
- gmx_fft_2d_real(fft->next,dir,data+i*ny*nzc,data+i*ny*nzc);
- }
-
- /* Transpose x & y slices to go from
- * (nx,ny,nzc) to (ny,nx,nzc).
- */
- rc=gmx_fft_transpose_2d_nelem(data,data,nx,ny,nzc,work_transp);
- if( rc != 0)
- {
- gmx_fatal(FARGS,"Cannot transpose X & Y/Z gmx_fft_3d_real().");
- return rc;
- }
-
- /* Then transpose from (ny,nx,nzc) to (ny,nzc,nx) */
- for(i=0;i<ny;i++)
- {
- gmx_fft_transpose_2d(data+i*nx*nzc,data+i*nx*nzc,nx,nzc);
- }
-
- /* Perform x transforms */
- for(i=0;i<ny*nzc;i++)
- {
- gmx_fft_1d(fft,GMX_FFT_FORWARD,data+i*nx,data+i*nx);
- }
-
- /* Transpose from (ny,nzc,nx) back to (ny,nx,nzc) */
- for(i=0;i<ny;i++)
- {
- gmx_fft_transpose_2d(data+i*nzc*nx,data+i*nzc*nx,nzc,nx);
- }
-
- /* Transpose back from (ny,nx,nzc) to (nx,ny,nz) */
- rc=gmx_fft_transpose_2d_nelem(data,data,ny,nx,nzc,work_transp);
- if( rc != 0)
- {
- gmx_fatal(FARGS,"Cannot transpose Y/Z & X in gmx_fft_3d_real().");
- return rc;
- }
-
- }
- else if(dir==GMX_FFT_COMPLEX_TO_REAL)
- {
- /* An in-place complex-to-real transform is straightforward,
- * since the output array must be large enough for the padding to fit.
- *
- * For out-of-place complex-to-real transforms we cannot just copy
- * data to the output array, since it is smaller than the input.
- * In this case there's nothing to do but employing temporary work data.
- */
- if(in_data != out_data)
- {
- memcpy(work_c2r,in_data,sizeof(t_complex)*nx*ny*nzc);
- data = (t_complex *)work_c2r;
- }
- else
- {
- /* in-place */
- data = (t_complex *)out_data;
- }
-
- /* Transpose x & y slices to go from
- * (nx,ny,nz) to (ny,nx,nz).
- */
- gmx_fft_transpose_2d_nelem(data,data,nx,ny,nzc,work_transp);
-
- /* Then go from (ny,nx,nzc) to (ny,nzc,nx) */
- for(i=0;i<ny;i++)
- {
- gmx_fft_transpose_2d(data+i*nx*nzc,data+i*nx*nzc,nx,nzc);
- }
-
-
- /* Perform x transforms */
- for(i=0;i<ny*nzc;i++)
- {
- gmx_fft_1d(fft,GMX_FFT_BACKWARD,data+i*nx,data+i*nx);
- }
-
- /* Transpose back from (ny,nzc,nx) to (ny,nx,nzc) */
- for(i=0;i<ny;i++)
- {
- gmx_fft_transpose_2d(data+i*nzc*nx,data+i*nzc*nx,nzc,nx);
- }
-
- /* Transpose back from (ny,nx,nzc) to (nx,ny,nz) */
- gmx_fft_transpose_2d_nelem(data,data,ny,nx,nzc,work_transp);
-
-
- /* Do 2D complex-to-real */
- for(i=0;i<nx;i++)
- {
- gmx_fft_2d_real(fft->next,dir,data+i*ny*nzc,data+i*ny*nzc);
- }
-
- if( in_data != out_data )
- {
- /* Output (pointed to by data) is now in padded format.
- * Pack it into out_data if we were doing an out-of-place transform.
- */
- p1 = (real *)data;
- p2 = (real *)out_data;
-
- for(i=0;i<nx;i++)
- {
- for(j=0;j<ny;j++)
- {
- for(k=0;k<nz;k++)
- {
- p2[(i*ny+j)*nz+k] = p1[(i*ny+j)*nzc*2+k];
- }
- }
- }
- }
-
- }
- else
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- return 0;
-}
-
-
-
-
void
gmx_fft_destroy(gmx_fft_t fft)
{
free(fft);
}
}
+
+void gmx_fft_cleanup()
+{
+}
#endif /* GMX_FFT_FFTPACK */
#include "gmx_fft.h"
#include "gmx_fatal.h"
-
#ifdef GMX_DOUBLE
#define FFTWPREFIX(name) fftw_ ## name
#else
#define FFTWPREFIX(name) fftwf_ ## name
#endif
-#ifdef GMX_THREAD_MPI
-#include "thread_mpi/threads.h"
-#endif
-
+#include "thread_mpi/mutex.h"
+#include "gromacs/utility/exceptions.h"
-#ifdef GMX_THREAD_MPI
/* none of the fftw3 calls, except execute(), are thread-safe, so
we need to serialize them with this mutex. */
-static tMPI_Thread_mutex_t big_fftw_mutex=TMPI_THREAD_MUTEX_INITIALIZER;
-static gmx_bool gmx_fft_threads_initialized=FALSE;
-#define FFTW_LOCK tMPI_Thread_mutex_lock(&big_fftw_mutex)
-#define FFTW_UNLOCK tMPI_Thread_mutex_unlock(&big_fftw_mutex)
-#else /* GMX_THREAD_MPI */
-#define FFTW_LOCK
-#define FFTW_UNLOCK
-#endif /* GMX_THREAD_MPI */
+static tMPI::mutex big_fftw_mutex;
+#define FFTW_LOCK try { big_fftw_mutex.lock(); } GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
+#define FFTW_UNLOCK try { big_fftw_mutex.unlock(); } GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
/* We assume here that aligned memory starts at multiple of 16 bytes and unaligned memory starts at multiple of 8 bytes. The later is guranteed for all malloc implementation.
Consequesences:
int ndim;
};
-
-
int
gmx_fft_init_1d(gmx_fft_t * pfft,
int nx,
- gmx_fft_flag flags)
+ gmx_fft_flag flags)
{
return gmx_fft_init_many_1d(pfft,nx,1,flags);
}
int
gmx_fft_init_many_1d(gmx_fft_t * pfft,
- int nx,
- int howmany,
- gmx_fft_flag flags)
+ int nx,
+ int howmany,
+ gmx_fft_flag flags)
{
gmx_fft_t fft;
FFTWPREFIX(complex) *p1,*p2,*up1,*up2;
size_t pc;
int i,j,k;
int fftw_flags;
-
+
#ifdef GMX_DISABLE_FFTW_MEASURE
flags |= GMX_FFT_FLAG_CONSERVATIVE;
#endif
-
- fftw_flags = (flags & GMX_FFT_FLAG_CONSERVATIVE) ? FFTW_ESTIMATE : FFTW_MEASURE;
+
+ fftw_flags = (flags & GMX_FFT_FLAG_CONSERVATIVE) ? FFTW_ESTIMATE : FFTW_MEASURE;
if(pfft==NULL)
{
return EINVAL;
}
*pfft = NULL;
-
+
FFTW_LOCK;
if( (fft = (gmx_fft_t)FFTWPREFIX(malloc)(sizeof(struct gmx_fft))) == NULL)
{
FFTW_UNLOCK;
return ENOMEM;
- }
-
+ }
+
/* allocate aligned, and extra memory to make it unaligned */
p1 = (FFTWPREFIX(complex) *) FFTWPREFIX(malloc)(sizeof(FFTWPREFIX(complex))*(nx+2)*howmany);
if(p1==NULL)
FFTW_UNLOCK;
return ENOMEM;
}
-
+
p2 = (FFTWPREFIX(complex) *) FFTWPREFIX(malloc)(sizeof(FFTWPREFIX(complex))*(nx+2)*howmany);
if(p2==NULL)
{
FFTW_UNLOCK;
return ENOMEM;
}
-
- /* make unaligned pointers.
+
+ /* make unaligned pointers.
* In double precision the actual complex datatype will be 16 bytes,
* so go to a char pointer and force an offset of 8 bytes instead.
*/
pc = (size_t)p1;
- pc += 8;
+ pc += 8;
up1 = (FFTWPREFIX(complex) *)pc;
-
+
pc = (size_t)p2;
- pc += 8;
+ pc += 8;
up2 = (FFTWPREFIX(complex) *)pc;
-
+
/* int rank, const int *n, int howmany,
fftw_complex *in, const int *inembed,
int istride, int idist,
fftw_complex *out, const int *onembed,
int ostride, int odist,
int sign, unsigned flags */
- fft->plan[0][0][0] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,up1,&nx,1,nx,up2,&nx,1,nx,FFTW_BACKWARD,fftw_flags);
- fft->plan[0][0][1] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,up1,&nx,1,nx,up2,&nx,1,nx,FFTW_FORWARD,fftw_flags);
- fft->plan[0][1][0] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,up1,&nx,1,nx,up1,&nx,1,nx,FFTW_BACKWARD,fftw_flags);
- fft->plan[0][1][1] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,up1,&nx,1,nx,up1,&nx,1,nx,FFTW_FORWARD,fftw_flags);
- fft->plan[1][0][0] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,p1,&nx,1,nx,p2,&nx,1,nx,FFTW_BACKWARD,fftw_flags);
- fft->plan[1][0][1] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,p1,&nx,1,nx,p2,&nx,1,nx,FFTW_FORWARD,fftw_flags);
- fft->plan[1][1][0] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,p1,&nx,1,nx,p1,&nx,1,nx,FFTW_BACKWARD,fftw_flags);
- fft->plan[1][1][1] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,p1,&nx,1,nx,p1,&nx,1,nx,FFTW_FORWARD,fftw_flags);
+ fft->plan[0][0][0] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,up1,&nx,1,nx,up2,&nx,1,nx,FFTW_BACKWARD,fftw_flags);
+ fft->plan[0][0][1] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,up1,&nx,1,nx,up2,&nx,1,nx,FFTW_FORWARD,fftw_flags);
+ fft->plan[0][1][0] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,up1,&nx,1,nx,up1,&nx,1,nx,FFTW_BACKWARD,fftw_flags);
+ fft->plan[0][1][1] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,up1,&nx,1,nx,up1,&nx,1,nx,FFTW_FORWARD,fftw_flags);
+ fft->plan[1][0][0] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,p1,&nx,1,nx,p2,&nx,1,nx,FFTW_BACKWARD,fftw_flags);
+ fft->plan[1][0][1] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,p1,&nx,1,nx,p2,&nx,1,nx,FFTW_FORWARD,fftw_flags);
+ fft->plan[1][1][0] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,p1,&nx,1,nx,p1,&nx,1,nx,FFTW_BACKWARD,fftw_flags);
+ fft->plan[1][1][1] = FFTWPREFIX(plan_many_dft)(1,&nx,howmany,p1,&nx,1,nx,p1,&nx,1,nx,FFTW_FORWARD,fftw_flags);
for(i=0;i<2;i++)
{
}
}
}
- }
-
+ }
+
FFTWPREFIX(free)(p1);
FFTWPREFIX(free)(p2);
-
+
fft->real_transform = 0;
fft->ndim = 1;
-
+
*pfft = fft;
FFTW_UNLOCK;
return 0;
}
-
int
gmx_fft_init_1d_real(gmx_fft_t * pfft,
int nx,
}
-
-int
-gmx_fft_init_2d(gmx_fft_t * pfft,
- int nx,
- int ny,
- gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- FFTWPREFIX(complex) *p1,*p2,*up1,*up2;
- size_t pc;
- int i,j,k;
- int fftw_flags;
-
-#ifdef GMX_DISABLE_FFTW_MEASURE
- flags |= GMX_FFT_FLAG_CONSERVATIVE;
-#endif
-
- fftw_flags = (flags & GMX_FFT_FLAG_CONSERVATIVE) ? FFTW_ESTIMATE : FFTW_MEASURE;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- FFTW_LOCK;
- if( (fft = (gmx_fft_t) FFTWPREFIX(malloc)(sizeof(struct gmx_fft))) == NULL)
- {
- FFTW_UNLOCK;
- return ENOMEM;
- }
-
- /* allocate aligned, and extra memory to make it unaligned */
- p1 = (FFTWPREFIX(complex) *) FFTWPREFIX(malloc)(sizeof(FFTWPREFIX(complex))*(nx*ny+2));
- if(p1==NULL)
- {
- FFTWPREFIX(free)(fft);
- FFTW_UNLOCK;
- return ENOMEM;
- }
-
- p2 = (FFTWPREFIX(complex) *) FFTWPREFIX(malloc)(sizeof(FFTWPREFIX(complex))*(nx*ny+2));
- if(p2==NULL)
- {
- FFTWPREFIX(free)(p1);
- FFTWPREFIX(free)(fft);
- FFTW_UNLOCK;
- return ENOMEM;
- }
-
- /* make unaligned pointers.
- * In double precision the actual complex datatype will be 16 bytes,
- * so go to a char pointer and force an offset of 8 bytes instead.
- */
- pc = (size_t)p1;
- pc += 8;
- up1 = (FFTWPREFIX(complex) *)pc;
-
- pc = (size_t)p2;
- pc += 8;
- up2 = (FFTWPREFIX(complex) *)pc;
-
-
- fft->plan[0][0][0] = FFTWPREFIX(plan_dft_2d)(nx,ny,up1,up2,FFTW_BACKWARD,fftw_flags);
- fft->plan[0][0][1] = FFTWPREFIX(plan_dft_2d)(nx,ny,up1,up2,FFTW_FORWARD,fftw_flags);
- fft->plan[0][1][0] = FFTWPREFIX(plan_dft_2d)(nx,ny,up1,up1,FFTW_BACKWARD,fftw_flags);
- fft->plan[0][1][1] = FFTWPREFIX(plan_dft_2d)(nx,ny,up1,up1,FFTW_FORWARD,fftw_flags);
-
- fft->plan[1][0][0] = FFTWPREFIX(plan_dft_2d)(nx,ny,p1,p2,FFTW_BACKWARD,fftw_flags);
- fft->plan[1][0][1] = FFTWPREFIX(plan_dft_2d)(nx,ny,p1,p2,FFTW_FORWARD,fftw_flags);
- fft->plan[1][1][0] = FFTWPREFIX(plan_dft_2d)(nx,ny,p1,p1,FFTW_BACKWARD,fftw_flags);
- fft->plan[1][1][1] = FFTWPREFIX(plan_dft_2d)(nx,ny,p1,p1,FFTW_FORWARD,fftw_flags);
-
-
- for(i=0;i<2;i++)
- {
- for(j=0;j<2;j++)
- {
- for(k=0;k<2;k++)
- {
- if(fft->plan[i][j][k] == NULL)
- {
- gmx_fatal(FARGS,"Error initializing FFTW3 plan.");
- FFTW_UNLOCK;
- gmx_fft_destroy(fft);
- FFTW_LOCK;
- FFTWPREFIX(free)(p1);
- FFTWPREFIX(free)(p2);
- FFTW_UNLOCK;
- return -1;
- }
- }
- }
- }
-
- FFTWPREFIX(free)(p1);
- FFTWPREFIX(free)(p2);
-
- fft->real_transform = 0;
- fft->ndim = 2;
-
- *pfft = fft;
- FFTW_UNLOCK;
- return 0;
-}
-
-
-
int
gmx_fft_init_2d_real(gmx_fft_t * pfft,
int nx,
return 0;
}
-
-
-int
-gmx_fft_init_3d(gmx_fft_t * pfft,
- int nx,
- int ny,
- int nz,
- gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- FFTWPREFIX(complex) *p1,*p2,*up1,*up2;
- size_t pc;
- int i,j,k;
- int fftw_flags;
-
-#ifdef GMX_DISABLE_FFTW_MEASURE
- flags |= GMX_FFT_FLAG_CONSERVATIVE;
-#endif
-
- fftw_flags = (flags & GMX_FFT_FLAG_CONSERVATIVE) ? FFTW_ESTIMATE : FFTW_MEASURE;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- FFTW_LOCK;
- if( (fft = (gmx_fft_t) FFTWPREFIX(malloc)(sizeof(struct gmx_fft))) == NULL)
- {
- FFTW_UNLOCK;
- return ENOMEM;
- }
-
- /* allocate aligned, and extra memory to make it unaligned */
- p1 = (FFTWPREFIX(complex) *) FFTWPREFIX(malloc)(sizeof(FFTWPREFIX(complex))*(nx*ny*nz+2));
- if(p1==NULL)
- {
- FFTWPREFIX(free)(fft);
- FFTW_UNLOCK;
- return ENOMEM;
- }
-
- p2 = (FFTWPREFIX(complex) *) FFTWPREFIX(malloc)(sizeof(FFTWPREFIX(complex))*(nx*ny*nz+2));
- if(p2==NULL)
- {
- FFTWPREFIX(free)(p1);
- FFTWPREFIX(free)(fft);
- FFTW_UNLOCK;
- return ENOMEM;
- }
-
- /* make unaligned pointers.
- * In double precision the actual complex datatype will be 16 bytes,
- * so go to a char pointer and force an offset of 8 bytes instead.
- */
- pc = (size_t)p1;
- pc += 8;
- up1 = (FFTWPREFIX(complex) *)pc;
-
- pc = (size_t)p2;
- pc += 8;
- up2 = (FFTWPREFIX(complex) *)pc;
-
-
- fft->plan[0][0][0] = FFTWPREFIX(plan_dft_3d)(nx,ny,nz,up1,up2,FFTW_BACKWARD,fftw_flags);
- fft->plan[0][0][1] = FFTWPREFIX(plan_dft_3d)(nx,ny,nz,up1,up2,FFTW_FORWARD,fftw_flags);
- fft->plan[0][1][0] = FFTWPREFIX(plan_dft_3d)(nx,ny,nz,up1,up1,FFTW_BACKWARD,fftw_flags);
- fft->plan[0][1][1] = FFTWPREFIX(plan_dft_3d)(nx,ny,nz,up1,up1,FFTW_FORWARD,fftw_flags);
-
- fft->plan[1][0][0] = FFTWPREFIX(plan_dft_3d)(nx,ny,nz,p1,p2,FFTW_BACKWARD,fftw_flags);
- fft->plan[1][0][1] = FFTWPREFIX(plan_dft_3d)(nx,ny,nz,p1,p2,FFTW_FORWARD,fftw_flags);
- fft->plan[1][1][0] = FFTWPREFIX(plan_dft_3d)(nx,ny,nz,p1,p1,FFTW_BACKWARD,fftw_flags);
- fft->plan[1][1][1] = FFTWPREFIX(plan_dft_3d)(nx,ny,nz,p1,p1,FFTW_FORWARD,fftw_flags);
-
-
- for(i=0;i<2;i++)
- {
- for(j=0;j<2;j++)
- {
- for(k=0;k<2;k++)
- {
- if(fft->plan[i][j][k] == NULL)
- {
- gmx_fatal(FARGS,"Error initializing FFTW3 plan.");
- FFTW_UNLOCK;
- gmx_fft_destroy(fft);
- FFTW_LOCK;
- FFTWPREFIX(free)(p1);
- FFTWPREFIX(free)(p2);
- FFTW_UNLOCK;
- return -1;
- }
- }
- }
- }
-
- FFTWPREFIX(free)(p1);
- FFTWPREFIX(free)(p2);
-
- fft->real_transform = 0;
- fft->ndim = 3;
-
- *pfft = fft;
- FFTW_UNLOCK;
- return 0;
-}
-
-
-
-int
-gmx_fft_init_3d_real(gmx_fft_t * pfft,
- int nx,
- int ny,
- int nz,
- gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- real *p1,*p2,*up1,*up2;
- size_t pc;
- int i,j,k;
- int fftw_flags;
-
-#ifdef GMX_DISABLE_FFTW_MEASURE
- flags |= GMX_FFT_FLAG_CONSERVATIVE;
-#endif
-
- fftw_flags = (flags & GMX_FFT_FLAG_CONSERVATIVE) ? FFTW_ESTIMATE : FFTW_MEASURE;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- FFTW_LOCK;
- if( (fft = (gmx_fft_t) FFTWPREFIX(malloc)(sizeof(struct gmx_fft))) == NULL)
- {
- FFTW_UNLOCK;
- return ENOMEM;
- }
-
- /* allocate aligned, and extra memory to make it unaligned */
- p1 = (real *) FFTWPREFIX(malloc)(sizeof(real)*( nx*ny*(nz/2+1)*2 + 2) );
- if(p1==NULL)
- {
- FFTWPREFIX(free)(fft);
- FFTW_UNLOCK;
- return ENOMEM;
- }
-
- p2 = (real *) FFTWPREFIX(malloc)(sizeof(real)*( nx*ny*(nz/2+1)*2 + 2) );
- if(p2==NULL)
- {
- FFTWPREFIX(free)(p1);
- FFTWPREFIX(free)(fft);
- FFTW_UNLOCK;
- return ENOMEM;
- }
-
- /* make unaligned pointers.
- * In double precision the actual complex datatype will be 16 bytes,
- * so go to a void pointer and force an offset of 8 bytes instead.
- */
- pc = (size_t)p1;
- pc += 8;
- up1 = (real *)pc;
-
- pc = (size_t)p2;
- pc += 8;
- up2 = (real *)pc;
-
-
- fft->plan[0][0][0] = FFTWPREFIX(plan_dft_c2r_3d)(nx,ny,nz,(FFTWPREFIX(complex) *)up1,up2,fftw_flags);
- fft->plan[0][0][1] = FFTWPREFIX(plan_dft_r2c_3d)(nx,ny,nz,up1,(FFTWPREFIX(complex) *)up2,fftw_flags);
- fft->plan[0][1][0] = FFTWPREFIX(plan_dft_c2r_3d)(nx,ny,nz,(FFTWPREFIX(complex) *)up1,up1,fftw_flags);
- fft->plan[0][1][1] = FFTWPREFIX(plan_dft_r2c_3d)(nx,ny,nz,up1,(FFTWPREFIX(complex) *)up1,fftw_flags);
-
- fft->plan[1][0][0] = FFTWPREFIX(plan_dft_c2r_3d)(nx,ny,nz,(FFTWPREFIX(complex) *)p1,p2,fftw_flags);
- fft->plan[1][0][1] = FFTWPREFIX(plan_dft_r2c_3d)(nx,ny,nz,p1,(FFTWPREFIX(complex) *)p2,fftw_flags);
- fft->plan[1][1][0] = FFTWPREFIX(plan_dft_c2r_3d)(nx,ny,nz,(FFTWPREFIX(complex) *)p1,p1,fftw_flags);
- fft->plan[1][1][1] = FFTWPREFIX(plan_dft_r2c_3d)(nx,ny,nz,p1,(FFTWPREFIX(complex) *)p1,fftw_flags);
-
-
- for(i=0;i<2;i++)
- {
- for(j=0;j<2;j++)
- {
- for(k=0;k<2;k++)
- {
- if(fft->plan[i][j][k] == NULL)
- {
- gmx_fatal(FARGS,"Error initializing FFTW3 plan.");
- FFTW_UNLOCK;
- gmx_fft_destroy(fft);
- FFTW_LOCK;
- FFTWPREFIX(free)(p1);
- FFTWPREFIX(free)(p2);
- FFTW_UNLOCK;
- return -1;
- }
- }
- }
- }
-
- FFTWPREFIX(free)(p1);
- FFTWPREFIX(free)(p2);
-
- fft->real_transform = 1;
- fft->ndim = 3;
-
- *pfft = fft;
- FFTW_UNLOCK;
- return 0;
-}
-
-
int
gmx_fft_1d (gmx_fft_t fft,
enum gmx_fft_direction dir,
int aligned = ((((size_t)in_data | (size_t)out_data) & 0xf)==0);
int inplace = (in_data == out_data);
int isforward = (dir == GMX_FFT_FORWARD);
-
+
/* Some checks */
if( (fft->real_transform == 1) || (fft->ndim != 1) ||
((dir != GMX_FFT_FORWARD) && (dir != GMX_FFT_BACKWARD)) )
{
gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
return EINVAL;
- }
+ }
FFTWPREFIX(execute_dft)(fft->plan[aligned][inplace][isforward],
(FFTWPREFIX(complex) *)in_data,
(FFTWPREFIX(complex) *)out_data);
-
+
return 0;
}
int
gmx_fft_many_1d (gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
+ enum gmx_fft_direction dir,
+ void * in_data,
+ void * out_data)
{
return gmx_fft_1d(fft,dir,in_data,out_data);
}
-int
+int
gmx_fft_1d_real (gmx_fft_t fft,
enum gmx_fft_direction dir,
void * in_data,
return gmx_fft_1d_real(fft,dir,in_data,out_data);
}
-int
-gmx_fft_2d (gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int aligned = ((((size_t)in_data | (size_t)out_data) & 0xf)==0);
- int inplace = (in_data == out_data);
- int isforward = (dir == GMX_FFT_FORWARD);
-
- /* Some checks */
- if( (fft->real_transform == 1) || (fft->ndim != 2) ||
- ((dir != GMX_FFT_FORWARD) && (dir != GMX_FFT_BACKWARD)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- FFTWPREFIX(execute_dft)(fft->plan[aligned][inplace][isforward],
- (FFTWPREFIX(complex) *)in_data,
- (FFTWPREFIX(complex) *)out_data);
-
- return 0;
-}
-
-
int
gmx_fft_2d_real (gmx_fft_t fft,
enum gmx_fft_direction dir,
return 0;
}
-
-int
-gmx_fft_3d (gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int aligned = ((((size_t)in_data | (size_t)out_data) & 0xf)==0);
- int inplace = (in_data == out_data);
- int isforward = (dir == GMX_FFT_FORWARD);
-
- /* Some checks */
- if( (fft->real_transform == 1) || (fft->ndim != 3) ||
- ((dir != GMX_FFT_FORWARD) && (dir != GMX_FFT_BACKWARD)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- FFTWPREFIX(execute_dft)(fft->plan[aligned][inplace][isforward],
- (FFTWPREFIX(complex) *)in_data,
- (FFTWPREFIX(complex) *)out_data);
-
- return 0;
-}
-
-
-int
-gmx_fft_3d_real (gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int aligned = ((((size_t)in_data | (size_t)out_data) & 0xf)==0);
- int inplace = (in_data == out_data);
- int isforward = (dir == GMX_FFT_REAL_TO_COMPLEX);
-
- /* Some checks */
- if( (fft->real_transform != 1) || (fft->ndim != 3) ||
- ((dir != GMX_FFT_REAL_TO_COMPLEX) && (dir != GMX_FFT_COMPLEX_TO_REAL)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- if(isforward)
- {
- FFTWPREFIX(execute_dft_r2c)(fft->plan[aligned][inplace][isforward],
- (real *)in_data,
- (FFTWPREFIX(complex) *)out_data);
- }
- else
- {
- FFTWPREFIX(execute_dft_c2r)(fft->plan[aligned][inplace][isforward],
- (FFTWPREFIX(complex) *)in_data,
- (real *)out_data);
- }
-
-
- return 0;
-}
-
-
void
gmx_fft_destroy(gmx_fft_t fft)
{
gmx_fft_destroy(fft);
}
-#else
-int
-gmx_fft_fftw3_empty;
+void gmx_fft_cleanup()
+{
+ FFTWPREFIX(cleanup)();
+}
#endif /* GMX_FFT_FFTW3 */
#include <stdlib.h>
#include <mkl_dfti.h>
+#include <mkl_service.h>
#include "gmx_fft.h"
}
-
-int
-gmx_fft_init_2d(gmx_fft_t * pfft,
- int nx,
- int ny,
- gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- int d;
- int status;
- MKL_LONG length[2];
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- if( (fft = (gmx_fft_t)malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- /* Mark all handles invalid */
- for(d=0;d<3;d++)
- {
- fft->inplace[d] = fft->ooplace[d] = NULL;
- }
- fft->ooplace[3] = NULL;
-
- length[0] = nx;
- length[1] = ny;
-
- status = DftiCreateDescriptor(&fft->inplace[0],GMX_DFTI_PREC,DFTI_COMPLEX,2,length);
-
- if( status == 0 )
- status = DftiSetValue(fft->inplace[0],DFTI_PLACEMENT,DFTI_INPLACE);
-
- if( status == 0 )
- status = DftiCommitDescriptor(fft->inplace[0]);
-
-
- if( status == 0 )
- status = DftiCreateDescriptor(&fft->ooplace[0],GMX_DFTI_PREC,DFTI_COMPLEX,2,length);
-
- if( status == 0 )
- status = DftiSetValue(fft->ooplace[0],DFTI_PLACEMENT,DFTI_NOT_INPLACE);
-
- if( status == 0 )
- status = DftiCommitDescriptor(fft->ooplace[0]);
-
-
- if( status != 0 )
- {
- gmx_fatal(FARGS,"Error initializing Intel MKL FFT; status=%d",status);
- gmx_fft_destroy(fft);
- return status;
- }
-
- fft->ndim = 2;
- fft->nx = nx;
- fft->ny = ny;
- fft->real_fft = 0;
- fft->work = NULL;
-
- *pfft = fft;
- return 0;
-}
-
-
int
gmx_fft_init_2d_real(gmx_fft_t * pfft,
return 0;
}
-
-
-int
-gmx_fft_init_3d(gmx_fft_t * pfft,
- int nx,
- int ny,
- int nz,
- gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- int d;
- MKL_LONG length[3];
- int status;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- if( (fft = (gmx_fft_t)malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- /* Mark all handles invalid */
- for(d=0;d<3;d++)
- {
- fft->inplace[d] = fft->ooplace[d] = NULL;
- }
- fft->ooplace[3] = NULL;
-
- length[0] = nx;
- length[1] = ny;
- length[2] = nz;
-
- status = DftiCreateDescriptor(&fft->inplace[0],GMX_DFTI_PREC,DFTI_COMPLEX,(MKL_LONG)3,length);
-
- if( status == 0 )
- status = DftiSetValue(fft->inplace[0],DFTI_PLACEMENT,DFTI_INPLACE);
-
- if( status == 0 )
- status = DftiCommitDescriptor(fft->inplace[0]);
-
-
- if( status == 0 )
- status = DftiCreateDescriptor(&fft->ooplace[0],GMX_DFTI_PREC,DFTI_COMPLEX,(MKL_LONG)3,length);
-
- if( status == 0 )
- status = DftiSetValue(fft->ooplace[0],DFTI_PLACEMENT,DFTI_NOT_INPLACE);
-
- if( status == 0 )
- status = DftiCommitDescriptor(fft->ooplace[0]);
-
-
- if( status != 0 )
- {
- gmx_fatal(FARGS,"Error initializing Intel MKL FFT; status=%d",status);
- gmx_fft_destroy(fft);
- return status;
- }
-
-
- fft->ndim = 3;
- fft->nx = nx;
- fft->ny = ny;
- fft->nz = nz;
- fft->real_fft = 0;
- fft->work = NULL;
-
- *pfft = fft;
- return 0;
-}
-
-
-
-
-int
-gmx_fft_init_3d_real(gmx_fft_t * pfft,
- int nx,
- int ny,
- int nz,
- gmx_fft_flag flags)
-{
- gmx_fft_t fft;
- int d;
- int status;
- MKL_LONG stride[2];
- int nzc;
-
- if(pfft==NULL)
- {
- gmx_fatal(FARGS,"Invalid opaque FFT datatype pointer.");
- return EINVAL;
- }
- *pfft = NULL;
-
- nzc = (nz/2 + 1);
-
- if( (fft = (gmx_fft_t)malloc(sizeof(struct gmx_fft))) == NULL)
- {
- return ENOMEM;
- }
-
- /* Mark all handles invalid */
- for(d=0;d<3;d++)
- {
- fft->inplace[d] = fft->ooplace[d] = NULL;
- }
- fft->ooplace[3] = NULL;
-
- /* Roll our own 3D real transform using multiple transforms in MKL,
- * since the current MKL versions does not support our storage format
- * or 3D real transforms.
- */
-
- /* In-place X FFT.
- * ny*nzc complex-to-complex transforms, length nx
- * transform distance: 1
- * element strides: ny*nzc
- */
- status = DftiCreateDescriptor(&fft->inplace[0],GMX_DFTI_PREC,DFTI_COMPLEX,1,(MKL_LONG)nx);
-
- if ( status == 0)
- {
- stride[0] = 0;
- stride[1] = ny*nzc;
-
- status =
- (DftiSetValue(fft->inplace[0],DFTI_PLACEMENT,DFTI_INPLACE) ||
- DftiSetValue(fft->inplace[0],DFTI_NUMBER_OF_TRANSFORMS,(MKL_LONG)ny*nzc) ||
- DftiSetValue(fft->inplace[0],DFTI_INPUT_DISTANCE,1) ||
- DftiSetValue(fft->inplace[0],DFTI_INPUT_STRIDES,stride) ||
- DftiSetValue(fft->inplace[0],DFTI_OUTPUT_DISTANCE,1) ||
- DftiSetValue(fft->inplace[0],DFTI_OUTPUT_STRIDES,stride) ||
- DftiCommitDescriptor(fft->inplace[0]));
- }
-
- /* Out-of-place X FFT:
- * ny*nzc complex-to-complex transforms, length nx
- * transform distance: 1
- * element strides: ny*nzc
- */
- if( status == 0 )
- status = DftiCreateDescriptor(&fft->ooplace[0],GMX_DFTI_PREC,DFTI_COMPLEX,1,(MKL_LONG)nx);
-
- if( status == 0 )
- {
- stride[0] = 0;
- stride[1] = ny*nzc;
-
- status =
- (DftiSetValue(fft->ooplace[0],DFTI_PLACEMENT,DFTI_NOT_INPLACE) ||
- DftiSetValue(fft->ooplace[0],DFTI_NUMBER_OF_TRANSFORMS,(MKL_LONG)ny*nzc) ||
- DftiSetValue(fft->ooplace[0],DFTI_INPUT_DISTANCE,1) ||
- DftiSetValue(fft->ooplace[0],DFTI_INPUT_STRIDES,stride) ||
- DftiSetValue(fft->ooplace[0],DFTI_OUTPUT_DISTANCE,1) ||
- DftiSetValue(fft->ooplace[0],DFTI_OUTPUT_STRIDES,stride) ||
- DftiCommitDescriptor(fft->ooplace[0]));
- }
-
-
- /* In-place Y FFT.
- * We cannot do all NX*NZC transforms at once, so define a handle to do
- * NZC transforms, and then execute it NX times.
- * nzc complex-to-complex transforms, length ny
- * transform distance: 1
- * element strides: nzc
- */
- if( status == 0 )
- status = DftiCreateDescriptor(&fft->inplace[1],GMX_DFTI_PREC,DFTI_COMPLEX,1,(MKL_LONG)ny);
-
- if( status == 0 )
- {
- stride[0] = 0;
- stride[1] = nzc;
-
- status =
- (DftiSetValue(fft->inplace[1],DFTI_PLACEMENT,DFTI_INPLACE) ||
- DftiSetValue(fft->inplace[1],DFTI_NUMBER_OF_TRANSFORMS,(MKL_LONG)nzc) ||
- DftiSetValue(fft->inplace[1],DFTI_INPUT_DISTANCE,1) ||
- DftiSetValue(fft->inplace[1],DFTI_INPUT_STRIDES,stride) ||
- DftiSetValue(fft->inplace[1],DFTI_OUTPUT_DISTANCE,1) ||
- DftiSetValue(fft->inplace[1],DFTI_OUTPUT_STRIDES,stride) ||
- DftiCommitDescriptor(fft->inplace[1]));
- }
-
-
- /* Out-of-place Y FFT:
- * We cannot do all NX*NZC transforms at once, so define a handle to do
- * NZC transforms, and then execute it NX times.
- * nzc complex-to-complex transforms, length ny
- * transform distance: 1
- * element strides: nzc
- */
- if( status == 0 )
- status = DftiCreateDescriptor(&fft->ooplace[1],GMX_DFTI_PREC,DFTI_COMPLEX,1,(MKL_LONG)ny);
-
- if( status == 0 )
- {
- stride[0] = 0;
- stride[1] = nzc;
-
- status =
- (DftiSetValue(fft->ooplace[1],DFTI_PLACEMENT,DFTI_NOT_INPLACE) ||
- DftiSetValue(fft->ooplace[1],DFTI_NUMBER_OF_TRANSFORMS,(MKL_LONG)nzc) ||
- DftiSetValue(fft->ooplace[1],DFTI_INPUT_DISTANCE,1) ||
- DftiSetValue(fft->ooplace[1],DFTI_INPUT_STRIDES,stride) ||
- DftiSetValue(fft->ooplace[1],DFTI_OUTPUT_DISTANCE,1) ||
- DftiSetValue(fft->ooplace[1],DFTI_OUTPUT_STRIDES,stride) ||
- DftiCommitDescriptor(fft->ooplace[1]));
- }
-
- /* In-place Z FFT:
- * nx*ny real-to-complex transforms, length nz
- * transform distance: nzc*2 -> nzc*2
- * element strides: 1
- */
- if( status == 0 )
- status = DftiCreateDescriptor(&fft->inplace[2],GMX_DFTI_PREC,DFTI_REAL,1,(MKL_LONG)nz);
-
- if( status == 0 )
- {
- stride[0] = 0;
- stride[1] = 1;
-
- status =
- (DftiSetValue(fft->inplace[2],DFTI_PLACEMENT,DFTI_INPLACE) ||
- DftiSetValue(fft->inplace[2],DFTI_NUMBER_OF_TRANSFORMS,(MKL_LONG)nx*ny) ||
- DftiSetValue(fft->inplace[2],DFTI_INPUT_DISTANCE,(MKL_LONG)nzc*2) ||
- DftiSetValue(fft->inplace[2],DFTI_INPUT_STRIDES,stride) ||
- DftiSetValue(fft->inplace[2],DFTI_OUTPUT_DISTANCE,(MKL_LONG)nzc*2) ||
- DftiSetValue(fft->inplace[2],DFTI_OUTPUT_STRIDES,stride) ||
- DftiCommitDescriptor(fft->inplace[2]));
- }
-
-
- /* Out-of-place real-to-complex (affects distance) Z FFT:
- * nx*ny real-to-complex transforms, length nz
- * transform distance: nz -> nzc*2
- * element STRIDES: 1
- */
- if( status == 0 )
- status = DftiCreateDescriptor(&fft->ooplace[2],GMX_DFTI_PREC,DFTI_REAL,1,(MKL_LONG)nz);
-
- if( status == 0 )
- {
- stride[0] = 0;
- stride[1] = 1;
-
- status =
- (DftiSetValue(fft->ooplace[2],DFTI_PLACEMENT,DFTI_NOT_INPLACE) ||
- DftiSetValue(fft->ooplace[2],DFTI_NUMBER_OF_TRANSFORMS,(MKL_LONG)nx*ny) ||
- DftiSetValue(fft->ooplace[2],DFTI_INPUT_DISTANCE,(MKL_LONG)nz) ||
- DftiSetValue(fft->ooplace[2],DFTI_INPUT_STRIDES,stride) ||
- DftiSetValue(fft->ooplace[2],DFTI_OUTPUT_DISTANCE,(MKL_LONG)nzc*2) ||
- DftiSetValue(fft->ooplace[2],DFTI_OUTPUT_STRIDES,stride) ||
- DftiCommitDescriptor(fft->ooplace[2]));
- }
-
-
- /* Out-of-place complex-to-real (affects distance) Z FFT:
- * nx*ny real-to-complex transforms, length nz
- * transform distance: nzc*2 -> nz
- * element STRIDES: 1
- */
- if( status == 0 )
- status = DftiCreateDescriptor(&fft->ooplace[3],GMX_DFTI_PREC,DFTI_REAL,1,(MKL_LONG)nz);
-
- if( status == 0 )
- {
- stride[0] = 0;
- stride[1] = 1;
-
- status =
- (DftiSetValue(fft->ooplace[3],DFTI_PLACEMENT,DFTI_NOT_INPLACE) ||
- DftiSetValue(fft->ooplace[3],DFTI_NUMBER_OF_TRANSFORMS,(MKL_LONG)nx*ny) ||
- DftiSetValue(fft->ooplace[3],DFTI_INPUT_DISTANCE,(MKL_LONG)nzc*2) ||
- DftiSetValue(fft->ooplace[3],DFTI_INPUT_STRIDES,stride) ||
- DftiSetValue(fft->ooplace[3],DFTI_OUTPUT_DISTANCE,(MKL_LONG)nz) ||
- DftiSetValue(fft->ooplace[3],DFTI_OUTPUT_STRIDES,stride) ||
- DftiCommitDescriptor(fft->ooplace[3]));
- }
-
-
- if ( status == 0 )
- {
- if ((fft->work = (t_complex *)malloc(sizeof(t_complex)*(nx*ny*(nz/2+1)))) == NULL)
- {
- status = ENOMEM;
- }
- }
-
-
- if( status != 0 )
- {
- gmx_fatal(FARGS,"Error initializing Intel MKL FFT; status=%d",status);
- gmx_fft_destroy(fft);
- return status;
- }
-
-
- fft->ndim = 3;
- fft->nx = nx;
- fft->ny = ny;
- fft->nz = nz;
- fft->real_fft = 1;
-
- *pfft = fft;
- return 0;
-}
-
-
-
-
int
gmx_fft_1d(gmx_fft_t fft,
enum gmx_fft_direction dir,
}
-int
-gmx_fft_2d(gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int inplace = (in_data == out_data);
- int status = 0;
-
- if( (fft->real_fft == 1) || (fft->ndim != 2) ||
- ((dir != GMX_FFT_FORWARD) && (dir != GMX_FFT_BACKWARD)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- if(dir==GMX_FFT_FORWARD)
- {
- if(inplace)
- {
- status = DftiComputeForward(fft->inplace[0],in_data);
- }
- else
- {
- status = DftiComputeForward(fft->ooplace[0],in_data,out_data);
- }
- }
- else
- {
- if(inplace)
- {
- status = DftiComputeBackward(fft->inplace[0],in_data);
- }
- else
- {
- status = DftiComputeBackward(fft->ooplace[0],in_data,out_data);
- }
- }
-
- if( status != 0 )
- {
- gmx_fatal(FARGS,"Error executing Intel MKL FFT.");
- status = -1;
- }
-
- return status;
-}
-
-
int
gmx_fft_2d_real(gmx_fft_t fft,
enum gmx_fft_direction dir,
}
else
{
- if(inplace)
- {
- /* complex-to-complex in X dimension, in-place */
- status = DftiComputeBackward(fft->inplace[0],in_data);
-
- /* complex-to-real in Y dimension, in-place */
- if ( status == 0 )
- status = DftiComputeBackward(fft->inplace[1],in_data);
-
- }
- else
- {
- /* complex-to-complex in X dimension, from in_data to work */
- status = DftiComputeBackward(fft->ooplace[0],in_data,fft->work);
-
- /* complex-to-real in Y dimension, from work to out_data */
- if ( status == 0 )
- status = DftiComputeBackward(fft->ooplace[1],fft->work,out_data);
-
- }
- }
-
- if( status != 0 )
- {
- gmx_fatal(FARGS,"Error executing Intel MKL FFT.");
- status = -1;
- }
-
- return status;
-}
-
-
-int
-gmx_fft_3d(gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int inplace = (in_data == out_data);
- int status = 0;
-
- if( (fft->real_fft == 1) || (fft->ndim != 3) ||
- ((dir != GMX_FFT_FORWARD) && (dir != GMX_FFT_BACKWARD)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- if(dir==GMX_FFT_FORWARD)
- {
- if(inplace)
- {
- status = DftiComputeForward(fft->inplace[0],in_data);
- }
- else
- {
- status = DftiComputeForward(fft->ooplace[0],in_data,out_data);
- }
- }
- else
- {
- if(inplace)
- {
- status = DftiComputeBackward(fft->inplace[0],in_data);
- }
- else
- {
- status = DftiComputeBackward(fft->ooplace[0],in_data,out_data);
- }
- }
-
- if( status != 0 )
- {
- gmx_fatal(FARGS,"Error executing Intel MKL FFT.");
- status = -1;
- }
-
- return status;
-}
-
-
-int
-gmx_fft_3d_real(gmx_fft_t fft,
- enum gmx_fft_direction dir,
- void * in_data,
- void * out_data)
-{
- int inplace = (in_data == out_data);
- int status = 0;
- int i;
- int nx,ny,nzc;
-
- nx = fft->nx;
- ny = fft->ny;
- nzc = fft->nz/2 + 1;
-
- if( (fft->real_fft != 1) || (fft->ndim != 3) ||
- ((dir != GMX_FFT_REAL_TO_COMPLEX) && (dir != GMX_FFT_COMPLEX_TO_REAL)) )
- {
- gmx_fatal(FARGS,"FFT plan mismatch - bad plan or direction.");
- return EINVAL;
- }
-
- if(dir==GMX_FFT_REAL_TO_COMPLEX)
- {
- if(inplace)
- {
- /* real-to-complex in Z dimension, in-place */
- status = DftiComputeForward(fft->inplace[2],in_data);
-
- /* complex-to-complex in Y dimension, in-place */
- for(i=0;i<nx;i++)
- {
- if ( status == 0 )
- status = DftiComputeForward(fft->inplace[1],(t_complex *)in_data+i*ny*nzc);
- }
-
- /* complex-to-complex in X dimension, in-place */
- if ( status == 0 )
- status = DftiComputeForward(fft->inplace[0],in_data);
- }
- else
- {
- /* real-to-complex in Z dimension, from in_data to out_data */
- status = DftiComputeForward(fft->ooplace[2],in_data,out_data);
-
- /* complex-to-complex in Y dimension, in-place */
- for(i=0;i<nx;i++)
- {
- if ( status == 0 )
- status = DftiComputeForward(fft->inplace[1],(t_complex *)out_data+i*ny*nzc);
- }
-
- /* complex-to-complex in X dimension, in-place */
- if ( status == 0 )
- status = DftiComputeForward(fft->inplace[0],out_data);
- }
- }
- else
- {
- if(inplace)
- {
- /* complex-to-complex in X dimension, in-place */
- status = DftiComputeBackward(fft->inplace[0],in_data);
-
- /* complex-to-complex in Y dimension, in-place */
- for(i=0;i<nx;i++)
- {
- if ( status == 0 )
- status = DftiComputeBackward(fft->inplace[1],(t_complex *)in_data+i*ny*nzc);
- }
-
- /* complex-to-real in Z dimension, in-place */
- if ( status == 0 )
- status = DftiComputeBackward(fft->inplace[2],in_data);
- }
- else
- {
- /* complex-to-complex in X dimension, from in_data to work */
- status = DftiComputeBackward(fft->ooplace[0],in_data,fft->work);
-
- /* complex-to-complex in Y dimension, in-place */
- for(i=0;i<nx;i++)
- {
- if ( status == 0 )
- status = DftiComputeBackward(fft->inplace[1],fft->work+i*ny*nzc);
- }
-
- /* complex-to-real in Z dimension, work to out_data */
- if ( status == 0 )
- status = DftiComputeBackward(fft->ooplace[2],fft->work,out_data);
- }
+ /* prior implementation was incorrect. See fft.cpp unit test */
+ gmx_incons("Complex -> Real is not supported by MKL.");
}
if( status != 0 )
return status;
}
-
-
void
gmx_fft_destroy(gmx_fft_t fft)
{
{
DftiFreeDescriptor(&fft->ooplace[3]);
}
+ if (fft->work != NULL)
+ {
+ free(fft->work);
+ }
free(fft);
}
}
+void gmx_fft_cleanup()
+{
+ mkl_free_buffers();
+}
+
#else
int
gmx_fft_mkl_empty;
ivec local_offset,
ivec local_size)
{
- int N1,M0,K0,K1,*coor;
- fft5d_local_size(p,&N1,&M0,&K0,&K1,&coor); /* M0=MG/P[0], K1=KG/P[1], NG,MG,KG global sizes */
-
local_offset[2]=0;
local_offset[1]=p->oM[0]; /*=p->coor[0]*p->MG/p->P[0]; */
local_offset[0]=p->oK[0]; /*=p->coor[1]*p->KG/p->P[1]; */
local_ndata[0]=p->pK[0];
if ((!(p->flags&FFT5D_BACKWARD)) && (p->flags&FFT5D_REALCOMPLEX)) {
+ //C is length in multiples of complex local_size in multiples of real
local_size[2]=p->C[0]*2;
} else {
local_size[2]=p->C[0];
int
gmx_parallel_3dfft_destroy(gmx_parallel_3dfft_t pfft_setup) {
- fft5d_destroy(pfft_setup->p2);
- fft5d_destroy(pfft_setup->p1);
- sfree(pfft_setup);
+ if (pfft_setup)
+ {
+ fft5d_destroy(pfft_setup->p2);
+ fft5d_destroy(pfft_setup->p1);
+ sfree(pfft_setup);
+ }
return 0;
}
// the declarations are still included for clarity.
virtual const char *typeString() const = 0;
virtual int valueCount() const { return static_cast<int>(values_->size()); }
- /*! \copydoc AbstractOptionStorage::formatValue()
+ /*! \copydoc gmx::AbstractOptionStorage::formatValue()
*
* OptionStorageTemplate implements handling of DefaultValueIfSetIndex
* in this method, as well as checking that \p i is a valid index.
#pragma warning(disable: 4065)
#endif
-/* These macros should be used at the beginning and end of each semantic action
+/*! \name Exception handling macros for actions
+ *
+ * These macros should be used at the beginning and end of each semantic action
* that may throw an exception. For robustness, it's best to wrap all actions
* that call functions declared outside parser.y should be wrapped.
* These macros take care to catch any exceptions, store the exception (or
* cleanly if necessary.
* The code calling the parser should use
* _gmx_sel_lexer_rethrow_exception_if_occurred() to rethrow any exceptions.
+ * \{
*/
+//! Starts an action that may throw exceptions.
#define BEGIN_ACTION \
try {
+//! Finishes an action that may throw exceptions.
#define END_ACTION \
} \
catch(const std::exception &ex) \
else \
YYABORT; \
}
+//!\}
/* Line 268 of yacc.c */
-#line 141 "parser.cpp"
+#line 147 "parser.cpp"
/* Enabling traces. */
#ifndef YYDEBUG
{
/* Line 293 of yacc.c */
-#line 95 "parser.y"
+#line 101 "parser.y"
int i;
real r;
/* Line 293 of yacc.c */
-#line 226 "parser.cpp"
+#line 232 "parser.cpp"
} YYSTYPE;
# define YYSTYPE_IS_TRIVIAL 1
# define yystype YYSTYPE /* obsolescent; will be withdrawn */
/* Line 343 of yacc.c */
-#line 262 "parser.cpp"
+#line 268 "parser.cpp"
#ifdef short
# undef short
/* YYRLINE[YYN] -- source line where rule number YYN was defined. */
static const yytype_uint16 yyrline[] =
{
- 0, 211, 211, 212, 223, 224, 246, 252, 253, 264,
- 276, 282, 288, 294, 300, 310, 318, 319, 329, 330,
- 337, 338, 352, 353, 357, 358, 361, 362, 365, 366,
- 374, 382, 390, 398, 402, 412, 420, 430, 431, 435,
- 442, 449, 459, 473, 482, 494, 501, 511, 517, 523,
- 529, 535, 541, 547, 554, 563, 577, 586, 590, 600,
- 613, 621, 629, 642, 644, 649, 650, 655, 664, 665,
- 666, 670, 671, 673, 678, 679, 683, 684, 686, 690,
- 696, 702, 708, 714, 718, 725, 732, 739, 743, 750,
- 757
+ 0, 217, 217, 218, 229, 230, 252, 258, 259, 270,
+ 282, 288, 294, 300, 306, 316, 324, 325, 335, 336,
+ 343, 344, 358, 359, 363, 364, 367, 368, 371, 372,
+ 380, 388, 396, 404, 408, 418, 426, 436, 437, 441,
+ 448, 455, 465, 479, 488, 500, 507, 517, 523, 529,
+ 535, 541, 547, 553, 560, 569, 583, 592, 596, 606,
+ 619, 627, 635, 648, 650, 655, 656, 661, 670, 671,
+ 672, 676, 677, 679, 684, 685, 689, 690, 692, 696,
+ 702, 708, 714, 720, 724, 731, 738, 745, 749, 756,
+ 763
};
#endif
case 5: /* "HELP_TOPIC" */
/* Line 1391 of yacc.c */
-#line 190 "parser.y"
+#line 196 "parser.y"
{ free((yyvaluep->str)); };
/* Line 1391 of yacc.c */
-#line 1368 "parser.cpp"
+#line 1374 "parser.cpp"
break;
case 8: /* "STR" */
/* Line 1391 of yacc.c */
-#line 190 "parser.y"
+#line 196 "parser.y"
{ free((yyvaluep->str)); };
/* Line 1391 of yacc.c */
-#line 1377 "parser.cpp"
+#line 1383 "parser.cpp"
break;
case 9: /* "IDENTIFIER" */
/* Line 1391 of yacc.c */
-#line 190 "parser.y"
+#line 196 "parser.y"
{ free((yyvaluep->str)); };
/* Line 1391 of yacc.c */
-#line 1386 "parser.cpp"
+#line 1392 "parser.cpp"
break;
case 25: /* "PARAM" */
/* Line 1391 of yacc.c */
-#line 191 "parser.y"
+#line 197 "parser.y"
{ if((yyvaluep->str)) free((yyvaluep->str)); };
/* Line 1391 of yacc.c */
-#line 1395 "parser.cpp"
+#line 1401 "parser.cpp"
break;
case 28: /* "CMP_OP" */
/* Line 1391 of yacc.c */
-#line 190 "parser.y"
+#line 196 "parser.y"
{ free((yyvaluep->str)); };
/* Line 1391 of yacc.c */
-#line 1404 "parser.cpp"
+#line 1410 "parser.cpp"
break;
case 51: /* "command" */
/* Line 1391 of yacc.c */
-#line 192 "parser.y"
+#line 198 "parser.y"
{ if((yyvaluep->sel)) _gmx_selelem_free((yyvaluep->sel)); };
/* Line 1391 of yacc.c */
-#line 1413 "parser.cpp"
+#line 1419 "parser.cpp"
break;
case 52: /* "cmd_plain" */
/* Line 1391 of yacc.c */
-#line 192 "parser.y"
+#line 198 "parser.y"
{ if((yyvaluep->sel)) _gmx_selelem_free((yyvaluep->sel)); };
/* Line 1391 of yacc.c */
-#line 1422 "parser.cpp"
+#line 1428 "parser.cpp"
break;
case 54: /* "help_topic" */
/* Line 1391 of yacc.c */
-#line 198 "parser.y"
+#line 204 "parser.y"
{ _gmx_selexpr_free_values((yyvaluep->val)); };
/* Line 1391 of yacc.c */
-#line 1431 "parser.cpp"
+#line 1437 "parser.cpp"
break;
case 55: /* "selection" */
/* Line 1391 of yacc.c */
-#line 193 "parser.y"
+#line 199 "parser.y"
{ _gmx_selelem_free_chain((yyvaluep->sel)); };
/* Line 1391 of yacc.c */
-#line 1440 "parser.cpp"
+#line 1446 "parser.cpp"
break;
case 59: /* "string" */
/* Line 1391 of yacc.c */
-#line 190 "parser.y"
+#line 196 "parser.y"
{ free((yyvaluep->str)); };
/* Line 1391 of yacc.c */
-#line 1449 "parser.cpp"
+#line 1455 "parser.cpp"
break;
case 60: /* "sel_expr" */
/* Line 1391 of yacc.c */
-#line 194 "parser.y"
+#line 200 "parser.y"
{ _gmx_selelem_free((yyvaluep->sel)); };
/* Line 1391 of yacc.c */
-#line 1458 "parser.cpp"
+#line 1464 "parser.cpp"
break;
case 62: /* "num_expr" */
/* Line 1391 of yacc.c */
-#line 194 "parser.y"
+#line 200 "parser.y"
{ _gmx_selelem_free((yyvaluep->sel)); };
/* Line 1391 of yacc.c */
-#line 1467 "parser.cpp"
+#line 1473 "parser.cpp"
break;
case 63: /* "str_expr" */
/* Line 1391 of yacc.c */
-#line 194 "parser.y"
+#line 200 "parser.y"
{ _gmx_selelem_free((yyvaluep->sel)); };
/* Line 1391 of yacc.c */
-#line 1476 "parser.cpp"
+#line 1482 "parser.cpp"
break;
case 64: /* "pos_expr" */
/* Line 1391 of yacc.c */
-#line 194 "parser.y"
+#line 200 "parser.y"
{ _gmx_selelem_free((yyvaluep->sel)); };
/* Line 1391 of yacc.c */
-#line 1485 "parser.cpp"
+#line 1491 "parser.cpp"
break;
case 65: /* "method_params" */
/* Line 1391 of yacc.c */
-#line 195 "parser.y"
+#line 201 "parser.y"
{ _gmx_selexpr_free_params((yyvaluep->param)); };
/* Line 1391 of yacc.c */
-#line 1494 "parser.cpp"
+#line 1500 "parser.cpp"
break;
case 66: /* "method_param_list" */
/* Line 1391 of yacc.c */
-#line 195 "parser.y"
+#line 201 "parser.y"
{ _gmx_selexpr_free_params((yyvaluep->param)); };
/* Line 1391 of yacc.c */
-#line 1503 "parser.cpp"
+#line 1509 "parser.cpp"
break;
case 67: /* "method_param" */
/* Line 1391 of yacc.c */
-#line 195 "parser.y"
+#line 201 "parser.y"
{ _gmx_selexpr_free_params((yyvaluep->param)); };
/* Line 1391 of yacc.c */
-#line 1512 "parser.cpp"
+#line 1518 "parser.cpp"
break;
case 68: /* "value_list" */
/* Line 1391 of yacc.c */
-#line 196 "parser.y"
+#line 202 "parser.y"
{ _gmx_selexpr_free_values((yyvaluep->val)); };
/* Line 1391 of yacc.c */
-#line 1521 "parser.cpp"
+#line 1527 "parser.cpp"
break;
case 69: /* "value_list_contents" */
/* Line 1391 of yacc.c */
-#line 196 "parser.y"
+#line 202 "parser.y"
{ _gmx_selexpr_free_values((yyvaluep->val)); };
/* Line 1391 of yacc.c */
-#line 1530 "parser.cpp"
+#line 1536 "parser.cpp"
break;
case 70: /* "basic_value_list" */
/* Line 1391 of yacc.c */
-#line 197 "parser.y"
+#line 203 "parser.y"
{ _gmx_selexpr_free_values((yyvaluep->val)); };
/* Line 1391 of yacc.c */
-#line 1539 "parser.cpp"
+#line 1545 "parser.cpp"
break;
case 71: /* "basic_value_list_contents" */
/* Line 1391 of yacc.c */
-#line 197 "parser.y"
+#line 203 "parser.y"
{ _gmx_selexpr_free_values((yyvaluep->val)); };
/* Line 1391 of yacc.c */
-#line 1548 "parser.cpp"
+#line 1554 "parser.cpp"
break;
case 72: /* "value_item" */
/* Line 1391 of yacc.c */
-#line 196 "parser.y"
+#line 202 "parser.y"
{ _gmx_selexpr_free_values((yyvaluep->val)); };
/* Line 1391 of yacc.c */
-#line 1557 "parser.cpp"
+#line 1563 "parser.cpp"
break;
case 73: /* "basic_value_item" */
/* Line 1391 of yacc.c */
-#line 197 "parser.y"
+#line 203 "parser.y"
{ _gmx_selexpr_free_values((yyvaluep->val)); };
/* Line 1391 of yacc.c */
-#line 1566 "parser.cpp"
+#line 1572 "parser.cpp"
break;
case 74: /* "value_item_range" */
/* Line 1391 of yacc.c */
-#line 196 "parser.y"
+#line 202 "parser.y"
{ _gmx_selexpr_free_values((yyvaluep->val)); };
/* Line 1391 of yacc.c */
-#line 1575 "parser.cpp"
+#line 1581 "parser.cpp"
break;
default:
case 2:
/* Line 1806 of yacc.c */
-#line 211 "parser.y"
+#line 217 "parser.y"
{ (yyval.sel) = NULL; }
break;
case 3:
/* Line 1806 of yacc.c */
-#line 213 "parser.y"
+#line 219 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_append_selection((yyvsp[(2) - (2)].sel), (yyvsp[(1) - (2)].sel), scanner);
case 4:
/* Line 1806 of yacc.c */
-#line 223 "parser.y"
+#line 229 "parser.y"
{ (yyval.sel) = (yyvsp[(1) - (2)].sel); }
break;
case 5:
/* Line 1806 of yacc.c */
-#line 225 "parser.y"
+#line 231 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = NULL;
case 6:
/* Line 1806 of yacc.c */
-#line 246 "parser.y"
+#line 252 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = NULL;
case 7:
/* Line 1806 of yacc.c */
-#line 252 "parser.y"
+#line 258 "parser.y"
{ (yyval.sel) = NULL; }
break;
case 8:
/* Line 1806 of yacc.c */
-#line 254 "parser.y"
+#line 260 "parser.y"
{
BEGIN_ACTION;
t_selelem *s, *p;
case 9:
/* Line 1806 of yacc.c */
-#line 265 "parser.y"
+#line 271 "parser.y"
{
BEGIN_ACTION;
t_selelem *s, *p;
case 10:
/* Line 1806 of yacc.c */
-#line 277 "parser.y"
+#line 283 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_selection(NULL, (yyvsp[(1) - (1)].sel), scanner);
case 11:
/* Line 1806 of yacc.c */
-#line 283 "parser.y"
+#line 289 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_selection((yyvsp[(1) - (2)].str), (yyvsp[(2) - (2)].sel), scanner);
case 12:
/* Line 1806 of yacc.c */
-#line 289 "parser.y"
+#line 295 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_assign_variable((yyvsp[(1) - (3)].str), (yyvsp[(3) - (3)].sel), scanner);
case 13:
/* Line 1806 of yacc.c */
-#line 295 "parser.y"
+#line 301 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_assign_variable((yyvsp[(1) - (3)].str), (yyvsp[(3) - (3)].sel), scanner);
case 14:
/* Line 1806 of yacc.c */
-#line 301 "parser.y"
+#line 307 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_assign_variable((yyvsp[(1) - (3)].str), (yyvsp[(3) - (3)].sel), scanner);
case 15:
/* Line 1806 of yacc.c */
-#line 311 "parser.y"
+#line 317 "parser.y"
{
BEGIN_ACTION;
_gmx_sel_handle_help_cmd(process_value_list((yyvsp[(2) - (2)].val), NULL), scanner);
case 16:
/* Line 1806 of yacc.c */
-#line 318 "parser.y"
+#line 324 "parser.y"
{ (yyval.val) = NULL; }
break;
case 17:
/* Line 1806 of yacc.c */
-#line 320 "parser.y"
+#line 326 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value(STR_VALUE);
case 18:
/* Line 1806 of yacc.c */
-#line 329 "parser.y"
+#line 335 "parser.y"
{ (yyval.sel) = (yyvsp[(1) - (1)].sel); }
break;
case 19:
/* Line 1806 of yacc.c */
-#line 331 "parser.y"
+#line 337 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_position((yyvsp[(1) - (1)].sel), NULL, scanner);
case 20:
/* Line 1806 of yacc.c */
-#line 337 "parser.y"
+#line 343 "parser.y"
{ (yyval.sel) = (yyvsp[(2) - (3)].sel); }
break;
case 21:
/* Line 1806 of yacc.c */
-#line 339 "parser.y"
+#line 345 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_modifier((yyvsp[(2) - (3)].meth), (yyvsp[(3) - (3)].param), (yyvsp[(1) - (3)].sel), scanner);
case 22:
/* Line 1806 of yacc.c */
-#line 352 "parser.y"
+#line 358 "parser.y"
{ (yyval.i) = (yyvsp[(1) - (1)].i); }
break;
case 23:
/* Line 1806 of yacc.c */
-#line 353 "parser.y"
+#line 359 "parser.y"
{ (yyval.i) = -(yyvsp[(2) - (2)].i); }
break;
case 24:
/* Line 1806 of yacc.c */
-#line 357 "parser.y"
+#line 363 "parser.y"
{ (yyval.r) = (yyvsp[(1) - (1)].r); }
break;
case 25:
/* Line 1806 of yacc.c */
-#line 358 "parser.y"
+#line 364 "parser.y"
{ (yyval.r) = -(yyvsp[(2) - (2)].r); }
break;
case 26:
/* Line 1806 of yacc.c */
-#line 361 "parser.y"
+#line 367 "parser.y"
{ (yyval.r) = (yyvsp[(1) - (1)].i); }
break;
case 27:
/* Line 1806 of yacc.c */
-#line 362 "parser.y"
+#line 368 "parser.y"
{ (yyval.r) = (yyvsp[(1) - (1)].r); }
break;
case 28:
/* Line 1806 of yacc.c */
-#line 365 "parser.y"
+#line 371 "parser.y"
{ (yyval.str) = (yyvsp[(1) - (1)].str); }
break;
case 29:
/* Line 1806 of yacc.c */
-#line 366 "parser.y"
+#line 372 "parser.y"
{ (yyval.str) = (yyvsp[(1) - (1)].str); }
break;
case 30:
/* Line 1806 of yacc.c */
-#line 375 "parser.y"
+#line 381 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_selelem_create(SEL_BOOLEAN);
case 31:
/* Line 1806 of yacc.c */
-#line 383 "parser.y"
+#line 389 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_selelem_create(SEL_BOOLEAN);
case 32:
/* Line 1806 of yacc.c */
-#line 391 "parser.y"
+#line 397 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_selelem_create(SEL_BOOLEAN);
case 33:
/* Line 1806 of yacc.c */
-#line 398 "parser.y"
+#line 404 "parser.y"
{ (yyval.sel) = (yyvsp[(2) - (3)].sel); }
break;
case 34:
/* Line 1806 of yacc.c */
-#line 403 "parser.y"
+#line 409 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_comparison((yyvsp[(1) - (3)].sel), (yyvsp[(3) - (3)].sel), (yyvsp[(2) - (3)].str), scanner);
case 35:
/* Line 1806 of yacc.c */
-#line 413 "parser.y"
+#line 419 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_group_by_name((yyvsp[(2) - (2)].str), scanner);
case 36:
/* Line 1806 of yacc.c */
-#line 421 "parser.y"
+#line 427 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_group_by_id((yyvsp[(2) - (2)].i), scanner);
case 37:
/* Line 1806 of yacc.c */
-#line 430 "parser.y"
+#line 436 "parser.y"
{ (yyval.str) = NULL; }
break;
case 38:
/* Line 1806 of yacc.c */
-#line 431 "parser.y"
+#line 437 "parser.y"
{ (yyval.str) = (yyvsp[(1) - (1)].str); }
break;
case 39:
/* Line 1806 of yacc.c */
-#line 436 "parser.y"
+#line 442 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_keyword((yyvsp[(2) - (2)].meth), NULL, (yyvsp[(1) - (2)].str), scanner);
case 40:
/* Line 1806 of yacc.c */
-#line 443 "parser.y"
+#line 449 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_keyword((yyvsp[(2) - (3)].meth), process_value_list((yyvsp[(3) - (3)].val), NULL), (yyvsp[(1) - (3)].str), scanner);
case 41:
/* Line 1806 of yacc.c */
-#line 450 "parser.y"
+#line 456 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_keyword((yyvsp[(2) - (3)].meth), process_value_list((yyvsp[(3) - (3)].val), NULL), (yyvsp[(1) - (3)].str), scanner);
case 42:
/* Line 1806 of yacc.c */
-#line 460 "parser.y"
+#line 466 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_method((yyvsp[(2) - (3)].meth), (yyvsp[(3) - (3)].param), (yyvsp[(1) - (3)].str), scanner);
case 43:
/* Line 1806 of yacc.c */
-#line 474 "parser.y"
+#line 480 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_selelem_create(SEL_CONST);
case 44:
/* Line 1806 of yacc.c */
-#line 483 "parser.y"
+#line 489 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_selelem_create(SEL_CONST);
case 45:
/* Line 1806 of yacc.c */
-#line 495 "parser.y"
+#line 501 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_keyword((yyvsp[(2) - (2)].meth), NULL, (yyvsp[(1) - (2)].str), scanner);
case 46:
/* Line 1806 of yacc.c */
-#line 502 "parser.y"
+#line 508 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_method((yyvsp[(2) - (3)].meth), (yyvsp[(3) - (3)].param), (yyvsp[(1) - (3)].str), scanner);
case 47:
/* Line 1806 of yacc.c */
-#line 512 "parser.y"
+#line 518 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_arithmetic((yyvsp[(1) - (3)].sel), (yyvsp[(3) - (3)].sel), '+', scanner);
case 48:
/* Line 1806 of yacc.c */
-#line 518 "parser.y"
+#line 524 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_arithmetic((yyvsp[(1) - (3)].sel), (yyvsp[(3) - (3)].sel), '-', scanner);
case 49:
/* Line 1806 of yacc.c */
-#line 524 "parser.y"
+#line 530 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_arithmetic((yyvsp[(1) - (3)].sel), (yyvsp[(3) - (3)].sel), '*', scanner);
case 50:
/* Line 1806 of yacc.c */
-#line 530 "parser.y"
+#line 536 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_arithmetic((yyvsp[(1) - (3)].sel), (yyvsp[(3) - (3)].sel), '/', scanner);
case 51:
/* Line 1806 of yacc.c */
-#line 536 "parser.y"
+#line 542 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_arithmetic((yyvsp[(2) - (2)].sel), NULL, '-', scanner);
case 52:
/* Line 1806 of yacc.c */
-#line 542 "parser.y"
+#line 548 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_arithmetic((yyvsp[(1) - (3)].sel), (yyvsp[(3) - (3)].sel), '^', scanner);
case 53:
/* Line 1806 of yacc.c */
-#line 547 "parser.y"
+#line 553 "parser.y"
{ (yyval.sel) = (yyvsp[(2) - (3)].sel); }
break;
case 54:
/* Line 1806 of yacc.c */
-#line 555 "parser.y"
+#line 561 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_selelem_create(SEL_CONST);
case 55:
/* Line 1806 of yacc.c */
-#line 564 "parser.y"
+#line 570 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_keyword((yyvsp[(2) - (2)].meth), NULL, (yyvsp[(1) - (2)].str), scanner);
case 56:
/* Line 1806 of yacc.c */
-#line 578 "parser.y"
+#line 584 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_const_position((yyvsp[(2) - (7)].r), (yyvsp[(4) - (7)].r), (yyvsp[(6) - (7)].r));
case 57:
/* Line 1806 of yacc.c */
-#line 586 "parser.y"
+#line 592 "parser.y"
{ (yyval.sel) = (yyvsp[(2) - (3)].sel); }
break;
case 58:
/* Line 1806 of yacc.c */
-#line 591 "parser.y"
+#line 597 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_method((yyvsp[(1) - (2)].meth), (yyvsp[(2) - (2)].param), NULL, scanner);
case 59:
/* Line 1806 of yacc.c */
-#line 601 "parser.y"
+#line 607 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_position((yyvsp[(3) - (3)].sel), (yyvsp[(1) - (3)].str), scanner);
case 60:
/* Line 1806 of yacc.c */
-#line 614 "parser.y"
+#line 620 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_variable_ref((yyvsp[(1) - (1)].sel));
case 61:
/* Line 1806 of yacc.c */
-#line 622 "parser.y"
+#line 628 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_variable_ref((yyvsp[(1) - (1)].sel));
case 62:
/* Line 1806 of yacc.c */
-#line 630 "parser.y"
+#line 636 "parser.y"
{
BEGIN_ACTION;
(yyval.sel) = _gmx_sel_init_variable_ref((yyvsp[(1) - (1)].sel));
case 63:
/* Line 1806 of yacc.c */
-#line 643 "parser.y"
+#line 649 "parser.y"
{ (yyval.param) = process_param_list((yyvsp[(1) - (1)].param)); }
break;
case 64:
/* Line 1806 of yacc.c */
-#line 645 "parser.y"
+#line 651 "parser.y"
{ (yyval.param) = process_param_list((yyvsp[(1) - (2)].param)); }
break;
case 65:
/* Line 1806 of yacc.c */
-#line 649 "parser.y"
+#line 655 "parser.y"
{ (yyval.param) = NULL; }
break;
case 66:
/* Line 1806 of yacc.c */
-#line 651 "parser.y"
+#line 657 "parser.y"
{ (yyvsp[(2) - (2)].param)->next = (yyvsp[(1) - (2)].param); (yyval.param) = (yyvsp[(2) - (2)].param); }
break;
case 67:
/* Line 1806 of yacc.c */
-#line 656 "parser.y"
+#line 662 "parser.y"
{
BEGIN_ACTION;
(yyval.param) = _gmx_selexpr_create_param((yyvsp[(1) - (2)].str));
case 68:
/* Line 1806 of yacc.c */
-#line 664 "parser.y"
+#line 670 "parser.y"
{ (yyval.val) = NULL; }
break;
case 69:
/* Line 1806 of yacc.c */
-#line 665 "parser.y"
+#line 671 "parser.y"
{ (yyval.val) = (yyvsp[(1) - (1)].val); }
break;
case 70:
/* Line 1806 of yacc.c */
-#line 666 "parser.y"
+#line 672 "parser.y"
{ (yyval.val) = (yyvsp[(2) - (3)].val); }
break;
case 71:
/* Line 1806 of yacc.c */
-#line 670 "parser.y"
+#line 676 "parser.y"
{ (yyval.val) = (yyvsp[(1) - (1)].val); }
break;
case 72:
/* Line 1806 of yacc.c */
-#line 672 "parser.y"
+#line 678 "parser.y"
{ (yyvsp[(2) - (2)].val)->next = (yyvsp[(1) - (2)].val); (yyval.val) = (yyvsp[(2) - (2)].val); }
break;
case 73:
/* Line 1806 of yacc.c */
-#line 674 "parser.y"
+#line 680 "parser.y"
{ (yyvsp[(3) - (3)].val)->next = (yyvsp[(1) - (3)].val); (yyval.val) = (yyvsp[(3) - (3)].val); }
break;
case 74:
/* Line 1806 of yacc.c */
-#line 678 "parser.y"
+#line 684 "parser.y"
{ (yyval.val) = (yyvsp[(1) - (1)].val); }
break;
case 75:
/* Line 1806 of yacc.c */
-#line 679 "parser.y"
+#line 685 "parser.y"
{ (yyval.val) = (yyvsp[(2) - (3)].val); }
break;
case 76:
/* Line 1806 of yacc.c */
-#line 683 "parser.y"
+#line 689 "parser.y"
{ (yyval.val) = (yyvsp[(1) - (1)].val); }
break;
case 77:
/* Line 1806 of yacc.c */
-#line 685 "parser.y"
+#line 691 "parser.y"
{ (yyvsp[(2) - (2)].val)->next = (yyvsp[(1) - (2)].val); (yyval.val) = (yyvsp[(2) - (2)].val); }
break;
case 78:
/* Line 1806 of yacc.c */
-#line 687 "parser.y"
+#line 693 "parser.y"
{ (yyvsp[(3) - (3)].val)->next = (yyvsp[(1) - (3)].val); (yyval.val) = (yyvsp[(3) - (3)].val); }
break;
case 79:
/* Line 1806 of yacc.c */
-#line 691 "parser.y"
+#line 697 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value_expr((yyvsp[(1) - (1)].sel));
case 80:
/* Line 1806 of yacc.c */
-#line 697 "parser.y"
+#line 703 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value_expr((yyvsp[(1) - (1)].sel));
case 81:
/* Line 1806 of yacc.c */
-#line 703 "parser.y"
+#line 709 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value_expr((yyvsp[(1) - (1)].sel));
case 82:
/* Line 1806 of yacc.c */
-#line 709 "parser.y"
+#line 715 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value_expr((yyvsp[(1) - (1)].sel));
case 83:
/* Line 1806 of yacc.c */
-#line 714 "parser.y"
+#line 720 "parser.y"
{ (yyval.val) = (yyvsp[(1) - (1)].val); }
break;
case 84:
/* Line 1806 of yacc.c */
-#line 719 "parser.y"
+#line 725 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value(INT_VALUE);
case 85:
/* Line 1806 of yacc.c */
-#line 726 "parser.y"
+#line 732 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value(REAL_VALUE);
case 86:
/* Line 1806 of yacc.c */
-#line 733 "parser.y"
+#line 739 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value(STR_VALUE);
case 87:
/* Line 1806 of yacc.c */
-#line 739 "parser.y"
+#line 745 "parser.y"
{ (yyval.val) = (yyvsp[(1) - (1)].val); }
break;
case 88:
/* Line 1806 of yacc.c */
-#line 744 "parser.y"
+#line 750 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value(INT_VALUE);
case 89:
/* Line 1806 of yacc.c */
-#line 751 "parser.y"
+#line 757 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value(REAL_VALUE);
case 90:
/* Line 1806 of yacc.c */
-#line 758 "parser.y"
+#line 764 "parser.y"
{
BEGIN_ACTION;
(yyval.val) = _gmx_selexpr_create_value(REAL_VALUE);
/* Line 1806 of yacc.c */
-#line 2823 "parser.cpp"
+#line 2829 "parser.cpp"
default: break;
}
/* User semantic actions sometimes alter yychar, and that requires
/* Line 2067 of yacc.c */
-#line 766 "parser.y"
+#line 772 "parser.y"
static t_selexpr_value *
{
/* Line 2068 of yacc.c */
-#line 95 "parser.y"
+#line 101 "parser.y"
int i;
real r;
#pragma warning(disable: 4065)
#endif
-/* These macros should be used at the beginning and end of each semantic action
+/*! \name Exception handling macros for actions
+ *
+ * These macros should be used at the beginning and end of each semantic action
* that may throw an exception. For robustness, it's best to wrap all actions
* that call functions declared outside parser.y should be wrapped.
* These macros take care to catch any exceptions, store the exception (or
* cleanly if necessary.
* The code calling the parser should use
* _gmx_sel_lexer_rethrow_exception_if_occurred() to rethrow any exceptions.
+ * \{
*/
+//! Starts an action that may throw exceptions.
#define BEGIN_ACTION \
try {
+//! Finishes an action that may throw exceptions.
#define END_ACTION \
} \
catch(const std::exception &ex) \
else \
YYABORT; \
}
+//!\}
%}
%union{
* \author Teemu Murtola <teemu.murtola@cbr.su.se>
* \ingroup module_selection
*/
-/*! \internal \file scanner_flex.h
+/*! \cond
+ * \internal \file scanner_flex.h
* \brief Generated (from scanner.l) header file by Flex.
*
* This file contains definitions of functions that are needed in
* scanner_internal.cpp.
*
* \ingroup module_selection
+ * \endcond
*/
#ifdef HAVE_CONFIG_H
#include <config.h>
root->registerSubTopic<SimpleHelpTopic<SyntaxHelpText> >();
return move(root);
}
-//! \cond
+//! \endcond
} // namespace gmx
* \param[in] checkItem Functor to check an individual item.
*
* This method creates a compound checker \c compound within which all
- * values of the sequence are checked. Calls checkItem(&compound, *i)
+ * values of the sequence are checked. Calls \c checkItem(&compound, *i)
* with that compound for each iterator \c i in the range [begin, end).
- * \p checkItem should use check*() methods in the passed checker to
- * check the each value.
+ * \p checkItem should use the various check methods in the passed
+ * checker to check each value.
*
* This method can be used to check a sequence made of compound types.
* Typically \p checkItem will create a compound within the passed
biod.pnpi.spb.ru / alexxy / gromacs.git / commitdiff