Commit 385956a8 authored by Ulrich Kemloh's avatar Ulrich Kemloh
Browse files

relocation building update in the simulation class

refactoring the base class for operative models
parent 01f2c6a9
......@@ -156,6 +156,7 @@ set ( source_files
math/GPU_GCFMModel.cpp
math/GCFMModel.cpp
math/GompertzModel.cpp
math/OperationalModel.cpp
mpi/LCGrid.cpp
......@@ -270,6 +271,7 @@ set ( header_files
math/GPU_GCFMModel.cpp
math/GCFMModel.h
math/GompertzModel.h
math/OperationalModel.h
poly2tri/poly2tri.h
poly2tri/common/shapes.h
......
......@@ -673,7 +673,7 @@ void TrajectoriesJPSV05::WriteGeometry(Building* building)
}
//write the grid
embed_geometry.append(building->GetGrid()->ToXML());
//embed_geometry.append(building->GetGrid()->ToXML());
embed_geometry.append("\t</geometry>\n");
_outputHandler->Write(embed_geometry);
......
......@@ -33,6 +33,15 @@
#include "math/GCFMModel.h"
#include "math/GompertzModel.h"
#include "geometry/Goal.h"
#ifdef _OPENMP
#include <omp.h>
#else
#define omp_get_thread_num() 0
#define omp_get_max_threads() 1
#endif
using namespace std;
......@@ -47,7 +56,7 @@ Simulation::Simulation()
_building = NULL;
_distribution = NULL;
_direction = NULL;
_model = NULL;
_operationalModel = NULL;
_solver = NULL;
_iod = new IODispatcher();
_fps=1;
......@@ -62,7 +71,7 @@ Simulation::~Simulation()
delete _building;
delete _distribution;
delete _direction;
delete _model;
delete _operationalModel;
delete _solver;
delete _iod;
delete _em;
......@@ -73,12 +82,6 @@ int Simulation::GetPedsNumber() const
return _nPeds;
}
Building * Simulation::GetBuilding() const
{
return _building;
}
void Simulation::InitArgs(ArgumentParser* args)
{
char tmp[CLENGTH];
......@@ -207,14 +210,14 @@ void Simulation::InitArgs(ArgumentParser* args)
break;
}
int model = args->GetModel();
if(model == 1)
if(model == MODEL_GFCM)
{ //GCFM
//if(args->GetHPCFlag()==0){
_model = new GCFMModel(_direction, args->GetNuPed(), args->GetNuWall(), args->GetDistEffMaxPed(),
_operationalModel = new GCFMModel(_direction, args->GetNuPed(), args->GetNuWall(), args->GetDistEffMaxPed(),
args->GetDistEffMaxWall(), args->GetIntPWidthPed(), args->GetIntPWidthWall(),
args->GetMaxFPed(), args->GetMaxFWall());
s.append("\tModel: GCFMModel\n");
s.append(_model->writeParameter());
s.append(_operationalModel->GetDescription());
//}
//else if(args->GetHPCFlag()==2){
// _model = new GPU_acc_GCFMModel(_direction, args->GetNuPed(), args->GetNuWall(), args->GetDistEffMaxPed(),
......@@ -231,28 +234,29 @@ void Simulation::InitArgs(ArgumentParser* args)
//s.append(_model->writeParameter());
//}
}
else if (model == 2)
else if (model == MODEL_GOMPERTZ)
{ //Gompertz
_model = new GompertzModel(_direction, args->GetNuPed(), args->GetaPed(), args->GetbPed(), args->GetcPed(),
_operationalModel = new GompertzModel(_direction, args->GetNuPed(), args->GetaPed(), args->GetbPed(), args->GetcPed(),
args->GetNuWall(), args->GetaWall(), args->GetbWall(), args->GetcWall() );
s.append("\tModel: GompertzModel\n");
s.append(_model->writeParameter());
s.append(_operationalModel->GetDescription());
}
// ODE solver
int solver = args->GetSolver();
sprintf(tmp, "\tODE Solver: %d\n", solver);
s.append(tmp);
switch (solver) {
case 1:
_solver = new EulerSolver(_model);
break;
//switch (solver) {
//case 1:
//_solver = new EulerSolver(_model);
//break;
//case 2:
// _solver = new VelocityVerletSolver(_model);
// break;
//case 3:
// _solver = new LeapfrogSolver(_model);
// break;
}
//}
sprintf(tmp, "\tnCPU: %d\n", args->GetMaxOpenMPThreads());
s.append(tmp);
_tmax = args->GetTmax();
......@@ -379,7 +383,7 @@ void Simulation::InitArgs(ArgumentParser* args)
//perform customs initialisation, like computing the phi for the gcfm
//this should be called after the routing engine has been initialised
// because a direction is needed for this initialisation.
_model->Init(_building);
_operationalModel->Init(_building);
//other initializations
const vector< Pedestrian* >& allPeds = _building->GetAllPedestrians();
......@@ -431,7 +435,7 @@ int Simulation::RunSimulation()
_iod->WriteFrame(0,_building);
//first initialisation needed by the linked-cells
Update();
UpdateRoutesAndLocations();
//needed to control the execution time PART 1
//in the case you want to run in no faster than realtime
......@@ -439,14 +443,31 @@ int Simulation::RunSimulation()
//time(&starttime);
// main program loop
for (t = 0; t < _tmax && _nPeds > 0; ++frameNr) {
for (t = 0; t < _tmax && _nPeds > 0; ++frameNr)
{
t = 0 + (frameNr - 1) * _deltaT;
// solve ODE
_solver->solveODE(t, t + _deltaT, _building);
// update and check if pedestrians change rooms
Update();
// update the positions
_operationalModel->ComputeNextTimeStep(t,_deltaT,_building);
// update the routes and locations
//Update();
UpdateRoutesAndLocations();
//update the events
_em->Update_Events(t,_deltaT);
// trajectories output
//other updates
//someone might have leave the building
_nPeds=_building->GetAllPedestrians().size();
//update the linked cells
_building->UpdateGrid();
// update the global time
Pedestrian::SetGlobalTime(t);
// write the trajectories
if (frameNr % writeInterval == 0) {
_iod->WriteFrame(frameNr / writeInterval, _building);
}
......@@ -486,18 +507,106 @@ int Simulation::RunSimulation()
}
// TODO: make the building class more independent by moving the update routing here.
void Simulation::Update()
void Simulation::UpdateRoutesAndLocations()
{
//_building->Update();
_building->Update();
//someone might have leave the building
_nPeds=_building->GetAllPedestrians().size();
// update the global time
Pedestrian::SetGlobalTime(Pedestrian::GetGlobalTime()+_deltaT);
//update the cells position
_building->UpdateGrid();
//pedestrians to be deleted
//you should better create this in the constructor and allocate it once.
vector<Pedestrian*> pedsToRemove;
pedsToRemove.reserve(100); //just reserve some space
// collect all pedestrians in the simulation.
const vector< Pedestrian* >& allPeds = _building->GetAllPedestrians();
const map<int, Goal*>& goals=_building->GetAllGoals();
unsigned int nSize = allPeds.size();
int nThreads = omp_get_max_threads();
int partSize = nSize / nThreads;
#pragma omp parallel default(shared) num_threads(nThreads)
{
const int threadID = omp_get_thread_num();
int start = threadID * partSize;
int end = (threadID + 1) * partSize - 1;
if ((threadID == nThreads - 1))
end = nSize - 1;
for (int p = start; p <= end; ++p)
{
Pedestrian* ped = allPeds[p];
Room* room = _building->GetRoom(ped->GetRoomID());
SubRoom* sub = room->GetSubRoom(ped->GetSubRoomID());
//set the new room if needed
if ((ped->GetFinalDestination() == FINAL_DEST_OUT)
&& (room->GetCaption() == "outside"))
{
#pragma omp critical
pedsToRemove.push_back(ped);
}
else if ((ped->GetFinalDestination() != FINAL_DEST_OUT)
&& (goals.at(ped->GetFinalDestination())->Contains(ped->GetPos())))
{
#pragma omp critical
pedsToRemove.push_back(ped);
}
// reposition in the case the pedestrians "accidently left the room" not via the intended exit.
// That may happen if the forces are too high for instance
// the ped is removed from the simulation, if it could not be reassigned
else if (!sub->IsInSubRoom(ped))
{
bool assigned = false;
const std::vector<Room*>& allRooms=_building->GetAllRooms();
for (Room* room : allRooms)
{
const vector<SubRoom*>& allSubs=room->GetAllSubRooms();
for(SubRoom* sub : allSubs)
{
SubRoom* old_sub= allRooms[ped->GetRoomID()]->GetSubRoom(ped->GetSubRoomID());
if ((sub->IsInSubRoom(ped->GetPos())) && (sub->IsDirectlyConnectedWith(old_sub))) {
ped->SetRoomID(room->GetID(), room->GetCaption());
ped->SetSubRoomID(sub->GetSubRoomID());
ped->ClearMentalMap(); // reset the destination
//ped->FindRoute();
assigned = true;
break;
}
}
if (assigned == true)
break; // stop the loop
}
if (assigned == false)
{
#pragma omp critical
pedsToRemove.push_back(ped);
}
}
//finally actualize the route
if (ped->FindRoute() == -1) {
//a destination could not be found for that pedestrian
Log->Write("\tINFO: \tCould not found a route for pedestrian %d",ped->GetID());
#pragma omp critical
pedsToRemove.push_back(ped);
}
}
}
// remove the pedestrians that have left the building
for (unsigned int p = 0; p < pedsToRemove.size(); p++)
{
_building->DeletePedestrian(pedsToRemove[p]);
}
// temporary fix for the safest path router
// if (dynamic_cast<SafestPathRouter*>(_routingEngine->GetRouter(1)))
// {
// SafestPathRouter* spr = dynamic_cast<SafestPathRouter*>(_routingEngine->GetRouter(1));
// spr->ComputeAndUpdateDestinations(_allPedestians);
// }
}
void Simulation::PrintStatistics()
......
......@@ -38,7 +38,7 @@
#include "geometry/Building.h"
#include "IO/OutputHandler.h"
#include "IO/IODispatcher.h"
#include "math/ForceModel.h"
#include "math/OperationalModel.h"
#include "math/ODESolver.h"
#include "routing/GlobalRouter.h"
#include "routing/QuickestPathRouter.h"
......@@ -74,7 +74,7 @@ private:
/// door crossing strategy for the pedestrians
DirectionStrategy* _direction;
/// Force model to use
ForceModel* _model;
OperationalModel* _operationalModel;
/// differential equation solver
ODESolver* _solver;
/// writing the trajectories to file
......@@ -97,26 +97,22 @@ public:
*/
int GetPedsNumber() const;
/**
* @return the building object containing all geometry elements
*/
Building* GetBuilding() const;
/**
* Read parameters from the argument parser class.
*/
void InitArgs(ArgumentParser *args);
/**
*
* @return the total simulated/evacuation time
*/
int RunSimulation();
/**
* Update the pedestrians states: positions, velocity, route
* Updathe route of the pedestrians and reassign rooms, in the case a room change happens
*/
void Update();
void UpdateRoutesAndLocations();
//void Update(double &b, double &p, double &t, double &g);
/**
......
......@@ -71,7 +71,6 @@
#define OBSTRUCTION 4
// Lenght of array
// Length of array
#define CLENGTH 1000
......@@ -113,7 +112,11 @@ enum RoutingStrategy {
enum OperativModels {
MODEL_GFCM=1,
MODEL_GOMPERTZ
MODEL_GOMPERTZ,
// MODEL_ORCA,
// MODEL_CFM,
// MODEL_VELO
// MODEL_GNM
};
//global functions for convenience
......
......@@ -370,9 +370,10 @@ void Building::LoadGeometry(const std::string &geometryfile)
Log->Write("ERROR:\tOnly the unit m (meters) is supported. \n\tYou supplied [%s]",xRootNode->Attribute("unit"));
exit(EXIT_FAILURE);
}
double version = xmltof(xRootNode->Attribute("version"), -1);
if (version != 0.5) { // @todo version number is hard coded
if (version != std::stod(JPS_VERSION)) {
Log->Write(" \tWrong geometry version!");
Log->Write(" \tOnly version >= %s supported",JPS_VERSION);
Log->Write(" \tPlease update the version of your geometry file to %s",JPS_VERSION);
......@@ -864,105 +865,8 @@ void Building::SanityCheck()
Log->Write("INFO: \t...Done!!!\n");
}
#ifdef _SIMULATOR
void Building::Update()
{
//pedestrians to be deleted
//you should better create this in the constructor and allocate it once.
vector<Pedestrian*> pedsToReposition;
pedsToReposition.reserve(100);
vector<Pedestrian*> pedsToRemove;
pedsToRemove.reserve(100);
for(unsigned int p=0;p<_allPedestians.size();p++)
{
Pedestrian* ped = _allPedestians[p];
Room* room = GetRoom(ped->GetRoomID());
SubRoom* sub = room->GetSubRoom(ped->GetSubRoomID());
//set the new room if needed
if ((ped->GetFinalDestination() == FINAL_DEST_OUT)
&& (GetRoom(ped->GetRoomID())->GetCaption() == "outside"))
{
pedsToRemove.push_back(ped);
} else if ((ped->GetFinalDestination() != FINAL_DEST_OUT)
&& (_goals[ped->GetFinalDestination()]->Contains(
ped->GetPos()))) {
pedsToRemove.push_back(ped);
} else if (!sub->IsInSubRoom(ped)) {
pedsToReposition.push_back(ped);
}
}
// reset that pedestrians who left their room not via the intended exit
for (unsigned int p = 0; p < pedsToReposition.size(); p++) {
Pedestrian* ped = pedsToReposition[p];
bool assigned = false;
for (int i = 0; i < GetNumberOfRooms(); i++) {
Room* room = GetRoom(i);
//if(room->GetCaption()=="outside") continue;
for (int j = 0; j < room->GetNumberOfSubRooms(); j++) {
SubRoom* sub = room->GetSubRoom(j);
SubRoom* old_sub= _rooms[ped->GetRoomID()]->GetSubRoom(ped->GetSubRoomID());
if ((sub->IsInSubRoom(ped->GetPos())) && (sub->IsDirectlyConnectedWith(old_sub))) {
ped->SetRoomID(room->GetID(), room->GetCaption());
ped->SetSubRoomID(sub->GetSubRoomID());
ped->ClearMentalMap(); // reset the destination
//ped->FindRoute();
assigned = true;
break;
}
}
if (assigned == true)
break; // stop the loop
}
if (assigned == false) {
pedsToRemove.push_back(ped);
}
}
// remove the pedestrians that have left the facility
for (unsigned int p = 0; p < pedsToRemove.size(); p++)
{
DeletePedestrian(pedsToRemove[p]);
}
// find the new goals, the parallel way
//FIXME temporary fix for the safest path router
if (dynamic_cast<SafestPathRouter*>(_routingEngine->GetRouter(1)))
{
SafestPathRouter* spr = dynamic_cast<SafestPathRouter*>(_routingEngine->GetRouter(1));
spr->ComputeAndUpdateDestinations(_allPedestians);
}
else
{
unsigned int nSize = _allPedestians.size();
int nThreads = omp_get_max_threads();
int partSize = nSize / nThreads;
#pragma omp parallel default(shared) num_threads(nThreads)
{
const int threadID = omp_get_thread_num();
int start = threadID * partSize;
int end = (threadID + 1) * partSize - 1;
if ((threadID == nThreads - 1))
end = nSize - 1;
for (int p = start; p <= end; ++p) {
if (_allPedestians[p]->FindRoute() == -1) {
//a destination could not be found for that pedestrian
Log->Write("\tINFO: \tCould not found a route for pedestrian %d",_allPedestians[p]->GetID());
DeletePedestrian(_allPedestians[p]);
//exit(EXIT_FAILURE);
}
}
}
}
}
#ifdef _SIMULATOR
void Building::UpdateGrid()
......
......@@ -169,7 +169,6 @@ public:
//void InitRoomsAndSubroomsMap();
void InitSavePedPathway(const std::string &filename);
void AddRoom(Room* room);
void Update();
void UpdateGrid();
void AddSurroundingRoom(); // add a final room (outside or world), that encompasses the complete geometry
......
......@@ -40,13 +40,6 @@
class Building;
/**
* @date Fri Apr 18 16:40:39 2014
*
* @brief The operative model. Definition of several force-based models
* for ped pedestrians dynamics
*
*/
class ForceModel {
public:
......@@ -55,7 +48,7 @@ public:
virtual ~ForceModel();
/**
* Solve the differential equations and update the positions and veloities
* Solve the differential equations and update the positions and velocities
* @param t the actual time
* @param tp the next timestep
* @param building the geometry object
......
......@@ -68,6 +68,7 @@ GCFMModel::~GCFMModel(void)
}
void GCFMModel::Init (Building* building) const
{
const vector< Pedestrian* >& allPeds = building->GetAllPedestrians();
......@@ -103,6 +104,132 @@ void GCFMModel::Init (Building* building) const
ped->SetEllipse(E);
}
}
void GCFMModel::ComputeNextTimeStep(double current, double deltaT, Building* building) const
{
double delta = 0.5;
//double h = tip1 - time;
//double h=deltaT;
// collect all pedestrians in the simulation.
const vector< Pedestrian* >& allPeds = building->GetAllPedestrians();
unsigned int nSize = allPeds.size();
int nThreads = omp_get_max_threads();
int partSize = nSize / nThreads;
int debugPed = -14;//10;
//building->GetGrid()->HighlightNeighborhood(debugPed, building);
#pragma omp parallel default(shared) num_threads(nThreads)
{
vector< Point > result_acc = vector<Point > ();
result_acc.reserve(2200);
const int threadID = omp_get_thread_num();
int start = threadID*partSize;
int end = (threadID + 1) * partSize - 1;
if ((threadID == nThreads - 1)) end = nSize - 1;
for (int p = start; p <= end; ++p) {
Pedestrian* ped = allPeds[p];
Room* room = building->GetRoom(ped->GetRoomID());
SubRoom* subroom = room->GetSubRoom(ped->GetSubRoomID());
//if(subroom->GetType()=="cellular") continue;
double normVi = ped->GetV().ScalarP(ped->GetV());
double tmp = (ped->GetV0Norm() + delta) * (ped->GetV0Norm() + delta);
if (normVi > tmp && ped->GetV0Norm() > 0) {
fprintf(stderr, "GCFMModel::calculateForce() WARNING: actual velocity (%f) of iped %d "
"is bigger than desired velocity (%f) at time: %fs\n",
sqrt(normVi), ped->GetID(), ped->GetV0Norm(), current);
// remove the pedestrian and abort
building->DeletePedestrian(ped);
Log->Write("\tERROR: one ped was removed due to high velocity");
exit(EXIT_FAILURE);
}
Point F_rep;
vector<Pedestrian*> neighbours;
building->GetGrid()->GetNeighbourhood(ped,neighbours);
//if(ped->GetID()==61) building->GetGrid()->HighlightNeighborhood(ped,building);
int nSize=neighbours.size();
for (int i = 0; i < nSize; i++) {
Pedestrian* ped1 = neighbours[i];
Point p1 = ped->GetPos();
Point p2 = ped1->GetPos();
bool isVisible = building->IsVisible(p1, p2, false);
if (!isVisible)
continue;
// if(debugPed == ped->GetID())
// {
// fprintf(stdout, "t=%f %f %f %f %f %d %d %d\n", time, p1.GetX(), p1.GetY(), p2.GetX(), p2.GetY(), isVisible, ped->GetID(), ped1->GetID());
// }
//if they are in the same subroom
if (ped->GetUniqueRoomID() == ped1->GetUniqueRoomID()) {
F_rep = F_rep + ForceRepPed(ped, ped1);
} else {
// or in neighbour subrooms
SubRoom* sb2=building->GetRoom(ped1->GetRoomID())->GetSubRoom(ped1->GetSubRoomID());
if(subroom->IsDirectlyConnectedWith(sb2)) {
F_rep = F_rep + ForceRepPed(ped, ped1);
}
}
}//for peds
//repulsive forces to the walls and transitions that are not my target
Point repwall = ForceRepRoom(allPeds[p], subroom);
Point fd = ForceDriv(ped, room);
// Point acc = (ForceDriv(ped, room) + F_rep + repwall) / ped->GetMass();
Point acc = (fd + F_rep + repwall) / ped->GetMass();
// if(ped->GetID() == debugPed){
// printf("%f %f %f %f %f %f\n", fd.GetX(), fd.GetY(), F_rep.GetX(), F_rep.GetY(), repwall.GetX(), repwall.GetY());
// }