introducing priority variables

src/RationalFunctionFields/RationalFunctionField.jl

@@ -278,6 +278,7 @@ Applies the following passes:
     rff::RationalFunctionField;
     discard_redundant = true,
     reversed_order = false,
+    priority_variables = [],
 )
     time_start = time_ns()
     fracs = dennums_to_fractions(rff.dennums)
@@ -290,7 +291,11 @@ Applies the following passes:
     end
     # Remove redundant pass
     if discard_redundant
-        sort!(fracs, lt = rational_function_cmp)
+        fracs_priority = filter(f -> issubset(vars(f), priority_variables), fracs)
+        fracs_rest = filter(f -> !(f in fracs_priority), fracs)
+        sort!(fracs_priority, lt = rational_function_cmp)
+        sort!(fracs_rest, lt = rational_function_cmp)
+        fracs = vcat(fracs_priority, fracs_rest)
         @debug "The pool of fractions:\n$(join(map(repr, fracs), ",\n"))"
         if reversed_order
             non_redundant = collect(1:length(fracs))
@@ -323,7 +328,7 @@ Applies the following passes:
         @debug "Out of $(length(fracs)) simplified generators there are $(length(non_redundant)) non redundant"
         fracs = fracs[non_redundant]
     end
-    sort!(fracs, lt = rational_function_cmp)
+    sort!(fracs, lt = (f, g) -> rational_function_cmp(f, g))
     spring_cleaning_pass!(fracs)
     _runtime_logger[:id_beautifulization] += (time_ns() - time_start) / 1e9
     return fracs
@@ -552,6 +557,7 @@ Result is correct (in the Monte-Carlo sense) with probability at least `prob_thr
     simplify = :standard,
     check_variables = false, # almost always slows down and thus turned off
     rational_interpolator = :VanDerHoevenLecerf,
+    priority_variables = [],
 )
     @info "Simplifying generating set. Simplification level: $simplify"
     _runtime_logger[:id_groebner_time] = 0.0
@@ -620,8 +626,10 @@ Result is correct (in the Monte-Carlo sense) with probability at least `prob_thr
     @debug """
 Final cleaning and simplification of generators. 
 Out of $(length(new_fracs)) fractions $(length(new_fracs_unique)) are syntactically unique."""
-    runtime =
-        @elapsed new_fracs = beautiful_generators(RationalFunctionField(new_fracs_unique))
+    runtime = @elapsed new_fracs = beautiful_generators(
+        RationalFunctionField(new_fracs_unique),
+        priority_variables = priority_variables,
+    )
     @debug "Checking inclusion with probability $prob_threshold"
     runtime =
         @elapsed result = issubfield(rff, RationalFunctionField(new_fracs), prob_threshold)

src/RationalFunctionFields/rankings.jl

@@ -37,7 +37,7 @@ function generating_set_rank(funcs)
 end
 
 """
-    rational_function_cmp(f, g; by=:naive)
+    rational_function_cmp(f, g; by=:naive, priority_variables=[])
 
 Returns `true` if `f < g`.
 
@@ -49,7 +49,6 @@ Provides keyword argument `by`, a sorting order. Possible options are:
 - `:naive`: Compare features one by one. Features in the order of importance:
     - Sum of total degrees
     - Total number of terms
-    - Higherst appearing variable
     - Total degree of the denominator
     - Leading terms of denominator and numerator
 """
@@ -61,14 +60,6 @@ function rational_function_cmp(f, g; by = :naive)
         flag = compare_rational_func_by(f, g, length, :additive)
         flag == 1 && return false
         flag == -1 && return true
-        lead_f = maximum(vars(f))
-        lead_g = maximum(vars(g))
-        if lead_f < lead_g
-            return true
-        end
-        if lead_g < lead_f
-            return false
-        end
         flag = compare_rational_func_by(f, g, total_degree, :denominator)
         flag == 1 && return false
         flag == -1 && return true

src/identifiable_functions.jl

@@ -123,6 +123,7 @@ function _find_identifiable_functions(
             seed = seed,
             simplify = simplify,
             rational_interpolator = rational_interpolator,
+            priority_variables = [parent_ring_change(p, bring) for p in ode.parameters],
         )
     else
         id_funcs_fracs = dennums_to_fractions(id_funcs)

test/identifiable_functions.jl

@@ -931,7 +931,7 @@ ode = StructuralIdentifiability.@ODEmodel(
     y2(t) = k6 * (x1 + x2 + 2x3),
     y3(t) = k7 * EpoR_A
 )
-ident_funcs = [k3, EpoR_A * k7, k1 // k7, k5 // k6, k2 // k6]
+ident_funcs = [k2 // k6, k3, EpoR_A * k7, EpoR_A * k1, k5 // k6]
 push!(test_cases, (ode = ode, ident_funcs = ident_funcs))
 
 ode = @ODEmodel(x1'(t) = x1, x2'(t) = x2, y(t) = x1 + x2(t))