|
Boost-Commit : |
Subject: [Boost-commit] svn:boost r64052 - in sandbox/odeint/branches/karsten: . boost/numeric boost/numeric/odeint/algebra boost/numeric/odeint/stepper boost/numeric/odeint/stepper/detail libs/numeric/odeint/test
From: karsten.ahnert_at_[hidden]
Date: 2010-07-15 12:58:15
Author: karsten
Date: 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
New Revision: 64052
URL: http://svn.boost.org/trac/boost/changeset/64052
Log:
* clean up and organizing test
* gsl_vector adaptor introduced
* vector space algebra refinded
* rk54_ck, rk4 implemented
* some base classes for error stepper
Added:
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/controlled_error_stepper.hpp (contents, props changed)
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_error_rk54_ck.hpp (contents, props changed)
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_rk4.hpp (contents, props changed)
sandbox/odeint/branches/karsten/libs/numeric/odeint/test/test_gsl_vector.cpp (contents, props changed)
sandbox/odeint/branches/karsten/libs/numeric/odeint/test/test_gsl_vector_with_euler.cpp (contents, props changed)
Removed:
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/runge_kutta_error_ck.hpp
Text files modified:
sandbox/odeint/branches/karsten/.cproject | 136 +++++++++++-----------
sandbox/odeint/branches/karsten/boost/numeric/odeint.hpp | 5
sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/gsl_vector_adaptor.hpp | 87 +++++++++++--
sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/standard_operations.hpp | 23 +++
sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/vector_space_algebra.hpp | 38 ++++++
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/detail/macros.hpp | 26 ++++
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_euler.hpp | 25 ----
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_stepper_base.hpp | 194 ++++++++++++++++++++++++++++++++
sandbox/odeint/branches/karsten/libs/numeric/odeint/test/Jamfile | 20 ++
sandbox/odeint/branches/karsten/libs/numeric/odeint/test/check_resize.cpp | 63 +++++----
sandbox/odeint/branches/karsten/libs/numeric/odeint/test/check_stepper_concepts.cpp | 236 ++++++++++++++++++++++++++++++----------
11 files changed, 645 insertions(+), 208 deletions(-)
Modified: sandbox/odeint/branches/karsten/.cproject
==============================================================================
--- sandbox/odeint/branches/karsten/.cproject (original)
+++ sandbox/odeint/branches/karsten/.cproject 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -5,64 +5,63 @@
<storageModule moduleId="org.eclipse.cdt.core.settings">
<cconfiguration id="cdt.managedbuild.config.gnu.exe.debug.860987696">
<storageModule buildSystemId="org.eclipse.cdt.managedbuilder.core.configurationDataProvider" id="cdt.managedbuild.config.gnu.exe.debug.860987696" moduleId="org.eclipse.cdt.core.settings" name="Debug">
- <externalSettings/>
- <extensions>
- <extension id="org.eclipse.cdt.core.ELF" point="org.eclipse.cdt.core.BinaryParser"/>
- <extension id="org.eclipse.cdt.core.GCCErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
- <extension id="org.eclipse.cdt.core.GASErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
- <extension id="org.eclipse.cdt.core.GLDErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
- <extension id="org.eclipse.cdt.core.GmakeErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
- <extension id="org.eclipse.cdt.core.CWDLocator" point="org.eclipse.cdt.core.ErrorParser"/>
- </extensions>
- </storageModule>
+<externalSettings/>
+<extensions>
+<extension id="org.eclipse.cdt.core.ELF" point="org.eclipse.cdt.core.BinaryParser"/>
+<extension id="org.eclipse.cdt.core.GCCErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
+<extension id="org.eclipse.cdt.core.GASErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
+<extension id="org.eclipse.cdt.core.GLDErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
+<extension id="org.eclipse.cdt.core.MakeErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
+</extensions>
+</storageModule>
<storageModule moduleId="cdtBuildSystem" version="4.0.0">
- <configuration artifactName="karsten" buildArtefactType="org.eclipse.cdt.build.core.buildArtefactType.exe" buildProperties="org.eclipse.cdt.build.core.buildType=org.eclipse.cdt.build.core.buildType.debug,org.eclipse.cdt.build.core.buildArtefactType=org.eclipse.cdt.build.core.buildArtefactType.exe" cleanCommand="rm -rf" description="" id="cdt.managedbuild.config.gnu.exe.debug.860987696" name="Debug" parent="cdt.managedbuild.config.gnu.exe.debug">
- <folderInfo id="cdt.managedbuild.config.gnu.exe.debug.860987696." name="/" resourcePath="">
- <toolChain id="cdt.managedbuild.toolchain.gnu.exe.debug.299425296" name="Linux GCC" superClass="cdt.managedbuild.toolchain.gnu.exe.debug">
- <targetPlatform id="cdt.managedbuild.target.gnu.platform.exe.debug.201156806" name="Debug Platform" superClass="cdt.managedbuild.target.gnu.platform.exe.debug"/>
- <builder buildPath="${workspace_loc:/karsten}" command="bjam" enabledIncrementalBuild="true" id="cdt.managedbuild.target.gnu.builder.exe.debug.1635191807" incrementalBuildTarget="debug" keepEnvironmentInBuildfile="false" managedBuildOn="false" name="Gnu Make Builder" superClass="cdt.managedbuild.target.gnu.builder.exe.debug"/>
- <tool id="cdt.managedbuild.tool.gnu.archiver.base.570769962" name="GCC Archiver" superClass="cdt.managedbuild.tool.gnu.archiver.base"/>
- <tool id="cdt.managedbuild.tool.gnu.cpp.compiler.exe.debug.1812732180" name="GCC C++ Compiler" superClass="cdt.managedbuild.tool.gnu.cpp.compiler.exe.debug">
- <option id="gnu.cpp.compiler.exe.debug.option.optimization.level.1153438008" name="Optimization Level" superClass="gnu.cpp.compiler.exe.debug.option.optimization.level" value="gnu.cpp.compiler.optimization.level.none" valueType="enumerated"/>
- <option id="gnu.cpp.compiler.exe.debug.option.debugging.level.1413060916" name="Debug Level" superClass="gnu.cpp.compiler.exe.debug.option.debugging.level" value="gnu.cpp.compiler.debugging.level.max" valueType="enumerated"/>
- <option id="gnu.cpp.compiler.option.include.paths.1220438677" name="Include paths (-I)" superClass="gnu.cpp.compiler.option.include.paths" valueType="includePath">
- <listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/ipp/em64t/include"/>
- <listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/mkl/include"/>
- <listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/tbb/include"/>
- <listOptionValue builtIn="false" value="/usr/local/cuda/include"/>
- </option>
- <inputType id="cdt.managedbuild.tool.gnu.cpp.compiler.input.1817643885" superClass="cdt.managedbuild.tool.gnu.cpp.compiler.input"/>
- </tool>
- <tool id="cdt.managedbuild.tool.gnu.c.compiler.exe.debug.394574166" name="GCC C Compiler" superClass="cdt.managedbuild.tool.gnu.c.compiler.exe.debug">
- <option defaultValue="gnu.c.optimization.level.none" id="gnu.c.compiler.exe.debug.option.optimization.level.868648107" name="Optimization Level" superClass="gnu.c.compiler.exe.debug.option.optimization.level" valueType="enumerated"/>
- <option id="gnu.c.compiler.exe.debug.option.debugging.level.589666025" name="Debug Level" superClass="gnu.c.compiler.exe.debug.option.debugging.level" value="gnu.c.debugging.level.max" valueType="enumerated"/>
- <option id="gnu.c.compiler.option.include.paths.2079736127" superClass="gnu.c.compiler.option.include.paths" valueType="includePath">
- <listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/ipp/em64t/include"/>
- <listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/mkl/include"/>
- <listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/tbb/include"/>
- <listOptionValue builtIn="false" value="/usr/local/cuda/include"/>
- </option>
- <inputType id="cdt.managedbuild.tool.gnu.c.compiler.input.1840035085" superClass="cdt.managedbuild.tool.gnu.c.compiler.input"/>
- </tool>
- <tool id="cdt.managedbuild.tool.gnu.c.linker.exe.debug.1392670668" name="GCC C Linker" superClass="cdt.managedbuild.tool.gnu.c.linker.exe.debug"/>
- <tool id="cdt.managedbuild.tool.gnu.cpp.linker.exe.debug.1920017020" name="GCC C++ Linker" superClass="cdt.managedbuild.tool.gnu.cpp.linker.exe.debug">
- <inputType id="cdt.managedbuild.tool.gnu.cpp.linker.input.2105648909" superClass="cdt.managedbuild.tool.gnu.cpp.linker.input">
- <additionalInput kind="additionalinputdependency" paths="$(USER_OBJS)"/>
- <additionalInput kind="additionalinput" paths="$(LIBS)"/>
- </inputType>
- </tool>
- <tool id="cdt.managedbuild.tool.gnu.assembler.exe.debug.1703113841" name="GCC Assembler" superClass="cdt.managedbuild.tool.gnu.assembler.exe.debug">
- <option id="gnu.both.asm.option.include.paths.1500689279" name="Include paths (-I)" superClass="gnu.both.asm.option.include.paths" valueType="includePath">
- <listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/ipp/em64t/include"/>
- <listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/mkl/include"/>
- <listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/tbb/include"/>
- </option>
- <inputType id="cdt.managedbuild.tool.gnu.assembler.input.1996174011" superClass="cdt.managedbuild.tool.gnu.assembler.input"/>
- </tool>
- </toolChain>
- </folderInfo>
- </configuration>
- </storageModule>
+<configuration artifactName="karsten" buildArtefactType="org.eclipse.cdt.build.core.buildArtefactType.exe" buildProperties="org.eclipse.cdt.build.core.buildType=org.eclipse.cdt.build.core.buildType.debug,org.eclipse.cdt.build.core.buildArtefactType=org.eclipse.cdt.build.core.buildArtefactType.exe" cleanCommand="rm -rf" description="" id="cdt.managedbuild.config.gnu.exe.debug.860987696" name="Debug" parent="cdt.managedbuild.config.gnu.exe.debug">
+<folderInfo id="cdt.managedbuild.config.gnu.exe.debug.860987696." name="/" resourcePath="">
+<toolChain id="cdt.managedbuild.toolchain.gnu.exe.debug.299425296" name="Linux GCC" superClass="cdt.managedbuild.toolchain.gnu.exe.debug">
+<targetPlatform id="cdt.managedbuild.target.gnu.platform.exe.debug.201156806" name="Debug Platform" superClass="cdt.managedbuild.target.gnu.platform.exe.debug"/>
+<builder buildPath="${workspace_loc:/karsten}" command="bjam" enabledIncrementalBuild="true" id="cdt.managedbuild.target.gnu.builder.exe.debug.1635191807" incrementalBuildTarget="debug" keepEnvironmentInBuildfile="false" managedBuildOn="false" name="Gnu Make Builder" superClass="cdt.managedbuild.target.gnu.builder.exe.debug"/>
+<tool id="cdt.managedbuild.tool.gnu.archiver.base.570769962" name="GCC Archiver" superClass="cdt.managedbuild.tool.gnu.archiver.base"/>
+<tool id="cdt.managedbuild.tool.gnu.cpp.compiler.exe.debug.1812732180" name="GCC C++ Compiler" superClass="cdt.managedbuild.tool.gnu.cpp.compiler.exe.debug">
+<option id="gnu.cpp.compiler.exe.debug.option.optimization.level.1153438008" name="Optimization Level" superClass="gnu.cpp.compiler.exe.debug.option.optimization.level" value="gnu.cpp.compiler.optimization.level.none" valueType="enumerated"/>
+<option id="gnu.cpp.compiler.exe.debug.option.debugging.level.1413060916" name="Debug Level" superClass="gnu.cpp.compiler.exe.debug.option.debugging.level" value="gnu.cpp.compiler.debugging.level.max" valueType="enumerated"/>
+<option id="gnu.cpp.compiler.option.include.paths.1220438677" name="Include paths (-I)" superClass="gnu.cpp.compiler.option.include.paths" valueType="includePath">
+<listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/ipp/em64t/include"/>
+<listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/mkl/include"/>
+<listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/tbb/include"/>
+<listOptionValue builtIn="false" value="/usr/local/cuda/include"/>
+</option>
+<inputType id="cdt.managedbuild.tool.gnu.cpp.compiler.input.1817643885" superClass="cdt.managedbuild.tool.gnu.cpp.compiler.input"/>
+</tool>
+<tool id="cdt.managedbuild.tool.gnu.c.compiler.exe.debug.394574166" name="GCC C Compiler" superClass="cdt.managedbuild.tool.gnu.c.compiler.exe.debug">
+<option defaultValue="gnu.c.optimization.level.none" id="gnu.c.compiler.exe.debug.option.optimization.level.868648107" name="Optimization Level" superClass="gnu.c.compiler.exe.debug.option.optimization.level" valueType="enumerated"/>
+<option id="gnu.c.compiler.exe.debug.option.debugging.level.589666025" name="Debug Level" superClass="gnu.c.compiler.exe.debug.option.debugging.level" value="gnu.c.debugging.level.max" valueType="enumerated"/>
+<option id="gnu.c.compiler.option.include.paths.2079736127" name="Include paths (-I)" superClass="gnu.c.compiler.option.include.paths" valueType="includePath">
+<listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/ipp/em64t/include"/>
+<listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/mkl/include"/>
+<listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/tbb/include"/>
+<listOptionValue builtIn="false" value="/usr/local/cuda/include"/>
+</option>
+<inputType id="cdt.managedbuild.tool.gnu.c.compiler.input.1840035085" superClass="cdt.managedbuild.tool.gnu.c.compiler.input"/>
+</tool>
+<tool id="cdt.managedbuild.tool.gnu.c.linker.exe.debug.1392670668" name="GCC C Linker" superClass="cdt.managedbuild.tool.gnu.c.linker.exe.debug"/>
+<tool id="cdt.managedbuild.tool.gnu.cpp.linker.exe.debug.1920017020" name="GCC C++ Linker" superClass="cdt.managedbuild.tool.gnu.cpp.linker.exe.debug">
+<inputType id="cdt.managedbuild.tool.gnu.cpp.linker.input.2105648909" superClass="cdt.managedbuild.tool.gnu.cpp.linker.input">
+<additionalInput kind="additionalinputdependency" paths="$(USER_OBJS)"/>
+<additionalInput kind="additionalinput" paths="$(LIBS)"/>
+</inputType>
+</tool>
+<tool id="cdt.managedbuild.tool.gnu.assembler.exe.debug.1703113841" name="GCC Assembler" superClass="cdt.managedbuild.tool.gnu.assembler.exe.debug">
+<option id="gnu.both.asm.option.include.paths.1500689279" name="Include paths (-I)" superClass="gnu.both.asm.option.include.paths" valueType="includePath">
+<listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/ipp/em64t/include"/>
+<listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/mkl/include"/>
+<listOptionValue builtIn="false" value="/opt/intel/Compiler/11.1/072/tbb/include"/>
+</option>
+<inputType id="cdt.managedbuild.tool.gnu.assembler.input.1996174011" superClass="cdt.managedbuild.tool.gnu.assembler.input"/>
+</tool>
+</toolChain>
+</folderInfo>
+</configuration>
+</storageModule>
<storageModule moduleId="org.eclipse.cdt.core.externalSettings"/>
<storageModule moduleId="org.eclipse.cdt.core.language.mapping"/>
<storageModule moduleId="scannerConfiguration">
@@ -484,16 +483,15 @@
</cconfiguration>
<cconfiguration id="cdt.managedbuild.config.gnu.exe.release.2045986849">
<storageModule buildSystemId="org.eclipse.cdt.managedbuilder.core.configurationDataProvider" id="cdt.managedbuild.config.gnu.exe.release.2045986849" moduleId="org.eclipse.cdt.core.settings" name="Release">
- <externalSettings/>
- <extensions>
- <extension id="org.eclipse.cdt.core.ELF" point="org.eclipse.cdt.core.BinaryParser"/>
- <extension id="org.eclipse.cdt.core.GCCErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
- <extension id="org.eclipse.cdt.core.GASErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
- <extension id="org.eclipse.cdt.core.GLDErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
- <extension id="org.eclipse.cdt.core.GmakeErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
- <extension id="org.eclipse.cdt.core.CWDLocator" point="org.eclipse.cdt.core.ErrorParser"/>
- </extensions>
- </storageModule>
+<externalSettings/>
+<extensions>
+<extension id="org.eclipse.cdt.core.ELF" point="org.eclipse.cdt.core.BinaryParser"/>
+<extension id="org.eclipse.cdt.core.GCCErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
+<extension id="org.eclipse.cdt.core.GASErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
+<extension id="org.eclipse.cdt.core.GLDErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
+<extension id="org.eclipse.cdt.core.MakeErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
+</extensions>
+</storageModule>
<storageModule moduleId="cdtBuildSystem" version="4.0.0">
<configuration artifactName="karsten" buildArtefactType="org.eclipse.cdt.build.core.buildArtefactType.exe" buildProperties="org.eclipse.cdt.build.core.buildType=org.eclipse.cdt.build.core.buildType.release,org.eclipse.cdt.build.core.buildArtefactType=org.eclipse.cdt.build.core.buildArtefactType.exe" cleanCommand="rm -rf" description="" id="cdt.managedbuild.config.gnu.exe.release.2045986849" name="Release" parent="cdt.managedbuild.config.gnu.exe.release">
<folderInfo id="cdt.managedbuild.config.gnu.exe.release.2045986849." name="/" resourcePath="">
@@ -947,6 +945,6 @@
</cconfiguration>
</storageModule>
<storageModule moduleId="cdtBuildSystem" version="4.0.0">
- <project id="karsten.cdt.managedbuild.target.gnu.exe.2055294126" name="Executable" projectType="cdt.managedbuild.target.gnu.exe"/>
- </storageModule>
+<project id="karsten.cdt.managedbuild.target.gnu.exe.2055294126" name="Executable" projectType="cdt.managedbuild.target.gnu.exe"/>
+</storageModule>
</cproject>
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -18,7 +18,10 @@
#include <boost/config.hpp>
#include <boost/numeric/odeint/stepper/explicit_euler.hpp>
+#include <boost/numeric/odeint/stepper/explicit_rk4.hpp>
-// #include <boost/numeric/odeint/stepper/runge_kutta_error_ck.hpp>
+#include <boost/numeric/odeint/stepper/explicit_error_rk54_ck.hpp>
+
+#include <boost/numeric/odeint/stepper/controlled_error_stepper.hpp>
#endif // BOOST_NUMERIC_ODEINT_HPP
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/gsl_vector_adaptor.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/gsl_vector_adaptor.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/gsl_vector_adaptor.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -13,14 +13,17 @@
#ifndef BOOST_NUMERIC_ODEINT_GSL_VECTOR_ADAPTOR_HPP_INCLUDED
#define BOOST_NUMERIC_ODEINT_GSL_VECTOR_ADAPTOR_HPP_INCLUDED
+#include <new>
+#include <iostream>
+
#include <gsl/gsl_vector.h>
-#include <boost/range.hpp>
+#include <boost/range.hpp>
#include <boost/iterator/iterator_facade.hpp>
-namespace boost {
-namespace numeric {
-namespace odeint {
+using namespace std;
+
+
class const_gsl_vector_iterator;
@@ -57,7 +60,7 @@
/*
* defines an const iterator for gsl_vector
*/
-class const_gsl_vector_iterator : public boost::iterator_facade< const_gsl_vector_iterator , double , boost::random_access_traversal_tag >
+class const_gsl_vector_iterator : public boost::iterator_facade< const_gsl_vector_iterator , const double , boost::random_access_traversal_tag >
{
public :
@@ -101,9 +104,8 @@
}
-} // namespace odeint
-} // namespace numeric
-} // namespace boost
+
+
namespace boost
@@ -111,51 +113,102 @@
template<>
struct range_mutable_iterator< gsl_vector >
{
- typedef boost::numeric::odeint::gsl_vector_iterator type;
+ typedef gsl_vector_iterator type;
};
template<>
struct range_const_iterator< gsl_vector >
{
- typedef boost::numeric::odeint::const_gsl_vector_iterator type;
+ typedef const_gsl_vector_iterator type;
};
} // namespace boost
-namespace boost {
-namespace numeric {
-namespace odeint {
-
+// template<>
inline gsl_vector_iterator range_begin( gsl_vector &x )
{
return gsl_vector_iterator( x );
}
+// template<>
inline const_gsl_vector_iterator range_begin( const gsl_vector &x )
{
return const_gsl_vector_iterator( x );
}
+// template<>
inline gsl_vector_iterator range_end( gsl_vector &x )
{
return end_iterator( x );
}
+// template<>
inline const_gsl_vector_iterator range_end( const gsl_vector &x )
{
return end_iterator( x );
}
-} // namespace odeint
-} // namespace numeric
-} // namespace boost
+namespace boost {
+namespace numeric {
+namespace odeint {
+
+
+template<>
+struct is_resizeable< gsl_vector >
+{
+ struct type : public boost::true_type { };
+ const static bool value = type::value;
+};
+
+
+template<>
+void resize( const gsl_vector &x , gsl_vector &dxdt )
+{
+ // ToDo : the vector could be uninitialized
+ if( dxdt.owner != 0 )
+ {
+ gsl_block_free( dxdt.block );
+ }
+ dxdt.size = 0;
+
+ if( x.size == 0 ) return;
+
+ gsl_block *block = gsl_block_alloc( x.size );
+ if( block == 0 ) throw std::bad_alloc( );
+
+ dxdt.data = block->data ;
+ dxdt.size = x.size;
+ dxdt.stride = 1;
+ dxdt.block = block;
+ dxdt.owner = 1;
+}
+
+template<>
+bool same_size( const gsl_vector &x1 , const gsl_vector &x2 )
+{
+ return x1.size == x2.size;
+}
+
+template<>
+void adjust_size( const gsl_vector &x1 , gsl_vector &x2 )
+{
+ if( !same_size( x1 , x2 ) ) resize( x1 , x2 );
+}
+
+
+} // odeint
+} // numeric
+} // boost
+
+
+
#endif // BOOST_NUMERIC_ODEINT_GSL_VECTOR_ADAPTOR_HPP_INCLUDED
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/standard_operations.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/standard_operations.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/standard_operations.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -26,11 +26,11 @@
{
typedef Time time_type;
- struct increment
+ struct increment1
{
time_type m_dt;
- increment( time_type dt ) : m_dt( dt ) { }
+ increment1( time_type dt ) : m_dt( dt ) { }
template< class T1 , class T2 >
void operator()( T1 &t1 , const T2 &t2 ) const
@@ -39,6 +39,22 @@
}
};
+
+ struct increment4
+ {
+ time_type m_alpha1 , m_alpha2 , m_alpha3 , m_alpha4;
+
+ increment4( time_type alpha1 , time_type alpha2 , time_type alpha3 , time_type alpha4 )
+ : m_alpha1( alpha1 ) , m_alpha2( alpha2 ) , m_alpha3( alpha3 ) , m_alpha4( alpha4 ) { }
+
+ template< class T1 , class T2 , class T3 , class T4 , class T5 >
+ void operator()( T1 &t1 , const T2 &t2 , const T3 &t3 , const T4 &t4 , const T5 &t5 )
+ {
+ t1 += m_alpha1 * t2 + m_alpha2 * t3 + m_alpha3 * t4 + m_alpha4 * t5;
+ }
+ };
+
+
struct scale_sum2
{
time_type m_alpha1;
@@ -107,7 +123,8 @@
{
t1 = m_alpha1 * t2 + m_alpha2 * t3 + m_alpha3 * t4 + m_alpha4 * t5 + m_alpha5 * t6 + m_alpha6 * t7;
}
- };
+ };
+
};
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/vector_space_algebra.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/vector_space_algebra.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/algebra/vector_space_algebra.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -21,11 +21,47 @@
struct vector_space_algebra
{
template< class StateType1 , class StateType2 , class Operation >
- void transform2( StateType1 &s1 , StateType2 &s2 , Operation op )
+ static void for_each2( StateType1 &s1 , StateType2 &s2 , Operation op )
{
// ToDo : build checks, that the +-*/ operators are well defined
op( s1 , s2 );
}
+
+
+ template< class StateType1 , class StateType2 , class StateType3 , class Operation >
+ static void for_each3( StateType1 &s1 , StateType2 &s2 , StateType3 &s3 , Operation op )
+ {
+ op( s1 , s2 , s3 );
+ }
+
+
+ template< class StateType1 , class StateType2 , class StateType3 , class StateType4 , class Operation >
+ static void for_each4( StateType1 &s1 , StateType2 &s2 , StateType3 &s3 , StateType4 &s4 , Operation op )
+ {
+ op( s1 , s2 , s3 , s4 );
+ }
+
+
+ template< class StateType1 , class StateType2 , class StateType3 , class StateType4 , class StateType5 , class Operation >
+ static void for_each5( StateType1 &s1 , StateType2 &s2 , StateType3 &s3 , StateType4 &s4 , StateType5 &s5 , Operation op )
+ {
+ op( s1 , s2 , s3 , s4 , s5 );
+ }
+
+
+ template< class StateType1 , class StateType2 , class StateType3 , class StateType4 , class StateType5 , class StateType6 , class Operation >
+ static void for_each6( StateType1 &s1 , StateType2 &s2 , StateType3 &s3 , StateType4 &s4 , StateType5 &s5 , StateType6 &s6 , Operation op )
+ {
+ op( s1 , s2 , s3 , s4 , s5 , s6 );
+ }
+
+
+ template< class StateType1 , class StateType2 , class StateType3 , class StateType4 , class StateType5 , class StateType6 ,class StateType7 , class Operation >
+ static void for_each7( StateType1 &s1 , StateType2 &s2 , StateType3 &s3 , StateType4 &s4 , StateType5 &s5 , StateType6 &s6 , StateType7 &s7 , Operation op )
+ {
+ op( s1 , s2 , s3 , s4 , s5 , s6 , s7 );
+ }
+
};
Added: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/controlled_error_stepper.hpp
==============================================================================
--- (empty file)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/controlled_error_stepper.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -0,0 +1,58 @@
+/*
+ boost header: NUMERIC_ODEINT_STEPPER/controlled_error_stepper.hpp
+
+ Copyright 2009 Karsten Ahnert
+ Copyright 2009 Mario Mulansky
+ Copyright 2009 Andre Bergner
+
+ Distributed under the Boost Software License, Version 1.0.
+ (See accompanying file LICENSE_1_0.txt or
+ copy at http://www.boost.org/LICENSE_1_0.txt)
+*/
+
+#ifndef BOOST_NUMERIC_ODEINT_STEPPER_CONTROLLED_ERROR_STEPPER_HPP_INCLUDED
+#define BOOST_NUMERIC_ODEINT_STEPPER_CONTROLLED_ERROR_STEPPER_HPP_INCLUDED
+
+#include <boost/numeric/odeint/algebra/standard_algebra.hpp>
+#include <boost/numeric/odeint/algebra/standard_operations.hpp>
+
+#include <boost/numeric/odeint/stepper/explicit_stepper_base.hpp>
+#include <boost/numeric/odeint/stepper/detail/macros.hpp>
+
+namespace boost {
+namespace numeric {
+namespace odeint {
+
+
+//template<
+// class State ,
+// class Time = double ,
+// class Algebra = standard_algebra< State > ,
+// class Operations = standard_operations< Time > ,
+// class AdjustSizePolicy = adjust_size_initially_tag
+// >
+//class explicit_euler
+//: public explicit_stepper_base<
+// explicit_euler< State , Time , Algebra , Operations , AdjustSizePolicy > ,
+// 1 , State , Time , Algebra , Operations , AdjustSizePolicy >
+//{
+//public :
+//
+// BOOST_ODEINT_EXPLICIT_STEPPERS_TYPEDEFS( explicit_euler , 1 );
+//
+// template< class System >
+// void do_step_impl( System system , state_type &x , const state_type &dxdt , time_type t , time_type dt )
+// {
+// algebra_type::for_each2( x , dxdt , typename operations_type::increment1( dt ) );
+// }
+//};
+
+
+
+
+} // odeint
+} // numeric
+} // boost
+
+
+#endif //BOOST_NUMERIC_ODEINT_STEPPER_CONTROLLED_ERROR_STEPPER_HPP_INCLUDED
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/detail/macros.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/detail/macros.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/detail/macros.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -22,7 +22,31 @@
typedef typename stepper_base_type::time_type time_type; \
typedef typename stepper_base_type::algebra_type algebra_type; \
typedef typename stepper_base_type::operations_type operations_type; \
-typedef typename stepper_base_type::adjust_size_policy adjust_size_policy
+typedef typename stepper_base_type::adjust_size_policy adjust_size_policy; \
+typedef typename stepper_base_type::stepper_type stepper_type
+#define BOOST_ODEINT_EXPLICIT_ERROR_STEPPERS_TYPEDEFS( STEPPER , STEPPER_ORDER , ERROR_ORDER ) \
+typedef explicit_error_stepper_base< \
+STEPPER< State , Time , Algebra , Operations , AdjustSizePolicy > , \
+STEPPER_ORDER , ERROR_ORDER , State , Time , Algebra , Operations , AdjustSizePolicy > stepper_base_type; \
+typedef typename stepper_base_type::state_type state_type; \
+typedef typename stepper_base_type::time_type time_type; \
+typedef typename stepper_base_type::algebra_type algebra_type; \
+typedef typename stepper_base_type::operations_type operations_type; \
+typedef typename stepper_base_type::adjust_size_policy adjust_size_policy; \
+typedef typename stepper_base_type::stepper_type stepper_type
+
+
+#define BOOST_ODEINT_EXPLICIT_STEPPERS_AND_ERROR_STEPPERS_TYPEDEFS( STEPPER , ORDER , STEPPER_ORDER , ERROR_ORDER ) \
+typedef explicit_stepper_and_error_stepper_base< \
+STEPPER< State , Time , Algebra , Operations , AdjustSizePolicy > , \
+ORDER , STEPPER_ORDER , ERROR_ORDER , State , Time , Algebra , Operations , AdjustSizePolicy > stepper_base_type; \
+typedef typename stepper_base_type::state_type state_type; \
+typedef typename stepper_base_type::time_type time_type; \
+typedef typename stepper_base_type::algebra_type algebra_type; \
+typedef typename stepper_base_type::operations_type operations_type; \
+typedef typename stepper_base_type::adjust_size_policy adjust_size_policy; \
+typedef typename stepper_base_type::stepper_type stepper_type
+
#endif //BOOST_BOOST_NUMERIC_ODEINT_DETAIL_MACROS_HPP_INCLUDED
Added: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_error_rk54_ck.hpp
==============================================================================
--- (empty file)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_error_rk54_ck.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -0,0 +1,186 @@
+/*
+ boost header: NUMERIC_ODEINT_STEPPER/explicit_error_rk54_ck.hpp
+
+ Copyright 2009 Karsten Ahnert
+ Copyright 2009 Mario Mulansky
+ Copyright 2009 Andre Bergner
+
+ Distributed under the Boost Software License, Version 1.0.
+ (See accompanying file LICENSE_1_0.txt or
+ copy at http://www.boost.org/LICENSE_1_0.txt)
+*/
+
+#ifndef BOOST_NUMERIC_ODEINT_STEPPER_EXPLICIT_ERROR_RK54_CK_HPP_INCLUDED
+#define BOOST_NUMERIC_ODEINT_STEPPER_EXPLICIT_ERROR_RK54_CK_HPP_INCLUDED
+
+
+#include <boost/numeric/odeint/algebra/standard_algebra.hpp>
+#include <boost/numeric/odeint/algebra/standard_operations.hpp>
+#include <boost/numeric/odeint/algebra/standard_resize.hpp>
+
+#include <boost/numeric/odeint/stepper/explicit_stepper_base.hpp>
+#include <boost/numeric/odeint/stepper/detail/macros.hpp>
+
+namespace boost {
+namespace numeric {
+namespace odeint {
+
+
+
+
+
+/*
+ * ToDo: Check orders rk_ckc
+ */
+template<
+ class State ,
+ class Time = double ,
+ class Algebra = standard_algebra< State > ,
+ class Operations = standard_operations< Time > ,
+ class AdjustSizePolicy = adjust_size_initially_tag
+ >
+class explicit_error_rk54_ck
+: public explicit_stepper_and_error_stepper_base<
+ explicit_error_rk54_ck< State , Time , Algebra , Operations , AdjustSizePolicy > ,
+ 5 , 5 , 4 , State , Time , Algebra , Operations , AdjustSizePolicy >
+{
+
+public :
+
+ BOOST_ODEINT_EXPLICIT_STEPPERS_AND_ERROR_STEPPERS_TYPEDEFS( explicit_error_rk54_ck , 5 , 5 , 4);
+
+ explicit_error_rk54_ck( void ) : m_size_adjuster( *this ) , m_x1() , m_x2() , m_x3() , m_x4() , m_x5() , m_x6()
+ { }
+
+
+ template< class System >
+ void do_step_impl( System system , state_type &x , const state_type &dxdt , time_type t , time_type dt , state_type &xerr )
+ {
+
+ const time_type c1 = static_cast<time_type> ( 37.0 ) / static_cast<time_type>( 378.0 );
+ const time_type c3 = static_cast<time_type> ( 250.0 ) / static_cast<time_type>( 621.0 );
+ const time_type c4 = static_cast<time_type> ( 125.0 ) / static_cast<time_type>( 594.0 );
+ const time_type c6 = static_cast<time_type> ( 512.0 ) / static_cast<time_type>( 1771.0 );
+
+ const time_type dc1 = c1 - static_cast<time_type> ( 2825.0 ) / static_cast<time_type>( 27648 );
+ const time_type dc3 = c3 - static_cast<time_type> ( 18575.0 ) / static_cast<time_type>( 48384.0 );
+ const time_type dc4 = c4 - static_cast<time_type> ( 13525.0 ) / static_cast<time_type>( 55296.0 );
+ const time_type dc5 = static_cast<time_type> ( -277.0 ) / static_cast<time_type>( 14336.0 );
+ const time_type dc6 = c6 - static_cast<time_type> ( 0.25 );
+
+ do_step_impl( system , x , dxdt , t , dt );
+
+ //error estimate
+ algebra_type::for_each6( xerr , dxdt , m_x3 , m_x4 , m_x5 , m_x6 ,
+ typename operations_type::scale_sum5( dt*dc1 , dt*dc3 , dt*dc4 , dt*dc5 , dt*dc6 ));
+
+ }
+
+
+
+ template< class System >
+ void do_step_impl( System system , state_type &x , const state_type &dxdt , time_type t , time_type dt )
+ {
+ /* ToDo: separate resize m_dxdt and m_x1..6 */
+ m_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() );
+
+ const time_type a2 = static_cast<time_type> ( 0.2 );
+ const time_type a3 = static_cast<time_type> ( 0.3 );
+ const time_type a4 = static_cast<time_type> ( 0.6 );
+ const time_type a5 = static_cast<time_type> ( 1.0 );
+ const time_type a6 = static_cast<time_type> ( 0.875 );
+
+ const time_type b21 = static_cast<time_type> ( 0.2 );
+ const time_type b31 = static_cast<time_type> ( 3.0 ) / static_cast<time_type>( 40.0 );
+ const time_type b32 = static_cast<time_type> ( 9.0 ) / static_cast<time_type>( 40.0 );
+ const time_type b41 = static_cast<time_type> ( 0.3 );
+ const time_type b42 = static_cast<time_type> ( -0.9 );
+ const time_type b43 = static_cast<time_type> ( 1.2 );
+ const time_type b51 = static_cast<time_type> ( -11.0 ) / static_cast<time_type>( 54.0 );
+ const time_type b52 = static_cast<time_type> ( 2.5 );
+ const time_type b53 = static_cast<time_type> ( -70.0 ) / static_cast<time_type>( 27.0 );
+ const time_type b54 = static_cast<time_type> ( 35.0 ) / static_cast<time_type>( 27.0 );
+ const time_type b61 = static_cast<time_type> ( 1631.0 ) / static_cast<time_type>( 55296.0 );
+ const time_type b62 = static_cast<time_type> ( 175.0 ) / static_cast<time_type>( 512.0 );
+ const time_type b63 = static_cast<time_type> ( 575.0 ) / static_cast<time_type>( 13824.0 );
+ const time_type b64 = static_cast<time_type> ( 44275.0 ) / static_cast<time_type>( 110592.0 );
+ const time_type b65 = static_cast<time_type> ( 253.0 ) / static_cast<time_type>( 4096.0 );
+
+ const time_type c1 = static_cast<time_type> ( 37.0 ) / static_cast<time_type>( 378.0 );
+ const time_type c3 = static_cast<time_type> ( 250.0 ) / static_cast<time_type>( 621.0 );
+ const time_type c4 = static_cast<time_type> ( 125.0 ) / static_cast<time_type>( 594.0 );
+ const time_type c6 = static_cast<time_type> ( 512.0 ) / static_cast<time_type>( 1771.0 );
+
+ //m_x1 = x + dt*b21*dxdt
+ algebra_type::for_each3( m_x1 , x , dxdt ,
+ typename operations_type::scale_sum2( 1.0 , dt*b21 ) );
+
+ system( m_x1 , m_x2 , t + dt*a2 );
+ // m_x1 = x + dt*b31*dxdt + dt*b32*m_x2
+ algebra_type::for_each4( m_x1 , x , dxdt , m_x2 ,
+ typename operations_type::scale_sum3( 1.0 , dt*b31 , dt*b32 ));
+
+ system( m_x1 , m_x3 , t + dt*a3 );
+ // m_x1 = x + dt * (b41*dxdt + b42*m_x2 + b43*m_x3)
+ algebra_type::for_each5( m_x1 , x , dxdt , m_x2 , m_x3 ,
+ typename operations_type::scale_sum4( 1.0 , dt*b41 , dt*b42 , dt*b43 ));
+
+ system( m_x1, m_x4 , t + dt*a4 );
+ algebra_type::for_each6( m_x1 , x , dxdt , m_x2 , m_x3 , m_x4 ,
+ typename operations_type::scale_sum5( 1.0 , dt*b51 , dt*b52 , dt*b53 , dt*b54 ));
+
+ system( m_x1 , m_x5 , t + dt*a5 );
+ algebra_type::for_each7( m_x1 , x , dxdt , m_x2 , m_x3 , m_x4 , m_x5 ,
+ typename operations_type::scale_sum6( 1.0 , dt*b61 , dt*b62 , dt*b63 , dt*b64 , dt*b65 ));
+
+ /*
+ * ToDo: use increment5
+ */
+ system( m_x1 , m_x6 , t + dt*a6 );
+ algebra_type::for_each6( x , x , dxdt , m_x3 , m_x4 , m_x6 ,
+ typename operations_type::scale_sum5( 1.0 , dt*c1 , dt*c3 , dt*c4 , dt*c6 ));
+
+ }
+
+
+ void adjust_size( const state_type &x )
+ {
+ m_size_adjuster.adjust_size( x );
+ stepper_base_type::adjust_size( x );
+ }
+
+
+private:
+
+ typedef size_adjuster< state_type , stepper_type > size_adjuster_type;
+ friend class size_adjuster< state_type , stepper_type >;
+
+ void adjust_size_impl( const state_type &x )
+ {
+ boost::numeric::odeint::adjust_size( x , m_x1 );
+ boost::numeric::odeint::adjust_size( x , m_x2 );
+ boost::numeric::odeint::adjust_size( x , m_x3 );
+ boost::numeric::odeint::adjust_size( x , m_x4 );
+ boost::numeric::odeint::adjust_size( x , m_x5 );
+ boost::numeric::odeint::adjust_size( x , m_x6 );
+ }
+
+ size_adjuster_type m_size_adjuster;
+ state_type m_x1, m_x2, m_x3, m_x4, m_x5, m_x6;
+
+};
+
+
+
+
+
+
+
+} // odeint
+} // numeric
+} // boost
+
+
+
+
+#endif //BOOST_NUMERIC_ODEINT_STEPPER_EXPLICIT_ERROR_RK54_CK_HPP_INCLUDED
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_euler.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_euler.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_euler.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -38,36 +38,13 @@
{
public :
-
BOOST_ODEINT_EXPLICIT_STEPPERS_TYPEDEFS( explicit_euler , 1 );
template< class System >
void do_step_impl( System system , state_type &x , const state_type &dxdt , time_type t , time_type dt )
{
- algebra_type::for_each2( x , dxdt , typename operations_type::increment( dt ) );
+ algebra_type::for_each2( x , dxdt , typename operations_type::increment1( dt ) );
}
-
-
-// explicit_euler( void ) : m_size_adjuster( *this ) { }
-//
-// void adjust_size( const state_type &x )
-// {
-// m_size_adjuster.adjust_size( x );
-// stepper_base_type::adjust_size( x );
-// }
-//
-//
-//private:
-//
-// typedef explicit_euler< State , Time , Algebra , Operations , AdjustSizePolicy > stepper_type;
-// typedef size_adjuster< state_type , stepper_type > size_adjuster_type;
-// friend class size_adjuster< state_type , stepper_type >;
-//
-// void adjust_size_impl( const state_type &x )
-// {
-// }
-//
-// size_adjuster_type m_size_adjuster;
};
Added: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_rk4.hpp
==============================================================================
--- (empty file)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_rk4.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -0,0 +1,123 @@
+/*
+ boost header: NUMERIC_ODEINT_STEPPER/explicit_rk4.hpp
+
+ Copyright 2009 Karsten Ahnert
+ Copyright 2009 Mario Mulansky
+ Copyright 2009 Andre Bergner
+
+ Distributed under the Boost Software License, Version 1.0.
+ (See accompanying file LICENSE_1_0.txt or
+ copy at http://www.boost.org/LICENSE_1_0.txt)
+*/
+
+#ifndef BOOST_NUMERIC_ODEINT_STEPPER_EXPLICIT_RK4_HPP_INCLUDED
+#define BOOST_NUMERIC_ODEINT_STEPPER_EXPLICIT_RK4_HPP_INCLUDED
+
+
+#include <boost/numeric/odeint/algebra/standard_algebra.hpp>
+#include <boost/numeric/odeint/algebra/standard_operations.hpp>
+#include <boost/numeric/odeint/algebra/standard_resize.hpp>
+
+#include <boost/numeric/odeint/stepper/explicit_stepper_base.hpp>
+#include <boost/numeric/odeint/stepper/detail/macros.hpp>
+
+namespace boost {
+namespace numeric {
+namespace odeint {
+
+
+template<
+ class State ,
+ class Time = double ,
+ class Algebra = standard_algebra< State > ,
+ class Operations = standard_operations< Time > ,
+ class AdjustSizePolicy = adjust_size_initially_tag
+ >
+class explicit_rk4
+: public explicit_stepper_base<
+ explicit_rk4< State , Time , Algebra , Operations , AdjustSizePolicy > ,
+ 4 , State , Time , Algebra , Operations , AdjustSizePolicy >
+{
+public :
+
+
+ BOOST_ODEINT_EXPLICIT_STEPPERS_TYPEDEFS( explicit_rk4 , 1 );
+
+
+ explicit_rk4( void ) : m_size_adjuster( *this ) { }
+
+ template< class System >
+ void do_step_impl( System system , state_type &x , const state_type &dxdt , time_type t , time_type dt )
+ {
+ m_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() );
+
+ const time_type val1 = static_cast<time_type>( 1.0 );
+
+ time_type dh = static_cast<time_type>( 0.5 ) * dt;
+ time_type th = t + dh;
+
+ // dt * dxdt = k1
+ // m_xt = x + dh*dxdt
+ algebra_type::for_each3( m_xt , x , dxdt , typename operations_type::scale_sum2( val1 , dh ) );
+
+
+ // dt * m_dxt = k2
+ system( m_xt , m_dxt , th );
+
+ // m_xt = x + dh*m_dxt
+ algebra_type::for_each3( m_xt , x , m_dxt , typename operations_type::scale_sum2( val1 , dh ) );
+
+
+ // dt * m_dxm = k3
+ system( m_xt , m_dxm , th );
+ //m_xt = x + dt*m_dxm
+ algebra_type::for_each3( m_xt , x , m_dxm , typename operations_type::scale_sum2( val1 , dt ) );
+
+
+ // dt * m_dxh = k4
+ system( m_xt , m_dxh , t + dt );
+ //x += dt/6 * ( m_dxdt + m_dxt + val2*m_dxm )
+ time_type dt6 = dt / static_cast<time_type>( 6.0 );
+ time_type dt3 = dt / static_cast<time_type>( 3.0 );
+ algebra_type::for_each5( x , dxdt , m_dxt , m_dxm , m_dxh , typename operations_type::increment4( dt6 , dt3 , dt3 , dt6 ) );
+ }
+
+
+ void adjust_size( const state_type &x )
+ {
+ m_size_adjuster.adjust_size( x );
+ stepper_base_type::adjust_size( x );
+ }
+
+
+private:
+
+ typedef size_adjuster< state_type , stepper_type > size_adjuster_type;
+ friend class size_adjuster< state_type , stepper_type >;
+
+ void adjust_size_impl( const state_type &x )
+ {
+ boost::numeric::odeint::adjust_size( x , m_dxt );
+ boost::numeric::odeint::adjust_size( x , m_dxm );
+ boost::numeric::odeint::adjust_size( x , m_dxh );
+ boost::numeric::odeint::adjust_size( x , m_xt );
+ }
+
+ size_adjuster_type m_size_adjuster;
+
+ state_type m_dxt;
+ state_type m_dxm;
+ state_type m_dxh;
+ state_type m_xt;
+
+};
+
+
+
+
+} // odeint
+} // numeric
+} // boost
+
+
+#endif //BOOST_NUMERIC_ODEINT_STEPPER_EXPLICIT_RK4_HPP_INCLUDED
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_stepper_base.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_stepper_base.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/explicit_stepper_base.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -24,6 +24,8 @@
+
+
/*
* Tags to specify resize behavior of steppers
*/
@@ -68,7 +70,6 @@
}
-
private:
@@ -92,8 +93,13 @@
+
+
+
+
/*
* base class for explicit steppers
+ * models the stepper concept
*/
template<
class Stepper ,
@@ -164,6 +170,192 @@
+
+
+
+
+
+
+
+/*
+ * base class for explicit error steppers
+ * models the error stepper concept
+ * ToDo : test
+ */
+template<
+ class ErrorStepper ,
+ unsigned short StepperOrder ,
+ unsigned short ErrorOrder ,
+ class State ,
+ class Time ,
+ class Algebra ,
+ class Operations ,
+ class AdjustSizePolicy
+>
+class explicit_error_stepper_base
+{
+public:
+
+
+ typedef State state_type;
+ typedef Time time_type;
+ typedef Algebra algebra_type;
+ typedef Operations operations_type;
+ typedef AdjustSizePolicy adjust_size_policy;
+ typedef ErrorStepper stepper_type;
+
+ typedef unsigned short order_type;
+ static const order_type stepper_order_value = StepperOrder;
+ static const order_type error_order_value = ErrorOrder;
+
+ order_type stepper_order( void ) const { return stepper_order_value; }
+
+ order_type error_order( void ) const { return error_order_value; }
+
+
+ explicit_error_stepper_base( void ) : m_size_adjuster( *this ) { }
+
+
+ stepper_type& stepper( void ) { return *static_cast< stepper_type* >( this ); }
+ const stepper_type& stepper( void ) const {return *static_cast< const stepper_type* >( this );}
+
+
+
+ template< class System >
+ void do_step( System system , state_type &x , time_type t , time_type dt , state_type &xerr )
+ {
+ m_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() );
+ system( x , m_dxdt ,t );
+ this->stepper().do_step_impl( system , x , m_dxdt , t , dt , xerr );
+ }
+
+
+ void adjust_size( const state_type &x )
+ {
+ m_size_adjuster.adjust_size( x );
+ }
+
+
+private:
+
+ typedef explicit_error_stepper_base< ErrorStepper , StepperOrder , ErrorOrder , State , Time , Algebra , Operations , AdjustSizePolicy > internal_stepper_base_type;
+ typedef size_adjuster< state_type , internal_stepper_base_type > base_size_adjuster_type;
+ friend class size_adjuster< state_type , internal_stepper_base_type >;
+
+ void adjust_size_impl( const state_type &x )
+ {
+ boost::numeric::odeint::adjust_size( x , m_dxdt );
+ }
+
+ state_type m_dxdt;
+ base_size_adjuster_type m_size_adjuster;
+};
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+/*
+ * base class for explicit stepper and error steppers
+ * models the stepper AND the error stepper concept
+ */
+template<
+ class Stepper ,
+ unsigned short Order ,
+ unsigned short StepperOrder ,
+ unsigned short ErrorOrder ,
+ class State ,
+ class Time ,
+ class Algebra ,
+ class Operations ,
+ class AdjustSizePolicy
+>
+class explicit_stepper_and_error_stepper_base
+{
+public:
+
+
+ typedef State state_type;
+ typedef Time time_type;
+ typedef Algebra algebra_type;
+ typedef Operations operations_type;
+ typedef AdjustSizePolicy adjust_size_policy;
+ typedef Stepper stepper_type;
+
+ typedef unsigned short order_type;
+ static const order_type order_value = Order;
+ static const order_type stepper_order_value = StepperOrder;
+ static const order_type error_order_value = ErrorOrder;
+
+ order_type order( void ) const { return order_value; }
+ order_type stepper_order( void ) const { return stepper_order_value; }
+ order_type error_order( void ) const { return error_order_value; }
+
+
+ explicit_stepper_and_error_stepper_base( void ) : m_size_adjuster( *this ) { }
+
+
+ stepper_type& stepper( void ) { return *static_cast< stepper_type* >( this ); }
+ const stepper_type& stepper( void ) const {return *static_cast< const stepper_type* >( this );}
+
+
+ template< class System >
+ void do_step( System system , state_type &x , time_type t , time_type dt )
+ {
+ m_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() );
+ system( x , m_dxdt ,t );
+ this->stepper().do_step_impl( system , x , m_dxdt , t , dt );
+ }
+
+
+ template< class System >
+ void do_step( System system , state_type &x , time_type t , time_type dt , state_type &xerr )
+ {
+ m_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() );
+ system( x , m_dxdt ,t );
+ this->stepper().do_step_impl( system , x , m_dxdt , t , dt , xerr );
+ }
+
+
+
+
+ void adjust_size( const state_type &x )
+ {
+ m_size_adjuster.adjust_size( x );
+ }
+
+
+private:
+
+ typedef explicit_stepper_and_error_stepper_base< Stepper , Order , StepperOrder , ErrorOrder , State , Time , Algebra , Operations , AdjustSizePolicy > internal_stepper_base_type;
+ typedef size_adjuster< state_type , internal_stepper_base_type > base_size_adjuster_type;
+ friend class size_adjuster< state_type , internal_stepper_base_type >;
+
+ void adjust_size_impl( const state_type &x )
+ {
+ boost::numeric::odeint::adjust_size( x , m_dxdt );
+ }
+
+ state_type m_dxdt;
+ base_size_adjuster_type m_size_adjuster;
+};
+
+
+
+
} // odeint
} // numeric
} // boost
Deleted: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/runge_kutta_error_ck.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/runge_kutta_error_ck.hpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
+++ (empty file)
@@ -1,159 +0,0 @@
-/*
- boost header: BOOST_NUMERIC_ODEINT/runge_kutta_error.hpp
-
- Copyright 2009 Karsten Ahnert
- Copyright 2009 Mario Mulansky
-
- Distributed under the Boost Software License, Version 1.0.
- (See accompanying file LICENSE_1_0.txt or
- copy at http://www.boost.org/LICENSE_1_0.txt)
-*/
-
-#ifndef BOOST_BOOST_NUMERIC_ODEINT_RUNGE_KUTTA_ERROR_HPP_INCLUDED
-#define BOOST_BOOST_NUMERIC_ODEINT_RUNGE_KUTTA_ERROR_HPP_INCLUDED
-
-#include <boost/numeric/odeint/algebra/standard_algebra.hpp>
-#include <boost/numeric/odeint/algebra/standard_operations.hpp>
-#include <boost/numeric/odeint/algebra/standard_resize.hpp>
-
-#include <boost/numeric/odeint/stepper/explicit_stepper_base.hpp>
-#include <boost/numeric/odeint/stepper/detail/macros.hpp>
-
-namespace boost {
-namespace numeric {
-namespace odeint {
-
-/*
- * ToDO: write error_stepper_base
- */
-
-template<
- class State ,
- class Time = double ,
- class Algebra = standard_algebra< State > ,
- class Operations = standard_operations< Time > ,
- class AdjustSizePolicy = adjust_size_initially_tag
- >
-class runge_kutta_error_ck
-: public explicit_stepper_base<
- runge_kutta_error_ck< State , Time , Algebra , Operations , AdjustSizePolicy > ,
- 5 , State , Time , Algebra , Operations , AdjustSizePolicy >
-{
-
-public :
-
- BOOST_ODEINT_EXPLICIT_STEPPERS_TYPEDEFS( runge_kutta_error_ck , 5 );
-
- friend class explicit_stepper_base< runge_kutta_error_ck< State , Time , Algebra , Operations , AdjustSizePolicy > ,
- 5 , State , Time , Algebra , Operations , AdjustSizePolicy >;
-
- runge_kutta_error_ck( void ) : m_dxdt() , m_x1() , m_x2() , m_x3() , m_x4() , m_x5() , m_x6()
- { }
-
-
- template< class System >
- void do_step( System system , state_type &x , time_type t , time_type dt , state_type &xerr)
- {
- this->adjust_size_by_policy( x , adjust_size_policy() );
- system( x , m_dxdt ,t );
- do_step( system , x , m_dxdt , t , dt , xerr);
- }
-
-
- template< class System >
- void do_step( System system , state_type &x , const state_type &dxdt , time_type t , time_type dt , state_type &xerr )
- {
- /* ToDo: separate resize m_dxdt and m_x1..6 */
- this->adjust_size_by_policy( x , adjust_size_policy() );
-
- const time_type a2 = static_cast<time_type> ( 0.2 );
- const time_type a3 = static_cast<time_type> ( 0.3 );
- const time_type a4 = static_cast<time_type> ( 0.6 );
- const time_type a5 = static_cast<time_type> ( 1.0 );
- const time_type a6 = static_cast<time_type> ( 0.875 );
-
- const time_type b21 = static_cast<time_type> ( 0.2 );
- const time_type b31 = static_cast<time_type> ( 3.0 ) / static_cast<time_type>( 40.0 );
- const time_type b32 = static_cast<time_type> ( 9.0 ) / static_cast<time_type>( 40.0 );
- const time_type b41 = static_cast<time_type> ( 0.3 );
- const time_type b42 = static_cast<time_type> ( -0.9 );
- const time_type b43 = static_cast<time_type> ( 1.2 );
- const time_type b51 = static_cast<time_type> ( -11.0 ) / static_cast<time_type>( 54.0 );
- const time_type b52 = static_cast<time_type> ( 2.5 );
- const time_type b53 = static_cast<time_type> ( -70.0 ) / static_cast<time_type>( 27.0 );
- const time_type b54 = static_cast<time_type> ( 35.0 ) / static_cast<time_type>( 27.0 );
- const time_type b61 = static_cast<time_type> ( 1631.0 ) / static_cast<time_type>( 55296.0 );
- const time_type b62 = static_cast<time_type> ( 175.0 ) / static_cast<time_type>( 512.0 );
- const time_type b63 = static_cast<time_type> ( 575.0 ) / static_cast<time_type>( 13824.0 );
- const time_type b64 = static_cast<time_type> ( 44275.0 ) / static_cast<time_type>( 110592.0 );
- const time_type b65 = static_cast<time_type> ( 253.0 ) / static_cast<time_type>( 4096.0 );
-
- const time_type c1 = static_cast<time_type> ( 37.0 ) / static_cast<time_type>( 378.0 );
- const time_type c3 = static_cast<time_type> ( 250.0 ) / static_cast<time_type>( 621.0 );
- const time_type c4 = static_cast<time_type> ( 125.0 ) / static_cast<time_type>( 594.0 );
- const time_type c6 = static_cast<time_type> ( 512.0 ) / static_cast<time_type>( 1771.0 );
-
- const time_type dc1 = c1 - static_cast<time_type> ( 2825.0 ) / static_cast<time_type>( 27648 );
- const time_type dc3 = c3 - static_cast<time_type> ( 18575.0 ) / static_cast<time_type>( 48384.0 );
- const time_type dc4 = c4 - static_cast<time_type> ( 13525.0 ) / static_cast<time_type>( 55296.0 );
- const time_type dc5 = static_cast<time_type> ( -277.0 ) / static_cast<time_type>( 14336.0 );
- const time_type dc6 = c6 - static_cast<time_type> ( 0.25 );
-
- //m_x1 = x + dt*b21*dxdt
- algebra_type::for_each3( m_x1 , x , dxdt ,
- typename operations_type::scale_sum2( 1.0 , dt*b21 ) );
-
- system( m_x1 , m_x2 , t + dt*a2 );
- // m_x1 = x + dt*b31*dxdt + dt*b32*m_x2
- algebra_type::for_each4( m_x1 , x , dxdt , m_x2 ,
- typename operations_type::scale_sum3( 1.0 , dt*b31 , dt*b32 ));
-
- system( m_x1 , m_x3 , t + dt*a3 );
- // m_x1 = x + dt * (b41*dxdt + b42*m_x2 + b43*m_x3)
- algebra_type::for_each5( m_x1 , x , dxdt , m_x2 , m_x3 ,
- typename operations_type::scale_sum4( 1.0 , dt*b41 , dt*b42 , dt*b43 ));
-
- system( m_x1, m_x4 , t + dt*a4 );
- algebra_type::for_each6( m_x1 , x , dxdt , m_x2 , m_x3 , m_x4 ,
- typename operations_type::scale_sum5( 1.0 , dt*b51 , dt*b52 , dt*b53 , dt*b54 ));
-
- system( m_x1 , m_x5 , t + dt*a5 );
- algebra_type::for_each7( m_x1 , x , dxdt , m_x2 , m_x3 , m_x4 , m_x5 ,
- typename operations_type::scale_sum6( 1.0 , dt*b61 , dt*b62 , dt*b63 , dt*b64 , dt*b65 ));
-
- system( m_x1 , m_x6 , t + dt*a6 );
- algebra_type::for_each6( x , x , dxdt , m_x3 , m_x4 , m_x6 ,
- typename operations_type::scale_sum5( 1.0 , dt*c1 , dt*c3 , dt*c4 , dt*c6 ));
-
- //error estimate
- algebra_type::for_each6( xerr , dxdt , m_x3 , m_x4 , m_x5 , m_x6 ,
- typename operations_type::scale_sum5( dt*dc1 , dt*dc3 , dt*dc4 , dt*dc5 , dt*dc6 ));
- }
-
-
-
-private:
-
- void adjust_size_impl( const state_type &x )
- {
- adjust_size( x , m_dxdt );
- adjust_size( x , m_x1 );
- adjust_size( x , m_x2 );
- adjust_size( x , m_x3 );
- adjust_size( x , m_x4 );
- adjust_size( x , m_x5 );
- adjust_size( x , m_x6 );
- }
-
- state_type m_dxdt;
- state_type m_x1, m_x2, m_x3, m_x4, m_x5, m_x6;
-
-
-};
-
-
-} // odeint
-} // numeric
-} // boost
-
-#endif //BOOST_BOOST_NUMERIC_ODEINT_RUNGE_KUTTA_ERROR_HPP_INCLUDED
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/test/Jamfile
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/test/Jamfile (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/test/Jamfile 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -3,8 +3,17 @@
# (See accompanying file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
# bring in rules for testing
+
import testing ;
+#project
+# : requirements
+# <include>../../../..
+# <include>$(BOOST_ROOT)
+# <define>BOOST_ALL_NO_LIB=1
+# : build-dir .
+# ;
+
project
: requirements
<library>/boost/test//boost_unit_test_framework
@@ -12,11 +21,16 @@
<include>../../../..
<include>$BOOST_ROOT
;
+
lib libgsl : : <name>gsl ;
lib libgslcblas : : <name>gslcblas ;
test-suite "odeint"
- : [ run check_stepper_concepts.cpp libgsl libgslcblas ]
- [ run check_resize.cpp ]
- ;
\ No newline at end of file
+ :
+ [ run check_stepper_concepts.cpp libgsl libgslcblas ]
+ [ run check_resize.cpp ]
+ [ run test_gsl_vector_with_euler.cpp libgsl libgslcblas ]
+ ;
+
+exe test_gsl_vector : test_gsl_vector.cpp libgsl libgslcblas ;
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/test/check_resize.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/test/check_resize.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/test/check_resize.cpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -14,6 +14,8 @@
#include <cmath>
#include <boost/array.hpp>
+#include <boost/bind.hpp>
+
#include <boost/test/unit_test.hpp>
#include <boost/numeric/odeint.hpp>
@@ -57,38 +59,19 @@
dxdt[0] = 1.0;
}
-
-void test_manual_resize( void )
+template< class Stepper >
+void test_resize( Stepper &stepper , size_t multiplicity , size_t resize_calls )
{
adjust_size_count = 0;
test_array_type x;
- explicit_euler< test_array_type , double , standard_algebra< test_array_type > , standard_operations< double > , adjust_size_manually_tag > euler;
- euler.do_step( constant_system , x , 0.0 , 0.1 );
+ stepper.do_step( constant_system , x , 0.0 , 0.1 );
+ stepper.do_step( constant_system , x , 0.0 , 0.1 );
+ stepper.do_step( constant_system , x , 0.0 , 0.1 );
- BOOST_CHECK( adjust_size_count == 0 );
-}
+ BOOST_TEST_MESSAGE( "adjust_size_count : " << adjust_size_count );
-void test_initially_resize( void )
-{
- adjust_size_count = 0;
- test_array_type x;
- explicit_euler< test_array_type , double , standard_algebra< test_array_type > , standard_operations< double > , adjust_size_initially_tag > euler;
- euler.do_step( constant_system , x , 0.0 , 0.1 );
- euler.do_step( constant_system , x , 0.0 , 0.1 );
- euler.do_step( constant_system , x , 0.0 , 0.1 );
- BOOST_CHECK( adjust_size_count == 1 );
-}
-
-void test_always_resize( void )
-{
- adjust_size_count = 0;
- test_array_type x;
- explicit_euler< test_array_type , double , standard_algebra< test_array_type > , standard_operations< double > , adjust_size_always_tag > euler;
- euler.do_step( constant_system , x , 0.0 , 0.1 );
- euler.do_step( constant_system , x , 0.0 , 0.1 );
- euler.do_step( constant_system , x , 0.0 , 0.1 );
- BOOST_CHECK( adjust_size_count == 3 );
+ BOOST_CHECK_MESSAGE( adjust_size_count == resize_calls * multiplicity , "adjust_size_count : " << adjust_size_count );
}
@@ -97,9 +80,31 @@
{
test_suite *test = BOOST_TEST_SUITE( "check resize functionality" );
- test->add( BOOST_TEST_CASE( &test_manual_resize ) );
- test->add( BOOST_TEST_CASE( &test_initially_resize ) );
- test->add( BOOST_TEST_CASE( &test_always_resize ) );
+ typedef explicit_euler< test_array_type , double , standard_algebra< test_array_type > , standard_operations< double > , adjust_size_manually_tag > euler_manual_type;
+ typedef explicit_euler< test_array_type , double , standard_algebra< test_array_type > , standard_operations< double > , adjust_size_initially_tag > euler_initially_type;
+ typedef explicit_euler< test_array_type , double , standard_algebra< test_array_type > , standard_operations< double > , adjust_size_always_tag > euler_always_type;
+
+ typedef explicit_rk4< test_array_type , double , standard_algebra< test_array_type > , standard_operations< double > , adjust_size_manually_tag > rk4_manual_type;
+ typedef explicit_rk4< test_array_type , double , standard_algebra< test_array_type > , standard_operations< double > , adjust_size_initially_tag > rk4_initially_type;
+ typedef explicit_rk4< test_array_type , double , standard_algebra< test_array_type > , standard_operations< double > , adjust_size_always_tag > rk4_always_type;
+
+
+ euler_manual_type euler_manual;
+ euler_initially_type euler_initially;
+ euler_always_type euler_always;
+
+ rk4_manual_type rk4_manual;
+ rk4_initially_type rk4_initially;
+ rk4_always_type rk4_always;
+
+
+ test->add( BOOST_TEST_CASE( boost::bind( &test_resize< euler_manual_type > , euler_manual , 1 , 0 ) ) );
+ test->add( BOOST_TEST_CASE( boost::bind( &test_resize< euler_initially_type > , euler_initially , 1 , 1 ) ) );
+ test->add( BOOST_TEST_CASE( boost::bind( &test_resize< euler_always_type > , euler_always , 1 , 3 ) ) );
+
+ test->add( BOOST_TEST_CASE( boost::bind( &test_resize< rk4_manual_type > , rk4_manual , 5 , 0 ) ) );
+ test->add( BOOST_TEST_CASE( boost::bind( &test_resize< rk4_initially_type > , rk4_initially , 5 , 1 ) ) );
+ test->add( BOOST_TEST_CASE( boost::bind( &test_resize< rk4_always_type > , rk4_always , 5 , 3 ) ) );
return test;
}
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/test/check_stepper_concepts.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/test/check_stepper_concepts.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/test/check_stepper_concepts.cpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -10,12 +10,23 @@
copy at http://www.boost.org/LICENSE_1_0.txt)
*/
+#define BOOST_TEST_MODULE test_stepper_concepts
+
#include <vector>
#include <cmath>
+#include <iostream>
#include <tr1/array>
#include <boost/test/unit_test.hpp>
+#include <boost/mpl/vector.hpp>
+#include <boost/mpl/for_each.hpp>
+#include <boost/mpl/insert_range.hpp>
+#include <boost/mpl/end.hpp>
+
+#include <boost/bind.hpp>
+#include <boost/utility.hpp>
+
#include <boost/numeric/odeint.hpp>
#include <boost/numeric/odeint/algebra/vector_space_algebra.hpp>
#include <boost/numeric/odeint/algebra/gsl_vector_adaptor.hpp>
@@ -23,42 +34,31 @@
#include "vector_space_1d.hpp"
-
-
using std::vector;
using namespace boost::unit_test;
using namespace boost::numeric::odeint;
+namespace mpl = boost::mpl;
-typedef std::vector< double > state_type1;
-typedef gsl_vector state_type2;
-typedef vector_space_1d< double > state_type3;
-typedef std::tr1::array< double , 1 > state_type4;
-void constant_system1( const state_type1 &x , state_type1 &dxdt , double t )
-{
- dxdt[0] = 1.0;
-}
+typedef std::vector< double > vector_type;
+typedef gsl_vector gsl_vector_type;
+typedef vector_space_1d< double > vector_space_type;
+typedef std::tr1::array< double , 1 > array_type;
-void constant_system2( const state_type2 &x , state_type2 &dxdt , double t )
-{
- gsl_vector_set( &dxdt , 0 , 1.0 );
-}
+template< class State > struct algebra_dispatcher { typedef standard_algebra< State > type; };
+template<> struct algebra_dispatcher< vector_space_type > { typedef vector_space_algebra type; };
-void constant_system3( const state_type3 &x , state_type3 &dxdt , double t )
-{
- dxdt.m_x = 1.0;
-}
-void constant_system4( const state_type4 &x , state_type4 &dxdt , double t )
-{
- dxdt[0] = 1.0;
-}
+void constant_system1( const vector_type &x , vector_type &dxdt , double t ) { dxdt[0] = 1.0; }
+void constant_system2( const gsl_vector_type &x , gsl_vector_type &dxdt , double t ) { gsl_vector_set( &dxdt , 0 , 1.0 ); }
+void constant_system3( const vector_space_type &x , vector_space_type &dxdt , double t ) { dxdt.m_x = 1.0; }
+void constant_system4( const array_type &x , array_type &dxdt , double t ) { dxdt[0] = 1.0; }
-const double eps = 1.0e-14;
+const double eps = 1.0e-14;
template< class Stepper , class System >
void check_stepper_concept( Stepper &stepper , System system , typename Stepper::state_type &x )
@@ -69,67 +69,185 @@
typedef typename stepper_type::time_type time_type;
stepper.do_step( system , x , 0.0 , 0.1 );
- double xval = * boost::begin( x );
- BOOST_CHECK_SMALL( fabs( xval - 0.1 ) , eps );
}
template< class Stepper , class System >
-void check_error_stepper_concept( Stepper &stepper , System system ,
- typename Stepper::state_type &x , typename Stepper::state_type &xerr )
+void check_error_stepper_concept( Stepper &stepper , System system , typename Stepper::state_type &x , typename Stepper::state_type &xerr )
{
typedef Stepper stepper_type;
typedef typename stepper_type::state_type container_type;
typedef typename stepper_type::order_type order_type;
typedef typename stepper_type::time_type time_type;
- stepper.do_step( system , x , 0.0 , 0.1 , xerr);
- double xval = * boost::begin( x );
- BOOST_CHECK_SMALL( fabs( xval - 0.1 ) , eps );
+ stepper.do_step( system , x , 0.0 , 0.1 , xerr );
}
-void test_euler_with_vector( void )
+
+
+template< class Stepper , class State > struct perform_stepper_test;
+
+template< class Stepper >
+struct perform_stepper_test< Stepper , vector_type >
{
- state_type1 x( 1 , 0.0 );
- explicit_euler< state_type1 > euler;
- check_stepper_concept( euler , constant_system1 , x );
-}
+ void operator()( void )
+ {
+ vector_type x( 1 , 0.0 );
+ Stepper stepper;
+ check_stepper_concept( stepper , constant_system1 , x );
+ BOOST_CHECK_SMALL( fabs( x[0] - 0.1 ) , eps );
+ }
+};
-void test_euler_with_gsl_vector( void )
+template< class Stepper >
+struct perform_stepper_test< Stepper , gsl_vector_type >
{
- state_type2 *x = gsl_vector_alloc( 1 );
- explicit_euler< state_type2 > euler;
-// check_stepper_concept( euler , constant_system2 , *x );
- gsl_vector_free( x );
-}
+ void operator()( void ) const
+ {
+ gsl_vector_type *x = gsl_vector_alloc( 1 );
+ gsl_vector_set( x , 0 , 0.0 );
+ Stepper stepper;
+// check_stepper_concept( stepper , constant_system2 , *x );
+// BOOST_CHECK_SMALL( fabs( gsl_vector_get( x , 0 ) - 0.1 ) , eps );
+ gsl_vector_free( x );
+ }
+};
+
+template< class Stepper >
+struct perform_stepper_test< Stepper , vector_space_type >
+{
+ void operator()( void ) const
+ {
+ vector_space_type x;
+ x.m_x = 0.0;
+ Stepper stepper;
+ check_stepper_concept( stepper , constant_system3 , x );
+ BOOST_CHECK_SMALL( fabs( x.m_x - 0.1 ) , eps );
+ }
+};
-void test_euler_with_array( void )
+template< class Stepper >
+struct perform_stepper_test< Stepper , array_type >
{
- state_type4 x;
- x[0] = 0.0;
- explicit_euler< state_type4 > euler;
- check_stepper_concept( euler , constant_system4 , x );
+ void operator()( void )
+ {
+ array_type x;
+ x[0] = 0.0;
+ Stepper stepper;
+ check_stepper_concept( stepper , constant_system4 , x );
+ BOOST_CHECK_SMALL( fabs( x[0] - 0.1 ) , eps );
+ }
+};
+
+
+
+
+template< class State > class stepper_methods : public mpl::vector<
+ explicit_euler< State , double , typename algebra_dispatcher< State >::type > ,
+ explicit_rk4< State , double , typename algebra_dispatcher< State >::type > ,
+ explicit_error_rk54_ck< State , double , typename algebra_dispatcher< State >::type >
+> { };
+
+
+typedef mpl::insert_range< mpl::vector0<> , mpl::end< mpl::vector0<> >::type , stepper_methods< vector_type > >::type first_stepper_type;
+typedef mpl::insert_range< first_stepper_type , mpl::end< first_stepper_type >::type , stepper_methods< gsl_vector_type > >::type second_stepper_type;
+typedef mpl::insert_range< second_stepper_type , mpl::end< second_stepper_type >::type , stepper_methods< vector_space_type > >::type third_stepper_type;
+typedef mpl::insert_range< third_stepper_type , mpl::end< third_stepper_type >::type , stepper_methods< array_type > >::type all_stepper_methods;
+
+
+
+BOOST_AUTO_TEST_CASE_TEMPLATE( stepper_test , Stepper, all_stepper_methods )
+{
+ perform_stepper_test< Stepper , typename Stepper::state_type > tester;
+ tester();
}
-//void test_runge_kutta_error_ck_with_vector( void )
-//{
-// state_type1 x( 1 , 0.0 );
-// state_type1 xerr( 1 , 0.0 );
-// runge_kutta_error_ck< state_type1 > rk_ck;
-// check_error_stepper_concept( rk_ck , constant_system1 , x , xerr );
-//}
-test_suite* init_unit_test_suite( int argc, char* argv[] )
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+template< class Stepper , class State > struct perform_error_stepper_test;
+
+template< class Stepper >
+struct perform_error_stepper_test< Stepper , vector_type >
{
- test_suite *test = BOOST_TEST_SUITE("check stepper concepts");
+ void operator()( void )
+ {
+ vector_type x( 1 , 0.0 ) , xerr( 1 );
+ Stepper stepper;
+ check_error_stepper_concept( stepper , constant_system1 , x , xerr );
+ BOOST_CHECK_SMALL( fabs( x[0] - 0.1 ) , eps );
+ }
+};
+template< class Stepper >
+struct perform_error_stepper_test< Stepper , gsl_vector_type >
+{
+ void operator()( void ) const
+ {
+ gsl_vector_type *x = gsl_vector_alloc( 1 ) , *xerr = gsl_vector_alloc( 1 );
+ gsl_vector_set( x , 0 , 0.0 );
+ Stepper stepper;
+ // check_error_stepper_concept( stepper , constant_system2 , *x , *xerr );
+ // BOOST_CHECK_SMALL( fabs( gsl_vector_get( x , 0 ) - 0.1 ) , eps );
+ gsl_vector_free( x ); gsl_vector_free( xerr );
+ }
+};
+template< class Stepper >
+struct perform_error_stepper_test< Stepper , vector_space_type >
+{
+ void operator()( void ) const
+ {
+ vector_space_type x , xerr;
+ x.m_x = 0.0;
+ Stepper stepper;
+ check_error_stepper_concept( stepper , constant_system3 , x , xerr );
+ BOOST_CHECK_SMALL( fabs( x.m_x - 0.1 ) , eps );
+ }
+};
- test->add( BOOST_TEST_CASE( &test_euler_with_vector ) );
- test->add( BOOST_TEST_CASE( &test_euler_with_array ) );
+template< class Stepper >
+struct perform_error_stepper_test< Stepper , array_type >
+{
+ void operator()( void )
+ {
+ array_type x , xerr;
+ x[0] = 0.0;
+ Stepper stepper;
+ check_error_stepper_concept( stepper , constant_system4 , x , xerr );
+ BOOST_CHECK_SMALL( fabs( x[0] - 0.1 ) , eps );
+ }
+};
-// test->add( BOOST_TEST_CASE( &test_euler_with_gsl_vector ) );
+template< class State > class error_stepper_methods : public mpl::vector<
+ explicit_error_rk54_ck< State , double , typename algebra_dispatcher< State >::type >
+> { };
- return test;
+typedef mpl::insert_range< mpl::vector0<> , mpl::end< mpl::vector0<> >::type , error_stepper_methods< vector_type > >::type first_error_stepper_type;
+typedef mpl::insert_range< first_error_stepper_type , mpl::end< first_error_stepper_type >::type , error_stepper_methods< gsl_vector_type > >::type second_error_stepper_type;
+typedef mpl::insert_range< second_error_stepper_type , mpl::end< second_error_stepper_type >::type , error_stepper_methods< vector_space_type > >::type third_error_stepper_type;
+typedef mpl::insert_range< third_error_stepper_type , mpl::end< third_error_stepper_type >::type , error_stepper_methods< array_type > >::type all_error_stepper_methods;
+
+
+BOOST_AUTO_TEST_CASE_TEMPLATE( error_stepper_test , Stepper , all_error_stepper_methods )
+{
+ perform_error_stepper_test< Stepper , typename Stepper::state_type > tester;
+ tester();
}
Added: sandbox/odeint/branches/karsten/libs/numeric/odeint/test/test_gsl_vector.cpp
==============================================================================
--- (empty file)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/test/test_gsl_vector.cpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -0,0 +1,40 @@
+/* Boost stepper_euler.cpp test file
+
+ Copyright 2009 Karsten Ahnert
+ Copyright 2009 Mario Mulansky
+
+ This file tests the use of the euler stepper
+
+ Distributed under the Boost Software License, Version 1.0.
+ (See accompanying file LICENSE_1_0.txt or
+ copy at http://www.boost.org/LICENSE_1_0.txt)
+*/
+
+#include <vector>
+#include <cmath>
+
+#include <boost/numeric/odeint.hpp>
+#include <boost/numeric/odeint/algebra/gsl_vector_adaptor.hpp>
+
+using namespace boost::numeric::odeint;
+
+
+void constant_system2( const gsl_vector &x , gsl_vector &dxdt , double t ) { gsl_vector_set( &dxdt , 0 , 1.0 ); }
+
+const double eps = 1.0e-14;
+
+int main( int argc , char **argv )
+{
+
+ gsl_vector *x_vec = gsl_vector_alloc( 1 );
+ gsl_vector &x = *x_vec;
+
+ explicit_euler< gsl_vector > euler;
+ euler.do_step( constant_system2 , x , 0.0 , 0.1 );
+
+ gsl_vector_free( x_vec );
+}
+
+
+
+
Added: sandbox/odeint/branches/karsten/libs/numeric/odeint/test/test_gsl_vector_with_euler.cpp
==============================================================================
--- (empty file)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/test/test_gsl_vector_with_euler.cpp 2010-07-15 12:58:13 EDT (Thu, 15 Jul 2010)
@@ -0,0 +1,44 @@
+/* Boost stepper_euler.cpp test file
+
+ Copyright 2009 Karsten Ahnert
+ Copyright 2009 Mario Mulansky
+
+ This file tests the use of the euler stepper
+
+ Distributed under the Boost Software License, Version 1.0.
+ (See accompanying file LICENSE_1_0.txt or
+ copy at http://www.boost.org/LICENSE_1_0.txt)
+*/
+
+#define BOOST_TEST_MODULE test_gsl_vector_with_euler
+
+#include <vector>
+#include <cmath>
+
+#include <boost/test/unit_test.hpp>
+
+#include <boost/numeric/odeint.hpp>
+#include <boost/numeric/odeint/algebra/gsl_vector_adaptor.hpp>
+
+using namespace boost::unit_test;
+using namespace boost::numeric::odeint;
+
+
+void constant_system2( const gsl_vector &x , gsl_vector &dxdt , double t ) { gsl_vector_set( &dxdt , 0 , 1.0 ); }
+
+const double eps = 1.0e-14;
+
+BOOST_AUTO_TEST_CASE( gsl_vector_test )
+{
+ gsl_vector *x_vec = gsl_vector_alloc( 1 );
+ gsl_vector &x = *x_vec;
+
+ explicit_euler< gsl_vector > euler;
+// euler.do_step( constant_system2 , x , 0.0 , 0.1 );
+
+ gsl_vector_free( x_vec );
+}
+
+
+
+
Boost-Commit list run by bdawes at acm.org, david.abrahams at rcn.com, gregod at cs.rpi.edu, cpdaniel at pacbell.net, john at johnmaddock.co.uk