diff --git a/.github/workflows/scala2.yml b/.github/workflows/scala2.yml
index 752ff222..e9d61862 100644
--- a/.github/workflows/scala2.yml
+++ b/.github/workflows/scala2.yml
@@ -19,12 +19,17 @@ jobs:
 
     steps:
     - uses: actions/checkout@v3
+
     - name: Set up JDK 17
       uses: actions/setup-java@v3
       with:
         java-version: '17'
         distribution: 'temurin'
         cache: 'sbt'
+
+    - name: Install Z3
+      run: sudo apt-get update && sudo apt-get install -y z3
+
     - name: Run tests
       working-directory: ./scala2
       run: sbt test
diff --git a/scala2/src/main/resources/2023/24-test.txt b/scala2/src/main/resources/2023/24-test.txt
new file mode 100644
index 00000000..cbe1492c
--- /dev/null
+++ b/scala2/src/main/resources/2023/24-test.txt
@@ -0,0 +1,5 @@
+19, 13, 30 @ -2,  1, -2
+18, 19, 22 @ -1, -1, -2
+20, 25, 34 @ -2, -2, -4
+12, 31, 28 @ -1, -2, -1
+20, 19, 15 @  1, -5, -3
\ No newline at end of file
diff --git a/scala2/src/main/resources/2023/24.txt b/scala2/src/main/resources/2023/24.txt
new file mode 100644
index 00000000..f54d52b4
--- /dev/null
+++ b/scala2/src/main/resources/2023/24.txt
@@ -0,0 +1,300 @@
+212542581053874, 357959731032403, 176793474286781 @ -88, -256, -240
+154677220587564, 207254130208265, 139183938188421 @ 184, 74, 235
+216869547613134, 38208083662943, 397740686492049 @ 109, 262, -66
+221241619250961, 303902532813154, 249144641113790 @ 48, 24, -112
+432610900189894, 347346225570463, 389169322099761 @ -166, -99, -81
+247078054674939, 279574769079583, 357168683293046 @ 68, -13, -42
+282963504691878, 53019767043895, 238787901458543 @ -17, 399, 62
+337163551253985, 323723171285276, 476752372944125 @ -8, -74, -103
+253985064168104, 111063839515573, 450754418443501 @ -10, 456, -513
+258074760233454, 457117599855759, 244857473415713 @ 49, -247, 92
+314400910480251, 261232162431814, 339725786848352 @ -14, 8, -14
+246590774392044, 268326724153423, 252283137292781 @ -76, 195, -225
+283441229604303, 368965998122557, 461958957382550 @ -8, -135, -282
+175273393192850, 325589084134993, 358406853832354 @ 167, -63, -156
+208105215242940, 346463859222197, 284004541154973 @ 25, -119, -783
+336054085870446, 539099420004527, 243053582150753 @ -125, -451, 41
+192157026866775, 221346796869463, 239193597650948 @ 138, 27, 138
+276769382377332, 404564336355307, 139075528478897 @ -118, -279, 221
+244171538000181, 413561793671125, 231547901494709 @ 6, -271, -7
+247208942851962, 296793553833683, 207661499518285 @ -14, 32, 36
+234192138385646, 334863379919503, 180826848665889 @ -28, -63, 52
+130437043915331, 177551493108966, 166468122354373 @ 338, 448, 129
+253559901772970, 311389465584945, 216687040440642 @ -9, -19, 42
+254322299286702, 321942519715615, 221793429844513 @ -17, -42, 21
+317611519474313, 509405643844478, 492629597734339 @ -24, -308, -221
+312434433255435, 526132170689517, 256217852328836 @ -50, -379, 50
+192520813054515, 351604734054654, 148333436370478 @ 130, -152, 132
+225639342892974, 304654108851055, 198053006284713 @ 19, 39, 12
+264606145094124, 239752767682501, 207210412874243 @ 29, 61, 131
+73119136879824, 45456753334198, 74136546109856 @ 253, 194, 300
+231632956914705, 457605979777570, 302698719514094 @ 26, -414, -234
+229858006692918, 295381637696215, 265115021539193 @ -34, 118, -354
+461793485716737, 253291754164921, 311510874631043 @ -158, 5, 47
+325336710582066, 310550712766927, 330361932845723 @ 5, -61, 46
+272391892781881, 310661548818291, 260226744483606 @ -34, -25, -30
+191540619797068, 314621834402429, 114137462106199 @ 138, -22, 288
+191146615936494, 272599606676084, 172343941415066 @ 139, 136, 110
+267901500435402, 344752144584613, 304856136274139 @ -39, -98, -158
+190703217190290, 247560977233459, 149406776287313 @ 142, 550, 121
+201153500629286, 339132186878755, 150323637135481 @ 35, -57, 45
+179357508767544, 21932393507743, 145226556713621 @ 154, 315, 227
+307987323598344, 301312391795996, 300888123246103 @ -11, -40, 27
+181344541654818, 250938115658739, 239246377640253 @ 196, 440, -384
+324550882422222, 211574061061843, 237485412093545 @ -29, 72, 111
+98537457698304, 238200330125587, 135289135303526 @ 249, 31, 240
+313205485078410, 121081864504891, 172268217697865 @ -104, 348, 163
+284584577159870, 144278803049487, 200681374736457 @ -57, 323, 99
+226456734486356, 282983479210501, 218425710839487 @ -108, 324, -366
+548706725108910, 277311257109343, 307850918038817 @ -217, -28, 69
+275949926976102, 282022314903463, 311454924607825 @ -113, 87, -291
+213263551904504, 318546624732143, 175742955114241 @ 36, 19, 37
+216681369963291, 296454359998620, 462583452318421 @ 82, 10, -495
+219979404537774, 211918469842543, 232459433192801 @ 77, 188, 27
+325980549629430, 547508205102316, 314702156232182 @ -13, -324, 38
+284385758172543, 352024785807952, 218382195747631 @ -128, -120, -5
+190048554539654, 289779355315179, 173904022587641 @ 149, 388, -145
+246873768466142, 272467557005199, 232828642561081 @ -39, 131, -80
+163218080831461, 148606013314170, 142402007207516 @ 176, 164, 230
+343986154477558, 348372021296886, 250152291857989 @ -205, -106, -23
+197767130311944, 180525916592167, 161798814812969 @ 114, 519, 129
+198370970228262, 298346938466104, 174124739099332 @ 115, 54, 102
+380807631586707, 369450765586399, 481309114237823 @ -200, -141, -381
+187175525195694, 309586416720283, 109735015379081 @ 171, 173, 417
+281478516067897, 366561714660494, 230365710762598 @ -204, -184, -132
+208535431591002, 325419351576169, 191516365446671 @ 57, -12, -40
+226624736718348, 308650992940298, 199188859982218 @ -23, 62, -66
+236728026197754, 172425510194479, 206288955041960 @ 79, 131, 146
+405157823215560, 285719256118206, 491621031364803 @ -83, -34, -129
+204644528670384, 339896525956405, 151222667113763 @ -56, -54, -46
+185135616546414, 293907013443535, 166394249527601 @ 179, 231, 10
+238567782335982, 328412098910479, 192473102467361 @ -85, -26, -45
+337502877434493, 307749379575088, 369033005648444 @ -8, -58, 6
+241635494672614, 304729449451093, 326150295954981 @ 51, -27, -95
+243264345725652, 264419513819446, 345867847848017 @ 73, 6, -27
+262661477328582, 274242124962679, 187580319707945 @ -42, 80, 102
+208059705739518, 217157359130755, 219168637501631 @ 115, 85, 120
+242872575971229, 425458083316108, 227500835165921 @ 36, -258, 53
+448263894001344, 301783842143683, 426560436553661 @ -176, -43, -117
+322545182274254, 452642407811057, 305182728508413 @ -21, -227, 33
+75201284514819, 160785119325533, 104951102720081 @ 352, 241, 293
+185323568951310, 233171350573267, 165775283064269 @ 163, 358, 102
+242854931589532, 319092328660893, 243043088450517 @ 18, -38, -17
+216818686888044, 223426989024383, 346281309297331 @ 58, 283, -434
+248624452315030, 297434951348619, 284011728475200 @ 83, -49, 96
+308774084417184, 254923461188383, 324690057597731 @ -22, 27, -20
+477086131906062, 293633862617983, 289717304778929 @ -153, -43, 83
+310556413418768, 388389186178424, 201218392651440 @ -67, -172, 124
+294335453543209, 239407270896468, 327371176340196 @ 6, 40, -8
+198704436992819, 423313637384734, 211494844150179 @ 114, -360, -21
+329848820747522, 171413007628813, 298751280762187 @ -19, 102, 54
+232791155115204, 389793252905721, 200914316887487 @ -56, -314, -82
+371238777825072, 364331369420770, 349467249711398 @ -35, -114, 34
+290194660125822, 286241186809255, 265989895196345 @ -93, 39, -71
+138500514403427, 150096241065480, 90891764345181 @ 270, 386, 345
+342432728761822, 169489883251125, 160900448427412 @ -69, 145, 204
+387284555873510, 268651056280939, 472485507519965 @ -235, 48, -406
+348183642749259, 358775680962680, 425167151051150 @ -26, -110, -64
+214564960481976, 332636582662001, 166341733116147 @ -35, -19, -17
+288352094438008, 172713636729509, 508999993126113 @ 32, 94, -171
+429522310493674, 350751066474876, 353772711184036 @ -241, -106, -110
+235317931818252, 321016228988431, 221782559196323 @ 8, -29, -24
+237395224369976, 271952629687205, 200198427527953 @ 48, 46, 109
+267689676493904, 197692176620362, 371002531720142 @ -151, 456, -664
+245022629101794, 336310573633898, 178670106406896 @ 79, -86, 192
+308555808260267, 235053234862080, 386617281992034 @ 20, 16, -14
+254939001555789, 196490967246088, 344406853138846 @ 46, 120, -66
+198067603029891, 228271655034477, 230280443116238 @ 112, 353, -142
+286182290367015, 324356938662982, 337470234197366 @ 40, -74, 30
+182772323746872, 243078357133489, 180935879286305 @ 182, 418, -11
+334555295661618, 283598866844113, 307163703216413 @ -19, -28, 51
+268032515539496, 179884088413109, 201707373095034 @ 53, 89, 166
+314152706857177, 231973181098005, 294537386515547 @ -10, 41, 47
+211145445957746, 256130578705417, 227544572919273 @ 105, 54, 81
+240281216930434, 308871905093706, 201208104074959 @ -81, 58, -69
+194990331403464, 353680311245143, 148152875882771 @ 96, -217, 54
+376506453102756, 282681817701967, 308014018133222 @ -59, -29, 56
+270559608594507, 293243401201593, 219465840278274 @ -38, 17, 48
+271434227107862, 301909819058233, 255037189953162 @ -157, 57, -212
+228819748864119, 35415178477933, 360982336746131 @ 100, 225, 7
+282571416535269, 237239909717928, 361818410574456 @ 34, 28, -20
+237878607070278, 316580360861695, 326197735422041 @ 42, -39, -160
+182576050978588, 187101572838317, 131079628110881 @ 153, 161, 245
+235076492312502, 318319071558547, 260428708551349 @ -203, 96, -762
+187583408959494, 313306547152243, 102716500573961 @ 189, 318, 643
+239870832918417, 26077551674877, 194534422319897 @ 10, 745, 77
+268344200395759, 336594882856173, 201457637971271 @ -144, -71, -13
+351502642713711, 343426968745552, 327993505552334 @ -54, -94, 8
+203247301060818, 304211148313729, 144668102389835 @ 17, 294, 108
+205287991254828, 326526363823258, 176717703001868 @ 73, -18, 32
+228237293076774, 492607291603871, 215976066898633 @ 49, -457, 39
+210744469193866, 242345551455857, 338364614486931 @ 87, 173, -305
+283795432748766, 324397233772201, 303141716476529 @ 27, -71, 37
+318328668468378, 290176478743420, 393177776909696 @ -9, -32, -60
+296886504660264, 330367685997713, 290768439652901 @ -8, -76, 23
+354925965270690, 204310714489876, 324988197369764 @ -73, 86, -6
+277940901110162, 351813200734015, 35068250706381 @ 26, -105, 370
+302119534267998, 319648005021051, 257122942759009 @ -185, -26, -123
+220643077695684, 355855985807773, 179518771733771 @ 30, -142, 66
+281128245154754, 261417029969318, 390265842706956 @ 47, -10, -20
+401214520647879, 209382315271573, 394433050114910 @ -66, 37, -12
+202419387388134, 370308338321798, 371095720268716 @ 115, -152, -266
+243933688132567, 284238541251243, 160523673406178 @ -60, 127, 134
+207949728449758, 342146025132292, 160034991638202 @ 27, -90, 52
+328107083354006, 506783656792851, 68374835799249 @ -27, -292, 321
+334371695727666, 295677547881937, 372669058678174 @ -209, 21, -342
+425555650358416, 141386636896383, 400700320063993 @ -94, 107, -23
+354111348691934, 360815768065903, 286452836048521 @ -183, -129, -62
+186640407106107, 302586436160656, 178053441373097 @ 172, 200, -99
+184504662208422, 279753007845271, 194263341780489 @ 217, 645, -497
+207498586652070, 122966697911623, 453309639270761 @ 88, 592, -760
+175420575327994, 240376844547933, 227794777853156 @ 159, 37, 122
+164379037426222, 314155556336019, 361113505933531 @ 171, -59, -30
+276158080997198, 255307550628353, 190284041278157 @ 27, 22, 167
+188538116291142, 332645165411655, 164893512402481 @ 166, 10, -105
+206485169181270, 324400373771635, 174538862304413 @ 37, 28, -44
+231705775027484, 354858179856273, 282940665637611 @ 53, -119, -77
+192935594972436, 232318329359044, 140280091724297 @ 125, 770, 173
+195865425610134, 156203126394403, 201075304936541 @ 121, 618, -22
+220291544685654, 60861797261303, 192839119984101 @ 97, 313, 156
+216421547979532, 166625355222409, 138232910764571 @ 33, 645, 215
+454745381031134, 359071031321598, 433358646779191 @ -229, -116, -177
+241211456940334, 274238344824863, 449761596808401 @ 87, -22, -86
+300489403237401, 150600505604740, 311021695401350 @ 12, 130, 36
+351714451571994, 313870670017933, 87623195530631 @ -79, -54, 304
+206304315814166, 322825195619263, 165678725657741 @ 47, 27, 35
+235260754952454, 192089045978143, 324835846533921 @ 88, 81, 21
+234012378392974, 220537666254223, 244419610199201 @ 21, 243, -67
+271646520250074, 328249590902743, 281718062915501 @ 38, -75, 56
+221872342890713, 311155829759331, 169665424751063 @ 96, -49, 191
+234965964033108, 289605268944487, 264575781614519 @ 53, 11, -17
+302656676541801, 314956862712124, 279759708917957 @ 22, -64, 89
+194778772231287, 358987377936095, 94944329588325 @ 100, -269, 633
+213962227655009, 360175350319088, 187183591353131 @ 78, -140, 95
+335872446070469, 310194803249308, 258956779512841 @ -196, -18, -51
+237632929340327, 282862400087895, 329308829827744 @ -54, 155, -578
+195165381881310, 334127994548035, 177008381765921 @ 83, 25, -395
+322431513433582, 335929297302159, 292621591828001 @ -117, -80, -70
+315117522492294, 188045711663619, 217859262699941 @ -11, 94, 140
+373704262775700, 207296292228736, 462531809606177 @ -143, 116, -267
+251431566932470, 262653021312867, 228387402272864 @ -45, 151, -52
+200968557727895, 335950123937744, 148782648944484 @ -28, 20, -56
+236191762393467, 21120063276068, 319730241477691 @ 22, 742, -245
+208913018909190, 375167847286459, 97238860543841 @ -29, -401, 565
+290435031159135, 56359235789263, 318028085602160 @ 28, 227, 36
+205348838719584, 129089359330730, 397134601580767 @ 109, 358, -317
+194757851501280, 310732265284483, 257897724121307 @ 122, 57, -352
+322105170096600, 281140697612407, 384217392018125 @ -40, -9, -101
+227263802927379, 293407939553692, 217416903679568 @ 34, 50, -6
+253024402494847, 472881498291876, 346563578286452 @ -60, -512, -448
+314479634908016, 119470812770471, 423291884678263 @ -102, 343, -326
+285064202741996, 260225765977046, 328922413430668 @ 17, 14, -12
+164764803937924, 131918495543174, 24421731030948 @ 209, 466, 528
+302401846001018, 47953136574059, 415330628001011 @ -43, 389, -220
+136889723461149, 102922183481659, 257049978466682 @ 254, 415, -22
+266977346640134, 399124471391351, 389593482782381 @ 29, -175, -130
+372790937915094, 329723676236983, 408552056960081 @ -242, -66, -337
+295046549007926, 307216389914204, 245279985581636 @ -93, -14, -10
+222962551504923, 300940192862044, 205863476663477 @ 42, 34, 17
+217936499850689, 313784346557728, 244304797861441 @ 86, -36, 21
+167754243520824, 183008145578501, 255319117162995 @ 184, 214, 6
+210177003182478, 429234450438847, 105538845228113 @ 101, -266, 296
+199347157290654, 270021317518867, 163266752925959 @ 99, 261, 88
+236595610312209, 357567541944595, 163695965251100 @ -116, -177, 62
+312426175948804, 126535043607303, 334910821408881 @ 20, 119, 45
+250185142335527, 339957515144567, 186490088641817 @ -40, -85, 77
+187486534782198, 153491915531223, 178050669591013 @ 145, 217, 168
+264452793637374, 364051575147343, 278407165734113 @ -66, -153, -167
+321652447174633, 234494711812364, 293511778578951 @ -94, 100, -45
+101178528277354, 45310253629503, 279722555547721 @ 231, 211, 93
+311483135370989, 117588178673003, 195876215498686 @ -30, 223, 154
+349074097657518, 71980788262470, 67744127629012 @ -53, 236, 322
+205665951950534, 317578175678991, 174079200152461 @ 9, 131, -140
+172224556223070, 207776433045037, 341587542422723 @ 163, 78, -22
+199777030731058, 159408594476631, 309534334024285 @ 108, 565, -396
+464121799742544, 439378219125328, 480156984371951 @ -191, -210, -177
+512533623338691, 269596699742212, 497024716185359 @ -200, -16, -141
+449327916519144, 205855641898298, 303200495165981 @ -131, 50, 65
+283724606928492, 240851979789205, 384065029237529 @ 38, 18, -31
+288013100270877, 208780552824757, 275880042837536 @ 42, 41, 100
+225018718395018, 305644724002975, 271186961007503 @ 73, -21, -28
+252580597893096, 289594241717095, 283259230473967 @ -14, 39, -134
+217696966439229, 236860502115098, 202812154030551 @ 82, 134, 91
+290383585895985, 172960262418583, 413523311841740 @ 22, 107, -88
+246394854383374, 404386901502983, 316452007110281 @ 63, -178, -10
+277440283560394, 298664786803743, 343675845620671 @ 29, -37, -26
+549218613325353, 365039367598781, 370814441177496 @ -212, -115, 10
+314422297180566, 290227116094987, 220524013374437 @ -127, 20, 52
+208288985823519, 332582248866208, 173680963176656 @ 38, -34, -6
+344954667419544, 271844404820980, 373070962444213 @ -11, -24, 9
+253707352980234, 336179622652543, 293095895839541 @ 22, -81, -58
+348619876072741, 233964822099034, 330658038016178 @ -48, 36, 8
+294065185957476, 303581225750131, 130406957718587 @ -14, -35, 246
+299497691666274, 143845614560983, 193269298796021 @ -30, 217, 148
+243761048715324, 305091461548017, 163187922986167 @ -56, 46, 126
+228852534457750, 321833988811099, 195275147618589 @ -19, -6, -24
+210326083482005, 405825122389483, 255008495328029 @ 48, -393, -343
+234154015900814, 320051906908983, 254552486072961 @ 15, -28, -111
+221069337881082, 268544927933243, 253793013055353 @ 40, 152, -161
+350843574391570, 546291208573433, 245148884150221 @ -57, -343, 105
+321884637993916, 151517460880797, 338566691457073 @ -108, 268, -147
+332310394683086, 296925474202815, 290407863681483 @ -133, -5, -62
+238868612579282, 442767280701488, 312738360614324 @ -73, -538, -562
+238005301025171, 315819717023939, 310603161892176 @ 62, -49, -50
+184500511341324, 354706787987715, 106858269141817 @ 184, -175, 409
+218719926572754, 248846852339899, 417842058727985 @ 109, 9, -67
+61527743340558, 41910468004471, 21054306023633 @ 311, 306, 391
+206891968805558, 359036109292941, 96347232076181 @ 71, -164, 395
+358569854122062, 275252683087759, 280731684589601 @ -219, 51, -75
+281134641567294, 255607683726743, 288258712879345 @ -11, 52, -17
+325364111053745, 294010589908703, 211245733791773 @ -82, -13, 113
+268578371026014, 236127445255093, 380797038274583 @ 59, 17, -13
+489485534335218, 226913670756739, 271116263277629 @ -155, 21, 107
+200091726654254, 345151854422543, 355346332534721 @ 125, -97, -106
+288773407987668, 371116294945099, 330355289788371 @ -128, -171, -300
+398345872685398, 211733420030191, 256489704564165 @ -127, 75, 84
+410444070061341, 232446702321753, 270268423194688 @ -80, 17, 106
+223446110503560, 342994867695369, 176936935212271 @ -23, -95, 15
+221456980555537, 346760891427937, 198873039052710 @ 8, -111, -48
+241138246161499, 378379591190555, 197910544306414 @ 44, -161, 118
+239309488603942, 346974551472951, 243277229211186 @ 51, -101, 40
+266720686719534, 249073195909171, 326627536261997 @ 59, 6, 38
+175265275057344, 233389584704983, 166232708483831 @ 228, 519, 48
+199647143220834, 333605166183208, 126338949433076 @ 87, -38, 274
+355059474706240, 112909300036023, 456817283755217 @ -18, 127, -67
+291234253229982, 225827198327447, 300458706607553 @ -17, 89, -19
+289865647604208, 451420237900663, 518182602681917 @ 12, -233, -253
+351557916167787, 226921027459900, 441676807044779 @ -112, 88, -241
+171724416708978, 120198249743035, 113845845697733 @ 210, 720, 308
+324992776605402, 546442272740313, 176464245808959 @ -35, -358, 186
+266166769986634, 268730110669199, 328440341073557 @ 9, 35, -97
+225525857274734, 233012276737543, 309959118198911 @ 75, 111, -88
+183405703797942, 309355721084695, 182533926443609 @ 173, 53, 19
+210515462842953, 291925001904295, 215164290510485 @ -72, 459, -671
+235537482505454, 310690173156743, 156326933471977 @ -21, 22, 154
+273084618253552, 382626888261135, 273689280998585 @ 8, -157, 17
+253722385917350, 315123819243295, 240205666004325 @ 57, -57, 102
+196398380835834, 355900579581511, 153312099518087 @ 109, -169, 118
+283444009699329, 366561045831077, 202002888581270 @ -146, -167, 26
+296532850164962, 253904723261403, 296984689916093 @ 38, -8, 86
+356774602967614, 373148843974943, 387561802202761 @ -67, -131, -74
+252770877244094, 318171858594683, 216752381636081 @ -25, -28, 17
+338919821368038, 232645806842808, 377368303830121 @ -29, 32, -35
+209566211606045, 338493165391956, 150634080385665 @ -72, -46, 21
+328388255345246, 295645021337031, 342359520358565 @ -13, -41, 11
+357923314873354, 213990132914548, 252051459029871 @ -145, 126, 39
+300794330404169, 291383403542233, 331991009111531 @ -28, -15, -61
+171452764273723, 106269888549931, 325980849738994 @ 168, 255, -42
+400651551538638, 227035846405175, 440348818761665 @ -113, 46, -127
+183528864416638, 337358904333282, 149647897440530 @ 155, -82, 206
+407834029671234, 363038317346083, 326507144649005 @ -126, -118, 6
+261232554625244, 272762277507558, 260158766990881 @ -139, 184, -267
+276291421501401, 231191690007043, 332403327250334 @ 32, 47, -8
+213146770430035, 320595954009642, 182342124018161 @ 36, 10, 5
+431072263340031, 549716908674601, 541461834661703 @ -234, -415, -393
+271067118142410, 252860807780751, 293724860670289 @ 25, 35, 13
+230424456187626, 346336855193767, 374228162998841 @ 76, -99, -145
+233683858918340, 451654359126852, 407119183631371 @ 45, -334, -365
+309530132675051, 360560922327322, 427729374335125 @ -118, -132, -399
+171253360815090, 226755112694191, 106314555408725 @ 286, 763, 428
+195185795409569, 365215513552403, 270027401984661 @ 134, -121, 73
\ No newline at end of file
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent03.scala b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent03.scala
index a08dd397..9904f70b 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent03.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent03.scala
@@ -11,7 +11,7 @@ import org.scalatest.matchers.should.Matchers._
 object Advent03 {
   type ClaimId = Int
 
-  final case class Claim(id: ClaimId, area: Area2D)
+  final case class Claim(id: ClaimId, area: Area2D[Int])
   object Claim {
     private val RegEx           = """#(\d+) @ (\d+),(\d+): (\d+)x(\d+)""".r
     def parse(s: String): Claim =
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent06.scala b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent06.scala
index eee8e2b0..d0965e49 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent06.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent06.scala
@@ -1,6 +1,7 @@
 package jurisk.adventofcode.y2018
 
 import cats.implicits._
+import jurisk.geometry.Coordinates2D
 import jurisk.geometry.Coords2D
 import jurisk.utils.CollectionOps.IterableOps
 import jurisk.utils.FileInput.parseFileLines
@@ -12,7 +13,7 @@ object Advent06 {
     parseFileLines(fileName, _.parseCoords2D)
 
   def part1(points: List[Coords2D]): Int = {
-    val boundingBox = Coords2D.boundingBoxInclusive(points)
+    val boundingBox = Coordinates2D.boundingBoxInclusive(points)
     println(s"Bounding box: $boundingBox")
 
     def onEdgeOfBoundingBox(point: Coords2D): Boolean =
@@ -47,7 +48,7 @@ object Advent06 {
   }
 
   def part2(points: List[Coords2D], limit: Int): Int =
-    Coords2D.boundingBoxInclusive(points).points.count { point =>
+    Coordinates2D.boundingBoxInclusive(points).points.count { point =>
       val distanceToAll = points.map(_.manhattanDistance(point)).sum
       distanceToAll < limit
     }
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent10.scala b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent10.scala
index 5431fab4..2d4a6e4c 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent10.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent10.scala
@@ -1,6 +1,7 @@
 package jurisk.adventofcode.y2018
 
 import cats.implicits._
+import jurisk.geometry.Coordinates2D
 import jurisk.geometry.Coords2D
 import jurisk.geometry.SparseBooleanField
 import jurisk.utils.FileInput._
@@ -34,7 +35,7 @@ object Advent10 {
 
     Simulation.runWithIterationCount(()) { case ((), time) =>
       val atTime      = field.map(x => x.atTime(time.toInt))
-      val boundingBox = Coords2D.boundingBoxInclusive(atTime)
+      val boundingBox = Coordinates2D.boundingBoxInclusive(atTime)
       if (boundingBox.height <= limit) {
         val fieldAtTime = SparseBooleanField(atTime, 100).toDebugRepresentation
         println(s"At time $time:")
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent17.scala b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent17.scala
index 6fd1348a..068216a5 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent17.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent17.scala
@@ -5,6 +5,7 @@ import jurisk.adventofcode.y2018.Advent17.Square.Clay
 import jurisk.adventofcode.y2018.Advent17.Square._
 import jurisk.adventofcode.y2018.Advent17.State.consecutiveRanges
 import jurisk.adventofcode.y2018.Advent17.State.isSolid
+import jurisk.geometry.Coordinates2D
 import jurisk.geometry.Coords2D
 import jurisk.geometry.Field2D
 import jurisk.utils.FileInput._
@@ -137,15 +138,15 @@ object Advent17 {
     }
   }
 
-  private def parsePoints(s: String): List[Coords2D] =
+  private def parsePoints(s: String): Seq[Coords2D] =
     s match {
       case s"x=$x, y=$y1..$y2" =>
-        Coords2D.allPointsInclusive(
+        Coordinates2D.allPointsInclusive(
           Coords2D.of(x.toInt, y1.toInt),
           Coords2D.of(x.toInt, y2.toInt),
         )
       case s"y=$y, x=$x1..$x2" =>
-        Coords2D.allPointsInclusive(
+        Coordinates2D.allPointsInclusive(
           Coords2D.of(x1.toInt, y.toInt),
           Coords2D.of(x2.toInt, y.toInt),
         )
@@ -157,7 +158,7 @@ object Advent17 {
   def parse(data: String): Field2D[Square] = {
     val points = data.parseLines(parsePoints).flatten
 
-    val boundingBox = Coords2D
+    val boundingBox = Coordinates2D
       .boundingBoxInclusive(SpringOfWater :: points)
       .expandInEachDirectionBy(1)
 
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent23.scala b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent23.scala
index b09a05c0..604c3f83 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent23.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2018/Advent23.scala
@@ -11,11 +11,13 @@ import jurisk.optimization.Optimizer
 import jurisk.utils.FileInput._
 import jurisk.utils.Parsing.StringOps
 
+import scala.annotation.nowarn
+
 object Advent23 {
   type Input = List[Nanobot]
 
   final case class Nanobot(
-    position: Coords3D,
+    position: Coords3D[Int],
     radius: Int,
   )
 
@@ -69,7 +71,7 @@ object Advent23 {
         nanobot.position.y,
         nanobot.position.z,
         nanobot.radius,
-      ).map(constant)
+      ).map(intConstant)
 
       val inRange = (x - nx).abs + (y - ny).abs + (z - nz).abs <= nr
 
@@ -83,19 +85,14 @@ object Advent23 {
       distanceFromOrigin === sum(abs(x), abs(y), abs(z)),
     )
 
-    // Objective - maximize nanobotsInRange and minimize distanceFromOrigin
-    val objective1 = maximize(nanobotsInRange)
-    val objective2 = minimize(distanceFromOrigin)
-
-    val model = checkAndGetModel()
+    // Objectives - maximize nanobotsInRange and minimize distanceFromOrigin
+    val _ = maximize(nanobotsInRange)
+    val _ = minimize(distanceFromOrigin)
 
-    println(model)
-
-    val List(xc, yc, zc, nir, dor) =
-      List(x, y, z, nanobotsInRange, distanceFromOrigin).map(extractInt)
+    val List(xc, yc, zc) = runExternal("x", "y", "z").map(resultToInt)
 
     val found = Coords3D(xc, yc, zc)
-    println(s"$found: nanobots in range: $nir, distance from origin: $dor")
+
     found.manhattanDistance(Coords3D.Zero)
   }
 
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent18.scala b/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent18.scala
index 63078044..78e55562 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent18.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent18.scala
@@ -7,22 +7,22 @@ import jurisk.utils.Parsing.StringOps
 import org.scalatest.matchers.should.Matchers._
 
 object Advent18 {
-  private def sidesFree(cube: Coords3D, cubes: Set[Coords3D]): Int =
+  private def sidesFree(cube: Coords3D[Int], cubes: Set[Coords3D[Int]]): Int =
     cube.adjacent6.count(n => !cubes.contains(n))
 
-  def parse(data: String): Set[Coords3D] =
-    data.parseLines(Coords3D.parse).toSet
+  def parse(data: String): Set[Coords3D[Int]] =
+    data.parseLines(Coords3D.parse[Int]).toSet
 
-  def part1(points: Set[Coords3D]): Int =
+  def part1(points: Set[Coords3D[Int]]): Int =
     points.toList.map(sidesFree(_, points)).sum
 
-  def part2(points: Set[Coords3D]): Int = {
+  def part2(points: Set[Coords3D[Int]]): Int = {
     val boundingBox         = Coords3D.boundingBoxInclusive(points)
     val expandedBoundingBox = boundingBox.expandInEachDirectionBy(1)
 
     val startingPointForSearch = expandedBoundingBox.min
     val reachable              = Bfs
-      .bfsReachable[Coords3D](
+      .bfsReachable[Coords3D[Int]](
         startingPointForSearch,
         _.adjacent6.filter { n =>
           !points.contains(n) && expandedBoundingBox.contains(n)
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent21.scala b/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent21.scala
index 2bf1920f..4537cb77 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent21.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent21.scala
@@ -184,10 +184,8 @@ object Advent21 {
       }
     }
 
-    val model = checkAndGetModel()
-    println(model)
-
-    extractLong(calculate.labeledInt)
+    val List(result) = runExternal(calculate)
+    resultToLong(result)
   }
 
   private def solvePart1Optimizer(
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent23.scala b/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent23.scala
index 8c6ab7b6..6ab6b36c 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent23.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent23.scala
@@ -1,5 +1,6 @@
 package jurisk.adventofcode.y2022
 
+import jurisk.geometry.Coordinates2D
 import jurisk.geometry.Coords2D
 import jurisk.geometry.Direction2D
 import jurisk.geometry.Direction2D._
@@ -81,7 +82,7 @@ object Advent23 {
     }
 
     def debugPrint(): Unit = {
-      val box       = Coords2D.boundingBoxInclusive(elves.map(_.location))
+      val box       = Coordinates2D.boundingBoxInclusive(elves.map(_.location))
       val field     = Field2D.forArea(box, '.')
       val resulting = elves.foldLeft(field) { case (acc, e) =>
         acc.updatedAtUnsafe(e.location, e.charRep)
@@ -107,7 +108,7 @@ object Advent23 {
         state.next(iteration.toInt)
     }
 
-    val box = Coords2D.boundingBoxInclusive(result.allCoords)
+    val box = Coordinates2D.boundingBoxInclusive(result.allCoords)
     (box.width * box.height) - result.allCoords.size
   }
 
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent24.scala b/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent24.scala
index 8857df79..4a077b93 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent24.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2022/Advent24.scala
@@ -24,7 +24,7 @@ object Advent24 {
     start: Coords2D,
     blizzards: Set[Blizzard],
     goal: Coords2D,
-    area: Area2D,
+    area: Area2D[Int],
   ) {
     def freeAtTimeNear(time: Int, location: Coords2D): List[Coords2D] = {
       val blizzardCoordinatesInNext = blizzardsAt(time)
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent11.scala b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent11.scala
index 86ff87ba..56a8588e 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent11.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent11.scala
@@ -1,7 +1,7 @@
 package jurisk.adventofcode.y2023
 
 import cats.implicits._
-import jurisk.geometry.BigIntCoords2D
+import jurisk.geometry.Coordinates2D
 import jurisk.geometry.Field2D
 import jurisk.utils.FileInput._
 import jurisk.utils.Parsing.StringOps
@@ -12,7 +12,8 @@ import scala.collection.immutable.ArraySeq
 import scala.collection.mutable
 
 object Advent11 {
-  private type Galaxies = ArraySeq[BigIntCoords2D]
+  private type BigIntCoords2D = Coordinates2D[BigInt]
+  private type Galaxies       = ArraySeq[BigIntCoords2D]
 
   def parse(input: String): Galaxies = {
     val field = Field2D.parseBooleanField(input)
@@ -20,7 +21,7 @@ object Advent11 {
     val coords = field.filterCoordsByValue(_ == true)
 
     val results = coords.map { c =>
-      BigIntCoords2D(c.x, c.y)
+      Coordinates2D[BigInt](c.x, c.y)
     }
 
     ArraySeq.from(results)
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent16.scala b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent16.scala
index f63e3f6a..2969d610 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent16.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent16.scala
@@ -20,6 +20,7 @@ import jurisk.utils.Parsing.StringOps
 import jurisk.utils.Simulation
 import mouse.all.booleanSyntaxMouse
 
+import scala.annotation.nowarn
 import scala.collection.immutable.ArraySeq
 
 object Advent16 {
@@ -250,9 +251,9 @@ object Advent16 {
       println(s"Upper:\n${objective.getUpper}")
     }
 
-    val model = checkAndGetModel()
-
     if (debug) {
+      @nowarn("cat=deprecation")
+      val model = checkAndGetModel()
       println(s"Model:\n$model")
 
       val debugField = field
@@ -261,7 +262,10 @@ object Advent16 {
             val incoming = incomingBool(c, d)
             val outgoing = outgoingBool(c, d)
 
+            @nowarn("cat=deprecation")
             val incm = extractBoolean(incoming).getOrElse("Unknown".fail)
+
+            @nowarn("cat=deprecation")
             val outg = extractBoolean(outgoing).getOrElse("Unknown".fail)
 
             (incm, outg) match {
@@ -288,7 +292,8 @@ object Advent16 {
       Field2D.printCharField(debugField)
     }
 
-    extractLong(energizedVar)
+    val List(result) = runExternal("energized")
+    resultToLong(result)
   }
 
   private def edgeIncomings(
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent21.scala b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent21.scala
index d302ffae..316532c8 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent21.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent21.scala
@@ -89,7 +89,10 @@ object Advent21 {
 
   // Working but slow
   def part2Simulation(data: Input, steps: Int): Long = {
-    def debugPrint(area: Area2D, map: mutable.HashMap[Coords2D, Long]): Long = {
+    def debugPrint(
+      area: Area2D[Int],
+      map: mutable.HashMap[Coords2D, Long],
+    ): Long = {
       var results = 0
 
       val stringField = Field2D.forArea(area, ' ').mapByCoords { c =>
@@ -411,6 +414,8 @@ object Advent21 {
     val realData: Input = parseFile(fileName(""))
 
     println(s"Part 1: ${part1(realData, 64)}")
+
+    // TODO:  Try solving Part 2 as extrapolating as a 2nd degree polynomial for 65 + n * 131 steps, where n = 0..4
     println(s"Part 2: ${part2FieldClassification(realData, 26501365)}")
   }
 }
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent22.scala b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent22.scala
index 1228a9ff..5acf52e3 100644
--- a/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent22.scala
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent22.scala
@@ -13,7 +13,7 @@ object Advent22 {
   type Input   = Vector[Brick]
   type BrickId = Int
 
-  final case class Brick(id: BrickId, blocks: Area3D) {
+  final case class Brick(id: BrickId, blocks: Area3D[Int]) {
     def name: String = if (id < 26) {
       ('A'.toInt + id).toChar.toString
     } else {
diff --git a/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent24.scala b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent24.scala
new file mode 100644
index 00000000..794867c8
--- /dev/null
+++ b/scala2/src/main/scala/jurisk/adventofcode/y2023/Advent24.scala
@@ -0,0 +1,414 @@
+package jurisk.adventofcode.y2023
+
+import cats.implicits._
+import jurisk.geometry.Area2D
+import jurisk.geometry.Coordinates2D
+import jurisk.geometry.Coords3D
+import jurisk.geometry.Coords3D.Axis
+import jurisk.geometry.lineLineIntersectionGivenTwoPointsOnEachLine
+import jurisk.math.GaussianElimination
+import jurisk.math.positiveAndNegativeDivisors
+import jurisk.optimization.ImplicitConversions.RichArithExprIntSort
+import jurisk.optimization.ImplicitConversions.RichExpr
+import jurisk.optimization.Optimizer
+import jurisk.utils.CollectionOps.IterableOps
+import jurisk.utils.FileInput._
+import jurisk.utils.Parsing.StringOps
+import org.scalatest.matchers.should.Matchers._
+
+object Advent24 {
+  private type InputPart2 = List[PositionAndVelocity3D]
+
+  final case class PositionAndVelocity2D(
+    position: Coordinates2D[Long],
+    velocity: Coordinates2D[Long],
+  )
+
+  final case class PositionAndVelocity3D(
+    position: Coords3D[Long],
+    velocity: Coords3D[Long],
+  ) {
+    def p: Coords3D[Long] = position
+    def v: Coords3D[Long] = velocity
+  }
+
+  def parse(input: String): InputPart2 = {
+    def parse3D(input: String): PositionAndVelocity3D =
+      input match {
+        case s"$position @ $velocity" =>
+          PositionAndVelocity3D(
+            Coords3D.parse[Long](position),
+            Coords3D.parse[Long](velocity),
+          )
+        case _                        => input.failedToParse
+      }
+
+    input.parseLines(parse3D)
+  }
+
+  // https://en.wikipedia.org/wiki/Chinese_remainder_theorem
+  // From https://github.com/ellentari/aoc2023/blob/main/src/main/scala/aoc/Day24.scala
+  // and https://github.com/ellentari/aoc2023/blob/main/src/main/scala/aoc/util/ChineseRemainderTheorem.scala
+  def solve2UsingChineseRemainderTheorem(
+    data: List[PositionAndVelocity3D]
+  ): Long = {
+    println("For getting to Chinese Remainder Theorem:")
+    println("v = vx + vy + vz")
+    println("p = px + py + pz")
+    println()
+
+    data.zipWithIndex foreach { case (r, id) =>
+      val idx = id + 1
+      val rp  = r.p.x + r.p.y + r.p.z
+      val rv  = r.v.x + r.v.y + r.v.z
+      println(s"p ${sgn(-rp)} = (${sgn(rv)} -v) * t$idx")
+      println(" ==> ")
+      // https://en.wikipedia.org/wiki/Modular_arithmetic#Congruence
+      println(s"p ≡ $rp (mod ($rv - v))")
+      println()
+    }
+
+    // We can assume that rv is rather small, e.g. -1000 to +1000.
+    // Then we can solve CRT for various assumed `coprime` (requirement for CRT!) `rv` values.
+    // Solution of the CRT is the `p` which is the solution we are looking for.
+
+    "Not implemented".fail
+  }
+
+  def solve2InferringVelocity(
+    data: List[PositionAndVelocity3D]
+  ): PositionAndVelocity3D = {
+    val debug = false
+
+    if (debug) printEquations(data)
+
+    val velocity = inferVelocity(data)
+
+    def solveAssumingV(
+      data: List[PositionAndVelocity3D],
+      v: Coords3D[Long],
+    ): Coords3D[Long] = {
+      if (debug) {
+        println(
+          s"Now that we know (vx, vy, vz) == $v it becomes a much simpler task"
+        )
+        data.zipWithIndex foreach { case (r, id) =>
+          val idx = id + 1
+          println(s"px = ${r.p.x} ${sgn(r.v.x - v.x)} * t$idx")
+          println(s"py = ${r.p.y} ${sgn(r.v.y - v.y)} * t$idx")
+          println(s"pz = ${r.p.z} ${sgn(r.v.z - v.z)} * t$idx")
+          println()
+        }
+      }
+
+      // Gaussian reduction (just the first 5 equations are enough, the rest are redundant)
+
+      val h = data.head
+      val g = data(1)
+
+      // Columns: px, py, pz, t1, t2
+      // Rows: Equations, 3 from h, 2 from g
+      val A = Array(
+        Array(1, 0, 0, v.x - h.v.x, 0),
+        Array(0, 1, 0, v.y - h.v.y, 0),
+        Array(0, 0, 1, v.z - h.v.z, 0),
+        Array(1, 0, 0, 0, v.x - g.v.x),
+        Array(0, 1, 0, 0, v.y - g.v.y),
+      ).map(_.map(_.toDouble))
+
+      val b = Array(h.p.x, h.p.y, h.p.z, g.p.x, g.p.y).map(_.toDouble)
+
+      val Array(px, py, pz, t1 @ _, t2 @ _) = GaussianElimination.solve(A, b)
+
+      Coords3D[Long](px.toLong, py.toLong, pz.toLong)
+    }
+
+    val position = solveAssumingV(data, velocity)
+
+    PositionAndVelocity3D(
+      position,
+      velocity,
+    )
+  }
+
+  // https://math.stackexchange.com/a/697278
+  def linesIntersect(
+    a: PositionAndVelocity3D,
+    b: PositionAndVelocity3D,
+  ): Boolean = {
+    val cp    = a.velocity crossProduct b.velocity
+    val pDiff = a.position - b.position
+
+    (cp dotProduct pDiff) == 0
+  }
+
+  private def vectorIntersection2D(
+    a: PositionAndVelocity2D,
+    b: PositionAndVelocity2D,
+  ): Option[Coordinates2D[BigDecimal]] = {
+
+    val result = lineLineIntersectionGivenTwoPointsOnEachLine(
+      a.position.map(BigDecimal(_)),
+      (a.position + a.velocity).map(BigDecimal(_)),
+      b.position.map(BigDecimal(_)),
+      (b.position + b.velocity).map(BigDecimal(_)),
+    )
+
+    result.flatMap { result =>
+      val px = result.x
+      val py = result.y
+
+      val Eps = 0.001
+      val tax = (px - a.position.x) / a.velocity.x
+      val tay = (py - a.position.y) / a.velocity.y
+      assert((tax - tay).abs <= Eps)
+      val ta  = (tax + tay) / 2
+
+      val tbx = (px - b.position.x) / b.velocity.x
+      val tby = (py - b.position.y) / b.velocity.y
+      assert((tbx - tby).abs <= Eps)
+      val tb  = (tbx + tby) / 2
+
+      if (ta > 0 && tb > 0) {
+        Coordinates2D[BigDecimal](px, py).some
+      } else {
+        none
+      }
+    }
+  }
+
+  def part1(data: InputPart2, minC: Long, maxC: Long): Long = {
+    val area =
+      Area2D[BigDecimal](
+        Coordinates2D[BigDecimal](minC, minC),
+        Coordinates2D[BigDecimal](maxC, maxC),
+      )
+
+    val input = data.map { c =>
+      PositionAndVelocity2D(
+        Coordinates2D[Long](
+          c.position.x,
+          c.position.y,
+        ),
+        Coordinates2D[Long](
+          c.velocity.x,
+          c.velocity.y,
+        ),
+      )
+    }
+
+    def solve1(
+      input: List[PositionAndVelocity2D],
+      area: Area2D[BigDecimal],
+    ): Long = {
+      def intersectWithinBounds2D(
+        a: PositionAndVelocity2D,
+        b: PositionAndVelocity2D,
+      ): Boolean = {
+        val debug = false
+
+        if (debug) println(s"a = $a, b= $b")
+        val result = vectorIntersection2D(a, b)
+        if (debug) println(s"result = $result\n")
+
+        result match {
+          case Some(c) =>
+            area.contains(c)
+          case None    => false
+        }
+      }
+
+      input.combinations(2).count { list =>
+        val List(a, b) = list
+        intersectWithinBounds2D(a, b)
+      }
+    }
+
+    solve1(input, area)
+  }
+
+  private def sgn(n: Long): String =
+    if (n < 0) {
+      s"- ${-n}"
+    } else {
+      s"+ $n"
+    }
+
+  private def toBasicEquation(
+    r: PositionAndVelocity3D,
+    idx: String,
+    axis: Axis,
+  ): String = {
+    val rp = r.p(axis)
+    val rv = r.v(axis)
+    val a  = axis.toChar
+
+    s"p$a ${sgn(-rp)} = t$idx * ($rv - v$a)"
+  }
+
+  private def inferVelocity(
+    data: List[PositionAndVelocity3D]
+  ): Coords3D[Long] = {
+    def deriveV(axis: Axis): Long = {
+      val debug                         = false
+      var candidates: Option[Set[Long]] = None
+
+      if (debug) println(s"Same r.v.${axis.toChar}: ")
+      data
+        .groupBy(_.v(axis))
+        .filter { case (_, equations) => equations.size >= 2 }
+        .foreach { case (n, list) =>
+          if (debug) {
+            list.zipWithIndex.foreach { case (r, idx) =>
+              println(toBasicEquation(r, ('a' + idx).toString, axis))
+            }
+          }
+
+          list.combinations(2) foreach { list2 =>
+            val List(a, b) = list2
+            val rpDiff     = (a.p(axis) - b.p(axis)).abs
+
+            assert(n == a.v(axis))
+            assert(n == b.v(axis))
+
+            val divisors    = positiveAndNegativeDivisors(rpDiff)
+            if (debug)
+              println(s"($n - v${axis.toChar}) is one of $divisors, thus...")
+            val validValues = divisors.map(n - _)
+            if (debug) println(s"v${axis.toChar} is one of $validValues")
+
+            candidates match {
+              case Some(filtered) =>
+                val n = filtered intersect validValues.toSet
+                candidates = n.some
+              case None           => candidates = validValues.toSet.some
+            }
+          }
+          if (debug) println()
+        }
+
+      if (debug) {
+        println(s"Outcome: Valid v${axis.toChar}-es: $candidates")
+        println()
+      }
+
+      candidates.get.toSeq.singleResultUnsafe
+    }
+
+    val vx = deriveV(Axis.X)
+    val vy = deriveV(Axis.Y)
+    val vz = deriveV(Axis.Z)
+
+    println(s"v = $vx, $vy, $vz")
+
+    Coords3D[Long](vx, vy, vz)
+  }
+
+  def printEquations(data: List[PositionAndVelocity3D]): Unit = {
+    println("Basic equations:")
+    data.zipWithIndex foreach { case (r, id) =>
+      val idx = id + 1
+      println(s"px + t$idx * vx = ${r.p.x} ${sgn(r.v.x)} * t$idx")
+      println(s"py + t$idx * vy = ${r.p.y} ${sgn(r.v.y)} * t$idx")
+      println(s"pz + t$idx * vz = ${r.p.z} ${sgn(r.v.z)} * t$idx")
+      println()
+    }
+
+    data.zipWithIndex foreach { case (r, id) =>
+      val idx = id + 1
+      Axis.All foreach { axes =>
+        println(toBasicEquation(r, idx.toString, axes))
+      }
+
+      println()
+    }
+
+    println()
+
+    println(s"${3 + 3 + data.length} variables")
+    println(s"${data.length * 3} equations")
+
+    println()
+    println("px + t_n * vx = rpx + rvx * t_n")
+    println("py + t_n * vy = rpy + rvy * t_n")
+    println("pz + t_n * vz = rpz + rvz * t_n")
+
+    println()
+  }
+
+  def solvePart2Optimizer(
+    data: List[PositionAndVelocity3D]
+  ): PositionAndVelocity3D = {
+    // Find a "result" PositionAndVelocity3D for which integer t exists where "position at t" for both
+    // "result" and all of "data" is identical
+
+    // t[n] - collision time with rock n, > 0
+
+    // find (px, py, pz) and (vx, vy, vz) and such t[] so that for all rocks:
+    //    px + t[n] * vx == px[n] + t[n] * vx[n]
+    //    py + t[n] * vy == py[n] + t[n] * vy[n]
+    //    pz + t[n] * vz == pz[n] + t[n] * vz[n]
+
+    implicit val o: Optimizer = Optimizer.z3()
+    import o._
+
+    val px = o.labeledInt("px")
+    val py = o.labeledInt("py")
+    val pz = o.labeledInt("pz")
+
+    val vx = o.labeledInt("vx")
+    val vy = o.labeledInt("vy")
+    val vz = o.labeledInt("vz")
+
+    data.zipWithIndex.foreach { case (rock, idx) =>
+      val t_n  = o.labeledInt(s"t_$idx")
+      val px_n = o.longConstant(rock.position.x)
+      val py_n = o.longConstant(rock.position.y)
+      val pz_n = o.longConstant(rock.position.z)
+
+      val vx_n = o.longConstant(rock.velocity.x)
+      val vy_n = o.longConstant(rock.velocity.y)
+      val vz_n = o.longConstant(rock.velocity.z)
+
+      o.addConstraints(
+        t_n >= Zero,
+        px + t_n * vx === px_n + t_n * vx_n,
+        py + t_n * vy === py_n + t_n * vy_n,
+        pz + t_n * vz === pz_n + t_n * vz_n,
+      )
+    }
+
+    val List(pxS, pyS, pzS, vxS, vyS, vzS) = o
+      .runExternal("px", "py", "pz", "vx", "vy", "vz")
+      .map(r => resultToLong(r))
+
+    val result = PositionAndVelocity3D(
+      Coords3D[Long](pxS, pyS, pzS),
+      Coords3D[Long](vxS, vyS, vzS),
+    )
+
+    println(result)
+
+    result
+  }
+
+  def part2(data: InputPart2): Long = {
+    // TODO:  Consider also trying Newton-Raphson (see https://github.com/rzikm/advent-of-code/blob/master/2023/24.fs,
+    //        or https://pastebin.com/s6nvy0jA) and/or gradient descent
+    val result = solve2InferringVelocity(data)
+    result.position.x + result.position.y + result.position.z
+  }
+
+  def parseFile(fileName: String): InputPart2 =
+    parse(readFileText(fileName))
+
+  def fileName(suffix: String): String =
+    s"2023/24$suffix.txt"
+
+  def main(args: Array[String]): Unit = {
+    val realData: InputPart2 = parseFile(fileName(""))
+
+    println(s"Part 1: ${part1(realData, 200000000000000L, 400000000000000L)}")
+    println(s"Part 2: ${part2(realData)}")
+  }
+}
diff --git a/scala2/src/main/scala/jurisk/geometry/Area2D.scala b/scala2/src/main/scala/jurisk/geometry/Area2D.scala
index 0431e7dc..924baa8d 100644
--- a/scala2/src/main/scala/jurisk/geometry/Area2D.scala
+++ b/scala2/src/main/scala/jurisk/geometry/Area2D.scala
@@ -1,55 +1,71 @@
 package jurisk.geometry
 
-final case class Area2D(min: Coords2D, max: Coords2D) {
-  def onInsideEdge(c: Coords2D): Boolean =
+import jurisk.math.Enumerated
+import jurisk.math.Enumerated.EnumeratedOps
+
+import scala.math.Numeric.Implicits.infixNumericOps
+import scala.math.Ordering.Implicits.infixOrderingOps
+
+// Inclusive
+final case class Area2D[N: Numeric](
+  min: Coordinates2D[N],
+  max: Coordinates2D[N],
+) {
+  private val One = implicitly[Numeric[N]].one
+
+  def onInsideEdge(c: Coordinates2D[N]): Boolean =
     (c.x == min.x) || (c.x == max.x) || (c.y == min.y) || (c.y == max.y)
 
-  def points: List[Coords2D] = {
+  def points(implicit enumerable: Enumerated[N]): Seq[Coordinates2D[N]] =
     (min.x to max.x) flatMap { x =>
       (min.y to max.y) map { y =>
-        Coords2D.of(x, y)
+        Coordinates2D.of[N](x, y)
       }
     }
-  }.toList
 
-  def shiftBy(c: Coords2D): Area2D =
-    Area2D(min + c, max + c)
+  def shiftBy(c: Coordinates2D[N]): Area2D[N] =
+    Area2D[N](min + c, max + c)
 
-  def +(c: Coords2D): Area2D =
+  def +(c: Coordinates2D[N]): Area2D[N] =
     shiftBy(c)
 
-  def topLeft: Coords2D     = min
-  def bottomRight: Coords2D = max
+  def topLeft: Coordinates2D[N]     = min
+  def bottomRight: Coordinates2D[N] = max
 
-  def height: Int = (max.y - min.y) + 1
-  def width: Int  = (max.x - min.x) + 1
+  def height: N = (max.y - min.y) + One
+  def width: N  = (max.x - min.x) + One
 
-  def left: Int   = min.x
-  def right: Int  = max.x
-  def top: Int    = min.y
-  def bottom: Int = max.y
+  def left: N   = min.x
+  def right: N  = max.x
+  def top: N    = min.y
+  def bottom: N = max.y
 
-  def expandInEachDirectionBy(n: Int): Area2D = Area2D(
-    min - Coords2D.of(n, n),
-    max + Coords2D.of(n, n),
+  def expandInEachDirectionBy(n: N): Area2D[N] = Area2D[N](
+    min - Coordinates2D.of[N](n, n),
+    max + Coordinates2D.of[N](n, n),
   )
 
-  def contains(c: Coords2D): Boolean =
+  def contains(c: Coordinates2D[N]): Boolean =
     c.x >= min.x && c.x <= max.x && c.y >= min.y && c.y <= max.y
 }
 
 object Area2D {
-  def fromLeftTopWidthHeight(
-    left: Int,
-    top: Int,
-    width: Int,
-    height: Int,
-  ): Area2D =
+  def fromLeftTopWidthHeight[N: Numeric](
+    left: N,
+    top: N,
+    width: N,
+    height: N,
+  ): Area2D[N] = {
+    val One = implicitly[Numeric[N]].one
+
     Area2D(
-      Coords2D.of(left, top),
-      Coords2D.of(left + width - 1, top + height - 1),
+      Coordinates2D.of(left, top),
+      Coordinates2D.of(left + width - One, top + height - One),
     )
+  }
 
-  def boundingBoxInclusive(coords: Seq[Coords2D]): Area2D =
-    Coords2D.boundingBoxInclusive(coords)
+  def boundingBoxInclusive[N: Numeric](
+    coords: Seq[Coordinates2D[N]]
+  ): Area2D[N] =
+    Coordinates2D.boundingBoxInclusive(coords)
 }
diff --git a/scala2/src/main/scala/jurisk/geometry/Area3D.scala b/scala2/src/main/scala/jurisk/geometry/Area3D.scala
index a97e958c..ac75de17 100644
--- a/scala2/src/main/scala/jurisk/geometry/Area3D.scala
+++ b/scala2/src/main/scala/jurisk/geometry/Area3D.scala
@@ -1,10 +1,15 @@
 package jurisk.geometry
 
-final case class Area3D(min: Coords3D, max: Coords3D) {
-  def contains(n: Coords3D): Boolean =
+import jurisk.math.Enumerated
+import jurisk.math.Enumerated.EnumeratedOps
+
+import scala.math.Ordering.Implicits.infixOrderingOps
+
+final case class Area3D[N: Numeric](min: Coords3D[N], max: Coords3D[N]) {
+  def contains(n: Coords3D[N]): Boolean =
     n.x >= min.x && n.x <= max.x && n.y >= min.y && n.y <= max.y && n.z >= min.z && n.z <= max.z
 
-  def points: Seq[Coords3D] =
+  def points(implicit enumerated: Enumerated[N]): Seq[Coords3D[N]] =
     (min.x to max.x) flatMap { x =>
       (min.y to max.y) flatMap { y =>
         (min.z to max.z) map { z =>
@@ -13,12 +18,12 @@ final case class Area3D(min: Coords3D, max: Coords3D) {
       }
     }
 
-  def expandInEachDirectionBy(n: Int): Area3D = Area3D(
+  def expandInEachDirectionBy(n: N): Area3D[N] = Area3D(
     min - Coords3D(n, n, n),
     max + Coords3D(n, n, n),
   )
 
-  def moveBy(c: Coords3D): Area3D =
+  def moveBy(c: Coords3D[N]): Area3D[N] =
     Area3D(
       min = min + c,
       max = max + c,
@@ -26,6 +31,8 @@ final case class Area3D(min: Coords3D, max: Coords3D) {
 }
 
 object Area3D {
-  def boundingBoxInclusive(points: Iterable[Coords3D]): Area3D =
+  def boundingBoxInclusive[N: Numeric](
+    points: Iterable[Coords3D[N]]
+  ): Area3D[N] =
     Coords3D.boundingBoxInclusive(points)
 }
diff --git a/scala2/src/main/scala/jurisk/geometry/BigIntCoords2D.scala b/scala2/src/main/scala/jurisk/geometry/BigIntCoords2D.scala
deleted file mode 100644
index 608baa37..00000000
--- a/scala2/src/main/scala/jurisk/geometry/BigIntCoords2D.scala
+++ /dev/null
@@ -1,6 +0,0 @@
-package jurisk.geometry
-
-object BigIntCoords2D {
-  def apply(x: BigInt, y: BigInt): BigIntCoords2D =
-    Coordinates2D.of[BigInt](x, y)
-}
diff --git a/scala2/src/main/scala/jurisk/geometry/Coordinates2D.scala b/scala2/src/main/scala/jurisk/geometry/Coordinates2D.scala
index 54b7d815..9495e6d1 100644
--- a/scala2/src/main/scala/jurisk/geometry/Coordinates2D.scala
+++ b/scala2/src/main/scala/jurisk/geometry/Coordinates2D.scala
@@ -1,12 +1,13 @@
 package jurisk.geometry
 
+import cats.Functor
 import cats.implicits._
 import jurisk.math.Enumerated
 import jurisk.utils.Parsing.StringOps
 
-import scala.math.Integral.Implicits.infixIntegralOps
+import scala.math.Numeric.Implicits.infixNumericOps
 
-final case class Coordinates2D[N: Integral](x: N, y: N) {
+final case class Coordinates2D[N: Numeric](x: N, y: N) {
   def +(other: Coordinates2D[N]): Coordinates2D[N] =
     Coordinates2D(x + other.x, y + other.y)
 
@@ -16,16 +17,25 @@ final case class Coordinates2D[N: Integral](x: N, y: N) {
   def *(n: N): Coordinates2D[N] =
     Coordinates2D(x * n, y * n)
 
-  def manhattanDistanceToOrigin: N =
+  def manhattanDistanceToOrigin(implicit integral: Integral[N]): N =
     x.abs + y.abs
 
-  def manhattanDistance(other: Coordinates2D[N]): N =
+  def manhattanDistance(other: Coordinates2D[N])(implicit
+    integral: Integral[N]
+  ): N =
     (this - other).manhattanDistanceToOrigin
 
-  def adjacent4: List[Coordinates2D[N]] = neighbours(includeDiagonal = false)
-  def adjacent8: List[Coordinates2D[N]] = neighbours(includeDiagonal = true)
+  def adjacent4(implicit integral: Integral[N]): List[Coordinates2D[N]] =
+    neighbours(includeDiagonal = false)
+  def adjacent8(implicit integral: Integral[N]): List[Coordinates2D[N]] =
+    neighbours(includeDiagonal = true)
+
+  def neighbours(
+    includeDiagonal: Boolean
+  )(implicit integral: Integral[N]): List[Coordinates2D[N]] = {
+    val _ =
+      integral.zero // This is a hack to ensure that the compiler does not complain about the unused Integral requirement
 
-  def neighbours(includeDiagonal: Boolean): List[Coordinates2D[N]] = {
     val directions = if (includeDiagonal) {
       Direction2D.AllDirections
     } else {
@@ -53,18 +63,21 @@ final case class Coordinates2D[N: Integral](x: N, y: N) {
   def NW: Coordinates2D[N] = this + Direction2D.NW
   def SE: Coordinates2D[N] = this + Direction2D.SE
   def SW: Coordinates2D[N] = this + Direction2D.SW
+
+  def map[M: Numeric](f: N => M): Coordinates2D[M] =
+    Coordinates2D(f(x), f(y))
 }
 
 object Coordinates2D {
   implicit def readingOrdering[N: Ordering]: Ordering[Coordinates2D[N]] =
     Ordering[(N, N)].contramap(c => (c.y, c.x))
 
-  def zero[N: Integral]: Coordinates2D[N] = {
+  def zero[N: Numeric]: Coordinates2D[N] = {
     val numeric = implicitly[Numeric[N]]
     Coordinates2D.of[N](numeric.zero, numeric.zero)
   }
 
-  def of[N: Integral](x: N, y: N): Coordinates2D[N] =
+  def of[N: Numeric](x: N, y: N): Coordinates2D[N] =
     Coordinates2D[N](x, y)
 
   def parse[N: Integral](s: String): Coordinates2D[N] = {
@@ -75,10 +88,22 @@ object Coordinates2D {
     Coordinates2D.of[N](an, bn)
   }
 
+  def boundingBoxInclusive[N: Numeric](
+    coords: Iterable[Coordinates2D[N]]
+  ): Area2D[N] = {
+    val xList = coords.map(_.x)
+    val minX  = xList.min
+    val maxX  = xList.max
+    val yList = coords.map(_.y)
+    val minY  = yList.min
+    val maxY  = yList.max
+    Area2D(Coordinates2D.of(minX, minY), Coordinates2D.of(maxX, maxY))
+  }
+
   def allPointsInclusive[N: Integral: Enumerated](
     a: Coordinates2D[N],
     b: Coordinates2D[N],
-  ): List[Coordinates2D[N]] = {
+  ): Seq[Coordinates2D[N]] = {
     import jurisk.math.Enumerated.EnumeratedOps
     val numeric          = implicitly[Numeric[N]]
     val min: (N, N) => N = numeric.min
@@ -102,5 +127,4 @@ object Coordinates2D {
     } else
       s"Expected $a and $b to have same x or y coordinates, but they do not".fail
   }
-
 }
diff --git a/scala2/src/main/scala/jurisk/geometry/Coords2D.scala b/scala2/src/main/scala/jurisk/geometry/Coords2D.scala
index 03556820..1a0ceec0 100644
--- a/scala2/src/main/scala/jurisk/geometry/Coords2D.scala
+++ b/scala2/src/main/scala/jurisk/geometry/Coords2D.scala
@@ -1,6 +1,7 @@
 package jurisk.geometry
 
 import cats.implicits._
+import jurisk.math.Matrix2x2
 
 object Coords2D {
   def apply(x: Int, y: Int): Coords2D = Coordinates2D[Int](x, y)
@@ -8,23 +9,12 @@ object Coords2D {
 
   val Zero: Coords2D = Coordinates2D.zero[Int]
 
-  // Note - This could be moved to `Coordinates2D` if `Area2D` is made generic
-  def boundingBoxInclusive(coords: Iterable[Coords2D]): Area2D = {
-    val xList = coords.map(_.x)
-    val minX  = xList.min
-    val maxX  = xList.max
-    val yList = coords.map(_.y)
-    val minY  = yList.min
-    val maxY  = yList.max
-    Area2D(Coordinates2D.of(minX, minY), Coordinates2D.of(maxX, maxY))
-  }
-
   def parse(s: String): Coords2D = Coordinates2D.parse[Int](s)
 
   implicit val readingOrdering: Ordering[Coords2D] =
     Coordinates2D.readingOrdering[Int]
 
-  def allPointsInclusive(a: Coords2D, b: Coords2D): List[Coords2D] =
+  def allPointsInclusive(a: Coords2D, b: Coords2D): Seq[Coords2D] =
     Coordinates2D.allPointsInclusive(a, b)
 
   private def adjacentCircularIndices[T](
@@ -47,10 +37,15 @@ object Coords2D {
   def areaOfSimplePolygon(seq: IndexedSeq[Coords2D]): Double = {
     val determinants = adjacentCircularValues(seq)
       .map { case (a, b) =>
-        a.x.toDouble * b.y.toDouble - a.y.toDouble * b.x.toDouble
+        Matrix2x2(
+          a.x.toDouble,
+          b.x.toDouble,
+          a.y.toDouble,
+          b.y.toDouble,
+        ).determinant
       }
 
-    (determinants.sum / 2.0).abs
+    determinants.sum / 2.0
   }
 
   def calculateBoundaryPoints(seq: IndexedSeq[Coords2D]): Long =
diff --git a/scala2/src/main/scala/jurisk/geometry/Coords3D.scala b/scala2/src/main/scala/jurisk/geometry/Coords3D.scala
index 360ddc1e..58f9aa64 100644
--- a/scala2/src/main/scala/jurisk/geometry/Coords3D.scala
+++ b/scala2/src/main/scala/jurisk/geometry/Coords3D.scala
@@ -1,44 +1,103 @@
 package jurisk.geometry
 
-import cats.implicits._
+import jurisk.geometry.Coords3D.Axis
 import jurisk.utils.Parsing.StringOps
 
-final case class Coords3D(x: Int, y: Int, z: Int) {
-  def adjacent6: List[Coords3D] = List(
-    Coords3D(x + 1, y, z),
-    Coords3D(x - 1, y, z),
-    Coords3D(x, y + 1, z),
-    Coords3D(x, y - 1, z),
-    Coords3D(x, y, z + 1),
-    Coords3D(x, y, z - 1),
+import scala.math.Numeric.Implicits.infixNumericOps
+
+final case class Coords3D[N: Numeric](x: N, y: N, z: N) {
+  private val One = implicitly[Numeric[N]].one
+
+  def adjacent6: List[Coords3D[N]] = List(
+    Coords3D(x + One, y, z),
+    Coords3D(x - One, y, z),
+    Coords3D(x, y + One, z),
+    Coords3D(x, y - One, z),
+    Coords3D(x, y, z + One),
+    Coords3D(x, y, z - One),
   )
 
-  def manhattanDistance(other: Coords3D): Int =
-    Math.abs(x - other.x) + Math.abs(y - other.y) + Math.abs(z - other.z)
+  def manhattanDistance(other: Coords3D[N]): N =
+    (x - other.x).abs + (y - other.y).abs + (z - other.z).abs
 
-  def +(other: Coords3D): Coords3D = Coords3D(
+  def +(other: Coords3D[N]): Coords3D[N] = Coords3D(
     x + other.x,
     y + other.y,
     z + other.z,
   )
 
-  def -(other: Coords3D): Coords3D = Coords3D(
+  def -(other: Coords3D[N]): Coords3D[N] = Coords3D(
     x - other.x,
     y - other.y,
     z - other.z,
   )
+
+  def apply(axis: Axis): N = axis match {
+    case Axis.X => x
+    case Axis.Y => y
+    case Axis.Z => z
+  }
+
+  // https://en.wikipedia.org/wiki/Cross_product
+  def crossProduct(b: Coords3D[N]): Coords3D[N] = {
+    val a = this
+    Coords3D[N](
+      a.y * b.z - a.z * b.y,
+      a.z * b.x - a.x * b.z,
+      a.x * b.y - a.y * b.x,
+    )
+  }
+
+  // https://en.wikipedia.org/wiki/Dot_product
+  def dotProduct(b: Coords3D[N]): N = {
+    val a = this
+    a.x * b.x + a.y * b.y + a.z * b.z
+  }
 }
 
 object Coords3D {
-  val Zero: Coords3D = Coords3D(0, 0, 0)
+  sealed trait Axis {
+    def toChar: Char
+  }
 
-  def parse(s: String): Coords3D =
+  object Axis {
+    case object X extends Axis {
+      override def toChar: Char = 'x'
+    }
+    case object Y extends Axis {
+      override def toChar: Char = 'y'
+    }
+    case object Z extends Axis {
+      override def toChar: Char = 'z'
+    }
+
+    val All: Set[Axis] = Set(X, Y, Z)
+  }
+
+  def Zero[N: Numeric]: Coords3D[N] = {
+    val Z = implicitly[Numeric[N]].zero
+    Coords3D(Z, Z, Z)
+  }
+
+  def parse[N: Numeric](s: String): Coords3D[N] = {
+    val numeric = implicitly[Numeric[N]]
     s match {
-      case s"$x,$y,$z" => Coords3D(x.trim.toInt, y.trim.toInt, z.trim.toInt)
-      case _           => s.failedToParse("Coords3D")
+      case s"$xs,$ys,$zs" =>
+        val x: N =
+          numeric.parseString(xs.trim) getOrElse s"Failed to parse $xs".fail
+        val y: N =
+          numeric.parseString(ys.trim) getOrElse s"Failed to parse $ys".fail
+        val z: N =
+          numeric.parseString(zs.trim) getOrElse s"Failed to parse $zs".fail
+
+        Coords3D[N](x, y, z)
+      case _              => s.failedToParse("Coords3D")
     }
+  }
 
-  def boundingBoxInclusive(points: Iterable[Coords3D]): Area3D = {
+  def boundingBoxInclusive[N: Numeric](
+    points: Iterable[Coords3D[N]]
+  ): Area3D[N] = {
     val minX = points.map(_.x).min
     val minY = points.map(_.y).min
     val minZ = points.map(_.z).min
@@ -47,4 +106,16 @@ object Coords3D {
     val maxZ = points.map(_.z).max
     Area3D(Coords3D(minX, minY, minZ), Coords3D(maxX, maxY, maxZ))
   }
+
+  def areVectorsParallel(a: Coords3D[Long], b: Coords3D[Long]): Boolean = {
+    val ax = BigDecimal(a.x)
+    val ay = BigDecimal(a.y)
+    val az = BigDecimal(a.z)
+
+    val bx = BigDecimal(b.x)
+    val by = BigDecimal(b.y)
+    val bz = BigDecimal(b.z)
+
+    List(ax / bx, ay / by, az / bz).distinct.size == 1
+  }
 }
diff --git a/scala2/src/main/scala/jurisk/geometry/Field2D.scala b/scala2/src/main/scala/jurisk/geometry/Field2D.scala
index 1e5291c0..43db6fed 100644
--- a/scala2/src/main/scala/jurisk/geometry/Field2D.scala
+++ b/scala2/src/main/scala/jurisk/geometry/Field2D.scala
@@ -336,7 +336,10 @@ object Field2D {
     }
   }
 
-  def forArea[T: ClassTag](boundingBox: Area2D, initialValue: T): Field2D[T] =
+  def forArea[T: ClassTag](
+    boundingBox: Area2D[Int],
+    initialValue: T,
+  ): Field2D[T] =
     ofSize(
       boundingBox.width,
       boundingBox.height,
diff --git a/scala2/src/main/scala/jurisk/geometry/SparseField.scala b/scala2/src/main/scala/jurisk/geometry/SparseField.scala
index 0b616f73..e7cf2529 100644
--- a/scala2/src/main/scala/jurisk/geometry/SparseField.scala
+++ b/scala2/src/main/scala/jurisk/geometry/SparseField.scala
@@ -7,7 +7,7 @@ final case class SparseField[T](points: Map[Coords2D, T]) {
     display: Option[T] => Char,
     limit: Int,
   ): String = {
-    val boundingBox = Coords2D.boundingBoxInclusive(points.keys)
+    val boundingBox = Coordinates2D.boundingBoxInclusive(points.keys)
     val xMin        = boundingBox.topLeft.x
     val yMin        = boundingBox.topLeft.y
     val xMax        = boundingBox.bottomRight.x
diff --git a/scala2/src/main/scala/jurisk/geometry/package.scala b/scala2/src/main/scala/jurisk/geometry/package.scala
index e633fca3..ba1e443a 100644
--- a/scala2/src/main/scala/jurisk/geometry/package.scala
+++ b/scala2/src/main/scala/jurisk/geometry/package.scala
@@ -1,5 +1,10 @@
 package jurisk
 
+import cats.implicits._
+import jurisk.math.Matrix2x2
+
+import scala.math.Fractional.Implicits.infixFractionalOps
+
 package object geometry {
   private[geometry] def mergeSeqSeqChar(data: Seq[Seq[Char]]): String =
     data.map(_.mkString).map(_ + '\n').mkString
@@ -9,4 +14,59 @@ package object geometry {
   type BigIntCoords2D = Coordinates2D[BigInt]
 
   def visualizeBoolean(b: Boolean): Char = if (b) '█' else '░'
+
+  // https://en.wikipedia.org/wiki/Line%E2%80%93line_intersection#Given_two_points_on_each_line
+  def lineLineIntersectionGivenTwoPointsOnEachLine[N: Fractional](
+    a1: Coordinates2D[N],
+    a2: Coordinates2D[N],
+    b1: Coordinates2D[N],
+    b2: Coordinates2D[N],
+  ): Option[Coordinates2D[N]] = {
+    val Zero = implicitly[Fractional[N]].zero
+    val One  = implicitly[Fractional[N]].one
+
+    val x1 = a1.x
+    val y1 = a1.y
+    val x2 = a2.x
+    val y2 = a2.y
+    val x3 = b1.x
+    val y3 = b1.y
+    val x4 = b2.x
+    val y4 = b2.y
+
+    val x1_x2 = Matrix2x2(x1, One, x2, One).determinant
+    val x3_x4 = Matrix2x2(x3, One, x4, One).determinant
+
+    val y1_y2 = Matrix2x2(y1, One, y2, One).determinant
+    val y3_y4 = Matrix2x2(y3, One, y4, One).determinant
+
+    val bottom = Matrix2x2(x1_x2, y1_y2, x3_x4, y3_y4).determinant
+
+    if (bottom == Zero) {
+      none
+    } else {
+      val a = Matrix2x2(x1, y1, x2, y2).determinant
+      val b = Matrix2x2(x3, y3, x4, y4).determinant
+
+      val pxTop = Matrix2x2(
+        a,
+        x1_x2,
+        b,
+        x3_x4,
+      ).determinant
+
+      val px = pxTop / bottom
+
+      val pyTop = Matrix2x2(
+        a,
+        y1_y2,
+        b,
+        y3_y4,
+      ).determinant
+
+      val py = pyTop / bottom
+
+      Coordinates2D(px, py).some
+    }
+  }
 }
diff --git a/scala2/src/main/scala/jurisk/math/DiscreteInterval.scala b/scala2/src/main/scala/jurisk/math/DiscreteInterval.scala
index 66527092..b70ff5f2 100644
--- a/scala2/src/main/scala/jurisk/math/DiscreteInterval.scala
+++ b/scala2/src/main/scala/jurisk/math/DiscreteInterval.scala
@@ -21,7 +21,9 @@ final case class DiscreteInterval[N: Numeric: Enumerated](
 
   def size: N = to - from + One
 
-  def toList: List[N] = from to to
+  def toSeq: Seq[N] = from to to
+
+  def toList: List[N] = toSeq.toList
 
   def union(
     other: DiscreteInterval[N]
diff --git a/scala2/src/main/scala/jurisk/math/Enumerated.scala b/scala2/src/main/scala/jurisk/math/Enumerated.scala
index 9472990c..fd50e769 100644
--- a/scala2/src/main/scala/jurisk/math/Enumerated.scala
+++ b/scala2/src/main/scala/jurisk/math/Enumerated.scala
@@ -1,23 +1,23 @@
 package jurisk.math
 
 trait Enumerated[N] {
-  def to(fromInclusive: N, toInclusive: N): List[N]
+  def to(fromInclusive: N, toInclusive: N): Seq[N]
 }
 
 object Enumerated {
   implicit class EnumeratedOps[N](val fromInclusive: N)(implicit
     ev: Enumerated[N]
   ) {
-    def to(toInclusive: N): List[N] = ev.to(fromInclusive, toInclusive)
+    def to(toInclusive: N): Seq[N] = ev.to(fromInclusive, toInclusive)
   }
 
   implicit object IntEnumerated extends Enumerated[Int] {
-    def to(fromInclusive: Int, toInclusive: Int): List[Int] =
-      (fromInclusive to toInclusive).toList
+    def to(fromInclusive: Int, toInclusive: Int): Seq[Int] =
+      fromInclusive to toInclusive
   }
 
   implicit object LongEnumerated extends Enumerated[Long] {
-    def to(fromInclusive: Long, toInclusive: Long): List[Long] =
-      (fromInclusive to toInclusive).toList
+    def to(fromInclusive: Long, toInclusive: Long): Seq[Long] =
+      fromInclusive to toInclusive
   }
 }
diff --git a/scala2/src/main/scala/jurisk/math/GaussianElimination.scala b/scala2/src/main/scala/jurisk/math/GaussianElimination.scala
new file mode 100644
index 00000000..af289fa8
--- /dev/null
+++ b/scala2/src/main/scala/jurisk/math/GaussianElimination.scala
@@ -0,0 +1,39 @@
+package jurisk.math
+
+// https://en.wikipedia.org/wiki/Gaussian_elimination
+object GaussianElimination {
+  def solve(A: Array[Array[Double]], b: Array[Double]): Array[Double] = {
+    val n         = A.length
+    val augmented = A.zip(b).map { case (row, bi) => row :+ bi }
+
+    for (col <- 0 until n) {
+      // Find pivot row
+      val pivotRow = (col until n).maxBy(row => math.abs(augmented(row)(col)))
+      val temp     = augmented(col)
+      augmented(col) = augmented(pivotRow)
+      augmented(pivotRow) = temp
+
+      // Make leading coefficient of pivot row 1
+      val pivotElement = augmented(col)(col)
+      for (j <- col until n + 1)
+        augmented(col)(j) /= pivotElement
+
+      // Eliminate below pivot
+      for (i <- col + 1 until n) {
+        val factor = augmented(i)(col)
+        for (j <- col until n + 1)
+          augmented(i)(j) -= factor * augmented(col)(j)
+      }
+    }
+
+    // Back substitution
+    val x = new Array[Double](n)
+    for (i <- n - 1 to 0 by -1) {
+      x(i) = augmented(i)(n)
+      for (j <- i + 1 until n)
+        x(i) -= augmented(i)(j) * x(j)
+    }
+
+    x
+  }
+}
diff --git a/scala2/src/main/scala/jurisk/math/Matrix2x2.scala b/scala2/src/main/scala/jurisk/math/Matrix2x2.scala
new file mode 100644
index 00000000..e9d425da
--- /dev/null
+++ b/scala2/src/main/scala/jurisk/math/Matrix2x2.scala
@@ -0,0 +1,13 @@
+package jurisk.math
+
+import scala.math.Numeric.Implicits.infixNumericOps
+
+final case class Matrix2x2[N: Numeric](a: N, b: N, c: N, d: N) {
+  override def toString: String =
+    s"""$a $b
+       |$c $d
+       |
+       |""".stripMargin
+
+  def determinant: N = a * d - b * c
+}
diff --git a/scala2/src/main/scala/jurisk/math/package.scala b/scala2/src/main/scala/jurisk/math/package.scala
index ce98db86..5f37c831 100644
--- a/scala2/src/main/scala/jurisk/math/package.scala
+++ b/scala2/src/main/scala/jurisk/math/package.scala
@@ -51,6 +51,25 @@ package object math {
   def pow(a: Int, b: Int): Long =
     Math.pow(a, b).toLong
 
+  def positiveDivisors(n: Long): Seq[Long] =
+    (1L to Math.sqrt(n.toDouble).toLong)
+      .flatMap(i =>
+        if (n % i == 0) {
+          if (i * i == n) {
+            // Avoid duplicates when n is a perfect square
+            Seq(i)
+          } else {
+            Seq(i, n / i)
+          }
+        } else Seq.empty[Long]
+      )
+      .sorted
+
+  def positiveAndNegativeDivisors(n: Long): Seq[Long] = {
+    val positive = positiveDivisors(n.abs)
+    positive.map(-_).reverse ++ positive
+  }
+
   def countLongsBetweenExclusive(low: Double, high: Double): Long = {
     def floorExclusive(x: Double): Long = {
       val result = x.floor
diff --git a/scala2/src/main/scala/jurisk/optimization/Optimizer.scala b/scala2/src/main/scala/jurisk/optimization/Optimizer.scala
index ab97d74a..8ea3dee7 100644
--- a/scala2/src/main/scala/jurisk/optimization/Optimizer.scala
+++ b/scala2/src/main/scala/jurisk/optimization/Optimizer.scala
@@ -12,11 +12,21 @@ import com.microsoft.z3.IntNum
 import com.microsoft.z3.IntSort
 import com.microsoft.z3.Model
 import com.microsoft.z3.Optimize
+import com.microsoft.z3.RatNum
+import com.microsoft.z3.RealExpr
 import com.microsoft.z3.Sort
 import com.microsoft.z3.Status
 import com.microsoft.z3.enumerations.Z3_lbool
+import jurisk.process.Runner
 import jurisk.utils.Parsing.StringOps
-
+import org.scalatest.matchers.should.Matchers.convertToAnyShouldWrapper
+
+// The https://github.com/tudo-aqua/z3-turnkey distribution did not work on task 2023-24 while the same SMT-LIB program
+// worked from the command line. Thus, some methods have ended up being deprecated, and this class is mostly
+// a way to generate SMT-LIB programs (see https://smtlib.cs.uiowa.edu/language.shtml).
+//
+// You could consider removing the Z3 dependency and just generating SMT-LIB format directly, but it's probably not
+// worth it.
 trait Optimizer {
   val context: Context
   val optimize: Optimize
@@ -33,17 +43,27 @@ trait Optimizer {
 
   def addConstraints(expressions: Expr[BoolSort]*): Unit
 
+  @deprecated("Use `runExternal` instead", "2023-12-24")
   def checkAndGetModel(): Model
 
+  // TODO: Make type-safe?
+  def runExternal(evaluate: String*): List[String]
+  def resultToInt(result: String): Int
+  def resultToLong(result: String): Long
+
   def debugPrint(): Unit
 
   def maximize[R <: Sort](expr: Expr[R]): Optimize.Handle[R]
   def minimize[R <: Sort](expr: Expr[R]): Optimize.Handle[R]
 
-  def constant(n: Int): IntNum
-  def constant(n: Long): IntNum
+  def realConstant(n: Int): RatNum
+  def intConstant(n: Int): IntNum
+  def longConstant(n: Long): IntNum
+
+  def intToReal(n: Expr[IntSort]): RealExpr
 
   def labeledBool(label: String): BoolExpr
+  def labeledReal(label: String): RealExpr
   def labeledInt(label: String): IntExpr
 
   def equal(a: Expr[_], b: Expr[_]): BoolExpr
@@ -51,6 +71,12 @@ trait Optimizer {
     a: Expr[A],
     b: Expr[B],
   ): BoolExpr
+
+  def greater[A <: ArithSort, B <: ArithSort](
+    a: Expr[A],
+    b: Expr[B],
+  ): BoolExpr
+
   def lessOrEqual[A <: ArithSort, B <: ArithSort](
     a: Expr[A],
     b: Expr[B],
@@ -67,26 +93,40 @@ trait Optimizer {
   def and(expressions: Expr[BoolSort]*): BoolExpr
   def or(expressions: Expr[BoolSort]*): BoolExpr
 
+  @deprecated("Use `runExternal` instead", "2023-12-24")
   def extractBoolean(b: BoolExpr): Option[Boolean]
+
+  @deprecated("Use `runExternal` instead", "2023-12-24")
   def extractInt(n: IntExpr): Int
+
+  @deprecated("Use `runExternal` instead", "2023-12-24")
   def extractLong(n: IntExpr): Long
 }
 
 private class Z3Optimizer(val context: Context, val optimize: Optimize)
     extends Optimizer {
-  val Zero: IntNum     = constant(0)
-  val One: IntNum      = constant(1)
-  val MinusOne: IntNum = constant(-1)
+  val Zero: IntNum     = intConstant(0)
+  val One: IntNum      = intConstant(1)
+  val MinusOne: IntNum = intConstant(-1)
 
   val False: BoolExpr = context.mkBool(false)
   val True: BoolExpr  = context.mkBool(true)
 
-  def constant(n: Int): IntNum =
+  def realConstant(n: Int): RatNum =
+    context.mkReal(n)
+
+  def intConstant(n: Int): IntNum =
     context.mkInt(n)
 
-  def constant(n: Long): IntNum =
+  def longConstant(n: Long): IntNum =
     context.mkInt(n)
 
+  def intToReal(n: Expr[IntSort]): RealExpr =
+    context.mkInt2Real(n)
+
+  def labeledReal(label: String): RealExpr =
+    context.mkRealConst(label)
+
   def labeledInt(label: String): IntExpr =
     context.mkIntConst(label)
 
@@ -96,6 +136,12 @@ private class Z3Optimizer(val context: Context, val optimize: Optimize)
   def equal(a: Expr[_], b: Expr[_]): BoolExpr =
     context.mkEq(a, b)
 
+  def greater[A <: ArithSort, B <: ArithSort](
+    a: Expr[A],
+    b: Expr[B],
+  ): BoolExpr =
+    context.mkGt(a, b)
+
   def greaterOrEqual[A <: ArithSort, B <: ArithSort](
     a: Expr[A],
     b: Expr[B],
@@ -119,10 +165,48 @@ private class Z3Optimizer(val context: Context, val optimize: Optimize)
 
   def checkAndGetModel(): Model = {
     val status = optimize.Check()
-    assert(status == Status.SATISFIABLE)
+    assert(status == Status.SATISFIABLE, "Model is not satisfiable")
     optimize.getModel
   }
 
+  def runExternal(evaluate: String*): List[String] = {
+    val debug = false
+
+    val programStart = optimize.toString
+
+    // Note - `get-value` is SMT-LIB spec, but `eval` works with Z3
+    val programEnd = evaluate
+      .map { what =>
+        s"(eval $what)"
+      }
+      .mkString("\n")
+
+    val program = s"$programStart\n$programEnd\n(exit)"
+
+    if (debug) println(program)
+
+    val results = Runner.runSync("z3", "-in")(program)
+
+    if (debug) println(results)
+
+    val lines = results.splitLines
+    lines.head.trim shouldEqual "sat"
+    lines.tail.size shouldEqual evaluate.length
+    lines.tail
+  }
+
+  def resultToInt(result: String): Int =
+    result.trim match {
+      case s"(- $n)" => -n.toInt
+      case other     => other.toInt
+    }
+
+  def resultToLong(result: String): Long =
+    result.trim match {
+      case s"(- $n)" => -n.toLong
+      case other     => other.toLong
+    }
+
   def debugPrint(): Unit =
     println(optimize)
 
@@ -197,12 +281,12 @@ object Optimizer {
 object ImplicitConversions {
   implicit class RichInt(val int: Int) {
     def constant(implicit optimizer: Optimizer): IntExpr =
-      optimizer.constant(int)
+      optimizer.intConstant(int)
   }
 
   implicit class RichLong(val long: Long) {
     def constant(implicit optimizer: Optimizer): IntExpr =
-      optimizer.constant(long)
+      optimizer.longConstant(long)
   }
 
   implicit class RichString(val string: String) {
@@ -225,6 +309,11 @@ object ImplicitConversions {
     def >=(other: Expr[B])(implicit
       optimizer: Optimizer
     ): BoolExpr = optimizer.greaterOrEqual(expr, other)
+
+    def >(other: Expr[B])(implicit
+      optimizer: Optimizer
+    ): BoolExpr = optimizer.greater(expr, other)
+
     def <=(other: Expr[B])(implicit
       optimizer: Optimizer
     ): BoolExpr = optimizer.lessOrEqual(expr, other)
diff --git a/scala2/src/main/scala/jurisk/process/Runner.scala b/scala2/src/main/scala/jurisk/process/Runner.scala
new file mode 100644
index 00000000..a00534d8
--- /dev/null
+++ b/scala2/src/main/scala/jurisk/process/Runner.scala
@@ -0,0 +1,28 @@
+package jurisk.process
+
+import org.scalatest.matchers.should.Matchers.convertToAnyShouldWrapper
+
+import scala.io.Source
+import scala.sys.process.ProcessIO
+import scala.sys.process.stringSeqToProcess
+
+object Runner {
+  def runSync(command: String*)(input: String): String = {
+    val outputBuilder = new StringBuilder
+    val io            = new ProcessIO(
+      in => {
+        in.write(input.getBytes)
+        in.close()
+      },
+      out => {
+        val output = Source.fromInputStream(out).mkString
+        outputBuilder.append(output)
+        out.close()
+      },
+      _.close(),
+    )
+    val process       = command.run(io)
+    process.exitValue() shouldEqual 0
+    outputBuilder.toString()
+  }
+}
diff --git a/scala2/src/test/scala/jurisk/adventofcode/y2023/Advent11Spec.scala b/scala2/src/test/scala/jurisk/adventofcode/y2023/Advent11Spec.scala
index 766ec3d4..e91a3719 100644
--- a/scala2/src/test/scala/jurisk/adventofcode/y2023/Advent11Spec.scala
+++ b/scala2/src/test/scala/jurisk/adventofcode/y2023/Advent11Spec.scala
@@ -2,7 +2,7 @@ package jurisk.adventofcode.y2023
 
 import org.scalatest.freespec.AnyFreeSpec
 import Advent11._
-import jurisk.geometry.BigIntCoords2D
+import jurisk.geometry.Coordinates2D
 import org.scalatest.matchers.should.Matchers._
 
 import scala.collection.immutable.ArraySeq
@@ -12,65 +12,65 @@ class Advent11Spec extends AnyFreeSpec {
     "test 1" in {
       expand(
         ArraySeq(
-          BigIntCoords2D(0, 0),
-          BigIntCoords2D(0, 1),
+          Coordinates2D[BigInt](0, 0),
+          Coordinates2D[BigInt](0, 1),
         ),
         2,
       ) shouldEqual ArraySeq(
-        BigIntCoords2D(0, 0),
-        BigIntCoords2D(0, 1),
+        Coordinates2D[BigInt](0, 0),
+        Coordinates2D[BigInt](0, 1),
       )
     }
 
     "test 2" in {
       expand(
         ArraySeq(
-          BigIntCoords2D(0, 0),
-          BigIntCoords2D(0, 2),
+          Coordinates2D[BigInt](0, 0),
+          Coordinates2D[BigInt](0, 2),
         ),
         2,
       ) shouldEqual ArraySeq(
-        BigIntCoords2D(0, 0),
-        BigIntCoords2D(0, 3),
+        Coordinates2D[BigInt](0, 0),
+        Coordinates2D[BigInt](0, 3),
       )
     }
 
     "test 3" in {
       expand(
         ArraySeq(
-          BigIntCoords2D(0, 0),
-          BigIntCoords2D(0, 3),
+          Coordinates2D[BigInt](0, 0),
+          Coordinates2D[BigInt](0, 3),
         ),
         2,
       ) shouldEqual ArraySeq(
-        BigIntCoords2D(0, 0),
-        BigIntCoords2D(0, 5),
+        Coordinates2D[BigInt](0, 0),
+        Coordinates2D[BigInt](0, 5),
       )
     }
 
     "test 4" in {
       expand(
         ArraySeq(
-          BigIntCoords2D(0, 0),
-          BigIntCoords2D(0, 1),
+          Coordinates2D[BigInt](0, 0),
+          Coordinates2D[BigInt](0, 1),
         ),
         100_000,
       ) shouldEqual ArraySeq(
-        BigIntCoords2D(0, 0),
-        BigIntCoords2D(0, 1),
+        Coordinates2D[BigInt](0, 0),
+        Coordinates2D[BigInt](0, 1),
       )
     }
 
     "test 5" in {
       expand(
         ArraySeq(
-          BigIntCoords2D(0, 0),
-          BigIntCoords2D(0, 2),
+          Coordinates2D[BigInt](0, 0),
+          Coordinates2D[BigInt](0, 2),
         ),
         100_000,
       ) shouldEqual ArraySeq(
-        BigIntCoords2D(0, 0),
-        BigIntCoords2D(0, 100_001),
+        Coordinates2D[BigInt](0, 0),
+        Coordinates2D[BigInt](0, 100_001),
       )
     }
 
diff --git a/scala2/src/test/scala/jurisk/adventofcode/y2023/Advent24Spec.scala b/scala2/src/test/scala/jurisk/adventofcode/y2023/Advent24Spec.scala
new file mode 100644
index 00000000..f64c6c32
--- /dev/null
+++ b/scala2/src/test/scala/jurisk/adventofcode/y2023/Advent24Spec.scala
@@ -0,0 +1,76 @@
+package jurisk.adventofcode.y2023
+
+import org.scalatest.freespec.AnyFreeSpec
+import Advent24._
+import jurisk.geometry.Coords3D
+import org.scalatest.matchers.should.Matchers._
+
+class Advent24Spec extends AnyFreeSpec {
+  // https://www.math3d.org/G5JHCaIly
+  private def testData = parseFile(fileName("-test"))
+
+  private def realData = parseFile(fileName(""))
+
+  val expectedTestAnswer: PositionAndVelocity3D = PositionAndVelocity3D(
+    Coords3D[Long](24, 13, 10),
+    Coords3D[Long](-3, 1, 2),
+  )
+
+  val expectedRealAnswer: PositionAndVelocity3D = PositionAndVelocity3D(
+    Coords3D[Long](191146615936494L, 342596108503183L, 131079628110881L),
+    Coords3D[Long](139, -93, 245),
+  )
+
+  "part 1" - {
+    "test" in {
+      part1(testData, 7, 27) shouldEqual 2
+    }
+
+    "real" in {
+      // Not 24198
+      part1(realData, 200000000000000L, 400000000000000L) shouldEqual 24192
+    }
+  }
+
+  "linesIntersect" in {
+    testData foreach { rock =>
+      linesIntersect(rock, expectedTestAnswer) shouldEqual true
+      linesIntersect(
+        rock,
+        PositionAndVelocity3D(
+          Coords3D[Long](-3, -4, -1),
+          Coords3D[Long](1, 2, 3),
+        ),
+      ) shouldEqual false
+    }
+  }
+
+  "part 2 Z3" - {
+    "test optimizer" in {
+      solvePart2Optimizer(testData) shouldEqual expectedTestAnswer
+    }
+
+    "real optimizer" in {
+      solvePart2Optimizer(realData) shouldEqual expectedRealAnswer
+    }
+  }
+
+  "part 2" - {
+    val expectedResult = 664822352550558L
+
+    // This works but is not fast enough to leave enabled
+    "real inferring velocity" ignore {
+      solve2InferringVelocity(realData) shouldEqual expectedRealAnswer
+    }
+
+    "real using CRT" ignore {
+      // Not actually implemented
+      solve2UsingChineseRemainderTheorem(realData) shouldEqual expectedResult
+    }
+
+    // This works but is not fast enough to leave enabled
+    "real" ignore {
+      part2(realData) shouldEqual expectedResult
+    }
+  }
+}
diff --git a/scala2/src/test/scala/jurisk/geometry/Coords3DSpec.scala b/scala2/src/test/scala/jurisk/geometry/Coords3DSpec.scala
new file mode 100644
index 00000000..0edbe7c4
--- /dev/null
+++ b/scala2/src/test/scala/jurisk/geometry/Coords3DSpec.scala
@@ -0,0 +1,27 @@
+package jurisk.geometry
+
+import org.scalatest.freespec.AnyFreeSpec
+import org.scalatest.matchers.should.Matchers._
+import jurisk.geometry.Coords3D.areVectorsParallel
+
+class Coords3DSpec extends AnyFreeSpec {
+  // https://www.mathsisfun.com/algebra/vectors-cross-product.html
+  "crossProduct" - {
+    (Coords3D[Long](2, 3, 4) crossProduct Coords3D[Long](
+      5,
+      6,
+      7,
+    )) shouldEqual Coords3D[Long](
+      -3,
+      6,
+      -3,
+    )
+  }
+
+  "areVectorsParallel" in {
+    val a = Coords3D[Long](5, 2, -1)
+    val b = Coords3D[Long](-10, -4, 2)
+    areVectorsParallel(a, b) shouldEqual true
+    areVectorsParallel(a, Coords3D[Long](-3, 1, 2)) shouldEqual false
+  }
+}
diff --git a/scala2/src/test/scala/jurisk/math/GaussianEliminationSpec.scala b/scala2/src/test/scala/jurisk/math/GaussianEliminationSpec.scala
new file mode 100644
index 00000000..ab2da0ec
--- /dev/null
+++ b/scala2/src/test/scala/jurisk/math/GaussianEliminationSpec.scala
@@ -0,0 +1,57 @@
+package jurisk.math
+
+import jurisk.math.GaussianElimination.solve
+import org.scalatest.freespec.AnyFreeSpec
+import org.scalatest.matchers.should.Matchers._
+
+class GaussianEliminationSpec extends AnyFreeSpec {
+  private def compare(
+    solution: Array[Double],
+    expected: Array[Double],
+  ): Unit = {
+    val Eps = 1e-5
+
+    solution.length shouldEqual expected.length
+    for ((v, i) <- solution.zipWithIndex)
+      assert(
+        math.abs(v - expected(i)) < Eps,
+        s"Value at index $i should be close to ${expected(i)}",
+      )
+  }
+
+  "GaussianElimination" - {
+    "test 1" in {
+      // 2x + y -z = 8
+      // -3x -y + 2z = -11
+      // -2x +y + 2z = -3
+      val A        = Array(
+        Array(2.0, 1.0, -1.0),
+        Array(-3.0, -1.0, 2.0),
+        Array(-2.0, 1.0, 2.0),
+      )
+      val b        = Array(8.0, -11.0, -3.0)
+      val expected = Array(2.0, 3.0, -1.0)
+      val solution = solve(A, b)
+
+      compare(solution, expected)
+    }
+
+    // https://en.wikipedia.org/wiki/Gaussian_elimination#Example_of_the_algorithm
+    "test from Wikipedia" in {
+      // 2x + y - z = 8
+      // -3x -y + 2z = -11
+      // -2x + y + 2z = -3
+      val A = Array(
+        Array(2.0, 1.0, -1.0),
+        Array(-3.0, -1.0, 2.0),
+        Array(-2.0, 1.0, 2.0),
+      )
+
+      val b        = Array(8.0, -11.0, -3.0)
+      val expected = Array(2.0, 3.0, -1.0)
+      val solution = solve(A, b)
+
+      compare(solution, expected)
+    }
+  }
+}