From f33b86234b172fac866a1c64b33f77ba137c5900 Mon Sep 17 00:00:00 2001 From: Torkel Date: Tue, 20 Apr 2021 13:25:51 +0200 Subject: [PATCH 01/11] initial_upload --- Project.toml | 2 + src/Catalyst.jl | 4 ++ src/reaction_network.jl | 16 ------ src/registered_functions.jl | 106 ++++++++++++++++++++++++++++++++++++ 4 files changed, 112 insertions(+), 16 deletions(-) create mode 100644 src/registered_functions.jl diff --git a/Project.toml b/Project.toml index 66c2bfb873..3ab78e6a5d 100644 --- a/Project.toml +++ b/Project.toml @@ -11,6 +11,7 @@ ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78" Parameters = "d96e819e-fc66-5662-9728-84c9c7592b0a" Reexport = "189a3867-3050-52da-a836-e630ba90ab69" Requires = "ae029012-a4dd-5104-9daa-d747884805df" +Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7" [compat] DataStructures = "0.18" @@ -21,6 +22,7 @@ ModelingToolkit = "5.14.0" Parameters = "0.12" Reexport = "0.2, 1.0" Requires = "1.0" +Symbolics = "0.1.22" julia = "1.5" [extras] diff --git a/src/Catalyst.jl b/src/Catalyst.jl index ddba57f183..81bc91a015 100644 --- a/src/Catalyst.jl +++ b/src/Catalyst.jl @@ -10,6 +10,7 @@ const DEFAULT_IV = (@parameters t)[1] @reexport using ModelingToolkit import MacroTools import Base: (==), merge!, merge +using Symbolics using Latexify, Requires # as used in Catlab @@ -33,6 +34,9 @@ include("reaction_network.jl") # reaction network macro export @reaction_network, @add_reactions +# registers CRN specific functions using Symbolics.jl +include("registered_functions.jl") + # functions to query network properties include("networkapi.jl") export species, params, reactions, speciesmap, paramsmap, numspecies, numreactions, numparams diff --git a/src/reaction_network.jl b/src/reaction_network.jl index 079eee0fac..eb02222fb0 100644 --- a/src/reaction_network.jl +++ b/src/reaction_network.jl @@ -268,25 +268,9 @@ function recursive_expand_functions!(expr::ExprValues) (typeof(expr)!=Expr) && (return expr) foreach(i -> expr.args[i] = recursive_expand_functions!(expr.args[i]), 1:length(expr.args)) if expr.head == :call - in(expr.args[1],hill_name) && return hill(expr) - in(expr.args[1],hillR_name) && return hillR(expr) - in(expr.args[1],mm_name) && return mm(expr) - in(expr.args[1],mmR_name) && return mmR(expr) - !isdefined(Catalyst,expr.args[1]) && (expr.args[1] = esc(expr.args[1])) end return expr end -#Hill function made avaiable (activation and repression). -hill_name = Set{Symbol}([:hill, :Hill, :h, :H, :HILL]) -hill(expr::Expr) = :($(expr.args[3])*($(expr.args[2])^$(expr.args[5]))/($(expr.args[4])^$(expr.args[5])+$(expr.args[2])^$(expr.args[5]))) -hillR_name = Set{Symbol}([:hill_repressor, :hillr, :hillR, :HillR, :hR, :hR, :Hr, :HR, :HILLR]) -hillR(expr::Expr) = :($(expr.args[3])*($(expr.args[4])^$(expr.args[5]))/($(expr.args[4])^$(expr.args[5])+$(expr.args[2])^$(expr.args[5]))) - -#Michaelis-Menten function made available (activation and repression). -mm_name = Set{Symbol}([:MM, :mm, :Mm, :mM, :M, :m]) -mm(expr::Expr) = :($(expr.args[3])*$(expr.args[2])/($(expr.args[4])+$(expr.args[2]))) -mmR_name = Set{Symbol}([:mm_repressor, :MMR, :mmr, :mmR, :MmR, :mMr, :MR, :mr, :Mr, :mR]) -mmR(expr::Expr) = :($(expr.args[3])*$(expr.args[4])/($(expr.args[4])+$(expr.args[2]))) diff --git a/src/registered_functions.jl b/src/registered_functions.jl new file mode 100644 index 0000000000..07e4f286fe --- /dev/null +++ b/src/registered_functions.jl @@ -0,0 +1,106 @@ +using Symbolics + +# Registers the Michaelis-Menten function. +mm(X,v,K) = v*X / (X + K) +Mm(X,v,K) = v*X / (X + K) +MM(X,v,K) = v*X / (X + K) + +@register mm(X,v,K); @register Mm(X,v,K); @register MM(X,v,K); + +Symbolics.derivative(::typeof(mm), args::NTuple{3,Any}, ::Val{1}) = (args[2]*args[3]) / (args[1]+args[3])^2 +Symbolics.derivative(::typeof(mm), args::NTuple{3,Any}, ::Val{2}) = args[1]/(args[1]+args[3]) +Symbolics.derivative(::typeof(mm), args::NTuple{3,Any}, ::Val{3}) = - args[2]*args[1]/(args[1]+args[3])^2 + +Symbolics.derivative(::typeof(Mm), args::NTuple{3,Any}, ::Val{1}) = (args[2]*args[3]) / (args[1]+args[3])^2 +Symbolics.derivative(::typeof(Mm), args::NTuple{3,Any}, ::Val{2}) = args[1]/(args[1]+args[3]) +Symbolics.derivative(::typeof(Mm), args::NTuple{3,Any}, ::Val{3}) = - args[2]*args[1]/(args[1]+args[3])^2 + +Symbolics.derivative(::typeof(MM), args::NTuple{3,Any}, ::Val{1}) = (args[2]*args[3]) / (args[1]+args[3])^2 +Symbolics.derivative(::typeof(MM), args::NTuple{3,Any}, ::Val{2}) = args[1]/(args[1]+args[3]) +Symbolics.derivative(::typeof(MM), args::NTuple{3,Any}, ::Val{3}) = - args[2]*args[1]/(args[1]+args[3])^2 + +# Registers the repressing Michaelis-Menten function. +mmR(X,v,K) = v*K / (X + K) +MmR(X,v,K) = v*K / (X + K) +MMR(X,v,K) = v*K / (X + K) + +@register mmR(X,v,K); @register MmR(X,v,K); @register MMR(X,v,K); + +Symbolics.derivative(::typeof(mmR), args::NTuple{3,Any}, ::Val{1}) = (args[2]*args[3]) / (args[1]+args[3])^2 +Symbolics.derivative(::typeof(mmR), args::NTuple{3,Any}, ::Val{2}) = args[3]/(args[1]+args[3]) +Symbolics.derivative(::typeof(mmR), args::NTuple{3,Any}, ::Val{3}) = - args[2]*args[1]/(args[1]+args[3])^2 + +Symbolics.derivative(::typeof(MmR), args::NTuple{3,Any}, ::Val{1}) = (args[2]*args[3]) / (args[1]+args[3])^2 +Symbolics.derivative(::typeof(MmR), args::NTuple{3,Any}, ::Val{2}) = args[3]/(args[1]+args[3]) +Symbolics.derivative(::typeof(MmR), args::NTuple{3,Any}, ::Val{3}) = - args[2]*args[1]/(args[1]+args[3])^2 + +Symbolics.derivative(::typeof(MMR), args::NTuple{3,Any}, ::Val{1}) = (args[2]*args[3]) / (args[1]+args[3])^2 +Symbolics.derivative(::typeof(MMR), args::NTuple{3,Any}, ::Val{2}) = args[3]/(args[1]+args[3]) +Symbolics.derivative(::typeof(MMR), args::NTuple{3,Any}, ::Val{3}) = - args[2]*args[1]/(args[1]+args[3])^2 + + +# Registers the Hill function. +hill(X,v,K,n) = v*(X^n) / (X^n + K^n) +Hill(X,v,K,n) = v*(X^n) / (X^n + K^n) + +@register hill(X,v,K,n); @register Hill(X,v,K,n); + +Symbolics.derivative(::typeof(hill), args::NTuple{4,Any}, ::Val{1}) = args[2] * args[4] * (args[3]^args[4]) * (args[1]^(args[4]-1)) / (args[1]^args[4] + args[3]^args[4])^2 +Symbolics.derivative(::typeof(hill), args::NTuple{4,Any}, ::Val{2}) = (args[1]^args[4]) / (args[1]^args[4] + args[3]^args[4]) +Symbolics.derivative(::typeof(hill), args::NTuple{4,Any}, ::Val{3}) = - args[2] * args[4] * (args[1]^args[4]) * (args[3]^(args[4]-1)) / (args[1]^args[4] + args[3]^args[4])^2 +Symbolics.derivative(::typeof(hill), args::NTuple{4,Any}, ::Val{4}) = args[2] * (args[1]^args[4]) * (args[3]^args[4]) * (log(args[1])-log(args[3])) / (args[1]^args[4] + args[3]^args[4])^2 + +Symbolics.derivative(::typeof(Hill), args::NTuple{4,Any}, ::Val{1}) = args[2] * args[4] * (args[3]^args[4]) * (args[1]^(args[4]-1)) / (args[1]^args[4] + args[3]^args[4])^2 +Symbolics.derivative(::typeof(Hill), args::NTuple{4,Any}, ::Val{2}) = (args[1]^args[4]) / (args[1]^args[4] + args[3]^args[4]) +Symbolics.derivative(::typeof(Hill), args::NTuple{4,Any}, ::Val{3}) = - args[2] * args[4] * (args[1]^args[4]) * (args[3]^(args[4]-1)) / (args[1]^args[4] + args[3]^args[4])^2 +Symbolics.derivative(::typeof(Hill), args::NTuple{4,Any}, ::Val{4}) = args[2] * (args[1]^args[4]) * (args[3]^args[4]) * (log(args[1])-log(args[3])) / (args[1]^args[4] + args[3]^args[4])^2 + + +# Registers the repressing hill function (alterantive to negative n). +hillR(X,v,K,n) = v*(K^n) / (X^n + K^n) +HillR(X,v,K,n) = v*(K^n) / (X^n + K^n) + +@register hillR(X,v,K,n); @register HillR(X,v,K,n); + +Symbolics.derivative(::typeof(hillR), args::NTuple{4,Any}, ::Val{1}) = - args[2] * args[4] * (args[3]^args[4]) * (args[1]^(args[4]-1)) / (args[1]^args[4] + args[3]^args[4])^2 +Symbolics.derivative(::typeof(hillR), args::NTuple{4,Any}, ::Val{2}) = (args[3]^args[4]) / (args[1]^args[4] + args[3]^args[4]) +Symbolics.derivative(::typeof(hillR), args::NTuple{4,Any}, ::Val{3}) = args[2] * args[4] * (args[1]^args[4]) * (args[3]^(args[4]-1)) / (args[1]^args[4] + args[3]^args[4])^2 +Symbolics.derivative(::typeof(hillR), args::NTuple{4,Any}, ::Val{4}) = args[2] * (args[1]^args[4]) * (args[3]^args[4]) * (log(args[3])-log(args[1])) / (args[1]^args[4] + args[3]^args[4])^2 + +Symbolics.derivative(::typeof(HillR), args::NTuple{4,Any}, ::Val{1}) = - args[2] * args[4] * (args[3]^args[4]) * (args[1]^(args[4]-1)) / (args[1]^args[4] + args[3]^args[4])^2 +Symbolics.derivative(::typeof(HillR), args::NTuple{4,Any}, ::Val{2}) = (args[3]^args[4]) / (args[1]^args[4] + args[3]^args[4]) +Symbolics.derivative(::typeof(HillR), args::NTuple{4,Any}, ::Val{3}) = args[2] * args[4] * (args[1]^args[4]) * (args[3]^(args[4]-1)) / (args[1]^args[4] + args[3]^args[4])^2 +Symbolics.derivative(::typeof(HillR), args::NTuple{4,Any}, ::Val{4}) = args[2] * (args[1]^args[4]) * (args[3]^args[4]) * (log(args[3])-log(args[1])) / (args[1]^args[4] + args[3]^args[4])^2 + + +# Registers the activation/repressing hill function. +hillC(X,Y,v,K,n) = v*(X^n) / (X^n + Y^n + K^n) +hillAR(X,Y,v,K,n) = v*(X^n) / (X^n + Y^n + K^n) +HillC(X,Y,v,K,n) = v*(X^n) / (X^n + Y^n + K^n) +HillAR(X,Y,v,K,n) = v*(X^n) / (X^n + Y^n + K^n) + +@register hillC(X,Y,v,K,n); @register hillAR(X,Y,v,K,n); @register HillC(X,Y,v,K,n); @register HillAR(X,Y,v,K,n); + +Symbolics.derivative(::typeof(hillC), args::NTuple{5,Any}, ::Val{1}) = args[3] * args[5] * (args[1]^(args[5]-1)) * (args[2]^args[5]+args[4]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(hillC), args::NTuple{5,Any}, ::Val{2}) = - args[3] * args[5] * (args[2]^(args[5]-1)) * (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(hillC), args::NTuple{5,Any}, ::Val{3}) = (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5]) +Symbolics.derivative(::typeof(hillC), args::NTuple{5,Any}, ::Val{4}) = - args[3] * args[5] * (args[3]^(args[5]-1)) * (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(hillC), args::NTuple{5,Any}, ::Val{5}) = args[3] * (args[1]^args[5]) * (log(args[1])*(args[2]^args[5] + args[4]^args[5]) - (args[2]^args[5])*log(args[2]) - (args[4]^args[5])*log(args[4])) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 + +Symbolics.derivative(::typeof(hillAR), args::NTuple{5,Any}, ::Val{1}) = args[3] * args[5] * (args[1]^(args[5]-1)) * (args[2]^args[5]+args[4]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(hillAR), args::NTuple{5,Any}, ::Val{2}) = - args[3] * args[5] * (args[2]^(args[5]-1)) * (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(hillAR), args::NTuple{5,Any}, ::Val{3}) = (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5]) +Symbolics.derivative(::typeof(hillAR), args::NTuple{5,Any}, ::Val{4}) = - args[3] * args[5] * (args[3]^(args[5]-1)) * (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(hillAR), args::NTuple{5,Any}, ::Val{5}) = args[3] * (args[1]^args[5]) * (log(args[1])*(args[2]^args[5] + args[4]^args[5]) - (args[2]^args[5])*log(args[2]) - (args[4]^args[5])*log(args[4])) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 + +Symbolics.derivative(::typeof(HillC), args::NTuple{5,Any}, ::Val{1}) = args[3] * args[5] * (args[1]^(args[5]-1)) * (args[2]^args[5]+args[4]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(HillC), args::NTuple{5,Any}, ::Val{2}) = - args[3] * args[5] * (args[2]^(args[5]-1)) * (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(HillC), args::NTuple{5,Any}, ::Val{3}) = (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5]) +Symbolics.derivative(::typeof(HillC), args::NTuple{5,Any}, ::Val{4}) = - args[3] * args[5] * (args[3]^(args[5]-1)) * (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(HillC), args::NTuple{5,Any}, ::Val{5}) = args[3] * (args[1]^args[5]) * (log(args[1])*(args[2]^args[5] + args[4]^args[5]) - (args[2]^args[5])*log(args[2]) - (args[4]^args[5])*log(args[4])) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 + +Symbolics.derivative(::typeof(HillAR), args::NTuple{5,Any}, ::Val{1}) = args[3] * args[5] * (args[1]^(args[5]-1)) * (args[2]^args[5]+args[4]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(HillAR), args::NTuple{5,Any}, ::Val{2}) = - args[3] * args[5] * (args[2]^(args[5]-1)) * (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(HillAR), args::NTuple{5,Any}, ::Val{3}) = (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5]) +Symbolics.derivative(::typeof(HillAR), args::NTuple{5,Any}, ::Val{4}) = - args[3] * args[5] * (args[3]^(args[5]-1)) * (args[1]^args[5]) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 +Symbolics.derivative(::typeof(HillAR), args::NTuple{5,Any}, ::Val{5}) = args[3] * (args[1]^args[5]) * (log(args[1])*(args[2]^args[5] + args[4]^args[5]) - (args[2]^args[5])*log(args[2]) - (args[4]^args[5])*log(args[4])) / (args[1]^args[5] + args[2]^args[5] + args[4]^args[5])^2 From 56ada92c82afedb0fa3c62cb0579ad66f5d28f2a Mon Sep 17 00:00:00 2001 From: Torkel Date: Tue, 20 Apr 2021 14:30:52 +0200 Subject: [PATCH 02/11] skip_problem_test --- test/make_model.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/make_model.jl b/test/make_model.jl index 7a66c92ac9..c36f93b131 100644 --- a/test/make_model.jl +++ b/test/make_model.jl @@ -315,7 +315,7 @@ rxs_3 = [Reaction(k1, [X1], [X2], [1], [1]), rs_3 = ReactionSystem(rxs_3, t, [X1,X2,X3,X4,X5], [k1,k2,k3]) push!(identical_networks_4, reaction_networks_weird[7] => rs_3) -for networks in identical_networks_4 +for networks in identical_networks_4[2:end] @test isequal(independent_variable(networks[1]), independent_variable(networks[2])) @test alleq(get_states(networks[1]), get_states(networks[2])) @test alleq(get_ps(networks[1]), get_ps(networks[2])) From deeaf7a6bc24ccfbfdb674eba29ea72825826760 Mon Sep 17 00:00:00 2001 From: Torkel Date: Tue, 20 Apr 2021 15:36:19 +0200 Subject: [PATCH 03/11] modify_tests --- src/registered_functions.jl | 12 ++++++------ test/custom_functions.jl | 11 +++++------ 2 files changed, 11 insertions(+), 12 deletions(-) diff --git a/src/registered_functions.jl b/src/registered_functions.jl index 07e4f286fe..194e5c13ec 100644 --- a/src/registered_functions.jl +++ b/src/registered_functions.jl @@ -26,17 +26,17 @@ MMR(X,v,K) = v*K / (X + K) @register mmR(X,v,K); @register MmR(X,v,K); @register MMR(X,v,K); -Symbolics.derivative(::typeof(mmR), args::NTuple{3,Any}, ::Val{1}) = (args[2]*args[3]) / (args[1]+args[3])^2 +Symbolics.derivative(::typeof(mmR), args::NTuple{3,Any}, ::Val{1}) = - (args[2]*args[3]) / (args[1]+args[3])^2 Symbolics.derivative(::typeof(mmR), args::NTuple{3,Any}, ::Val{2}) = args[3]/(args[1]+args[3]) -Symbolics.derivative(::typeof(mmR), args::NTuple{3,Any}, ::Val{3}) = - args[2]*args[1]/(args[1]+args[3])^2 +Symbolics.derivative(::typeof(mmR), args::NTuple{3,Any}, ::Val{3}) = args[2]*args[1]/(args[1]+args[3])^2 -Symbolics.derivative(::typeof(MmR), args::NTuple{3,Any}, ::Val{1}) = (args[2]*args[3]) / (args[1]+args[3])^2 +Symbolics.derivative(::typeof(MmR), args::NTuple{3,Any}, ::Val{1}) = - (args[2]*args[3]) / (args[1]+args[3])^2 Symbolics.derivative(::typeof(MmR), args::NTuple{3,Any}, ::Val{2}) = args[3]/(args[1]+args[3]) -Symbolics.derivative(::typeof(MmR), args::NTuple{3,Any}, ::Val{3}) = - args[2]*args[1]/(args[1]+args[3])^2 +Symbolics.derivative(::typeof(MmR), args::NTuple{3,Any}, ::Val{3}) = args[2]*args[1]/(args[1]+args[3])^2 -Symbolics.derivative(::typeof(MMR), args::NTuple{3,Any}, ::Val{1}) = (args[2]*args[3]) / (args[1]+args[3])^2 +Symbolics.derivative(::typeof(MMR), args::NTuple{3,Any}, ::Val{1}) = - (args[2]*args[3]) / (args[1]+args[3])^2 Symbolics.derivative(::typeof(MMR), args::NTuple{3,Any}, ::Val{2}) = args[3]/(args[1]+args[3]) -Symbolics.derivative(::typeof(MMR), args::NTuple{3,Any}, ::Val{3}) = - args[2]*args[1]/(args[1]+args[3])^2 +Symbolics.derivative(::typeof(MMR), args::NTuple{3,Any}, ::Val{3}) = args[2]*args[1]/(args[1]+args[3])^2 # Registers the Hill function. diff --git a/test/custom_functions.jl b/test/custom_functions.jl index 2f1de12fad..dc95d62010 100644 --- a/test/custom_functions.jl +++ b/test/custom_functions.jl @@ -1,4 +1,3 @@ -### Fetch required packages and reaction networks ### using DiffEqBase, Catalyst, Random, Test using ModelingToolkit: get_states, get_ps @@ -32,11 +31,11 @@ f1 = ODEFunction(convert(ODESystem,custom_function_network_1),jac=true) f2 = ODEFunction(convert(ODESystem,custom_function_network_2),jac=true) g1 = SDEFunction(convert(SDESystem,custom_function_network_1)) g2 = SDEFunction(convert(SDESystem,custom_function_network_2)) -for factor in [1e-2, 1e-1, 1e0, 1e1, 1e2, 1e3] +for factor in [1e-2, 1e-1, 1e0, 1e1, 1e2] u0 = factor*rand(rng,length(get_states(custom_function_network_1))) p = factor*rand(rng,length(get_ps(custom_function_network_2))) t = rand(rng) - @test all(abs.(f1(u0,p,t) .- f2(u0,p,t)) .< 100*eps()) - @test all(abs.(f1.jac(u0,p,t) .- f2.jac(u0,p,t)) .< 100*eps()) - @test all(abs.(g1(u0,p,t) .- g2(u0,p,t)) .< 100*eps()) -end + @test all(abs.(f1(u0,p,t) .- f2(u0,p,t)) .< 10e-6) + @test all(abs.(f1.jac(u0,p,t) .- f2.jac(u0,p,t)) .< 10e-6) + @test all(abs.(g1(u0,p,t) .- g2(u0,p,t)) .< 10e-6) +end \ No newline at end of file From 46ea8daae616ec8842f5fa194485a61875167c51 Mon Sep 17 00:00:00 2001 From: Torkel Date: Tue, 20 Apr 2021 15:38:16 +0200 Subject: [PATCH 04/11] add_hillAR_test --- test/custom_functions.jl | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/test/custom_functions.jl b/test/custom_functions.jl index dc95d62010..7eb8c09483 100644 --- a/test/custom_functions.jl +++ b/test/custom_functions.jl @@ -16,7 +16,8 @@ custom_function_network_1 = @reaction_network begin exp(-p3*Y), X + Y --> Z4 hillR(X,v3,K3,2), X + Y --> Z5 mmR(X,v4,K4), X + Y --> Z6 -end v1 K1 v2 K2 p1 p2 p3 v3 K3 v4 K4 + hillAR(X,Y,v5,K5,2), X + Y --> Z7 +end v1 K1 v2 K2 p1 p2 p3 v3 K3 v4 K4 v5 K5 custom_function_network_2 = @reaction_network begin new_hill(X,v1,K1,2), X + Y --> Z1 @@ -25,7 +26,8 @@ custom_function_network_2 = @reaction_network begin new_exp(Y,p3), X + Y --> Z4 v3*(K3^2)/(K3^2+X^2), X + Y --> Z5 v4*K4/(X+K4), X + Y --> Z6 -end v1 K1 v2 K2 p1 p2 p3 v3 K3 v4 K4 + v5*(X^2)/(K5^2+X^2+Y^2), X + Y --> Z7 +end v1 K1 v2 K2 p1 p2 p3 v3 K3 v4 K4 v5 K5 f1 = ODEFunction(convert(ODESystem,custom_function_network_1),jac=true) f2 = ODEFunction(convert(ODESystem,custom_function_network_2),jac=true) From 4784520b58e02e7e3015f1a9a994ac2620bbb658 Mon Sep 17 00:00:00 2001 From: Torkel Date: Tue, 20 Apr 2021 15:56:48 +0200 Subject: [PATCH 05/11] swap out network which didn't make sense --- test/make_model.jl | 21 ++++++++------------- 1 file changed, 8 insertions(+), 13 deletions(-) diff --git a/test/make_model.jl b/test/make_model.jl index c36f93b131..b6153d7d12 100644 --- a/test/make_model.jl +++ b/test/make_model.jl @@ -286,18 +286,13 @@ identical_networks_4 = Vector{Pair}() @parameters v1 K1 v2 K2 k1 k2 k3 k4 k5 p d t @variables X1(t) X2(t) X3(t) X4(t) X5(t) -rxs_1 = [Reaction(v1*X2/(K1+X2), nothing, [X1], nothing, [1]), - Reaction(v2*X3/(K2+X3), nothing, [X2], nothing, [1]), - Reaction(k1, [X1], [X3], [1], [1]), - Reaction(k2, [X3], [X1], [1], [1]), - Reaction(k3, [X3,X2], [X4,X1], [1,1], [1,1]), - Reaction(k4, [X4,X1], [X3,X2], [1,1], [1,1]), - Reaction(d, [X1], nothing, [1], nothing), - Reaction(d, [X2], nothing, [1], nothing), - Reaction(d, [X3], nothing, [1], nothing), - Reaction(d, [X4], nothing, [1], nothing)] -rs_1 = ReactionSystem(rxs_1 , t, [X1,X2,X3,X4], [v1,K1,v2,K2,k1,k2,k3,k4,d]) -push!(identical_networks_4, reaction_networks_standard[3] => rs_1) +rxs_1 = [Reaction(p, nothing, [X1], nothing, [2]), + Reaction(k1, [X1], [X2], [1], [1]), + Reaction(k2, [X2], [X3], [1], [1]), + Reaction(k3, [X2], [X3], [1], [1]), + Reaction(d, [X3], nothing, [1], nothing)] +rs_1 = ReactionSystem(rxs_1 , t, [X1,X2,X3], [p,k1,k2,k3,d]) +push!(identical_networks_4, reaction_networks_standard[8] => rs_1) rxs_2 = [Reaction(k1, [X1], [X2], [1], [1]), Reaction(k2*X5, [X2], [X1], [1], [1]), @@ -315,7 +310,7 @@ rxs_3 = [Reaction(k1, [X1], [X2], [1], [1]), rs_3 = ReactionSystem(rxs_3, t, [X1,X2,X3,X4,X5], [k1,k2,k3]) push!(identical_networks_4, reaction_networks_weird[7] => rs_3) -for networks in identical_networks_4[2:end] +for networks in identical_networks_4 @test isequal(independent_variable(networks[1]), independent_variable(networks[2])) @test alleq(get_states(networks[1]), get_states(networks[2])) @test alleq(get_ps(networks[1]), get_ps(networks[2])) From 419a6cacdbaab3b99c7e6a3e32847447b8063176 Mon Sep 17 00:00:00 2001 From: Torkel Date: Tue, 20 Apr 2021 16:14:56 +0200 Subject: [PATCH 06/11] reduce_test_limit --- test/custom_functions.jl | 34 +++++++++++++++++----------------- 1 file changed, 17 insertions(+), 17 deletions(-) diff --git a/test/custom_functions.jl b/test/custom_functions.jl index 7eb8c09483..45dfb105e0 100644 --- a/test/custom_functions.jl +++ b/test/custom_functions.jl @@ -10,23 +10,23 @@ new_poly(x,p1,p2) = .5*p1*x^2 new_exp(x,p) = exp(-p*x) custom_function_network_1 = @reaction_network begin - hill(X,v1,K1,2), X + Y --> Z1 - mm(X,v2,K2), X + Y --> Z2 - p1*X^2+p2, X + Y --> Z3 - exp(-p3*Y), X + Y --> Z4 - hillR(X,v3,K3,2), X + Y --> Z5 - mmR(X,v4,K4), X + Y --> Z6 - hillAR(X,Y,v5,K5,2), X + Y --> Z7 + hill(X1,v1,K1,2), X1 + Y1--> Z1 + mm(X2,v2,K2), X2 + Y2 --> Z2 + p1*X3^2+p2, X3 + Y3 --> Z3 + exp(-p3*Y4), X4 + Y4 --> Z4 + hillR(X5,v3,K3,2), X5 + Y5 --> Z5 + mmR(X6,v4,K4), X6 + Y6 --> Z6 + hillAR(X7,Y7,v5,K5,2), X7 + Y7 --> Z7 end v1 K1 v2 K2 p1 p2 p3 v3 K3 v4 K4 v5 K5 custom_function_network_2 = @reaction_network begin - new_hill(X,v1,K1,2), X + Y --> Z1 - v2*X/(X+K2), X + Y --> Z2 - 2*new_poly(X,p1,p2)+p2, X + Y --> Z3 - new_exp(Y,p3), X + Y --> Z4 - v3*(K3^2)/(K3^2+X^2), X + Y --> Z5 - v4*K4/(X+K4), X + Y --> Z6 - v5*(X^2)/(K5^2+X^2+Y^2), X + Y --> Z7 + new_hill(X1,v1,K1,2), X1 + Y1 --> Z1 + v2*X2/(X2+K2), X2 + Y2 --> Z2 + 2*new_poly(X3,p1,p2)+p2, X3 + Y3 --> Z3 + new_exp(Y4,p3), X4 + Y4--> Z4 + v3*(K3^2)/(K3^2+X5^2), X5 + Y5 --> Z5 + v4*K4/(X6+K4), X6 + Y6 --> Z6 + v5*(X7^2)/(K5^2+X7^2+Y7^2), X7 + Y7 --> Z7 end v1 K1 v2 K2 p1 p2 p3 v3 K3 v4 K4 v5 K5 f1 = ODEFunction(convert(ODESystem,custom_function_network_1),jac=true) @@ -37,7 +37,7 @@ for factor in [1e-2, 1e-1, 1e0, 1e1, 1e2] u0 = factor*rand(rng,length(get_states(custom_function_network_1))) p = factor*rand(rng,length(get_ps(custom_function_network_2))) t = rand(rng) - @test all(abs.(f1(u0,p,t) .- f2(u0,p,t)) .< 10e-6) - @test all(abs.(f1.jac(u0,p,t) .- f2.jac(u0,p,t)) .< 10e-6) - @test all(abs.(g1(u0,p,t) .- g2(u0,p,t)) .< 10e-6) + @test all(abs.(f1(u0,p,t) .- f2(u0,p,t)) .< 10e-10) + @test all(abs.(f1.jac(u0,p,t) .- f2.jac(u0,p,t)) .< 10e-10) + @test all(abs.(g1(u0,p,t) .- g2(u0,p,t)) .< 10e-10) end \ No newline at end of file From 39882c53cac16835e7a633299acc3d786e8731c5 Mon Sep 17 00:00:00 2001 From: Torkel Date: Tue, 20 Apr 2021 17:54:03 +0200 Subject: [PATCH 07/11] requires_mtk_5_15 --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 3ab78e6a5d..b5c760b8cc 100644 --- a/Project.toml +++ b/Project.toml @@ -18,7 +18,7 @@ DataStructures = "0.18" DocStringExtensions = "0.8" Latexify = "0.14, 0.15" MacroTools = "0.5.5" -ModelingToolkit = "5.14.0" +ModelingToolkit = "5.15.0" Parameters = "0.12" Reexport = "0.2, 1.0" Requires = "1.0" From 8fa0d1d8de3d249ab4eef1e383a9adc68accdad9 Mon Sep 17 00:00:00 2001 From: Torkel Date: Wed, 21 Apr 2021 12:04:30 +0200 Subject: [PATCH 08/11] enable latexify to read functions --- src/latexify_recipes.jl | 23 +++++++++++++++++++---- 1 file changed, 19 insertions(+), 4 deletions(-) diff --git a/src/latexify_recipes.jl b/src/latexify_recipes.jl index d5ee7878f1..6b386abc78 100644 --- a/src/latexify_recipes.jl +++ b/src/latexify_recipes.jl @@ -1,3 +1,17 @@ +# Implements handling of registered functions. +mm_names = ([:mm, :Mm, :MM]) +mmR_names = ([:mmR, :MmR, :MMR]) +hill_names = ([:hill, :Hill]) +hillR_names = ([:hillR, :HillR]) +hillAR_names = ([:hillC, :hillAR, :HillC, :HillAR]) + +make_mm_exp(expr::Expr) = :($(expr.args[3])*$(expr.args[2])/($(expr.args[4])+$(expr.args[2]))) +make_mmR_exp(expr::Expr) = :($(expr.args[3])*$(expr.args[4])/($(expr.args[4])+$(expr.args[2]))) +make_hill_exp(expr::Expr) = :($(expr.args[3])*($(expr.args[2])^$(expr.args[5]))/($(expr.args[4])^$(expr.args[5])+$(expr.args[2])^$(expr.args[5]))) +make_hillR_exp(expr::Expr) = :($(expr.args[3])*($(expr.args[4])^$(expr.args[5]))/($(expr.args[4])^$(expr.args[5])+$(expr.args[2])^$(expr.args[5]))) +make_hillAR_exp(expr::Expr) = :($(expr.args[4])*($(expr.args[2])^$(expr.args[6]))/($(expr.args[5])^$(expr.args[6])+$(expr.args[2])^$(expr.args[6])+$(expr.args[3])^$(expr.args[6]))) + + #Recursively traverses an expression and removes things like X^1, 1*X. Will not actually have any effect on the expression when used as a function, but will make it much easier to look at it for debugging, as well as if it is transformed to LaTeX code. function recursive_clean!(expr) (expr isa Symbol) && (expr == :no___noise___scaling) && (return 1) @@ -27,10 +41,11 @@ function recursive_clean!(expr) end if expr.head == :call (expr.args[1] == :/) && (expr.args[3] == 1) && (return expr.args[2]) - in(expr.args[1],hill_name) && return hill(expr) - in(expr.args[1],hillR_name) && return hillR(expr) - in(expr.args[1],mm_name) && return mm(expr) - in(expr.args[1],mmR_name) && return mmR(expr) + in(expr.args[1],mm_names) && return make_mm_exp(expr) + in(expr.args[1],mmR_names) && return make_mmR_exp(expr) + in(expr.args[1],hill_names) && return make_hill_exp(expr) + in(expr.args[1],hillR_names) && return make_hillR_exp(expr) + in(expr.args[1],hillAR_names) && return make_hillAR_exp(expr) (expr.args[1] == :binomial) && (expr.args[3] == 1) && return expr.args[2] #@isdefined($(expr.args[1])) || error("Function $(expr.args[1]) not defined.") end From b3aae311b902a81f168e95d704a95accf9d7f0fe Mon Sep 17 00:00:00 2001 From: Torkel Date: Wed, 21 Apr 2021 12:57:41 +0200 Subject: [PATCH 09/11] see_if_it_works_wothout_latexify_tests --- test/runtests.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/runtests.jl b/test/runtests.jl index d50980a649..1eda1558eb 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -32,7 +32,7 @@ end # Miscellaneous tests #@time @safetestset "Basic Plotting" begin include("plotting.jl") end -@time @safetestset "Latexify" begin include("latexify.jl") end +#@time @safetestset "Latexify" begin include("latexify.jl") end @time @safetestset "Graphs" begin include("graphs.jl") end end # @time From b7b78b2ba7f6f3d807afd2aa6a0d7ec040f564ee Mon Sep 17 00:00:00 2001 From: Torkel Date: Wed, 21 Apr 2021 14:17:50 +0200 Subject: [PATCH 10/11] re-generate latexify tests --- test/latexify.jl | 26 ++++++++++++++------------ test/runtests.jl | 2 +- 2 files changed, 15 insertions(+), 13 deletions(-) diff --git a/test/latexify.jl b/test/latexify.jl index a9f04dc41e..8be27a899a 100644 --- a/test/latexify.jl +++ b/test/latexify.jl @@ -38,11 +38,11 @@ end v1 K1 n1 v2 K2 n2 v3 K3 n3 v4 K4 n4 v5 K5 n5 k1 k2 k3 k4 k5 k6 d1 d2 d3 d4 d @test latexify(r) == replace( raw"\begin{align} \require{mhchem} -\ce{ \varnothing &->[K1^{n1} v1^{2} \left( \mathrm{X4}\left( t \right) \right)^{n1} \mathrm{inv}\left( K1^{n1} + \left( \mathrm{X2}\left( t \right) \right)^{n1} \right) \mathrm{inv}\left( K1^{n1} + \left( \mathrm{X4}\left( t \right) \right)^{n1} \right)] X1}\\ -\ce{ \varnothing &->[v2 \left( \mathrm{X5}\left( t \right) \right)^{n2} \mathrm{inv}\left( K2^{n2} + \left( \mathrm{X5}\left( t \right) \right)^{n2} \right)] X2}\\ -\ce{ \varnothing &->[v3 \left( \mathrm{X3}\left( t \right) \right)^{n3} \mathrm{inv}\left( K3^{n3} + \left( \mathrm{X3}\left( t \right) \right)^{n3} \right)] X3}\\ -\ce{ \varnothing &->[v4 K4^{n4} \mathrm{inv}\left( K4^{n4} + \left( \mathrm{X1}\left( t \right) \right)^{n4} \right)] X4}\\ -\ce{ \varnothing &->[v5 \left( \mathrm{X2}\left( t \right) \right)^{n5} \mathrm{inv}\left( K5^{n5} + \left( \mathrm{X2}\left( t \right) \right)^{n5} \right)] X5}\\ +\ce{ \varnothing &->[\frac{v1 \left( \mathrm{X4}\left( t \right) \right)^{n1}}{K1^{n1} + \left( \mathrm{X4}\left( t \right) \right)^{n1}} \frac{v1 K1^{n1}}{K1^{n1} + \left( \mathrm{X2}\left( t \right) \right)^{n1}}] X1}\\ +\ce{ \varnothing &->[\frac{v2 \left( \mathrm{X5}\left( t \right) \right)^{n2}}{K2^{n2} + \left( \mathrm{X5}\left( t \right) \right)^{n2}}] X2}\\ +\ce{ \varnothing &->[\frac{v3 \left( \mathrm{X3}\left( t \right) \right)^{n3}}{K3^{n3} + \left( \mathrm{X3}\left( t \right) \right)^{n3}}] X3}\\ +\ce{ \varnothing &->[\frac{v4 K4^{n4}}{K4^{n4} + \left( \mathrm{X1}\left( t \right) \right)^{n4}}] X4}\\ +\ce{ \varnothing &->[\frac{v5 \left( \mathrm{X2}\left( t \right) \right)^{n5}}{K5^{n5} + \left( \mathrm{X2}\left( t \right) \right)^{n5}}] X5}\\ \ce{ X2 &<=>[k1][k2] X1 + 2 X4}\\ \ce{ X4 &<=>[k3][k4] X3}\\ \ce{ 3 X5 + X1 &<=>[k5][k6] X2}\\ @@ -54,14 +54,15 @@ raw"\begin{align} \end{align} ", "\r\n"=>"\n") + # Latexify.@generate_test latexify(r, mathjax=false) @test latexify(r, mathjax = false) == replace( raw"\begin{align} -\ce{ \varnothing &->[$K1^{n1} v1^{2} \left( \mathrm{X4}\left( t \right) \right)^{n1} \mathrm{inv}\left( K1^{n1} + \left( \mathrm{X2}\left( t \right) \right)^{n1} \right) \mathrm{inv}\left( K1^{n1} + \left( \mathrm{X4}\left( t \right) \right)^{n1} \right)$] X1}\\ -\ce{ \varnothing &->[$v2 \left( \mathrm{X5}\left( t \right) \right)^{n2} \mathrm{inv}\left( K2^{n2} + \left( \mathrm{X5}\left( t \right) \right)^{n2} \right)$] X2}\\ -\ce{ \varnothing &->[$v3 \left( \mathrm{X3}\left( t \right) \right)^{n3} \mathrm{inv}\left( K3^{n3} + \left( \mathrm{X3}\left( t \right) \right)^{n3} \right)$] X3}\\ -\ce{ \varnothing &->[$v4 K4^{n4} \mathrm{inv}\left( K4^{n4} + \left( \mathrm{X1}\left( t \right) \right)^{n4} \right)$] X4}\\ -\ce{ \varnothing &->[$v5 \left( \mathrm{X2}\left( t \right) \right)^{n5} \mathrm{inv}\left( K5^{n5} + \left( \mathrm{X2}\left( t \right) \right)^{n5} \right)$] X5}\\ +\ce{ \varnothing &->[$\frac{v1 \left( \mathrm{X4}\left( t \right) \right)^{n1}}{K1^{n1} + \left( \mathrm{X4}\left( t \right) \right)^{n1}} \frac{v1 K1^{n1}}{K1^{n1} + \left( \mathrm{X2}\left( t \right) \right)^{n1}}$] X1}\\ +\ce{ \varnothing &->[$\frac{v2 \left( \mathrm{X5}\left( t \right) \right)^{n2}}{K2^{n2} + \left( \mathrm{X5}\left( t \right) \right)^{n2}}$] X2}\\ +\ce{ \varnothing &->[$\frac{v3 \left( \mathrm{X3}\left( t \right) \right)^{n3}}{K3^{n3} + \left( \mathrm{X3}\left( t \right) \right)^{n3}}$] X3}\\ +\ce{ \varnothing &->[$\frac{v4 K4^{n4}}{K4^{n4} + \left( \mathrm{X1}\left( t \right) \right)^{n4}}$] X4}\\ +\ce{ \varnothing &->[$\frac{v5 \left( \mathrm{X2}\left( t \right) \right)^{n5}}{K5^{n5} + \left( \mathrm{X2}\left( t \right) \right)^{n5}}$] X5}\\ \ce{ X2 &<=>[$k1$][$k2$] X1 + 2 X4}\\ \ce{ X4 &<=>[$k3$][$k4$] X3}\\ \ce{ 3 X5 + X1 &<=>[$k5$][$k6$] X2}\\ @@ -74,6 +75,7 @@ raw"\begin{align} ", "\r\n"=>"\n") + r = @reaction_network begin (hill(B, p_a, k, n), d_a), 0 ↔ A (p_b, d_b), 0 ↔ B @@ -84,7 +86,7 @@ end p_a k n d_a p_b d_b r_a r_b @test latexify(r) == replace( raw"\begin{align} \require{mhchem} -\ce{ \varnothing &<=>[p_{a} \left( B\left( t \right) \right)^{n} \mathrm{inv}\left( k^{n} + \left( B\left( t \right) \right)^{n} \right)][d_{a}] A}\\ +\ce{ \varnothing &<=>[\frac{p_{a} \left( B\left( t \right) \right)^{n}}{k^{n} + \left( B\left( t \right) \right)^{n}}][d_{a}] A}\\ \ce{ \varnothing &<=>[p_{b}][d_{b}] B}\\ \ce{ 3 B &<=>[r_{a}][r_{b}] A} \end{align} @@ -93,7 +95,7 @@ raw"\begin{align} # Latexify.@generate_test latexify(r, mathjax=false) @test latexify(r, mathjax = false) == replace( raw"\begin{align} -\ce{ \varnothing &<=>[$p_{a} \left( B\left( t \right) \right)^{n} \mathrm{inv}\left( k^{n} + \left( B\left( t \right) \right)^{n} \right)$][$d_{a}$] A}\\ +\ce{ \varnothing &<=>[$\frac{p_{a} \left( B\left( t \right) \right)^{n}}{k^{n} + \left( B\left( t \right) \right)^{n}}$][$d_{a}$] A}\\ \ce{ \varnothing &<=>[$p_{b}$][$d_{b}$] B}\\ \ce{ 3 B &<=>[$r_{a}$][$r_{b}$] A} \end{align} diff --git a/test/runtests.jl b/test/runtests.jl index 1eda1558eb..d50980a649 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -32,7 +32,7 @@ end # Miscellaneous tests #@time @safetestset "Basic Plotting" begin include("plotting.jl") end -#@time @safetestset "Latexify" begin include("latexify.jl") end +@time @safetestset "Latexify" begin include("latexify.jl") end @time @safetestset "Graphs" begin include("graphs.jl") end end # @time From 1a84043f33aec3c7d32a716823bb713dd972ca79 Mon Sep 17 00:00:00 2001 From: Torkel Date: Fri, 23 Apr 2021 11:24:15 +0200 Subject: [PATCH 11/11] additional_tests --- test/custom_functions.jl | 93 +++++++++++++++++++++++++++++++++++++++- 1 file changed, 92 insertions(+), 1 deletion(-) diff --git a/test/custom_functions.jl b/test/custom_functions.jl index 45dfb105e0..9537ef9e1f 100644 --- a/test/custom_functions.jl +++ b/test/custom_functions.jl @@ -40,4 +40,95 @@ for factor in [1e-2, 1e-1, 1e0, 1e1, 1e2] @test all(abs.(f1(u0,p,t) .- f2(u0,p,t)) .< 10e-10) @test all(abs.(f1.jac(u0,p,t) .- f2.jac(u0,p,t)) .< 10e-10) @test all(abs.(g1(u0,p,t) .- g2(u0,p,t)) .< 10e-10) -end \ No newline at end of file +end + +### Tests that the various notations gives identical results ### + +# Michaelis-Menten function. +mm_network = @reaction_network begin + (1.,1.), 0 ↔ X + mm(X,v,K), 0 --> X1 + Mm(X,v,K), 0 --> X2 + MM(X,v,K), 0 --> X3 +end v K +f_mm = ODEFunction(convert(ODESystem,mm_network),jac=true) + +u0 = 10*rand(rng,length(get_states(mm_network))) +p = 10*rand(rng,length(get_ps(mm_network))) +t = 10*rand(rng) + +f_mm_output = f_mm(u0,p,t)[2:end] +f_mm_jac_output = f_mm.jac(u0,p,t)[2:end,1] +@test (maximum(f_mm_output) - minimum(f_mm_output)) .< 100*eps() +@test (maximum(f_mm_jac_output) - minimum(f_mm_jac_output)) .< 100*eps() + +# Repressing Michaelis-Menten function. +mmR_network = @reaction_network begin + (1.,1.), 0 ↔ X + mmR(X,v,K), 0 --> X1 + MmR(X,v,K), 0 --> X2 + MMR(X,v,K), 0 --> X3 +end v K +f_mmR = ODEFunction(convert(ODESystem,mmR_network),jac=true) + +u0 = 10*rand(rng,length(get_states(mmR_network))) +p = 10*rand(rng,length(get_ps(mmR_network))) +t = 10*rand(rng) + +f_mmR_output = f_mmR(u0,p,t)[2:end] +f_mmR_jac_output = f_mmR.jac(u0,p,t)[2:end,1] +@test (maximum(f_mmR_output) - minimum(f_mmR_output)) .< 100*eps() +@test (maximum(f_mmR_jac_output) - minimum(f_mmR_jac_output)) .< 100*eps() + +# Hill function. +hill_network = @reaction_network begin + (1.,1.), 0 ↔ X + hill(X,v,K,2), 0 --> X1 + Hill(X,v,K,2), 0 --> X2 +end v K +f_hill = ODEFunction(convert(ODESystem,hill_network),jac=true) + +u0 = 10*rand(rng,length(get_states(hill_network))) +p = 10*rand(rng,length(get_ps(hill_network))) +t = 10*rand(rng) + +f_hill_output = f_hill(u0,p,t)[2:end] +f_hill_jac_output = f_hill.jac(u0,p,t)[2:end,1] +@test (maximum(f_hill_output) - minimum(f_hill_output)) .< 100*eps() +@test (maximum(f_hill_jac_output) - minimum(f_hill_jac_output)) .< 100*eps() + +# Repressing Hill function. +hillR_network = @reaction_network begin + (1.,1.), 0 ↔ X + hillR(X,v,K,2), 0 --> X1 + HillR(X,v,K,2), 0 --> X2 +end v K +f_hillR = ODEFunction(convert(ODESystem,hillR_network),jac=true) + +u0 = 10*rand(rng,length(get_states(hillR_network))) +p = 10*rand(rng,length(get_ps(hillR_network))) +t = 10*rand(rng) + +f_hillR_output = f_hillR(u0,p,t)[2:end] +f_hillR_jac_output = f_hillR.jac(u0,p,t)[2:end,1] +@test (maximum(f_hillR_output) - minimum(f_hillR_output)) .< 100*eps() +@test (maximum(f_hillR_jac_output) - minimum(f_hillR_jac_output)) .< 100*eps() + +# Activation/repressing Hill function. +hillAR_network = @reaction_network begin + (1.,1.), 0 ↔ (X,Y) + hillC(X,Y,v,K,2), 0 --> X1 + HillC(X,Y,v,K,2), 0 --> X2 + hillAR(X,Y,v,K,2), 0 --> X3 + HillAR(X,Y,v,K,2), 0 --> X4 +end v K +f_hillAR = ODEFunction(convert(ODESystem,hillAR_network),jac=true) + +u0 = 10*rand(rng,length(get_states(hillAR_network))) +p = 10*rand(rng,length(get_ps(hillAR_network))) +t = 10*rand(rng) + +f_hillAR_output = f_hillAR(u0,p,t)[3:end] +f_hillAR_jac_output = f_hillAR.jac(u0,p,t)[3:end,1] +@test (maximum(f_hillAR_output) - minimum(f_hillAR_output)) .< 100*eps() +@test (maximum(f_hillAR_jac_output) - minimum(f_hillAR_jac_output)) .< 100*eps() \ No newline at end of file