`remake` doesn't recompute conservation law constants

[jonathanfischer97]

Made a comment on @TorkelE 's issue here that I'm making a seperate issue about:

"Was just about to make an issue related to this, it's a particular problem when using remake to update the initial conditions of an ODEProblem(rs, ...; remove_conserved = true).

• This is because the Γ conservation law constants are defined as parameters, which have symbolic dependencies on the initial conditions.

• Perhaps an overload is necessary for the internal SciMLBase.remake function here:

  function _updated_u0_p_symmap(prob, u0, ::Val{true}, p, ::Val{true})
  isu0dep = any(symbolic_type(v) !== NotSymbolic() for (_, v) in u0)
  ispdep = any(symbolic_type(v) !== NotSymbolic() for (_, v) in p)
  
  if !isu0dep && !ispdep
      return remake_buffer(prob, state_values(prob), u0),
      remake_buffer(prob, parameter_values(prob), p)
  end
  if !isu0dep
      u0 = remake_buffer(prob, state_values(prob), u0)
      return _updated_u0_p_symmap(prob, u0, Val(false), p, Val(true))
  end
  if !ispdep
      p = remake_buffer(prob, parameter_values(prob), p)
      return _updated_u0_p_symmap(prob, u0, Val(true), p, Val(false))
  end
  
  varmap = merge(u0, p)
  u0 = anydict(k => symbolic_type(v) === NotSymbolic() ? v : symbolic_evaluate(v, varmap)
  for (k, v) in u0)
  p = anydict(k => symbolic_type(v) === NotSymbolic() ? v : symbolic_evaluate(v, varmap)
  for (k, v) in p)
  return remake_buffer(prob, state_values(prob), u0),
  remake_buffer(prob, parameter_values(prob), p)
  end

  • Here the main issue I've found are that the ispdep check doesn't find that Γ in the parameters is a symbolic type, because it is evaluated to numerical values in the ODEProblem. The only place where the symbolic expressions for Γ are exposed are in defaults field of the ODESystem.
  • This obstructs the ostensible allowance for interdependencies between u0 and p, seen with the varmap = merge(u0, p) line.
  • Currently my workaround is calling symbolic_evaluate prior to calling remake using the symbolic expressions from the defaults map:

    
    function compute_Γ(tunable_var_value_map, constraint_syms, sys_defs::T) where T

  * Merge the default system definitions with the new symbolic map

  merged_defs = merge(sys_defs, Dict(tunable_var_value_map))

  NTuple{4, eltype(tunable_var_value_map)}(constraint => symbolic_evaluate(merged_defs[constraint], merged_defs) for constraint in constraint_syms) end

I can make a separate issue for this in this repo (not sure if it pertains most here, or ModelingToolkit.jl or SciMLBase.jl) if this is truly a problem and I'm not missing anything."

Originally posted by @jonathanfischer97 in https://github.com/SciML/Catalyst.jl/issues/806#issuecomment-2059826318

[TorkelE]

I can confirm that this is an issue coming from MTK, I have raised it there as well and will push for them sorting it out asap.