diff --git a/src/fft/projection_base.cc b/src/fft/projection_base.cc index 63ec35c..a5d4c62 100644 --- a/src/fft/projection_base.cc +++ b/src/fft/projection_base.cc @@ -1,57 +1,56 @@ /** * file projection_base.cc * * @author Till Junge * * @date 06 Dec 2017 * * @brief implementation of base class for projections * * @section LICENCE * * Copyright (C) 2017 Till Junge * * µSpectre 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, or (at * your option) any later version. * * µSpectre 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 GNU Emacs; see the file COPYING. If not, write to the * Free Software Foundation, Inc., 59 Temple Place - Suite 330, * Boston, MA 02111-1307, USA. */ #include "fft/projection_base.hh" namespace muSpectre { /* ---------------------------------------------------------------------- */ template ProjectionBase::ProjectionBase(FFT_Engine & engine) - : fft_engine{engine}, projection_container{engine.get_field_collection()}, - sizes{engine.get_sizes()} + : fft_engine{engine}, projection_container{engine.get_field_collection()} { static_assert((DimS == FFT_Engine::sdim), "spatial dimensions are incompatible"); static_assert((DimM == FFT_Engine::mdim), "material dimensions are incompatible"); } /* ---------------------------------------------------------------------- */ template void ProjectionBase:: initialise(FFT_PlanFlags flags) { fft_engine.initialise(flags); } template class ProjectionBase; template class ProjectionBase; template class ProjectionBase; } // muSpectre diff --git a/src/fft/projection_base.hh b/src/fft/projection_base.hh index 78d43cc..8fc8735 100644 --- a/src/fft/projection_base.hh +++ b/src/fft/projection_base.hh @@ -1,94 +1,93 @@ /** * file projection_base.hh * * @author Till Junge * * @date 03 Dec 2017 * * @brief Base class for Projection operators * * @section LICENCE * * Copyright (C) 2017 Till Junge * * µSpectre 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, or (at * your option) any later version. * * µSpectre 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 GNU Emacs; see the file COPYING. If not, write to the * Free Software Foundation, Inc., 59 Temple Place - Suite 330, * Boston, MA 02111-1307, USA. */ #ifndef PROJECTION_BASE_H #define PROJECTION_BASE_H #include "common/common.hh" #include "common/field_collection.hh" #include "common/field.hh" #include "fft/fft_engine_base.hh" namespace muSpectre { template struct Projection_traits { }; template class ProjectionBase { public: using FFT_Engine = FFT_Engine_base; using Ccoord = typename FFT_Engine::Ccoord; using GFieldCollection_t = typename FFT_Engine::GFieldCollection_t; using LFieldCollection_t = typename FFT_Engine::LFieldCollection_t; using Field_t = typename FFT_Engine::Field_t; using iterator = typename FFT_Engine::iterator; //! Default constructor ProjectionBase() = delete; //! Constructor with system sizes ProjectionBase(FFT_Engine & engine); //! Copy constructor ProjectionBase(const ProjectionBase &other) = delete; //! Move constructor ProjectionBase(ProjectionBase &&other) noexcept = default; //! Destructor virtual ~ProjectionBase() noexcept = default; //! Copy assignment operator ProjectionBase& operator=(const ProjectionBase &other) = delete; //! Move assignment operator ProjectionBase& operator=(ProjectionBase &&other) noexcept = default; //! initialises the fft engine (plan the transform) void initialise(FFT_PlanFlags flags = FFT_PlanFlags::estimate); //! apply the projection operator to a field void apply_projection(Field_t & field) const; protected: FFT_Engine & fft_engine; LFieldCollection_t & projection_container{}; - Ccoord sizes; private: }; } // muSpectre #endif /* PROJECTION_BASE_H */ diff --git a/src/fft/projection_finite_strain.cc b/src/fft/projection_finite_strain.cc index f6487c1..b6ba58d 100644 --- a/src/fft/projection_finite_strain.cc +++ b/src/fft/projection_finite_strain.cc @@ -1,83 +1,87 @@ /** * file projection_finite_strain.cc * * @author Till Junge * * @date 05 Dec 2017 * * @brief implementation of standard finite strain projection operator * * @section LICENCE * * Copyright (C) 2017 Till Junge * * µSpectre 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, or (at * your option) any later version. * * µSpectre 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 GNU Emacs; see the file COPYING. If not, write to the * Free Software Foundation, Inc., 59 Temple Place - Suite 330, * Boston, MA 02111-1307, USA. */ #include +#include "Eigen/Dense" #include "fft/projection_finite_strain.hh" #include "fft/fftw_engine.hh" #include "fft/fft_utils.hh" #include "common/field_map.hh" +#include "common/tensor_algebra.hh" namespace muSpectre { /* ---------------------------------------------------------------------- */ template ProjectionFiniteStrain:: ProjectionFiniteStrain(FFT_Engine & engine) :Parent{engine}, Ghat{make_field("Projection Operator", this->projection_container)} {} /* ---------------------------------------------------------------------- */ template void ProjectionFiniteStrain:: initialise(FFT_PlanFlags flags) { Parent::initialise(flags); - FFT_freqs fft_freqs(this->sizes); + FFT_freqs fft_freqs(this->fft_engine.get_sizes()); for (auto && tup: boost::combine(this->fft_engine, this->Ghat)) { const auto & ccoord = boost::get<0> (tup); auto & G = boost::get<1>(tup); auto xi = fft_freqs.get_unit_xi(ccoord); - for (Dim_t im = 0; im < DimS; ++im) { - for (Dim_t j = 0; j < DimS; ++j) { - for (Dim_t l = 0; l < DimS; ++l) { - get(G, im, j, l, im) = xi(j)*xi(l); - } - } - } + //! this is simplyfiable usinc Curnier's Méthodes numériques, 6.69(c) + G = Matrices::outer_over(Matrices::I2(), xi*xi.transpose()); + // for (Dim_t im = 0; im < DimS; ++im) { + // for (Dim_t j = 0; j < DimS; ++j) { + // for (Dim_t l = 0; l < DimS; ++l) { + // get(G, im, j, l, im) = xi(j)*xi(l); + // } + // } + // } } + this->Ghat[0].setZero(); } /* ---------------------------------------------------------------------- */ template void ProjectionFiniteStrain::apply_projection(Field_t & field) { Vector_map field_map{this->fft_engine.fft(field)}; Real factor = this->fft_engine.normalisation(); for (auto && tup: boost::combine(this->Ghat, field_map)) { auto & G{boost::get<0>(tup)}; auto & f{boost::get<1>(tup)}; f = factor * (G*f).eval(); } } template class ProjectionFiniteStrain; - template class ProjectionFiniteStrain; template class ProjectionFiniteStrain; } // muSpectre diff --git a/tests/test_fftw_engine.cc b/tests/test_fftw_engine.cc index d7b1219..95118d7 100644 --- a/tests/test_fftw_engine.cc +++ b/tests/test_fftw_engine.cc @@ -1,110 +1,110 @@ /** * file test_fftw_engine.cc * * @author Till Junge * * @date 05 Dec 2017 * * @brief tests for the fftw fft engine implementation * * @section LICENCE * * Copyright (C) 2017 Till Junge * * µSpectre 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, or (at * your option) any later version. * * µSpectre 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 GNU Emacs; see the file COPYING. If not, write to the * Free Software Foundation, Inc., 59 Temple Place - Suite 330, * Boston, MA 02111-1307, USA. */ #include -#include "boost/range/combine.hpp" +#include #include "tests.hh" #include "fft/fftw_engine.hh" #include "common/ccoord_operations.hh" #include "common/field_collection.hh" #include "common/field_map.hh" namespace muSpectre { BOOST_AUTO_TEST_SUITE(fftw_engine); /* ---------------------------------------------------------------------- */ template struct FFTW_fixture { constexpr static Dim_t box_size{3}; constexpr static Dim_t sdim{DimS}; constexpr static Dim_t mdim{DimM}; FFTW_fixture() :engine(CcoordOps::get_cube(box_size)){} FFTW_Engine engine; }; using fixlist = boost::mpl::list, FFTW_fixture, FFTW_fixture>; /* ---------------------------------------------------------------------- */ BOOST_FIXTURE_TEST_CASE_TEMPLATE(Constructor_test, Fix, fixlist, Fix) { BOOST_CHECK_NO_THROW(Fix::engine.initialise(FFT_PlanFlags::estimate)); BOOST_CHECK_EQUAL(Fix::engine.size(), ipow(Fix::box_size, Fix::sdim)); } /* ---------------------------------------------------------------------- */ BOOST_FIXTURE_TEST_CASE_TEMPLATE(fft_test, Fix, fixlist, Fix) { Fix::engine.initialise(FFT_PlanFlags::estimate); constexpr Dim_t order{2}; using FC_t = FieldCollection; FC_t fc; auto & input{make_field>("input", fc)}; auto & ref {make_field>("reference", fc)}; auto & result{make_field>("result", fc)}; fc.initialise(CcoordOps::get_cube(Fix::box_size)); using map_t = MatrixFieldMap; auto inmap{map_t{input}}; auto refmap{map_t{ref}}; auto resultmap{map_t{result}}; for (auto tup: boost::combine(inmap, refmap)) { boost::get<0>(tup).setRandom(); boost::get<1>(tup) = boost::get<0>(tup); } auto & complex_field = Fix::engine.fft(input); using cmap_t = MatrixFieldMap, Complex, Fix::mdim, Fix::mdim>; cmap_t complex_map(complex_field); Real error = complex_map[0].imag().norm(); BOOST_CHECK_LT(error, tol); /* make sure, the engine has not modified input (which is unfortunately const-casted internally, hence this test) */ for (auto && tup: boost::combine(inmap, refmap)) { Real error{(boost::get<0>(tup) - boost::get<1>(tup)).norm()}; BOOST_CHECK_LT(error, tol); } /* make sure that the ifft of fft returns the original*/ Fix::engine.ifft(result); for (auto && tup: boost::combine(resultmap, refmap)) { Real error{(boost::get<0>(tup)*Fix::engine.normalisation() - boost::get<1>(tup)).norm()}; BOOST_CHECK_LT(error, tol); if (error > tol) { std::cout << boost::get<0>(tup).array()/boost::get<1>(tup).array() << std::endl << std::endl; } } } BOOST_AUTO_TEST_SUITE_END(); } // muSpectre diff --git a/tests/test_projection_finite.cc b/tests/test_projection_finite.cc new file mode 100644 index 0000000..20b2f75 --- /dev/null +++ b/tests/test_projection_finite.cc @@ -0,0 +1,82 @@ +/** + * file test_projection_finite.cc + * + * @author Till Junge + * + * @date 07 Dec 2017 + * + * @brief tests for standard finite strain projection operator + * + * @section LICENCE + * + * Copyright (C) 2017 Till Junge + * + * µSpectre 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, or (at + * your option) any later version. + * + * µSpectre 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 GNU Emacs; see the file COPYING. If not, write to the + * Free Software Foundation, Inc., 59 Temple Place - Suite 330, + * Boston, MA 02111-1307, USA. + */ + +#include + +#include "tests.hh" +#include "fft/fftw_engine.hh" +#include "fft/projection_finite_strain.hh" +#include "common/common.hh" +#include "common/field_collection.hh" + +namespace muSpectre { + + BOOST_AUTO_TEST_SUITE(projection_finite_strain); + + + /* ---------------------------------------------------------------------- */ + template + struct Sizes { + }; + template<> + struct Sizes { + constexpr static Ccoord_t get_value() { + return Ccoord_t{3, 5};} + }; + template<> + struct Sizes { + constexpr static Ccoord_t get_value() { + return Ccoord_t{3, 5, 7};} + }; + /* ---------------------------------------------------------------------- */ + template + struct ProjectionFixture { + using Engine = FFTW_Engine; + using Parent = ProjectionFiniteStrain; + ProjectionFixture(): engine{Sizes::get_value()}, + parent(engine){} + Engine engine; + Parent parent; + }; + + /* ---------------------------------------------------------------------- */ + using fixlist = boost::mpl::list, + ProjectionFixture>; + + /* ---------------------------------------------------------------------- */ + BOOST_FIXTURE_TEST_CASE_TEMPLATE(constructor_test, fix, fixlist, fix) { + BOOST_CHECK_NO_THROW(fix::parent.initialise(FFT_PlanFlags::estimate)); + } + + + /* ---------------------------------------------------------------------- */ + BOOST_FIXTURE_TEST_CASE_TEMPLATE(Ctest_name, type_name, TL, F) + BOOST_AUTO_TEST_SUITE_END(); + +} // muSpectre