text:
/* matrix/gsl_matrix_double.h
 *
 * Copyright (C) 1996, 1997, 1998, 1999, 2000, 2007 Gerard Jungman, Brian Gough
 *
 * 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 3 of the License, or (at
 * your option) any later version.
 *
 * This program is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
 */

#ifndef __GSL_MATRIX_DOUBLE_H__
#define __GSL_MATRIX_DOUBLE_H__

#include <stdlib.h>
#include <gsl/gsl_types.h>
#include <gsl/gsl_errno.h>
#include <gsl/gsl_inline.h>
#include <gsl/gsl_check_range.h>
#include <gsl/gsl_vector_double.h>
#include <gsl/gsl_blas_types.h>

#undef __BEGIN_DECLS
#undef __END_DECLS
#ifdef __cplusplus
# define __BEGIN_DECLS extern "C" {
# define __END_DECLS }
#else
# define __BEGIN_DECLS /* empty */
# define __END_DECLS /* empty */
#endif

__BEGIN_DECLS

typedef struct
{
  size_t size1;
  size_t size2;
  size_t tda;
  double * data;
  gsl_block * block;
  int owner;
} gsl_matrix;

typedef struct
{
  gsl_matrix matrix;
} _gsl_matrix_view;

typedef _gsl_matrix_view gsl_matrix_view;

typedef struct
{
  gsl_matrix matrix;
} _gsl_matrix_const_view;

typedef const _gsl_matrix_const_view gsl_matrix_const_view;

/* Allocation */

gsl_matrix *
gsl_matrix_alloc (const size_t n1, const size_t n2);

gsl_matrix *
gsl_matrix_calloc (const size_t n1, const size_t n2);

gsl_matrix *
gsl_matrix_alloc_from_block (gsl_block * b,
                                   const size_t offset,
                                   const size_t n1,
                                   const size_t n2,
                                   const size_t d2);

gsl_matrix *
gsl_matrix_alloc_from_matrix (gsl_matrix * m,
                                    const size_t k1,
                                    const size_t k2,
                                    const size_t n1,
                                    const size_t n2);

gsl_vector *
gsl_vector_alloc_row_from_matrix (gsl_matrix * m,
                                        const size_t i);

gsl_vector *
gsl_vector_alloc_col_from_matrix (gsl_matrix * m,
                                        const size_t j);

void gsl_matrix_free (gsl_matrix * m);

/* Views */

_gsl_matrix_view
gsl_matrix_submatrix (gsl_matrix * m,
                            const size_t i, const size_t j,
                            const size_t n1, const size_t n2);

_gsl_vector_view
gsl_matrix_row (gsl_matrix * m, const size_t i);

_gsl_vector_view
gsl_matrix_column (gsl_matrix * m, const size_t j);

_gsl_vector_view
gsl_matrix_diagonal (gsl_matrix * m);

_gsl_vector_view
gsl_matrix_subdiagonal (gsl_matrix * m, const size_t k);

_gsl_vector_view
gsl_matrix_superdiagonal (gsl_matrix * m, const size_t k);

_gsl_vector_view
gsl_matrix_subrow (gsl_matrix * m, const size_t i,
                         const size_t offset, const size_t n);

_gsl_vector_view
gsl_matrix_subcolumn (gsl_matrix * m, const size_t j,
                            const size_t offset, const size_t n);

_gsl_matrix_view
gsl_matrix_view_array (double * base,
                             const size_t n1,
                             const size_t n2);

_gsl_matrix_view
gsl_matrix_view_array_with_tda (double * base,
                                      const size_t n1,
                                      const size_t n2,
                                      const size_t tda);

_gsl_matrix_view
gsl_matrix_view_vector (gsl_vector * v,
                              const size_t n1,
                              const size_t n2);

_gsl_matrix_view
gsl_matrix_view_vector_with_tda (gsl_vector * v,
                                       const size_t n1,
                                       const size_t n2,
                                       const size_t tda);

_gsl_matrix_const_view
gsl_matrix_const_submatrix (const gsl_matrix * m,
                                  const size_t i, const size_t j,
                                  const size_t n1, const size_t n2);

_gsl_vector_const_view
gsl_matrix_const_row (const gsl_matrix * m,
                            const size_t i);

_gsl_vector_const_view
gsl_matrix_const_column (const gsl_matrix * m,
                               const size_t j);

_gsl_vector_const_view
gsl_matrix_const_diagonal (const gsl_matrix * m);

_gsl_vector_const_view
gsl_matrix_const_subdiagonal (const gsl_matrix * m,
                                    const size_t k);

_gsl_vector_const_view
gsl_matrix_const_superdiagonal (const gsl_matrix * m,
                                      const size_t k);

_gsl_vector_const_view
gsl_matrix_const_subrow (const gsl_matrix * m, const size_t i,
                               const size_t offset, const size_t n);

_gsl_vector_const_view
gsl_matrix_const_subcolumn (const gsl_matrix * m, const size_t j,
                                  const size_t offset, const size_t n);

_gsl_matrix_const_view
gsl_matrix_const_view_array (const double * base,
                                   const size_t n1,
                                   const size_t n2);

_gsl_matrix_const_view
gsl_matrix_const_view_array_with_tda (const double * base,
                                            const size_t n1,
                                            const size_t n2,
                                            const size_t tda);

_gsl_matrix_const_view
gsl_matrix_const_view_vector (const gsl_vector * v,
                                    const size_t n1,
                                    const size_t n2);

_gsl_matrix_const_view
gsl_matrix_const_view_vector_with_tda (const gsl_vector * v,
                                             const size_t n1,
                                             const size_t n2,
                                             const size_t tda);

/* Operations */

void gsl_matrix_set_zero (gsl_matrix * m);
void gsl_matrix_set_identity (gsl_matrix * m);
void gsl_matrix_set_all (gsl_matrix * m, double x);

int gsl_matrix_fread (FILE * stream, gsl_matrix * m) ;
int gsl_matrix_fwrite (FILE * stream, const gsl_matrix * m) ;
int gsl_matrix_fscanf (FILE * stream, gsl_matrix * m);
int gsl_matrix_fprintf (FILE * stream, const gsl_matrix * m, const char * format);

int gsl_matrix_memcpy(gsl_matrix * dest, const gsl_matrix * src);
int gsl_matrix_swap(gsl_matrix * m1, gsl_matrix * m2);
int gsl_matrix_tricpy(CBLAS_UPLO_t Uplo, CBLAS_DIAG_t Diag, gsl_matrix * dest, const gsl_matrix * src);

int gsl_matrix_swap_rows(gsl_matrix * m, const size_t i, const size_t j);
int gsl_matrix_swap_columns(gsl_matrix * m, const size_t i, const size_t j);
int gsl_matrix_swap_rowcol(gsl_matrix * m, const size_t i, const size_t j);
int gsl_matrix_transpose (gsl_matrix * m);
int gsl_matrix_transpose_memcpy (gsl_matrix * dest, const gsl_matrix * src);
int gsl_matrix_transpose_tricpy (CBLAS_UPLO_t Uplo_src, CBLAS_DIAG_t Diag, gsl_matrix * dest, const gsl_matrix * src);

double gsl_matrix_max (const gsl_matrix * m);
double gsl_matrix_min (const gsl_matrix * m);
void gsl_matrix_minmax (const gsl_matrix * m, double * min_out, double * max_out);

void gsl_matrix_max_index (const gsl_matrix * m, size_t * imax, size_t *jmax);
void gsl_matrix_min_index (const gsl_matrix * m, size_t * imin, size_t *jmin);
void gsl_matrix_minmax_index (const gsl_matrix * m, size_t * imin, size_t * jmin, size_t * imax, size_t * jmax);

int gsl_matrix_equal (const gsl_matrix * a, const gsl_matrix * b);

int gsl_matrix_isnull (const gsl_matrix * m);
int gsl_matrix_ispos (const gsl_matrix * m);
int gsl_matrix_isneg (const gsl_matrix * m);
int gsl_matrix_isnonneg (const gsl_matrix * m);

double gsl_matrix_norm1 (const gsl_matrix * m);

int gsl_matrix_add (gsl_matrix * a, const gsl_matrix * b);
int gsl_matrix_sub (gsl_matrix * a, const gsl_matrix * b);
int gsl_matrix_mul_elements (gsl_matrix * a, const gsl_matrix * b);
int gsl_matrix_div_elements (gsl_matrix * a, const gsl_matrix * b);
int gsl_matrix_scale (gsl_matrix * a, const double x);
int gsl_matrix_scale_rows (gsl_matrix * a, const gsl_vector * x);
int gsl_matrix_scale_columns (gsl_matrix * a, const gsl_vector * x);
int gsl_matrix_add_constant (gsl_matrix * a, const double x);
int gsl_matrix_add_diagonal (gsl_matrix * a, const double x);

/***********************************************************************/
/* The functions below are obsolete                                    */
/***********************************************************************/
int gsl_matrix_get_row(gsl_vector * v, const gsl_matrix * m, const size_t i);
int gsl_matrix_get_col(gsl_vector * v, const gsl_matrix * m, const size_t j);
int gsl_matrix_set_row(gsl_matrix * m, const size_t i, const gsl_vector * v);
int gsl_matrix_set_col(gsl_matrix * m, const size_t j, const gsl_vector * v);
/***********************************************************************/

/* inline functions if you are using GCC */

INLINE_DECL double   gsl_matrix_get(const gsl_matrix * m, const size_t i, const size_t j);
INLINE_DECL void    gsl_matrix_set(gsl_matrix * m, const size_t i, const size_t j, const double x);
INLINE_DECL double * gsl_matrix_ptr(gsl_matrix * m, const size_t i, const size_t j);
INLINE_DECL const double * gsl_matrix_const_ptr(const gsl_matrix * m, const size_t i, const size_t j);

#ifdef HAVE_INLINE
INLINE_FUN
double
gsl_matrix_get(const gsl_matrix * m, const size_t i, const size_t j)
{
#if GSL_RANGE_CHECK
  if (GSL_RANGE_COND(1))
    {
      if (i >= m->size1)
        {
          GSL_ERROR_VAL("first index out of range", GSL_EINVAL, 0) ;
        }
      else if (j >= m->size2)
        {
          GSL_ERROR_VAL("second index out of range", GSL_EINVAL, 0) ;
        }
    }
#endif
  return m->data[i * m->tda + j] ;
}

INLINE_FUN
void
gsl_matrix_set(gsl_matrix * m, const size_t i, const size_t j, const double x)
{
#if GSL_RANGE_CHECK
  if (GSL_RANGE_COND(1))
    {
      if (i >= m->size1)
        {
          GSL_ERROR_VOID("first index out of range", GSL_EINVAL) ;
        }
      else if (j >= m->size2)
        {
          GSL_ERROR_VOID("second index out of range", GSL_EINVAL) ;
        }
    }
#endif
  m->data[i * m->tda + j] = x ;
}

INLINE_FUN
double *
gsl_matrix_ptr(gsl_matrix * m, const size_t i, const size_t j)
{
#if GSL_RANGE_CHECK
  if (GSL_RANGE_COND(1))
    {
      if (i >= m->size1)
        {
          GSL_ERROR_NULL("first index out of range", GSL_EINVAL) ;
        }
      else if (j >= m->size2)
        {
          GSL_ERROR_NULL("second index out of range", GSL_EINVAL) ;
        }
    }
#endif
  return (double *) (m->data + (i * m->tda + j)) ;
}

INLINE_FUN
const double *
gsl_matrix_const_ptr(const gsl_matrix * m, const size_t i, const size_t j)
{
#if GSL_RANGE_CHECK
  if (GSL_RANGE_COND(1))
    {
      if (i >= m->size1)
        {
          GSL_ERROR_NULL("first index out of range", GSL_EINVAL) ;
        }
      else if (j >= m->size2)
        {
          GSL_ERROR_NULL("second index out of range", GSL_EINVAL) ;
        }
    }
#endif
  return (const double *) (m->data + (i * m->tda + j)) ;
}

#endif

__END_DECLS

#endif /* __GSL_MATRIX_DOUBLE_H__ */


text:
//
// Created by peixinho on 19/12/16.
//

#include <ml_matrix.h>
#include <cblas.h>
#include <ml_memory.h>

MlMatrix* ml_create_matrix(size_t nrow, size_t ncol) {
    MlMatrix* m = NULL;
    m = ml_alloc(MlMatrix);
    m->nrows = nrow;
    m->ncols = ncol;
    m->n = nrow*ncol;
    m->data = ml_alloc_array(m->n, float);
    return m;
}

void mlMultMatrixIp(MlMatrix* A, MlMatrix* B, MlMatrix* M) {
    double alpha = 1.0, beta = 0.0;

    cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, A->nrows, B->ncols, \
        A->ncols, alpha, A->data, A->ncols, B->data, B->ncols, beta, \
        M->data, B->ncols);
}

MlMatrix* mlMultMatrix(MlMatrix* A, MlMatrix* B) {
    MlMatrix* M = ml_create_matrix(B->ncols, A->nrows);
    mlMultMatrixIp(A, B, M);
    return M;
}

float mlMeanMatrix(MlMatrix* m) {
    float mean = 0.0f;
    for (int i = 0; i < m->n; ++i) {
        mean += m->data[i] / m->n;
    }
    return mean;
}

MlMatrix* mlMatrixSumCol(MlMatrix* m) {
    MlMatrix* sum = ml_create_matrix(1, m->ncols);

    for (int i = 0; i < m->nrows; ++i) {
        for (int j = 0; j < m->ncols; ++j) {
            ml_mat_elem(sum, 0, j) += ml_mat_elem(m, i, j);
        }
    }

    return sum;
}

