Units don't seem to work correctly

[ctessum]

Hello,

Here is a small example with units:

using ModelingToolkit, MethodOfLines, OrdinaryDiffEq, DomainSets
using Unitful

function ModelingToolkit.get_unit(a::Pair{Num, Float64})
    ModelingToolkit.get_unit(a[1])
end

@parameters x [unit=u"m"] 
@parameters t [unit=u"s"]
@variables u(..) [unit=u"kg"]
Dt = Differential(t)
Dx = Differential(x)

x_min = t_min = 0.0
x_max = t_max = 1.0

@parameters α = 10. [unit=u"m/s"]

eq = [Dt(u(x,t)) ~ α * Dx(u(x,t))]

@assert ModelingToolkit.get_unit(eq[1].lhs) == u"kg/s"
@assert ModelingToolkit.get_unit(eq[1].rhs) == u"kg/s"

domains = [x ∈ Interval(x_min, x_max),
              t ∈ Interval(t_min, t_max)]

bcs = [u(x,t_min) ~ 0.0,
       u(x_min,t) ~ 0.0,
       u(x_max,t) ~ sin(1.0),
] 

@named pdesys = PDESystem(eq,bcs,domains,[x,t],[u(x,t)], [α => 10.])

discretization = MOLFiniteDifference([x=>6], t)
@time prob = discretize(pdesys,discretization)

It works if there aren't any units, and as you can see the units on the two sides of the equation match. However, running the code to the end gives this error: ERROR: ModelingToolkit.ValidationError("Some equations had invalid units. See warnings for details.")

Here are the warnings in question:

┌ Warning:  in eq. #1left, in comparison >, units [m s⁻¹] and [] do not match.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/rrbUl/src/systems/validation.jl:152
┌ Warning:  in eq. #2left, in comparison >, units [m s⁻¹] and [] do not match.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/rrbUl/src/systems/validation.jl:152
┌ Warning:  in eq. #3left, in comparison >, units [m s⁻¹] and [] do not match.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/rrbUl/src/systems/validation.jl:152
┌ Warning:  in eq. #4left, in comparison >, units [m s⁻¹] and [] do not match.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/rrbUl/src/systems/validation.jl:152

It looks like what is happening is that the upwind differencing is comparing the velocity to zero, but the zero doesn't have any units, which causes the problem. I also tried advection_scheme=WENOScheme(), but it gave a similar error.

[aditya-sengupta]

Running into the same issue - unit validation passes on the PDESystem but not on the discretized system. I can post an MWE for my case if it'd help

[xtalax]

What would be the successful case, that the result array has the right units?

[aditya-sengupta]

Potentially it'd be running unit validation on the PDESystem and trusting that that correctness is preserved on the discretization. The units in the result array being correct is ensured by the user declaring it as a variable with those units and it passing the PDESystem validation - it seems like it'd run just as correctly and with fewer hurdles if the discretized system didn't have to be unit checked (or we figure out autopromotion rules for constants to unit aware constants, but that seems harder)

[ctessum]

I guess the immediate successful case would be that the example above would run without error. As far as I can tell it's just an issue of the constants that are added being units (e.g. when the upwind scheme checks for u > 0, the 0 is unitless but u is not).