System.h

Go to the documentation of this file.

#ifndef SimTK_SimTKCOMMON_SYSTEM_H_
#define SimTK_SimTKCOMMON_SYSTEM_H_

/* -------------------------------------------------------------------------- *
 *                      SimTK Core: SimTKcommon                               *
 * -------------------------------------------------------------------------- *
 * This is part of the SimTK Core biosimulation toolkit originating from      *
 * Simbios, the NIH National Center for Physics-Based Simulation of           *
 * Biological Structures at Stanford, funded under the NIH Roadmap for        *
 * Medical Research, grant U54 GM072970. See https://simtk.org.               *
 *                                                                            *
 * Portions copyright (c) 2006-10 Stanford University and the Authors.        *
 * Authors: Michael Sherman                                                   *
 * Contributors: Peter Eastman                                                *
 *                                                                            *
 * Permission is hereby granted, free of charge, to any person obtaining a    *
 * copy of this software and associated documentation files (the "Software"), *
 * to deal in the Software without restriction, including without limitation  *
 * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
 * and/or sell copies of the Software, and to permit persons to whom the      *
 * Software is furnished to do so, subject to the following conditions:       *
 *                                                                            *
 * The above copyright notice and this permission notice shall be included in *
 * all copies or substantial portions of the Software.                        *
 *                                                                            *
 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
 * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
 * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
 * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
 * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
 * -------------------------------------------------------------------------- */

#include "SimTKcommon/basics.h"
#include "SimTKcommon/Simmatrix.h"
#include "SimTKcommon/internal/State.h"
#include "SimTKcommon/internal/Subsystem.h"

#include <cassert>

namespace SimTK {

class DecorativeGeometry;
class DefaultSystemSubsystem;
class ScheduledEventHandler;
class ScheduledEventReporter;
class TriggeredEventHandler;
class TriggeredEventReporter;

class SimTK_SimTKCOMMON_EXPORT System {
public:
    class Guts; // local; name is System::Guts
    friend class Guts;
private:
    // This is the only data member in this class. Also, any class derived from
    // System must have *NO* data members at all (data goes in the Guts class).
    Guts* guts;
public:
    System() : guts(0) { }
    System(const System&);
    System& operator=(const System&);
    ~System();

    const String& getName()    const;
    const String& getVersion() const;

    inline void addEventHandler(ScheduledEventHandler* handler);
    inline void addEventHandler(TriggeredEventHandler* handler);
    inline void addEventReporter(ScheduledEventReporter* handler) const;
    inline void addEventReporter(TriggeredEventReporter* handler) const;

    System& setUpDirection(const CoordinateDirection& up);
    CoordinateDirection getUpDirection() const;

    System& setUseUniformBackground(bool useUniformBackground);
    bool getUseUniformBackground() const;


    class ProjectOptions {
        unsigned int optionSet;
        explicit ProjectOptions(unsigned o) : optionSet(o) { }
    public:

        enum Option {
            None   = 0x00,

            Q      = 0x01,
            U      = 0x02,
            QError = 0x04,
            UError = 0x08,

            PositionOnly = (Q|QError),
            VelocityOnly = (U|UError),
            All          = (PositionOnly|VelocityOnly),

            // This option says we expect the state to be close to a 
            // solution already and restricts projection to move downhill
            // in the local vicinity. This should be used during any
            // continuous integration to prevent erroneous jumps in the state.
            LocalOnly = 0x1000
        };

        ProjectOptions() : optionSet(0) { }

        // This is an implicit conversion
        ProjectOptions(Option opt) : optionSet((unsigned)opt) { }

        // Implicit conversion to bool when needed
        operator bool() const {return optionSet != 0;}
        bool hasAnyPositionOptions() const {return (optionSet&(unsigned)PositionOnly) != 0;}
        bool hasAnyVelocityOptions() const {return (optionSet&(unsigned)VelocityOnly) != 0;}
        bool isEmpty() const {return optionSet==0;}

        bool isOptionSet(Option opt) const {return (optionSet&(unsigned)opt) != 0;}
        void clear() {optionSet=0;}
        void clearOption(Option opt) {optionSet &= ~(unsigned)opt;}
        void setOption  (Option opt) {optionSet |= (unsigned)opt;}