void ml_set_matrix(MlMatrix *m, float val) {
    for (int i = 0; i < m->n; ++i)
        m->data[i] = val;
}

void ml_free_matrix(MlMatrix *m) {
    ml_free(m->data);
    ml_free(m);
}

void ml_print_matrix(MlMatrix *m) {
    printf("[");
    for (int r = 0; r < m->nrows; ++r) {
        printf("[");
        for (int c = 0; c < m->ncols; ++c) {
            printf("%4.2f", ml_mat_elem(m, r, c));
            if(c < m->ncols-1)
                printf(", ");
        }
        printf("]");
        if(r<m->nrows-1)
            printf("\n");
    }
    printf("]\n");
}

ext: c

max_stars_repo_name: alanpeixinho/lib-machine-learning

max_issues_repo_name: alanpeixinho/lib-machine-learning

max_forks_repo_name: alanpeixinho/lib-machine-learning


text:
//
// Created by david on 2018-11-30.
//

#ifndef DMRG_CLASS_XDMRG_FULL_FUNCTOR_H
#define DMRG_CLASS_XDMRG_FULL_FUNCTOR_H
#ifdef OpenMP_AVAILABLE
#include <omp.h>
#endif

#ifdef OpenBLAS_AVAILABLE
#include <cblas.h>
#endif

#include <Eigen/Core>
#include <unsupported/Eigen/CXX11/Tensor>
#include <general/nmspc_tensor_extra.h>
#include <general/class_tic_toc.h>
template<typename Scalar>
class class_xDMRG_full_functor {
private:
    double variance;
    double energy  ;
    double energy_lower_bound;
    double energy_upper_bound;
    double energy_target;
    double energy_window;
public:
    template <typename T>
    int sgn(const T val) const {
        return (T(0) < val) - (val < T(0));
    }

    using MatrixType_ = Eigen::Matrix<Scalar,Eigen::Dynamic, Eigen::Dynamic>;
    using VectorType_ = Eigen::Matrix<Scalar,Eigen::Dynamic, 1>;
    size_t   counter = 0;
//    const size_t shape;
    void set_energy_bounds(double E_lower, double E_upper);
    bool have_bounds_on_energy = false;
    double get_variance(){return variance;}
    double get_energy  (){return energy  ;}
    size_t get_count   (){return counter;}
    Eigen::Tensor<double,4> HA_MPO;
    Eigen::Tensor<double,4> HB_MPO;
    Eigen::Tensor<double,3> Lblock;
    Eigen::Tensor<double,3> Rblock;
    Eigen::Tensor<double,4> Lblock2;
    Eigen::Tensor<double,4> Rblock2;

    Eigen::Tensor<double,6> HAHB;
    Eigen::Tensor<double,8> HAHB2;
    Eigen::DSizes<long,4>   dsizes;
    class_tic_toc t_lbfgs;
    class_xDMRG_full_functor(
            const Eigen::Tensor<Scalar,4> &HA_MPO_,
            const Eigen::Tensor<Scalar,4> &HB_MPO_,
            const Eigen::Tensor<Scalar,3> &Lblock_,
            const Eigen::Tensor<Scalar,3> &Rblock_,
            const Eigen::Tensor<Scalar,4> &Lblock2_,
            const Eigen::Tensor<Scalar,4> &Rblock2_,
            const Eigen::DSizes<long,4>   &dsizes_
            );

    double get_vH2v(const Eigen::Matrix<double,Eigen::Dynamic,1> &v);
    double get_vHv(const Eigen::Matrix<double,Eigen::Dynamic,1> &v);
    Eigen::VectorXd get_vH2 (const Eigen::Matrix<double,Eigen::Dynamic,1> &v);
    Eigen::VectorXd get_vH (const Eigen::Matrix<double,Eigen::Dynamic,1> &v);
    double operator()(const Eigen::Matrix<double,Eigen::Dynamic,1> &v, Eigen::Matrix<double,Eigen::Dynamic,1> &grad);

};

#endif //DMRG_CLASS_XDMRG_FULL_FUNCTOR_H

ext: h

max_stars_repo_name: DavidAce/DMRG

max_issues_repo_name: DavidAce/DMRG

max_forks_repo_name: DavidAce/DMRG


text:
#pragma once

#include <iostream>
#include <map>
#include <memory>
#include <string>
#include <vector>

#include <gsl/gsl_complex.h>
#include <gsl/gsl_math.h>

#include "../Scanner/scanner.h"
#include "../Utils/Types.h"

class Expression;

typedef std::shared_ptr<Expression> expression;
typedef std::map<std::string, expression> Variables;

extern const Variables emptyVars;

typedef expression(*TransformerFunction)(expression);

class Expression: public std::enable_shared_from_this<Expression> {

    public:
        const Scanner::Type kind;

        Expression(Scanner::Type kind);
        virtual ~Expression();

        virtual expression simplify();
        virtual expression derivative(const std::string& var);
        virtual expression integrate(const std::string& var);

        virtual bool isComplex() const = 0;
        virtual bool isEvaluable(const Variables& vars = emptyVars) const = 0;
        virtual bool isNumber() const { return false; }

        virtual expression eval(const Variables& vars = emptyVars) = 0;
        virtual double value(const Variables& vars = emptyVars) const = 0;
        virtual gsl_complex complex(const Variables& vars = emptyVars) const;

        virtual bool equals(expression, double precision=1e-15) const = 0;

        std::vector<double> array();
        virtual expression at(const int index);
        virtual double get(const int index);
        virtual size_t shape(const int axis) const;
        virtual size_t size() const;

        virtual expression call(expression e);
        inline expression operator()(expression e){ return call(e); }

        virtual expression apply(TransformerFunction f);
        inline expression operator()(TransformerFunction f){ return apply(f); }

        inline bool is(Scanner::Type kind) const { return this->kind == kind; }
        virtual std::string repr() const;
        virtual int id() const;

        expression copy();

        virtual std::ostream& print(std::ostream&, const bool pretty=false) const = 0;
        virtual std::ostream& postfix(std::ostream&) const = 0;

    private:
        struct gsl_expression_struct {
            Expression* e;
            const std::string var;
            Variables vars;

            gsl_expression_struct(Expression*, const std::string var = "x");
        };
        static double gsl_expression_function(double x, void* p);
        std::unique_ptr<gsl_expression_struct> gsl_params;

    public:
        virtual gsl_function function(const std::string& var = "x");

    private:
        struct ExpressionIterator {
            expression e;
            size_t index;

            ExpressionIterator(expression e, size_t index);

            expression operator*();

            bool operator==(const ExpressionIterator& other) const;
            bool operator!=(const ExpressionIterator& other) const;
            bool operator<(const ExpressionIterator& other) const;
            bool operator<=(const ExpressionIterator& other) const;
            bool operator>(const ExpressionIterator& other) const;
            bool operator>=(const ExpressionIterator& other) const;

            ExpressionIterator& operator++();
            ExpressionIterator& operator--();
        };
    public:
        ExpressionIterator begin();
        ExpressionIterator end();
};

std::ostream& operator<<(std::ostream&, const expression);

#define EXPRESSION_PARTIAL_OVERRIDES                                                \
    bool isComplex() const override;                                                \
    bool isEvaluable(const Variables& vars = emptyVars) const override;             \
    bool equals(expression, double precision) const override;                       \
    std::ostream& print(std::ostream&, const bool pretty = false) const override;   \
    std::ostream& postfix(std::ostream&) const override;

#define EXPRESSION_OVERRIDES                                                        \
    EXPRESSION_PARTIAL_OVERRIDES                                                    \
    expression eval(const Variables& vars = emptyVars) override;                    \
    double value(const Variables& vars = emptyVars) const override;

max_stars_repo_name: antoniojkim/CalcPlusPlus

max_issues_repo_name: antoniojkim/CalcPlusPlus

max_forks_repo_name: antoniojkim/CalcPlusPlus


text:
#pragma once
#include <Hypo/System/Exports.h>
#include <Hypo/Config.h>
#include <ostream>
#include <gsl/span>
#include <Hypo/System/Exports.h>

#include <Hypo/System/Buffer/Buffer.h>
#include <Hypo/System/Streams/MemoryStream.h>

namespace Hypo
{
	class HYPO_SYSTEM_API BinaryWriter
	{
	public:
		BinaryWriter(std::ostream& stream);

		BinaryWriter& operator <<(uInt8 value);
		BinaryWriter& operator <<(uInt16 value);
		BinaryWriter& operator <<(uInt32 value);
		BinaryWriter& operator <<(uInt64 value);

		BinaryWriter& operator <<(Int8 value);
		BinaryWriter& operator <<(Int16 value);
		BinaryWriter& operator <<(Int32 value);
		BinaryWriter& operator <<(Int64 value);

		BinaryWriter& operator <<(const std::string& value);
		BinaryWriter& operator <<(const char* value);

		template<typename T>
		BinaryWriter& operator <<(const std::vector<T>& value)
		{
			const uInt32 size(static_cast<uInt32>(value.size()));
			*this << size;
			for (auto& element : value)
			{
				*this << element;
			}
			return *this;
		}

		void WriteRaw(const std::string& data);
		void WriteRaw(const gsl::span<uInt8>& data);
		void WriteRaw(const void* data, int size);

		void WriteBOM();

		void Write7BitEncoded(uInt32 value);
		void Write7BitEncoded(uInt64 value);

		void Flush()const { m_Stream.flush(); }
		bool Good()const { return m_Stream.good(); }
		bool Bad()const { return m_Stream.bad(); }
		bool Fail() const { return m_Stream.fail(); }

		std::ostream& GetStream() const {return m_Stream;}
	private:
		std::ostream& m_Stream;
	};

	template <typename T>
	class BasicMemoryBinaryWriter : public BinaryWriter
		/// A convenient wrapper for using Buffer and MemoryStream with BinarWriter.
	{
	public:
		BasicMemoryBinaryWriter(Buffer<T>& data) :
			BinaryWriter(_ostr),
			_data(data),
			_ostr(data.begin(), data.Capacity())
		{
		}

		~BasicMemoryBinaryWriter()
		{
			try
			{
				Flush();
			}
			catch (...)
			{
			}
		}

		Buffer<T>& data()
		{
			return _data;
		}

		const Buffer<T>& data() const
		{
			return _data;
		}

		const MemoryOutputStream& stream() const
		{
			return _ostr;
		}

		MemoryOutputStream& stream()
		{
			return _ostr;
		}

	private:
		Buffer<T>& _data;
		MemoryOutputStream _ostr;
	};

	typedef BasicMemoryBinaryWriter<char> MemoryBinaryWriter;

}

max_stars_repo_name: TheodorLindberg/Hypo

max_issues_repo_name: TheodorLindberg/Hypo

max_forks_repo_name: TheodorLindberg/Hypo


text:
/************************************************************************
 Tuple.h.h - Copyright marcello

 Here you can write a license for your code, some comments or any other
 information you want to have in your generated code. To to this simply
 configure the "headings" directory in uml to point to a directory
 where you have your heading files.

 or you can just replace the contents of this file with your own.
 If you want to do this, this file is located at

 /usr/share/apps/umbrello/headings/heading.h

 -->Code Generators searches for heading files based on the file extension
 i.e. it will look for a file name ending in ".h" to include in C++ header
 files, and for a file name ending in ".java" to include in all generated
 java code.
 If you name the file "heading.<extension>", Code Generator will always
 choose this file even if there are other files with the same extension in the
 directory. If you name the file something else, it must be the only one with that
 extension in the directory to guarantee that Code Generator will choose it.

 you can use variables in your heading files which are replaced at generation
 time. possible variables are : author, date, time, filename and filepath.
 just write %variable_name%

 This file was generated on Sat Nov 10 2007 at 15:05:38
 The original location of this file is /home/marcello/Projects/fitted/Developing/Tuple.h
 **************************************************************************/

#ifndef TUPLE_H
#define TUPLE_H

#include <fstream>
#include <vector>
//#include <gsl/gsl_vector.h>

//<< Sam
typedef enum
{
  L1,
  L2,
  Linfty
} NormType;

using namespace std;

namespace PoliFitted
{

/**
 * Tuples are sequence containers representing arrays with additional functionalities.
 *
 * Just like arrays, tuples use contiguous storage locations for their elements, which
 * means that their elements can also be accessed using offsets on regular pointers to its elements,
 * and just as efficiently as in arrays.
 *
 * Like arrays, internally, tuples do not store information about the dimension that must be handled
 * manually from outside the class.
 */
class Tuple
{
  public:

    /**
     * Converts N vectors to a single tuple obtained by
     * concatenation of the vectors (double) in the given order
     * Vectors are received as sequence of const_iterator:
     * v1.begin, v1.end, v2.begin, v2.end, v3.begin ...
     *
     * @param count    final number of elements of the tuple.
     Sum of size of each vector
     * @param ...      sequence of const_iterator (begin,end,begin,end,...)
     * @return the reference of the new tuple
     */
    //static Tuple* CONVERT(unsigned int count, ...);

    // Constructors/Destructors
    //

    /**
     * Empty Constructor
     */
    Tuple();

    /**
     * Construct a tuple of specified size
     * @param size the size of the tuple
     */
    Tuple(unsigned int size);

    /**
     * Empty Destructor
     */
    ~Tuple();

