Changeset 1422 for applications
- Timestamp:
- 01/26/12 12:23:52 (13 years ago)
- Location:
- applications/doprava
- Files:
-
- 8 added
- 2 modified
Legend:
- Unmodified
- Added
- Removed
-
applications/doprava/QuadraticMinimalizator.h
r1420 r1422 25 25 e << Q << endl; 26 26 it_warning(e.str()); 27 return false; 27 28 } 28 29 if ( !isDiagonal(R) ) { … … 31 32 e << R << endl; 32 33 it_warning(e.str()); 34 return false; 33 35 } 34 36 if ( !isPositiveDefinite(Q) ) { … … 37 39 e << Q << endl; 38 40 it_warning(e.str()); 41 return false; 39 42 } 40 43 if ( !isPositiveDefinite(R) ) { … … 43 46 e << R << endl; 44 47 it_warning(e.str()); 48 return false; 45 49 } 50 return true; 46 51 } 47 52 … … 56 61 } 57 62 58 59 /**60 * x(t+1) = Ax(t) + Bu(t) + C61 * @param horizont62 * @param C63 * @return L: u(t) = L [ x(t) 1 ]64 */65 mat L( const int horizont, const mat & C ) {66 it_warning("Nijak se tam C neprojevuje!!!");67 int xdim = A.cols();68 int udim = B.cols();69 mat sqQ = sqrt(Q);70 mat sqR = sqrt(R);71 mat M0, M;72 mat qrQ;73 mat Lall, Lc1, Lc2, Lc3, Lx1, Lx2, Lu, L;74 M0 = concat_vertical(75 concat_vertical(76 concat_horizontal( sqQ * B, concat_horizontal( sqQ * A, zeros(xdim, udim) ) ),77 concat_horizontal( sqR, concat_horizontal( zeros( udim, xdim ), zeros(udim, udim) ) )78 ),79 concat_horizontal( zeros( xdim, xdim + udim ) , sqQ * C)80 );81 82 M = M0;83 for (int h = horizont; h >= 0; h --) {84 qr(M, Lall);85 86 Lu = Lall(0, udim-1, 0, udim-1 );87 Lx1 = Lall(0, udim-1, udim, udim+xdim-1);88 Lc1 = Lall(0, udim-1, udim+xdim, udim+xdim );89 Lx2 = Lall(udim, udim+xdim-1, udim, udim+xdim-1);90 Lc2 = Lall(udim, udim+xdim-1, udim+xdim, udim+xdim);91 Lc3 = Lall(udim+xdim, udim+xdim+xdim-1, udim+xdim, udim+xdim);92 M = concat_vertical(93 M0,94 concat_vertical(95 concat_horizontal( Lx2 * B, concat_horizontal( Lx2 * A, Lc2 ) ),96 concat_horizontal( zeros(xdim, xdim + udim), Lc3 )97 )98 );99 }100 L = concat_horizontal( - inv(Lu) * Lx1, - inv(Lu) * Lc1);101 return L;102 }103 63 /** 104 64 * x(t+1) = Ax(t) + Bu(t) … … 134 94 } 135 95 96 /* 136 97 double maxDif( const mat & M1, const mat & M2 ) { 137 98 return max(max(abs(M1 - M2), 1)); … … 150 111 return false; 151 112 } 152 113 */ 153 114 bool isDiagonal( const mat M ) { 154 115 mat D = abs(M); … … 162 123 return true; 163 124 } 125 164 126 165 127 bool isPositiveDefinite( const mat M ) { -
applications/doprava/traffic_agent_LQ.h
r1421 r1422 1 1 2 #include "BaseTrafficAgentCt.h" 2 3 #include "QuadraticMinimalizator.h" 3 4 4 class TrafficAgentLQ : public BaseTrafficAgent {5 class TrafficAgentLQ : public BaseTrafficAgentCt { 5 6 protected: 6 mat A, B, Q, R;7 static const int T = 90;8 int L;9 int nOfBroadcast;10 7 QuadraticMinimalizator * minimizer; 11 //map<string, int> queueIndex;8 12 9 public: 13 10 // FUNKCE VOLANE V main_loop NA ZACATKU 14 11 12 void from_setting(const Setting& set) { 13 BaseTrafficAgentCt::from_setting(set); 14 } 15 15 16 void validate (){ 16 L = 90; 17 BaseTrafficAgent::validate(); 18 19 rv_action = RV("Tc",1); 20 rv_action.add( RV( stage_names, ones_i(stage_names.length()) ) ); 21 22 for ( int i=0; i<lanehs.length(); i ++ ) { 23 unsigned int index = find_index( green_names, name + "_" + lanehs(i)->getSG() ); 24 lanehs(i)->green_time_ratio = green_times( index ); 25 } 26 27 Q = diag( ones( queues.length()) ); 28 R = "1"; 29 A = diag( ones( queues.length()) ); 30 B = ones( queues.length(), 1); 31 32 minimizer = new QuadraticMinimalizator(A,B,Q,R); 17 BaseTrafficAgentCt::validate(); 33 18 34 19 cout << "HELLO! " << name << " IS HERE!" << endl; 35 /*36 cout << "inputs " << inputs.length() << endl;37 cout << "queues " << queues.length() << endl;38 cout << "rv_inputs" << rv_inputs.length() << endl << rv_inputs << endl;39 cout << "rv_queues" << rv_queues.length() << endl << rv_queues << endl;40 cout << "lanehs " << lanehs.length() << endl << endl;41 cout << "A B Q R " << endl << A << endl << B << endl << Q << endl << R << endl << endl;42 */43 20 } 44 21 45 22 // FUNKCE VOLANE V main_loop V KAZDEM CYKLU 46 void adapt (const vec &glob_dt) { 47 for ( int i=0; i<lanehs.length(); i ++ ) {48 //lanehs(i)->inputs( 0 ) = inputs( 2*find_index(rv_inputs, lanehs(i)->rv_inputs.name(0) ) ); 49 //lanehs(i)->inputs( 1 ) = inputs( 2*find_index(rv_inputs, lanehs(i)->rv_inputs.name(0) ) + 1);50 lanehs(i)->queue = queues( find_index( rv_queues, lanehs(i)->getQueueName() ) );51 lanehs(i)->countAvgs();52 }53 23 void adapt (const vec &glob_dt) { 24 BaseTrafficAgentCt::adapt(glob_dt); 25 26 Q = 100 * diag( ones( queues.length()) ) + diag(computedQueues); 27 R = "0.0000001"; 28 A = diag( ones( queues.length()) ); 29 B = zeros( queues.length(), 1); 30 54 31 // jeden clen matice B 55 32 for ( int i=0; i<queues.length(); i ++ ) { 56 B(i,0) = 0.5* lanehs(i)->green_time_ratio;33 B(i,0) = - s_flow(i) * lanehs(i)->green_time_ratio; 57 34 } 58 35 59 36 60 61 BaseTrafficAgent::adapt(glob_dt);62 37 } 63 38 … … 71 46 stringstream coefStr; 72 47 coefStr << "coef" << lanehs(i)->rv_outputs.name(j); 73 double val = 0.5* lanehs(i)->green_time_ratio * lanehs(i)->getLane().alpha(j);48 double val = s_flow(i) * lanehs(i)->green_time_ratio * lanehs(i)->getLane().alpha(j); 74 49 sendToAll<double>(set, coefStr.str(), val); 75 50 } … … 77 52 } 78 53 } 54 55 56 if ( nOfBroadcast == 1 ) { 57 broadcastTc(set); 58 } 59 79 60 nOfBroadcast ++; 61 } 62 63 double votingWeight() { 64 //vec computedQueues = getComputedQueues(); 65 double qsum = 0; 66 for ( int i = 0; i < computedQueues.length(); i++ ) { 67 qsum += computedQueues(i); 68 } 69 if ( qsum > 1 ) 70 return qsum; 71 else 72 return 1; 73 } 74 75 void broadcastTc( Setting & set) { 76 // posli idealni tc 77 Tc_computed = computeTc(); 78 stringstream tcstr; 79 tcstr << "timecycle" << name; 80 81 // posli vahu hlasu (suma front) 82 //vec computedQueues = getComputedQueues(); 83 84 vec tc_vec = "0.0 0.0"; 85 tc_vec(0) = Tc_computed; 86 tc_vec(1) = votingWeight(); 87 // posilam vektor [tc, w] 88 sendToAll<vec>(set, tcstr.str(), tc_vec); 80 89 } 81 90 … … 88 97 double val; 89 98 UI::get(val, msg, "value"); 90 B(i,0) = B(i,0) - val; 91 //cout << name << " received coefs " << inputName << " " << what << endl << val <<endl; 99 B(i,0) = B(i,0) + val; 92 100 } 101 102 if ( what.substr(0,9) == "timecycle" ) { 103 vec val; 104 UI::get(val, msg, "value"); 105 cout << "receiving Tc " << val << endl; 106 Tc_sum += val(0) * val(1); 107 Tc_w_sum += val(1); 108 } 109 93 110 } 94 111 95 void act (vec &glob_ut){ 96 //cout <<"A act"<<endl<<A<<endl; 97 B = B * T * L; 98 minimizer = new QuadraticMinimalizator(A, B, Q, R); 99 vec computedQueues = getComputedQueues(); 100 vec res = minimizer->minimize(computedQueues, 100); 101 //cout << "minimalization res " << name << endl << res << endl << endl; 102 delete minimizer; 103 int Tc = 80; 104 vec action; 105 action.set_size(rv_action._dsize()); 106 action(find_index(rv_action, "Tc")) = Tc; 107 int st; 108 int stage_time_sum = 0; // soucet delky fazi 109 action(find_index(rv_action, "Tc")) = Tc; 110 for ( int i =0; i < stage_names.length(); i ++) { 111 if ( (i+1) < stage_names.length() ) { 112 st = round(((double)(stage_times(i)*Tc))/80.0); 113 stage_time_sum += st; 114 action(find_index(rv_action, stage_names(i))) = st; 115 } 116 else { // dopocitani posledni faze - oprava zaokrouhlovaci chyby 117 action(find_index(rv_action, stage_names(i))) = Tc - stage_time_sum; 118 } 119 } 120 action2ds.filldown(action,glob_ut); 112 void act (vec &glob_ut){ 113 //vec computedQueues = getComputedQueues(); 114 cout << endl << name << " ACT" << endl; 115 cout << "Tc computed = " << Tc_computed << endl; 116 cout << "queues = " << round( computedQueues )<< endl; 117 double w = votingWeight(); 118 Tc_sum += Tc_computed * w; 119 Tc_w_sum += w; 120 int Tc = (int)round( Tc_sum / Tc_w_sum ); 121 setCycleTime(Tc, glob_ut); 122 looper ++; 123 } 124 // VYPOCETNI FUNKCE 125 126 // saturovany tok pruhu i 127 double s_flow(int i) { 128 double min_flow = 0.3; 129 double max_flow = 0.5; 130 double input = 0; 131 if ( 2*i < inputs.length() ) 132 input = inputs(2*i); 133 134 vec q = getComputedQueues(); 135 double flow = q(i); 136 if ( input > 0 ) 137 flow += input; 138 flow = flow/T; 139 if ( flow > max_flow ) 140 return max_flow; 141 if ( flow < min_flow ) 142 return min_flow; 143 else 144 return flow; 121 145 } 122 146 123 void setCycleTime( int Tc, vec &glob_ut ) { 124 147 double computeTc() { 148 B = B * T; 149 //vec computedQueues = getComputedQueues(); 150 minimizer = new QuadraticMinimalizator(A, B, Q, R); 151 mat L_mat = minimizer->L( 100); 152 vec u_vec = L_mat * computedQueues; 153 double u = u_vec(0); 154 double Tc_computed = L / (1 - u); 155 delete minimizer; 156 if ( u >= 1 || Tc_computed > Tc_max ) { 157 cout << endl << "Tc error : " << Tc_computed << " setting to max" << endl; 158 Tc_computed = Tc_max; 159 } 160 if ( Tc_computed < Tc_min ) { 161 cout << endl << "Tc error : " << Tc_computed << " setting to min" << endl; 162 Tc_computed = Tc_min; 163 } 164 return Tc_computed; 125 165 } 166 167 168 126 169 127 170 128 171 129 // POMOCNE FCE130 vec getComputedQueues() {131 vec q = zeros( queues.length() );132 for ( int i = 0; i < lanehs.length(); i++ ) {133 q(i) = lanehs(i)->queue_avg;134 }135 return q;136 }137 172 138 template <class T> void sendToAll( Setting &set, string msgCode, T value) {139 for ( int i = 0; i < neighbours.length(); i ++ ) {140 Setting & msg = set.add(Setting::TypeGroup);141 UI::save( neighbours(i), msg, "to" );142 UI::save (name,msg,"from");143 UI::save( msgCode, msg, "what" );144 UI::save( value, msg, "value" );145 }146 }147 173 148 void printVector ( RV rv_vector, vec vector, string description ) { 149 cout << endl << description << " " << name << endl; 150 int k = 0; 151 for ( int i = 0; i < rv_vector.length(); i ++ ) { 152 cout << rv_vector.name(i) << " : "; 153 for ( int j = 0; j < rv_vector.size(i); j ++ ) { 154 cout << vector(k) << " "; 155 k ++; 156 } 157 cout << endl; 158 } 159 cout << endl; 160 } 161 162 unsigned int find_index ( RV rv_vector, string index_string ) { 163 for ( unsigned int i = 0; i < rv_vector.length(); i ++) { 164 if ( rv_vector.name(i) == index_string ) 165 return i; 166 } 167 return -1; 168 } 169 170 unsigned int find_index ( Array <string> arr, string index_string ) { 171 for ( unsigned int i = 0; i < arr.length(); i ++) { 172 if ( arr(i) == index_string ) 173 return i; 174 } 175 return -1; 174 // KONEC 175 ~TrafficAgentLQ() { 176 176 177 } 177 178