Skip to content

Commit 5486e7f

Browse files
fit-alessandro-bertifit-alessandro-berti
committed
further bug fix
1 parent 73bdaa7 commit 5486e7f

File tree

1 file changed

+120
-75
lines changed

1 file changed

+120
-75
lines changed

‎pm4py/objects/petri_net/utils/murata.py‎

Lines changed: 120 additions & 75 deletions
Original file line numberDiff line numberDiff line change
@@ -11,10 +11,18 @@ def apply_reduction(
1111
net: PetriNet, im: Marking, fm: Marking
1212
) -> Tuple[PetriNet, Marking, Marking]:
1313
"""
14-
Apply the Murata reduction to an accepting Petri net, removing the structurally redundant places.
14+
Apply the Murata reduction to an accepting Petri net, removing structurally redundant (implicit) places.
1515
16-
The implementation follows the Berthelot algorithm as in:
17-
https://svn.win.tue.nl/repos/prom/Packages/Murata/Trunk/src/org/processmining/algorithms/BerthelotAlgorithm.java
16+
This implementation follows the (Berthelot) implicit-place check by searching, via ILP,
17+
for a place-flow f = a_p * p - Σ a_q * q (a_p >= 1, a_q >= 0) such that:
18+
19+
(1) f is a P-invariant (flow): a_p*C(p,t) - Σ a_q*C(q,t) = 0 for all transitions t
20+
where C(p,t) = Post(p,t) - Pre(p,t)
21+
22+
(2) a_p*Pre(p,t) - Σ a_q*Pre(q,t) <= k for all transitions t
23+
with k = a_p*M0(p) - Σ a_q*M0(q), and k >= 0
24+
25+
If such a solution exists, the candidate place is implicit and can be removed without changing behavior.
1826
1927
Parameters
2028
---------------
@@ -28,101 +36,137 @@ def apply_reduction(
2836
Returns
2937
--------------
3038
net
31-
Petri net
39+
Reduced Petri net
3240
im
33-
Initial marking
41+
Initial marking (unchanged)
3442
fm
35-
Final marking
43+
Final marking (unchanged)
3644
"""
3745
places = sorted(list(net.places), key=lambda x: x.name)
38-
place_index = {p: i for i, p in enumerate(places)}
39-
n_places = len(places)
40-
4146
redundant = set()
4247

48+
# Choose solver backend
49+
proposed_solver = solver.SCIPY
50+
if importlib.util.find_spec("pulp"):
51+
proposed_solver = solver.PULP
52+
else:
53+
if constants.SHOW_INTERNAL_WARNINGS:
54+
warnings.warn(
55+
"solution from scipy may be unstable. Please install PuLP (pip install pulp) for fully reliable results."
56+
)
57+
4358
for place in places:
59+
# Skip already proven redundant places
60+
if place in redundant:
61+
continue
62+
4463
# Skip places in the initial or final markings
4564
if place in im or place in fm:
4665
continue
4766

67+
# Work on the "current" net as if previously found redundant places were removed:
68+
active_places = [p for p in places if p not in redundant]
69+
active_places = sorted(active_places, key=lambda x: x.name)
70+
place_index = {p: i for i, p in enumerate(active_places)}
71+
72+
n_places = len(active_places)
73+
n_vars = n_places + 1 # last variable is k
74+
k_idx = n_vars - 1
75+
4876
Aeq = []
4977
Aub = []
5078
beq = []
5179
bub = []
5280

53-
# first constraint
54-
constraint = [0] * (n_places + 1)
55-
for p2 in im:
56-
if p2 not in redundant:
57-
if p2 == place:
58-
constraint[place_index[p2]] = im[p2]
59-
else:
60-
constraint[place_index[p2]] = -im[p2]
61-
constraint[-1] = -1
62-
Aeq.append(constraint)
81+
# ---------------------------------------------------------------------
82+
# (A) k = a_p*M0(p) - Σ a_q*M0(q)
83+
# => a_p*M0(p) - Σ a_q*M0(q) - k = 0
84+
# ---------------------------------------------------------------------
85+
eq = [0] * n_vars
86+
for p2 in active_places:
87+
m0 = im[p2] if p2 in im else 0
88+
if m0 == 0:
89+
continue
90+
if p2 == place:
91+
eq[place_index[p2]] += m0
92+
else:
93+
eq[place_index[p2]] -= m0
94+
eq[k_idx] = -1
95+
Aeq.append(eq)
6396
beq.append(0)
6497

65-
# second constraints
98+
# ---------------------------------------------------------------------
99+
# (B) Flow constraints: for all transitions t,
100+
# a_p*C(p,t) - Σ a_q*C(q,t) = 0
101+
# where C(p,t)=Post(p,t)-Pre(p,t)
102+
# ---------------------------------------------------------------------
66103
for trans in net.transitions:
67-
constraint = [0] * (n_places + 1)
104+
pre = {}
105+
post = {}
68106