        // Set operators: not, or, and, set difference
        ProjectOptions operator~() const {return ProjectOptions( (~optionSet) & (unsigned)All );}
        ProjectOptions& operator|=(ProjectOptions opts) {optionSet |= opts.optionSet; return *this;}
        ProjectOptions& operator&=(ProjectOptions opts) {optionSet &= opts.optionSet; return *this;}
        ProjectOptions& operator-=(ProjectOptions opts) {optionSet &= ~opts.optionSet; return *this;}

        ProjectOptions& operator|=(Option opt) {setOption(opt); return *this;}
        ProjectOptions& operator-=(Option opt) {clearOption(opt); return *this;}
    };

        // STATISTICS //

    void resetAllCountersToZero();

        // Realization

    int getNumRealizationsOfThisStage(Stage) const;

    int getNumRealizeCalls() const;

        // Prescribed motion

    int getNumPrescribeCalls() const;

        // Projection

    int getNumQProjections() const;
    int getNumUProjections() const;
    int getNumQErrorEstimateProjections() const;
    int getNumUErrorEstimateProjections() const;

    int getNumProjectCalls() const;

        // Event handling and reporting

    int getNumHandlerCallsThatChangedStage(Stage) const;

    int getNumHandleEventCalls() const;

    int getNumReportEventCalls() const;


        // REALIZATION //

    const State& realizeTopology() const;

    const State& getDefaultState() const;
    State&       updDefaultState();

    void realizeModel(State&) const;

    void realize(const State& s, Stage g = Stage::HighestRuntime) const;

    void calcDecorativeGeometryAndAppend(const State&, Stage, 
                                         Array_<DecorativeGeometry>&) const;

        // THE CONTINUOUS SYSTEM //

            // UNCONSTRAINED

    Real calcTimescale(const State&) const;

    void calcYUnitWeights(const State&, Vector& weights) const;

            // CONSTRAINED

    void prescribe(State&, Stage) const;

    void project(State&, Real consAccuracy, const Vector& yWeights,
                 const Vector& cWeights, Vector& yerrest, 
                 ProjectOptions=ProjectOptions::All) const;

    void calcYErrUnitTolerances(const State&, Vector& tolerances) const;

            // FAST VARIABLES

    void relax(State&, Stage, Real accuracy, 
               const Vector& yWeights, const Vector& cWeights) const;


        // THE DISCRETE (SLOW) SYSTEM //

    void setHasTimeAdvancedEvents(bool); // default=false
    bool hasTimeAdvancedEvents() const;

    void handleEvents
        (State&, Event::Cause, const Array_<EventId>& eventIds,
        Real accuracy, const Vector& yWeights, const Vector& cWeights,
        Stage& lowestModified, bool& shouldTerminate) const;
    
    void reportEvents(const State& s, Event::Cause cause, 
                      const Array_<EventId>& eventIds) const;

    void calcEventTriggerInfo(const State&, Array_<EventTriggerInfo>&) const;

    void calcTimeOfNextScheduledEvent(const State&, Real& tNextEvent, 
                                      Array_<EventId>& eventIds, bool includeCurrentTime) const;

    void calcTimeOfNextScheduledReport(const State&, Real& tNextEvent, 
                                       Array_<EventId>& eventIds, bool includeCurrentTime) const;
    
    //TODO: these operators should be provided by the Vector class where they
    //can be performed more efficiently.

    static Real calcWeightedRMSNorm(const Vector& values, const Vector& weights) {
        assert(weights.size() == values.size());
        if (values.size()==0) return 0;
        Real sumsq = 0;
        for (int i=0; i<values.size(); ++i) {
            const Real wv = weights[i]*values[i];
            sumsq += wv*wv;
        }
        return std::sqrt(sumsq/weights.size());
    }

    static Real calcWeightedInfinityNorm(const Vector& values, const Vector& weights) {
        assert(weights.size() == values.size());
        if (values.size()==0) return 0;
        Real maxval = 0;
        for (int i=0; i<values.size(); ++i) {
            const Real wv = std::abs(weights[i]*values[i]);
            if (wv > maxval) maxval=wv;
        }
        return maxval;
    }

    SubsystemIndex adoptSubsystem(Subsystem& child);

