C++ TimeStep类代码示例

void GnuplotExportModule::outputBoundaryCondition(PrescribedGradientBCNeumann &iBC, TimeStep *tStep)
{
    FloatArray stress;
    iBC.computeField(stress, tStep);

    printf("Mean stress computed in Gnuplot export module: "); stress.printYourself();

    double time = 0.0;

    TimeStep *ts = emodel->giveCurrentStep();
    if ( ts != NULL ) {
        time = ts->giveTargetTime();
    }

    int bcIndex = iBC.giveNumber();

    std :: stringstream strMeanStress;
    strMeanStress << "PrescribedGradientGnuplotMeanStress" << bcIndex << "Time" << time << ".dat";
    std :: string nameMeanStress = strMeanStress.str();
    std::vector<double> componentArray, stressArray;

    for(int i = 1; i <= stress.giveSize(); i++) {
        componentArray.push_back(i);
        stressArray.push_back(stress.at(i));
    }

    XFEMDebugTools::WriteArrayToGnuplot(nameMeanStress, componentArray, stressArray);

    // Homogenized strain
    
    FloatArray grad;
    iBC.giveGradientVoigt(grad);
    outputGradient(iBC.giveNumber(), *iBC.giveDomain(), grad, tStep);
}

void FE2FluidMaterialStatus :: setTimeStep(TimeStep *tStep)
{
    TimeStep *rveTStep = this->rve->giveCurrentStep(); // Should i create a new one if it is empty?
    rveTStep->setNumber( tStep->giveNumber() );
    rveTStep->setTime( tStep->giveTargetTime() );
    rveTStep->setTimeIncrement( tStep->giveTimeIncrement() );
}

void main()
{
	LOG("Burger Bash running on the Sifteo SDK\n");
	
	EventHandler eventHandler;
	eventHandler.init();
	
	// Plate cube initialize
	pCube.init(gNumCubes - 1);
	
	// Init all other cubes
	for (int i = 0; i < gNumCubes - 1; i++)
	{
		cubes[i].init(i);
	}
	
	// Game loop
	TimeStep ts;
    while (1) {
		for (int i = 0; i < gNumCubes; i++)
		{
			cubes[i].update();
		}
		
        System::paint();
        ts.next();
    }
    
    //AudioTracker::play(Music);
}

void GnuplotExportModule::outputXFEMGeometry(const std::vector< std::vector<FloatArray> > &iEnrItemPoints)
{
    double time = 0.0;

    TimeStep *ts = emodel->giveCurrentStep();
    if ( ts != NULL ) {
        time = ts->giveTargetTime();
    }

    std :: stringstream strCracks;
    strCracks << "CracksTime" << time << ".dat";
    std :: string nameCracks = strCracks.str();
    WritePointsToGnuplot(nameCracks, iEnrItemPoints);
}

contextIOResultType
PrimaryField :: restoreContext(DataStream *stream, ContextMode mode)
{
    int i, class_id;
    contextIOResultType iores;
    // read class header
    if ( !stream->read(& class_id, 1) ) {
        return CIO_IOERR;
    }

    if ( class_id != PrimaryFieldClass ) {
        return CIO_BADVERSION;
    }

    if ( !stream->read(& actualStepNumber, 1) ) {
        THROW_CIOERR(CIO_IOERR);
    }

    if ( !stream->read(& actualStepIndx, 1) ) {
        THROW_CIOERR(CIO_IOERR);
    }

    for ( i = 0; i <= nHistVectors; i++ ) {
        if ( ( iores = solutionVectors.at(i + 1)->restoreYourself(stream, mode) ) != CIO_OK ) {
            THROW_CIOERR(iores);
        }
    }

    int flag;
    TimeStep *iStep;
    for ( i = 0; i <= nHistVectors; i++ ) {
        if ( !stream->read(& flag, 1) ) {
            THROW_CIOERR(CIO_IOERR);
        }

        if ( flag ) {
            iStep = new TimeStep(emodel);
            if ( ( iores = iStep->restoreContext(stream, mode) ) != CIO_OK ) {
                THROW_CIOERR(iores);
            }
        } else {
            iStep = NULL;
        }

        solStepList.put(i + 1, iStep);
    }

    return CIO_OK;
}