    /**
     * Returns a reference to the element at position \f$pos\f$ in the tuple.
     *
     * Programs should never call this function with an argument \f$pos\f$ that is
     * out of range, since this causes undefined behavior.
     * @brief Access element
     * @param pos Position of an element in the container. Notice that the first element has a position of 0 (not 1).
     * @return The element at the specified position in the tuple.
     */
    float& operator[](const unsigned int pos);

    /**
     * Returns a reference to the element at position \f$pos\f$ in the tuple.
     *
     * The function does not automatically check whether \f$pos\f$ is within the
     * bounds of valid elements in the vector.
     * Programs should never call this function with an argument \f$pos\f$ that is
     * out of range, since this causes undefined behavior.
     *
     * @brief Access element
     * @param pos Position of an element in the container. Notice that the first element has a position of 0 (not 1).
     * @return The element at the specified position in the container.
     */
    float& at(const unsigned int pos);

    /**
     * Compute a distance between the current and a given tuple. In other words,
     * it computes the p-norm between tuples.
     * \a p must be one of: "L1", "L2", "Linfty". "L1" is the absolute distance,
     * "L2" is the Euclidean distance, and "Linfty" is the maximum norm.
     * The norm is computed taking into account the first \a size elements of
     * the current and given tuples.
     *
     * Programs should never call this function with an argument \f$size\f$ that is
     * out of range for one of the two tuples, since this causes undefined behavior.
     *
     * @brief Distance computation
     * @param t An other tuple.
     * @param size the number of elements to be considered.
     * @param type the norm type
     * @return the \a p-norm
     */
    float Distance(Tuple* t, unsigned int size, NormType type = L2);

    /**
     * Assigns new contents to the tuple, replacing its current contents, without modifying its size.
     *
     * Since the dimension of the tuple is not changed accordingly to the given array,
     * the behavior is undefined when \a size is greater than the current tuple dimension.
     * In this case it is necessary to call the complementary function Tuple::Resize.
     *
     * The tuple is filled starting from the first elements.
     *
     * @brief Assign elements
     * @param values Array storing the new elements.
     * @param size The size of the given array
     */
    void SetValues(float* values, unsigned int size);

    /**
     * Returns a direct pointer to the memory array used internally by the tuple to store its owned elements.
     *
     * Because elements in the tuple are guaranteed to be stored in contiguous storage locations in the same
     * order as represented by the tuple, the pointer retrieved can be offset to access any element in the array.
     * @brief Access data
     * @return A pointer to the first element in the array used internally by the tuple.
     */
    float* GetValues();

    /**
     * Returns a clone of the tuple by copying the first \a input_size elements.
     *
     * @brief Clone object
     * @param input_size the number of elements to be copied
     * @return A pointer to the new tuple.
     */
    Tuple* Clone(unsigned int input_size);

    /**
     * Resizes the container so that it contains \a new_size elements.
     *
     * If \a new_size is smaller than the current container size, the content is reduced to its first n elements,
     * removing those beyond (and destroying them).
     * If \a new_size is greater than the current container size, the content is expanded by inserting at the end
     * as many elements as needed to reach a size of \a new_size.
     * The value of the newly allocated portion is indeterminate.
     *
     * Notice that this function changes the actual content of the container by inserting or erasing elements
     * from it.
     *
     * @brief Change size
     * @param new_size New container size, expressed in number of elements.
     */
    void Resize(unsigned int new_size);

    template<class T>
    T* GetValues(unsigned int input_size);

    //gsl_vector* GetGSLVector(unsigned int input_size);

    std::vector<double> ToVector(unsigned int input_size);

  protected:

  private:

    float* mValues;

};

inline float& Tuple::at(const unsigned int pos)
{
  return mValues[pos];
}

inline float& Tuple::operator[](const unsigned int pos)
{
  return mValues[pos];
}

inline void Tuple::SetValues(float* values, unsigned int size)
{
  for (unsigned int i = 0; i < size; i++)
  {
    mValues[i] = values[i];
  }
}

inline float* Tuple::GetValues()
{
  return mValues;
}

template<class T>
inline T* Tuple::GetValues(unsigned int input_size)
{
  T* values = new T[input_size];
  for (unsigned int i = 0; i < input_size; i++)
  {
    values[i] = (T) mValues[i];
  }
  return values;
}

/*inline gsl_vector* Tuple::GetGSLVector(unsigned int input_size)
 {
 gsl_vector* values = gsl_vector_alloc(input_size);
 for (unsigned int i = 0; i < input_size; i++)
 {
 gsl_vector_set(values, i, mValues[i]);
 }
 return values;
 }*/

inline std::vector<double> Tuple::ToVector(unsigned int input_size)
{
  std::vector<double> values;
  for (unsigned int i = 0; i < input_size; i++)
  {
    values.push_back(mValues[i]);
  }
  return values;
}

}

#endif // TUPLE_H

max_stars_repo_name: akangasr/sdirl

max_issues_repo_name: akangasr/sdirl

max_forks_repo_name: akangasr/sdirl


text:
#pragma once
#include "nkg_point.h"
#include "nkg_rect.h"
#include "nk_essential.h"
#include "utility/types.h"
#include "utility/c++std/string_view.h"
#include <nuklear.h>
#include <gsl/gsl>

namespace cws80 {

// drawing
void im_draw_point(nk_command_buffer *out, f32 x, f32 y, nk_color color);
void im_fill_polygon(nk_command_buffer *out, gsl::span<const im_pointf> points, nk_color color);
void im_draw_image(nk_command_buffer *out, im_rectf rect, const im_texture &tex,
                   nk_color color);
void im_draw_tiled_image(nk_command_buffer *out, im_rectf rect,
                         const im_texture &tex, nk_color color);
void im_draw_title_box(nk_context *ctx, im_rectf bounds, cxx::string_view title,
                       im_rectf *content_bounds = nullptr);
void im_draw_title_box_alt(nk_context *ctx, im_rectf bounds, cxx::string_view title);
void im_draw_led_screen(nk_context *ctx, im_rectf bounds, cxx::string_view text,
                        im_pointf textoff);

// drawing in screen coordinates
void ims_stroke_line(nk_context *ctx, f32 x0, f32 y0, f32 x1, f32 y1,
                     f32 line_thickness, nk_color color);
void ims_stroke_rect(nk_context *ctx, im_rectf rect, f32 rounding,
                     f32 line_thickness, nk_color color);
void ims_fill_rect(nk_context *ctx, im_rectf rect, f32 rounding, nk_color color);
void ims_fill_polygon(nk_context *ctx, gsl::span<const im_pointf> points, nk_color color);
void ims_draw_image(nk_context *ctx, im_rectf rect, const im_texture &tex, nk_color color);
void ims_draw_tiled_image(nk_context *ctx, im_rectf rect, const im_texture &tex,
                          nk_color color);
void ims_draw_text(nk_context *ctx, im_rectf rect, const char *text, int len,
                   const nk_user_font *fnt, nk_color bg, nk_color fg);
void ims_push_scissor(nk_context *ctx, im_rectf rect);
//
void ims_draw_title_box(nk_context *ctx, im_rectf bounds, cxx::string_view title,
                        im_rectf *content_bounds = nullptr);
void ims_draw_title_box_alt(nk_context *ctx, im_rectf bounds, cxx::string_view title);
void ims_draw_led_screen(nk_context *ctx, im_rectf bounds,
                         cxx::string_view text, im_pointf textoff);

// text
f32 im_text_width(const nk_user_font *fnt, cxx::string_view text);

}  // namespace cws80


text:
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include <complex.h>
#include <gsl/gsl_linalg.h>
#include <gsl/gsl_matrix.h>

#define _USE_MATH_DEFINES

// global variables
int n = 10;

double incident(double x);
double green_function(double r1[n][3], double r2[n][3]);

int main()
{
	int i, j, k;
	double x, y, z;
	double scatterers[n][3];

	gsl_matrix * mmat = gsl_matrix_alloc (n,n);

	srand((unsigned)time(NULL));
	for (i = 1; i < n; i++)
	{
		x = 20*(double)rand()/RAND_MAX;
		y = 50*(double)rand()/RAND_MAX;
		z = 30*(double)rand()/RAND_MAX;

		scatterers[i][0] = x;
		scatterers[i][1] = y;
		scatterers[i][2] = z;

		printf("Scatterer %d position: %lf %lf %lf \n", i, scatterers[i][0], scatterers[i][1], scatterers[i][2]);

	//double field_scatterer = incident(scatterers[i][0]);

	//printf("Incident field at scatterer %d: %lf\n", i, field_scatterer);

	// fixed below
		double field_scatterer = incident(scatterers[i][0]);

	// solving the matrix equation

		gsl_matrix_set (mmat, i,i, field_scatterer);
		// continue rewriting this part

	//	gsl_vector_view b
	//		= gsl_vector_view_array (field_scatterer, n);
	// gsl vector set instead

	//	gsl_vector *total_field = gsl_vector_alloc (n);
	}
	double green = green_function(scatterers, scatterers);
	int s;
	//gsl_permutation * p = gsl_permutation_alloc (n);
	//gsl_linalg_LU_decomp (&green.matrix, p, &s);
	//gsl_linalg_LU_solve (&green.matrix, p, &b.vector, total_field);

	//printf("Total field = \n");
	//gsl_vector_fprintf (stdout, total_field, "%g");

	//gsl_permutation_free (p);

	// defining surface of sphere
	double phi[150];
	double theta[150];
	int r = 500;

	for (i=(M_PI-2); i < (M_PI+2); i++)
	{
		phi[i] =(double)(i)/((M_PI-2) - (M_PI+2));

	}
	for (i=(M_PI/2 - 2); i < (M_PI/2 +2); i++)
	{
		theta[i] = (double)(i)/((M_PI/2 - 2) - (M_PI/2 +2));
	}

	return 0;
}

double incident(double x)
{
	double E_0 = 1.0;

	return E_0*cexp(1*I*2*M_PI*x);
}

double green_function(double r1[n][3], double r2[n][3])
{
	int i, j;
	gsl_matrix * m = gsl_matrix_alloc (n,n);

	double r[n], g[n];
	for (i = 0; i < n; i++)
	{
		for (j = 0; j < n; j++)
		{
			r[n] = sqrt((r1[i][1]-r2[j][1])*(r1[i][1]-r2[j][1]) + (r1[i][2]-r2[j][2])*(r1[i][2]-r2[j][2]) + (r1[i][3] - r2[j][3])*(r1[i][3] - r2[j][3]));

			if (r[n] == 0.0)
			{
				g[n] = 1.0;
			}
			else
			{
				g[n] = exp(1*I*2*M_PI*r[n])/(4*M_PI*r[n]);
			}

			gsl_matrix_set (m, i, j, g[n]);
		}
	}

	for (i = 0; i < 100; i++)
	{
		for (j = 0; j < n; j++)
		{
			printf ("m(%d,%d) = %g\n", i, j,
					gsl_matrix_get (m, i, j));
		}
	}
	gsl_matrix_free (m);
}


text:
/*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*
**                                                                                **
** This file forms part of the Underworld geophysics modelling application.         **
**                                                                                **
** For full license and copyright information, please refer to the LICENSE.md file  **
** located at the project root, or contact the authors.                             **
**                                                                                **
**~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*~*/
#include <stdlib.h>
#include <string.h>
#include <assert.h>
#include <mpi.h>
#include <petsc.h>
#include <petscvec.h>
#include <petscmat.h>
#include <petscksp.h>
#include <petscpc.h>
#include <petscsnes.h>

#include <petscversion.h>
#if ( (PETSC_VERSION_MAJOR >= 3) && (PETSC_VERSION_MINOR >=3) )
  #if (PETSC_VERSION_MINOR >=6)
     #include "petsc/private/kspimpl.h"
  #else
     #include "petsc-private/kspimpl.h"   /*I "petscksp.h" I*/
  #endif
#else
  #include "private/kspimpl.h"   /*I "petscksp.h" I*/
#endif

#include <StGermain/libStGermain/src/StGermain.h>
#include <StgDomain/libStgDomain/src/StgDomain.h>
#include <StgFEM/libStgFEM/src/StgFEM.h>
#include <PICellerator/libPICellerator/src/PICellerator.h>
#include <Underworld/libUnderworld/src/Underworld.h>

#include "Solvers/KSPSolvers/src/KSPSolvers.h"
#include "BSSCR/petsccompat.h"

#include "Test/TestKSP.h"
#include "BSSCR/BSSCR.h"

/************************************************************************/
/*** This function is called from _StokesBlockKSPInterface_Initialise **************/
/************************************************************************/
#undef __FUNCT__
#define __FUNCT__ "KSPRegisterAllKSP"
PetscErrorCode KSPRegisterAllKSP(const char path[])
{

    PetscFunctionBegin;
    KSPRegisterTEST(path);/* not sure if the path matters much here. everything still worked even though I had it wrong */
    KSPRegisterBSSCR (path);
    PetscFunctionReturn(0);
}


text:
// SPDX-License-Identifier: MIT
// The MIT License (MIT)
//
// Copyright (c) 2014-2018, Institute for Software & Systems Engineering
// Copyright (c) 2018-2019, Johannes Leupolz
//
// Permission is hereby granted, free of charge, to any person obtaining a copy
// of this software and associated documentation files (the "Software"), to deal
// in the Software without restriction, including without limitation the rights
// to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
// copies of the Software, and to permit persons to whom the Software is
// furnished to do so, subject to the following conditions:
//
// The above copyright notice and this permission notice shall be included in
// all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
// OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
// THE SOFTWARE.

#ifndef PEMC_LMC_LMC_H_
#define PEMC_LMC_LMC_H_

#include <atomic>
#include <functional>
#include <gsl/span>
#include <string>
#include <vector>

