Automatic Uncertainty Quantification, Arbitrary Precision, and Unit Checking in ODE Solutions using Julia's Type System

One of the nice things about DifferentialEquations.jl is that it is designed with Julia's type system in mind. What this means is, if you have properly defined a Number type, you can use this number type in DifferentialEquations.jl's algorithms! There's more than a few useful/interesting types that can be used:

and on and on. That's only a subset of types people have effectively used on the SciML tools.

We will look into the BigFloat, Measurement, and Unitful cases to demonstrate the utility of alternative numerical types.

How Type Support Works in DifferentialEquations.jl / SciML

DifferentialEquations.jl determines the numbers to use in its solvers via the types that are designated by tspan and the initial condition u0 of the problem. It will keep the time values in the same type as tspan, and the solution values in the same type as the initial condition.

Let's use this feature in some cool ways!

Arbitrary Precision: Rationals and BigFloats

Let's solve the linear ODE. First, define an easy way to get ODEProblems for the linear ODE:

using DifferentialEquations
f(u, p, t) = p * u
prob_ode_linear = ODEProblem(f, 1 / 2, (0.0, 1.0), 1.01);

ODEProblem with uType Float64 and tType Float64. In-place: false
timespan: (0.0, 1.0)
u0: 0.5

Next, let's solve it using Float64s. To do so, we just need to set u0 to a Float64 (which is done by the default) and dt should be a float as well.

sol = solve(prob_ode_linear, Tsit5())

retcode: Success
Interpolation: specialized 4th order "free" interpolation
t: 5-element Vector{Float64}:
 0.0
 0.09964258706516003
 0.3457024214672858
 0.6776921708765304
 1.0
u: 5-element Vector{Float64}:
 0.5
 0.552938681151017
 0.7089376222328616
 0.991359430285871
 1.3728004409033048

Notice that both the times and the solutions were saved as Float64. Let's change the state to use BigFloat values. We do this by changing the u0 to use BigFloats like:

prob_ode_linear_bigu = ODEProblem(f, big(1 / 2), (0.0, 1.0), 1.01);
sol = solve(prob_ode_linear_bigu, Tsit5())

retcode: Success
Interpolation: specialized 4th order "free" interpolation
t: 5-element Vector{Float64}:
 0.0
 0.09964258706516002
 0.34570242146728575
 0.6776921708765303
 1.0
u: 5-element Vector{BigFloat}:
 0.5
 0.5529386811510169868158546583396554107562766958646715721373574023383501046221058
 0.7089376222328617224502897748356594658455684696955511956359006091387463746837088
 0.9913594302858714240183997742098005064281191747328421549570415223008330473386488
 1.372800440903305954330299558342953663171497171911604970344459251505789637218533

Now we see that u is in arbitrary precision BigFloats, while t is in Float64. We can then change t to be arbitrary precision BigFloats by changing the types of the tspan like:

prob_ode_linear_big = ODEProblem(f, big(1 / 2), (big(0.0), big(1.0)), 1.01);
sol = solve(prob_ode_linear_big, Tsit5())

retcode: Success
Interpolation: specialized 4th order "free" interpolation
t: 5-element Vector{BigFloat}:
 0.0
 0.09964258706516002757241113786493079821949560824572094442904017056206815655027755
 0.3457024214672857657475098820237325621051648241581980706974757091926545753310105
 0.677692170876530371514288101002609220871911443860015471896386211177424067229719
 1.0
u: 5-element Vector{BigFloat}:
 0.5
 0.5529386811510169905121549222942681803517017137036250612133811084046740072664042
 0.7089376222328617354344942756503376924689703213410589145373354448957336352383165
 0.9913594302858714646480239922983656519745236607948619810174815468913895034199903
 1.372800440903305954330308294683133150311893252665247472635661549462376339096003

Now let's send it into the bizarre territory. Let's use rational values for everything. Let's start by making the time type Rational. Rationals are incompatible with adaptive time stepping since they do not have an L2 norm (this can be worked around by defining internalnorm, but we will skip that in this tutorial). To account for this, let's turn off adaptivity as well. Thus the following is a valid use of rational time (and parameter):

