@@ -78,4 +78,140 @@ const D_unitful = MT.Differential(t_unitful)
7878
7979#  Extension loaded - all Unitful-specific functionality is now available
8080
81- end  #  module
81+ end  #  module
82+ 
83+ #  Create the UnitfulUnitCheck module inside ModelingToolkit for backward compatibility
84+ @eval  ModelingToolkit module  UnitfulUnitCheck
85+ 
86+ using  ModelingToolkit, Symbolics, SciMLBase, Unitful, RecursiveArrayTools
87+ using  ModelingToolkit:  ValidationError, Connection, instream, JumpType, VariableUnit,
88+                          get_systems, Conditional, Comparison, Integral, Differential
89+ using  JumpProcesses:  MassActionJump, ConstantRateJump, VariableRateJump
90+ using  Symbolics:  Symbolic, value, issym, isadd, ismul, ispow, iscall, operation, arguments, getmetadata
91+ 
92+ const  MT =  ModelingToolkit
93+ 
94+ Base.:* (x:: Union{MT.Num, Symbolic} , y:: Unitful.AbstractQuantity ) =  x *  y
95+ Base.:/ (x:: Union{MT.Num, Symbolic} , y:: Unitful.AbstractQuantity ) =  x /  y
96+ 
97+ """ 
98+ Throw exception on invalid unit types, otherwise return argument. 
99+ """ 
100+ function  screen_unit (result)
101+     result isa  Unitful. Unitlike || 
102+         throw (ValidationError (" Unit must be a subtype of Unitful.Unitlike, not $(typeof (result)) ." 
103+     result isa  Unitful. ScalarUnits || 
104+         throw (ValidationError (" Non-scalar units such as $result  are not supported. Use a scalar unit instead." 
105+     result ==  Unitful. u "" && 
106+         throw (ValidationError (" Degrees are not supported. Use radians instead." 
107+     result
108+ end 
109+ 
110+ """ 
111+ Test unit equivalence. 
112+ """ 
113+ equivalent (x, y) =  isequal (1  *  x, 1  *  y)
114+ const  unitless =  Unitful. unit (1 )
115+ 
116+ """ 
117+ Find the unit of a symbolic item. 
118+ """ 
119+ get_unit (x:: Real ) =  unitless
120+ get_unit (x:: Unitful.Quantity ) =  screen_unit (Unitful. unit (x))
121+ get_unit (x:: AbstractArray ) =  map (get_unit, x)
122+ get_unit (x:: MT.Num ) =  get_unit (value (x))
123+ function  get_unit (x:: Union{Symbolics.ArrayOp, Symbolics.Arr, Symbolics.CallWithMetadata} )
124+     get_literal_unit (x)
125+ end 
126+ get_unit (op:: Differential , args) =  get_unit (args[1 ]) /  get_unit (op. x)
127+ get_unit (op:: typeof (getindex), args) =  get_unit (args[1 ])
128+ get_unit (x:: SciMLBase.NullParameters ) =  unitless
129+ get_unit (op:: typeof (instream), args) =  get_unit (args[1 ])
130+ 
131+ get_literal_unit (x) =  screen_unit (getmetadata (x, VariableUnit, unitless))
132+ 
133+ function  get_unit (op, args) #  Fallback
134+     result =  op (1  .*  get_unit .(args)... )
135+     try 
136+         Unitful. unit (result)
137+     catch 
138+         throw (ValidationError (" Unable to get unit for operation $op  with arguments $args ." 
139+     end 
140+ end 
141+ 
142+ function  get_unit (op:: Integral , args)
143+     unit =  1 
144+     if  op. domain. variables isa  Vector
145+         for  u in  op. domain. variables
146+             unit *=  get_unit (u)
147+         end 
148+     else 
149+         unit *=  get_unit (op. domain. variables)
150+     end 
151+     return  get_unit (args[1 ]) *  unit
152+ end 
153+ 
154+ function  get_unit (op:: Conditional , args)
155+     terms =  get_unit .(args)
156+     terms[1 ] ==  unitless || 
157+         throw (ValidationError (" , in $op , [$(terms[1 ]) ] is not dimensionless." 
158+     equivalent (terms[2 ], terms[3 ]) || 
159+         throw (ValidationError (" , in $op , units [$(terms[2 ]) ] and [$(terms[3 ]) ] do not match." 
160+     return  terms[2 ]
161+ end 
162+ 
163+ function  get_unit (op:: typeof (Symbolics. _mapreduce), args)
164+     if  args[2 ] ==  + 
165+         get_unit (args[3 ])
166+     else 
167+         throw (ValidationError (" Unsupported array operation $op " 
168+     end 
169+ end 
170+ 
171+ function  get_unit (op:: Comparison , args)
172+     terms =  get_unit .(args)
173+     equivalent (terms[1 ], terms[2 ]) || 
174+         throw (ValidationError (" , in comparison $op , units [$(terms[1 ]) ] and [$(terms[2 ]) ] do not match." 
175+     return  unitless
176+ end 
177+ 
178+ function  get_unit (x:: Symbolic )
179+     if  issym (x)
180+         get_literal_unit (x)
181+     elseif  isadd (x)
182+         terms =  get_unit .(arguments (x))
183+         firstunit =  terms[1 ]
184+         for  other in  terms[2 : end ]
185+             termlist =  join (map (repr, terms), " , " 
186+             equivalent (other, firstunit) || 
187+                 throw (ValidationError (" , in sum $x , units [$termlist ] do not match." 
188+         end 
189+         return  firstunit
190+     elseif  ispow (x)
191+         pargs =  arguments (x)
192+         base, expon =  get_unit .(pargs)
193+         @assert  expon isa  Unitful. DimensionlessUnits
194+         if  base ==  unitless
195+             unitless
196+         else 
197+             pargs[2 ] isa  Number ?  base^ pargs[2 ] :  (1  *  base)^ pargs[2 ]
198+         end 
199+     elseif  iscall (x)
200+         op =  operation (x)
201+         if  issym (op) ||  (iscall (op) &&  iscall (operation (op))) #  Dependent variables, not function calls
202+             return  screen_unit (getmetadata (x, VariableUnit, unitless)) #  Like x(t) or x[i]
203+         elseif  iscall (op) &&  ! iscall (operation (op))
204+             gp =  getmetadata (x, Symbolics. GetindexParent, nothing ) #  Like x[1](t)
205+             return  screen_unit (getmetadata (gp, VariableUnit, unitless))
206+         end   #  Actual function calls:
207+         args =  arguments (x)
208+         return  get_unit (op, args)
209+     else  #  This function should only be reached by Terms, for which `iscall` is true
210+         throw (ArgumentError (" Unsupported value $x ." 
211+     end 
212+ end 
213+ 
214+ #  Re-use validation functions from main package
215+ using  ModelingToolkit:  safe_get_unit, _validate, validate
216+ 
217+ end  #  module UnitfulUnitCheck
0 commit comments