Merge pull request #1240 from SciML/permit_function_name_interpolation

Enables DSL interpolation of a function name only

Author: TorkelE

docs/src/model_creation/functional_parameters.md

@@ -14,11 +14,10 @@ spline = LinearInterpolation((2 .+ ts) ./ (1 .+ ts), ts)
 @parameters (pIn::typeof(spline))(..)
 plot(spline)
 ```
-Next, we create our model, [interpolating](@ref dsl_advanced_options_symbolics_and_DSL_interpolation) the input parameter into it (making it a function of `t`).
+Next, we create our model, [interpolating](@ref dsl_advanced_options_symbolics_and_DSL_interpolation) the input parameter into the `@reaction_network` declaration.
 ```@example functional_parameters_basic_example
-input = pIn(default_t())
 bd_model = @reaction_network begin
-    $input, 0 --> X
+    $pIn(t), 0 --> X
     d, X --> 0
 end
 ```
@@ -76,11 +75,10 @@ plot(sol)
 ```
 
 ### [Interpolating the input into the DSL](@id functional_parameters_circ_rhythm_dsl)
-It is possible to use time-dependent inputs when creating models [through the DSL](@ref dsl_description) as well. However, it can still be convenient to declare the input parameter programmatically as above. Using it, we form an expression of it as a function of time, and then [interpolate](@ref dsl_advanced_options_symbolics_and_DSL_interpolation) it into our DSL-declaration:
+It is possible to use time-dependent inputs when creating models [through the DSL](@ref dsl_description) as well. However, it can still be convenient to declare the input parameter programmatically as above. Next, we can [interpolate](@ref dsl_advanced_options_symbolics_and_DSL_interpolation) it into our DSL-declaration (ensuring to also make it a function of `t`):
 ```@example functional_parameters_circ_rhythm
-input = light_in(t)
 rs_dsl = @reaction_network rs begin
-    (kA*$input, kD), Pᵢ <--> Pₐ
+    (kA*$light_in(t), kD), Pᵢ <--> Pₐ
 end
 ```
 We can confirm that this model is identical to our programmatic one (and should we wish to, we can simulate it using identical syntax).
@@ -112,14 +110,12 @@ I_rate = LinearInterpolation(I_measured, I_grid)
 plot(I_rate; label = "Measured infection rate")
 plot!(I_grid, I_grid; label = "Normal SIR infection rate")
 ```
-Next, we create our model (using the DSL approach). As `I_rate` will be a function of $I$ we will need to declare this species first as well.
+Next, we create our model (using the DSL approach).
 ```@example functional_parameters_sir
 using Catalyst
 @parameters (inf_rate::typeof(I_rate))(..)
-@species I(default_t())
-inf_rate_in = inf_rate(I)
 sir = @reaction_network rs begin
-    k1*$inf_rate_in, S --> I
+    k1*$inf_rate(I), S --> I
     k2, I --> R
 end
 nothing # hide

src/dsl.jl

@@ -984,7 +984,8 @@ function recursive_escape_functions!(expr::ExprValues, syms_skip = [])
     (typeof(expr) != Expr) && (return expr)
     foreach(i -> expr.args[i] = recursive_escape_functions!(expr.args[i], syms_skip),
         1:length(expr.args))
-    if (expr.head == :call) && !isdefined(Catalyst, expr.args[1]) && expr.args[1] ∉ syms_skip
+    if (expr.head == :call) && (expr.args[1] isa Symbol) &&!isdefined(Catalyst, expr.args[1]) && 
+            expr.args[1] ∉ syms_skip
         expr.args[1] = esc(expr.args[1])
     end
     expr

test/reactionsystem_core/functional_parameters.jl

@@ -41,9 +41,8 @@ let
     rs_pIn = complete(rs_pIn)
 
     # Defines a `ReactionSystem` with the input parameter (DSL).
-    input = pIn(t)
     rs_pIn_dsl = @reaction_network rs_pIn begin
-        ($input,d), 0 <--> X
+        ($pIn(t),d), 0 <--> X
         (k1*X,k2), Y1 <--> Y2
     end
 
@@ -123,14 +122,13 @@ let
     spline1d = LinearInterpolation(100.0*Is ./ (100.0 .+ Is), Is)
     @parameters (i_rate::typeof(spline1d))(..)
 
-    # Decalres the models.
+    # Declares the models.
     sir = @reaction_network begin
         k1*100*I/(100 + I), S --> I
         k2, I --> R
     end
-    input = i_rate(sir.I)
     sir_funcp = @reaction_network rs begin
-        k1*$(input), S --> I
+        k1*$i_rate(I), S --> I
         k2, I --> R
     end
 
@@ -169,3 +167,29 @@ let
     @named rs = ReactionSystem(rxs, t)
     @test issetequal([Catalyst.get_unit(rx.rate) for rx in reactions(rs)], [u"mol/(s*m^3)", u"1/s"])
 end
+
+# Tests that a functional parameter can be interpolated as the function only, as a function of a 
+# symbolic variable, or interpolated with the functional parameter and argument separately.
+let 
+    # Prepares the functional parameter.
+    ts = collect(0.0:0.1:1.0)
+    spline = LinearInterpolation(log.(ts), ts)
+    t_var = default_t()
+    @parameters (pIn::typeof(spline))(..)
+    pIn_func = pIn(t_var)
+
+    # Creates models using different approaches, check that all yield the same model.
+    rn1 = @reaction_network rn begin
+        p, 0 --> X
+        $pIn_func, 0 --> X
+    end
+    rn2 = @reaction_network rn begin
+        p, 0 --> X
+        $pIn(t), 0 --> X
+    end
+    rn3 = @reaction_network rn begin
+        p, 0 --> X
+        $pIn($t_var), 0 --> X
+    end
+    @test rn1 == rn2 == rn3
+end