Go to the documentation of this file.00001 #ifndef OPENTISSUE_DYNAMICS_MBD_UTIL_SIMULATORS_MBD_FIX_POINT_STEP_SIMULATOR_H
00002 #define OPENTISSUE_DYNAMICS_MBD_UTIL_SIMULATORS_MBD_FIX_POINT_STEP_SIMULATOR_H
00003
00004
00005
00006
00007
00008
00009
00010 #include <OpenTissue/configuration.h>
00011
00012 #include <OpenTissue/dynamics/mbd/interfaces/mbd_simulator_interface.h>
00013
00014 #include <boost/cast.hpp>
00015
00016 namespace OpenTissue
00017 {
00018 namespace mbd
00019 {
00025 template< typename mbd_types >
00026 class FixPointStepSimulator
00027 : public SimulatorInterface<mbd_types>
00028 {
00029 protected:
00030
00031 typedef typename mbd_types::math_policy math_policy;
00032 typedef typename math_policy::real_type real_type;
00033 typedef typename math_policy::value_traits value_traits;
00034 typedef typename math_policy::vector_type vector_type;
00035 typedef typename mbd_types::group_type group_type;
00036 typedef typename mbd_types::group_ptr_container group_ptr_container;
00037
00038 public:
00039
00040 class node_traits{};
00041 class edge_traits{};
00042 class constraint_traits{};
00043
00044 protected:
00045
00046 vector_type m_st;
00047 vector_type m_s;
00048 vector_type m_ss;
00049 vector_type m_u;
00050 group_type * m_all;
00051 group_ptr_container m_groups;
00052
00053
00054 public:
00055
00056 FixPointStepSimulator(){}
00057
00058 virtual ~FixPointStepSimulator(){}
00059
00060 public:
00061
00062 void run(real_type const & time_step)
00063 {
00064 real_type m_epsilon_fix = boost::numeric_cast<real_type>(10e-4);
00065
00066 m_all = this->get_configuration()->get_all_body_group();
00067
00068 mbd::get_position_vector(*m_all, m_s);
00069 mbd::get_velocity_vector(*m_all, m_u);
00070 mbd::compute_position_update(*m_all,m_st,m_u,time_step,m_s);
00071 mbd::set_position_vector(*m_all,m_s);
00072
00073 mbd::compute_scripted_motions(*m_all,this->time() + time_step);
00074
00075 do
00076 {
00077 mbd::get_position_vector(*m_all, m_s);
00078
00079 this->get_collision_detection()->run( m_groups );
00080
00081 for(typename group_ptr_container::iterator tmp=m_groups.begin();tmp!=m_groups.end();++tmp)
00082 {
00083 group_type * group = (*tmp);
00084
00085 this->get_sleepy()->evaluate(group->body_begin(),group->body_end());
00086
00087 if(!mbd::is_all_bodies_sleepy(*group))
00088 this->get_stepper()->run(*group,time_step);
00089 }
00090 mbd::get_velocity_vector(*m_all, m_u);
00091 mbd::compute_position_update(*m_all,m_st,m_u,time_step,m_ss);
00092 mbd::set_position_vector(*m_all,m_ss);
00093 }
00094
00095 while(norm(m_s,m_ss)>=m_epsilon_fix);
00096
00097 update_time(time_step);
00098 }
00099
00100 protected:
00101
00110 real_type norm(vector_type & A,vector_type & B)
00111 {
00112 using std::fabs;
00113
00114 real_type L_inf = value_traits::zero();
00115
00116 typename vector_type::iterator end = A.end();
00117 typename vector_type::iterator a = A.begin();
00118 typename vector_type::iterator b = B.begin();
00119
00120 for(;a!=end;++a,++b)
00121 {
00122 real_type diff = fabs( (*a)-(*b) );
00123 if(diff > L_inf)
00124 L_inf = diff;
00125 }
00126 return L_inf;
00127 }
00128
00129 };
00130
00131 }
00132 }
00133
00134 #endif