double NLTransientTransportProblem :: giveUnknownComponent(ValueModeType mode, TimeStep *tStep, Domain *d, Dof *dof)
// returns unknown quantity like displacement, velocity of equation
// This function translates this request to numerical method language
{
    if ( this->requiresUnknownsDictionaryUpdate() ) {
        if ( mode == VM_Incremental ) { //get difference between current and previous time variable
            return dof->giveUnknowns()->at(0) - dof->giveUnknowns()->at(1);
        }
        int hash = this->giveUnknownDictHashIndx(mode, tStep);
        if ( dof->giveUnknowns()->includes(hash) ) {
            return dof->giveUnknowns()->at(hash);
        } else {
            OOFEM_ERROR("Dof unknowns dictionary does not contain unknown of value mode (%s)", __ValueModeTypeToString(mode) );
        }
    }

    double t = tStep->giveTargetTime();
    TimeStep *previousStep = this->givePreviousStep(), *currentStep = this->giveCurrentStep();

    if ( dof->__giveEquationNumber() == 0 ) {
        OOFEM_ERROR("invalid equation number on DoF %d", dof->giveDofID() );
    }

    if ( ( t >= previousStep->giveTargetTime() ) && ( t <= currentStep->giveTargetTime() ) ) {
        ///@todo Shouldn't it be enough to just run this?
        //UnknownsField->giveUnknownValue(dof, mode, currentStep);
        double rtdt = UnknownsField->giveUnknownValue(dof, VM_Total, currentStep);
        double rt   = UnknownsField->giveUnknownValue(dof, VM_Total, previousStep);
        double psi = ( t - previousStep->giveTargetTime() ) / currentStep->giveTimeIncrement();
        if ( mode == VM_Velocity ) {
            return ( rtdt - rt ) / currentStep->giveTimeIncrement();
        } else if ( mode == VM_Total ) {
            return psi * rtdt + ( 1. - psi ) * rt;
        } else if ( mode == VM_Incremental ) {
            if ( previousStep->isIcApply() ) {
                return 0;
            } else {
                return ( rtdt - rt );
            }
        } else {
            OOFEM_ERROR("Unknown mode %s is undefined for this problem", __ValueModeTypeToString(mode) );
        }
    } else {
        OOFEM_ERROR("time value %f not within bounds %f and %f", t, previousStep->giveTargetTime(), currentStep->giveTargetTime() );
    }

    return 0.; // to make compiler happy;
}

void PlayerComponent::OnUpdate(TimeStep timeStep)
{
	Application& app = Application::GetApplication();
	float speed = (float)(0.1 * timeStep.GetMills());
	vec3 pos = m_Transform->m_Transform.GetPosition();

	if (app.IsKeyPressed(GLFW_KEY_UP))
	{
		pos.y += speed;
	}
	else if (app.IsKeyPressed(GLFW_KEY_DOWN))
	{
		pos.y -= speed;
	}

	if (app.IsKeyPressed(GLFW_KEY_LEFT))
	{
		pos.x -= speed;
	}
	else if (app.IsKeyPressed(GLFW_KEY_RIGHT))
	{
		pos.x += speed;
	}

	m_Transform->m_Transform.SetPosition(pos);
}

void main()
{
    assetConfig.append(MainSlot, BetterflowAssets);
    loader.init();

    for (CubeID cube : CubeSet::connected())
    {
        onConnect(NULL, cube);
    }
    Events::cubeConnect.set(&onConnect);

    gVideo[0].setMode(BG0);

    MainMenu gameMenu(&gVideo[0]);

    int colour;
    Colormap colourList;
    colourList.setEGA();

    

    while(1)
    {
        colour = gameMenu.runMenu();

        gVideo[0].setMode(SOLID_MODE);

        gVideo[0].colormap[0].set(colourList[colour].get());
        System::paint();

        TimeStep ts;
        TimeTicker ticker = TimeTicker(1);
        int tickIncr;
        int tickCount  = 0;
        do
        {
            tickIncr = ticker.tick( ts.delta() );
            ts.next();
            tickCount += tickIncr;
        } while (tickCount < 3);

        gVideo[0].setMode(BG0);
        gVideo[0].bg0.erase(StripeTile);
    }

}