prob = ODEProblem(f, 1 / 2, (0 // 1, 1 // 1), 101 // 100);
sol = solve(prob, RK4(), dt = 1 // 2^(6), adaptive = false)

retcode: Success
Interpolation: 3rd order Hermite
t: 65-element Vector{Rational{Int64}}:
   0
  1//64
  1//32
  3//64
  1//16
  5//64
  3//32
  7//64
  1//8
  9//64
   ⋮
  7//8
 57//64
 29//32
 59//64
 15//16
 61//64
 31//32
 63//64
   1
u: 65-element Vector{Float64}:
 0.5
 0.5079532157789419
 0.5160329388403366
 0.5242411814636141
 0.5325799879363893
 0.5410514350635981
 0.549657632684732
 0.5584007241993002
 0.5672828871006491
 0.5763063335182743
 ⋮
 1.2099787760366485
 1.2292252206241676
 1.2487778074652507
 1.268641406190701
 1.288820963889771
 1.3093215063422494
 1.3301481392701477
 1.3513060496092948
 1.3728005068011595

Now let's change the state to use Rational{BigInt}. You will see that we will need to use the arbitrary-sized integers because... well... there's a reason people use floating-point numbers with ODE solvers:

prob = ODEProblem(f, BigInt(1) // BigInt(2), (0 // 1, 1 // 1), 101 // 100);
sol = solve(prob, RK4(), dt = 1 // 2^(6), adaptive = false)

retcode: Success
Interpolation: 3rd order Hermite
t: 65-element Vector{Rational{Int64}}:
   0
  1//64
  1//32
  3//64
  1//16
  5//64
  3//32
  7//64
  1//8
  9//64
   ⋮
  7//8
 57//64
 29//32
 59//64
 15//16
 61//64
 31//32
 63//64
   1
u: 65-element Vector{Rational{BigInt}}:
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                         1//2
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                         13635265310428667//26843545600000000
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                         185920460085779372649467295396889//360287970189639680000000000000000
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                         2535074799906565070064434989240252460736488216963//4835703278458516698824704000000000000000000000000
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                         34566417478507880850358857402427148787055026137565875648030878321//64903710731685345356631204115251200000000000000000000000000000000
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                         471322273150493660801770967151032391984466213543795788952681156195576797027228107//871122859317602466466238995025326621327360000000000000000000000000000000000000000
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                         6426604241121310927369897121054534260287845822775576065191752910027218849129163851685338960943369//11692013098647223345629478661730264157247460343808000000000000000000000000000000000000000000000000
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                         87628453872815159549919436744286633306700458799675030242287144817323241473991649398590738650552963396600101159123//156927543384667019095894735580191660402558886111600862822400000000000000000000000000000000000000000000000000000000
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                        1194837217298495123489429765132896187014406474671976544708139096975464385034479106227218239981906730284477718437811007486707779041//2106245833371143733958360553673408646377901908010982225086219550720000000000000000000000000000000000000000000000000000000000000000
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                        16291922460639289757947281837298870307759032284133099989470765273134679312684297502599368900453650828381789478033888032849237996811196931228168347//28269553036454149273332760011886696253239742350009903329945699220681916416000000000000000000000000000000000000000000000000000000000000000000000000
                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                          ⋮
                                                                                                                                  3476651764144701299310823013868414852914252425511969801823326943470139627183050321356269729758757892616986280081068426854388132587192659833281759109495586510217757379507472058821725242337421483829232094724350000335051444864163716743131954661196447426549808629937028209252581525400046451480100419001530953781616768025169970806686494551090052818656160455734295279711817930495039244376027618813690810243490836308240922572885036008357232995165091370374018505314343843108268301984782218425989406930319958551532270117232742567343778042005440462725562350111368484320102066864357344195777652146963879004877113342262696777549666508001903654324815937791481538838675335919013649404331120990596974956064711019198691285927797569234244700072843046601685002725165218418237056611786017220421490343908940716757353909791039304274005391221502042135010059552204316515216866746639520694954996904124971579780497757666750446881//2873316320086756454205652005706207886829487380916701497268289374190875758638147354923588066229040681914626081324682241425350209449285219103061691508162945093459145235278159647906758878132917284580285004541144156524655039141064404379679426100372681576010961105039478156124917692281743798746817307433826768839799128901030530861680784239247749388313984379590631715911498407693680085160154180294124335267854340370884940503485168233544639735700169786044364881920000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
                                                                                                                  47405069196082873328585408620421318517516390291920914673056124180124925970261893410339173429287494238135595611089781966147422284876218622543131269354235054173066292325370254052047479196233718967762509800342189006216218284815334686783879790109432868299224476127915815042082361972609813617003873110744660057301821702132621624893170611073772099444022185151482507338829550979383225592740578974181160289252871652171292405008199989659405158877156438081779650473274842639963837184419397179015646801970478919138945275747462931559006763638612060474604554509550795966812356642042917363452960338457115020031120789490076363823059477978169077204763046262491762113637668511018482223361774121845706177451044161738746050766877434757946831223565527239248794315463305271779949208058882442095280247508627968155300497343796852202862275891120750390256408617409094008210261398024974884224628684715716894544417457681835995541069229911923137627//38564998830736521417281865696453025806593491967131023221754800625044118265468851210705360385717536794615180260494208076605798671660719333199513807806252394423283413430106003596332513246682903994829528690198205120921557533726473585751382193953592127439965050261476810842071573684505878854588706623484573925925903505747545471088867712185004135201289273405614415899438276535626346098904241020877974002916168099951885406379295536200413493190419727789712076165429067776000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
                                                                                                  646380695547819379379462110595350483323442009248798219582259024446464544132551698473809581856967370658487548204659626702292564992976158496529375784071763787592253430041646066314009418251081347658033077807894727674506660673057920454240886643468816490881054458605117293952007109054634081637211268136338650821833722847677488664849766203163787324732106652618777185431293544445474668250279073481131086314724562094446770240396090985922895943040057962860047030005371781373835035785147257040477418516657611117881834875218299885716266834779541157085670534199549342497110048855301431141761423177209449504715692545645531288355844442793377192971285135291013238861059535547536381876669229982900292241340160011803760056601165836520905723932492682377301319454029065863067972071145347687240835433388266357742534589282848560823758438409092468036851229133683325521477136527739510675196517661031060688353455828921203165442555541627063437450443175007153209//517610652338411247125091194937906278248634591141457661513916951298640145335453006427092274839221143832965212887398076520128024933571973574771371398838386057405297005067351460185957906049968291062814308810921896100805671839805065813156369093210647278367856544550122347590861343331896786250613857746265095538581506488867049470623751042202917069732186894727238934606968195284747923033659051938664223591447345684161896856058575310936125171656160723031361663708684162205273292800000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
                                                                                  8813572275333935103416906132960621777114439444141647438438393753862988179987938832721184273969585102967467754980010070157566337119239854562632863958739529252701552386911176964405034482000877845003064456723043889756774450867335661101187923620690769719662220354335695798311755154012564019508697637433705198641352232164934271731689599505104756286149010104803971173588887614751600240293596938071050020381715583108625023759753393738598032521536448316590020179966102744430869710510004750450586729772303515586404571938281315619339180632083313992191945702244726898267871642409305486953268724327679433424951407071104066109102285627362450306825448652531631455682741766958888491690416428140178416796287551293015562399117395448453688862961061347809034853245381286200796207330060769671545355757344311231491079771497941834167729303081908176522671807233326506847483489756135442309818267518863626272611798974800987095037233771816523002200031070842864277033671134642403//6947252574545944471877627197737088174346755392521537393659097358399012979651430037341372277527182721482180237678288966210172977571138121768185158859627001581181653918534439996364172765366419849849364181449231848010219624385221589633232036842015330311191725363944820361567168306503513658105903059201890160889134614375299267401072260772210764407187169348166518951608784413937917928229658547784150744532806296935881538019052501315074027366329988944810663524927964212077525097869475840000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
                                                                 120175396306816661567962374621260628614816898908130603456791258506661752158816764138355946608581111648754811275356266553194146487289970011521206302938350516439692698830565485671063124957581001985813711933633645259078108062496048585468892270573397374250419348298467866026065051327496126625966211919529492809946883329115587011078666703280156236855586771535787690208072551405979966710615670765008341976837702247734425622133716186858699168908901364438717264342619451907519562821993757807949284311972035986557650598491323972524603931519064752851585706054312058594719949638576121660783232383108957295631067768960287620832338037163277354938280475722726742662910468041915790311921524235141338322754706846323314939020948898862794933449368601120381306723818943912893036079998982040701877643402517997677027213445418199720233720320010396665109155653989028150312257462380364161291864607152883841218256979741199145521849881542278617665978410806969940456720443825181089585213684966801//93244445639770729862957501651127871609648939295598894004396565397511723957132517156291414749188392311819502398767993997219818753659311707791315470545820707967860914422798962606432752916695795974088280267869564598517299876599783853712075728174759257553768635074864690629768384949250922720794309199432719070009410782760880699063638360466965433587592874066215113935387298888055887737745183249936814734070168263882738171006284743922624810491863481443760464848830272021035557849901868779569152000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
                                                 1638623412429354568616301751922853165062178143570237809868003552419139731525207008014840385057094808550364022265598937908854560007761359201302534206737402807692400464255340884583551113600457074478398193829260671903408623257256918638344550520063772668911716725034710930551991674894972611899808216730842062998272110294468577188489521117566036367472120519200995130413089843947574203037759878714438894957215775888827327087684601532223092366221531321381581156578870580679772475683937879311488116796619715835741276270959196507639605482261723621713433086664848878328456065693446846127336706287037161239359880305409587491274511777885262373075163575581721009527514233124879843185501554931387437408536169975692089151929986328911859307334752043996902855680228690265685371832496623752200243917392455804140656929928494581593961635191153288675520454159469305085294704574693513943256409827056994632177252792094052051470442150934455756986752540284228302747677207523056397782475349694314971666973684267//1251505764238953380310790723217063156492278350986520395258292901865744044288949584763845459354175586006508115805319015342448238968794450346355026058791097531879778895383050816179236228126428300279367584898068436076222415813190335413631917039956177449583366403740944868375040179411785414948859053607735844243893605388086697470737326815552111955066160454615751415166882205681178758918321833875017541716382217695604557573133901205035650986064847896555168936823386256828536078182639385054790440032665600000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
                                 22343064972354224957479728113350949780887535331771171875613387648500150040113223878247360980060744795749139600967047698773976023825716749550816790048959630214024291450253822775597752506859958652993156369432205621106046955897352932772730437075711261897922223752076981370603845666571148563739751397868267477530534170057941586126867460170518844864494697419135689305804862898660459973467138685238383474453740580946286484621665119314190873607079738136863819984192230557369652571594272717804907202344994813159997266200695798241334728503709447203221423207356469728025768552535361151641594813114163751201865741331096235472064289773307514336130797690384014695762583913450057967476452114645805950514182594981422860829708262919876496710557833096659829409940579069682311025892212895827818852274344280241244125649458665933112767576600946408242226001592726350162502705017014799678405577891527233732317745485132196617377835361746048692485877388072155967761600004670740471374233078485232398111732065728284486083682089//16797426025505597180323426475367106769690204981299732607723004644694712665399418877354675409763379446686211251695080855446054459199265402449658467899173553735531149248066277025612272731441890321839848252186585707866380747311767925077561613338740335684791398530164032080655694478936221481797309836746050823857727758891755480154571048997687953221061814590197326795648367049015650538387689411549829040931804163870543930391122605951634723435058354534121380135661044065795059797970422730537112837730462261575680000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
                 304653618746195407227163730322053916014053723110422517841663212509580133414870563836575120110686536695480977952772624365535866244332873061241508781399832428387630674513930691330101904718500293023050050607427115063669681175960942195506424032272471651202997578108749930421242453943184047368369010644749823498643293649551273090398280784531128183631032000069615743032947667239062468509040079935236177208481442271959808690696140664325113131314265183037989326932456072811588119452547655042672272627027134583718064095094039463733284097409615865832861209712359926308642772539377944587409962967122237590175479982785779905543006714406871097094026357618247456813205608858618077379847671759450329112666484299745038339155180488295039857022873030587080343463926426043822286460042109518151839354796204791318240797287518363486427798004869303012233120170679276644118435222844774404863091415468587316888561285613163612304906504998288703317237135069355519646441949776015247310544320628566775036316038856240879188327554671006066340045363//225451235739143130482621660669882103656123857644433259761609693674437158756273421423885518367590980693612040313054990119426585607799410158579879393738800746006890372530436409579627705493048468302253321227325456578709735148712759856718454352234022223757002566466115585318472006322496350419162828280010585455071241504098321847789561501104635433430119850185353627418544411222894525170437245598787501266228575781564409246128662754026769046536188789298827011803275711261889543970772915009224750475981532113842946965504000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
 4154032919386558883432944248380343483762044089219885824293861963690668280133800624271545564442460641100421478069957127705133139131053171319939289915624722195403241736871340745589519387833493153871994750550507166424767604170338332253959630697516305444248796250106488696552824425774652891031781638156634640665726706553562695794716367646798636566490125595141712720380867485868916531456648814528917577693417533965049279568879801863167212171389128029079788394889712773514836798543384276326561054294342851708282050876790968869065128360584151770000714515194551497614161342119347668187950856166437783338125107242946094385126468080818490755092469614835748767521966870937090173768929887202086899128132689201712566935821453568568851761907310360889009454819233203019261511646422045122043461427963067831419822632761257565485308244276118163333934078610669354885645888806741789229076806586507072844471249752898840782835318816592414922484506856439857852070928805249944302969170900303083044962139908567605824428891872081720287044135359380045755621121//3025955263570019164018502277869853398058543745963126397283707470775892712704232437030043920740033026198847216426264951289188498307633591122471111874163926157374989814610878574225506571713008520940845805558579429855707382314196875257835647882856218717417250856125102284683546912020709544155188247379716859572950811281937944702307676679453365814328593305957854274867553594143460475201489987084725797475032257007739929467758191052369579260681352907875927458926484892315482757871323905647524505025315981027903769053444125491200000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000

Yeah...

sol[end]

4154032919386558883432944248380343483762044089219885824293861963690668280133800624271545564442460641100421478069957127705133139131053171319939289915624722195403241736871340745589519387833493153871994750550507166424767604170338332253959630697516305444248796250106488696552824425774652891031781638156634640665726706553562695794716367646798636566490125595141712720380867485868916531456648814528917577693417533965049279568879801863167212171389128029079788394889712773514836798543384276326561054294342851708282050876790968869065128360584151770000714515194551497614161342119347668187950856166437783338125107242946094385126468080818490755092469614835748767521966870937090173768929887202086899128132689201712566935821453568568851761907310360889009454819233203019261511646422045122043461427963067831419822632761257565485308244276118163333934078610669354885645888806741789229076806586507072844471249752898840782835318816592414922484506856439857852070928805249944302969170900303083044962139908567605824428891872081720287044135359380045755621121//3025955263570019164018502277869853398058543745963126397283707470775892712704232437030043920740033026198847216426264951289188498307633591122471111874163926157374989814610878574225506571713008520940845805558579429855707382314196875257835647882856218717417250856125102284683546912020709544155188247379716859572950811281937944702307676679453365814328593305957854274867553594143460475201489987084725797475032257007739929467758191052369579260681352907875927458926484892315482757871323905647524505025315981027903769053444125491200000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000

That's one huge fraction! 0 floating-point error ODE solve achieved.

Unit Checked Arithmetic via Unitful.jl

Units and dimensional analysis are standard tools across the sciences for checking the correctness of your equation. However, most ODE solvers only allow for the equation to be in dimensionless form, leaving it up to the user to both convert the equation to a dimensionless form, punch in the equations, and hopefully not make an error along the way.

DifferentialEquations.jl allows for one to use Unitful.jl to have unit-checked arithmetic natively in the solvers. Given the dispatch implementation of the Unitful, this has little overhead because the unit checks occur at compile-time and not at runtime, and thus it does not have a runtime effect unless conversions are required (i.e. converting cm to m), which automatically adds a floating-point operation for the multiplication.

Let's see this in action.

Using Unitful

To use Unitful, you need to have the package installed. Then you can add units to your variables. For example:

using Unitful
t = 1.0u"s"

1.0 s

Notice that t is a variable with units in seconds. If we make another value with seconds, they can add:

t2 = 1.02u"s"
t + t2

2.02 s

and they can multiply:

t * t2

1.02 s^2

You can even do rational roots:

sqrt(t)

1.0 s^1/2

Many operations work. These operations will check to make sure units are correct, and will throw an error for incorrect operations:

t + sqrt(t)

Using Unitful with DifferentialEquations.jl

Just like with other number systems, you can choose the units for your numbers by simply specifying the units of the initial condition and the timespan. For example, to solve the linear ODE where the variable has units of Newton's and t is in seconds, we would use:

using DifferentialEquations
f(u, p, t) = 0.5 * u
u0 = 1.5u"N"
prob = ODEProblem(f, u0, (0.0u"s", 1.0u"s"))
#sol = solve(prob,Tsit5())

ODEProblem with uType Unitful.Quantity{Float64, 𝐋 𝐌 𝐓^-2, Unitful.FreeUnits{(N,), 𝐋 𝐌 𝐓^-2, nothing}} and tType Unitful.Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}. In-place: false
timespan: (0.0 s, 1.0 s)
u0: 1.5 N

Notice that we received a unit mismatch error. This is correctly so! Remember that for an ODE:

\[\frac{dy}{dt} = f(t,y)\]

we must have that f is a rate, i.e. f is a change in y per unit time. So, we need to fix the units of f in our example to be N/s. Notice that we then do not receive an error if we do the following:

f(y, p, t) = 0.5 * y / 3.0u"s"
prob = ODEProblem(f, u0, (0.0u"s", 1.0u"s"))
sol = solve(prob, Tsit5())

retcode: Success
Interpolation: specialized 4th order "free" interpolation
t: 3-element Vector{Unitful.Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}}:
                 0.0 s
 0.14311598261241779 s
                 1.0 s
u: 3-element Vector{Unitful.Quantity{Float64, 𝐋 𝐌 𝐓^-2, Unitful.FreeUnits{(N,), 𝐋 𝐌 𝐓^-2, nothing}}}:
                1.5 N
 1.5362091208988309 N
 1.7720406194871123 N

This gives a normal solution object. Notice that the values are all with the correct units:

print(sol[:])

Unitful.Quantity{Float64, 𝐋 𝐌 𝐓^-2, Unitful.FreeUnits{(N,), 𝐋 𝐌 𝐓^-2, nothing}}[1.5 N, 1.5362091208988309 N, 1.7720406194871123 N]

And when we plot the solution, it automatically adds the units:

using Plots
gr()
plot(sol, lw = 3)

Measurements.jl: Numbers with Linear Uncertainty Propagation

The result of a measurement should be given as a number with an attached uncertainty, besides the physical unit, and all operations performed involving the result of the measurement should propagate the uncertainty, taking care of correlation between quantities.

There is a Julia package for dealing with numbers with uncertainties: Measurements.jl. Thanks to Julia's features, DifferentialEquations.jl easily works together with Measurements.jl out-of-the-box.

Let's try to automate uncertainty propagation through number types on some classical physics examples!

Warning about Measurement type

Before going on with the tutorial, we must point up a subtlety of Measurements.jl that you should be aware of:

using Measurements
5.23 ± 0.14 === 5.23 ± 0.14

false

(5.23 ± 0.14) - (5.23 ± 0.14)

\[0.0 \pm 0.2\]

(5.23 ± 0.14) / (5.23 ± 0.14)

\[1.0 \pm 0.038\]

The two numbers above, even though have the same nominal value and the same uncertainties, are actually two different measurements that only by chance share the same figures and their difference and their ratio have a non-zero uncertainty. It is common in physics to get very similar, or even equal, results for a repeated measurement, but the two measurements are not the same thing.

Instead, if you have one measurement and want to perform some operations involving it, you have to assign it to a variable:

x = 5.23 ± 0.14
x === x

true

x - x

\[0.0 \pm 0.0\]

x / x

\[1.0 \pm 0.0\]

With that in mind, let's start using Measurements.jl for realsies.

Automated UQ on an ODE: Radioactive Decay of Carbon-14

The rate of decay of carbon-14 is governed by a first order linear ordinary differential equation:

\[\frac{\mathrm{d}u(t)}{\mathrm{d}t} = -\frac{u(t)}{\tau}\]

where $\tau$ is the mean lifetime of carbon-14, which is related to the half-life $t_{1/2} = (5730 \pm 40)$ years by the relation $\tau = t_{1/2}/\ln(2)$. Writing this in DifferentialEquations.jl syntax, this looks like:

# Half-life and mean lifetime of radiocarbon, in years
t_12 = 5730 ± 40
τ = t_12 / log(2)

#Setup
u₀ = 1 ± 0
tspan = (0.0, 10000.0)

#Define the problem
radioactivedecay(u, p, t) = -u / τ

#Pass to solver
prob = ODEProblem(radioactivedecay, u₀, tspan)
sol = solve(prob, Tsit5(), reltol = 1e-8)

retcode: Success
Interpolation: specialized 4th order "free" interpolation
t: 15-element Vector{Float64}:
     0.0
     0.15287276015505386
     1.6816003617055924
    16.968876377210975
   169.84163653226483
   630.9509989706019
  1273.5406363650245
  2036.8208634612938
  2960.4190551032307
  4014.4725998532945
  5207.541840774114
  6524.063040171728
  7962.81793637761
  9515.824059252558
 10000.0
u: 15-element Vector{Measurements.Measurement{Float64}}:
        1.0 ± 0.0
 0.99998151 ± 1.3e-7
  0.9997966 ± 1.4e-6
   0.997949 ± 1.4e-5
    0.97966 ± 0.00014
    0.92652 ± 0.00049
    0.85722 ± 0.00092
     0.7816 ± 0.0013
      0.699 ± 0.0017
     0.6153 ± 0.0021
     0.5326 ± 0.0023
     0.4542 ± 0.0025
     0.3817 ± 0.0026
     0.3163 ± 0.0025
     0.2983 ± 0.0025

And bingo: numbers with uncertainty went in, so numbers with uncertainty came out. But can we trust those values for the uncertainty?

We can check the uncertainty quantification by evaluating an analytical solution to the ODE. Since it's a linear ODE, the analytical solution is simply given by the exponential:

u = exp.(-sol.t / τ)

15-element Vector{Measurements.Measurement{Float64}}:
        1.0 ± 0.0
 0.99998151 ± 1.3e-7
  0.9997966 ± 1.4e-6
   0.997949 ± 1.4e-5
    0.97966 ± 0.00014
    0.92652 ± 0.00049
    0.85722 ± 0.00092
     0.7816 ± 0.0013
      0.699 ± 0.0017
     0.6153 ± 0.0021
     0.5326 ± 0.0023
     0.4542 ± 0.0025
     0.3817 ± 0.0026
     0.3163 ± 0.0025
     0.2983 ± 0.0025

How do the two solutions compare?

plot(sol.t, sol.u, label = "Numerical", xlabel = "Years", ylabel = "Fraction of Carbon-14")
plot!(sol.t, u, label = "Analytic")

The two curves are perfectly superimposed, indicating that the numerical solution matches the analytic one. We can check that also the uncertainties are correctly propagated in the numerical solution:

println("Quantity of carbon-14 after ", sol.t[11], " years:")
println("Numerical: ", sol[11])
println("Analytic:  ", u[11])

Quantity of carbon-14 after 5207.541840774114 years:
Numerical: 0.5326 ± 0.0023
Analytic:  0.5326 ± 0.0023

Bullseye. Both the value of the numerical solution and its uncertainty match the analytic solution within the requested tolerance. We can also note that close to 5730 years after the beginning of the decay (half-life of the radioisotope), the fraction of carbon-14 that survived is about 0.5.

Simple pendulum: Small angles approximation

The next problem we are going to study is the simple pendulum in the approximation of small angles. We address this simplified case because there exists an easy analytic solution to compare.

The differential equation we want to solve is:

\[\ddot{\theta} + \frac{g}{L} \theta = 0\]

where $g = (9.79 \pm 0.02)~\mathrm{m}/\mathrm{s}^2$ is the gravitational acceleration measured where the experiment is carried out, and $L = (1.00 \pm 0.01)~\mathrm{m}$ is the length of the pendulum.

When you set up the problem for DifferentialEquations.jl remember to define the measurements as variables, as seen above.

using DifferentialEquations, Measurements, Plots

g = 9.79 ± 0.02; # Gravitational constants
L = 1.00 ± 0.01; # Length of the pendulum

#Initial Conditions
u₀ = [0 ± 0, π / 60 ± 0.01] # Initial speed and initial angle
tspan = (0.0, 6.3)

#Define the problem
function simplependulum(du, u, p, t)
    θ = u[1]
    dθ = u[2]
    du[1] = dθ
    du[2] = -(g / L) * θ
end

#Pass to solvers
prob = ODEProblem(simplependulum, u₀, tspan)
sol = solve(prob, Tsit5(), reltol = 1e-6)

retcode: Success
Interpolation: specialized 4th order "free" interpolation
t: 68-element Vector{Float64}:
 0.0
 0.07248248631909444
 0.11977474469092585
 0.178067672884797
 0.23890971597314792
 0.3037539593300621
 0.3717638701111954
 0.44269637457333977
 0.516266922273873
 0.5923146402591175
 ⋮
 5.525954381107212
 5.627441267149779
 5.72857561469214
 5.829200078041837
 5.929630055307474
 6.030620433405114
 6.133250712998771
 6.239645388649417
 6.3
u: 68-element Vector{Vector{Measurements.Measurement{Float64}}}:
 [0.0 ± 0.0, 0.052 ± 0.01]
 [0.00376 ± 0.00072, 0.051 ± 0.0097]
 [0.0061 ± 0.0012, 0.0487 ± 0.0093]
 [0.0088 ± 0.0017, 0.0444 ± 0.0085]
 [0.0114 ± 0.0022, 0.0384 ± 0.0073]
 [0.0136 ± 0.0026, 0.0304 ± 0.0058]
 [0.0154 ± 0.0029, 0.0208 ± 0.004]
 [0.0164 ± 0.0031, 0.0097 ± 0.0019]
 [0.0167 ± 0.0032, -0.00233 ± 0.00062]
 [0.0161 ± 0.0031, -0.0146 ± 0.0028]
 ⋮
 [-0.0167 ± 0.0032, 0.0006 ± 0.0046]
 [-0.0158 ± 0.0031, 0.0169 ± 0.0055]
 [-0.0134 ± 0.0027, 0.0315 ± 0.0071]
 [-0.0096 ± 0.0023, 0.0429 ± 0.0087]
 [-0.0049 ± 0.0018, 0.0501 ± 0.0097]
 [0.00033 ± 0.0016, 0.0523 ± 0.01]
 [0.0056 ± 0.0019, 0.0493 ± 0.0096]
 [0.0104 ± 0.0024, 0.0409 ± 0.0085]
 [0.0127 ± 0.0026, 0.0341 ± 0.0076]

And that's it! What about comparing it this time to the analytical solution?

u = u₀[2] .* cos.(sqrt(g / L) .* sol.t)

plot(sol.t, getindex.(sol.u, 2), label = "Numerical")
plot!(sol.t, u, label = "Analytic")

Bingo. Also in this case there is a perfect superimposition between the two curves, including their uncertainties.

We can also have a look at the difference between the two solutions:

plot(sol.t, getindex.(sol.u, 2) .- u, label = "")

Tiny difference on the order of the chosen 1e-6 tolerance.

Simple pendulum: Arbitrary amplitude

Now that we know how to solve differential equations involving numbers with uncertainties, we can solve the simple pendulum problem without any approximation. This time, the differential equation to solve is the following:

\[\ddot{\theta} + \frac{g}{L} \sin(\theta) = 0\]

That would be done via:

g = 9.79 ± 0.02; # Gravitational constants
L = 1.00 ± 0.01; # Length of the pendulum

#Initial Conditions
u₀ = [0 ± 0, π / 3 ± 0.02] # Initial speed and initial angle
tspan = (0.0, 6.3)

#Define the problem
function simplependulum(du, u, p, t)
    θ = u[1]
    dθ = u[2]
    du[1] = dθ
    du[2] = -(g / L) * sin(θ)
end

#Pass to solvers
prob = ODEProblem(simplependulum, u₀, tspan)
sol = solve(prob, Tsit5(), reltol = 1e-6)

plot(sol.t, getindex.(sol.u, 2), label = "Numerical")

Warning about Linear Uncertainty Propagation

Measurements.jl uses linear uncertainty propagation, which has an error associated with it. MonteCarloMeasurements.jl has a page which showcases where this method can lead to incorrect uncertainty measurements. Thus for more nonlinear use cases, it's suggested that one uses one of the more powerful UQ methods, such as:

• MonteCarloMeasurements.jl
• PolyChaos.jl
• SciMLExpectations.jl
• The ProbInts Uncertainty Quantification callbacks

Basically, types can make the algorithm you want to run exceedingly simple to do, but make sure it's the correct algorithm!