Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

non-SI units not working #1037

Open
isaacsas opened this issue Sep 1, 2024 · 8 comments
Open

non-SI units not working #1037

isaacsas opened this issue Sep 1, 2024 · 8 comments
Labels

Comments

@isaacsas
Copy link
Member

isaacsas commented Sep 1, 2024

SI units works, i.e.

    @test_nowarn @reaction_network begin
        @ivs t [unit=u"s"]
        @species begin
            X1(t), [unit=u"mol/m^3"]
            Z1(t), [unit=u"mol/m^3"]
            X2(t), [unit=u"mol/m^3"]
            Z2(t), [unit=u"mol/m^3"]
            X3(t), [unit=u"mol/m^3"]
            Y3(t), [unit=u"mol/m^3"]
            Z3(t), [unit=u"mol/m^3"]
        end
        @parameters begin
            k1, [unit=u"(m^6)/(s*mol^2)"]
            v2, [unit=u"(m^6)/(s*mol^2)"]
            K2, [unit=u"mol/m^3"]
            v3, [unit=u"(m^3)/(s*mol)"]
            K3, [unit=u"mol/m^3"]
            n3
        end
        k1*X1, 2X1 --> Z1
        mm(X2, v2, K2), 3X2 --> Z2
        hill(X3, v3, K3, n3), X3 + Y3--> Z3
    end

passes, but non-SI units still have issues right now:

    @test_nowarn @reaction_network begin
        @ivs t [unit=u"s"]
        @species begin
            X1(t), [unit=u"μM"]
            Z1(t), [unit=u"μM"]
            X2(t), [unit=u"μM"]
            Z2(t), [unit=u"μM"]
            X3(t), [unit=u"μM"]
            Y3(t), [unit=u"μM"]
            Z3(t), [unit=u"μM"]
        end
        @parameters begin
            k1, [unit=u"1/(s*(μM)^2)"]
            v2, [unit=u"1/(s*(μM)^2)"]
            K2, [unit=u"μM"]
            v3, [unit=u"1/(s*μM)"]
            K3, [unit=u"μM"]
            n3
        end
        k1*X1, 2X1 --> Z1
        mm(X2, v2, K2), 3X2 --> Z2
        hill(X3, v3, K3, n3), X3 + Y3--> Z3
    end

gives

┌ Warning: Reaction rate laws are expected to have units of 0.0009999999999999998 m⁻³ s⁻¹ mol however, k1*X1(t), 2*X1 --> Z1 has units of 9.999999999999995e-7 m⁻³ s⁻¹ mol.
└ @ Catalyst ~/.julia/dev/Catalyst/src/reactionsystem.jl:1492
┌ Warning: Reaction rate laws are expected to have units of 0.0009999999999999998 m⁻³ s⁻¹ mol however, Catalyst.mm(X2(t), v2, K2), 3*X2 --> Z2 has units of 9.999999999999992e-10 m⁻³ s⁻¹ mol.
└ @ Catalyst ~/.julia/dev/Catalyst/src/reactionsystem.jl:1492
┌ Warning: Reaction rate laws are expected to have units of 0.0009999999999999998 m⁻³ s⁻¹ mol however, Catalyst.hill(X3(t), v3, K3, n3), X3 + Y3 --> Z3 has units of 9.999999999999995e-7 m⁻³ s⁻¹ mol.
└ @ Catalyst ~/.julia/dev/Catalyst/src/reactionsystem.jl:1492
Test Failed at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-R17H3W25T9.0/build/default-honeycrisp-R17H3W25T9-0/julialang/julia-release-1-dot-10/usr/share/julia/stdlib/v1.10/Test/src/Test.jl:903
  Expression: isempty(stderr_content)
   Evaluated: isempty("┌ Warning: Reaction rate laws are expected to have units of 0.0009999999999999998 m⁻³ s⁻¹ mol however, k1*X1(t), 2*X1 --> Z1 has units of 9.999999999999995e-7 m⁻³ s⁻¹ mol.\n└ @ Catalyst ~/.julia/dev/Catalyst/src/reactionsystem.jl:1492\n┌ Warning: Reaction rate laws are expected to have units of 0.0009999999999999998 m⁻³ s⁻¹ mol however, Catalyst.mm(X2(t), v2, K2), 3*X2 --> Z2 has units of 9.999999999999992e-10 m⁻³ s⁻¹ mol.\n└ @ Catalyst ~/.julia/dev/Catalyst/src/reactionsystem.jl:1492\n┌ Warning: Reaction rate laws are expected to have units of 0.0009999999999999998 m⁻³ s⁻¹ mol however, Catalyst.hill(X3(t), v3, K3, n3), X3 + Y3 --> Z3 has units of 9.999999999999995e-7 m⁻³ s⁻¹ mol.\n└ @ Catalyst ~/.julia/dev/Catalyst/src/reactionsystem.jl:1492\n")

ERROR: There was an error during testing

So it seems there are still issues with unit conversion to SI going on.

@isaacsas isaacsas added the bug label Sep 1, 2024
@isaacsas
Copy link
Member Author

isaacsas commented Sep 1, 2024

This really requires someone to dig into MTK and DynamicQuantites to understand how they are handling numeric coefficients that arise when converting and comparing units, and ensure that only numeric components arising from unit conversions are getting compared (i.e. numeric components that are part of expressions like X^2/2 should get dropped probably before conversion, including when tracing through registered functions).

@ctessum
Copy link

ctessum commented Sep 3, 2024

MTK version 9.34 supports non-SI units as described here: SciML/ModelingToolkit.jl#2621 . However Catalyst is currently pinned to MTK v9.32. So it may just be as simple as updating that dependency.

@isaacsas
Copy link
Member Author

isaacsas commented Sep 3, 2024

Catalyst isn't pinned to 9.32; are you sure this isn't something you did locally? I currently have ModelingToolkit 9.35.1 with Catalyst 14.4 locally.

Note that the line

ModelingToolkit = "9.32"

in the Project.toml means support is declared for any 9.X release with X >= 32 (see https://pkgdocs.julialang.org/v1/compatibility/).

@isaacsas
Copy link
Member Author

isaacsas commented Sep 3, 2024

I suspect the issue here is that ModelingToolkit converts non-SI units to SI internally it seems, and that can lead to expressions where numeric coefficients in a symbolic expression get merged with numeric coefficients generated from the unit conversions (like in converting micromolar to mol/m^3). There is then no effective way to unit check after this point. I think expressions need to have numeric coefficients removed prior to unit conversion to avoid this issue.

@isaacsas
Copy link
Member Author

isaacsas commented Sep 3, 2024

Or perhaps a way is needed for users to indicate how registered functions transform units to avoid having to try to trace through them.

@ctessum
Copy link

ctessum commented Sep 3, 2024

Catalyst isn't pinned to 9.32; are you sure this isn't something you did locally? I currently have ModelingToolkit 9.35.1 with Catalyst 14.4 locally.

Note that the line

ModelingToolkit = "9.32"

in the Project.toml means support is declared for any 9.X release with X >= 32 (see https://pkgdocs.julialang.org/v1/compatibility/).

Oh, sorry, I guess I should have read the manual more carefully :)

@isaacsas
Copy link
Member Author

isaacsas commented Sep 3, 2024

No worries! I still get messed up on aspects of semver in Julia (and finding out where to even read about it can take a bit of effort).

@ctessum
Copy link

ctessum commented Sep 4, 2024

This may be related to SciML/ModelingToolkit.jl#3017

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

No branches or pull requests

2 participants