int  TimeIntegrationScheme :: requiresNewLhsAt (TimeStep* stepN)
   // Returns True if the time increment deltaT of stepN differs from that
   // of the previous step, else returns False. This is the default imple-
   // mentation of the method.
{
   TimeStep* pastStep ;
   int       n,answer ;

   n = stepN->giveNumber() ;
   if (n <= 1)
      answer = TRUE ;
   else {
      pastStep = new TimeStep(n-1,this) ;
      answer = (pastStep->giveTimeIncrement() != stepN->giveTimeIncrement());
      delete pastStep ;}

   return answer ;
}

void GnuplotExportModule::outputMesh(Domain &iDomain)
{
    std::vector< std::vector<FloatArray> > pointArray;

    if(iDomain.giveNumberOfSpatialDimensions() == 2) {
        // Write all element edges to gnuplot
        for ( auto &el : iDomain.giveElements() ) {
            int numEdges = el->giveNumberOfNodes();


            for ( int edgeIndex = 1; edgeIndex <= numEdges; edgeIndex++ ) {
                std::vector<FloatArray> points;

                IntArray bNodes;
                el->giveInterpolation()->boundaryGiveNodes(bNodes, edgeIndex);

                int niLoc = bNodes.at(1);
                const FloatArray &xS = *(el->giveNode(niLoc)->giveCoordinates() );
                points.push_back(xS);

                int njLoc = bNodes.at( bNodes.giveSize() );
                const FloatArray &xE = *(el->giveNode(njLoc)->giveCoordinates() );
                points.push_back(xE);

                pointArray.push_back(points);
            }

        }


        double time = 0.0;

        TimeStep *ts = emodel->giveCurrentStep();
        if ( ts != NULL ) {
            time = ts->giveTargetTime();
        }

        std :: stringstream strMesh;
        strMesh << "MeshGnuplotTime" << time << ".dat";
        std :: string nameMesh = strMesh.str();

        WritePointsToGnuplot(nameMesh, pointArray);
    }
}

void
StaticFracture :: solveYourself()
{
    MetaStep *activeMStep;
    FILE *out = this->giveOutputStream();
    this->timer.startTimer(EngngModelTimer :: EMTT_AnalysisTimer);
    this->giveNumberOfSlaveProblems();
    
    this->timer.startTimer(EngngModelTimer :: EMTT_SolutionStepTimer);
    this->timer.initTimer(EngngModelTimer :: EMTT_NetComputationalStepTimer);


    int numMetaSteps = this->giveNumberOfMetaSteps();
    for (int imstep = 1; imstep <= numMetaSteps; imstep++) { // don't know what will happen if we have several meta steps?
        activeMStep = this->giveMetaStep(imstep);

        int nTimeSteps = activeMStep->giveNumberOfSteps();
        for ( int tStepNum = 1; tStepNum <= nTimeSteps; tStepNum++ ) { //loop over time steps in opt analysis

            for ( int subProb = 1; subProb <= this->giveNumberOfSlaveProblems(); subProb++ ) {
        
                EngngModel *sp = this->giveSlaveProblem(subProb); 
                sp->solveYourself();

                //this->updateYourself( this->giveCurrentStep()); // not neccessary
            
                // optimization
                this->optimize( this->giveCurrentStep() );    

                // Resetting the time step number for each sp after each optimization time step
                TimeStep *tStep =  sp->giveCurrentStep();
                tStep->setNumber(tStepNum); 
                sp->giveExportModuleManager()->doOutput(tStep); // turn off export during regular analysis
            
                tStep->setNumber(0); // otherwise the anlysis wont restart at time 0
                }
        }
    }
    

}

void TimeSteps::clear()
{
    Tools::ClassManager *cm = Tools::ClassManager::getInstance();
    vector<TimeStep *>::iterator tsIterator = timeSteps.begin();

    while (tsIterator != timeSteps.end())
    {
        TimeStep *timeStep = *tsIterator;

        if (timeStep != NULL)
        {
            timeStep->clear();
            cm->deleteObject(timeStep->getID());

            timeStep = NULL;
        }

        ++tsIterator;
    }

    timeSteps.clear();
}