#include "pemc/basic/dll_defines.h"
#include "pemc/basic/label.h"
#include "pemc/basic/model_capacity.h"
#include "pemc/basic/probability.h"
#include "pemc/basic/tsc_index.h"
#include "pemc/formula/formula.h"

namespace pemc {

struct LmcStateEntry {
  TransitionIndex from;
  int32_t elements;
};

// Transition = Target + Probability
struct LmcTransitionEntry {
  Probability probability;
  Label label;
  StateIndex state;

  LmcTransitionEntry() = default;
  LmcTransitionEntry(Probability _probability, Label _label, StateIndex _state)
      : probability(_probability), label(_label), state(_state) {}
};

class Lmc {
 private:
  TransitionIndex maxNumberOfTransitions = 0;
  std::atomic<TransitionIndex> transitionCount{0};
  std::vector<LmcTransitionEntry> transitions;
  TransitionIndex initialTransitionFrom =
      -1;  // is uninitialized at first, but may be something else than 0
  int32_t initialTransitionElements = 0;

  StateIndex maxNumberOfStates = 0;
  StateIndex stateCount = 0;
  std::vector<LmcStateEntry> states;

  std::vector<std::string> labelIdentifier;

  TransitionIndex getPlaceForNewTransitionEntries(NoOfElements number);

 public:
  Lmc();

  gsl::span<LmcStateEntry> getStates();

  gsl::span<std::string> getLabelIdentifier();
  void setLabelIdentifier(const std::vector<std::string>& _labelIdentifier);
  std::function<bool(TransitionIndex)> createLabelBasedFormulaEvaluator(
      Formula* formula);

  gsl::span<LmcTransitionEntry> getTransitions();
  gsl::span<LmcTransitionEntry> getInitialTransitions();
  std::tuple<TransitionIndex, TransitionIndex> getInitialTransitionIndexes();
  gsl::span<LmcTransitionEntry> getTransitionsOfState(StateIndex state);
  std::tuple<TransitionIndex, TransitionIndex> getTransitionIndexesOfState(
      StateIndex state);

  TransitionIndex getPlaceForNewTransitionEntriesOfState(StateIndex stateIndex,
                                                         NoOfElements number);
  TransitionIndex getPlaceForNewInitialTransitionEntries(NoOfElements number);
  void setLmcTransitionEntry(TransitionIndex index,
                             const LmcTransitionEntry& entry);
  void createStutteringState(StateIndex stutteringStateIndex);

  void initialize(ModelCapacity& modelCapacity);
  void finishCreation(StateIndex _stateCount);
  void validate();
};

}  // namespace pemc

#endif  // PEMC_LMC_LMC_H_

max_stars_repo_name: joleuger/pemc

max_issues_repo_name: joleuger/pemc

max_forks_repo_name: joleuger/pemc


text:
/*
 * Copyright (c) 1997-1999 Massachusetts Institute of Technology
 *
 * 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.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 *
 */

/* This file was automatically generated --- DO NOT EDIT */
/* Generated on Mon Mar  8 17:46:56 EST 1999 */

#include <fftw-int.h>
#include <fftw.h>

/* Generated by: ./genfft -magic-alignment-check -magic-twiddle-load-all -magic-variables 4 -magic-loopi -twiddleinv 32 */

/*
 * This function contains 434 FP additions, 208 FP multiplications,
 * (or, 340 additions, 114 multiplications, 94 fused multiply/add),
 * 90 stack variables, and 128 memory accesses
 */
static const fftw_real K555570233 = FFTW_KONST(+0.555570233019602224742830813948532874374937191);
static const fftw_real K831469612 = FFTW_KONST(+0.831469612302545237078788377617905756738560812);
static const fftw_real K980785280 = FFTW_KONST(+0.980785280403230449126182236134239036973933731);
static const fftw_real K195090322 = FFTW_KONST(+0.195090322016128267848284868477022240927691618);
static const fftw_real K923879532 = FFTW_KONST(+0.923879532511286756128183189396788286822416626);
static const fftw_real K382683432 = FFTW_KONST(+0.382683432365089771728459984030398866761344562);
static const fftw_real K707106781 = FFTW_KONST(+0.707106781186547524400844362104849039284835938);

/*
 * Generator Id's :
 * $Id: exprdag.ml,v 1.36 1999/02/19 17:22:11 athena Exp $
 * $Id: fft.ml,v 1.41 1999/02/19 17:22:13 athena Exp $
 * $Id: to_c.ml,v 1.24 1999/02/19 17:22:17 athena Exp $
 */

