<> <> DIRECTORY PriorityQueue USING[Ref, Predict, Insert, Empty, Remove], Real USING[FixI, RoundI], RealVec USING[Handle, Object]; RealVecImpl: MONITOR -- protects sortee IMPORTS PriorityQueue, Real EXPORTS RealVec = { OPEN RealVec; <> LengthFault: PUBLIC ERROR[vec: Handle, expectedLength: NAT] = CODE; <> <> <> Smooth: PUBLIC PROC[h: Handle, nElements: NAT] RETURNS[ans: Handle] = { IF nElements = h.length THEN RETURN[h]; <> ans _ NEW[Object[nElements]]; FOR j: NAT IN [0..nElements) DO x: REAL = j; startT: REAL _ x / nElements; endT: REAL _ (x + 1.0) / nElements; a: NAT _ Real.FixI[startT * h.length]; b: NAT _ Real.RoundI[endT * h.length]; -- may be too small by 1 dtOld: REAL _ 1.0/h.length; dtNew: REAL _ endT - startT; IF b*dtOld < endT THEN b _ b + 1; IF a = b THEN ERROR; ans.elements[j] _ 0.0; <> <> FOR z: NAT IN [a..b) -- left edge index (from 0) identifies the bucket <> DO --add in the value in the zth bucket times the fraction of the time <> ans.elements[j] _ ans.elements[j] + h.elements[z] -- the value in the zth bucket * (MIN[dtOld, MIN[endT, (z + 1)*dtOld] - MAX[startT, z*dtOld]])/dtNew ENDLOOP; ENDLOOP}; <> All: PUBLIC PROC[length: NAT, value: REAL _ 0.0] RETURNS[ans: Handle] = {ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ value ENDLOOP}; <> IndexVec: PUBLIC PROC[length: NAT] RETURNS[ans: Handle] = {ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ i ENDLOOP}; <> SubVec: PUBLIC PROC[h: Handle, firstIndex, lastIndex: NAT] RETURNS[ans: Handle] = {length: NAT _ (IF h = NIL THEN 0 ELSE h.length); IF lastIndex < firstIndex OR h = NIL OR lastIndex >= h.length THEN ERROR LengthFault[h, length]; length _ lastIndex - firstIndex + 1; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ h.elements[firstIndex + i] ENDLOOP}; ConCat: PUBLIC PROC[h1, h2: Handle] RETURNS[ans: Handle] = { length: NAT; IF h1 = NIL OR h2 = NIL THEN ERROR LengthFault[NIL, 0]; IF LOOPHOLE[LOOPHOLE[h1.length, CARDINAL] + LOOPHOLE[h2.length, CARDINAL], INTEGER] < 0 THEN ERROR LengthFault[NIL, 0]; -- overflow check length _ h1.length + h2.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..h1.length) DO ans.elements[i] _ h1.elements[i] ENDLOOP; FOR i: NAT IN [0..h2.length) DO ans.elements[h1.length + i] _ h2.elements[i] ENDLOOP}; <> <> Permute: PUBLIC PROC[h: Handle, perms: Handle] RETURNS[ans: Handle] = { length: NAT _ perms.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO x: NAT = Real.FixI[perms.elements[i]]; IF x >= h.length THEN ERROR LengthFault[h, x+1]; ans.elements[i] _ h.elements[x]; ENDLOOP; }; sortee: Handle; sortPred: SAFE PROC[x, y: REF ANY, data: REF ANY] RETURNS[BOOLEAN] = TRUSTED {i: NAT = LOOPHOLE--NARROW--[x, REF NAT]^; j: NAT = LOOPHOLE--NARROW--[y, REF NAT]^; RETURN[sortee.elements[i] < sortee.elements[j]]}; <> <> SortOrder: PUBLIC ENTRY PROC[h: Handle] RETURNS[perms: Handle] = {ENABLE UNWIND => NULL; length: NAT _ h.length; q: PriorityQueue.Ref; sortee _ h; q _ PriorityQueue.Predict[pred: sortPred, size: h.length]; perms _ NEW[Object[length]]; FOR i: NAT IN [0..h.length) DO PriorityQueue.Insert[q, NEW[NAT _ i]] ENDLOOP; FOR i: NAT IN [0..h.length) DO perms.elements[i] _ NARROW[PriorityQueue.Remove[q], REF NAT]^ ENDLOOP; IF NOT PriorityQueue.Empty[q] THEN ERROR; }; <> Add: PUBLIC PROC[h1, h2: Handle] RETURNS[ans: Handle] = { length: NAT; IF h1 = NIL OR h2 = NIL THEN ERROR LengthFault[NIL, 0]; IF h1.length # h2.length THEN ERROR LengthFault[h2, h1.length]; length _ h1.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ h1.elements[i] + h2.elements[i] ENDLOOP}; Subtract: PUBLIC PROC[h1, h2: Handle] RETURNS[ans: Handle] = { length: NAT; IF h1 = NIL OR h2 = NIL THEN ERROR LengthFault[NIL, 0]; IF h1.length # h2.length THEN ERROR LengthFault[h2, h1.length]; length _ h1.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ h1.elements[i] - h2.elements[i] ENDLOOP}; Multiply: PUBLIC PROC[h1, h2: Handle] RETURNS[ans: Handle] = { length: NAT; IF h1 = NIL OR h2 = NIL THEN ERROR LengthFault[NIL, 0]; IF h1.length # h2.length THEN ERROR LengthFault[h2, h1.length]; length _ h1.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ h1.elements[i] * h2.elements[i] ENDLOOP}; Divide: PUBLIC PROC[h1, h2: Handle] RETURNS[ans: Handle] = { length: NAT; IF h1 = NIL OR h2 = NIL THEN ERROR LengthFault[NIL, 0]; IF h1.length # h2.length THEN ERROR LengthFault[h2, h1.length]; length _ h1.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ h1.elements[i] / h2.elements[i] ENDLOOP}; DotProduct: PUBLIC PROC[h1, h2: Handle] RETURNS[ans: REAL] = { length: NAT; IF h1 = NIL OR h2 = NIL THEN ERROR LengthFault[NIL, 0]; IF h1.length # h2.length THEN ERROR LengthFault[h2, h1.length]; length _ h1.length; ans _ 0.0; FOR i: NAT IN [0..length) DO ans _ ans + h1.elements[i] * h2.elements[i] ENDLOOP}; <> ScalarAdd: PUBLIC PROC[h: Handle, s: REAL] RETURNS[ans: Handle] = { length: NAT; IF h = NIL THEN RETURN[NIL]; length _ h.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ h.elements[i] + s ENDLOOP}; ScalarSubtract: PUBLIC PROC[h: Handle, s: REAL] RETURNS[ans: Handle] = { length: NAT; IF h = NIL THEN RETURN[NIL]; length _ h.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ h.elements[i] - s ENDLOOP}; ScalarMultiply: PUBLIC PROC[h: Handle, s: REAL] RETURNS[ans: Handle] = { length: NAT; IF h = NIL THEN RETURN[NIL]; length _ h.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ h.elements[i] * s ENDLOOP}; ScalarDivide: PUBLIC PROC[h: Handle, s: REAL] RETURNS[ans: Handle] = { length: NAT; IF h = NIL THEN RETURN[NIL]; length _ h.length; ans _ NEW[Object[length]]; FOR i: NAT IN [0..length) DO ans.elements[i] _ h.elements[i] / s ENDLOOP}; Equal: PUBLIC PROC[h1, h2: Handle] RETURNS[BOOLEAN] = { IF h1 = NIL OR h2 = NIL THEN RETURN[h1 = h2]; FOR i: NAT IN [0..h1.length) DO IF NOT almost[h1.elements[i], h2.elements[i]] THEN RETURN[FALSE] ENDLOOP; RETURN[TRUE]}; almost: PROC[p, q: REAL] RETURNS[BOOLEAN] = INLINE {RETURN[MAX[p, q] = 0.0 OR ABS[p - q]/MAX[p, q] < 0.00001]}; }.