C ---------------------------------------------------------------------- INTEGER FUNCTION PENALTYCAL(sim,net) C ---------------------------------------------------------------------- C - Top level routine to update the link to link movement penalties. C - INCLUDED FILES: #include "dyna.inc" #include "sim.inc" #include "network.inc" C - UNMODIFIED ARGUMENTS: RECORD /Sim_Data/ sim C - MODIFIED ARGUMENTS: RECORD /Road_Network/ net C - MODIFIED GLOBAL DATA: ! NONE C - LOCAL VARIABLES: INTEGER i C - FUNCTIONS CALLED: INTEGER LINK_PENALTYCAL !function C - RETURN VALUE: ! Simulation status C ---------------------------------------------------------------------- i = 1 DO WHILE (i.LE.net.nlinks + .AND. + LINK_PENALTYCAL(sim,net,i).EQ.SIM_IN_PROGRESS) i = i + 1 ENDDO PENALTYCAL = sim.status RETURN END C ---------------------------------------------------------------------- INTEGER FUNCTION LINK_PENALTYCAL(sim,net,i) C ---------------------------------------------------------------------- C - Determines the average movement delay associated with each C - movement eminating from link i. The delays are calculated using a C - rolling horizon average of the delays experienced by vehicles C - making the movement during the preceding NU_DS timesteps. The C - average delay is given in minutes. C - INCLUDED FILES: #include "dyna.inc" #include "sim.inc" #include "network.inc" C - UNMODIFIED ARGUMENTS: RECORD /Sim_Data/ sim C - MODIFIED ARGUMENTS: RECORD /Road_Network/ net C - MODIFIED GLOBAL DATA: ! none C - LOCAL VARIABLES: REAL turnveh(MAX_MOVE_TYPE) INTEGER i,im,is,isc,nodeup,ip,il,ili,itmp,mtmp,kk,itmp2,j + ,mv REAL c,nl,s,v,vs,gc,vc REAL dsum REAL tmp C - FUNCTIONS CALLED: ! NONE C - RETURN VALUE: ! Simulation status C ---------------------------------------------------------------------- Cr - Loop over all movements feeding from this link DO im = 1,net.link(i).numuslinks mv = net.link(i).usmoveptr(im) Cr - Loop over horizon to sum recent delays dsum = 0.D0 ! sum of delays for this movement isc = 0 ! # of periods with vehicle delay stored DO is = 1,NU_DS IF (net.movement(mv).delaytime(is).GT.0.D0) THEN tmp = MOD(INT(net.movement(mv).delaytime(is)),10000) dsum = dsum + tmp isc = isc + INT(net.movement(mv).delaytime(is)/10000) ENDIF ENDDO IF (isc.EQ.0) THEN Cr - No vehicles exited during the rolling horizon so no delay Cr - estimates are available from simulated vehicles. Use Cr - expected waiting time for control type cr turnveh(im) = net.movement(mv).queued*sim.timestep nodeup=net.link(i).iunod ip =net.node(nodeup).curplan IF ((net.node(nodeup).plan(ip).type.EQ.CTL_PRETIMED + .OR. + net.node(nodeup).plan(ip).type.EQ.CTL_ACTUATED) + .AND.net.movement(mv).type.NE.CONNECTOR) THEN c = net.node(nodeup).plan(ip).cycle nl = net.movement(mv).nlanes s = nl*1800 ! total sat flow = 1800*nlanes: HARDCODED IF (c.EQ.0.OR.s.EQ.0) THEN turnveh(im) = INFINITY/2 ELSE v = net.movement(mv).totmoved/ + (sim.time.minutes+sim.timestep)*60.0 vs = v/s gc = net.movement(mv).cycgreen/c IF (gc.EQ.0) THEN turnveh(im) = INFINITY/2.0 ELSE vc = vs/gc IF (vc.GE.1) THEN ! Overflow turnveh(im) = INFINITY/2 ELSE turnveh(im) = + c*(1-gc**2)/(2*(1-vs)) IF (v.GT.0) turnveh(im) = turnveh(im) + + vc**2/(2*v*(1-vc)) turnveh(im) = turnveh(im)*0.90/60.0 ENDIF ENDIF ENDIF ELSE IF (net.node(nodeup).plan(ip).type.EQ.CTL_STOP + .AND. + net.movement(mv).type.NE.CONNECTOR) THEN CR - Make the *extremely weak* assumption that people will CR - assess their stop-controlled intersection delay as CR - 1/2 of a minute crxxxx turnveh(im) = 0.50 CR - Use the 1994 HCM method to estimate delay based on CR - historical demand (weak approximation) v = net.movement(mv).totmoved/ + (sim.time.minutes+sim.timestep)*60.0 c = net.movement(mv).totcap/ + (sim.time.minutes+sim.timestep)*60.0 IF (c.LE.0) THEN turnveh(im) = 999.9 ELSE turnveh(im) = EXP(3.8*v/c)/60.D0 ENDIF ELSE IF (net.node(nodeup).plan(ip).type.EQ.CTL_NONE) THEN turnveh(im) = 0.D0 ENDIF ELSE turnveh(im) = dsum / isc ! average waiting time for ! movement index im ENDIF ENDDO c -- convert the delay to penalty c -- all the delay is in turnveh(i,j), which j=1,2,3,4,5 CR - Loop over all links in the network nodeup = net.link(i).iunod ip = net.link(i).backindex CR - Loop over all links that *feed* this link DO j = 1,net.link(i).numuslinks il=net.movement(net.link(i).usmoveptr(j)).fromlink ili = j CR - Find the upstream link's relative backpointr index which CR - is used to assign penalties in the penalty array itmp=net.link(il).iunod mtmp=net.spd.bsd.backpointr(nodeup+1) - + net.spd.bsd.backpointr(nodeup) do kk=1,mtmp itmp2=net.spd.bsd.backstr1((kk-1) + + net.spd.bsd.backpointr(nodeup)) IF (itmp.eq.itmp2) THEN CR - This backpointr points to the upstream link associated CR - with the current feeding link (the upstream nodes are CR - identical) IF (net.movement(net.link(i).usmoveptr(j)).type + .EQ.0) THEN CR - This movement is *not* allowed according to the CR - movement file. Penalize it heavilty so it won't CR - belong to a shortest path net.spd.pn.penalty(ip,kk) = INFINITY net.movement(net.link(i).usmoveptr(j)).penalty = + INFINITY ELSE CR - This movement *is* allowed by the movement file. CR - Penalize the movement according to the expected CR - delay as given by current conditions in the CR - network. net.spd.pn.penalty(ip,kk)=turnveh(ili) net.movement(net.link(i).usmoveptr(j)).penalty = + turnveh(ili) ENDIF go to 601 endif ENDDO WRITE(ostr,'(A,3I5,A)') 'ERROR : ',i,nodeup,itmp,CHAR(0) CALL DYNA_ERROR(ostr + ,DYNA_FATAL_ERROR + ,DYNA_LOGIC_BUG + ,DYNA_UNKNOWN_ERROR) 601 CONTINUE ENDDO Cr - Verify that all movements have been properly penalized and that Cr - all movements which are not allowed (don't exist) have infinite Cr - penalty ip = net.link(i).backindex itmp = 0 DO j = 1,7 IF (net.spd.pn.penalty(ip,j).LT.INFINITY-1) itmp = itmp+1 ENDDO IF (itmp.NE.net.link(i).numuslinks) THEN WRITE(ostr,'(A,I5,I5,I5,I5,I5,A)') + 'PENALTY ERROR! i,un,dn:', + i,net.link(i).iunod,net.link(i).idnod + ,itmp,net.link(i).numuslinks,CHAR(0) CALL DYNA_ERROR(ostr + ,DYNA_FATAL_ERROR + ,DYNA_LOGIC_BUG + ,DYNA_UNKNOWN_ERROR) ENDIF LINK_PENALTYCAL = sim.status RETURN END