// integrator_euler.cc // // Atomic class that defines Euler integrator defined as CLASSTYPE INTEGRATOR. // // - The initial condition for the dynamic variable x is returned for a static step at initial time // - The integration scheme is invoked for a dynamic step // - If VariableTimeStep ( 1 ()) in the *.run file, then the time step is adapted in order // to satisfy the integration tolerance specified in the GlobalSettings section of the // *.prf file with the key Tolerance. // /////////////////////////////////////////////////////////////////////////////// #ifdef spark_parser CLASSTYPE INTEGRATOR ; PORT x; PORT xdot; PORT dt; EQUATIONS { x = integrator_euler( xdot ); bad_inverses = xdot, dt ; } FUNCTIONS { x = EVALUATE integrator_euler__evaluate( xdot, dt ) PREDICT integrator_euler__predict( xdot, dt ) STATIC_CONSTRUCT integrator_euler__static_construct() CONSTRUCT integrator_euler__construct( x, xdot ) STATIC_PREPARE_STEP integrator_euler__static_prepare_step() PREPARE_STEP integrator_euler__prepare_step( x, xdot, dt ) CHECK_INTEGRATION_STEP integrator_euler__check_integration_step( x ) STATIC_CHECK_INTEGRATION_STEP integrator_euler__static_check_integration_step() COMMIT integrator_euler__commit( x, xdot ) STATIC_COMMIT integrator_euler__static_commit() ROLLBACK integrator_euler__rollback( x, xdot ) STATIC_ROLLBACK integrator_euler__static_rollback() DESTRUCT integrator_euler__destruct() STATIC_DESTRUCT integrator_euler__static_destruct() ; } #endif // spark_parser // The debugger can't handle symbols more than 255 characters long. // STL often creates symbols longer than that. // When symbols are longer than 255 characters, the warning is disabled. #pragma warning(disable:4786) #include "spark.h" ////////////////// // PRIVATE DATA // ////////////////// /////////////////////////////////////////////////////////////////////////////// // Include implementation files with support for the integration methods /////////////////////////////////////////////////////////////////////////////// #include "integrators/euler.h" /////////////////////////////////////////////////////////////////////////////// // Activate the namespace of the desired integration method /////////////////////////////////////////////////////////////////////////////// using namespace SPARK::Euler; /////////////////////////////////////////////////////////////////////////////// // Define shortcut types /////////////////////////////////////////////////////////////////////////////// typedef TIntegrationMethod::integration_controller_type TIntegrationController; typedef TIntegrationMethod::integrator_type TIntegrator; /////////////// // CALLBACKS // /////////////// /////////////////////////////////////////////////////////////////////////////// STATIC_CONSTRUCT( integrator_euler__static_construct ) { TIntegrationController* MyController = TIntegrationMethod::make_integration_controller( THIS->GetNumObjects(), // Number of instances ACTIVE_PROBLEM->GetName(), // Problem name false // Generate verbose integration diagnostic? ); // Abort simulation if memory could not be allocated! if ( !MyController ) { REQUEST__ABORT( "Could not allocate memory for the integration controller!" ); } // Store pointer to private static data within this inverse SET_DATA( THIS, TIntegrationController, MyController ); } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// CONSTRUCT( integrator_euler__construct ) { ARGUMENT( 0, x ); ARGUMENT( 1, xdot ); GET_DATA( ACTIVE_INVERSE, TIntegrationController, MyController ); TIntegrator* MyIntegrator = TIntegrationMethod::make_integrator( THIS->GetInstanceHandle(), // Unique instance handle x, // Dynamic variable's name xdot, // Time-derivative's name MyController // Integration controller ); // Abort simulation if memory could not be allocated! if ( !MyIntegrator ) { REQUEST__ABORT( "Could not allocate memory for the integrator!" ); } // Store pointer to private data within this object SET_DATA( THIS, TIntegrator, MyIntegrator ); } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// STATIC_PREPARE_STEP( integrator_euler__static_prepare_step ) { // Retrieve integration step if dynamic step if ( !ACTIVE_PROBLEM->IsStaticStep() ) { GET_DATA( THIS, TIntegrationController, MyController); MyController->PrepareIntegrationStep( ACTIVE_PROBLEM->GetGlobalTimeStep() ); } } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// PREPARE_STEP( integrator_euler__prepare_step ) { // Prepare integrator private data at beginning of integration step if ( !ACTIVE_PROBLEM->IsStaticStep() ) { ARGUMENT( 0, x ); ARGUMENT( 1, xdot ); ARGUMENT( 2, dt ); GET_DATA( THIS, TIntegrator, MyIntegrator); MyIntegrator->PrepareIntegrationStep( x, xdot, dt ); } } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// PREDICT( integrator_euler__predict ) { TARGET( 0, x ); ARGUMENT( 0, xdot ); ARGUMENT( 1, dt ); double result; if ( ACTIVE_PROBLEM->IsStaticStep() ) { // Static calculation result = ( ACTIVE_PROBLEM->IsInitialTime() ? x.GetInit() // Initial time solution special case : x[1] // Past value for restart after initial time solution ); } else { // Compute prediction using past value of time-derivative GET_DATA( THIS, TIntegrator, MyIntegrator); result = MyIntegrator->Predict( x, xdot, dt ); } // Update dynamic variable with result x = result; } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// EVALUATE( integrator_euler__evaluate ) { TARGET( 0, x ); ARGUMENT( 0, xdot); ARGUMENT( 1, dt); double result; if ( ACTIVE_PROBLEM->IsStaticStep() ) { // Static calculation result = ( ACTIVE_PROBLEM->IsInitialTime() ? x.GetInit() // Initial time solution : x[1] // Past value for warm restart ); } else { // Perform the actual integration using current value of time-derivative GET_DATA( THIS, TIntegrator, MyIntegrator); result = MyIntegrator->Integrate( x, xdot, dt ); } // Update dynamic variable with result x = result; } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// CHECK_INTEGRATION_STEP( integrator_euler__check_integration_step ) { // Estimate integration error only if: // - variable time step mode // - this is not a static step if ( ACTIVE_PROBLEM->IsStaticStep() ) { return true; } else { ARGUMENT( 0, x); GET_DATA( THIS, TIntegrator, MyIntegrator); MyIntegrator->EstimateLocalError( x ); // Nothing to do. Decision is postponed to static check integration step callback return true; } } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// STATIC_CHECK_INTEGRATION_STEP( integrator_euler__static_check_integration_step ) { // Check integration error only if: // - variable time step mode // - this is not a static step if ( ACTIVE_PROBLEM->IsStaticStep() ) { return true; } else { if ( ACTIVE_PROBLEM->IsTimeStepVariable() ) { GET_DATA( THIS, TIntegrationController, MyController); return MyController->CheckIntegrationStep( THIS->GetTolerance() ); } else { return true; } } } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// COMMIT( integrator_euler__commit ) { GET_DATA( THIS, TIntegrator, MyIntegrator); if ( ACTIVE_PROBLEM->IsStaticStep() ) { MyIntegrator->Restart(); } else { ARGUMENT( 0, x); ARGUMENT( 1, xdot); MyIntegrator->Commit( x, xdot ); } } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// STATIC_COMMIT( integrator_euler__static_commit ) { // Predict next time step only if: // - variable time step mode // - this is not a static step GET_DATA( THIS, TIntegrationController, MyController); if ( ACTIVE_PROBLEM->IsStaticStep() ) { MyController->Restart(); } else { if ( ACTIVE_PROBLEM->IsTimeStepVariable() ) { // Requests time step for next step to satisfy the user-requested // tolerance for the next dynamic step. REQUEST__SET_TIME_STEP( MyController->GetName(), MyController->PredictTimeStep( THIS->GetTolerance() ) ); } MyController->Commit(); } } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// ROLLBACK( integrator_euler__rollback ) { GET_DATA( THIS, TIntegrator, MyIntegrator); if ( ACTIVE_PROBLEM->IsStaticStep() ) { MyIntegrator->Restart(); } else { ARGUMENT( 0, x); ARGUMENT( 1, xdot); MyIntegrator->Rollback( x, xdot ); } } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// STATIC_ROLLBACK( integrator_euler__static_rollback ) { // Correct current time step only if: // - variable time step mode // - this is not a static step GET_DATA( THIS, TIntegrationController, MyController); if ( ACTIVE_PROBLEM->IsStaticStep() ) { MyController->Restart(); } else{ if ( ACTIVE_PROBLEM->IsTimeStepVariable() ) { // Requests time step for current step after rejection to satisfy the user-requested // tolerance for the next dynamic step. REQUEST__SET_TIME_STEP( MyController->GetName(), MyController->CorrectTimeStep( THIS->GetTolerance() ) ); } MyController->Rollback(); } } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// DESTRUCT( integrator_euler__destruct ) { // Deallocate instance private data DELETE_DATA( THIS, TIntegrator, false ); } /////////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////////// STATIC_DESTRUCT( integrator_euler__static_destruct ) { // Deallocate static private data DELETE_DATA( THIS, TIntegrationController, false ); } ///////////////////////////////////////////////////////////////////////////////