SM_Beam.h

Go to the documentation of this file.

#ifndef SM_Beam_H
#define SM_Beam_H
 
//define the std::valarray object
#include <valarray>
//define the std::array object
#include <array>
 
//base class
#include "SM_Object.h"
 
//core run for class factoru
#include "CORE_Run.h"
 
//system class
#include "SM_System.h"
#include "SM_StratonovichNormalizedSystem.h"
#include "SM_StratonovichSystem.h"
#include "SM_ItoSystem.h"
 
//constant noise rate
#include "SM_ConstantNoiseRateFunction.h"
 
//stochastic data
#include "SM_StochasticFunction.h"
 
//network data
#include "SM_Network.h"
 
//Time stepper data
#include "SM_TimeStepper.h"
 
//operator interface
#include "SM_Operator.h"
 
class SM_Beam : public SM_Object {
    
    //attributes
private :
    
    
    //number of simulations
    tIndex mBeamSize;
    
    //mu array of size mBeamSize x nTimeStepsNumber x mParticlesNumber x mDimension
    std::valarray<tReal> mMu;
 
    //mu array of size mBeamSize x nTimeStepsNumber x (nOperatorsNumber+1) 
    std::valarray<tReal> mEs;
 
    //templated method called
    tBoolean mIsTemplatedSimulateVersion;
    
    //associations
    //============
    
    //stochastic function
    CORE_UniquePointer<SM_StochasticFunction> mStochasticFunction;
 
    //system
    CORE_UniquePointer<SM_System> mSystem;
 
    //deterministic system
    CORE_UniquePointer<SM_System> mDeterministicSystem;
    std::valarray<tReal> mDMu;
    std::valarray<tReal> mDEs;
 
protected:
    // CONSTRUCTORS 
    SM_Beam(void) {
        mBeamSize=1000;
        mIsTemplatedSimulateVersion=false;
    }
    
    // DESTRUCTORS 
    virtual ~SM_Beam(void) {
    }
    
public : 
    // CREATE class as a pointer
    
    inline static CORE_UniquePointer<SM_Beam> New() {
        return CORE_UniquePointer<SM_Beam>(new SM_Beam(),
                                           SM_Beam::Delete());
    }
 
 
    
    //MEMORY
    
    virtual tMemSize getMemorySize() const {
        return sizeof(*this)+getContentsMemorySize();
    }
 
    virtual tMemSize getContentsMemorySize() const {
        tMemSize mem=SM_Object::getContentsMemorySize();
        mem+=mMu.size()*sizeof(tReal);
        mem+=mEs.size()*sizeof(tReal);
        mem+=mDMu.size()*sizeof(tReal);
        mem+=mDEs.size()*sizeof(tReal);
        mem+=(mStochasticFunction.get()==null)?0:mStochasticFunction->getMemorySize();
        mem+=(mSystem.get()==null)?0:mSystem->getMemorySize();
        mem+=(mDeterministicSystem.get()==null)?0:mDeterministicSystem->getMemorySize();
        
        return mem;
    }
 
    
    //SET & GET methods
 
 
    inline void setIsTemplatedSimulateVersion(const tBoolean& v) {
        mIsTemplatedSimulateVersion=v;
    }
    
    inline void setBeamSize(const tIndex& n) {
        mBeamSize=n;
    }
    inline const tIndex& getBeamSize() const {
        return mBeamSize;
    }
    
    inline void setStochasticFunction(CORE_UniquePointer<SM_StochasticFunction>& f) {
        mStochasticFunction=std::move(f);
    }
    inline const SM_StochasticFunction&  getStochasticFunction() const {
        return *mStochasticFunction.get();
    }
    inline void setSystem(CORE_UniquePointer<SM_System>& sys) {
        mSystem=std::move(sys);
    }
    virtual const SM_System& getSystem() const {
        return *mSystem.get();
    }
    