void fftwi_twiddle_32(fftw_complex *A, const fftw_complex *W, int iostride, int m, int dist)
{
     int i;
     fftw_complex *inout;
     inout = A;
     for (i = m; i > 0; i = i - 1, inout = inout + dist, W = W + 31) {
	  fftw_real tmp19;
	  fftw_real tmp387;
	  fftw_real tmp472;
	  fftw_real tmp486;
	  fftw_real tmp442;
	  fftw_real tmp456;
	  fftw_real tmp191;
	  fftw_real tmp303;
	  fftw_real tmp161;
	  fftw_real tmp403;
	  fftw_real tmp276;
	  fftw_real tmp316;
	  fftw_real tmp372;
	  fftw_real tmp400;
	  fftw_real tmp259;
	  fftw_real tmp319;
	  fftw_real tmp42;
	  fftw_real tmp455;
	  fftw_real tmp201;
	  fftw_real tmp304;
	  fftw_real tmp390;
	  fftw_real tmp437;
	  fftw_real tmp196;
	  fftw_real tmp305;
	  fftw_real tmp184;
	  fftw_real tmp401;
	  fftw_real tmp375;
	  fftw_real tmp404;
	  fftw_real tmp270;
	  fftw_real tmp317;
	  fftw_real tmp279;
	  fftw_real tmp320;
	  fftw_real tmp66;
	  fftw_real tmp395;
	  fftw_real tmp224;
	  fftw_real tmp312;
	  fftw_real tmp357;
	  fftw_real tmp396;
	  fftw_real tmp219;
	  fftw_real tmp311;
	  fftw_real tmp114;
	  fftw_real tmp410;
	  fftw_real tmp249;
	  fftw_real tmp323;
	  fftw_real tmp363;
	  fftw_real tmp407;
	  fftw_real tmp232;
	  fftw_real tmp326;
	  fftw_real tmp89;
	  fftw_real tmp393;
	  fftw_real tmp213;
	  fftw_real tmp309;
	  fftw_real tmp354;
	  fftw_real tmp392;
	  fftw_real tmp208;
	  fftw_real tmp308;
	  fftw_real tmp137;
	  fftw_real tmp408;
	  fftw_real tmp366;
	  fftw_real tmp411;
	  fftw_real tmp243;
	  fftw_real tmp324;
	  fftw_real tmp252;
	  fftw_real tmp327;
	  ASSERT_ALIGNED_DOUBLE();
	  {
	       fftw_real tmp1;
	       fftw_real tmp440;
	       fftw_real tmp6;
	       fftw_real tmp439;
	       fftw_real tmp12;
	       fftw_real tmp188;
	       fftw_real tmp17;
	       fftw_real tmp189;
	       ASSERT_ALIGNED_DOUBLE();
	       tmp1 = c_re(inout[0]);
	       tmp440 = c_im(inout[0]);
	       {
		    fftw_real tmp3;
		    fftw_real tmp5;
		    fftw_real tmp2;
		    fftw_real tmp4;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp3 = c_re(inout[16 * iostride]);
		    tmp5 = c_im(inout[16 * iostride]);
		    tmp2 = c_re(W[15]);
		    tmp4 = c_im(W[15]);
		    tmp6 = (tmp2 * tmp3) + (tmp4 * tmp5);
		    tmp439 = (tmp2 * tmp5) - (tmp4 * tmp3);
	       }
	       {
		    fftw_real tmp9;
		    fftw_real tmp11;
		    fftw_real tmp8;
		    fftw_real tmp10;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp9 = c_re(inout[8 * iostride]);
		    tmp11 = c_im(inout[8 * iostride]);
		    tmp8 = c_re(W[7]);
		    tmp10 = c_im(W[7]);
		    tmp12 = (tmp8 * tmp9) + (tmp10 * tmp11);
		    tmp188 = (tmp8 * tmp11) - (tmp10 * tmp9);
	       }
	       {
		    fftw_real tmp14;
		    fftw_real tmp16;
		    fftw_real tmp13;
		    fftw_real tmp15;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp14 = c_re(inout[24 * iostride]);
		    tmp16 = c_im(inout[24 * iostride]);
		    tmp13 = c_re(W[23]);
		    tmp15 = c_im(W[23]);
		    tmp17 = (tmp13 * tmp14) + (tmp15 * tmp16);
		    tmp189 = (tmp13 * tmp16) - (tmp15 * tmp14);
	       }
	       {
		    fftw_real tmp7;
		    fftw_real tmp18;
		    fftw_real tmp470;
		    fftw_real tmp471;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp7 = tmp1 + tmp6;
		    tmp18 = tmp12 + tmp17;
		    tmp19 = tmp7 + tmp18;
		    tmp387 = tmp7 - tmp18;
		    tmp470 = tmp12 - tmp17;
		    tmp471 = tmp440 - tmp439;
		    tmp472 = tmp470 + tmp471;
		    tmp486 = tmp471 - tmp470;
	       }
	       {
		    fftw_real tmp438;
		    fftw_real tmp441;
		    fftw_real tmp187;
		    fftw_real tmp190;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp438 = tmp188 + tmp189;
		    tmp441 = tmp439 + tmp440;
		    tmp442 = tmp438 + tmp441;
		    tmp456 = tmp441 - tmp438;
		    tmp187 = tmp1 - tmp6;
		    tmp190 = tmp188 - tmp189;
		    tmp191 = tmp187 - tmp190;
		    tmp303 = tmp187 + tmp190;
	       }
	  }
	  {
	       fftw_real tmp143;
	       fftw_real tmp272;
	       fftw_real tmp159;
	       fftw_real tmp257;
	       fftw_real tmp148;
	       fftw_real tmp273;
	       fftw_real tmp154;
	       fftw_real tmp256;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp140;
		    fftw_real tmp142;
		    fftw_real tmp139;
		    fftw_real tmp141;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp140 = c_re(inout[31 * iostride]);
		    tmp142 = c_im(inout[31 * iostride]);
		    tmp139 = c_re(W[30]);
		    tmp141 = c_im(W[30]);
		    tmp143 = (tmp139 * tmp140) + (tmp141 * tmp142);
		    tmp272 = (tmp139 * tmp142) - (tmp141 * tmp140);
	       }
	       {
		    fftw_real tmp156;
		    fftw_real tmp158;
		    fftw_real tmp155;
		    fftw_real tmp157;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp156 = c_re(inout[23 * iostride]);
		    tmp158 = c_im(inout[23 * iostride]);
		    tmp155 = c_re(W[22]);
		    tmp157 = c_im(W[22]);
		    tmp159 = (tmp155 * tmp156) + (tmp157 * tmp158);
		    tmp257 = (tmp155 * tmp158) - (tmp157 * tmp156);
	       }
	       {
		    fftw_real tmp145;
		    fftw_real tmp147;
		    fftw_real tmp144;
		    fftw_real tmp146;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp145 = c_re(inout[15 * iostride]);
		    tmp147 = c_im(inout[15 * iostride]);
		    tmp144 = c_re(W[14]);
		    tmp146 = c_im(W[14]);
		    tmp148 = (tmp144 * tmp145) + (tmp146 * tmp147);
		    tmp273 = (tmp144 * tmp147) - (tmp146 * tmp145);
	       }
	       {
		    fftw_real tmp151;
		    fftw_real tmp153;
		    fftw_real tmp150;
		    fftw_real tmp152;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp151 = c_re(inout[7 * iostride]);
		    tmp153 = c_im(inout[7 * iostride]);
		    tmp150 = c_re(W[6]);
		    tmp152 = c_im(W[6]);
		    tmp154 = (tmp150 * tmp151) + (tmp152 * tmp153);
		    tmp256 = (tmp150 * tmp153) - (tmp152 * tmp151);
	       }
	       {
		    fftw_real tmp149;
		    fftw_real tmp160;
		    fftw_real tmp274;
		    fftw_real tmp275;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp149 = tmp143 + tmp148;
		    tmp160 = tmp154 + tmp159;
		    tmp161 = tmp149 + tmp160;
		    tmp403 = tmp149 - tmp160;
		    tmp274 = tmp272 - tmp273;
		    tmp275 = tmp154 - tmp159;
		    tmp276 = tmp274 + tmp275;
		    tmp316 = tmp274 - tmp275;
	       }
	       {
		    fftw_real tmp370;
		    fftw_real tmp371;
		    fftw_real tmp255;
		    fftw_real tmp258;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp370 = tmp272 + tmp273;
		    tmp371 = tmp256 + tmp257;
		    tmp372 = tmp370 + tmp371;
		    tmp400 = tmp370 - tmp371;
		    tmp255 = tmp143 - tmp148;
		    tmp258 = tmp256 - tmp257;
		    tmp259 = tmp255 - tmp258;
		    tmp319 = tmp255 + tmp258;
	       }
	  }
	  {
	       fftw_real tmp24;
	       fftw_real tmp193;
	       fftw_real tmp40;
	       fftw_real tmp199;
	       fftw_real tmp29;
	       fftw_real tmp194;
	       fftw_real tmp35;
	       fftw_real tmp198;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp21;
		    fftw_real tmp23;
		    fftw_real tmp20;
		    fftw_real tmp22;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp21 = c_re(inout[4 * iostride]);
		    tmp23 = c_im(inout[4 * iostride]);
		    tmp20 = c_re(W[3]);
		    tmp22 = c_im(W[3]);
		    tmp24 = (tmp20 * tmp21) + (tmp22 * tmp23);
		    tmp193 = (tmp20 * tmp23) - (tmp22 * tmp21);
	       }
	       {
		    fftw_real tmp37;
		    fftw_real tmp39;
		    fftw_real tmp36;
		    fftw_real tmp38;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp37 = c_re(inout[12 * iostride]);
		    tmp39 = c_im(inout[12 * iostride]);
		    tmp36 = c_re(W[11]);
		    tmp38 = c_im(W[11]);
		    tmp40 = (tmp36 * tmp37) + (tmp38 * tmp39);
		    tmp199 = (tmp36 * tmp39) - (tmp38 * tmp37);
	       }
	       {
		    fftw_real tmp26;
		    fftw_real tmp28;
		    fftw_real tmp25;
		    fftw_real tmp27;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp26 = c_re(inout[20 * iostride]);
		    tmp28 = c_im(inout[20 * iostride]);
		    tmp25 = c_re(W[19]);
		    tmp27 = c_im(W[19]);
		    tmp29 = (tmp25 * tmp26) + (tmp27 * tmp28);
		    tmp194 = (tmp25 * tmp28) - (tmp27 * tmp26);
	       }
	       {
		    fftw_real tmp32;
		    fftw_real tmp34;
		    fftw_real tmp31;
		    fftw_real tmp33;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp32 = c_re(inout[28 * iostride]);
		    tmp34 = c_im(inout[28 * iostride]);
		    tmp31 = c_re(W[27]);
		    tmp33 = c_im(W[27]);
		    tmp35 = (tmp31 * tmp32) + (tmp33 * tmp34);
		    tmp198 = (tmp31 * tmp34) - (tmp33 * tmp32);
	       }
	       {
		    fftw_real tmp30;
		    fftw_real tmp41;
		    fftw_real tmp197;
		    fftw_real tmp200;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp30 = tmp24 + tmp29;
		    tmp41 = tmp35 + tmp40;
		    tmp42 = tmp30 + tmp41;
		    tmp455 = tmp30 - tmp41;
		    tmp197 = tmp35 - tmp40;
		    tmp200 = tmp198 - tmp199;
		    tmp201 = tmp197 + tmp200;
		    tmp304 = tmp200 - tmp197;
	       }
	       {
		    fftw_real tmp388;
		    fftw_real tmp389;
		    fftw_real tmp192;
		    fftw_real tmp195;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp388 = tmp198 + tmp199;
		    tmp389 = tmp193 + tmp194;
		    tmp390 = tmp388 - tmp389;
		    tmp437 = tmp389 + tmp388;
		    tmp192 = tmp24 - tmp29;
		    tmp195 = tmp193 - tmp194;
		    tmp196 = tmp192 - tmp195;
		    tmp305 = tmp192 + tmp195;
	       }
	  }
	  {
	       fftw_real tmp166;
	       fftw_real tmp261;
	       fftw_real tmp171;
	       fftw_real tmp262;
	       fftw_real tmp260;
	       fftw_real tmp263;
	       fftw_real tmp177;
	       fftw_real tmp266;
	       fftw_real tmp182;
	       fftw_real tmp267;
	       fftw_real tmp265;
	       fftw_real tmp268;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp163;
		    fftw_real tmp165;
		    fftw_real tmp162;
		    fftw_real tmp164;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp163 = c_re(inout[3 * iostride]);
		    tmp165 = c_im(inout[3 * iostride]);
		    tmp162 = c_re(W[2]);
		    tmp164 = c_im(W[2]);
		    tmp166 = (tmp162 * tmp163) + (tmp164 * tmp165);
		    tmp261 = (tmp162 * tmp165) - (tmp164 * tmp163);
	       }
	       {
		    fftw_real tmp168;
		    fftw_real tmp170;
		    fftw_real tmp167;
		    fftw_real tmp169;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp168 = c_re(inout[19 * iostride]);
		    tmp170 = c_im(inout[19 * iostride]);
		    tmp167 = c_re(W[18]);
		    tmp169 = c_im(W[18]);
		    tmp171 = (tmp167 * tmp168) + (tmp169 * tmp170);
		    tmp262 = (tmp167 * tmp170) - (tmp169 * tmp168);
	       }
	       tmp260 = tmp166 - tmp171;
	       tmp263 = tmp261 - tmp262;
	       {
		    fftw_real tmp174;
		    fftw_real tmp176;
		    fftw_real tmp173;
		    fftw_real tmp175;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp174 = c_re(inout[27 * iostride]);
		    tmp176 = c_im(inout[27 * iostride]);
		    tmp173 = c_re(W[26]);
		    tmp175 = c_im(W[26]);
		    tmp177 = (tmp173 * tmp174) + (tmp175 * tmp176);
		    tmp266 = (tmp173 * tmp176) - (tmp175 * tmp174);
	       }
	       {
		    fftw_real tmp179;
		    fftw_real tmp181;
		    fftw_real tmp178;
		    fftw_real tmp180;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp179 = c_re(inout[11 * iostride]);
		    tmp181 = c_im(inout[11 * iostride]);
		    tmp178 = c_re(W[10]);
		    tmp180 = c_im(W[10]);
		    tmp182 = (tmp178 * tmp179) + (tmp180 * tmp181);
		    tmp267 = (tmp178 * tmp181) - (tmp180 * tmp179);
	       }
	       tmp265 = tmp177 - tmp182;
	       tmp268 = tmp266 - tmp267;
	       {
		    fftw_real tmp172;
		    fftw_real tmp183;
		    fftw_real tmp373;
		    fftw_real tmp374;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp172 = tmp166 + tmp171;
		    tmp183 = tmp177 + tmp182;
		    tmp184 = tmp172 + tmp183;
		    tmp401 = tmp172 - tmp183;
		    tmp373 = tmp261 + tmp262;
		    tmp374 = tmp266 + tmp267;
		    tmp375 = tmp373 + tmp374;
		    tmp404 = tmp374 - tmp373;
	       }
	       {
		    fftw_real tmp264;
		    fftw_real tmp269;
		    fftw_real tmp277;
		    fftw_real tmp278;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp264 = tmp260 - tmp263;
		    tmp269 = tmp265 + tmp268;
		    tmp270 = K707106781 * (tmp264 + tmp269);
		    tmp317 = K707106781 * (tmp264 - tmp269);
		    tmp277 = tmp260 + tmp263;
		    tmp278 = tmp268 - tmp265;
		    tmp279 = K707106781 * (tmp277 + tmp278);
		    tmp320 = K707106781 * (tmp278 - tmp277);
	       }
	  }
	  {
	       fftw_real tmp48;
	       fftw_real tmp215;
	       fftw_real tmp64;
	       fftw_real tmp222;
	       fftw_real tmp53;
	       fftw_real tmp216;
	       fftw_real tmp59;
	       fftw_real tmp221;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp45;
		    fftw_real tmp47;
		    fftw_real tmp44;
		    fftw_real tmp46;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp45 = c_re(inout[2 * iostride]);
		    tmp47 = c_im(inout[2 * iostride]);
		    tmp44 = c_re(W[1]);
		    tmp46 = c_im(W[1]);
		    tmp48 = (tmp44 * tmp45) + (tmp46 * tmp47);
		    tmp215 = (tmp44 * tmp47) - (tmp46 * tmp45);
	       }
	       {
		    fftw_real tmp61;
		    fftw_real tmp63;
		    fftw_real tmp60;
		    fftw_real tmp62;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp61 = c_re(inout[26 * iostride]);
		    tmp63 = c_im(inout[26 * iostride]);
		    tmp60 = c_re(W[25]);
		    tmp62 = c_im(W[25]);
		    tmp64 = (tmp60 * tmp61) + (tmp62 * tmp63);
		    tmp222 = (tmp60 * tmp63) - (tmp62 * tmp61);
	       }
	       {
		    fftw_real tmp50;
		    fftw_real tmp52;
		    fftw_real tmp49;
		    fftw_real tmp51;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp50 = c_re(inout[18 * iostride]);
		    tmp52 = c_im(inout[18 * iostride]);
		    tmp49 = c_re(W[17]);
		    tmp51 = c_im(W[17]);
		    tmp53 = (tmp49 * tmp50) + (tmp51 * tmp52);
		    tmp216 = (tmp49 * tmp52) - (tmp51 * tmp50);
	       }
	       {
		    fftw_real tmp56;
		    fftw_real tmp58;
		    fftw_real tmp55;
		    fftw_real tmp57;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp56 = c_re(inout[10 * iostride]);
		    tmp58 = c_im(inout[10 * iostride]);
		    tmp55 = c_re(W[9]);
		    tmp57 = c_im(W[9]);
		    tmp59 = (tmp55 * tmp56) + (tmp57 * tmp58);
		    tmp221 = (tmp55 * tmp58) - (tmp57 * tmp56);
	       }
	       {
		    fftw_real tmp54;
		    fftw_real tmp65;
		    fftw_real tmp220;
		    fftw_real tmp223;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp54 = tmp48 + tmp53;
		    tmp65 = tmp59 + tmp64;
		    tmp66 = tmp54 + tmp65;
		    tmp395 = tmp54 - tmp65;
		    tmp220 = tmp48 - tmp53;
		    tmp223 = tmp221 - tmp222;
		    tmp224 = tmp220 - tmp223;
		    tmp312 = tmp220 + tmp223;
	       }
	       {
		    fftw_real tmp355;
		    fftw_real tmp356;
		    fftw_real tmp217;
		    fftw_real tmp218;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp355 = tmp215 + tmp216;
		    tmp356 = tmp221 + tmp222;
		    tmp357 = tmp355 + tmp356;
		    tmp396 = tmp355 - tmp356;
		    tmp217 = tmp215 - tmp216;
		    tmp218 = tmp59 - tmp64;
		    tmp219 = tmp217 + tmp218;
		    tmp311 = tmp217 - tmp218;
	       }
	  }
	  {
	       fftw_real tmp96;
	       fftw_real tmp245;
	       fftw_real tmp112;
	       fftw_real tmp230;
	       fftw_real tmp101;
	       fftw_real tmp246;
	       fftw_real tmp107;
	       fftw_real tmp229;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp93;
		    fftw_real tmp95;
		    fftw_real tmp92;
		    fftw_real tmp94;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp93 = c_re(inout[iostride]);
		    tmp95 = c_im(inout[iostride]);
		    tmp92 = c_re(W[0]);
		    tmp94 = c_im(W[0]);
		    tmp96 = (tmp92 * tmp93) + (tmp94 * tmp95);
		    tmp245 = (tmp92 * tmp95) - (tmp94 * tmp93);
	       }
	       {
		    fftw_real tmp109;
		    fftw_real tmp111;
		    fftw_real tmp108;
		    fftw_real tmp110;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp109 = c_re(inout[25 * iostride]);
		    tmp111 = c_im(inout[25 * iostride]);
		    tmp108 = c_re(W[24]);
		    tmp110 = c_im(W[24]);
		    tmp112 = (tmp108 * tmp109) + (tmp110 * tmp111);
		    tmp230 = (tmp108 * tmp111) - (tmp110 * tmp109);
	       }
	       {
		    fftw_real tmp98;
		    fftw_real tmp100;
		    fftw_real tmp97;
		    fftw_real tmp99;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp98 = c_re(inout[17 * iostride]);
		    tmp100 = c_im(inout[17 * iostride]);
		    tmp97 = c_re(W[16]);
		    tmp99 = c_im(W[16]);
		    tmp101 = (tmp97 * tmp98) + (tmp99 * tmp100);
		    tmp246 = (tmp97 * tmp100) - (tmp99 * tmp98);
	       }
	       {
		    fftw_real tmp104;
		    fftw_real tmp106;
		    fftw_real tmp103;
		    fftw_real tmp105;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp104 = c_re(inout[9 * iostride]);
		    tmp106 = c_im(inout[9 * iostride]);
		    tmp103 = c_re(W[8]);
		    tmp105 = c_im(W[8]);
		    tmp107 = (tmp103 * tmp104) + (tmp105 * tmp106);
		    tmp229 = (tmp103 * tmp106) - (tmp105 * tmp104);
	       }
	       {
		    fftw_real tmp102;
		    fftw_real tmp113;
		    fftw_real tmp247;
		    fftw_real tmp248;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp102 = tmp96 + tmp101;
		    tmp113 = tmp107 + tmp112;
		    tmp114 = tmp102 + tmp113;
		    tmp410 = tmp102 - tmp113;
		    tmp247 = tmp245 - tmp246;
		    tmp248 = tmp107 - tmp112;
		    tmp249 = tmp247 + tmp248;
		    tmp323 = tmp247 - tmp248;
	       }
	       {
		    fftw_real tmp361;
		    fftw_real tmp362;
		    fftw_real tmp228;
		    fftw_real tmp231;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp361 = tmp245 + tmp246;
		    tmp362 = tmp229 + tmp230;
		    tmp363 = tmp361 + tmp362;
		    tmp407 = tmp361 - tmp362;
		    tmp228 = tmp96 - tmp101;
		    tmp231 = tmp229 - tmp230;
		    tmp232 = tmp228 - tmp231;
		    tmp326 = tmp228 + tmp231;
	       }
	  }
	  {
	       fftw_real tmp71;
	       fftw_real tmp204;
	       fftw_real tmp87;
	       fftw_real tmp211;
	       fftw_real tmp76;
	       fftw_real tmp205;
	       fftw_real tmp82;
	       fftw_real tmp210;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp68;
		    fftw_real tmp70;
		    fftw_real tmp67;
		    fftw_real tmp69;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp68 = c_re(inout[30 * iostride]);
		    tmp70 = c_im(inout[30 * iostride]);
		    tmp67 = c_re(W[29]);
		    tmp69 = c_im(W[29]);
		    tmp71 = (tmp67 * tmp68) + (tmp69 * tmp70);
		    tmp204 = (tmp67 * tmp70) - (tmp69 * tmp68);
	       }
	       {
		    fftw_real tmp84;
		    fftw_real tmp86;
		    fftw_real tmp83;
		    fftw_real tmp85;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp84 = c_re(inout[22 * iostride]);
		    tmp86 = c_im(inout[22 * iostride]);
		    tmp83 = c_re(W[21]);
		    tmp85 = c_im(W[21]);
		    tmp87 = (tmp83 * tmp84) + (tmp85 * tmp86);
		    tmp211 = (tmp83 * tmp86) - (tmp85 * tmp84);
	       }
	       {
		    fftw_real tmp73;
		    fftw_real tmp75;
		    fftw_real tmp72;
		    fftw_real tmp74;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp73 = c_re(inout[14 * iostride]);
		    tmp75 = c_im(inout[14 * iostride]);
		    tmp72 = c_re(W[13]);
		    tmp74 = c_im(W[13]);
		    tmp76 = (tmp72 * tmp73) + (tmp74 * tmp75);
		    tmp205 = (tmp72 * tmp75) - (tmp74 * tmp73);
	       }
	       {
		    fftw_real tmp79;
		    fftw_real tmp81;
		    fftw_real tmp78;
		    fftw_real tmp80;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp79 = c_re(inout[6 * iostride]);
		    tmp81 = c_im(inout[6 * iostride]);
		    tmp78 = c_re(W[5]);
		    tmp80 = c_im(W[5]);
		    tmp82 = (tmp78 * tmp79) + (tmp80 * tmp81);
		    tmp210 = (tmp78 * tmp81) - (tmp80 * tmp79);
	       }
	       {
		    fftw_real tmp77;
		    fftw_real tmp88;
		    fftw_real tmp209;
		    fftw_real tmp212;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp77 = tmp71 + tmp76;
		    tmp88 = tmp82 + tmp87;
		    tmp89 = tmp77 + tmp88;
		    tmp393 = tmp77 - tmp88;
		    tmp209 = tmp71 - tmp76;
		    tmp212 = tmp210 - tmp211;
		    tmp213 = tmp209 - tmp212;
		    tmp309 = tmp209 + tmp212;
	       }
	       {
		    fftw_real tmp352;
		    fftw_real tmp353;
		    fftw_real tmp206;
		    fftw_real tmp207;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp352 = tmp204 + tmp205;
		    tmp353 = tmp210 + tmp211;
		    tmp354 = tmp352 + tmp353;
		    tmp392 = tmp352 - tmp353;
		    tmp206 = tmp204 - tmp205;
		    tmp207 = tmp82 - tmp87;
		    tmp208 = tmp206 + tmp207;
		    tmp308 = tmp206 - tmp207;
	       }
	  }
	  {
	       fftw_real tmp119;
	       fftw_real tmp234;
	       fftw_real tmp124;
	       fftw_real tmp235;
	       fftw_real tmp233;
	       fftw_real tmp236;
	       fftw_real tmp130;
	       fftw_real tmp239;
	       fftw_real tmp135;
	       fftw_real tmp240;
	       fftw_real tmp238;
	       fftw_real tmp241;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp116;
		    fftw_real tmp118;
		    fftw_real tmp115;
		    fftw_real tmp117;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp116 = c_re(inout[5 * iostride]);
		    tmp118 = c_im(inout[5 * iostride]);
		    tmp115 = c_re(W[4]);
		    tmp117 = c_im(W[4]);
		    tmp119 = (tmp115 * tmp116) + (tmp117 * tmp118);
		    tmp234 = (tmp115 * tmp118) - (tmp117 * tmp116);
	       }
	       {
		    fftw_real tmp121;
		    fftw_real tmp123;
		    fftw_real tmp120;
		    fftw_real tmp122;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp121 = c_re(inout[21 * iostride]);
		    tmp123 = c_im(inout[21 * iostride]);
		    tmp120 = c_re(W[20]);
		    tmp122 = c_im(W[20]);
		    tmp124 = (tmp120 * tmp121) + (tmp122 * tmp123);
		    tmp235 = (tmp120 * tmp123) - (tmp122 * tmp121);
	       }
	       tmp233 = tmp119 - tmp124;
	       tmp236 = tmp234 - tmp235;
	       {
		    fftw_real tmp127;
		    fftw_real tmp129;
		    fftw_real tmp126;
		    fftw_real tmp128;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp127 = c_re(inout[29 * iostride]);
		    tmp129 = c_im(inout[29 * iostride]);
		    tmp126 = c_re(W[28]);
		    tmp128 = c_im(W[28]);
		    tmp130 = (tmp126 * tmp127) + (tmp128 * tmp129);
		    tmp239 = (tmp126 * tmp129) - (tmp128 * tmp127);
	       }
	       {
		    fftw_real tmp132;
		    fftw_real tmp134;
		    fftw_real tmp131;
		    fftw_real tmp133;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp132 = c_re(inout[13 * iostride]);
		    tmp134 = c_im(inout[13 * iostride]);
		    tmp131 = c_re(W[12]);
		    tmp133 = c_im(W[12]);
		    tmp135 = (tmp131 * tmp132) + (tmp133 * tmp134);
		    tmp240 = (tmp131 * tmp134) - (tmp133 * tmp132);
	       }
	       tmp238 = tmp130 - tmp135;
	       tmp241 = tmp239 - tmp240;
	       {
		    fftw_real tmp125;
		    fftw_real tmp136;
		    fftw_real tmp364;
		    fftw_real tmp365;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp125 = tmp119 + tmp124;
		    tmp136 = tmp130 + tmp135;
		    tmp137 = tmp125 + tmp136;
		    tmp408 = tmp125 - tmp136;
		    tmp364 = tmp234 + tmp235;
		    tmp365 = tmp239 + tmp240;
		    tmp366 = tmp364 + tmp365;
		    tmp411 = tmp365 - tmp364;
	       }
	       {
		    fftw_real tmp237;
		    fftw_real tmp242;
		    fftw_real tmp250;
		    fftw_real tmp251;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp237 = tmp233 - tmp236;
		    tmp242 = tmp238 + tmp241;
		    tmp243 = K707106781 * (tmp237 + tmp242);
		    tmp324 = K707106781 * (tmp237 - tmp242);
		    tmp250 = tmp233 + tmp236;
		    tmp251 = tmp241 - tmp238;
		    tmp252 = K707106781 * (tmp250 + tmp251);
		    tmp327 = K707106781 * (tmp251 - tmp250);
	       }
	  }
	  {
	       fftw_real tmp91;
	       fftw_real tmp383;
	       fftw_real tmp444;
	       fftw_real tmp446;
	       fftw_real tmp186;
	       fftw_real tmp445;
	       fftw_real tmp386;
	       fftw_real tmp435;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp43;
		    fftw_real tmp90;
		    fftw_real tmp436;
		    fftw_real tmp443;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp43 = tmp19 + tmp42;
		    tmp90 = tmp66 + tmp89;
		    tmp91 = tmp43 + tmp90;
		    tmp383 = tmp43 - tmp90;
		    tmp436 = tmp357 + tmp354;
		    tmp443 = tmp437 + tmp442;
		    tmp444 = tmp436 + tmp443;
		    tmp446 = tmp443 - tmp436;
	       }
	       {
		    fftw_real tmp138;
		    fftw_real tmp185;
		    fftw_real tmp384;
		    fftw_real tmp385;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp138 = tmp114 + tmp137;
		    tmp185 = tmp161 + tmp184;
		    tmp186 = tmp138 + tmp185;
		    tmp445 = tmp138 - tmp185;
		    tmp384 = tmp372 + tmp375;
		    tmp385 = tmp363 + tmp366;
		    tmp386 = tmp384 - tmp385;
		    tmp435 = tmp385 + tmp384;
	       }
	       c_re(inout[16 * iostride]) = tmp91 - tmp186;
	       c_re(inout[0]) = tmp91 + tmp186;
	       c_re(inout[24 * iostride]) = tmp383 - tmp386;
	       c_re(inout[8 * iostride]) = tmp383 + tmp386;
	       c_im(inout[0]) = tmp435 + tmp444;
	       c_im(inout[16 * iostride]) = tmp444 - tmp435;
	       c_im(inout[8 * iostride]) = tmp445 + tmp446;
	       c_im(inout[24 * iostride]) = tmp446 - tmp445;
	  }
	  {
	       fftw_real tmp359;
	       fftw_real tmp379;
	       fftw_real tmp450;
	       fftw_real tmp452;
	       fftw_real tmp368;
	       fftw_real tmp381;
	       fftw_real tmp377;
	       fftw_real tmp380;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp351;
		    fftw_real tmp358;
		    fftw_real tmp448;
		    fftw_real tmp449;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp351 = tmp19 - tmp42;
		    tmp358 = tmp354 - tmp357;
		    tmp359 = tmp351 + tmp358;
		    tmp379 = tmp351 - tmp358;
		    tmp448 = tmp66 - tmp89;
		    tmp449 = tmp442 - tmp437;
		    tmp450 = tmp448 + tmp449;
		    tmp452 = tmp449 - tmp448;
	       }
	       {
		    fftw_real tmp360;
		    fftw_real tmp367;
		    fftw_real tmp369;
		    fftw_real tmp376;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp360 = tmp114 - tmp137;
		    tmp367 = tmp363 - tmp366;
		    tmp368 = tmp360 - tmp367;
		    tmp381 = tmp360 + tmp367;
		    tmp369 = tmp161 - tmp184;
		    tmp376 = tmp372 - tmp375;
		    tmp377 = tmp369 + tmp376;
		    tmp380 = tmp376 - tmp369;
	       }
	       {
		    fftw_real tmp378;
		    fftw_real tmp451;
		    fftw_real tmp382;
		    fftw_real tmp447;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp378 = K707106781 * (tmp368 + tmp377);
		    c_re(inout[20 * iostride]) = tmp359 - tmp378;
		    c_re(inout[4 * iostride]) = tmp359 + tmp378;
		    tmp451 = K707106781 * (tmp368 - tmp377);
		    c_im(inout[12 * iostride]) = tmp451 + tmp452;
		    c_im(inout[28 * iostride]) = tmp452 - tmp451;
		    tmp382 = K707106781 * (tmp380 - tmp381);
		    c_re(inout[28 * iostride]) = tmp379 - tmp382;
		    c_re(inout[12 * iostride]) = tmp379 + tmp382;
		    tmp447 = K707106781 * (tmp381 + tmp380);
		    c_im(inout[4 * iostride]) = tmp447 + tmp450;
		    c_im(inout[20 * iostride]) = tmp450 - tmp447;
	       }
	  }
	  {
	       fftw_real tmp391;
	       fftw_real tmp419;
	       fftw_real tmp398;
	       fftw_real tmp454;
	       fftw_real tmp422;
	       fftw_real tmp462;
	       fftw_real tmp406;
	       fftw_real tmp417;
	       fftw_real tmp457;
	       fftw_real tmp463;
	       fftw_real tmp426;
	       fftw_real tmp433;
	       fftw_real tmp413;
	       fftw_real tmp416;
	       fftw_real tmp429;
	       fftw_real tmp432;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp394;
		    fftw_real tmp397;
		    fftw_real tmp424;
		    fftw_real tmp425;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp391 = tmp387 - tmp390;
		    tmp419 = tmp387 + tmp390;
		    tmp394 = tmp392 - tmp393;
		    tmp397 = tmp395 + tmp396;
		    tmp398 = K707106781 * (tmp394 - tmp397);
		    tmp454 = K707106781 * (tmp397 + tmp394);
		    {
			 fftw_real tmp420;
			 fftw_real tmp421;
			 fftw_real tmp402;
			 fftw_real tmp405;
			 ASSERT_ALIGNED_DOUBLE();
			 tmp420 = tmp395 - tmp396;
			 tmp421 = tmp393 + tmp392;
			 tmp422 = K707106781 * (tmp420 + tmp421);
			 tmp462 = K707106781 * (tmp420 - tmp421);
			 tmp402 = tmp400 - tmp401;
			 tmp405 = tmp403 - tmp404;
			 tmp406 = (K382683432 * tmp402) - (K923879532 * tmp405);
			 tmp417 = (K923879532 * tmp402) + (K382683432 * tmp405);
		    }
		    tmp457 = tmp455 + tmp456;
		    tmp463 = tmp456 - tmp455;
		    tmp424 = tmp400 + tmp401;
		    tmp425 = tmp403 + tmp404;
		    tmp426 = (K923879532 * tmp424) - (K382683432 * tmp425);
		    tmp433 = (K382683432 * tmp424) + (K923879532 * tmp425);
		    {
			 fftw_real tmp409;
			 fftw_real tmp412;
			 fftw_real tmp427;
			 fftw_real tmp428;
			 ASSERT_ALIGNED_DOUBLE();
			 tmp409 = tmp407 - tmp408;
			 tmp412 = tmp410 - tmp411;
			 tmp413 = (K382683432 * tmp409) + (K923879532 * tmp412);
			 tmp416 = (K382683432 * tmp412) - (K923879532 * tmp409);
			 tmp427 = tmp407 + tmp408;
			 tmp428 = tmp410 + tmp411;
			 tmp429 = (K923879532 * tmp427) + (K382683432 * tmp428);
			 tmp432 = (K923879532 * tmp428) - (K382683432 * tmp427);
		    }
	       }
	       {
		    fftw_real tmp399;
		    fftw_real tmp414;
		    fftw_real tmp415;
		    fftw_real tmp418;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp399 = tmp391 - tmp398;
		    tmp414 = tmp406 - tmp413;
		    c_re(inout[30 * iostride]) = tmp399 - tmp414;
		    c_re(inout[14 * iostride]) = tmp399 + tmp414;
		    tmp415 = tmp391 + tmp398;
		    tmp418 = tmp416 + tmp417;
		    c_re(inout[22 * iostride]) = tmp415 - tmp418;
		    c_re(inout[6 * iostride]) = tmp415 + tmp418;
	       }
	       {
		    fftw_real tmp465;
		    fftw_real tmp466;
		    fftw_real tmp461;
		    fftw_real tmp464;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp465 = tmp416 - tmp417;
		    tmp466 = tmp463 - tmp462;
		    c_im(inout[14 * iostride]) = tmp465 + tmp466;
		    c_im(inout[30 * iostride]) = tmp466 - tmp465;
		    tmp461 = tmp413 + tmp406;
		    tmp464 = tmp462 + tmp463;
		    c_im(inout[6 * iostride]) = tmp461 + tmp464;
		    c_im(inout[22 * iostride]) = tmp464 - tmp461;
	       }
	       {
		    fftw_real tmp423;
		    fftw_real tmp430;
		    fftw_real tmp431;
		    fftw_real tmp434;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp423 = tmp419 - tmp422;
		    tmp430 = tmp426 - tmp429;
		    c_re(inout[26 * iostride]) = tmp423 - tmp430;
		    c_re(inout[10 * iostride]) = tmp423 + tmp430;
		    tmp431 = tmp419 + tmp422;
		    tmp434 = tmp432 + tmp433;
		    c_re(inout[18 * iostride]) = tmp431 - tmp434;
		    c_re(inout[2 * iostride]) = tmp431 + tmp434;
	       }
	       {
		    fftw_real tmp459;
		    fftw_real tmp460;
		    fftw_real tmp453;
		    fftw_real tmp458;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp459 = tmp432 - tmp433;
		    tmp460 = tmp457 - tmp454;
		    c_im(inout[10 * iostride]) = tmp459 + tmp460;
		    c_im(inout[26 * iostride]) = tmp460 - tmp459;
		    tmp453 = tmp429 + tmp426;
		    tmp458 = tmp454 + tmp457;
		    c_im(inout[2 * iostride]) = tmp453 + tmp458;
		    c_im(inout[18 * iostride]) = tmp458 - tmp453;
	       }
	  }
	  {
	       fftw_real tmp307;
	       fftw_real tmp335;
	       fftw_real tmp338;
	       fftw_real tmp492;
	       fftw_real tmp487;
	       fftw_real tmp493;
	       fftw_real tmp314;
	       fftw_real tmp484;
	       fftw_real tmp322;
	       fftw_real tmp333;
	       fftw_real tmp342;
	       fftw_real tmp349;
	       fftw_real tmp329;
	       fftw_real tmp332;
	       fftw_real tmp345;
	       fftw_real tmp348;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp306;
		    fftw_real tmp336;
		    fftw_real tmp337;
		    fftw_real tmp485;
		    fftw_real tmp310;
		    fftw_real tmp313;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp306 = K707106781 * (tmp304 - tmp305);
		    tmp307 = tmp303 - tmp306;
		    tmp335 = tmp303 + tmp306;
		    tmp336 = (K382683432 * tmp312) - (K923879532 * tmp311);
		    tmp337 = (K923879532 * tmp308) + (K382683432 * tmp309);
		    tmp338 = tmp336 + tmp337;
		    tmp492 = tmp336 - tmp337;
		    tmp485 = K707106781 * (tmp196 - tmp201);
		    tmp487 = tmp485 + tmp486;
		    tmp493 = tmp486 - tmp485;
		    tmp310 = (K382683432 * tmp308) - (K923879532 * tmp309);
		    tmp313 = (K382683432 * tmp311) + (K923879532 * tmp312);
		    tmp314 = tmp310 - tmp313;
		    tmp484 = tmp313 + tmp310;
	       }
	       {
		    fftw_real tmp318;
		    fftw_real tmp321;
		    fftw_real tmp340;
		    fftw_real tmp341;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp318 = tmp316 - tmp317;
		    tmp321 = tmp319 - tmp320;
		    tmp322 = (K195090322 * tmp318) - (K980785280 * tmp321);
		    tmp333 = (K980785280 * tmp318) + (K195090322 * tmp321);
		    tmp340 = tmp316 + tmp317;
		    tmp341 = tmp319 + tmp320;
		    tmp342 = (K831469612 * tmp340) - (K555570233 * tmp341);
		    tmp349 = (K555570233 * tmp340) + (K831469612 * tmp341);
	       }
	       {
		    fftw_real tmp325;
		    fftw_real tmp328;
		    fftw_real tmp343;
		    fftw_real tmp344;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp325 = tmp323 - tmp324;
		    tmp328 = tmp326 - tmp327;
		    tmp329 = (K195090322 * tmp325) + (K980785280 * tmp328);
		    tmp332 = (K195090322 * tmp328) - (K980785280 * tmp325);
		    tmp343 = tmp323 + tmp324;
		    tmp344 = tmp326 + tmp327;
		    tmp345 = (K831469612 * tmp343) + (K555570233 * tmp344);
		    tmp348 = (K831469612 * tmp344) - (K555570233 * tmp343);
	       }
	       {
		    fftw_real tmp315;
		    fftw_real tmp330;
		    fftw_real tmp331;
		    fftw_real tmp334;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp315 = tmp307 - tmp314;
		    tmp330 = tmp322 - tmp329;
		    c_re(inout[31 * iostride]) = tmp315 - tmp330;
		    c_re(inout[15 * iostride]) = tmp315 + tmp330;
		    tmp331 = tmp307 + tmp314;
		    tmp334 = tmp332 + tmp333;
		    c_re(inout[23 * iostride]) = tmp331 - tmp334;
		    c_re(inout[7 * iostride]) = tmp331 + tmp334;
	       }
	       {
		    fftw_real tmp495;
		    fftw_real tmp496;
		    fftw_real tmp491;
		    fftw_real tmp494;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp495 = tmp332 - tmp333;
		    tmp496 = tmp493 - tmp492;
		    c_im(inout[15 * iostride]) = tmp495 + tmp496;
		    c_im(inout[31 * iostride]) = tmp496 - tmp495;
		    tmp491 = tmp329 + tmp322;
		    tmp494 = tmp492 + tmp493;
		    c_im(inout[7 * iostride]) = tmp491 + tmp494;
		    c_im(inout[23 * iostride]) = tmp494 - tmp491;
	       }
	       {
		    fftw_real tmp339;
		    fftw_real tmp346;
		    fftw_real tmp347;
		    fftw_real tmp350;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp339 = tmp335 - tmp338;
		    tmp346 = tmp342 - tmp345;
		    c_re(inout[27 * iostride]) = tmp339 - tmp346;
		    c_re(inout[11 * iostride]) = tmp339 + tmp346;
		    tmp347 = tmp335 + tmp338;
		    tmp350 = tmp348 + tmp349;
		    c_re(inout[19 * iostride]) = tmp347 - tmp350;
		    c_re(inout[3 * iostride]) = tmp347 + tmp350;
	       }
	       {
		    fftw_real tmp489;
		    fftw_real tmp490;
		    fftw_real tmp483;
		    fftw_real tmp488;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp489 = tmp348 - tmp349;
		    tmp490 = tmp487 - tmp484;
		    c_im(inout[11 * iostride]) = tmp489 + tmp490;
		    c_im(inout[27 * iostride]) = tmp490 - tmp489;
		    tmp483 = tmp345 + tmp342;
		    tmp488 = tmp484 + tmp487;
		    c_im(inout[3 * iostride]) = tmp483 + tmp488;
		    c_im(inout[19 * iostride]) = tmp488 - tmp483;
	       }
	  }
	  {
	       fftw_real tmp203;
	       fftw_real tmp287;
	       fftw_real tmp290;
	       fftw_real tmp478;
	       fftw_real tmp473;
	       fftw_real tmp479;
	       fftw_real tmp226;
	       fftw_real tmp468;
	       fftw_real tmp254;
	       fftw_real tmp285;
	       fftw_real tmp294;
	       fftw_real tmp301;
	       fftw_real tmp281;
	       fftw_real tmp284;
	       fftw_real tmp297;
	       fftw_real tmp300;
	       ASSERT_ALIGNED_DOUBLE();
	       {
		    fftw_real tmp202;
		    fftw_real tmp288;
		    fftw_real tmp289;
		    fftw_real tmp469;
		    fftw_real tmp214;
		    fftw_real tmp225;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp202 = K707106781 * (tmp196 + tmp201);
		    tmp203 = tmp191 - tmp202;
		    tmp287 = tmp191 + tmp202;
		    tmp288 = (K923879532 * tmp224) - (K382683432 * tmp219);
		    tmp289 = (K382683432 * tmp208) + (K923879532 * tmp213);
		    tmp290 = tmp288 + tmp289;
		    tmp478 = tmp288 - tmp289;
		    tmp469 = K707106781 * (tmp305 + tmp304);
		    tmp473 = tmp469 + tmp472;
		    tmp479 = tmp472 - tmp469;
		    tmp214 = (K923879532 * tmp208) - (K382683432 * tmp213);
		    tmp225 = (K923879532 * tmp219) + (K382683432 * tmp224);
		    tmp226 = tmp214 - tmp225;
		    tmp468 = tmp225 + tmp214;
	       }
	       {
		    fftw_real tmp244;
		    fftw_real tmp253;
		    fftw_real tmp292;
		    fftw_real tmp293;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp244 = tmp232 - tmp243;
		    tmp253 = tmp249 - tmp252;
		    tmp254 = (K555570233 * tmp244) - (K831469612 * tmp253);
		    tmp285 = (K831469612 * tmp244) + (K555570233 * tmp253);
		    tmp292 = tmp232 + tmp243;
		    tmp293 = tmp249 + tmp252;
		    tmp294 = (K980785280 * tmp292) - (K195090322 * tmp293);
		    tmp301 = (K195090322 * tmp292) + (K980785280 * tmp293);
	       }
	       {
		    fftw_real tmp271;
		    fftw_real tmp280;
		    fftw_real tmp295;
		    fftw_real tmp296;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp271 = tmp259 - tmp270;
		    tmp280 = tmp276 - tmp279;
		    tmp281 = (K555570233 * tmp271) + (K831469612 * tmp280);
		    tmp284 = (K555570233 * tmp280) - (K831469612 * tmp271);
		    tmp295 = tmp259 + tmp270;
		    tmp296 = tmp276 + tmp279;
		    tmp297 = (K980785280 * tmp295) + (K195090322 * tmp296);
		    tmp300 = (K980785280 * tmp296) - (K195090322 * tmp295);
	       }
	       {
		    fftw_real tmp227;
		    fftw_real tmp282;
		    fftw_real tmp283;
		    fftw_real tmp286;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp227 = tmp203 + tmp226;
		    tmp282 = tmp254 + tmp281;
		    c_re(inout[21 * iostride]) = tmp227 - tmp282;
		    c_re(inout[5 * iostride]) = tmp227 + tmp282;
		    tmp283 = tmp203 - tmp226;
		    tmp286 = tmp284 - tmp285;
		    c_re(inout[29 * iostride]) = tmp283 - tmp286;
		    c_re(inout[13 * iostride]) = tmp283 + tmp286;
	       }
	       {
		    fftw_real tmp477;
		    fftw_real tmp480;
		    fftw_real tmp481;
		    fftw_real tmp482;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp477 = tmp285 + tmp284;
		    tmp480 = tmp478 + tmp479;
		    c_im(inout[5 * iostride]) = tmp477 + tmp480;
		    c_im(inout[21 * iostride]) = tmp480 - tmp477;
		    tmp481 = tmp254 - tmp281;
		    tmp482 = tmp479 - tmp478;
		    c_im(inout[13 * iostride]) = tmp481 + tmp482;
		    c_im(inout[29 * iostride]) = tmp482 - tmp481;
	       }
	       {
		    fftw_real tmp291;
		    fftw_real tmp298;
		    fftw_real tmp299;
		    fftw_real tmp302;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp291 = tmp287 + tmp290;
		    tmp298 = tmp294 + tmp297;
		    c_re(inout[17 * iostride]) = tmp291 - tmp298;
		    c_re(inout[iostride]) = tmp291 + tmp298;
		    tmp299 = tmp287 - tmp290;
		    tmp302 = tmp300 - tmp301;
		    c_re(inout[25 * iostride]) = tmp299 - tmp302;
		    c_re(inout[9 * iostride]) = tmp299 + tmp302;
	       }
	       {
		    fftw_real tmp467;
		    fftw_real tmp474;
		    fftw_real tmp475;
		    fftw_real tmp476;
		    ASSERT_ALIGNED_DOUBLE();
		    tmp467 = tmp301 + tmp300;
		    tmp474 = tmp468 + tmp473;
		    c_im(inout[iostride]) = tmp467 + tmp474;
		    c_im(inout[17 * iostride]) = tmp474 - tmp467;
		    tmp475 = tmp294 - tmp297;
		    tmp476 = tmp473 - tmp468;
		    c_im(inout[9 * iostride]) = tmp475 + tmp476;
		    c_im(inout[25 * iostride]) = tmp476 - tmp475;
	       }
	  }
     }
}

