diff --git a/examples/building_models/01-building-models.py b/examples/building_models/01-building-models.py
index 68ee92d..2ecccbf 100644
--- a/examples/building_models/01-building-models.py
+++ b/examples/building_models/01-building-models.py
@@ -5,7 +5,7 @@ import stormpy.examples
 import stormpy.examples.files
 
 
-def example_building_models():
+def example_building_models_01():
     path = stormpy.examples.files.drn_ctmc_dft
     model = stormpy.build_model_from_drn(path)
     print(model.model_type)
@@ -24,6 +24,5 @@ def example_building_models():
     print("Number of states: {}".format(model.nr_states))
 
 
-
 if __name__ == '__main__':
-    example_building_models()
\ No newline at end of file
+    example_building_models_01()
\ No newline at end of file
diff --git a/examples/building_models/02-building-models.py b/examples/building_models/02-building-models.py
new file mode 100644
index 0000000..12f76ac
--- /dev/null
+++ b/examples/building_models/02-building-models.py
@@ -0,0 +1,60 @@
+import stormpy
+import stormpy.core
+
+import stormpy.examples
+import stormpy.examples.files
+
+import stormpy.info
+import pycarl
+
+def example_building_models_02():
+
+
+    import stormpy.pars
+    if stormpy.info.storm_ratfunc_use_cln():
+        import pycarl.cln as pc
+    else:
+        import pycarl.gmp as pc
+
+    def make_factorized_rf(var, cache):
+        num = pc.FactorizedPolynomial(pc.Polynomial(var), cache)
+        denom = pc.FactorizedPolynomial(pc.Rational(1))
+        return pc.FactorizedRationalFunction(num, denom)
+
+        # And the parametric
+    path = stormpy.examples.files.drn_pdtmc_die
+    model = stormpy.build_parametric_model_from_drn(path)
+    print(model.model_type)
+    print("Number of states: {}".format(model.nr_states))
+
+
+
+    parameters = model.collect_probability_parameters()
+    bar_parameters = dict()
+
+    for p in parameters:
+        # Ensure that variables with that name are not recognized by pycarl.
+        assert pycarl.variable_with_name(p.name + "_bar").is_no_variable
+        bar_parameters[p] = pycarl.Variable(p.name + "_bar")
+        print(pc.Polynomial(1) - p)
+    substitutions = dict([[pc.Polynomial(1) - p, bar_parameters[p]] for p in parameters])
+    print(substitutions)
+
+    matrix = model.transition_matrix
+    for e in matrix:
+        val = e.value()
+        if val.is_constant():
+            continue
+        val_pol = val.numerator.polynomial()
+        cache = val.numerator.cache()
+        for sub, repl in substitutions.items():
+            print(sub)
+            print(val_pol - sub)
+            if val_pol - sub == 0:
+                print("Found substitution")
+                e.set_value(make_factorized_rf(repl, cache))
+                break # Assume only one substitution per entry
+    print(matrix)
+
+if __name__ == '__main__':
+    example_building_models_02()
\ No newline at end of file