Line data Source code
1 : /*******************************************************************************
2 : * *
3 : * Cosmos:(C)oncept et (O)utils (S)tatistique pour les (Mo)deles *
4 : * (S)tochastiques *
5 : * *
6 : * Copyright (C) 2009-2012 LSV & LACL *
7 : * Authors: Paolo Ballarini, Benoit Barbot & Hilal Djafri *
8 : * Website: http://www.lsv.ens-cachan.fr/Software/cosmos *
9 : * *
10 : * This program is free software; you can redistribute it and/or modify *
11 : * it under the terms of the GNU General Public License as published by *
12 : * the Free Software Foundation; either version 3 of the License, or *
13 : * (at your option) any later version. *
14 : * *
15 : * This program is distributed in the hope that it will be useful, *
16 : * but WITHOUT ANY WARRANTY; without even the implied warranty of *
17 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
18 : * GNU General Public License for more details. *
19 : * *
20 : * You should have received a copy of the GNU General Public License along *
21 : * with this program; if not, write to the Free Software Foundation, Inc., *
22 : * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. *
23 : *******************************************************************************
24 : */
25 :
26 : #include <iostream>
27 : #include <set>
28 : #include <unistd.h>
29 : #include <math.h>
30 : #include <float.h>
31 : #include <time.h>
32 : #include <ratio>
33 : #include <chrono>
34 : #include <iomanip>
35 :
36 : #include "Simulator.hpp"
37 : #include "marking.hpp"
38 :
39 : using namespace std;
40 :
41 : /**
42 : * Constructor for the Simulator initialize the event queue
43 : * but don't fill it.
44 : */
45 : template <class S, class EQT,class DEDS>
46 26 : SimulatorBase<S,EQT,DEDS>::SimulatorBase(DEDS& spn,LHA_orig<typename DEDS::stateType>& automate):N(spn),A(automate){
47 26 : EQ = new EQT(N.getNbTransition()); //initialization of the event queue
48 26 : logResult=false;
49 26 : sampleTrace = 0.0;
50 26 : Result.quantR.resize(A.FormulaVal.size());
51 26 : Result.qualR.resize(A.FormulaValQual.size());
52 26 : BatchSize = 1000;
53 26 : minInteractiveTime = 0.0;
54 26 : waitForTransition = -1;
55 26 : }
56 :
57 :
58 : /*Simulator::Simulator():verbose(0)){
59 : EQ = new EventsQueue(N); //initialization of the event queue
60 : logResult=false;
61 : sampleTrace = 0.0;
62 : Result.second.resize(A.FormulaVal.size());
63 : BatchSize = 1000;
64 : minInteractiveTime = 0.0;
65 : }*/
66 :
67 : template <class S, class EQT,class DEDS>
68 26 : SimulatorBase<S,EQT,DEDS>::~SimulatorBase() {
69 26 : delete EQ;
70 26 : }
71 :
72 :
73 : template <class S, class EQT,class DEDS>
74 26 : void SimulatorBase<S,EQT,DEDS>::logValue(const string path){
75 26 : logvalue.open(path,fstream::out);
76 26 : logvalue.precision(15);
77 26 : logResult=true;
78 26 : }
79 :
80 : template <class S, class EQT,class DEDS>
81 :
82 26 : void SimulatorBase<S,EQT,DEDS>::logTrace(const string path,double sample){
83 26 : if(!logtrace.is_open()){
84 26 : sampleTrace = sample;
85 26 : logtrace.open(path,fstream::out);
86 : //logtrace << "# sampling at:" << sample << endl;
87 26 : logtrace.precision(15);
88 26 : logtrace << "Time ";
89 26 : N.printHeader(logtrace);
90 26 : A.printHeader(logtrace);
91 26 : logtrace << endl;
92 : }
93 26 : }
94 :
95 : template <class S, class EQT,class DEDS>
96 104358764 : void SimulatorBase<S,EQT,DEDS>::printLog(double eTime,size_t t){
97 104358764 : if(logtrace.is_open())
98 272 : if((A.CurrentTime - lastSampled) >= sampleTrace){
99 272 : lastSampled = A.CurrentTime;
100 272 : logtrace <<setw(9)<<left<<setprecision(8)<< A.CurrentTime << " ";
101 272 : logtrace << right;
102 272 : N.printState(logtrace,eTime);
103 272 : A.printState(logtrace);
104 272 : if(t!=string::npos)logtrace << " ->" <<N.getTransitionLabel(t);
105 272 : logtrace << endl;
106 : }
107 104358764 : }
108 :
109 : template <class S, class EQT,class DEDS>
110 52068230 : void SimulatorBase<S,EQT,DEDS>::printLog(double eTime){
111 52068230 : printLog(eTime, string::npos);
112 52068230 : }
113 :
114 : template <class S, class EQT,class DEDS>
115 26 : void SimulatorBase<S,EQT,DEDS>::SetBatchSize(const size_t RI) {
116 26 : BatchSize = RI;
117 26 : }
118 :
119 : /**
120 : * Reset the SPN, The LHA and the Event Queue to the initial state.
121 : */
122 : template <class S, class EQT,class DEDS>
123 216102 : void SimulatorBase<S,EQT,DEDS>::reset() {
124 216102 : timeGen::reset();
125 216102 : N.reset();
126 216102 : A.reset(N.getState());
127 216102 : EQ->reset();
128 216102 : }
129 :
130 : /**
131 : * This function is called when an accepting state is reached in
132 : * the automaton. It update the automaton variable before updating the
133 : * Hasl formula.
134 : */
135 : template <class S, class EQT,class DEDS>
136 143094 : void SimulatorBase<S,EQT,DEDS>::returnResultTrue(){
137 143094 : A.getFinalValues(N.getState(),Result.quantR,Result.qualR);
138 143094 : Result.accept = true;
139 143094 : }
140 :
141 : /**
142 : * This function update the likelihood in the rare event context.
143 : * Do nothing when not in rare event context.
144 : * @param i Number of the transition of the SPN
145 : */
146 : template <class S, class EQT,class DEDS>
147 51652480 : void SimulatorBase<S,EQT,DEDS>::updateLikelihood(size_t){
148 51652480 : return;
149 : }
150 :
151 : /**
152 : * This function return true if the transition is the sink transition.
153 : * always return true when not in rare event context.
154 : * @param i Number of the transition of the SPN
155 : */
156 : template <class S, class EQT,class DEDS>
157 51652480 : bool SimulatorBase<S,EQT,DEDS>::transitionSink(size_t ){
158 51652480 : return false;
159 : }
160 :
161 :
162 : /**
163 : * Simulate one step of simulation
164 : * @return true if the simulation did not reach an accepting are refusing state.
165 : */
166 : template <class S, class EQT,class DEDS>
167 52384327 : bool SimulatorBase<S,EQT,DEDS>::SimulateOneStep(){
168 :
169 52384327 : AutEdge AE = A.GetEnabled_A_Edges( N.getState());
170 :
171 : //If there is no enabled transition in the Petri net
172 : //try to reach an accepting state by using autonomous edge of
173 : //the automata refuse the simulation otherwise.
174 52384327 : if (EQ->isEmpty()) {
175 70000 : while (AE.Index>-1) {
176 70000 : if(P.verbose>3){
177 0 : cerr << "Autonomous transition with Empty Queue:";
178 0 : cerr << AE.Index << endl;
179 0 : A.printState(cerr);
180 0 : cerr << endl;
181 : }
182 70000 : A.updateLHA( AE.FiringTime - A.CurrentTime, N.getState() );
183 70000 : A.fireAutonomous(AE.Index,N.getState());
184 70000 : if (A.isFinal()) {
185 70000 : static_cast<S*>(this)->returnResultTrue();
186 70000 : return false;
187 0 : } else AE = A.GetEnabled_A_Edges( N.getState());
188 : }
189 0 : Result.accept=false;
190 0 : return false;
191 : } else {
192 : //Take the first event in the queue
193 52314327 : const Event &E1 = EQ->InPosition(0);
194 :
195 : //If this transition is the sink transition refuse the simulation
196 : //Only usefull for Rare Event handling.
197 52314327 : if(static_cast<S*>(this)->transitionSink(E1.transition)){
198 653 : if(P.verbose>3)cerr << "\033[1;33mFiring:\033[0m" << "Transition Sink\n";
199 653 : Result.accept=false;
200 653 : return false;
201 : }
202 :
203 : //Hook for rare event simulation
204 52313674 : static_cast<S*>(this)->updateLikelihood(E1.transition);
205 :
206 : //Take all autonomous edge in the automata before the fire time
207 : //of the transition of the Petri net.
208 52613144 : while (E1.time >= AE.FiringTime) {
209 : //cerr << "looping on autonomous edge";
210 172875 : double eTime = AE.FiringTime - A.CurrentTime;
211 172875 : A.updateLHA(eTime , N.getState());
212 172875 : static_cast<S*>(this)->printLog(eTime);
213 172875 : A.fireAutonomous(AE.Index,N.getState());
214 172875 : if(P.verbose>3){
215 0 : cerr << "Autonomous transition:" << AE.Index << " "<<E1.time << ">=" << AE.FiringTime<< endl;
216 0 : A.printState(cerr);
217 0 : cerr << endl;
218 : }
219 172875 : static_cast<S*>(this)->printLog(eTime);
220 172875 : if (A.isFinal()) {
221 23140 : static_cast<S*>(this)->returnResultTrue();
222 23140 : return false;
223 149735 : } else AE = A.GetEnabled_A_Edges( N.getState());
224 : }
225 52290534 : if(P.verbose>3){
226 0 : cerr << "\033[1;33mFiring:\033[0m" << N.getTransitionLabel(E1.transition) ;
227 0 : E1.binding.print();
228 : }
229 :
230 : //Make time elapse in the LHA
231 52290534 : double eTime = E1.time - A.CurrentTime;
232 52290534 : A.updateLHA( eTime, N.getState() );
233 52290534 : if(eTime>0.0){ current_weight = 0.0; }
234 : else {
235 10742378 : current_weight = E1.weight;
236 : }
237 : //Print the state of the system after the time elapse and the transition name
238 52290534 : static_cast<S*>(this)->printLog(eTime,E1.transition);
239 :
240 : //Fire the transition in the SPN
241 52290534 : N.fire(E1.transition, E1.binding, A.CurrentTime);
242 :
243 52290534 : if(E1.transition >= A.NbTrans ){
244 0 : if(P.verbose>3)cerr << " Transition filtered " << endl;
245 0 : N.update(A.CurrentTime, E1.transition, E1.binding, *EQ,*this);
246 0 : return true;
247 : }
248 :
249 : //Check if there exist a valid Synchronisation in the automata and fire it.
250 52290534 : int SE = A.synchroniseWith(E1.transition, N.getState(), E1.binding);
251 :
252 : //If no synchronisation is possible the trajectory is rejected
253 52290534 : if (SE < 0) {
254 73008 : if(P.verbose>3)cerr << " No synchronisation" << endl;
255 73008 : Result.accept=false;
256 73008 : return false;
257 : } else {
258 52217526 : if(P.verbose>3)cerr << " Synchronisation with " << SE << endl;
259 : //If synchronisation is possible check if the
260 : // reached state is final. Then update the SPN.
261 52217526 : if (A.isFinal()) {
262 55301 : static_cast<S*>(this)->returnResultTrue();
263 55301 : return false;
264 : } else {
265 52162225 : N.update(A.CurrentTime, E1.transition, E1.binding, *EQ,*this);
266 : }
267 : }
268 : }
269 52162225 : return true;
270 : }
271 :
272 : /**
273 : * Interactive mode stop the simulation until the user choose a transition.
274 : */
275 : template <class S, class EQT,class DEDS>
276 0 : void SimulatorBase<S,EQT,DEDS>::interactiveSimulation(){
277 0 : string input_line;
278 0 : if((waitForTransition >0
279 0 : && (size_t)waitForTransition != N.lastTransition))return;
280 0 : if(A.CurrentTime < minInteractiveTime
281 0 : && (waitForTransition < 0
282 0 : || (size_t)waitForTransition != N.lastTransition))return;
283 0 : bool continueLoop = true;
284 0 : waitForTransition = -1;
285 0 : while(continueLoop){
286 0 : cerr << "\033[1;31mCosmosSimulator>\033[0m";
287 0 : auto trLabels = N.getTransitionLabels();
288 0 : if (cin.good()) {
289 0 : getline(cin, input_line);
290 0 : if(input_line.substr(0,5).compare("fire ")==0){
291 0 : string trans = input_line.substr(5,input_line.size()-5);
292 : size_t traid;
293 0 : for(traid=0; traid < trLabels.size() && trLabels[traid]!= trans ; traid++) ;
294 0 : if(traid == trLabels.size())cerr << "Unknown transition: "<< trans << endl;
295 : else{
296 0 : if(EQ->isScheduled(traid, 0)){
297 0 : Event &F = EQ->getEvent(traid,0);
298 0 : F.setTime(A.CurrentTime);
299 0 : F.setPriority((size_t)-1);
300 0 : F.setWeight(0.0);
301 0 : EQ->replace(F);
302 0 : continueLoop = false;
303 0 : }else cerr << "Transition: "<< trans << " is not enabled" << endl;
304 : }
305 :
306 :
307 0 : }else if(input_line.compare("step")==0)continueLoop=false;
308 0 : else if(input_line.compare("stop")==0)exit(EXIT_SUCCESS);
309 0 : else if(input_line == "display" || input_line == "d"){
310 0 : if(!dotFile.empty()){
311 0 : stringstream ss;
312 0 : ss << "sed ";
313 0 : N.printSedCmd(ss);
314 0 : EQ->printSedCmd(N.getTransitionLabels(), ss);
315 0 : ss << P.tmpPath << "/templatePetriNet.dot > " << P.tmpPath << "/PetriNet.dot; ";
316 0 : ss << "dot "<< P.tmpPath <<"/PetriNet.dot -Tpdf -o " << dotFile << endl;
317 0 : const auto retval = system(ss.str().c_str());
318 0 : if(retval >0)cerr << "Fail to lauch graphviz with error code" << retval;
319 0 : } else cerr << "No dot output specified!" << endl;
320 0 : } else if(input_line.substr(0,5)=="wait "){
321 0 : const auto arg = input_line.substr(5,input_line.length()-5);
322 0 : const auto transitionsLabels = N.getTransitionLabels();
323 0 : const auto trid = find_if(transitionsLabels.begin(), transitionsLabels.end(),
324 0 : [&] (string tr){return (tr == arg);});
325 0 : if(trid != transitionsLabels.end()){
326 0 : waitForTransition = trid - transitionsLabels.begin();
327 0 : continueLoop = false;
328 : }else try {
329 0 : minInteractiveTime= A.CurrentTime + stod(arg);
330 0 : continueLoop = false;
331 0 : } catch (const invalid_argument& ia) {
332 0 : cerr << "Fail to parse time!" << endl;
333 0 : } catch (const out_of_range& ia) {
334 0 : cerr << "Fail to parse time!" << endl;
335 : }
336 :
337 0 : } else if(input_line.compare("exit")==0){
338 0 : exit(EXIT_SUCCESS);
339 0 : } else if(input_line.compare("help")==0 || input_line.compare("h")==0){
340 0 : cerr << "Available command:\n\thelp:\tdisplay this message"<<endl;
341 0 : cerr << "\ts, step:\tmake one step of simulation" << endl;
342 0 : cerr << "\td, draw:\tdraw the GSPN with dot" << endl;
343 0 : cerr << "\tfire tr:\tfire transition tr" << endl;
344 0 : cerr << "\twait tr:\twait until transition tr occurs" << endl;
345 0 : cerr << "\twait t:\twait t times unit" << endl;
346 0 : cerr << "\texit t:\texit simulator" << endl;
347 :
348 0 : } else if (input_line.compare("s")==0 || input_line.compare("s")==0)continueLoop=false;
349 0 : else if (input_line.compare("")==0);
350 : else {
351 0 : cerr << "Command not found:" << input_line << endl;
352 : }
353 :
354 : }else {
355 0 : if(cin.eof())exit(EXIT_SUCCESS);
356 0 : cerr << "error on the input stream\n";
357 0 : exit(EXIT_FAILURE);
358 : }
359 :
360 : }
361 :
362 :
363 : }
364 :
365 : /**
366 : * Simulate a whole trajectory in the system. Result is store in SimOutput
367 : */
368 : template <class S, class EQT,class DEDS>
369 216102 : void SimulatorBase<S,EQT,DEDS>::SimulateSinglePath() {
370 :
371 216102 : static_cast<S*>(this)->reset();
372 216102 : N.initialEventsQueue(*EQ,*this);
373 216102 : minInteractiveTime=0.0;
374 216102 : waitForTransition= -1;
375 :
376 216102 : if(logtrace.is_open())logtrace << "New Path"<< endl;
377 : //cerr << "start path"<< endl;
378 :
379 216102 : bool continueb = true;
380 216102 : lastSampled = -sampleTrace;
381 103661062 : while (continueb) {
382 : //cerr << "continue path"<< endl;
383 51722480 : static_cast<S*>(this)->printLog(0.0);
384 51722480 : if(P.verbose>3){
385 : //Print marking and location of the automata
386 : //Usefull to track a simulation
387 0 : N.printHeader(cerr);
388 0 : A.printHeader(cerr);
389 0 : cerr << endl;
390 0 : N.printState(cerr,0.0);
391 0 : A.printState(cerr);
392 0 : cerr << endl;
393 0 : if(P.verbose>4)EQ->view(N.getTransitionLabels());
394 0 : if(P.verbose==6)static_cast<S*>(this)->interactiveSimulation();
395 : }
396 :
397 51722480 : continueb = static_cast<S*>(this)->SimulateOneStep();
398 : }
399 216102 : if(P.verbose>3){
400 : //Print marking and location of the automata
401 : //Usefull to track a simulation
402 0 : N.printHeader(cerr);
403 0 : A.printHeader(cerr);
404 0 : cerr << endl;
405 0 : N.printState(cerr,0.0);
406 0 : A.printState(cerr);
407 0 : cerr << endl;
408 0 : if(P.verbose>4)EQ->view(N.getTransitionLabels());
409 : }
410 : //cerr << "finish path"<< endl;
411 216102 : }
412 :
413 : template <class S, class EQT,class DEDS>
414 253 : BatchR SimulatorBase<S,EQT,DEDS>::RunBatch(){
415 253 : auto starttime = chrono::steady_clock::now();
416 253 : auto currenttime = chrono::steady_clock::now();
417 253 : chrono::duration<double> timesize(0.03);
418 253 : BatchR batchResult(A.FormulaVal.size(),A.FormulaValQual.size());
419 444457 : while ((batchResult.I < BatchSize && BatchSize!=0) || (currenttime-starttime < timesize && BatchSize==0) ) {
420 222102 : static_cast<S*>(this)->SimulateSinglePath();
421 222102 : batchResult.addSim(Result);
422 222102 : if(P.verbose>3)cerr << batchResult;
423 :
424 222102 : if (Result.accept && logResult){
425 1237487 : for(size_t i=0; i<Result.quantR.size();i++){
426 1089046 : if (i>0)logvalue << "\t";
427 1089046 : logvalue << Result.quantR[i];
428 : }
429 148441 : logvalue << endl;
430 : }
431 222102 : currenttime=chrono::steady_clock::now();
432 : }
433 253 : batchResult.simTime = (currenttime - starttime).count();
434 253 : return batchResult;
435 : }
|