void GnuplotExportModule::outputXFEM(Crack &iCrack, TimeStep *tStep)
{
    const std::vector<GaussPoint*> &czGaussPoints = iCrack.giveCohesiveZoneGaussPoints();

    std::vector<double> arcLengthPositions, normalJumps, tangJumps, normalTractions;

    const BasicGeometry *bg = iCrack.giveGeometry();

    for( GaussPoint *gp: czGaussPoints ) {

        StructuralInterfaceMaterialStatus *matStat = dynamic_cast<StructuralInterfaceMaterialStatus*> ( gp->giveMaterialStatus() );
        if(matStat != NULL) {

            // Compute arc length position of the Gauss point
            const FloatArray &coord = (gp->giveGlobalCoordinates());
            double tangDist = 0.0, arcPos = 0.0;
            bg->computeTangentialSignDist(tangDist, coord, arcPos);
            arcLengthPositions.push_back(arcPos);

            // Compute displacement jump in normal and tangential direction
            // Local numbering: (tang_z, tang, normal)
            const FloatArray &jumpLoc = matStat->giveJump();

            double normalJump = jumpLoc.at(3);
            normalJumps.push_back(normalJump);


            tangJumps.push_back( jumpLoc.at(2) );


            const FloatArray &trac = matStat->giveFirstPKTraction();
            normalTractions.push_back(trac.at(3));
        }
    }



    Domain *domain = emodel->giveDomain(1);
    XfemManager *xMan = domain->giveXfemManager();
    if ( xMan != NULL ) {
        double time = 0.0;

        TimeStep *ts = emodel->giveCurrentStep();
        if ( ts != NULL ) {
            time = ts->giveTargetTime();
        }

        int eiIndex = iCrack.giveNumber();

        std :: stringstream strNormalJump;
        strNormalJump << "NormalJumpGnuplotEI" << eiIndex << "Time" << time << ".dat";
        std :: string nameNormalJump = strNormalJump.str();
        XFEMDebugTools::WriteArrayToGnuplot(nameNormalJump, arcLengthPositions, normalJumps);

        std :: stringstream strTangJump;
        strTangJump << "TangJumpGnuplotEI" << eiIndex << "Time" << time << ".dat";
        std :: string nameTangJump = strTangJump.str();
        XFEMDebugTools::WriteArrayToGnuplot(nameTangJump, arcLengthPositions, tangJumps);

        std :: stringstream strNormalTrac;
        strNormalTrac << "NormalTracGnuplotEI" << eiIndex << "Time" << time << ".dat";
        std :: string nameNormalTrac = strNormalTrac.str();
        XFEMDebugTools::WriteArrayToGnuplot(nameNormalTrac, arcLengthPositions, normalTractions);


        std::vector<FloatArray> matForcesStart, matForcesEnd;
        std::vector<double> radii;

        // Material forces
        for(double matForceRadius : mMatForceRadii) {

            radii.push_back(matForceRadius);

            EnrichmentFront *efStart = iCrack.giveEnrichmentFrontStart();
            const TipInfo &tipInfoStart = efStart->giveTipInfo();

            FloatArray matForceStart;
            mpMatForceEvaluator->computeMaterialForce(matForceStart, *domain, tipInfoStart, tStep, matForceRadius);

            if(matForceStart.giveSize() > 0) {
                matForcesStart.push_back(matForceStart);
            }
            else {
                matForcesStart.push_back({0.0,0.0});
            }


            EnrichmentFront *efEnd = iCrack.giveEnrichmentFrontEnd();
            const TipInfo &tipInfoEnd = efEnd->giveTipInfo();

            FloatArray matForceEnd;
            mpMatForceEvaluator->computeMaterialForce(matForceEnd, *domain, tipInfoEnd, tStep, matForceRadius);

            if(matForceEnd.giveSize() > 0) {
                matForcesEnd.push_back(matForceEnd);
            }
            else {
                matForcesEnd.push_back({0.0,0.0});
            }

//.........这里部分代码省略.........

//.........这里部分代码省略.........
    for(size_t i = 0; i < mTriangles.size(); i++) {
    	if( mTriangles[i].getArea() > triTol ) {
    		triToKeep.push_back(i);
    	}
    }

    int nPointsTot = nPoints * triToKeep.size();
    FloatArray coords_xi1, coords_xi2, weights;
    this->giveTriCoordsAndWeights(nPoints, coords_xi1, coords_xi2, weights);
    this->gaussPointArray = new GaussPoint * [ nPointsTot ];
    ////////////////////////////////////////////


    std :: vector< FloatArray >newGPCoord;

    double parentArea = this->elem->computeArea();

    // Loop over triangles
    for ( int i = 0; i < int( triToKeep.size() ); i++ ) {
        // TODO: Probably unnecessary to allocate here
        const FloatArray **coords = new const FloatArray * [ mTriangles [ triToKeep[i] ].giveNrVertices() ];
        // this we should put into the function before
        for ( int k = 0; k < mTriangles [ triToKeep[i] ].giveNrVertices(); k++ ) {
            coords [ k ] = new FloatArray( ( mTriangles [ triToKeep[i] ].giveVertex(k + 1) ) );
        }

        // Can not be used because it writes to the start of the array instead of appending.
        //		int nPointsTri = GaussIntegrationRule :: SetUpPointsOnTriangle(nPoints, mode);

        for ( int j = 0; j < nPoints; j++ ) {
            FloatArray global;
            GaussPoint * &gp = this->gaussPointArray [ pointsPassed ];

            FloatArray *coord = new FloatArray(2);
            coord->at(1) = coords_xi1.at(j + 1);
            coord->at(2) = coords_xi2.at(j + 1);
            gp = new GaussPoint(this, pointsPassed + 1, coord, weights.at(j + 1), mode);



            mTriInterp.local2global( global, * gp->giveCoordinates(),
                                     FEIVertexListGeometryWrapper(mTriangles [ triToKeep[i] ].giveNrVertices(), coords) );

            newGPCoord.push_back(global);


            FloatArray local;
            this->elem->computeLocalCoordinates(local, global);

            gp->setCoordinates(local);




            double refElArea = this->elem->giveParentElSize();

            gp->setWeight(2.0 * refElArea * gp->giveWeight() * mTriangles [ triToKeep[i] ].getArea() / parentArea); // update integration weight


            pointsPassed++;
        }


        for ( int k = 0; k < mTriangles [ triToKeep[i] ].giveNrVertices(); k++ ) {
            delete coords [ k ];
        }

        delete [] coords;
    }

#if PATCH_INT_DEBUG > 0

    double time = 0.0;

    Element *el = this->elem;
    if(el != NULL) {
    	Domain *dom = el->giveDomain();
    	if(dom != NULL) {
    		EngngModel *em = dom->giveEngngModel();
    		if(em != NULL) {
    			TimeStep *ts = em->giveCurrentStep();
    			if(ts != NULL) {
    				time = ts->giveTargetTime();
    			}
    		}
    	}
    }
    int elIndex = this->elem->giveGlobalNumber();
    std :: stringstream str;
    str << "GaussPointsTime" << time << "El" << elIndex << ".vtk";
    std :: string name = str.str();

    XFEMDebugTools :: WritePointsToVTK(name, newGPCoord);
#endif

    numberOfIntegrationPoints = pointsPassed;


    return numberOfIntegrationPoints;
}

void main()
{
    GameState state = StateResults;
	Random randomGen;
  	randomGen.seed();

    for (unsigned i = 0; i < arraysize(cubeVideo); i++)
    {
        InitCube(i, i);
    }

    EventHandler handler;
    handler.init();

    TimeStep ts;
    float Delay;
    bool isDone;

    while (1)
    {
        switch (state)
        {
            case StateInit:
				CurrentGame->init();
				Delay = 3;
				playSfx (CountSound);
				LOG ("Init\n");
				state = StateStart;
            break;

            case StateStart:
                Delay -= float(ts.delta());

                if (Delay <= 0)
                {
                    CurrentGame->start();
                    for (unsigned i = 0; i < gNumCubes; i++)
                    {
                        CurrentGameCube[i]->start();
                    }
                    LOG ("Start Done\n");
		            state = StatePlay;
                }
            break;

            case StatePlay:
				isDone = CurrentGame->update(ts.delta());

				if (isDone)
				{
					LOG ("Game Done\n");
					state = StateResults;
					Delay = 3;
					playSfx (CheersSound);
				}
				else
				{
					for (unsigned i = 0; i < gNumCubes; i++)
					{
						CurrentGameCube[i]->update(ts.delta());
					}
				}
            break;

            case StateResults:
                Delay -= float(ts.delta());

                if (Delay <= 0)
                {
					state = StateInit;
					LOG ("Results Done\n");
					
					int rand = randomGen.randint(0, 3);
					switch (rand)
					{
					    case 0:
						case 2:
						case 3:
						LOG ("Playing Flip It\n");
						CurrentGame = &flipItGame;
					    break;
						
					    case 1:
						LOG ("Playing Shake\n");
						CurrentGame = &shakeGame;
						break;
					} 
                }
            break;

        }

        System::paint();
        ts.next();
    }
}