69107
for arc in trans.in_arcs:
70-
p2 = arc.source
71-
if p2 not in redundant:
72-
if p2 == place:
73-
constraint[place_index[p2]] = arc.weight
74-
else:
75-
constraint[place_index[p2]] = -arc.weight
76-
77-
constraint[-1] = -1
78-
Aub.append(constraint)
79-
bub.append(0)
80-
81-
# third constraints (FIXED SIGNS)
82-
for trans in net.transitions:
83-
constraint = [0] * (n_places + 1)
108+
src = arc.source
109+
if src in place_index:
110+
pre[src] = pre.get(src, 0) + arc.weight
84111

85112
for arc in trans.out_arcs:
86-
p2 = arc.target
87-
if p2 not in redundant:
88-
# FIX: candidate place must have +weight, other places -weight
89-
# (the previous implementation inverted these signs)
90-
if p2 == place:
91-
constraint[place_index[p2]] = arc.weight
92-
else:
93-
constraint[place_index[p2]] = -arc.weight
94-
95-
Aub.append(constraint)
96-
bub.append(0)
97-
98-
# fourth constraint
99-
for p2 in places:
100-
if p2 not in redundant:
101-
constraint = [0] * (n_places + 1)
102-
constraint[place_index[p2]] = -1
103-
Aub.append(constraint)
113+
tgt = arc.target
114+
if tgt in place_index:
115+
post[tgt] = post.get(tgt, 0) + arc.weight
116+
117+
eq = [0] * n_vars
118+
for p2 in active_places:
119+
c_val = post.get(p2, 0) - pre.get(p2, 0)
120+
if c_val == 0:
121+
continue
104122
if p2 == place:
105-
bub.append(-1)
123+
eq[place_index[p2]] += c_val
106124
else:
107-
bub.append(0)
108-
109-
# fifth constraint
110-
constraint = [0] * (n_places + 1)
111-
constraint[-1] = -1
112-
Aub.append(constraint)
113-
bub.append(0)
114-
115-
c = [1] * (n_places + 1)
116-
integrality = [1] * (n_places + 1)
117-
118-
proposed_solver = solver.SCIPY
119-
if importlib.util.find_spec("pulp"):
120-
proposed_solver = solver.PULP
121-
else:
122-
if constants.SHOW_INTERNAL_WARNINGS:
123-
warnings.warn(
124-
"solution from scipy may be unstable. Please install PuLP (pip install pulp) for fully reliable results."
125-
)
125+
eq[place_index[p2]] -= c_val
126+
# k has coefficient 0 here
127+
Aeq.append(eq)
128+
beq.append(0)
129+
130+
# ---------------------------------------------------------------------
131+
# (C) Implicitness inequalities: for all transitions t,
132+
# a_p*Pre(p,t) - Σ a_q*Pre(q,t) - k <= 0
133+
# ---------------------------------------------------------------------
134+
for trans in net.transitions:
135+
ineq = [0] * n_vars
136+
for arc in trans.in_arcs:
137+
p2 = arc.source
138+
if p2 not in place_index:
139+
continue
140+
if p2 == place:
141+
ineq[place_index[p2]] += arc.weight
142+
else:
143+
ineq[place_index[p2]] -= arc.weight
144+
ineq[k_idx] = -1
145+
Aub.append(ineq)
146+
bub.append(0)
147+
148+
# ---------------------------------------------------------------------
149+
# (D) Variable lower bounds:
150+
# a_p >= 1, a_q >= 0, k >= 0
151+
# Implemented as -a <= b
152+
# ---------------------------------------------------------------------
153+
for p2 in active_places:
154+
ineq = [0] * n_vars
155+
ineq[place_index[p2]] = -1
156+
Aub.append(ineq)
157+
if p2 == place:
158+
bub.append(-1) # -a_p <= -1 -> a_p >= 1
159+
else:
160+
bub.append(0) # -a_q <= 0 -> a_q >= 0
161+
162+
ineq = [0] * n_vars
163+
ineq[k_idx] = -1
164+
Aub.append(ineq)
165+
bub.append(0) # -k <= 0 -> k >= 0
166+
167+
# Objective (any feasible solution is enough; minimize sum to keep it small)
168+
c = [1] * n_vars
169+
integrality = [1] * n_vars
126170

127171
xx = solver.apply(
128172
c,
@@ -139,7 +183,8 @@ def apply_reduction(
139183
):
140184
redundant.add(place)
141185

142-
for place in redundant:
143-
net = remove_place(net, place)
186+
# Remove redundant places (deterministic order)
187+
for pl in sorted(list(redundant), key=lambda x: x.name):
188+
net = remove_place(net, pl)
144189

145190
return net, im, fm

0 commit comments

Comments
 (0)