static const int twiddle_order[] =
{1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31};
fftw_codelet_desc fftwi_twiddle_32_desc =
{
     "fftwi_twiddle_32",
     (void (*)()) fftwi_twiddle_32,
     32,
     FFTW_BACKWARD,
     FFTW_TWIDDLE,
     520,
     31,
     twiddle_order,
};

max_stars_repo_name: jlprieur/jlplib

max_issues_repo_name: jlprieur/jlplib

max_forks_repo_name: jlprieur/jlplib


text:
#pragma once

#include <memory>
#include <gsl/gsl>

/// Holds code related to computer memory.
namespace MEMORY
{
    /// A type alias to improve readability for non-null raw pointers.
    /// @tparam Type - The type of the underyling item being pointed to.
    template <typename Type>
    using NonNullRawPointer = gsl::strict_not_null<Type*>;

    /// A type alias to improve readability for non-null unique pointers.
    /// @tparam Type - The type of the underlying item in the pointer.
    template <typename Type>
    using NonNullUniquePointer = gsl::strict_not_null<std::unique_ptr<Type>>;

    /// A type alias to improve readability for non-null shared pointers.
    /// @tparam Type - The type of the underlying item in the pointer.
    template <typename Type>
    using NonNullSharedPointer = gsl::strict_not_null<std::shared_ptr<Type>>;
}

