LINE Solver
MATLAB API documentation
Loading...
Searching...
No Matches
afterEventStation.m
1function [outspace, outrate, outprob, eventCache] = afterEventStation(sn, ind, inspace, event, class, isSimulation, eventCache, ...
2 M, R, S, phasessz, phaseshift, pie, isf, ismkvmod, ismkvmodclass, lldscaling, lldlimit, cdscaling, ...
3 hasOnlyExp, ist, K, Ks, mu, phi, proc, capacity, classcap, V, space_buf, space_srv, space_var, key)
4outspace = [];
5outrate = [];
6outprob = 1;
7switch event
8 case EventType.ARV %% passive
9 % return if there is no space to accept the arrival
10 [ni,nir] = State.toMarginalAggr(sn,ind,inspace,K,Ks,space_buf,space_srv,space_var);
11 % otherwise check scheduling strategy
12 pentry = pie{ist}{class};
13 % For Place nodes (INF scheduling with NaN service), use uniform entry probability
14 if all(isnan(pentry))
15 pentry = ones(size(pentry)) / length(pentry);
16 end
17 outprob = [];
18 outprob_k = [];
19 for kentry = 1:K(class)
20 space_var_k = space_var;
21 space_srv_k = space_srv;
22 space_buf_k = space_buf;
23 switch sn.sched(ist)
24 case SchedStrategy.EXT % source, can receive any "virtual" arrival from the sink as long as it is from an open class
25 if isinf(sn.njobs(class))
26 outspace = inspace;
27 outrate = -1*zeros(size(outspace,1)); % passive action, rate is unspecified
28 outprob = ones(size(outspace,1));
29 break
30 end
31 case {SchedStrategy.PS, SchedStrategy.INF, SchedStrategy.DPS, SchedStrategy.GPS, SchedStrategy.PSPRIO, SchedStrategy.DPSPRIO, SchedStrategy.GPSPRIO, SchedStrategy.LPS}
32 % job enters service immediately
33 if space_srv_k(:,Ks(class)+kentry) < classcap(ist,class)
34 space_srv_k(:,Ks(class)+kentry) = space_srv_k(:,Ks(class)+kentry) + 1;
35 outprob_k = pentry(kentry)*ones(size(space_srv_k,1));
36 else
37 outprob_k = pentry(kentry)*zeros(size(space_srv_k,1));
38 end
39 case {SchedStrategy.SIRO, SchedStrategy.SEPT, SchedStrategy.LEPT}
40 if ni<S(ist)
41 space_srv_k(:,Ks(class)+kentry) = space_srv_k(:,Ks(class)+kentry) + 1;
42 outprob_k = pentry(kentry)*ones(size(space_srv_k,1));
43 else
44 space_buf_k(:,class) = space_buf_k(:,class) + 1;
45 outprob_k = pentry(kentry)*ones(size(space_srv_k,1));
46 end
47 case {SchedStrategy.FCFS, SchedStrategy.HOL, SchedStrategy.LCFS, SchedStrategy.LCFSPRIO}
48 % find states with all servers busy - this
49 % needs not to be moved
50
51 % if MAP service, when empty restart from the phase
52 % stored in space_var for this class
53 if ~ismkvmodclass(class) || (ismkvmodclass(class) && kentry == space_var(sum(sn.nvars(ind,1:class))))
54 if ismkvmodclass(class)
55 pentry = zeros(size(pentry));
56 pentry(kentry) = 1.0;
57 end
58 all_busy_srv = sum(space_srv_k,2) >= S(ist);
59
60 % find and modify states with an idle server
61 idle_srv = sum(space_srv_k,2) < S(ist);
62 space_srv_k(idle_srv, end-sum(K)+Ks(class)+kentry) = space_srv_k(idle_srv,end-sum(K)+Ks(class)+kentry) + 1; % job enters service
63
64 % this section dynamically grows the number of
65 % elements in the buffer
66
67 if any(ni < capacity(ist))
68 if any(nir(:,class) < classcap(ist,class)) % if there is room
69 if ~any(space_buf_k(:)==0) % but the buffer has no empty slots
70 % append job slot
71 space_buf_k = [zeros(size(space_buf_k,1),1),space_buf_k];
72 end
73 end
74 end
75 %end
76 %get position of first empty slot
77 empty_slots = -1*ones(size(all_busy_srv,1),1);
78 if size(space_buf_k,2) == 0
79 empty_slots(all_busy_srv) = false;
80 elseif size(space_buf_k,2) == 1
81 empty_slots(all_busy_srv) = space_buf_k(all_busy_srv,:)==0;
82 else
83 empty_slots(all_busy_srv) = max(bsxfun(@times, space_buf_k(all_busy_srv,:)==0, [1:size(space_buf_k,2)]),[],2);
84 end
85
86 % ignore states where the buffer has no empty slots
87 wbuf_empty = empty_slots>0;
88 if any(wbuf_empty)
89 space_srv_k = space_srv_k(wbuf_empty,:);
90 space_buf_k = space_buf_k(wbuf_empty,:);
91 space_var_k = space_var_k(wbuf_empty,:);
92 empty_slots = empty_slots(wbuf_empty);
93 space_buf_k(sub2ind(size(space_buf_k),1:size(space_buf_k,1),empty_slots')) = class;
94 %outspace(all_busy_srv(wbuf_empty),:) = [space_buf, space_srv, space_var];
95 end
96 outprob_k = pentry(kentry)*ones(size(space_srv_k,1),1);
97 else
98 outprob_k = 0*ones(size(space_srv_k,1),1); % zero probability event
99 end
100 case {SchedStrategy.FCFSPR,SchedStrategy.FCFSPI,SchedStrategy.FCFSPRPRIO,SchedStrategy.FCFSPIPRIO,SchedStrategy.LCFSPR,SchedStrategy.LCFSPI,SchedStrategy.LCFSPRPRIO,SchedStrategy.LCFSPIPRIO}
101 % find states with all servers busy - this
102 % must not be moved
103 all_busy_srv = sum(space_srv_k,2) >= S(ist);
104 % find states with an idle server
105 idle_srv = sum(space_srv_k,2) < S(ist);
106
107 % reorder states so that idle ones come first
108 space_buf_k_reord = space_buf_k(idle_srv,:);
109 space_srv_k_reord = space_srv_k(idle_srv,:);
110 space_var_k_reord = space_var_k(idle_srv,:);
111
112 % if idle, the job enters service in phase kentry
113 if any(idle_srv)
114 space_srv_k_reord(:, end-sum(K)+Ks(class)+kentry) = space_srv_k_reord(:,end-sum(K)+Ks(class)+kentry) + 1;
115 outprob_k = pentry(kentry);
116 else
117 % if all busy, expand output states for all possible choices of job class to preempt
118 psentry = ones(size(space_buf_k_reord,1),1); % probability scaling due to preemption
119 for classpreempt = 1:R
120 % For priority variants, only higher-priority (lower class number) jobs can preempt
121 if (sn.sched(ist) == SchedStrategy.FCFSPRPRIO || sn.sched(ist) == SchedStrategy.FCFSPIPRIO || sn.sched(ist) == SchedStrategy.LCFSPRPRIO || sn.sched(ist) == SchedStrategy.LCFSPIPRIO)
122 if class >= classpreempt % arriving job has same or lower priority, cannot preempt
123 continue;
124 end
125 end
126 for phasepreempt = 1:K(classpreempt) % phase of job to preempt
127 si_preempt = space_srv_k(:, (end-sum(K)+Ks(classpreempt)+phasepreempt));
128 busy_preempt = si_preempt > 0; % states where there is at least on class-r job in execution
129 if any(busy_preempt)
130 psentry = [psentry; si_preempt(busy_preempt) ./ sum(space_srv_k,2)];
131 space_srv_k_preempt = space_srv_k(busy_preempt,:);
132 space_buf_k_preempt = space_buf_k(busy_preempt,:);
133 space_var_k_preempt = space_var_k(busy_preempt,:);
134 space_srv_k_preempt(:, end-sum(K)+Ks(classpreempt)+phasepreempt) = space_srv_k_preempt(:,end-sum(K)+Ks(classpreempt)+phasepreempt) - 1; % remove preempted job
135 space_srv_k_preempt(:, end-sum(K)+Ks(class)+kentry) = space_srv_k_preempt(:,end-sum(K)+Ks(class)+kentry) + 1;
136
137 % dynamically grow buffer lenght in
138 % simulation
139 if isSimulation
140 if ni < capacity(ist) && nir(class) < classcap(ist,class) % if there is room
141 if ~any(space_buf_k_preempt(:)==0) % but the buffer has no empty slots
142 % append job slot
143 space_buf_k_preempt = [zeros(size(space_buf_k_preempt,1),2),space_buf_k]; % append two columns for (class, preempt-phase)
144 end
145 end
146 end
147
148 %get position of first empty slot
149 empty_slots = -1*ones(sum(busy_preempt),1);
150 if size(space_buf_k_preempt,2) == 0
151 empty_slots(busy_preempt) = false;
152 elseif size(space_buf_k_preempt,2) == 2 % 2 due to (class, preempt-phase) pairs
153 empty_slots(busy_preempt) = space_buf_k_preempt(busy_preempt,1:2:end)==0;
154 else
155 empty_slots(busy_preempt) = max(bsxfun(@times, space_buf_k_preempt(busy_preempt,:)==0, [1:size(space_buf_k_preempt,2)]),[],2)-1; %-1 due to (class, preempt-phase) pairs
156 end
157
158 % ignore states where the buffer has no empty slots
159 wbuf_empty = empty_slots>0;
160 if any(wbuf_empty)
161 space_srv_k_preempt = space_srv_k_preempt(wbuf_empty,:);
162 space_buf_k_preempt = space_buf_k_preempt(wbuf_empty,:);
163 space_var_k_preempt = space_var_k_preempt(wbuf_empty,:);
164 empty_slots = empty_slots(wbuf_empty);
165 if sn.sched(ist) == SchedStrategy.LCFSPR || sn.sched(ist) == SchedStrategy.LCFSPRPRIO || sn.sched(ist) == SchedStrategy.FCFSPR || sn.sched(ist) == SchedStrategy.FCFSPRPRIO % preempt-resume
166 space_buf_k_preempt(sub2ind(size(space_buf_k_preempt),1:size(space_buf_k_preempt,1),empty_slots')+1) = phasepreempt;
167 elseif sn.sched(ist) == SchedStrategy.LCFSPI || sn.sched(ist) == SchedStrategy.LCFSPIPRIO || sn.sched(ist) == SchedStrategy.FCFSPI || sn.sched(ist) == SchedStrategy.FCFSPIPRIO % preempt-independent
168 space_buf_k_preempt(sub2ind(size(space_buf_k_preempt),1:size(space_buf_k_preempt,1),empty_slots')+1) = 1;
169 end
170 space_buf_k_preempt(sub2ind(size(space_buf_k_preempt),1:size(space_buf_k_preempt,1),empty_slots')) = classpreempt;
171 %outspace(all_busy_srv(wbuf_empty),:) = [space_buf, space_srv, space_var];
172 end
173 space_srv_k_reord = [space_srv_k_reord; space_srv_k_preempt];
174 space_buf_k_reord = [space_buf_k_reord; space_buf_k_preempt];
175 space_var_k_reord = [space_var_k_reord; space_var_k_preempt];
176 end
177 end
178 end
179 outprob_k = pentry(kentry) * psentry .* ones(size(space_srv_k_reord,1),1);
180 end
181 space_buf_k = space_buf_k_reord; % save reordered output states
182 space_srv_k = space_srv_k_reord; % save reordered output states
183 space_var_k = space_var_k_reord; % save reordered output states
184 end
185 % form the new state
186 outspace_k = [space_buf_k, space_srv_k, space_var_k];
187 % remove states where new arrival violates capacity or cutoff constraints
188 [oi,oir] = State.toMarginalAggr(sn,ind,outspace_k,K,Ks,space_buf_k,space_srv_k,space_var_k);
189 en_o = classcap(ist,class)>= oir(:,class) | capacity(ist)*ones(size(oi,1),1) >= oi;
190
191 if size(outspace,2)>size(outspace_k(en_o,:),2)
192 outspace = [outspace; zeros(1,size(outspace,2)-size(outspace_k(en_o,:),2)),outspace_k(en_o,:)];
193 elseif size(outspace,2)<size(outspace_k(en_o,:),2)
194 outspace = [zeros(size(outspace,1),size(outspace_k(en_o,:),2)-size(outspace,2)), outspace; outspace_k(en_o,:)];
195 else
196 outspace = [outspace; outspace_k(en_o,:)];
197 end
198 outrate = [outrate; -1*ones(size(outspace_k(en_o,:),1),1)]; % passive action, rate is unspecified
199 outprob = [outprob; outprob_k(en_o,:)];
200 end
201 if isSimulation
202 if size(outprob,1) > 1
203 cum_prob = cumsum(outprob) / sum(outprob);
204 firing_ctr = 1 + max([0,find( rand > cum_prob' )]); % select action
205 outspace = outspace(firing_ctr,:);
206 outrate = -1;
207 outprob = 1;
208 end
209 end
210 case EventType.DEP
211 if any(any(space_srv(:,(Ks(class)+1):(Ks(class)+K(class))))) % something is busy
212 if hasOnlyExp && (sn.sched(ist) == SchedStrategy.PS || sn.sched(ist) == SchedStrategy.DPS || sn.sched(ist) == SchedStrategy.GPS || sn.sched(ist) == SchedStrategy.INF || sn.sched(ist) == SchedStrategy.PSPRIO || sn.sched(ist) == SchedStrategy.DPSPRIO || sn.sched(ist) == SchedStrategy.GPSPRIO || sn.sched(ist) == SchedStrategy.LPS)
213 nir = space_srv;
214 ni = sum(nir,2);
215 sir = nir;
216 kir = sir;
217 else
218 [ni,nir,sir,kir] = State.toMarginal(sn,ind,inspace,K,Ks,space_buf,space_srv,space_var);
219 end
220 switch sn.routing(ind,class)
221 case RoutingStrategy.RROBIN
222 idx = find(space_var(sum(sn.nvars(ind,1:(R+class)))) == sn.nodeparam{ind}{class}.outlinks);
223 if idx < length(sn.nodeparam{ind}{class}.outlinks)
224 space_var(sum(sn.nvars(ind,1:(R+class)))) = sn.nodeparam{ind}{class}.outlinks(idx+1);
225 else
226 space_var(sum(sn.nvars(ind,1:(R+class)))) = sn.nodeparam{ind}{class}.outlinks(1);
227 end
228 end
229 if sir(class)>0 % is a job of class is in service
230 outprob = [];
231 for k=1:K(class)
232 space_srv = inspace(:,(end-sum(K)-V+1):(end-V)); % server state
233 space_buf = inspace(:,1:(end-sum(K)-V)); % buffer state
234 rate = zeros(size(space_srv,1),1);
235 en = space_srv(:,Ks(class)+k) > 0;
236 if any(en)
237 switch sn.sched(ist)
238 case SchedStrategy.EXT % source, can produce an arrival from phase-k as long as it is from an open class
239 if isinf(sn.njobs(class))
240 pentry = pie{ist}{class};
241 for kentry = 1:K(class)
242 space_srv = inspace(:,(end-sum(K)-V+1):(end-V)); % server state
243 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
244 space_srv(en,Ks(class)+kentry) = space_srv(en,Ks(class)+kentry) + 1; % new job
245 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
246 if isinf(ni) % hit limited load-dependence
247 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*pentry(kentry)*mu{ist}{class}(k)*phi{ist}{class}(k)*ones(size(inspace(en,:),1),1)];
248 else
249 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*pentry(kentry)*mu{ist}{class}(k)*phi{ist}{class}(k)*ones(size(inspace(en,:),1),1)];
250 end
251 outprob = [outprob; ones(size(space_buf(en,:),1),1)];
252 end
253 end
254 case SchedStrategy.INF % move first job in service
255 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
256 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(en,class,k); % assume active
257 % if state is unchanged, still add with rate 0
258 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
259 if isinf(ni) % hit limited load-dependence
260 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
261 else
262 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
263 end
264 outprob = [outprob; ones(size(rate(en,:),1),1)];
265 case SchedStrategy.PS % move first job in service
266 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
267 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*(kir(en,class,k)./ni(en)).*min(ni(en),S(ist)); % assume active
268 % if state is unchanged, still add with rate 0
269 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
270 if isinf(ni) % hit limited load-dependence
271 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
272 else
273 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
274 end
275 outprob = [outprob; ones(size(rate(en,:),1),1)];
276 % end
277 case SchedStrategy.PSPRIO
278 % unclear if LD scaling should be with
279 % ni or with niprio, for now left as ni
280 % for consistency with HOL multiserver
281 if all(ni(en) <= S(ist))
282 % n <= c: all jobs get service, priority doesn't matter
283 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
284 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*(kir(en,class,k)./ni(en)).*min(ni(en),S(ist)); % assume active
285 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
286 if isinf(ni) % hit limited load-dependence
287 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
288 else
289 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
290 end
291 outprob = [outprob; ones(size(rate(en,:),1),1)];
292 elseif sn.classprio(class) == min(sn.classprio(nir>0)) % if this class is in the most urgent priority group (lower value = higher priority in LINE)
293 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
294 niprio(en) = sum(nir(sn.classprio==sn.classprio(class))); % jobs at the same priority class
295 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*(kir(en,class,k)./niprio(en)).*min(niprio(en),S(ist)); % assume active
296 % if state is unchanged, still add with rate 0
297 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
298 if isinf(ni) % hit limited load-dependence
299 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
300 else
301 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(niprio(en),lldlimit)).*rate(en,:)];
302 end
303 outprob = [outprob; ones(size(rate(en,:),1),1)];
304 else % n > c and not highest priority: set rate to zero
305 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
306 outrate = [outrate; zeros(size(rate(en,:),1),1)];
307 outprob = [outprob; ones(size(rate(en,:),1),1)];
308 end
309 case SchedStrategy.DPS
310 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
311 if S(ist) > 1
312 line_error(mfilename,'Multi-server DPS stations are not supported yet.');
313 end
314 % in GPS, the scheduling parameter are the weights
315 w_i = sn.schedparam(ist,:);
316 w_i = w_i / sum(w_i);
317 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k))*(kir(en,class,k)/nir(class))*w_i(class)*nir(class)./(sum(repmat(w_i,sum(en),1)*nir',2));
318 % if state is unchanged, still add with rate 0
319 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
320 if isinf(ni) % hit limited load-dependence
321 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
322 else
323 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
324 end
325 outprob = [outprob; ones(size(rate(en,:),1),1)];
326 case SchedStrategy.DPSPRIO
327 if all(ni(en) <= S(ist))
328 % n <= c: all jobs get service, behave like regular DPS
329 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
330 if S(ist) > 1
331 line_error(mfilename,'Multi-server DPS stations are not supported yet.');
332 end
333 w_i = sn.schedparam(ist,:);
334 w_i = w_i / sum(w_i);
335 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k))*(kir(en,class,k)/nir(class))*w_i(class)*nir(class)./(sum(repmat(w_i,sum(en),1)*nir',2));
336 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
337 if isinf(ni) % hit limited load-dependence
338 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
339 else
340 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
341 end
342 outprob = [outprob; ones(size(rate(en,:),1),1)];
343 elseif sn.classprio(class) == min(sn.classprio(nir>0)) % if this class is in the most urgent priority group (lower value = higher priority in LINE)
344 nirprio = nir;
345 nirprio(sn.classprio~=sn.classprio(class)) = 0; % ignore jobs of lower priority
346 niprio = sum(nirprio);
347 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
348 if S(ist) > 1
349 line_error(mfilename,'Multi-server DPS stations are not supported yet.');
350 end
351 % in GPS, the scheduling parameter are the weights
352 w_i = sn.schedparam(ist,:);
353 w_i = w_i / sum(w_i);
354 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k))*(kir(en,class,k)/nirprio(class))*w_i(class)*nirprio(class)./(sum(repmat(w_i,sum(en),1)*nirprio',2));
355 % if state is unchanged, still add with rate 0
356 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
357 if isinf(ni) % hit limited load-dependence
358 outrate = [outrate; cdscaling{ist}(nirprio).*lldscaling(ist,end).*rate(en,:)];
359 else
360 outrate = [outrate; cdscaling{ist}(nirprio).*lldscaling(ist,min(niprio(en),lldlimit)).*rate(en,:)];
361 end
362 outprob = [outprob; ones(size(rate(en,:),1),1)];
363 else % n > c and not most urgent priority: set rate to zero
364 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
365 outrate = [outrate; zeros(size(rate(en,:),1),1)];
366 outprob = [outprob; ones(size(rate(en,:),1),1)];
367 end
368 case SchedStrategy.GPS
369 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
370 if S(ist) > 1
371 line_error(mfilename,'Multi-server GPS stations are not supported yet.');
372 end
373 % in GPS, the scheduling parameter are the weights
374 w_i = sn.schedparam(ist,:);
375 w_i = w_i / sum(w_i);
376 cir = min(nir,ones(size(nir)));
377 rate = mu{ist}{class}(k)*(phi{ist}{class}(k))*(kir(en,class,k)/nir(class))*w_i(class)/(w_i*cir(:)); % assume active
378 % if state is unchanged, still add with rate 0
379 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
380 if isinf(ni) % hit limited load-dependence
381 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
382 else
383 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
384 end
385 outprob = [outprob; ones(size(rate(en,:),1),1)];
386 case SchedStrategy.GPSPRIO
387 if all(ni(en) <= S(ist))
388 % n <= c: all jobs get service, behave like regular GPS
389 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
390 if S(ist) > 1
391 line_error(mfilename,'Multi-server GPS stations are not supported yet.');
392 end
393 w_i = sn.schedparam(ist,:);
394 w_i = w_i / sum(w_i);
395 cir = min(nir,ones(size(nir)));
396 rate = mu{ist}{class}(k)*(phi{ist}{class}(k))*(kir(en,class,k)/nir(class))*w_i(class)/(w_i*cir(:)); % assume active
397 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
398 if isinf(ni) % hit limited load-dependence
399 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
400 else
401 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
402 end
403 outprob = [outprob; ones(size(rate(en,:),1),1)];
404 elseif sn.classprio(class) == min(sn.classprio(nir>0)) % if this class is in the most urgent priority group (lower value = higher priority in LINE)
405 nirprio = nir;
406 nirprio(sn.classprio~=sn.classprio(class)) = 0; % ignore jobs of lower priority
407 niprio = sum(nirprio);
408 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
409 if S(ist) > 1
410 line_error(mfilename,'Multi-server DPS stations are not supported yet.');
411 end
412 % in GPS, the scheduling parameter are the weights
413 w_i = sn.schedparam(ist,:);
414 w_i = w_i / sum(w_i);
415 cir = min(nirprio,ones(size(nirprio)));
416 rate = mu{ist}{class}(k)*(phi{ist}{class}(k))*(kir(en,class,k)/nirprio(class))*w_i(class)/(w_i*cir(:)); % assume active
417 % if state is unchanged, still add with rate 0
418 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
419 if isinf(ni) % hit limited load-dependence
420 outrate = [outrate; cdscaling{ist}(nirprio).*lldscaling(ist,end).*rate(en,:)];
421 else
422 outrate = [outrate; cdscaling{ist}(nirprio).*lldscaling(ist,min(niprio(en),lldlimit)).*rate(en,:)];
423 end
424 outprob = [outprob; ones(size(rate(en,:),1),1)];
425 else % n > c and not most urgent priority: set rate to zero
426 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
427 outrate = [outrate; zeros(size(rate(en,:),1),1)];
428 outprob = [outprob; ones(size(rate(en,:),1),1)];
429 end
430 case SchedStrategy.FCFS % move first job in service
431 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
432 en_wbuf = en & ni>S(ist); %states with jobs in buffer
433 for kdest=1:K(class) % new phase
434 space_buf_kd = space_buf;
435 space_var_kd = space_var;
436 if ismkvmodclass(class)
437 space_var_kd(en,sum(sn.nvars(ind,1:class))) = kdest;
438 end
439 rate_kd = rate;
440 rate_kd(en) = proc{ist}{class}{2}(k,kdest).*kir(en,class,k); % assume active
441 % first process states without jobs in buffer
442 en_wobuf = ~en_wbuf;
443 if any(en_wobuf) %any state without jobs in buffer
444 outspace = [outspace; space_buf_kd(en_wobuf,:), space_srv(en_wobuf,:), space_var_kd(en_wobuf,:)];
445 if isinf(ni) % hit limited load-dependence
446 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_kd(en_wobuf,:)];
447 else
448 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_kd(en_wobuf,:)];
449 end
450 end
451 % now process states with jobs in buffer
452 outprob = [outprob; ones(size(rate_kd(en_wobuf,:),1),1)];
453 if any(en_wbuf) %any state with jobs in buffer
454 % get class of job at head
455 start_svc_class = space_buf_kd(en_wbuf,end);
456 if start_svc_class > 0 % redunant if?
457 % update input buffer
458 space_buf_kd(en_wbuf,:) = [zeros(sum(en_wbuf),1),space_buf_kd(en_wbuf,1:end-1)];
459 % probability vector for the next job of starting in phase kentry
460 if ismkvmodclass(start_svc_class) % if markov-modulated
461 if start_svc_class==class % if successive service from the same class
462 kentry_range = kdest; % new job enters in phase left by departing job
463 else % resume phase from local variables
464 kentry_range = space_var_kd(en,sum(sn.nvars(ind,1:start_svc_class)));
465 end
466 pentry_svc_class = 0*pie{ist}{start_svc_class};
467 pentry_svc_class(kentry_range) = 1.0;
468 else % if i.i.d.
469 pentry_svc_class = pie{ist}{start_svc_class};
470 kentry_range = 1:K(start_svc_class);
471 end
472 for kentry = kentry_range
473 space_srv(en_wbuf,Ks(start_svc_class)+kentry) = space_srv(en_wbuf,Ks(start_svc_class)+kentry) + 1;
474 outspace = [outspace; space_buf_kd(en,:), space_srv(en,:), space_var_kd(en,:)];
475 rate_k = rate_kd;
476 rate_k(en_wbuf,:) = rate_kd(en_wbuf,:)*pentry_svc_class(kentry);
477 if isinf(ni) % use limited load-dependence at the latest user-provided level
478 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en,:)];
479 else
480 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en,:)];
481 end
482 outprob_cur = ones(size(rate_kd(en,:),1),1);
483 outprob_cur(outrate==0.0) = 0;
484 outprob = [outprob; outprob_cur'];
485 space_srv(en_wbuf,Ks(start_svc_class)+kentry) = space_srv(en_wbuf,Ks(start_svc_class)+kentry) - 1;
486 end
487 end
488 end
489 end
490 % if state is unchanged, still add with rate 0
491 case SchedStrategy.HOL % FCFS priority
492 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
493 en_wbuf = en & ni>S(ist); %states with jobs in buffer
494 en_wobuf = ~en_wbuf;
495 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
496 priogroup = [Inf,sn.classprio]; % Inf for empty positions (lower value = higher priority)
497 space_buf_groupg = arrayfun(@(x) priogroup(1+x), space_buf);
498 start_classprio = min(space_buf_groupg(en_wbuf,:),[],2); % min finds highest priority
499 isrowmax = space_buf_groupg == repmat(start_classprio, 1, size(space_buf_groupg,2));
500 [~,rightmostMaxPosFlipped]=max(fliplr(isrowmax),[],2);
501 rightmostMaxPos = size(isrowmax,2) - rightmostMaxPosFlipped + 1;
502 start_svc_class = space_buf(en_wbuf, rightmostMaxPos);
503 outspace = [outspace; space_buf(en_wobuf,:), space_srv(en_wobuf,:), space_var(en_wobuf,:)];
504 if isinf(ni) % hit limited load-dependence
505 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en_wobuf,:)];
506 else
507 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en_wobuf,:)];
508 end
509 outprob = [outprob; ones(size(rate(en_wobuf,:),1),1)];
510 if start_svc_class > 0
511 pentry_svc_class = pie{ist}{start_svc_class};
512 for kentry = 1:K(start_svc_class)
513 space_srv_k = space_srv;
514 space_buf_k = space_buf;
515 space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) = space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) + 1;
516 for j=find(en_wbuf)'
517 space_buf_k(j,:) = [0, space_buf_k(j,1:rightmostMaxPos(j)-1), space_buf_k(j,(rightmostMaxPos(j)+1):end)];
518 end
519 % if state is unchanged, still add with rate 0
520 outspace = [outspace; space_buf_k(en_wbuf,:), space_srv_k(en_wbuf,:), space_var(en_wbuf,:)];
521 rate_k = rate;
522 rate_k(en_wbuf,:) = rate(en_wbuf,:) * pentry_svc_class(kentry);
523 if isinf(ni) % hit limited load-dependence
524 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en_wbuf,:)];
525 else
526 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en_wbuf,:)];
527 end
528 outprob = [outprob; ones(size(rate_k(en_wbuf,:),1),1)];
529 end
530 end
531 case SchedStrategy.LCFSPRIO % LCFS priority - like HOL but LCFS order within priority groups
532 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
533 en_wbuf = en & ni>S(ist); %states with jobs in buffer
534 en_wobuf = ~en_wbuf;
535 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
536 priogroup = [Inf,sn.classprio]; % Inf for empty positions (lower value = higher priority)
537 space_buf_groupg = arrayfun(@(x) priogroup(1+x), space_buf);
538 start_classprio = min(space_buf_groupg(en_wbuf,:),[],2); % min finds highest priority
539 isrowmax = space_buf_groupg == repmat(start_classprio, 1, size(space_buf_groupg,2));
540 % LCFS: Find leftmost (first) position instead of rightmost for LCFS order
541 [~,leftmostMaxPos]=max(isrowmax,[],2);
542 start_svc_class = space_buf(en_wbuf, leftmostMaxPos);
543 outspace = [outspace; space_buf(en_wobuf,:), space_srv(en_wobuf,:), space_var(en_wobuf,:)];
544 if isinf(ni) % hit limited load-dependence
545 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en_wobuf,:)];
546 else
547 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en_wobuf,:)];
548 end
549 outprob = [outprob; ones(size(rate(en_wobuf,:),1),1)];
550 if start_svc_class > 0
551 pentry_svc_class = pie{ist}{start_svc_class};
552 for kentry = 1:K(start_svc_class)
553 space_srv_k = space_srv;
554 space_buf_k = space_buf;
555 space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) = space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) + 1;
556 for j=find(en_wbuf)'
557 % LCFS: Remove from leftmost position instead of rightmost
558 space_buf_k(j,:) = [0, space_buf_k(j,1:leftmostMaxPos(j)-1), space_buf_k(j,(leftmostMaxPos(j)+1):end)];
559 end
560 % if state is unchanged, still add with rate 0
561 outspace = [outspace; space_buf_k(en_wbuf,:), space_srv_k(en_wbuf,:), space_var(en_wbuf,:)];
562 rate_k = rate;
563 rate_k(en_wbuf,:) = rate_k(en_wbuf,:)*pentry_svc_class(kentry);
564 if isinf(ni) % hit limited load-dependence
565 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en_wbuf,:)];
566 else
567 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en_wbuf,:)];
568 end
569 outprob = [outprob; ones(size(rate_k(en_wbuf,:),1),1)];
570 space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) = space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) - 1;
571 end
572 end
573 case SchedStrategy.LCFS % move last job in service
574 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
575 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
576 en_wbuf = en & ni>S(ist); %states with jobs in buffer
577 [~, colfirstnnz] = max( space_buf(en_wbuf,:) ~=0, [], 2 ); % find first nnz column
578 start_svc_class = space_buf(en_wbuf,colfirstnnz); % job entering service
579 space_buf(en_wbuf,colfirstnnz)=0;
580 if isempty(start_svc_class)
581 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
582 if isinf(ni) % hit limited load-dependence
583 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
584 else
585 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
586 end
587 outprob = [outprob; ones(size(rate(en,:),1),1)];
588 if isSimulation && nargin>=7 && isobject(eventCache)
589 eventCache(key) = {outprob, outspace,outrate};
590 end
591 return
592 end
593 for kentry = 1:K(start_svc_class)
594 pentry_svc_class = pie{ist}{start_svc_class};
595 space_srv(en_wbuf,Ks(start_svc_class)+kentry) = space_srv(en_wbuf,Ks(start_svc_class)+kentry) + 1;
596 % if state is unchanged, still add with rate 0
597 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
598 rate_k = rate;
599 rate_k(en_wbuf,:) = rate(en_wbuf,:)*pentry_svc_class(kentry);
600 if isinf(ni) % hit limited load-dependence
601 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en,:)];
602 else
603 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en,:)];
604 end
605 outprob = [outprob; ones(size(rate(en,:),1),1)];
606 space_srv(en_wbuf,Ks(start_svc_class)+kentry) = space_srv(en_wbuf,Ks(start_svc_class)+kentry) - 1;
607 end
608 case SchedStrategy.LCFSPR % move last job in service (preempt-resume)
609 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
610 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
611 en_wbuf = en & ni>S(ist); %states with jobs in buffer
612 [~, colfirstnnz] = max( space_buf(en_wbuf,:) ~=0, [], 2 ); % find first nnz column
613 start_svc_class = space_buf(en_wbuf,colfirstnnz); % job entering service
614 kentry = space_buf(en_wbuf,colfirstnnz+1); % entry phase of job resuming service
615 space_buf(en_wbuf,colfirstnnz)=0;% zero popped job
616 space_buf(en_wbuf,colfirstnnz+1)=0; % zero popped phase
617 if isempty(start_svc_class)
618 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
619 if isinf(ni) % hit limited load-dependence
620 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
621 else
622 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
623 end
624 outprob = [outprob; ones(size(rate(en,:),1),1)];
625 if isSimulation && nargin>=7 && isobject(eventCache)
626 eventCache(key) = {outprob, outspace,outrate};
627 end
628 return
629 end
630 space_srv(en_wbuf,Ks(start_svc_class)+kentry) = space_srv(en_wbuf,Ks(start_svc_class)+kentry) + 1;
631 % if state is unchanged, still add with rate 0
632 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
633 if isinf(ni) % hit limited load-dependence
634 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
635 else
636 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
637 end
638 outprob = [outprob; ones(size(rate(en,:),1),1)];
639 case SchedStrategy.LCFSPI % move last job in service (preempt-independent)
640 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
641 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
642 en_wbuf = en & ni>S(ist); %states with jobs in buffer
643 [~, colfirstnnz] = max( space_buf(en_wbuf,:) ~=0, [], 2 ); % find first nnz column
644 start_svc_class = space_buf(en_wbuf,colfirstnnz); % job entering service
645 space_buf(en_wbuf,colfirstnnz)=0;% zero popped job
646 space_buf(en_wbuf,colfirstnnz+1)=0; % zero popped phase (ignored for LCFSPI)
647 if isempty(start_svc_class)
648 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
649 if isinf(ni) % hit limited load-dependence
650 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
651 else
652 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
653 end
654 outprob = [outprob; ones(size(rate(en,:),1),1)];
655 if isSimulation && nargin>=7 && isobject(eventCache)
656 eventCache(key) = {outprob, outspace,outrate};
657 end
658 return
659 end
660 % For LCFSPI, jobs restart from pie distribution instead of stored phase
661 pentry_svc_class = pie{ist}{start_svc_class};
662 for kentry = 1:K(start_svc_class)
663 space_srv_k = space_srv;
664 space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) = space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) + 1;
665 % if state is unchanged, still add with rate 0
666 outspace = [outspace; space_buf(en,:), space_srv_k(en,:), space_var(en,:)];
667 rate_k = rate;
668 rate_k(en_wbuf,:) = rate_k(en_wbuf,:) * pentry_svc_class(kentry);
669 if isinf(ni) % hit limited load-dependence
670 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en,:)];
671 else
672 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en,:)];
673 end
674 outprob = [outprob; ones(size(rate_k(en,:),1),1)];
675 end
676 case SchedStrategy.FCFSPR % FCFS preempt-resume (no priority)
677 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
678 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
679 en_wbuf = en & ni>S(ist); %states with jobs in buffer
680 % FCFS: Find rightmost (last) non-zero position (longest waiting job)
681 [~, colLastNnz] = max(fliplr(space_buf(en_wbuf,:) ~= 0), [], 2);
682 colLastNnz = size(space_buf,2) - colLastNnz; % convert to actual position (odd index for class)
683 start_svc_class = space_buf(en_wbuf, colLastNnz); % job entering service
684 kentry = space_buf(en_wbuf, colLastNnz+1); % entry phase of job resuming service
685 space_buf(en_wbuf, colLastNnz) = 0; % zero popped job
686 space_buf(en_wbuf, colLastNnz+1) = 0; % zero popped phase
687 if isempty(start_svc_class)
688 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
689 if isinf(ni)
690 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
691 else
692 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
693 end
694 outprob = [outprob; ones(size(rate(en,:),1),1)];
695 if isSimulation && nargin>=7 && isobject(eventCache)
696 eventCache(key) = {outprob, outspace, outrate};
697 end
698 return
699 end
700 space_srv(en_wbuf, Ks(start_svc_class)+kentry) = space_srv(en_wbuf, Ks(start_svc_class)+kentry) + 1;
701 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
702 if isinf(ni)
703 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
704 else
705 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
706 end
707 outprob = [outprob; ones(size(rate(en,:),1),1)];
708 case SchedStrategy.FCFSPI % FCFS preempt-independent (no priority)
709 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
710 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
711 en_wbuf = en & ni>S(ist); %states with jobs in buffer
712 % FCFS: Find rightmost (last) non-zero position (longest waiting job)
713 [~, colLastNnz] = max(fliplr(space_buf(en_wbuf,:) ~= 0), [], 2);
714 colLastNnz = size(space_buf,2) - colLastNnz; % convert to actual position (odd index for class)
715 start_svc_class = space_buf(en_wbuf, colLastNnz); % job entering service
716 space_buf(en_wbuf, colLastNnz) = 0; % zero popped job
717 space_buf(en_wbuf, colLastNnz+1) = 0; % zero popped phase (ignored for FCFSPI)
718 if isempty(start_svc_class)
719 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
720 if isinf(ni)
721 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en,:)];
722 else
723 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en,:)];
724 end
725 outprob = [outprob; ones(size(rate(en,:),1),1)];
726 if isSimulation && nargin>=7 && isobject(eventCache)
727 eventCache(key) = {outprob, outspace, outrate};
728 end
729 return
730 end
731 % For FCFSPI, jobs restart from pie distribution instead of stored phase
732 pentry_svc_class = pie{ist}{start_svc_class};
733 for kentry = 1:K(start_svc_class)
734 space_srv_k = space_srv;
735 space_srv_k(en_wbuf, Ks(start_svc_class)+kentry) = space_srv_k(en_wbuf, Ks(start_svc_class)+kentry) + 1;
736 outspace = [outspace; space_buf(en,:), space_srv_k(en,:), space_var(en,:)];
737 rate_k = rate;
738 rate_k(en_wbuf,:) = rate_k(en_wbuf,:) * pentry_svc_class(kentry);
739 if isinf(ni)
740 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en,:)];
741 else
742 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en,:)];
743 end
744 outprob = [outprob; ones(size(rate_k(en,:),1),1)];
745 end
746 case SchedStrategy.FCFSPRPRIO % FCFS preempt-resume with priority groups
747 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
748 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
749 en_wbuf = en & ni>S(ist); %states with jobs in buffer
750 en_wobuf = ~en_wbuf;
751 priogroup = [Inf,sn.classprio]; % Inf for empty positions (lower value = higher priority)
752 space_buf_groupg = arrayfun(@(x) priogroup(1+x), space_buf);
753 start_classprio = min(space_buf_groupg(en_wbuf,:),[],2); % min finds highest priority
754 isrowmax = space_buf_groupg == repmat(start_classprio, 1, size(space_buf_groupg,2));
755 % FCFS: Find rightmost (last) position for highest priority class (longest waiting)
756 [~,rightmostMaxPos]=max(fliplr(isrowmax),[],2);
757 rightmostMaxPos = size(space_buf_groupg,2) - rightmostMaxPos + 1;
758 start_svc_class = space_buf(en_wbuf, rightmostMaxPos); % job entering service
759 kentry = space_buf(en_wbuf, rightmostMaxPos+1); % entry phase of job resuming service (preempt-resume)
760
761 % Handle states without buffer jobs
762 outspace = [outspace; space_buf(en_wobuf,:), space_srv(en_wobuf,:), space_var(en_wobuf,:)];
763 if isinf(ni)
764 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en_wobuf,:)];
765 else
766 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en_wobuf,:)];
767 end
768 outprob = [outprob; ones(size(rate(en_wobuf,:),1),1)];
769
770 % Handle states with buffer jobs
771 if any(en_wbuf) && start_svc_class > 0
772 space_srv_k = space_srv;
773 space_buf_k = space_buf;
774 space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) = space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) + 1;
775 for j=find(en_wbuf)'
776 % Remove both class and phase from rightmost position (preempt-resume)
777 space_buf_k(j,:) = [0, space_buf_k(j,1:rightmostMaxPos(j)-1), space_buf_k(j,(rightmostMaxPos(j)+2):end), 0];
778 end
779 outspace = [outspace; space_buf_k(en_wbuf,:), space_srv_k(en_wbuf,:), space_var(en_wbuf,:)];
780 if isinf(ni)
781 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en_wbuf,:)];
782 else
783 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en_wbuf,:)];
784 end
785 outprob = [outprob; ones(size(rate(en_wbuf,:),1),1)];
786 end
787 case SchedStrategy.LCFSPRPRIO % LCFS preempt-resume with priority groups
788 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
789 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
790 en_wbuf = en & ni>S(ist); %states with jobs in buffer
791 en_wobuf = ~en_wbuf;
792 priogroup = [Inf,sn.classprio]; % Inf for empty positions (lower value = higher priority)
793 space_buf_groupg = arrayfun(@(x) priogroup(1+x), space_buf);
794 start_classprio = min(space_buf_groupg(en_wbuf,:),[],2); % min finds highest priority
795 isrowmax = space_buf_groupg == repmat(start_classprio, 1, size(space_buf_groupg,2));
796 % LCFS: Find leftmost (first) position for highest priority class (most recently preempted)
797 [~,leftmostMaxPos]=max(isrowmax,[],2);
798 start_svc_class = space_buf(en_wbuf, leftmostMaxPos); % job entering service
799 kentry = space_buf(en_wbuf, leftmostMaxPos+1); % entry phase of job resuming service (preempt-resume)
800
801 % Handle states without buffer jobs
802 outspace = [outspace; space_buf(en_wobuf,:), space_srv(en_wobuf,:), space_var(en_wobuf,:)];
803 if isinf(ni)
804 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en_wobuf,:)];
805 else
806 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en_wobuf,:)];
807 end
808 outprob = [outprob; ones(size(rate(en_wobuf,:),1),1)];
809
810 % Handle states with buffer jobs
811 if any(en_wbuf) && start_svc_class > 0
812 space_srv_k = space_srv;
813 space_buf_k = space_buf;
814 space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) = space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) + 1;
815 for j=find(en_wbuf)'
816 % Remove both class and phase from leftmost position (preempt-resume)
817 space_buf_k(j,:) = [0, space_buf_k(j,1:leftmostMaxPos(j)-1), space_buf_k(j,(leftmostMaxPos(j)+2):end), 0];
818 end
819 outspace = [outspace; space_buf_k(en_wbuf,:), space_srv_k(en_wbuf,:), space_var(en_wbuf,:)];
820 if isinf(ni)
821 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en_wbuf,:)];
822 else
823 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en_wbuf,:)];
824 end
825 outprob = [outprob; ones(size(rate(en_wbuf,:),1),1)];
826 end
827 case SchedStrategy.FCFSPIPRIO % FCFS preempt-independent with priority groups
828 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
829 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
830 en_wbuf = en & ni>S(ist); %states with jobs in buffer
831 en_wobuf = ~en_wbuf;
832 priogroup = [Inf,sn.classprio]; % Inf for empty positions (lower value = higher priority)
833 space_buf_groupg = arrayfun(@(x) priogroup(1+x), space_buf);
834 start_classprio = min(space_buf_groupg(en_wbuf,:),[],2); % min finds highest priority
835 isrowmax = space_buf_groupg == repmat(start_classprio, 1, size(space_buf_groupg,2));
836 % FCFS: Find rightmost (last) position for highest priority class
837 [~,rightmostMaxPos]=max(fliplr(isrowmax),[],2);
838 rightmostMaxPos = size(space_buf_groupg,2) - rightmostMaxPos + 1;
839 start_svc_class = space_buf(en_wbuf, rightmostMaxPos); % job entering service
840
841 % Handle states without buffer jobs
842 outspace = [outspace; space_buf(en_wobuf,:), space_srv(en_wobuf,:), space_var(en_wobuf,:)];
843 if isinf(ni)
844 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en_wobuf,:)];
845 else
846 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en_wobuf,:)];
847 end
848 outprob = [outprob; ones(size(rate(en_wobuf,:),1),1)];
849
850 % Handle states with buffer jobs
851 if any(en_wbuf) && start_svc_class > 0
852 % For FCFSPIPRIO, jobs restart from pie distribution instead of stored phase
853 pentry_svc_class = pie{ist}{start_svc_class};
854 for kentry = 1:K(start_svc_class)
855 space_srv_k = space_srv;
856 space_buf_k = space_buf;
857 space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) = space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) + 1;
858 for j=find(en_wbuf)'
859 % Remove both class and phase from rightmost position (preempt-independent ignores stored phase)
860 space_buf_k(j,:) = [0, space_buf_k(j,1:rightmostMaxPos(j)-1), space_buf_k(j,(rightmostMaxPos(j)+2):end), 0];
861 end
862 outspace = [outspace; space_buf_k(en_wbuf,:), space_srv_k(en_wbuf,:), space_var(en_wbuf,:)];
863 rate_k = rate;
864 rate_k(en_wbuf,:) = rate_k(en_wbuf,:) * pentry_svc_class(kentry);
865 if isinf(ni)
866 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en_wbuf,:)];
867 else
868 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en_wbuf,:)];
869 end
870 outprob = [outprob; ones(size(rate_k(en_wbuf,:),1),1)];
871 end
872 end
873 case SchedStrategy.LCFSPIPRIO % LCFS preempt-independent with priority groups
874 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
875 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % assume active
876 en_wbuf = en & ni>S(ist); %states with jobs in buffer
877 en_wobuf = ~en_wbuf;
878 priogroup = [Inf,sn.classprio]; % Inf for empty positions (lower value = higher priority)
879 space_buf_groupg = arrayfun(@(x) priogroup(1+x), space_buf);
880 start_classprio = min(space_buf_groupg(en_wbuf,:),[],2); % min finds highest priority
881 isrowmax = space_buf_groupg == repmat(start_classprio, 1, size(space_buf_groupg,2));
882 % LCFS: Find leftmost (first) position for highest priority class
883 [~,leftmostMaxPos]=max(isrowmax,[],2);
884 start_svc_class = space_buf(en_wbuf, leftmostMaxPos); % job entering service
885
886 % Handle states without buffer jobs
887 outspace = [outspace; space_buf(en_wobuf,:), space_srv(en_wobuf,:), space_var(en_wobuf,:)];
888 if isinf(ni)
889 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en_wobuf,:)];
890 else
891 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en_wobuf,:)];
892 end
893 outprob = [outprob; ones(size(rate(en_wobuf,:),1),1)];
894
895 % Handle states with buffer jobs
896 if any(en_wbuf) && start_svc_class > 0
897 % For LCFSPIPRIO, jobs restart from pie distribution instead of stored phase
898 pentry_svc_class = pie{ist}{start_svc_class};
899 for kentry = 1:K(start_svc_class)
900 space_srv_k = space_srv;
901 space_buf_k = space_buf;
902 space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) = space_srv_k(en_wbuf,Ks(start_svc_class)+kentry) + 1;
903 for j=find(en_wbuf)'
904 % Remove both class and phase from leftmost position (preempt-independent ignores stored phase)
905 space_buf_k(j,:) = [0, space_buf_k(j,1:leftmostMaxPos(j)-1), space_buf_k(j,(leftmostMaxPos(j)+2):end), 0];
906 end
907 outspace = [outspace; space_buf_k(en_wbuf,:), space_srv_k(en_wbuf,:), space_var(en_wbuf,:)];
908 rate_k = rate;
909 rate_k(en_wbuf,:) = rate_k(en_wbuf,:) * pentry_svc_class(kentry);
910 if isinf(ni)
911 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en_wbuf,:)];
912 else
913 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en_wbuf,:)];
914 end
915 outprob = [outprob; ones(size(rate_k(en_wbuf,:),1),1)];
916 end
917 end
918 case SchedStrategy.SIRO
919 rate = zeros(size(space_srv,1),1);
920 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % this is for states not in en_buf
921 space_srv = inspace(:,end-sum(K)+1:end); % server state
922 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
923 % first record departure in states where the buffer is empty
924 en_wobuf = en & sum(space_buf(en,:),2) == 0;
925 outspace = [outspace; space_buf(en_wobuf,:), space_srv(en_wobuf,:), space_var(en_wobuf,:)];
926 if isinf(ni) % hit limited load-dependence
927 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate(en_wobuf,:)];
928 else
929 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate(en_wobuf,:)];
930 end
931 outprob = [outprob; ones(size(rate(en_wobuf,:),1),1)];
932 % let's go now to states where the buffer is non-empty
933 for r=1:R % pick a job of a random class
934 rate_r = rate;
935 space_buf = inspace(:,1:(end-sum(K))); % buffer state
936 en_wbuf = en & space_buf(en,r) > 0; % states where the buffer is non-empty
937 space_buf(en_wbuf,r) = space_buf(en_wbuf,r) - 1; % remove from buffer
938 space_srv_r = space_srv;
939 pentry_svc_class = pie{ist}{r};
940 pick_prob = (nir(r)-sir(r)) / (ni-sum(sir));
941 if pick_prob >= 0
942 rate_r(en_wbuf,:) = rate_r(en_wbuf,:) * pick_prob;
943 end
944 for kentry=1:K(r)
945 space_srv_r(en_wbuf,Ks(r)+kentry) = space_srv_r(en_wbuf,Ks(r)+kentry) + 1; % bring job in service
946 outspace = [outspace; space_buf(en_wbuf,:), space_srv_r(en_wbuf,:), space_var(en_wbuf,:)];
947 rate_k = rate_r;
948 rate_k(en_wbuf,:) = rate_k(en_wbuf,:) * pentry_svc_class(kentry);
949 if isinf(ni) % hit limited load-dependence
950 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en_wbuf,:)];
951 else
952 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en_wbuf,:)];
953 end
954 outprob = [outprob; ones(size(rate(en_wbuf,:),1),1)];
955 space_srv_r(en_wbuf,Ks(r)+kentry) = space_srv_r(en_wbuf,Ks(r)+kentry) - 1; % bring job in service
956 end
957 end
958 case {SchedStrategy.SEPT,SchedStrategy.LEPT} % move last job in service
959 rate = zeros(size(space_srv,1),1);
960 rate(en) = mu{ist}{class}(k)*(phi{ist}{class}(k)).*kir(:,class,k); % this is for states not in en_buf
961 space_srv = inspace(:,end-sum(K)+1:end); % server state
962 space_srv(en,Ks(class)+k) = space_srv(en,Ks(class)+k) - 1; % record departure
963 space_buf = inspace(:,1:(end-sum(K))); % buffer state
964 % in SEPT, the scheduling parameter is the priority order of the class means
965 % en_wbuf: states where the buffer is non-empty
966 % sept_class: class to pick in service
967 [en_wbuf, first_class_inrow] = max(space_buf(:,sn.schedparam(ist,:))~=0, [], 2);
968 sept_class = sn.schedparam(ist,first_class_inrow); % this is different for sept and lept
969
970 space_buf(en_wbuf,sept_class) = space_buf(en_wbuf,sept_class) - 1; % remove from buffer
971 pentry = pie{ist}{sept_class};
972 for kentry=1:K(sept_class)
973 space_srv(en_wbuf,Ks(sept_class)+kentry) = space_srv(en_wbuf,Ks(sept_class)+kentry) + 1; % bring job in service
974 if isSimulation
975 % break the tie
976 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
977 rate_k = rate;
978 rate_k(en,:) = rate_k(en,:) * pentry(kentry);
979 if isinf(ni) % hit limited load-dependence
980 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en,:)];
981 else
982 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en,:)];
983 end
984 outprob = [outprob; ones(size(rate(en,:),1),1)];
985 else
986 outspace = [outspace; space_buf(en,:), space_srv(en,:), space_var(en,:)];
987 rate_k = rate;
988 rate_k(en,:) = rate_k(en,:) * pentry(kentry);
989 if isinf(ni) % hit limited load-dependence
990 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate_k(en,:)];
991 else
992 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate_k(en,:)];
993 end
994 outprob = [outprob; ones(size(rate(en,:),1),1)];
995 end
996 space_srv(en_wbuf,Ks(sept_class)+kentry) = space_srv(en_wbuf,Ks(sept_class)+kentry) - 1; % bring job in service
997 end
998 otherwise
999 line_error(mfilename,sprintf('Scheduling strategy %s is not supported.', SchedStrategy.toText(sn.sched(ist))));
1000 end
1001 end
1002 end
1003 if isSimulation
1004 if nargin>=7 && isobject(eventCache)
1005 eventCache(key) = {outprob, outspace,outrate};
1006 end
1007
1008 if size(outspace,1) > 1
1009 tot_rate = sum(outrate);
1010 cum_rate = cumsum(outrate) / tot_rate;
1011 firing_ctr = 1 + max([0,find( rand > cum_rate' )]); % select action
1012 outspace = outspace(firing_ctr,:);
1013 outrate = sum(outrate);
1014 outprob = outprob(firing_ctr,:);
1015 end
1016 end
1017 end
1018 end
1019 case EventType.PHASE
1020 outspace = [];
1021 outrate = [];
1022 outprob = [];
1023 [ni,nir,~,kir] = State.toMarginal(sn,ind,inspace,K,Ks,space_buf,space_srv,space_var);
1024 if nir(class)>0
1025 for k=1:K(class)
1026 en = space_srv(:,Ks(class)+k) > 0;
1027 if any(en)
1028 for kdest=setdiff(1:K(class),k) % new phase
1029 rate = 0;
1030 space_srv_k = space_srv(en,:);
1031 space_buf_k = space_buf(en,:);
1032 space_var_k = space_var(en,:);
1033 if ismkvmodclass(class)
1034 space_var_k(sum(sn.nvars(ind,1:class))) = kdest;
1035 end
1036 space_srv_k(:,Ks(class)+k) = space_srv_k(:,Ks(class)+k) - 1;
1037 space_srv_k(:,Ks(class)+kdest) = space_srv_k(:,Ks(class)+kdest) + 1;
1038 switch sn.sched(ist)
1039 case SchedStrategy.EXT
1040 rate = proc{ist}{class}{1}(k,kdest); % move next job forward
1041 case SchedStrategy.INF
1042 rate = proc{ist}{class}{1}(k,kdest)*kir(:,class,k); % assume active
1043 case {SchedStrategy.PS, SchedStrategy.LPS}
1044 rate = proc{ist}{class}{1}(k,kdest)*kir(:,class,k)./ni(:).*min(ni(:),S(ist)); % assume active
1045 case SchedStrategy.PSPRIO
1046 if all(ni <= S(ist)) || sn.classprio(class) == min(sn.classprio(nir>0))
1047 rate = proc{ist}{class}{1}(k,kdest)*kir(:,class,k)./ni(:).*min(ni(:),S(ist)); % assume active
1048 else
1049 rate = 0; % not in most urgent priority group
1050 end
1051 case SchedStrategy.DPSPRIO
1052 if all(ni <= S(ist)) || sn.classprio(class) == min(sn.classprio(nir>0))
1053 w_i = sn.schedparam(ist,:);
1054 w_i = w_i / sum(w_i);
1055 nirprio = nir;
1056 if ~all(ni <= S(ist))
1057 nirprio(sn.classprio~=sn.classprio(class)) = 0;
1058 end
1059 rate = proc{ist}{class}{1}(k,kdest)*kir(:,class,k)*w_i(class)./(sum(repmat(w_i,size(nirprio,1),1)*nirprio',2)); % assume active
1060 else
1061 rate = 0; % not in most urgent priority group
1062 end
1063 case SchedStrategy.GPSPRIO
1064 if all(ni <= S(ist)) || sn.classprio(class) == min(sn.classprio(nir>0))
1065 w_i = sn.schedparam(ist,:);
1066 w_i = w_i / sum(w_i);
1067 nirprio = nir;
1068 if ~all(ni <= S(ist))
1069 nirprio(sn.classprio~=sn.classprio(class)) = 0;
1070 end
1071 cir = min(nirprio,ones(size(nirprio)));
1072 rate = proc{ist}{class}{1}(k,kdest)*kir(:,class,k)/nirprio(class)*w_i(class)/(w_i*cir(:)); % assume active
1073 else
1074 rate = 0; % not in most urgent priority group
1075 end
1076 case SchedStrategy.DPS
1077 if S(ist) > 1
1078 line_error(mfilename,'Multi-server DPS not supported yet');
1079 end
1080 w_i = sn.schedparam(ist,:);
1081 w_i = w_i / sum(w_i);
1082 rate = proc{ist}{class}{1}(k,kdest)*kir(:,class,k)*w_i(class)./(sum(repmat(w_i,size(nir,1),1)*nir',2)); % assume active
1083 case SchedStrategy.GPS
1084 if S(ist) > 1
1085 line_error(mfilename,'Multi-server GPS not supported yet');
1086 end
1087 cir = min(nir,ones(size(nir)));
1088 w_i = sn.schedparam(ist,:); w_i = w_i / sum(w_i);
1089 rate = proc{ist}{class}{1}(k,kdest)*kir(:,class,k)/nir(class)*w_i(class)/(w_i*cir(:)); % assume active
1090
1091 case {SchedStrategy.FCFS, SchedStrategy.HOL, SchedStrategy.LCFS, SchedStrategy.LCFSPR, SchedStrategy.LCFSPI, SchedStrategy.LCFSPRIO, SchedStrategy.SIRO, SchedStrategy.SEPT, SchedStrategy.LEPT}
1092 rate = proc{ist}{class}{1}(k,kdest)*kir(:,class,k); % assume active
1093 end
1094 % if the class cannot be served locally,
1095 % then rate = NaN since mu{i,class}=NaN
1096 if isinf(ni) % hit limited load-dependence
1097 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,end).*rate];
1098 else
1099 outrate = [outrate; cdscaling{ist}(nir).*lldscaling(ist,min(ni(en),lldlimit)).*rate];
1100 end
1101 outspace = [outspace; space_buf_k, space_srv_k, space_var_k];
1102 outprob = [outprob; ones(size(rate,1),1)];
1103 end
1104 end
1105 end
1106 if isSimulation
1107 if nargin>=7 && isobject(eventCache)
1108 eventCache(key) = {outprob, outspace,outrate};
1109 end
1110
1111 if size(outspace,1) > 1
1112 tot_rate = sum(outrate);
1113 cum_rate = cumsum(outrate) / tot_rate;
1114 firing_ctr = 1 + max([0,find( rand > cum_rate' )]); % select action
1115 outspace = outspace(firing_ctr,:);
1116 outrate = sum(outrate);
1117 outprob = outprob(firing_ctr,:);
1118 end
1119 end
1120 end
1121end
1122
1123end
Definition mmt.m:92