Subject: [Boost-commit] svn:boost r68573 - in sandbox/odeint/branches/karsten: . boost/numeric/odeint/external boost/numeric/odeint/stepper boost/numeric/odeint/stepper/base boost/numeric/odeint/util libs/numeric/odeint/doc libs/numeric/odeint/ideas/rosenbrock4 libs/numeric/odeint/performance libs/numeric/odeint/test
From: karsten.ahnert_at_[hidden]
Date: 2011-01-30 14:18:02
Author: karsten
Date: 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
New Revision: 68573
URL: http://svn.boost.org/trac/boost/changeset/68573
Log:
* dense output for units
* some refactoring
* testing explicit steppers with boost::range
* starting migration of rosenbrock
Added:
sandbox/odeint/branches/karsten/boost/numeric/odeint/external/
sandbox/odeint/branches/karsten/boost/numeric/odeint/util/
sandbox/odeint/branches/karsten/libs/numeric/odeint/performance/
Text files modified:
sandbox/odeint/branches/karsten/.cproject | 83 ++++++++++++++++++++-------------------
sandbox/odeint/branches/karsten/TODO | 71 ++++++++++++++++-----------------
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/base/explicit_stepper_base.hpp | 22 +++++++++-
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/controlled_error_stepper.hpp | 36 +++++++++-------
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp | 4
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_explicit.hpp | 2
sandbox/odeint/branches/karsten/libs/numeric/odeint/doc/tutorial.qbk | 6 ++
sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4.cpp | 55 ++++++++++++++-----------
sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4.hpp | 79 ++++++++++++++++++++-----------------
sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4_stepper.cpp | 20 +++++---
sandbox/odeint/branches/karsten/libs/numeric/odeint/test/Jamfile | 3
sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_ranges.cpp | 27 +++++++-----
sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_units.cpp | 38 +++++++++++++++++
13 files changed, 263 insertions(+), 183 deletions(-)
Modified: sandbox/odeint/branches/karsten/.cproject
==============================================================================
--- sandbox/odeint/branches/karsten/.cproject (original)
+++ sandbox/odeint/branches/karsten/.cproject 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -5,44 +5,45 @@
<storageModule moduleId="org.eclipse.cdt.core.settings">
<cconfiguration id="0.1427786045">
<storageModule buildSystemId="org.eclipse.cdt.managedbuilder.core.configurationDataProvider" id="0.1427786045" moduleId="org.eclipse.cdt.core.settings" name="Default">
-<externalSettings/>
-<extensions>
-<extension id="org.eclipse.cdt.core.VCErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
-<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>
+ <externalSettings/>
+ <extensions>
+ <extension id="org.eclipse.cdt.core.VCErrorParser" point="org.eclipse.cdt.core.ErrorParser"/>
+ <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>
<storageModule moduleId="cdtBuildSystem" version="4.0.0">
-<configuration buildProperties="" description="" id="0.1427786045" name="Default" parent="org.eclipse.cdt.build.core.prefbase.cfg">
-<folderInfo id="0.1427786045." name="/" resourcePath="">
-<toolChain id="org.eclipse.cdt.build.core.prefbase.toolchain.1701227309" name="No ToolChain" resourceTypeBasedDiscovery="false" superClass="org.eclipse.cdt.build.core.prefbase.toolchain">
-<targetPlatform id="org.eclipse.cdt.build.core.prefbase.toolchain.1701227309.1923743366" name=""/>
-<builder cleanBuildTarget="--clean" command="bjam" id="org.eclipse.cdt.build.core.settings.default.builder.1866451942" incrementalBuildTarget="" keepEnvironmentInBuildfile="false" managedBuildOn="false" name="Gnu Make Builder" superClass="org.eclipse.cdt.build.core.settings.default.builder"/>
-<tool id="org.eclipse.cdt.build.core.settings.holder.libs.295828857" name="holder for library settings" superClass="org.eclipse.cdt.build.core.settings.holder.libs"/>
-<tool id="org.eclipse.cdt.build.core.settings.holder.1728088817" name="Assembly" superClass="org.eclipse.cdt.build.core.settings.holder">
-<option id="org.eclipse.cdt.build.core.settings.holder.incpaths.558519945" superClass="org.eclipse.cdt.build.core.settings.holder.incpaths" valueType="includePath">
-<listOptionValue builtIn="false" value=""${BOOST_ROOT}""/>
-</option>
-<inputType id="org.eclipse.cdt.build.core.settings.holder.inType.1379843062" languageId="org.eclipse.cdt.core.assembly" languageName="Assembly" sourceContentType="org.eclipse.cdt.core.asmSource" superClass="org.eclipse.cdt.build.core.settings.holder.inType"/>
-</tool>
-<tool id="org.eclipse.cdt.build.core.settings.holder.780742962" name="GNU C++" superClass="org.eclipse.cdt.build.core.settings.holder">
-<option id="org.eclipse.cdt.build.core.settings.holder.incpaths.1009256344" superClass="org.eclipse.cdt.build.core.settings.holder.incpaths" valueType="includePath">
-<listOptionValue builtIn="false" value=""${BOOST_ROOT}""/>
-</option>
-<inputType id="org.eclipse.cdt.build.core.settings.holder.inType.157834134" languageId="org.eclipse.cdt.core.g++" languageName="GNU C++" sourceContentType="org.eclipse.cdt.core.cxxSource,org.eclipse.cdt.core.cxxHeader" superClass="org.eclipse.cdt.build.core.settings.holder.inType"/>
-</tool>
-<tool id="org.eclipse.cdt.build.core.settings.holder.421502995" name="GNU C" superClass="org.eclipse.cdt.build.core.settings.holder">
-<option id="org.eclipse.cdt.build.core.settings.holder.incpaths.388283580" superClass="org.eclipse.cdt.build.core.settings.holder.incpaths" valueType="includePath">
-<listOptionValue builtIn="false" value=""${BOOST_ROOT}""/>
-</option>
-<inputType id="org.eclipse.cdt.build.core.settings.holder.inType.1344941182" languageId="org.eclipse.cdt.core.gcc" languageName="GNU C" sourceContentType="org.eclipse.cdt.core.cSource,org.eclipse.cdt.core.cHeader" superClass="org.eclipse.cdt.build.core.settings.holder.inType"/>
-</tool>
-</toolChain>
-</folderInfo>
-</configuration>
-</storageModule>
+ <configuration buildProperties="" description="" id="0.1427786045" name="Default" parent="org.eclipse.cdt.build.core.prefbase.cfg">
+ <folderInfo id="0.1427786045." name="/" resourcePath="">
+ <toolChain id="org.eclipse.cdt.build.core.prefbase.toolchain.1701227309" name="No ToolChain" resourceTypeBasedDiscovery="false" superClass="org.eclipse.cdt.build.core.prefbase.toolchain">
+ <targetPlatform id="org.eclipse.cdt.build.core.prefbase.toolchain.1701227309.1923743366" name=""/>
+ <builder cleanBuildTarget="--clean" command="bjam" id="org.eclipse.cdt.build.core.settings.default.builder.1866451942" incrementalBuildTarget="" keepEnvironmentInBuildfile="false" managedBuildOn="false" name="Gnu Make Builder" superClass="org.eclipse.cdt.build.core.settings.default.builder"/>
+ <tool id="org.eclipse.cdt.build.core.settings.holder.libs.295828857" name="holder for library settings" superClass="org.eclipse.cdt.build.core.settings.holder.libs"/>
+ <tool id="org.eclipse.cdt.build.core.settings.holder.1728088817" name="Assembly" superClass="org.eclipse.cdt.build.core.settings.holder">
+ <option id="org.eclipse.cdt.build.core.settings.holder.incpaths.558519945" name="Include Paths" superClass="org.eclipse.cdt.build.core.settings.holder.incpaths" valueType="includePath">
+ <listOptionValue builtIn="false" value=""${BOOST_ROOT}""/>
+ </option>
+ <inputType id="org.eclipse.cdt.build.core.settings.holder.inType.1379843062" languageId="org.eclipse.cdt.core.assembly" languageName="Assembly" sourceContentType="org.eclipse.cdt.core.asmSource" superClass="org.eclipse.cdt.build.core.settings.holder.inType"/>
+ </tool>
+ <tool id="org.eclipse.cdt.build.core.settings.holder.780742962" name="GNU C++" superClass="org.eclipse.cdt.build.core.settings.holder">
+ <option id="org.eclipse.cdt.build.core.settings.holder.incpaths.1009256344" name="Include Paths" superClass="org.eclipse.cdt.build.core.settings.holder.incpaths" valueType="includePath">
+ <listOptionValue builtIn="false" value=""${BOOST_ROOT}""/>
+ </option>
+ <inputType id="org.eclipse.cdt.build.core.settings.holder.inType.157834134" languageId="org.eclipse.cdt.core.g++" languageName="GNU C++" sourceContentType="org.eclipse.cdt.core.cxxSource,org.eclipse.cdt.core.cxxHeader" superClass="org.eclipse.cdt.build.core.settings.holder.inType"/>
+ </tool>
+ <tool id="org.eclipse.cdt.build.core.settings.holder.421502995" name="GNU C" superClass="org.eclipse.cdt.build.core.settings.holder">
+ <option id="org.eclipse.cdt.build.core.settings.holder.incpaths.388283580" name="Include Paths" superClass="org.eclipse.cdt.build.core.settings.holder.incpaths" valueType="includePath">
+ <listOptionValue builtIn="false" value=""${BOOST_ROOT}""/>
+ </option>
+ <inputType id="org.eclipse.cdt.build.core.settings.holder.inType.1344941182" languageId="org.eclipse.cdt.core.gcc" languageName="GNU C" sourceContentType="org.eclipse.cdt.core.cSource,org.eclipse.cdt.core.cHeader" superClass="org.eclipse.cdt.build.core.settings.holder.inType"/>
+ </tool>
+ </toolChain>
+ </folderInfo>
+ </configuration>
+ </storageModule>
<storageModule moduleId="scannerConfiguration">
<autodiscovery enabled="true" problemReportingEnabled="true" selectedProfileId="org.eclipse.cdt.make.core.GCCStandardMakePerProjectProfile"/>
<profile id="org.eclipse.cdt.make.core.GCCStandardMakePerProjectProfile">
@@ -543,10 +544,10 @@
</storageModule>
<storageModule moduleId="org.eclipse.cdt.core.externalSettings"/>
<storageModule moduleId="org.eclipse.cdt.core.language.mapping"/>
- <storageModule moduleId="org.eclipse.cdt.internal.ui.text.commentOwnerProjectMappings"/>
-</cconfiguration>
+ <storageModule moduleId="org.eclipse.cdt.internal.ui.text.commentOwnerProjectMappings"/>
+ </cconfiguration>
</storageModule>
<storageModule moduleId="cdtBuildSystem" version="4.0.0">
-<project id="odeint.null.1351650296" name="odeint"/>
-</storageModule>
+ <project id="odeint.null.1351650296" name="odeint"/>
+ </storageModule>
</cproject>
Modified: sandbox/odeint/branches/karsten/TODO
==============================================================================
--- sandbox/odeint/branches/karsten/TODO (original)
+++ sandbox/odeint/branches/karsten/TODO 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -1,61 +1,58 @@
* Unit test
NEARLY DONE * test operations
- OK * test standard_algebra
- OK * test fusion_algebra
* test vector_space_algebra, maybe with some proto lib
NEARLY DONE * test copying
* include controlled_error_stepper_fsal
* include dense_output_explicit
* inlcude dense_output_controlled_explicit_fsal
- OK * test, if copy construct of stepper_base is called when explicit_euler is used
* test gsl
- * test explicit stepper with ranges
- * test units with dense output
+ DIFFICULT * test explicit stepper with ranges
+ * discuss
+ * split check_concepts into check_stepper_concept, check_error_stepper_concept, check_controlled_stepper_concept
+ * include test/thrust in jam system, use system from
+* include rosenbrock4 in trunk
+DIFFICULT * finishing change of controlled_stepper to units
+ * check if rosenbrock controller and controlled_stepper can both be used with the explicit steppers
+ OK * move error_checker into controlled_stepper
+ * rename controlled_stepper to a more specific name
+* file cleanup
+ * split resizing and copy/destruct/construct in different files
+ * subfolder algebra, operations, util
+* general:
+ * check if everywhere static_cast< value_type > is used
+ * check header guards
+ * check copyright note
+ * documente every file in the preamble
+ * check once more, if all contructor, destructors and assign-operators are present
+
+* Integrate functions
+* skript for setting the include defines according to the position in file system an writing a general copyright comment at the beginning
+
+* Documentation
+* Tutorials
+
+
+
+DONE:
+
OK * change standard_operations::rel_error in order to word with units and test it
OK * include implicit euler
OK * call via std::pair< deriv , jacobi >
OK * resizing
-* include rosenbrock4 in trunk
OK * operations that fit units
OK * operations that fit result_of
-* include test/thrust in jam system, use system from
OK * change stepper to stepper_units
OK * change error_stepper to error_stepper_units
-DIFFICULT * finishing change of controlled_stepper to units
- * check if rosenbrock controller and controlled_stepper can both be used with the explicit steppers
- OK * move error_checker into controlled_stepper
- * rename controlled_stepper to a more specific name
OK * change dense_output to units
OK * roll out dense_output_explicit_euler::calc_state() to explicit_euler::dense_output()
OK * roll out dense_output_dopri5::calc_state() to explicit_error_dopri5::dense_output()
OK * create dense_output_explicit
OK * create dense_output_explicit_controlled_fsal
-* split check_concepts into check_stepper_concept, check_error_stepper_concept, check_controlled_stepper_concept
-* file cleanup
- * split resizing and copy/destruct/construct in different files
- * subfolder algebra, operations, util
OK * change resizing concept, in order to word within the implicit steppers
OK * in all tests and regression test do not include odeint.hpp, only include the headers which are really needed
OK * start new doc or cleanup the old project
-* check header guards
-* check copyright note
-* documente every file in the preamble
-* check once more, if all contructor, destructors and assign-operators are present
-
-* Integrate functions
-* skript for setting the include defines according to the position in file system an writing a general copyright comment at the beginning
-
-* Documentation
-
-
-
-
-* Adaptoren
- * GMP
- * gsl_complex, gsl_complex_vector
-
-
-* Examples
- * Lyapunov exponents
- * Thrust example with a parameter study
- * A test with boost.graph, phase oscillators on a complex network
\ No newline at end of file
+OK * unit test
+ OK * test standard_algebra
+ OK * test fusion_algebra
+ OK * test, if copy construct of stepper_base is called when explicit_euler is used
+ OK * test units with dense output
\ No newline at end of file
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/base/explicit_stepper_base.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/base/explicit_stepper_base.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/base/explicit_stepper_base.hpp 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -89,19 +89,35 @@
}
-
-
-
+ /*
+ * test if boost::range works
+ *
+ * is ok, but takes many implementations of do_step() with various const
+ */
// do_step( sys , x , t , dt )
template< class System , class StateInOut >
void do_step( System system , StateInOut &x , const time_type &t , const time_type &dt )
{
+ do_step_caller( system , x , t , dt );
+ }
+
+ template< class System , class StateInOut >
+ void do_step( System system , const StateInOut &x , const time_type &t , const time_type &dt )
+ {
+ do_step_caller( system , x , t , dt );
+ }
+
+ template< class System , class StateInOut >
+ void do_step_caller( System system , StateInOut &x , const time_type &t , const time_type &dt )
+ {
typename boost::unwrap_reference< System >::type &sys = system;
m_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() );
sys( x , m_dxdt ,t );
this->stepper().do_step_impl( system , x , m_dxdt , t , x , dt );
}
+
+
// do_step( sys , x , dxdt , t , dt )
template< class System , class StateInOut , class DerivIn >
void do_step( System system , StateInOut &x , const DerivIn &dxdt , const time_type &t , const time_type &dt )
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/controlled_error_stepper.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/controlled_error_stepper.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/controlled_error_stepper.hpp 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -318,8 +318,8 @@
boost::numeric::odeint::construct( m_dxdtnew );
m_dxdt_size_adjuster.register_state( 0 , m_dxdt );
m_xerr_size_adjuster.register_state( 0 , m_xerr );
- m_new_size_adjuster.register_state( 0 , m_xnew );
- m_new_size_adjuster.register_state( 1 , m_dxdtnew );
+ m_x_new_size_adjuster.register_state( 0 , m_xnew );
+ m_dxdt_new_size_adjuster.register_state( 0 , m_dxdtnew );
}
void copy( const controlled_error_stepper &stepper )
@@ -347,7 +347,7 @@
const error_checker_type &error_checker = error_checker_type()
)
: m_stepper( stepper ) , m_error_checker( error_checker ) ,
- m_dxdt_size_adjuster() , m_xerr_size_adjuster() , m_new_size_adjuster() ,
+ m_dxdt_size_adjuster() , m_xerr_size_adjuster() , m_x_new_size_adjuster() , m_dxdt_new_size_adjuster() ,
m_dxdt() , m_xerr() , m_xnew() , m_dxdtnew() ,
m_first_call( true )
{
@@ -356,7 +356,7 @@
controlled_error_stepper( const controlled_error_stepper &stepper )
: m_stepper( stepper.m_stepper ) , m_error_checker( stepper.m_error_checker ) ,
- m_dxdt_size_adjuster() , m_xerr_size_adjuster() , m_new_size_adjuster() ,
+ m_dxdt_size_adjuster() , m_xerr_size_adjuster() , m_x_new_size_adjuster() , m_dxdt_new_size_adjuster() ,
m_dxdt() , m_xerr() , m_xnew() , m_dxdtnew() ,
m_first_call( true )
{
@@ -383,8 +383,8 @@
// try_step( sys , x , t , dt )
- template< class System >
- controlled_step_result try_step( System system , state_type &x , time_type &t , time_type &dt )
+ template< class System , class StateInOut >
+ controlled_step_result try_step( System system , StateInOut &x , time_type &t , time_type &dt )
{
if( m_dxdt_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() ) || m_first_call )
{
@@ -396,8 +396,8 @@
}
// try_step( sys , in , t , out , dt );
- template< class System >
- controlled_step_result try_step( System system , const state_type &in , time_type &t , state_type &out , time_type &dt )
+ template< class System , class StateIn , class StateOut >
+ controlled_step_result try_step( System system , const StateIn &in , time_type &t , StateOut &out , time_type &dt )
{
if( m_dxdt_size_adjuster.adjust_size_by_policy( in , adjust_size_policy() ) || m_first_call )
{
@@ -409,10 +409,11 @@
}
// try_step( sys , x , dxdt , t , dt )
- template< class System >
- controlled_step_result try_step( System system , state_type &x , state_type &dxdt , time_type &t , time_type &dt )
+ template< class System , class StateInOut , class DerivInOut >
+ controlled_step_result try_step( System system , StateInOut &x , DerivInOut &dxdt , time_type &t , time_type &dt )
{
- m_new_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() );
+ m_x_new_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() );
+ m_dxdt_new_size_adjuster.adjust_size_by_policy( x , adjust_size_policy() );
controlled_step_result res = try_step( system , x , dxdt , t , m_xnew , m_dxdtnew , dt );
if( ( res == success_step_size_increased ) || ( res == success_step_size_unchanged) )
{
@@ -423,9 +424,9 @@
}
// try_step( sys , in , dxdt , t , out , dt )
- template< class System >
- controlled_step_result try_step( System system , const state_type &in , const state_type &dxdt_in , time_type &t ,
- state_type &out , state_type &dxdt_out , time_type &dt )
+ template< class System , class StateIn , class DerivIn , class StateOut , class DerivOut >
+ controlled_step_result try_step( System system , const StateIn &in , const DerivIn &dxdt_in , time_type &t ,
+ StateOut &out , DerivOut &dxdt_out , time_type &dt )
{
using std::max;
using std::min;
@@ -438,7 +439,7 @@
m_stepper.do_step( system , in , dxdt_in , t , out , dxdt_out , dt , m_xerr );
// this potentially overwrites m_x_err! (standard_error_checker does, at least)
- time_type max_rel_err = m_error_checker.error( in , dxdt_in , m_xerr , dt );
+ value_type max_rel_err = m_error_checker.error( in , dxdt_in , m_xerr , dt );
if( max_rel_err > 1.1 )
{
@@ -472,6 +473,8 @@
bool changed = false;
changed |= m_dxdt_size_adjuster.adjust_size( x );
changed |= m_xerr_size_adjuster.adjust_size( x );
+ changed |= m_x_new_size_adjuster.adjust_size( x );
+ changed |= m_dxdt_new_size_adjuster.adjust_size( x );
changed |= m_stepper.adjust_size( x );
if( changed )
m_first_call = true;
@@ -498,7 +501,8 @@
size_adjuster< deriv_type , 1 > m_dxdt_size_adjuster;
size_adjuster< state_type , 1 > m_xerr_size_adjuster;
- size_adjuster< state_type , 2 > m_new_size_adjuster;
+ size_adjuster< state_type , 1 > m_x_new_size_adjuster;
+ size_adjuster< deriv_type , 1 > m_dxdt_new_size_adjuster;
deriv_type m_dxdt;
state_type m_xerr;
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -94,7 +94,7 @@
m_state_adjuster() , m_deriv_adjuster() ,
m_x1() , m_x2() , m_current_state( &m_x1 ) , m_old_state( &m_x2 ) ,
m_dxdt1() , m_dxdt2() , m_current_deriv( &m_dxdt1 ) , m_old_deriv( &m_dxdt2 ) ,
- m_t( 0.0 ) , m_t_old( 0.0 ) , m_dt( 1.0 ) , m_is_deriv_initialized( false )
+ m_t() , m_t_old() , m_dt() , m_is_deriv_initialized( false )
{
initialize_variables();
}
@@ -210,7 +210,7 @@
state_type m_x1 , m_x2;
state_type *m_current_state , *m_old_state;
deriv_type m_dxdt1 , m_dxdt2;
- state_type *m_current_deriv , *m_old_deriv;
+ deriv_type *m_current_deriv , *m_old_deriv;
time_type m_t , m_t_old , m_dt;
bool m_is_deriv_initialized;
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_explicit.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_explicit.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_explicit.hpp 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -71,7 +71,7 @@
dense_output_explicit( const stepper_type &stepper = stepper_type() )
: m_stepper( stepper ) , m_size_adjuster() ,
m_x1() , m_x2() , m_current_state( &m_x1 ) , m_old_state( &m_x2 ) ,
- m_t( 0.0 ) , m_t_old( 0.0 ) , m_dt( 1.0 )
+ m_t() , m_t_old() , m_dt()
{
initialize_variables();
}
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/doc/tutorial.qbk
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/doc/tutorial.qbk (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/doc/tutorial.qbk 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -93,6 +93,12 @@
[endsect]
+[section Using arbitrary precision floating point types]
+
+GMP
+
+[endsect]
+
[endsect]
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4.cpp 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -7,6 +7,7 @@
#include <iostream>
#include <fstream>
+#include <utility>
#include <tr1/array>
#include "rosenbrock4.hpp"
@@ -52,30 +53,36 @@
const time_type R = 28.0;
const time_type b = 8.0 / 3.0;
-template< class StateType >
-void system( const StateType &x , StateType &dxdt , time_type t )
+struct lorenz
{
- dxdt[0] = sigma * ( x[1] - x[0] );
- dxdt[1] = R * x[0] - x[1] - x[0] * x[2];
- dxdt[2] = x[0] * x[1] - b * x[2];
-}
+ template< class StateType >
+ void operator()( const StateType &x , StateType &dxdt , time_type t )
+ {
+ dxdt[0] = sigma * ( x[1] - x[0] );
+ dxdt[1] = R * x[0] - x[1] - x[0] * x[2];
+ dxdt[2] = x[0] * x[1] - b * x[2];
+ }
+};
-void jacobi( const state_type &x , matrix_type &J , time_type t , state_type &dfdt )
+struct jacobi
{
- J( 0 , 0 ) = -sigma;
- J( 0 , 1 ) = sigma;
- J( 0 , 2 ) = 0.0;
- J( 1 , 0 ) = R - x[2];
- J( 1 , 1 ) = -1.0;
- J( 1 , 2 ) = -x[0];
- J( 2 , 0 ) = x[1];
- J( 2 , 1 ) = x[0];
- J( 2 , 2 ) = -b;
-
- dfdt[0] = 0.0;
- dfdt[1] = 0.0;
- dfdt[2] = 0.0;
-}
+ void operator()( const state_type &x , matrix_type &J , time_type t , state_type &dfdt )
+ {
+ J( 0 , 0 ) = -sigma;
+ J( 0 , 1 ) = sigma;
+ J( 0 , 2 ) = 0.0;
+ J( 1 , 0 ) = R - x[2];
+ J( 1 , 1 ) = -1.0;
+ J( 1 , 2 ) = -x[0];
+ J( 2 , 0 ) = x[1];
+ J( 2 , 1 ) = x[0];
+ J( 2 , 2 ) = -b;
+
+ dfdt[0] = 0.0;
+ dfdt[1] = 0.0;
+ dfdt[2] = 0.0;
+ }
+};
@@ -87,7 +94,7 @@
time_type t = 0.0 , dt = 0.00001;
stepper_type stepper;
- stepper.do_step( system< state_type > , jacobi , x , t , dt , xerr );
+ stepper.do_step( make_pair( lorenz() , jacobi() ) , x , t , dt , xerr );
controlled_stepper_type controlled_stepper;
x[0] = 1.0 ; x[1] = 1.0 ; x[2] = 0.0;
@@ -105,7 +112,7 @@
size_t trials = 0;
while( trials < 100 )
{
- if( controlled_stepper.try_step( system< state_type > , jacobi , x , t , dt ) != step_size_decreased )
+ if( controlled_stepper.try_step( make_pair( lorenz() , jacobi() ) , x , t , dt ) != step_size_decreased )
break;
// clog.precision( 14 );
// clog << dt << "\n";
@@ -146,7 +153,7 @@
size_t trials = 0;
while( trials < 100 )
{
- if( stepper.try_step( system< state_type2 > , x , t , dt ) != step_size_decreased )
+ if( stepper.try_step( lorenz() , x , t , dt ) != step_size_decreased )
break;
++trials;
}
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4.hpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4.hpp 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -8,6 +8,8 @@
#ifndef ROSENBROCK4_HPP_
#define ROSENBROCK4_HPP_
+#include <boost/ref.hpp>
+
#include <boost/numeric/ublas/vector.hpp>
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/lu.hpp>
@@ -35,8 +37,10 @@
{
public:
- typedef Value time_type;
- typedef boost::numeric::ublas::vector< time_type > state_type;
+ typedef Value value_type;
+ typedef boost::numeric::ublas::vector< value_type > state_type;
+ typedef state_type deriv_type;
+ typedef value_type time_type;
typedef boost::numeric::ublas::matrix< time_type > matrix_type;
typedef boost::numeric::ublas::permutation_matrix< size_t > pmatrix_type;
@@ -44,33 +48,36 @@
{
}
-// void do_step( System &system , Jacobi &jacobi , ... );
-// void do_step( pair< System& , Jacobi& > system , ... );
-// void do_step( pair< System , Jacobi > & , ... );
-
-
- template< class System , class Jacobi >
- void do_step( System &system , Jacobi &jacobi , const state_type &x , time_type t , state_type &xout , time_type dt , state_type &xerr )
- {
- const double gamma = 0.25;
- const double d1 = 0.25 , d2 = -0.1043 , d3 = 0.1035 , d4 = 0.3620000000000023e-01;
- const double c2 = 0.386 , c3 = 0.21 , c4 = 0.63;
- const double c21 = -0.5668800000000000e+01;
- const double a21 = 0.1544000000000000e+01;
- const double c31 = -0.2430093356833875e+01 , c32 = -0.2063599157091915e+00;
- const double a31 = 0.9466785280815826e+00 , a32 = 0.2557011698983284e+00;
- const double c41 = -0.1073529058151375e+00 , c42 = -0.9594562251023355e+01 , c43 = -0.2047028614809616e+02;
- const double a41 = 0.3314825187068521e+01 , a42 = 0.2896124015972201e+01 , a43 = 0.9986419139977817e+00;
- const double c51 = 0.7496443313967647e+01 , c52 = -0.1024680431464352e+02 , c53 = -0.3399990352819905e+02 , c54 = 0.1170890893206160e+02;
- const double a51 = 0.1221224509226641e+01 , a52 = 0.6019134481288629e+01 , a53 = 0.1253708332932087e+02 , a54 = -0.6878860361058950e+00 ;
- const double c61 = 0.8083246795921522e+01 , c62 = -0.7981132988064893e+01 , c63 = -0.3152159432874371e+02 , c64 = 0.1631930543123136e+02 , c65 = -0.6058818238834054e+01;
+ template< class System >
+ void do_step( System system , const state_type &x , time_type t , state_type &xout , time_type dt , state_type &xerr )
+ {
+ const value_type gamma = 0.25;
+ const value_type d1 = 0.25 , d2 = -0.1043 , d3 = 0.1035 , d4 = 0.3620000000000023e-01;
+ const value_type c2 = 0.386 , c3 = 0.21 , c4 = 0.63;
+ const value_type c21 = -0.5668800000000000e+01;
+ const value_type a21 = 0.1544000000000000e+01;
+ const value_type c31 = -0.2430093356833875e+01 , c32 = -0.2063599157091915e+00;
+ const value_type a31 = 0.9466785280815826e+00 , a32 = 0.2557011698983284e+00;
+ const value_type c41 = -0.1073529058151375e+00 , c42 = -0.9594562251023355e+01 , c43 = -0.2047028614809616e+02;
+ const value_type a41 = 0.3314825187068521e+01 , a42 = 0.2896124015972201e+01 , a43 = 0.9986419139977817e+00;
+ const value_type c51 = 0.7496443313967647e+01 , c52 = -0.1024680431464352e+02 , c53 = -0.3399990352819905e+02 , c54 = 0.1170890893206160e+02;
+ const value_type a51 = 0.1221224509226641e+01 , a52 = 0.6019134481288629e+01 , a53 = 0.1253708332932087e+02 , a54 = -0.6878860361058950e+00 ;
+ const value_type c61 = 0.8083246795921522e+01 , c62 = -0.7981132988064893e+01 , c63 = -0.3152159432874371e+02 , c64 = 0.1631930543123136e+02 , c65 = -0.6058818238834054e+01;
+
+ typedef typename boost::unwrap_reference< System >::type system_type;
+ typedef typename boost::unwrap_reference< typename system_type::first_type >::type deriv_func_type;
+ typedef typename boost::unwrap_reference< typename system_type::second_type >::type jacobi_func_type;
+ system_type &sys = system;
+ deriv_func_type &deriv_func = sys.first;
+ jacobi_func_type &jacobi_func = sys.second;
+
const size_t n = x.size();
matrix_type jac( n , n );
pmatrix_type pm( n );
state_type dfdt( n ) , dxdt( n );
- system( x , dxdt , t );
- jacobi( x , jac , t , dfdt );
+ deriv_func( x , dxdt , t );
+ jacobi_func( x , jac , t , dfdt );
state_type g1( n ) , g2( n ) , g3( n ) , g4( n ) , g5( n );
state_type xtmp( n ) , dxdtnew( n );
@@ -87,7 +94,7 @@
for( size_t i=0 ; i<n ; ++i )
xtmp[i] = x[i] + a21 * g1[i];
- system( xtmp , dxdtnew , t + c2 * dt );
+ deriv_func( xtmp , dxdtnew , t + c2 * dt );
for( size_t i=0 ; i<n ; ++i )
g2[i] = dxdtnew[i] + dt * d2 * dfdt[i] + c21 * g1[i] / dt;
boost::numeric::ublas::lu_substitute( jac , pm , g2 );
@@ -95,7 +102,7 @@
for( size_t i=0 ; i<n ; ++i )
xtmp[i] = x[i] + a31 * g1[i] + a32 * g2[i];
- system( xtmp , dxdtnew , t + c3 * dt );
+ deriv_func( xtmp , dxdtnew , t + c3 * dt );
for( size_t i=0 ; i<n ; ++i )
g3[i] = dxdtnew[i] + dt * d3 * dfdt[i] + ( c31 * g1[i] + c32 * g2[i] ) / dt;
boost::numeric::ublas::lu_substitute( jac , pm , g3 );
@@ -103,7 +110,7 @@
for( size_t i=0 ; i<n ; ++i )
xtmp[i] = x[i] + a41 * g1[i] + a42 * g2[i] + a43 * g3[i];
- system( xtmp , dxdtnew , t + c4 * dt );
+ deriv_func( xtmp , dxdtnew , t + c4 * dt );
for( size_t i=0 ; i<n ; ++i )
g4[i] = dxdtnew[i] + dt * d4 * dfdt[i] + ( c41 * g1[i] + c42 * g2[i] + c43 * g3[i] ) / dt;
boost::numeric::ublas::lu_substitute( jac , pm , g4 );
@@ -111,14 +118,14 @@
for( size_t i=0 ; i<n ; ++i )
xtmp[i] = x[i] + a51 * g1[i] + a52 * g2[i] + a53 * g3[i] + a54 * g4[i];
- system( xtmp , dxdtnew , t + dt );
+ deriv_func( xtmp , dxdtnew , t + dt );
for( size_t i=0 ; i<n ; ++i )
g5[i] = dxdtnew[i] + ( c51 * g1[i] + c52 * g2[i] + c53 * g3[i] + c54 * g4[i] ) / dt;
boost::numeric::ublas::lu_substitute( jac , pm , g5 );
for( size_t i=0 ; i<n ; ++i )
xtmp[i] += g5[i];
- system( xtmp , dxdtnew , t + dt );
+ deriv_func( xtmp , dxdtnew , t + dt );
for( size_t i=0 ; i<n ; ++i )
xerr[i] = dxdtnew[i] + ( c61 * g1[i] + c62 * g2[i] + c63 * g3[i] + c64 * g4[i] + c65 * g5[i] ) / dt;
boost::numeric::ublas::lu_substitute( jac , pm , xerr );
@@ -127,12 +134,10 @@
xout[i] = xtmp[i] + xerr[i];
}
- template< class System , class Jacobi >
- void do_step( System &system , Jacobi &jacobi , state_type &x , time_type t , time_type dt , state_type &xerr )
+ template< class System >
+ void do_step( System system , state_type &x , time_type t , time_type dt , state_type &xerr )
{
- state_type out( x.size() );
- do_step( system , jacobi , x , t , out , dt , xerr );
- x = out;
+ do_step( system , x , t , x , dt , xerr );
}
@@ -177,15 +182,15 @@
}
- template< class System , class Jacobi >
+ template< class System >
boost::numeric::odeint::controlled_step_result
- try_step( System &sys , Jacobi &jacobi , state_type &x , time_type &t , time_type &dt )
+ try_step( System sys , state_type &x , time_type &t , time_type &dt )
{
static const time_type safe = 0.9 , fac1 = 5.0 , fac2 = 1.0 / 6.0;
const size_t n = x.size();
state_type xnew( n ) , xerr( n );
- m_rb4.do_step( sys , jacobi , x , t , xnew , dt , xerr );
+ m_rb4.do_step( sys , x , t , xnew , dt , xerr );
time_type err = error( xnew , x , xerr );
time_type fac = std::max( fac2 ,std::min( fac1 , std::pow( err , 0.25 ) / safe ) );
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4_stepper.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4_stepper.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/rosenbrock4/rosenbrock4_stepper.cpp 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -7,6 +7,7 @@
#include <iostream>
#include <fstream>
+#include <utility>
#include <tr1/array>
#include "rosenbrock4.hpp"
@@ -52,13 +53,16 @@
const time_type R = 28.0;
const time_type b = 8.0 / 3.0;
-template< class StateType >
-void system( const StateType &x , StateType &dxdt , time_type t )
+struct lorenz
{
- dxdt[0] = sigma * ( x[1] - x[0] );
- dxdt[1] = R * x[0] - x[1] - x[0] * x[2];
- dxdt[2] = x[0] * x[1] - b * x[2];
-}
+ template< class StateType >
+ void operator()( const StateType &x , StateType &dxdt , time_type t )
+ {
+ dxdt[0] = sigma * ( x[1] - x[0] );
+ dxdt[1] = R * x[0] - x[1] - x[0] * x[2];
+ dxdt[2] = x[0] * x[1] - b * x[2];
+ }
+};
void jacobi( const state_type &x , matrix_type &J , time_type t , state_type &dfdt )
{
@@ -102,7 +106,7 @@
fout << xerr[0] << "\t" << xerr[1] << "\t" << xerr[2] << "\t";
fout <<std::endl;
- stepper.do_step( system< state_type > , jacobi , x , t , dt , xerr );
+ stepper.do_step( make_pair( lorenz() , jacobi ) , x , t , dt , xerr );
++count;
t += dt;
}
@@ -130,7 +134,7 @@
fout << xerr[0] << "\t" << xerr[1] << "\t" << xerr[2] << "\t";
fout <<std::endl;
- rk_stepper.do_step( system< state_type2 > , x , t , dt , xerr );
+ rk_stepper.do_step( lorenz() , x , t , dt , xerr );
++count;
t += dt;
}
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 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -27,4 +27,5 @@
[ run stepper_with_units.cpp ]
[ run stepper_copying.cpp ]
[ run stepper_with_ranges.cpp ]
- ;
+ : <testing.launcher>valgrind
+ ;
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_ranges.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_ranges.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_ranges.cpp 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -13,6 +13,8 @@
#include <tr1/array>
#include <utility>
+#include <boost/range.hpp>
+
#include <boost/numeric/odeint/stepper/explicit_euler.hpp>
typedef std::vector< double > state_type;
@@ -21,15 +23,14 @@
struct lorenz
{
template< class State , class Deriv >
- void operator()( const State &x , Deriv &dxdt , double t )
+ void operator()( const State &x_ , Deriv &dxdt_ , double t )
{
- const double sigma = 10.0;
- const double R = 28.0;
- const double b = 8.0 / 3.0;
-
-// dxdt[0] = sigma * ( x[1] - x[0] );
-// dxdt[1] = R * x[0] - x[1] - x[0] * x[2];
-// dxdt[2] = x[0]*x[1] - b * x[2];
+ typename boost::range_iterator< const State >::type x = boost::begin( x_ );
+ typename boost::range_iterator< Deriv >::type dxdt = boost::begin( dxdt_ );
+
+ dxdt[0] = x[0];
+ dxdt[1] = 2.0;
+ dxdt[2] = 3.0;
}
};
@@ -38,10 +39,14 @@
BOOST_AUTO_TEST_CASE( explicit_euler_with_range )
{
std::vector< double > x( 3 * 2 );
+ x[0] = 1.0;
+ x[1] = 1.0;
+ x[2] = 1.0;
boost::numeric::odeint::explicit_euler< state_type > euler;
- std::pair< std::vector< double >::iterator , std::vector< double >::iterator > r( x.begin() , x.begin() + 3 );
- euler.do_step( lorenz() , r , 0.0 , 0.1 );
- euler.do_step( lorenz() , std::make_pair( x.begin() + 3 , x.begin() + 6 ) , 0.1 , 0.1 );
+ euler.do_step( lorenz() , std::make_pair( x.begin() , x.begin() + 3 ) , 0.1 , 0.1 );
+ BOOST_CHECK_CLOSE( x[0] , 1.1 , 1.0e-10 );
+ BOOST_CHECK_CLOSE( x[1] , 1.2 , 1.0e-10 );
+ BOOST_CHECK_CLOSE( x[2] , 1.3 , 1.0e-10 );
}
BOOST_AUTO_TEST_CASE( explicit_euler_with_array )
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_units.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_units.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_units.cpp 2011-01-30 14:17:56 EST (Sun, 30 Jan 2011)
@@ -32,6 +32,8 @@
#include <boost/numeric/odeint/stepper/explicit_error_rk54_ck.hpp>
#include <boost/numeric/odeint/stepper/explicit_error_dopri5.hpp>
#include <boost/numeric/odeint/stepper/controlled_error_stepper.hpp>
+#include <boost/numeric/odeint/stepper/dense_output_explicit.hpp>
+#include <boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp>
#include <boost/numeric/odeint/algebra/fusion_algebra.hpp>
@@ -53,7 +55,6 @@
void oscillator( const state_type &x , deriv_type &dxdt , time_type t )
{
-
}
template< class Stepper >
@@ -202,6 +203,26 @@
}
+template< class Stepper >
+void check_dense_output_stepper( Stepper &stepper )
+{
+ typedef Stepper stepper_type;
+ typedef typename stepper_type::state_type state_type;
+ typedef typename stepper_type::value_type value_type;
+ typedef typename stepper_type::deriv_type deriv_type;
+ typedef typename stepper_type::time_type time_type;
+// typedef typename stepper_type::order_type order_type;
+
+ time_type t( 0.0 * si::second );
+ time_type dt( 0.1 * si::second );
+ state_type x( 1.0 * si::meter , 0.0 * si::meter_per_second ) , x2;
+
+ stepper.initialize( x , t , dt );
+ stepper.do_step( oscillator );
+ stepper.calc_state( dt / 2.0 , x2 );
+}
+
+
@@ -230,9 +251,16 @@
class controlled_stepper_types : public mpl::vector
<
- controlled_error_stepper< explicit_error_rk54_ck< state_type , value_type , deriv_type , time_type , fusion_algebra > >
+ controlled_error_stepper< explicit_error_rk54_ck< state_type , value_type , deriv_type , time_type , fusion_algebra > > ,
+ controlled_error_stepper< explicit_error_dopri5< state_type , value_type , deriv_type , time_type , fusion_algebra > >
> { };
+class dense_output_stepper_types : public mpl::vector
+<
+ dense_output_explicit< explicit_euler< state_type , value_type , deriv_type , time_type , fusion_algebra > > ,
+ dense_output_controlled_explicit_fsal<
+ controlled_error_stepper< explicit_error_dopri5< state_type , value_type , deriv_type , time_type , fusion_algebra > > >
+> { };
@@ -269,5 +297,11 @@
check_controlled_stepper( stepper );
}
+BOOST_AUTO_TEST_CASE_TEMPLATE( dense_ouput_test , Stepper , dense_output_stepper_types )
+{
+ Stepper stepper;
+ check_dense_output_stepper( stepper );
+}
+
BOOST_AUTO_TEST_SUITE_END()