max_stars_repo_name: jpike/CppLibraries

max_issues_repo_name: jpike/CppLibraries

max_forks_repo_name: jpike/CppLibraries


text:
#include <stdio.h>
#include <stdlib.h>
#include <stdbool.h>
#include <math.h>
#include <string.h>

#ifdef ACCELERATE
#include <Accelerate/Accelerate.h>
#else
#include <cblas.h>
#endif

#include "neural.h"

void
neural_feed_forward(struct neural_net_layer *head, float *x, int batch_size)
{
	/* point "current" layer to head and set input as feature array x */
	struct neural_net_layer *curr = head;
	struct neural_net_layer *next;
	curr->act = x;

	int in_nodes;
	int out_nodes;

	while ((next = curr->next) != NULL) {
		in_nodes = curr->in_nodes;
		out_nodes = curr->out_nodes;

		/* act = X * Theta^T + bias */
		for (int i = 0; i < batch_size; i++) {
			memcpy(next->act + i * out_nodes, curr->bias, out_nodes * sizeof(float));
		}
		cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasTrans, batch_size, out_nodes, in_nodes,
			1.0, curr->act, in_nodes, curr->theta, in_nodes, 1.0, next->act, out_nodes);

		/* apply sigmoid function */
		for (int i = 0; i < out_nodes * batch_size; i++) {
			next->act[i] = 1.0 / (1 + expf(-next->act[i]));
		}

		/* transition to next layer */
		curr = next;
	}
}