    int getNumSubsystems() const;
    const Subsystem& getSubsystem(SubsystemIndex)   const;
    Subsystem&       updSubsystem(SubsystemIndex);
    const DefaultSystemSubsystem& getDefaultSubsystem() const;
    DefaultSystemSubsystem& updDefaultSubsystem();

    inline operator const Subsystem&() const; // implemented below
    inline operator Subsystem&();

    // Internal use only
    bool isOwnerHandle() const;
    bool isEmptyHandle() const;

    // There can be multiple handles on the same System.
    bool isSameSystem(const System& otherSystem) const;

    bool systemTopologyHasBeenRealized() const;

    // dynamic_cast the returned reference to a reference to your concrete Guts
    // class.
    const System::Guts& getSystemGuts() const {assert(guts); return *guts;}
    System::Guts&       updSystemGuts()       {assert(guts); return *guts;}

    // Put new *unowned* Guts into this *empty* handle and take over ownership.
    // If this handle is already in use, or if Guts is already owned this
    // routine will throw an exception.
    void adoptSystemGuts(System::Guts* g);

    explicit System(System::Guts* g) : guts(g) { }
    bool hasGuts() const {return guts!=0;}

private:
    // OBSOLETE
    // Part of our ongoing crusade to turn getN's into getNums for API consistency.
    int getNSubsystems() const {return getNumSubsystems();}
};


class SimTK_SimTKCOMMON_EXPORT DefaultSystemSubsystem : public Subsystem {
public:
    explicit DefaultSystemSubsystem(System& sys);
    void addEventHandler(ScheduledEventHandler* handler);
    void addEventHandler(TriggeredEventHandler* handler);
    void addEventReporter(ScheduledEventReporter* handler) const;
    void addEventReporter(TriggeredEventReporter* handler) const;
    EventId createEventId(SubsystemIndex subsys, const State& state) const;
    void findSubsystemEventIds
       (SubsystemIndex subsys, const State& state, 
        const Array_<EventId>& allEvents, 
        Array_<EventId>& eventsForSubsystem) const;
  // don't let doxygen see this private class
    class Guts;
private:
    const Guts& getGuts() const;
    Guts& updGuts();
};

inline void System::addEventHandler(ScheduledEventHandler* handler)
{   updDefaultSubsystem().addEventHandler(handler); }
inline void System::addEventHandler(TriggeredEventHandler* handler)
{   updDefaultSubsystem().addEventHandler(handler); }
inline void System::addEventReporter(ScheduledEventReporter* handler) const
{   getDefaultSubsystem().addEventReporter(handler); }
inline void System::addEventReporter(TriggeredEventReporter* handler) const
{   getDefaultSubsystem().addEventReporter(handler); }

inline System::operator const Subsystem&() const {return getDefaultSubsystem();}
inline System::operator Subsystem&() {return updDefaultSubsystem();}

inline static System::ProjectOptions operator|(System::ProjectOptions::Option  o1,    System::ProjectOptions::Option o2)    {return System::ProjectOptions(o1) |= o2;}
inline static System::ProjectOptions operator|(System::ProjectOptions          opts,  System::ProjectOptions::Option o)     {return opts |= o;}
inline static System::ProjectOptions operator|(System::ProjectOptions::Option  o,     System::ProjectOptions         opts)  {return opts |= o;}
inline static System::ProjectOptions operator&(System::ProjectOptions::Option  o1,    System::ProjectOptions::Option o2)    {return System::ProjectOptions(o1) &= o2;}
inline static System::ProjectOptions operator&(System::ProjectOptions          opts,  System::ProjectOptions::Option o)     {return opts &= o;}
inline static System::ProjectOptions operator&(System::ProjectOptions::Option  o,     System::ProjectOptions         opts)  {return opts &= o;}
inline static System::ProjectOptions operator~(System::ProjectOptions::Option  o)                                           {return ~System::ProjectOptions(o);}
inline static System::ProjectOptions operator-(System::ProjectOptions          opts,  System::ProjectOptions::Option o)     {return opts -= o;}
inline static System::ProjectOptions operator-(System::ProjectOptions          opts1, System::ProjectOptions         opts2) {return opts1 -= opts2;}

} // namespace SimTK

#endif // SimTK_SimTKCOMMON_SYSTEM_H_