    virtual SM_System& getSystem() {
        return *mSystem.get();
    }
    
    
    inline const std::valarray<tReal>& getMagneticMoment() const {
        return mMu;
    }
    inline std::valarray<tReal>& getMagneticMoment()  {
        return mMu;
    }
    inline const std::valarray<tReal>& getEnergies() const {
        return mEs;
    }
    inline std::valarray<tReal>& getEnergies()  {
        return mEs;
    }
    // MAIN Method
    
public:
    inline void discretize() {
        //get system
        const SM_System& system=getSystem();
        
        //get network
        const SM_Network& network=system.getNetwork();
        
        //get time stepper
        const SM_TimeStepper& stepper=system.getTimeStepper();
        
        //number of time steps
        const tIndex& nTimeSteps=stepper.getMaximumTimeStepsNumber();
        
        //the mu by simulation
        tIndex muSize=network.getDimension()*network.getParticlesNumber()*(nTimeSteps+1);
 
        //the E by simulation
        tIndex ESize=(system.getOperatorsNumber()+1)*nTimeSteps;
        
        //size of mu
        try { 
            mMu.resize(mBeamSize * muSize);
        } catch (std::exception& e) {
            CORE_Run::Out().printError("not enough memory to alocate "+std::to_string(mBeamSize*muSize*sizeof(tReal)/1000000)+" Mo of memory for mMu");
            CORE_Run::Out()<<e.what()<<"\n";
        }
        //size of E
        try { 
            mEs.resize(mBeamSize * ESize);
        } catch (std::exception& e) {
            CORE_Run::Out().printError("not enough memory to alocate "+std::to_string(mBeamSize*ESize*sizeof(tReal)/1000000)+" Mo of memory for mEs");
            CORE_Run::Out()<<e.what()<<"\n";
        }
 
        //discretize the system
        mSystem->discretize();
 
        //define the determinitic system
        mDeterministicSystem=SM_ItoSystem::New();//create a deterministic system
        mDeterministicSystem->getNetwork().copy(network);//same network
        mDeterministicSystem->getTimeStepper().copy(stepper);//same stepper
        mDeterministicSystem->getLandauLifschitzFunction().copy(system.getLandauLifschitzFunction());//same LL
        CORE_UniquePointer<SM_NoiseRateFunction> noise=SM_ConstantNoiseRateFunction::New();
        mDeterministicSystem->setNoiseRateFunction(noise);//no noise
        mDeterministicSystem->getNoiseRateFunction().setNoiseRate(0);
        mDeterministicSystem->copyOperators(system);//same operators
        mDeterministicSystem->setInitialMagneticMoment(system.getInitialMagneticMoment());//same mu0
        mDMu.resize(muSize);
        mDEs.resize(ESize);
        mDeterministicSystem->discretize();//discretize
        
    }
 
    inline const std::valarray<tReal>& getDeterministicEnergies() const {
        return mDEs;
    }
 
    inline tIndex runSimulations() {
        if (mDeterministicSystem.get()!=null) {
            mDeterministicSystem->makeRelaxation(null,
                                                 &mDMu[0],
                                                 &mDEs[0]);
        }
        if (mIsTemplatedSimulateVersion) return templatedSimulate();
        return simulate();
    }
private:
    tIndex  simulate();
 
    inline tIndex  templatedSimulate() {
        //get the system
        SM_System &system=getSystem();
        
        //get the noise of the system
        const SM_NoiseRateFunction *nF=&system.getNoiseRateFunction();
 
        //cast to templated classes
        
        SM_ItoSystem* sys1=dynamic_cast<SM_ItoSystem*>(&getSystem());
        if (sys1!=null) {
            if (nF->isInstanceOf<SM_ConstantNoiseRateFunction>()) {
                return templatedSimulate(*sys1,*dynamic_cast<const SM_ConstantNoiseRateFunction*>(nF));
            } else if (nF->isInstanceOf<SM_InverseNoiseRateFunction>()) { 
                return templatedSimulate(*sys1,*dynamic_cast<const SM_InverseNoiseRateFunction*>(nF));
            } 
        }
        SM_StratonovichSystem* sys2=dynamic_cast<SM_StratonovichSystem*>(&getSystem());
        if (sys2!=null) {
            if (nF->isInstanceOf<SM_ConstantNoiseRateFunction>()) {
                return templatedSimulate(*sys2,*dynamic_cast<const SM_ConstantNoiseRateFunction*>(nF));
            } else if (nF->isInstanceOf<SM_InverseNoiseRateFunction>()) { 
                return templatedSimulate(*sys2,*dynamic_cast<const SM_InverseNoiseRateFunction*>(nF));
            } 
        }
        SM_StratonovichNormalizedSystem* sys3=dynamic_cast<SM_StratonovichNormalizedSystem*>(&getSystem());
        if (sys3!=null) {
            if (nF->isInstanceOf<SM_ConstantNoiseRateFunction>()) {
                return templatedSimulate(*sys3,*dynamic_cast<const SM_ConstantNoiseRateFunction*>(nF));
            } else if (nF->isInstanceOf<SM_InverseNoiseRateFunction>()) { 
                return templatedSimulate(*sys3,*dynamic_cast<const SM_InverseNoiseRateFunction*>(nF));
            } 
        }
 
        //default simulation
        return simulate();
    }
 
public:
    virtual tString toString() const override {
        std::stringstream ret;
        //print id
        ret<<SM_Object::toString()<<"\n";
        //print beam size
        ret<<"beam size: "<<mBeamSize<<"\n";
       
        ret<<"stochastic Function:";
        if (mStochasticFunction.get()!=null) 
            ret<<mStochasticFunction->toString()<<"\n";
        else
            ret<<"null\n";
        ret<<"system:";
        if (mSystem.get()!=null) 
            ret<<mSystem->toString()<<"\n";
        else
            ret<<"null\n";
        
        return ret.str();
    }
private:
    template<class SystemImpl,class NoiseImpl>
    tIndex  templatedSimulate(SystemImpl& system,const NoiseImpl& noise);
    
public:
    void computeL2Convergence(std::valarray<tReal>& L2conv) const;
    
    
};
 
#include "SM_Beam.hpp"
 
#endif