void
neural_sgd_feed_forward(struct neural_net_layer *head, float *x)
{
	/* point "current" layer to head and set input as feature array x */
	struct neural_net_layer *curr = head;
	struct neural_net_layer *next;
	curr->act = x;

	int rows;
	int cols;

	while ((next = curr->next) != NULL) {
		rows = curr->out_nodes;
		cols = curr->in_nodes;

		/* act = theta * x + act */
		memcpy(next->act, curr->bias, rows * sizeof(float));
		cblas_sgemv(CblasRowMajor, CblasNoTrans, rows, cols, 1.0, curr->theta,
			cols, curr->act, 1, 1.0, next->act, 1);

		/* apply sigmoid function */
		for (int i = 0; i < rows; i++) {
			next->act[i] = 1.0 / (1 + expf(-next->act[i]));
		}

		/* transition to next layer */
		curr = next;
	}
}

void
neural_back_prop(struct neural_net_layer *tail, float *y, const int batch_size,
		const float alpha, const float lambda)
{
	/* get delta from final predictions vs. actual; delta = a - y */
	cblas_scopy(tail->in_nodes * batch_size, tail->act, 1, tail->delta, 1);
	cblas_saxpy(tail->in_nodes * batch_size, -1.0, y, 1, tail->delta, 1);

	struct neural_net_layer *curr;
	struct neural_net_layer *prev;
	int in_nodes;
	int out_nodes;

	for (curr = tail; curr->prev != NULL; curr = prev) {
		prev = curr->prev;
		in_nodes = prev->in_nodes;
		out_nodes = prev->out_nodes;

		/* Delta^(n-1) = Delta^(n) * Theta^(n-1) .* A^(n-1) (1 - A^(n-1)) */
		cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, batch_size, in_nodes, out_nodes,
				1.0, curr->delta, out_nodes, prev->theta, in_nodes, 0.0, prev->delta, in_nodes);
		for (int j = 0; j < in_nodes * batch_size; j++) {
			prev->delta[j] *= (prev->act[j] * (1 - prev->act[j]));
		}
	}

	for (curr = tail; curr->prev != NULL; curr = prev) {
		prev = curr->prev;
		in_nodes = prev->in_nodes;
		out_nodes = prev->out_nodes;

		/* bias -= alpha * delta */
		for (int i = 0; i < batch_size; i++) {
			cblas_saxpy(out_nodes, -alpha, curr->delta + i * out_nodes, 1, prev->bias, 1);
		}

		/* Theta -= alpha * (Delta^T * Act + lambda * Theta) */
		cblas_sgemm(CblasRowMajor, CblasTrans, CblasNoTrans, out_nodes, in_nodes, batch_size,
				-alpha, curr->delta, out_nodes, prev->act, in_nodes, 1 - alpha * lambda,
				prev->theta, in_nodes);
	}
}

void
neural_sgd_back_prop(struct neural_net_layer *tail, float *y, const float alpha,
		const float lambda)
{
	/* get delta from final predictions vs. actual; delta = a - y */
	cblas_scopy(tail->in_nodes, tail->act, 1, tail->delta, 1);
	cblas_saxpy(tail->in_nodes, -1.0, y, 1, tail->delta, 1);

	struct neural_net_layer *curr;
	struct neural_net_layer *prev;
	int rows;
	int cols;

	for (curr = tail; curr->prev != NULL; curr = prev) {
		prev = curr->prev;
		rows = prev->out_nodes;
		cols = prev->in_nodes;

		/* delta^(n-1) = theta^(n-1)T * delta^(n) .* a^(n-1) (1 - a^(n-1)) */
		cblas_sgemv(CblasRowMajor, CblasTrans, rows, cols, 1.0, prev->theta,
				cols, curr->delta, 1, 0.0, prev->delta, 1);
		for (int j = 0; j < cols; j++) {
			prev->delta[j] *= (prev->act[j] * (1 - prev->act[j]));
		}
	}

	for (curr = tail; curr->prev != NULL; curr = prev) {
		prev = curr->prev;
		rows = prev->out_nodes;
		cols = prev->in_nodes;

		/* bias -= alpha * delta */
		cblas_saxpy(rows, -alpha, curr->delta, 1, prev->bias, 1);

		/* theta -= alpha * (delta * actT + lambda * theta) */
		cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans, rows, cols, 1, -alpha,
				curr->delta, 1, prev->act, cols, 1 - alpha * lambda, prev->theta, cols);
	}
}

void
neural_train_iteration(struct neural_net_layer *head, struct neural_net_layer *tail,
		float *features, float *labels, const int n_samples, const int batch_size,
		const float alpha, const float lambda)
{
	int n_features = head->in_nodes;
	int n_labels = tail->in_nodes;
	int i;

	if (batch_size == 1) {
		for (i = 0; i < n_samples; i++) {
			neural_sgd_feed_forward(head, features + i * n_features);
			neural_sgd_back_prop(tail, labels + i * n_labels, alpha, lambda);
		}
	} else {
		for (i = 0; i < n_samples / batch_size; i++) {
			neural_feed_forward(head, features + i * n_features * batch_size, batch_size);
			neural_back_prop(tail, labels + i * n_labels * batch_size, batch_size, alpha, lambda);
		}

		/* take care of remaining examples */
		int remainder = n_samples % batch_size;
		if (remainder > 0) {
			neural_feed_forward(head, features + i * n_features * batch_size, remainder);
			neural_back_prop(tail, labels + i * n_labels * batch_size, remainder, alpha, lambda);
		}
	}
}

void
neural_predict_prob(struct neural_net_layer *head, struct neural_net_layer *tail,
		float *features, float *preds, const int n_samples, const int batch_size)
{
	int n_features = head->in_nodes;
	int n_labels = tail->in_nodes;
	int i;

	if (batch_size == 1) {
		for (i = 0; i < n_samples; i++) {
			neural_sgd_feed_forward(head, features + i * n_features);
			memcpy(preds + i * n_labels, tail->act, n_labels * sizeof(float));
		}
	} else {
		for (i = 0; i < n_samples / batch_size; i++) {
			neural_feed_forward(head, features + i * n_features * batch_size, batch_size);
			memcpy(preds + i * n_labels * batch_size, tail->act,
					n_labels * batch_size * sizeof(float));
		}

		/* take care of remaining examples */
		int remainder = n_samples % batch_size;
		if (remainder > 0) {
			neural_feed_forward(head, features + i * n_features * batch_size, remainder);
			memcpy(preds + i * n_labels * batch_size, tail->act,
					n_labels * remainder * sizeof(float));
		}
	}
}

max_stars_repo_name: jwilk-forks/pyneural

max_issues_repo_name: jwilk-forks/pyneural

max_forks_repo_name: jwilk-forks/pyneural


text:
#ifndef BATCH_EXTRACTOR_H_IVH3CLLQ
#define BATCH_EXTRACTOR_H_IVH3CLLQ

#include <gsl/gsl>
#include <optional>
#include <sens_loc/util/correctness_util.h>
#include <string_view>

namespace sens_loc::apps {

/// \ingroup feature-plotter-driver
/// Defines a strongly typed switch for the color to use. Will be mapped
/// to the actual OpenCV-color using \c color_to_bgr.
enum class feature_color {
    green,
    blue,
    red,
    orange,
    purple,  ///< Nice purple.
    all      ///< Draws all keypoints with different colors.
};

/// Convert a color name to the corresponding \c feature_color.
/// The string must match with the enum-value, e.g. "green".
/// \ingroup feature-plotter-driver
inline feature_color str_to_color(std::string_view color_string) noexcept {

#define COLOR_SWITCH(COLOR)                                                    \
    if (color_string == #COLOR)                                                \
        return feature_color::COLOR;
    COLOR_SWITCH(all)
    COLOR_SWITCH(green)
    COLOR_SWITCH(blue)
    COLOR_SWITCH(red)
    COLOR_SWITCH(purple)
    COLOR_SWITCH(orange)

    UNREACHABLE("Unexpected color to convert");  // LCOV_EXCL_LINE
}

/// Provide a conversion from the \¢ feature_color to actual OpenCV colors.
/// \ingroup feature-plotter-driver
struct color_to_bgr {
    // Color Space in BGR.
    static cv::Scalar convert(feature_color c) {
        using cv::Scalar;
        switch (c) {
        case feature_color::all: return Scalar::all(-1);
        case feature_color::green: return Scalar(5, 117, 65);
        case feature_color::blue: return Scalar(226, 144, 74);
        case feature_color::red: return Scalar(27, 2, 208);
        case feature_color::purple: return Scalar(254, 19, 144);
        case feature_color::orange: return Scalar(35, 166, 245);
        }
        UNREACHABLE("Invalid enum-value!");  // LCOV_EXCL_LINE
    }
};

/// Helper class that visits a list of images plots keypoints onto them.
/// \ingroup feature-plotter-driver
class batch_plotter {
  public:
    batch_plotter(std::string_view                feature_file_pattern,
                  std::string_view                output_file_pattern,
                  feature_color                   color,
                  std::optional<std::string_view> target_image_file_pattern)
        : _feature_file_pattern{feature_file_pattern}
        , _ouput_file_pattern{output_file_pattern}
        , _color{color}
        , _target_image_file_pattern{target_image_file_pattern} {
        Expects(!_feature_file_pattern.empty());
        Expects(!_ouput_file_pattern.empty());
        if (_target_image_file_pattern)
            Expects(!_target_image_file_pattern->empty());
    }

    /// Process a whole batch of files in the range [start, end].
    [[nodiscard]] bool process_batch(int start, int end) const noexcept;

  private:
    /// Process a single feature file. Called in parallel from \c process_batch.
    [[nodiscard]] bool process_index(int idx) const noexcept;

    /// The feature-detector creates files with the keypoints. This file-pattern
    /// needs to be provided in order to plot those keypoints.
    std::string_view _feature_file_pattern;

    /// File pattern for the output image to be written to.
    /// These images are 8-bit RGB images!
    std::string_view _ouput_file_pattern;

    /// Color to draw the keypoints in.
    feature_color _color;

    /// Even though the feature-files should have a path to the original image,
    /// the features were detected on, this path might not be the desired target
    /// for plotting. This is the case for plotting multiple feature keypoints
    /// or if the path to the file is incorrect.
    std::optional<std::string_view> _target_image_file_pattern;
};
}  // namespace sens_loc::apps

#endif /* end of include guard: BATCH_EXTRACTOR_H_IVH3CLLQ */

max_stars_repo_name: JonasToth/depth-conversions

max_issues_repo_name: JonasToth/depth-conversions

max_forks_repo_name: JonasToth/depth-conversions