text:
Gm, Gn = 100, 100;
Gcolptr = [1, 17, 33, 50, 65, 74, 87, 101, 113, 124, 136, 149, 170, 182, 195, 211, 224, 234, 251, 271, 284, 299, 317, 329, 345, 362, 377, 390, 402, 415, 429, 450, 464, 478, 496, 506, 525, 543, 567, 580, 590, 602, 612, 633, 649, 662, 673, 686, 699, 721, 736, 747, 769, 786, 797, 810, 830, 839, 851, 867, 883, 898, 919, 938, 952, 962, 978, 994, 999, 1016, 1027, 1039, 1049, 1067, 1079, 1093, 1109, 1129, 1147, 1168, 1186, 1198, 1213, 1222, 1232, 1246, 1261, 1275, 1288, 1297, 1309, 1319, 1337, 1358, 1373, 1387, 1401, 1414, 1432, 1442, 1461];
Growval = [3, 6, 8, 16, 18, 22, 36, 37, 42, 51, 59, 66, 87, 90, 93, 97, 4, 9, 10, 15, 26, 35, 40, 42, 46, 47, 52, 77, 80, 85, 92, 94, 1, 8, 15, 16, 19, 22, 32, 39, 42, 49, 50, 56, 58, 61, 80, 82, 98, 2, 6, 15, 19, 25, 29, 36, 47, 53, 73, 85, 87, 88, 89, 100, 11, 39, 43, 44, 47, 77, 82, 93, 98, 1, 4, 18, 22, 24, 57, 63, 64, 75, 78, 87, 91, 100, 13, 15, 19, 21, 34, 36, 38, 56, 59, 60, 67, 70, 78, 90, 1, 3, 20, 33, 35, 45, 63, 70, 74, 80, 82, 91, 2, 14, 22, 25, 36, 40, 44, 50, 56, 59, 93, 2, 11, 35, 38, 41, 47, 61, 62, 63, 81, 83, 92, 5, 10, 13, 23, 27, 29, 38, 47, 52, 55, 76, 89, 93, 20, 24, 28, 32, 39, 47, 49, 50, 57, 63, 64, 69, 72, 73, 78, 86, 91, 93, 95, 98, 100, 7, 11, 22, 28, 34, 41, 43, 49, 54, 55, 78, 79, 9, 23, 27, 32, 36, 37, 38, 57, 65, 67, 77, 88, 90, 2, 3, 4, 7, 20, 27, 30, 31, 38, 53, 62, 75, 89, 95, 98, 100, 1, 3, 19, 31, 32, 38, 40, 42, 55, 61, 78, 94, 95, 17, 26, 39, 61, 62, 80, 81, 83, 87, 92, 1, 6, 19, 26, 34, 37, 43, 46, 51, 53, 64, 65, 79, 82, 91, 94, 96, 3, 4, 7, 16, 18, 33, 36, 43, 48, 51, 58, 60, 62, 67, 69, 73, 93, 94, 95, 100, 8, 12, 15, 32, 56, 61, 63, 69, 70, 73, 85, 86, 98, 7, 29, 31, 32, 35, 43, 45, 52, 53, 55, 60, 66, 74, 77, 81, 1, 3, 6, 9, 13, 23, 31, 38, 48, 50, 53, 54, 63, 75, 76, 85, 89, 100, 14, 22, 34, 38, 52, 76, 78, 84, 86, 91, 92, 95, 6, 12, 26, 31, 36, 38, 39, 50, 52, 53, 57, 74, 82, 89, 97, 98, 4, 9, 30, 33, 34, 41, 43, 45, 46, 49, 56, 62, 63, 66, 69, 77, 87, 2, 17, 18, 24, 29, 38, 46, 51, 58, 69, 71, 73, 78, 79, 80, 11, 14, 15, 31, 41, 43, 56, 60, 69, 71, 90, 93, 96, 12, 13, 31, 41, 42, 49, 52, 60, 62, 66, 75, 84, 4, 11, 21, 26, 31, 33, 46, 61, 72, 73, 79, 88, 96, 15, 25, 34, 39, 41, 44, 49, 51, 62, 63, 66, 82, 83, 94, 15, 16, 21, 22, 24, 27, 29, 32, 36, 45, 48, 52, 53, 60, 71, 75, 77, 82, 86, 96, 100, 3, 12, 14, 16, 20, 21, 31, 49, 64, 65, 69, 73, 76, 79, 8, 19, 25, 29, 40, 49, 67, 71, 74, 79, 84, 85, 97, 98, 7, 13, 18, 23, 25, 30, 50, 51, 58, 62, 63, 65, 66, 71, 73, 80, 88, 99, 2, 8, 10, 21, 50, 60, 75, 80, 93, 96, 1, 4, 7, 9, 14, 19, 21, 24, 31, 37, 53, 58, 61, 73, 74, 78, 83, 93, 100, 1, 14, 18, 24, 36, 38, 48, 58, 59, 62, 66, 73, 79, 86, 87, 89, 92, 99, 7, 10, 11, 14, 16, 22, 23, 24, 26, 37, 39, 42, 49, 59, 69, 71, 77, 79, 82, 84, 87, 95, 96, 98, 3, 12, 17, 24, 30, 38, 45, 48, 59, 74, 78, 80, 92, 2, 9, 16, 33, 41, 49, 77, 78, 88, 93, 10, 13, 25, 27, 28, 30, 40, 54, 67, 71, 77, 80, 1, 3, 16, 38, 55, 61, 64, 67, 76, 88, 5, 10, 13, 18, 21, 25, 27, 43, 44, 46, 49, 53, 54, 55, 56, 59, 71, 75, 82, 85, 96, 5, 9, 30, 43, 47, 49, 56, 64, 73, 79, 88, 90, 93, 94, 98, 100, 8, 21, 25, 31, 39, 48, 52, 62, 63, 72, 81, 93, 97, 2, 18, 25, 26, 29, 43, 54, 72, 92, 94, 97, 2, 4, 5, 10, 11, 12, 44, 50, 71, 83, 90, 93, 97, 19, 22, 31, 37, 39, 45, 52, 56, 62, 75, 76, 77, 86, 3, 12, 13, 25, 28, 30, 32, 33, 38, 40, 43, 44, 52, 53, 59, 61, 75, 76, 84, 92, 93, 97, 3, 6, 9, 12, 22, 24, 31, 34, 35, 47, 57, 59, 83, 92, 94, 1, 18, 19, 26, 30, 34, 56, 59, 66, 96, 100, 1, 2, 11, 17, 18, 21, 23, 24, 28, 31, 45, 48, 49, 57, 76, 77, 78, 85, 88, 95, 98, 100, 4, 15, 18, 21, 22, 24, 31, 43, 61, 64, 66, 73, 74, 75, 76, 79, 100, 13, 22, 41, 43, 46, 59, 62, 64, 77, 81, 88, 11, 13, 16, 21, 42, 43, 55, 59, 62, 67, 78, 88, 89, 3, 7, 9, 20, 25, 27, 43, 44, 48, 51, 60, 62, 65, 69, 73, 78, 79, 81, 86, 87, 6, 12, 14, 24, 52, 64, 66, 88, 93, 3, 19, 26, 34, 37, 63, 70, 71, 85, 86, 91, 95, 1, 7, 9, 37, 38, 39, 43, 49, 50, 51, 54, 55, 63, 68, 74, 90, 7, 19, 21, 27, 28, 31, 35, 56, 61, 63, 76, 77, 80, 92, 94, 100, 3, 10, 16, 17, 20, 29, 42, 49, 53, 60, 72, 77, 78, 90, 100, 10, 15, 17, 28, 30, 34, 37, 45, 48, 54, 55, 56, 70, 78, 79, 80, 81, 82, 86, 97, 100, 6, 8, 10, 12, 20, 22, 25, 34, 45, 58, 59, 60, 63, 67, 84, 86, 91, 93, 96, 6, 12, 18, 32, 42, 44, 53, 54, 57, 78, 84, 85, 95, 98, 14, 18, 32, 34, 56, 65, 86, 91, 93, 95, 1, 21, 25, 28, 30, 34, 37, 43, 51, 53, 57, 70, 72, 75, 82, 85, 7, 14, 19, 33, 41, 42, 55, 63, 67, 68, 70, 74, 78, 87, 92, 94, 59, 67, 69, 70, 94, 8, 12, 14, 19, 20, 25, 26, 27, 32, 38, 56, 68, 71, 92, 96, 98, 100, 7, 8, 20, 58, 62, 66, 67, 68, 73, 81, 89, 26, 27, 31, 33, 34, 38, 41, 43, 47, 58, 71, 93, 12, 29, 45, 46, 61, 66, 74, 80, 97, 100, 4, 8, 12, 19, 20, 26, 29, 32, 33, 34, 37, 44, 53, 56, 57, 70, 82, 87, 8, 21, 24, 33, 36, 39, 53, 59, 67, 72, 80, 97, 6, 15, 22, 28, 31, 35, 43, 48, 49, 53, 66, 87, 92, 100, 11, 22, 23, 32, 40, 42, 48, 49, 52, 53, 60, 62, 77, 91, 94, 96, 2, 5, 14, 21, 25, 31, 38, 40, 41, 48, 52, 54, 60, 61, 76, 77, 86, 88, 94, 96, 6, 7, 12, 13, 16, 23, 26, 39, 40, 52, 55, 56, 61, 62, 64, 67, 78, 83, 13, 18, 26, 29, 32, 33, 37, 38, 44, 53, 56, 62, 70, 72, 82, 90, 92, 94, 97, 98, 99, 3, 8, 17, 26, 34, 35, 39, 41, 60, 62, 72, 74, 80, 84, 85, 95, 96, 99, 10, 17, 21, 32, 45, 54, 56, 59, 62, 70, 81, 99, 3, 5, 8, 18, 24, 30, 31, 38, 43, 62, 66, 73, 79, 90, 92, 10, 17, 24, 30, 36, 47, 50, 78, 96, 23, 28, 33, 38, 49, 63, 64, 80, 90, 94, 2, 4, 20, 22, 33, 43, 52, 58, 64, 66, 80, 95, 98, 99, 12, 20, 23, 31, 37, 48, 56, 58, 62, 63, 65, 77, 89, 94, 98, 1, 4, 6, 10, 17, 25, 37, 38, 56, 67, 73, 75, 87, 99, 4, 29, 34, 40, 42, 44, 52, 54, 55, 57, 77, 94, 96, 4, 15, 22, 24, 37, 55, 70, 86, 90, 1, 7, 14, 27, 44, 47, 59, 61, 82, 84, 89, 91, 6, 8, 12, 18, 23, 58, 63, 65, 76, 90, 2, 10, 17, 23, 37, 39, 46, 49, 50, 60, 67, 69, 75, 79, 82, 93, 95, 98, 1, 5, 9, 11, 12, 19, 27, 35, 36, 40, 44, 45, 47, 49, 57, 59, 63, 65, 71, 92, 98, 16, 18, 19, 30, 44, 46, 50, 67, 68, 77, 79, 84, 86, 88, 97, 12, 15, 16, 19, 23, 38, 52, 58, 64, 65, 80, 85, 92, 100, 18, 27, 29, 31, 35, 38, 43, 51, 63, 69, 77, 80, 83, 88, 1, 24, 33, 45, 46, 47, 49, 62, 72, 74, 79, 94, 97, 3, 5, 12, 15, 20, 24, 33, 38, 44, 52, 64, 69, 79, 85, 86, 92, 93, 98, 15, 23, 34, 37, 79, 80, 81, 85, 87, 100, 4, 6, 12, 15, 19, 22, 31, 44, 51, 52, 53, 60, 61, 62, 69, 72, 75, 95, 99];
Gnzval = [0.3086314079058456, 0.9986607788380012, 0.5610064890329525, 0.832224131008404, 0.3294462924259989, 0.4727573932124432, 0.5229458380085813, 0.27087256414111605, 0.5599319687443616, 0.5475404884489812, 0.5640382080194491, 0.8584089047683475, 0.553751100697396, 0.8421209009050953, 0.28010922081880957, 0.4894927072798642, 0.5202896174970701, 0.6858392093301448, 0.786405109732176, 0.6097607678487019, 1.0, 0.6345427724624766, 0.7779020535541951, 0.25997986371934084, 0.5685732722128215, 0.5820882261494038, 1.0, 0.6691488365587362, 0.35646490992886976, 0.6804455112852272, 0.8616653682395536, 0.3844247395512977, 0.31061552646480806, 0.5557190259501611, 0.8344029925675192, 1.0, 0.7358803644161728, 0.8944100283447786, 0.9498171618334101, 0.3549166996682608, 0.9456124456334171, 0.41963854845994036, 0.2932962208715996, 0.9366584806486538, 0.4948237860641257, 0.6916280365570763, 0.890500050923668, 0.3407142556667635, 0.6368968698231599, 0.462159532657063, 0.46642958841536775, 1.0, 0.7621858334934764, 0.4498556705758161, 0.7381551854211765, 0.2977304585846391, 0.6309825101084697, 0.7946260974635482, 1.0, 0.7657409178932815, 0.4205882539525326, 1.0, 0.43161617848205563, 0.5048452061495825, 0.3059623993159188, 0.3085328961710782, 0.4790575785932825, 1.0, 0.872322982002817, 1.0, 0.5786347998173016, 0.5071192436353391, 0.4049779585811116, 0.920366678728016, 0.3081132295673938, 0.6059122748103212, 0.6827108182507844, 0.802145824129018, 0.6732235741265936, 0.2722809797775942, 0.48633048354281944, 0.9130987251694483, 0.16792928811770605, 0.34465045261583066, 0.47065602553857, 0.4152126837429248, 0.5174992252770355, 0.35010018108319557, 0.708378714834061, 0.9235904721119229, 0.7261615415705166, 0.47566022237080685, 0.6872735540303211, 0.5599559422337107, 0.7685025598243236, 1.0, 0.8874895387402504, 0.6068787255265603, 0.41977628238218145, 0.9542714496902132, 0.5574670888347338, 0.548685613998101, 0.8908416680709261, 0.6504234380203574, 0.8941542269678509, 0.29802766178781015, 0.47812860970903437, 1.0, 0.9886239285166574, 0.9267754556210938, 0.7789307726436544, 0.3768980043566626, 0.5341665029010937, 1.0, 0.3349242596057137, 0.47693120172909714, 0.2790158466233321, 0.5763517627725415, 0.5137867593451328, 0.40790826601308333, 0.6612067783335817, 0.9999999999999999, 0.31913420867592435, 0.4312155324016517, 0.37405891327138613, 0.3007648865513607, 0.3750386866611699, 0.309152831000324, 0.6018677494973111, 0.5313504822646872, 0.8855888387317028, 0.6593257046490747, 0.2734026971173483, 0.6256687678067726, 0.722787256262097, 0.355041268121685, 0.7804479145324065, 0.9296632825397264, 0.36783724094285836, 0.8902546936560085, 0.48038465212718906, 0.7086538024223994, 0.8022704681788595, 0.5001305302545863, 1.0, 1.0, 0.41837665694554327, 0.4075814930363359, 0.7929646349176189, 0.4916207805477212, 0.31387063055087305, 0.657301102842506, 0.9657214191689005, 0.8693893582752115, 0.4054034595807778, 0.31936383368753135, 0.8748287382049544, 0.8303891613523123, 0.4006804870166255, 0.4916055177188628, 0.3103717560842533, 0.3447168546007089, 0.6290276366365739, 0.3711216629292957, 1.0, 0.26831572897486544, 0.27176949563348846, 0.20486921466022615, 0.9535941533408456, 0.37801038581718915, 0.47912326991080945, 0.847614122252833, 0.3747480665523656, 0.9214176036453711, 0.6404184434802861, 0.48441445048853726, 0.6167570644152066, 0.59675201571596, 0.2762906868028233, 0.4817333511756801, 0.4581126124044482, 1.0, 0.9650416483958933, 0.4009033228853521, 0.3803680611928475, 0.48180188014863484, 0.5950729876118684, 0.350326371357497, 0.7574416672442622, 0.40523520126285845, 0.9208640751666125, 0.6080920934486442, 0.2629408708546492, 0.9134491924862159, 0.44727467968541407, 0.955167486928693, 0.8257870713518156, 0.31814684997020104, 0.4772647277943498, 0.7371304308743143, 0.2760468041615327, 0.8321000281285578, 0.27239515144078236, 0.6747324130335252, 0.7009431652813863, 0.5474741686501907, 1.0, 0.43827760510730873, 0.6154592493116096, 0.40112166434289687, 0.837574303742781, 1.0, 0.7621509137437208, 0.8155251268045675, 0.3904904672349383, 0.6366371934445823, 0.410935953054335, 0.4839843513641727, 0.3101664281028919, 0.3393953336740818, 0.3298485572204908, 0.6757876973942154, 0.91993060313145, 0.2707960676353462, 0.8622383089816502, 0.5457270194362043, 0.5486697671725652, 0.952958041849733, 0.26428354744482635, 0.27508440313548116, 0.3816366169968296, 0.3662130020276344, 0.8889719791354472, 0.32444472687347065, 0.6474749239281873, 0.6411402419432995, 0.6245097323233761, 0.6510442541338313, 0.4702533398908854, 0.37598154112230303, 0.9999999999999999, 0.26138834570880637, 0.7945072898658906, 0.7074072722650838, 0.4229775883253915, 0.4400490072405359, 0.9172833432509847, 0.8118319673613743, 0.9561356611899414, 0.8347696202928084, 0.8594870728673177, 0.6421809346782505, 0.6567941132491838, 0.8901703179109118, 0.7652654670164794, 0.4453368303254978, 0.28916246417069635, 0.17795419188577394, 0.6702677239346234, 0.7098350156700307, 1.0, 0.47468107139365234, 0.2941245487007703, 0.7257573879434771, 0.6507952217976634, 0.43954956743878454, 0.45952066000717945, 0.9337439682708898, 0.6854476125940218, 0.9487761238241454, 0.8573318833509846, 0.752112912544934, 0.39616164044455604, 0.3166397285430773, 0.5841486462005989, 0.9319630547364821, 0.9999999999999999, 0.5383042760589416, 0.8899049990468427, 0.46938107865435424, 0.25567160455017157, 0.49823507223720315, 0.5471428244808886, 0.7962371481875565, 0.3238765253529882, 0.6733703195618738, 0.7192527233194029, 0.4131090917796851, 0.2691562739573891, 0.4256821317592017, 0.25439192925453463, 0.6358771028367011, 0.657956321892106, 0.6563570883838363, 1.0, 0.9215848533463602, 0.5343279983717063, 0.558369331848037, 0.46990112530424144, 0.8833275486211594, 0.7363122865303411, 0.2721375921489712, 0.9097475297136448, 0.7613312874547992, 0.7196003470106229, 1.0, 0.7232728494264885, 0.6663223366241441, 0.9274687279300606, 0.8998834927757241, 0.935561005534593, 0.2604182043087865, 0.3537976264100425, 0.978984603815119, 0.4777341524354235, 0.9999999999999999, 0.7962171337091215, 0.7730671304326604, 0.49940274218768654, 0.44622816545327215, 0.3321449022834124, 0.36044979596884097, 0.619905267027834, 0.34917121244250726, 0.3727204870081016, 0.26346749345638143, 1.0, 0.706003368851335, 0.8818591720737212, 0.49402446840048186, 0.8093452273883495, 0.9564689751182814, 0.4634300374278596, 0.8884708791703964, 0.5161402332720586, 0.49456436270377074, 0.5451831220279059, 0.65658680859289, 1.0, 0.2729025348554793, 0.9999999999999999, 0.6456532480586678, 0.43318781552554086, 0.3558695906054176, 0.46743240537946834, 0.5651890889543015, 0.41760134213989253, 0.4402829666798088, 0.6491186685078362, 0.43902813194792023, 0.29904093618817956, 0.35262270256777767, 0.5641866517325925, 0.6772352917025517, 0.9953380728771916, 0.3601756090150912, 0.8555236022009636, 0.7794569214792634, 1.0, 0.995124718310694, 0.3262465932019833, 0.6419430379775529, 0.9380675084389861, 0.6393688768607003, 0.8054032624000077, 0.43996064278903346, 0.3406747075460717, 0.42102092250636064, 0.7914618877981048, 0.701233635327357, 0.2915836611573606, 0.9246839136724336, 0.41152481895822207, 0.24350864496574667, 0.4745312733810382, 0.6502358407148993, 0.5032435291507724, 0.33289036749328604, 0.6498983261647092, 0.5165878402860212, 0.7023031233844951, 0.5166004956360988, 0.8256138341278212, 0.4202378828379737, 0.6962437161766372, 1.0, 0.44222921829094464, 0.4825606599276818, 0.6710159392890304, 0.3948938915853573, 0.5133271763205812, 0.2760479916618351, 0.256138376060717, 0.2716985357630993, 1.0, 0.382603844013551, 0.4457420861707507, 0.6825692906575722, 0.9155577380561006, 0.7504881703906944, 0.8037299095931034, 0.5032892997033424, 0.25430442424902366, 0.28187315392847645, 0.41508220263838863, 0.31678737850158967, 0.31237847802715524, 0.6816426142634668, 0.8719434864674187, 0.96049199319049, 0.3445584277408837, 0.5934597136558226, 0.9152142090460307, 0.8089271612484954, 0.3307484272690676, 0.397648557429565, 0.9358925992417119, 0.4387671258559621, 1.0, 0.38985552603408025, 0.56703699678337, 0.4245550391872948, 0.9651729998467905, 0.2757142763685669, 0.9030755117610851, 0.7868025600441677, 0.7588565341217742, 1.0, 0.7103207213826934, 0.7969296254820437, 0.605905950761958, 0.7120145595120534, 0.9284269767348193, 0.5001750782201424, 0.32752531609727864, 0.6043747336116498, 0.32058585445291765, 0.27944254355910797, 0.37972112745119935, 0.46455015919201526, 0.9092388253466036, 0.37486297610865993, 0.9682285164575469, 0.30448420765603024, 0.33869965540444813, 0.43329737156662806, 0.3423161645026107, 0.3764650105620389, 0.3466546013277709, 0.9556331620185714, 0.660127805661213, 0.31485606200044924, 0.3928815512466846, 0.3352739753306991, 0.6663483661513072, 0.6222642836138808, 0.4041153943935262, 0.5910918679880185, 0.5095198592266572, 0.3558946261198822, 0.5220346224267965, 0.5635213451455362, 0.6987583928737265, 0.6645644847059751, 0.38464999639940506, 0.3083613375179511, 0.3334146948854633, 0.5959152372266439, 0.6003738110945085, 0.5652784799033094, 0.9412320687395962, 0.8079874122865662, 0.2763575225860902, 0.5396467695173596, 1.0, 0.368479973077396, 0.2816393654593703, 0.575666248642928, 1.0, 0.664421020978018, 0.49731113816055605, 0.4500375914942832, 0.6814173383136499, 0.6204514992751057, 0.6726548942295332, 0.7050416596431076, 0.33530711344026143, 0.7928331180670295, 0.3745990555617443, 0.30373549925994575, 0.7131305473045673, 0.6772938448089398, 0.6957950769768124, 0.5689165058260527, 0.6521585509066127, 0.43432805729507806, 0.9673640701399818, 0.37543611465128474, 0.4063302938734346, 0.5529924818602802, 0.8479640505437054, 0.9085971209937749, 0.7906136901088885, 0.5048875477328656, 0.3144689801327111, 1.0, 0.40808244821349815, 0.7174451386121848, 0.43616010107985687, 0.7531570672033683, 0.470402944104025, 0.3539756103048682, 0.8746587920488481, 0.6233356683623907, 0.638523285358317, 0.26270118331093906, 0.32774133150760454, 0.39650630788815117, 0.2958412589817853, 0.562806186357827, 0.16259753569594224, 0.6815888356646476, 1.0, 0.4270055875901325, 0.3774614371920048, 0.6778780679457582, 0.6654008683901576, 0.27335668426687665, 0.46530853711820047, 0.4780562808535975, 0.5161740382970629, 0.45006200238276983, 0.6653385391099522, 0.4725501975801607, 0.6149559497605293, 0.4248326399201499, 1.0, 0.7312773380497218, 0.35611638969255266, 0.7116244324017098, 0.3732539229128152, 0.5336007784682106, 0.5300159183438783, 0.42576192761818965, 0.38558573578782324, 0.2799995373982361, 0.6184714744422988, 0.9836544277628815, 0.43227089957532366, 0.8573637205350025, 0.3303560174760293, 0.33541475634223655, 0.2625603418307162, 0.4633122526436575, 0.6292661870103966, 0.50922717426517, 0.5759212361504547, 0.45774336535498544, 0.7689315294108274, 0.646473779995219, 0.6948505219645864, 0.9181620895622001, 0.40719764137238385, 0.9003724549354264, 0.2844102491472594, 0.31203638686983115, 0.35293618726384607, 0.9585926284621957, 0.5893409990333635, 0.5098374194529022, 0.3181156932723296, 0.268764075537574, 0.3771560898385238, 0.8817076507939354, 0.29810970460827185, 0.8488702976787862, 0.6243842619725074, 0.9999999999999999, 0.5951096342970024, 0.4995266470079993, 0.4751897103353425, 0.4118098451032538, 0.59196744338709, 0.6060179864746851, 0.7950699054467367, 0.9352305920789038, 0.6044734153049013, 0.863839470165202, 0.26223246335194483, 0.49551300578501883, 0.9545972057666796, 0.418051585862971, 0.954417922471762, 0.2787251280067237, 1.0, 0.8294336886407009, 0.8923195456499379, 0.8562728084784927, 0.48927689243879846, 0.379322895569401, 0.8109889491857846, 0.5002124492066917, 0.8393621305194641, 0.42960320390770124, 0.42333318201520653, 0.28241246430358097, 0.885745422761262, 0.8340683401116982, 0.4159050448592875, 0.3079139213859449, 0.6510028341534729, 0.41956104882869427, 0.26097208470242467, 0.7094016167007122, 0.5177629458976718, 0.33397243888563666, 0.27976084104334004, 0.6897434111277928, 1.0, 0.3110070287870447, 0.7219541011878216, 0.5906416205330725, 0.6898321262488695, 0.49662505161958015, 0.9728819281023432, 0.8947404416299449, 0.7919896888927374, 0.7177397494822383, 0.586329571628228, 0.5472241445342726, 0.47162143588976024, 0.9185224145802212, 1.0, 0.6503838900646768, 0.4373357728884155, 0.35510988103436963, 0.5683753942155829, 0.726733339796927, 0.7308786207942628, 0.5752443280775552, 0.360548827340001, 0.5902174708408113, 0.4464281978411166, 0.42024594629214274, 0.2932707155439862, 0.8303523088010327, 0.6384493609125711, 0.43815335447254533, 0.5292358085638685, 0.6872503930540391, 0.42922993884158495, 0.5013844422390095, 0.6649857252368991, 0.8446241372253419, 0.5914884853579543, 0.8455681009206062, 0.5034210237150176, 0.537775784839607, 1.0, 0.46420158928294764, 0.42903749170617383, 0.35650912671471297, 0.9999999999999999, 0.3832786650310891, 0.41123625677243164, 0.7057065811413229, 0.5242580948189572, 0.6207834017328003, 0.3982905935174588, 0.8207774117716758, 0.475114392877417, 0.1268731426805701, 0.3573921202088831, 0.43614080776382563, 0.5432406558466367, 0.6739212669548067, 0.4305513727575639, 0.8315182521822679, 0.35395676634800266, 0.6377237827095186, 0.722298302870812, 0.5008986387584945, 0.6282585917003262, 0.5042951048521527, 0.40026274255190925, 0.5755171748885688, 0.7344448054655102, 0.3896944784487999, 0.37386408220000367, 0.8891949140469476, 0.4831084689966445, 0.2717348251611834, 0.37341651641506474, 1.0, 1.0, 0.31734136182613915, 0.4187860885834938, 0.9871813920561348, 0.5289232361522495, 0.5100552082080594, 0.8448293416552926, 0.5926190874142127, 1.0, 0.5210866223691706, 0.5062574945063679, 0.7501182929945124, 0.8164394847855628, 1.0, 0.5524637911750025, 0.5953362180436003, 0.2503510185363501, 0.265811413755178, 0.8300460366148369, 0.6735198647618312, 0.5968020573288837, 0.2809643864511108, 1.0, 0.890158281328718, 0.4065368174188575, 0.4971706487338218, 0.3880328795736865, 0.3493685044571082, 0.4439983768042092, 0.42175501773505975, 0.8926091130399741, 0.6499487554948054, 0.31597055428096626, 0.8219851827529027, 0.6849130796240683, 0.8805411815130135, 0.2707916979388709, 0.35687181808142787, 1.0, 0.9920958768190388, 0.8655009729484989, 0.5522637121054796, 0.26770741047727087, 0.6090774386394321, 0.7917569759030904, 0.3102147801632941, 0.6750464000766107, 1.0, 0.9361868602187147, 0.8593890908716064, 0.4628928940688722, 0.41810300560863434, 0.30949397584156285, 1.0, 0.437150904133923, 0.3727629942198941, 0.30443189438182205, 0.3655617893857382, 0.5168776920142614, 0.8451338241709254, 0.47371910676371953, 0.7418888274741012, 0.5264706731755406, 0.4228555426127055, 0.30997090785691833, 0.8370248886762044, 0.31318586341287574, 0.43933468860559705, 0.5614597418122041, 0.3904823129453313, 0.9656134001340161, 0.8112187383597694, 0.5107940177527805, 0.6485127619073818, 0.5949893160048346, 0.61180051373391, 0.8427655474932756, 0.6077482980825794, 0.9747564147376784, 0.6769683008077197, 0.9657779569121266, 0.5313587786723944, 0.995716796091152, 0.3167861405850656, 0.8818715290441921, 0.7487663920226688, 0.963161249698199, 0.8321440682391955, 0.4006707546156044, 0.5956955905096392, 0.34638512756376066, 0.6212987802844527, 0.8729114509246172, 0.6232620507342075, 0.353069453814275, 0.60776552554974, 0.38806220343038333, 0.4113715208165741, 0.604878889294329, 0.588385940719837, 0.5463966990493874, 0.4857324216313917, 0.5361529630532237, 0.2874308445997272, 0.2969096601458405, 0.5799562404473391, 0.3045630463902186, 0.6911179573109365, 0.0055489581662255, 0.3742846590654724, 0.6754741922167704, 0.5267559776328825, 0.5815667917760697, 0.31013558221962545, 0.37714352193939715, 0.9387036308394823, 0.4454846045686635, 0.5451817907955155, 0.616103735901633, 0.6925268689954808, 0.8198217453201901, 0.6654367741400747, 0.476060815083481, 0.3731630509178315, 1.0, 0.5476724061835141, 0.5191554112563257, 0.25204994140158093, 0.5173583087877003, 0.25618511202342203, 0.5337079799039907, 0.6649651860595284, 0.4841365236203405, 0.7055831565890588, 0.5490601665861466, 0.7401252563328183, 0.8791060496907026, 0.6273526989784678, 1.0, 0.5399521696711336, 0.8168511531918854, 0.4743820263760467, 0.7650555655609802, 0.5075101733080543, 0.5779399315903697, 1.0, 0.8175409915018916, 0.8245955842547197, 0.3796469828317284, 0.568733675479175, 0.34148786539829484, 0.4042871067762432, 0.5546589955753632, 0.7837315487553226, 0.9077828376860282, 0.6203465827944645, 0.5718458931294086, 0.6145601194871151, 0.8287895345230815, 0.4840671287983366, 0.528045500711047, 0.3019908611605266, 0.5920651127359726, 0.6633952407714433, 0.2717300292368672, 0.9410284529321048, 0.40805599852975677, 0.3938763723945892, 0.9037975178233014, 1.0, 0.5388776790134677, 1.0, 1.0, 0.5119625169087315, 0.761337441382025, 0.5245161076572342, 0.4454867145122507, 0.4832349190296145, 1.0, 0.6884018327340987, 0.3871598593733973, 0.3865059449003242, 0.26122092667909996, 0.7855607355995992, 0.6516425101448562, 0.35665517685207765, 0.31580867556595144, 0.8365977186538159, 0.6931063296916856, 0.6222593658724235, 0.34012076087652204, 0.597117429789852, 0.9828997190250446, 0.9245065488241969, 0.9999999999999999, 0.5355011576186107, 0.9613491931397905, 0.30799893344724044, 0.3079065101263343, 0.7768325679696754, 0.6292300682309973, 0.3717428479750089, 0.3002061501054825, 1.0, 0.5916759340195396, 1.0, 0.5232801401529532, 0.6799722663465934, 0.3244454108881641, 0.3297829488833795, 0.4660069727460682, 0.9635242769667345, 0.4149253599615462, 0.6818423921591732, 0.5005586085332286, 0.9246463436000913, 0.4710487008377107, 0.6789851756097529, 0.7769049960546438, 0.33626663293279196, 0.7657830326946345, 0.5247364295376176, 0.29218225593917235, 0.5982371312851179, 0.2825362733693641, 0.4601433858987366, 0.7321159285921941, 0.7876094130670995, 0.9999999999999999, 0.9001269342348065, 0.5643747028510756, 0.7543842879730561, 0.27282687921370147, 0.46269368290715934, 0.8525119254935838, 0.29154429851624847, 0.1161422805173723, 0.6562489734034517, 0.3240681846384767, 0.9479385606640501, 0.3443287124336051, 0.4723741732929433, 0.6752244409691538, 0.5257340972627579, 0.3959179240012579, 0.7156340067832128, 0.5815650375368426, 0.9301649048604338, 0.37009805229260256, 0.3021194851381369, 0.6250199963416462, 0.8032147397667786, 0.5408567025051199, 0.5115622600724877, 0.4647381777537374, 1.0, 0.6497068798606109, 0.6308866558106152, 0.8179592396604881, 0.9582907081902652, 0.709819833549875, 0.4759120856110369, 1.0, 0.9999999999999999, 0.6001652248810057, 0.6436732277002749, 0.49681115932125774, 0.9382285239917659, 1.0, 0.587681850630745, 0.7888598353808676, 0.7150504777222365, 0.3303094674571879, 0.4823373083417599, 0.18485675924952247, 0.9312665112083346, 0.9226877553289299, 0.5117452750618938, 0.65025953301384, 0.2849671475257664, 0.4256005344639421, 0.6237606745961042, 0.4004271686784926, 0.6840336186018563, 0.7355885109020938, 0.6033119411934473, 0.5440462097668168, 0.6061442655352145, 0.45492582011326077, 0.534272134395567, 0.6705550951622102, 0.24142709993973696, 0.8701092014902211, 0.6639111253182824, 0.6562934695976097, 0.6096278271054203, 0.3119267619039908, 0.32509226481750836, 0.6901453653872811, 0.897108699389913, 0.8098856073025985, 0.4048330922320003, 0.31038287943024184, 0.8569562337423862, 0.2622487276188536, 0.7068226063177803, 0.2693915446946901, 0.8025894340462817, 0.8160039970034775, 1.0, 0.41923647935172204, 0.9936598666956731, 0.5095085161402483, 0.8539189335290018, 0.7435760908492823, 0.6452450962670923, 0.2604390372412532, 0.5764258660575072, 1.0, 0.4618715523367149, 0.7783458069547369, 0.4756656993517585, 1.0, 0.9246879407737996, 0.28867651206318046, 0.9447984818058941, 0.5911853428542229, 0.675142836384889, 0.33910450212711973, 0.43142775263620464, 0.4482318863350352, 0.6861781409627846, 0.9315640541142023, 1.0, 1.0, 0.94165587054529, 0.7207221792416874, 0.47555324321101683, 0.360031681791678, 0.5375690500819532, 0.7608102402090884, 0.3493514201113307, 0.6201489244881461, 0.4971765668189029, 0.30486427216258427, 1.0, 0.30182288173199856, 0.4005335194670301, 0.3023454716745363, 0.2582113271409211, 1.0, 0.7235086590193208, 0.259939650705488, 0.3926648744291177, 1.0, 0.25047768362308626, 0.476991437272595, 0.521688569885495, 0.6847695138066671, 0.5694378031289298, 0.5053506805477829, 0.7192563527214755, 0.3377759883175256, 0.9401709299052567, 0.42533286886280486, 0.3248401000939976, 1.0, 0.35435801338098605, 0.3063192690129409, 0.6057248203153163, 0.9546280021385912, 0.3823472271689831, 0.9853140648933334, 0.8350564271062805, 0.2667556826563614, 0.7881469762006607, 0.29148115795757684, 0.8400080624127885, 0.9131222634141362, 0.5263984785296858, 0.6651618994668821, 1.0, 0.8934031464060519, 0.3565293574876279, 0.9999999999999999, 0.4796252855677836, 0.264440726006657, 0.5971297331572046, 0.31248146579271413, 0.8530365619109169, 0.4363614905033652, 1.0, 0.4748485949551296, 0.8040049005900616, 0.3930184180999957, 0.6188572556680967, 0.30373764724900304, 0.5840815104273308, 1.0, 0.6396032008936419, 0.6417368337983786, 1.0, 0.252668600178402, 0.6193468906279678, 0.5649700084172461, 0.3979924494680678, 0.35830604593451, 0.6892020027166835, 0.37023750331147076, 0.4372363945145301, 0.5916136810967754, 0.8773890959484252, 0.42399570101959305, 0.4577618191608436, 0.39915576019146554, 0.6050321687886541, 0.46573051065875765, 0.45802066267638863, 0.6431389525820794, 0.3226566691103885, 0.7368320693541283, 0.8687696402043263, 0.25745351856915527, 0.3887394781672973, 0.26206870882522254, 0.5717442203332083, 0.580523625669574, 0.7216201744650723, 0.30846500185735987, 0.5624080576531992, 0.6733443892930497, 0.7443481744529276, 0.39270141124067137, 0.8447454167102542, 0.2920682859588468, 0.5307222784679029, 1.0, 0.8338004163345377, 0.525330234091529, 0.9844972068755765, 0.3871211108577018, 1.0, 0.96066601564995, 0.5729907171267085, 0.9675821950492496, 0.8481531429916626, 0.9999999999999999, 0.7763880463896716, 0.8701728294291401, 1.0, 1.0, 0.2620019649708491, 0.6893414036275599, 0.6226622263225624, 0.6558492598977669, 0.8693810135665435, 0.9616130915578651, 0.37930325280590643, 1.0, 0.9213990890089069, 0.5311369799416545, 0.7946336885971379, 0.9187052754109976, 0.9729265262948814, 0.28757710492059424, 0.2799776003965103, 0.4943666390701472, 0.2902415029547936, 0.6130617256671121, 0.9569342536384492, 0.36573418091371834, 0.5151027489635106, 0.3711926992925339, 0.8265041888057472, 0.8711989155681541, 0.26638065409653683, 0.602416249699157, 0.6466757166086051, 0.6869570959923448, 0.933343172675266, 0.5487914095314416, 0.4801225813581666, 0.7462326501829734, 0.6824062117185221, 0.8426490119727528, 0.25429037885178524, 0.2666761369603077, 0.567934992398271, 0.7324004284235474, 0.3728118355292142, 0.5956068278236786, 0.5601790155751898, 0.5237793396323196, 0.9014696980491468, 0.3669263636585938, 1.0, 1.0, 0.8469964163576164, 0.3197456008459634, 0.6101546861765202, 0.795967607052283, 0.1459631645904522, 0.5066807948390265, 0.32793372231199736, 0.3619680783182243, 0.7593226021470355, 0.9022486187099299, 0.9560707063589859, 0.5521318441185364, 0.4267752277451563, 0.6884295368562907, 0.3908620690286871, 0.7627755003817481, 0.5528182694087352, 0.6628180766433622, 0.3436532060135707, 0.6890523459698933, 0.454235115753991, 0.9673187897912604, 0.5483193092160413, 0.3188771931224632, 0.731986371680395, 0.4745962247245148, 0.6916189910884466, 1.0, 0.21764505314015709, 0.8457034776232927, 0.33839532807430034, 0.6197370106730804, 0.5198941149219323, 0.9185573227536804, 0.771076207661388, 0.3314623530222701, 0.5022484042378285, 0.27924642258507754, 0.34465587607879683, 0.4599681423428597, 0.4398592848874709, 0.5572569697697632, 0.2714407146118095, 0.9999999999999999, 1.0, 0.3530736372438787, 0.2524188463069196, 0.47907463771865294, 0.4736124378752765, 0.37196443071576957, 0.6173007019791437, 0.3786432796181109, 0.5345124256882472, 0.41738612052990115, 0.3587837000504286, 0.5078336426237405, 0.35996357178756705, 0.3392725438906634, 0.5115187039292588, 0.9533640638562828, 0.6818752230587695, 0.649888913318365, 0.5259450491762055, 0.9662007334889993, 0.5118191982704826, 0.7005395548747968, 0.6979182939464459, 0.5589875232081593, 0.21801593735190072, 0.8533955208041076, 0.4645947060642481, 0.5763306501547387, 0.38226799611567625, 0.6089868548732689, 0.34683066175747307, 0.3938402950846586, 0.6418452083401203, 0.8389294213001933, 0.49195203646211433, 0.4223524930896001, 0.9449331541538666, 0.8222675473725303, 0.9999999999999999, 0.7241809393586591, 0.9720751070542792, 0.5719068703503046, 0.8653526879724195, 0.736256897161814, 0.2912757011591064, 0.6211509155566186, 0.5089917708413657, 0.32524059413213474, 0.7629040242348248, 0.8026945227856959, 0.34083731188745275, 0.45451676093570975, 0.29734364156820403, 0.834243557607054, 0.4619529169978442, 0.6124774701495476, 0.49082068853093513, 0.37624752587278676, 0.6636897504604966, 0.6881723017397747, 0.971674701438198, 0.8109461350055639, 0.724833180995697, 0.7190775568629929, 0.4192539842987171, 0.2613412738267231, 0.6273760806162776, 0.9498326181335748, 0.3549358507517958, 0.6850006180482371, 0.3393170131238329, 0.9735513812279997, 0.8268943689023313, 0.5044913887188193, 0.37468431847196104, 0.9797661990473747, 0.8143349610797487, 0.2578744894612354, 0.7628863080887179, 0.3556857492040174, 0.3113102739984362, 0.48234900945966425, 0.848253187187421, 0.553704994347155, 0.8470366854960708, 0.9681399923268599, 1.0, 0.27033887162375575, 0.6001445566013698, 0.8411645981505789, 0.7082370001715657, 0.7520041385280705, 0.2839283717189735, 0.33577835064041733, 0.40132007946391207, 0.30372822635538443, 0.4440259163499051, 0.47327501432046654, 0.25588507074098815, 0.8363490778364355, 0.4227938589182584, 0.5356371548638525, 0.6262392572965, 0.49887348836512835, 0.4337443948965278, 0.5088840587507965, 0.7475377966636774, 0.9084633394631072, 0.8980677097667479, 1.0, 0.4971919835249001, 0.36215867844515115, 0.5030439358256541, 0.9999999999999999, 0.6128340129322238, 0.957853827520242, 0.26173711293857144, 0.9445601800195532, 0.6463062168253, 0.9359790834644321, 0.48481556150407, 0.4606726776681997, 0.6523234407592361, 0.9431899510005763, 0.7716997840253819, 0.7037752041364878, 0.5209368152030175, 0.4686815261968328, 0.631477970451007, 0.6861367630758495, 0.9153324016497332, 0.40719894945531665, 0.6671296984803738, 0.45684185233732216, 0.5571227393940452, 0.3277307156780763, 0.40664590299519765, 0.2685391927099636, 0.37900078258929576, 0.9106459733208746, 0.5809471676873356, 0.696903479360243, 0.9174640138932246, 0.514569020188006, 0.525771165640361, 1.0, 0.25989533628093325, 0.7811455881507257, 0.37154849348521657, 0.8255691233488658, 0.34372999822715966, 0.8150702376429599, 0.6352278571558805, 0.8177211895058523, 0.368081040181522, 0.28127649415584643, 0.2797171091416061, 0.6512589663473924, 0.44130082628674205, 0.67894227049179, 0.663315640939975, 0.8669993705189127, 0.3837779650490343, 0.610741612892118, 0.4325822240062867, 0.9500731897169573, 0.8298474562556176, 0.7044212228512335, 0.47191315927145916, 0.2500486066687155, 0.8056157758952187, 0.7607045387568007, 0.583145403367163, 0.5094346200273572, 0.6987329453869764, 0.5094115849729446, 0.8124403673268045, 1.0, 0.8766623126853672, 0.36229685050944144, 0.2814831423962957, 0.7400462241954264, 0.7925874112259261, 0.48657337318632843, 0.38427959671173495, 0.4421206466111112, 0.571162424441088, 0.450526334456363, 0.9034126825248033, 0.4413704155428708, 0.5120738051974955, 0.7293093387206996, 0.3098471684828302, 0.5803142665004766, 0.3337301388087098, 1.0, 0.25747612749506243, 0.3074992708743665, 1.0, 0.5901258524797363, 0.4105136591687168, 0.982655463899435, 0.4425131750701062, 0.6335508324331612, 0.9999999999999999, 0.629231912878056, 0.3576480481693306, 0.44529106472878766, 0.9347796875685308, 0.34411663066839193, 0.7977463652562109, 0.9943702129577762, 0.34837546804964287, 0.27965932932363335, 0.40970313100628736, 0.40809653305001736, 0.3277257124827053, 0.6522010767748746, 0.294188391303247, 0.6079975216577458, 0.5306700711172905, 0.3918712505166603, 0.5719585597208344, 0.3410013923755241, 0.40354432676449276];
G = SparseMatrixCSC(Gm, Gn, Gcolptr, Growval, Gnzval);


text:
#= Fornberg algorithm
This implements the Fornberg (1988) algorithm (https://doi.org/10.1090/S0025-5718-1988-0935077-0)
and hermite-based finite difference Fornberg(2020) algorithm (https://doi.org/10.1093/imanum/draa006)
to obtain Finite Difference weights over arbitrary points to arbitrary order.
Inputs:
        order: The derivative order for which we need the coefficients
        x0   : The point in the array 'x' for which we need the coefficients
        x    : A dummy array with relative coordinates, e.g., central differences
               need coordinates centred at 0 while those at boundaries need
               coordinates starting from 0 to the end point
        dfdx : optional argument to consider weights of the first-derivative of function or not
                if    dfdx == false (default kwarg), implies Fornberg(1988)
                      dfdx == true,                  implies Fornberg(2020)
    Outputs:
        if dfdx == false (default kwarg),   _C : weights to approximate derivative of required order using function values only.
                                 else,   _D,_E : weights to approximate derivative of required order using function and its first-
                                                 derivative values respectively.
=#

function calculate_weights(order::Int, x0::T, x::AbstractVector; dfdx::Bool = false) where T<:Real
    N = length(x)
    @assert order < N "Not enough points for the requested order."
    M = order
    c1 = one(T)
    c4 = x[1] - x0
    C = zeros(T, N, M+1)
    C[1,1] = 1
    @inbounds for i in 1 : N-1
        i1 = i + 1
        mn = min(i, M)
        c2 = one(T)
        c5 = c4
        c4 = x[i1] - x0
        for j in 0 : i-1
            j1 = j + 1
            c3 = x[i1] - x[j1]
            c2 *= c3
            if j == i-1
                for s in mn : -1 : 1
                    s1 = s + 1
                    C[i1,s1] = c1*(s*C[i,s] - c5*C[i,s1]) / c2
                end
                C[i1,1] = -c1*c5*C[i,1] / c2
           end
            for s in mn : -1 : 1
                s1 = s + 1
                C[j1,s1] = (c4*C[j1,s1] - s*C[j1,s]) / c3
            end
            C[j1,1] = c4 * C[j1,1] / c3
        end
        c1 = c2
    end
    #=
        This is to fix the problem of numerical instability which occurs when the sum of the stencil_coefficients is not
        exactly 0.
        https://scicomp.stackexchange.com/questions/11249/numerical-derivative-and-finite-difference-coefficients-any-update-of-the-fornb
        Stack Overflow answer on this issue.
        http://epubs.siam.org/doi/pdf/10.1137/S0036144596322507 - Modified Fornberg Algorithm
    =#
    _C = C[:,end]
    if order != 0
        _C[div(N,2)+1] -= sum(_C)
    end
    if dfdx == false
        return _C
    else
        A = x .- x';
        s = sum(1 ./ (A + I(N)), dims = 1) .- 1;
        cp = factorial.(0:M);
        cc = C./cp'
        c̃ = zeros(N, M+2);
        for k in 1:M+1
           c̃[:,k+1] = sum(cc[:,1:k].*cc[:,k:-1:1], dims = 2);
        end
        E = c̃[:,1:M+1] - (x .- x0).*c̃[:,2:M+2];
        D = c̃[:,2:M+2] + 2*E.*s';
        D = D.*cp';
        E = E.*cp';

        _D = D[:,end];   _E = E[:,end]
        return _D, _E
    end
end

ext: jl


text:
# This file was generated, do not modify it. # hide
sum(coefs .≈ 0) / length(coefs)


text:
module ArrayInterface

using Requires

function ismutable end

"""
    ismutable(x::DataType)

Query whether a type is mutable or not, see
https://github.com/JuliaDiffEq/RecursiveArrayTools.jl/issues/19.
"""
Base.@pure ismutable(x::DataType) = x.mutable
ismutable(x) = ismutable(typeof(x))

ismutable(::Type{Array}) = true
ismutable(::Type{<:Number}) = false

function __init__()

  @require StaticArrays="90137ffa-7385-5640-81b9-e52037218182" begin
    ismutable(::Type{StaticArrays.StaticArray}) = false
    ismutable(::Type{StaticArrays.MArray}) = true
  end

  @require LabelledArrays="2ee39098-c373-598a-b85f-a56591580800" begin
    ismutable(::Type{LArray{T,N,Syms}}) where {T,N,Syms} = ismutable(T)
  end

  @require Flux="587475ba-b771-5e3f-ad9e-33799f191a9c" begin
    ismutable(x::Flux.Tracker.TrackedArray) = false
  end
end

end


text:
function LCR(model; args...)
    Linesearch(model, CRlin, profil=true)
end


text:
# Autogenerated wrapper script for coreutils_jll for armv6l-linux-gnueabihf
export b2sum, base32, base64, basename_bin, basenc, cat_bin, chcon, chgrp, chmod_bin, chown_bin, chroot, cksum, comm, cp_bin, csplit, cut, date, dd, df, dir, dircolors, dirname_bin, du, echo, env, expand, expr, factor, false_bin, fmt, fold, groups, head, hostid, id, install, join_bin, kill_bin, left_bracket, libstdbuf, link, ln, logname, ls, md5sum, mkdir_bin, mkfifo, mknod, mktemp_bin, mv_bin, nice, nl, nohup, nproc, numfmt, od, paste, pathchk, pinky, pr, printenv, printf, ptx, pwd_bin, readlink_bin, realpath_bin, rm_bin, rmdir, runcon, seq, sha1sum, sha224sum, sha256sum, sha384sum, sha512sum, shred, shuf, sleep_bin, sort_bin, split_bin, stat_bin, stdbuf, stty, sum_bin, sync, tac, tail_bin, tee, test, timeout, touch_bin, tr, true_bin, truncate_bin, tsort, tty, uname, unexpand, uniq, unlink_bin, uptime, users, vdir, wc, who, whoami, yes

JLLWrappers.@generate_wrapper_header("coreutils")
JLLWrappers.@declare_library_product(libstdbuf, "libstdbuf.so")
JLLWrappers.@declare_executable_product(b2sum)
JLLWrappers.@declare_executable_product(base32)
JLLWrappers.@declare_executable_product(base64)
JLLWrappers.@declare_executable_product(basename_bin)
JLLWrappers.@declare_executable_product(basenc)
JLLWrappers.@declare_executable_product(cat_bin)
JLLWrappers.@declare_executable_product(chcon)
JLLWrappers.@declare_executable_product(chgrp)
JLLWrappers.@declare_executable_product(chmod_bin)
JLLWrappers.@declare_executable_product(chown_bin)
JLLWrappers.@declare_executable_product(chroot)
JLLWrappers.@declare_executable_product(cksum)
JLLWrappers.@declare_executable_product(comm)
JLLWrappers.@declare_executable_product(cp_bin)
JLLWrappers.@declare_executable_product(csplit)
JLLWrappers.@declare_executable_product(cut)
JLLWrappers.@declare_executable_product(date)
JLLWrappers.@declare_executable_product(dd)
JLLWrappers.@declare_executable_product(df)
JLLWrappers.@declare_executable_product(dir)
JLLWrappers.@declare_executable_product(dircolors)
JLLWrappers.@declare_executable_product(dirname_bin)
JLLWrappers.@declare_executable_product(du)
JLLWrappers.@declare_executable_product(echo)
JLLWrappers.@declare_executable_product(env)
JLLWrappers.@declare_executable_product(expand)
JLLWrappers.@declare_executable_product(expr)
JLLWrappers.@declare_executable_product(factor)
JLLWrappers.@declare_executable_product(false_bin)
JLLWrappers.@declare_executable_product(fmt)
JLLWrappers.@declare_executable_product(fold)
JLLWrappers.@declare_executable_product(groups)
JLLWrappers.@declare_executable_product(head)
JLLWrappers.@declare_executable_product(hostid)
JLLWrappers.@declare_executable_product(id)
JLLWrappers.@declare_executable_product(install)
JLLWrappers.@declare_executable_product(join_bin)
JLLWrappers.@declare_executable_product(kill_bin)
JLLWrappers.@declare_executable_product(left_bracket)
JLLWrappers.@declare_executable_product(link)
JLLWrappers.@declare_executable_product(ln)
JLLWrappers.@declare_executable_product(logname)
JLLWrappers.@declare_executable_product(ls)
JLLWrappers.@declare_executable_product(md5sum)
JLLWrappers.@declare_executable_product(mkdir_bin)
JLLWrappers.@declare_executable_product(mkfifo)
JLLWrappers.@declare_executable_product(mknod)
JLLWrappers.@declare_executable_product(mktemp_bin)
JLLWrappers.@declare_executable_product(mv_bin)
JLLWrappers.@declare_executable_product(nice)
JLLWrappers.@declare_executable_product(nl)
JLLWrappers.@declare_executable_product(nohup)
JLLWrappers.@declare_executable_product(nproc)
JLLWrappers.@declare_executable_product(numfmt)
JLLWrappers.@declare_executable_product(od)
JLLWrappers.@declare_executable_product(paste)
JLLWrappers.@declare_executable_product(pathchk)
JLLWrappers.@declare_executable_product(pinky)
JLLWrappers.@declare_executable_product(pr)
JLLWrappers.@declare_executable_product(printenv)
JLLWrappers.@declare_executable_product(printf)
JLLWrappers.@declare_executable_product(ptx)
JLLWrappers.@declare_executable_product(pwd_bin)
JLLWrappers.@declare_executable_product(readlink_bin)
JLLWrappers.@declare_executable_product(realpath_bin)
JLLWrappers.@declare_executable_product(rm_bin)
JLLWrappers.@declare_executable_product(rmdir)
JLLWrappers.@declare_executable_product(runcon)
JLLWrappers.@declare_executable_product(seq)
JLLWrappers.@declare_executable_product(sha1sum)
JLLWrappers.@declare_executable_product(sha224sum)
JLLWrappers.@declare_executable_product(sha256sum)
JLLWrappers.@declare_executable_product(sha384sum)
JLLWrappers.@declare_executable_product(sha512sum)
JLLWrappers.@declare_executable_product(shred)
JLLWrappers.@declare_executable_product(shuf)
JLLWrappers.@declare_executable_product(sleep_bin)
JLLWrappers.@declare_executable_product(sort_bin)
JLLWrappers.@declare_executable_product(split_bin)
JLLWrappers.@declare_executable_product(stat_bin)
JLLWrappers.@declare_executable_product(stdbuf)
JLLWrappers.@declare_executable_product(stty)
JLLWrappers.@declare_executable_product(sum_bin)
JLLWrappers.@declare_executable_product(sync)
JLLWrappers.@declare_executable_product(tac)
JLLWrappers.@declare_executable_product(tail_bin)
JLLWrappers.@declare_executable_product(tee)
JLLWrappers.@declare_executable_product(test)
JLLWrappers.@declare_executable_product(timeout)
JLLWrappers.@declare_executable_product(touch_bin)
JLLWrappers.@declare_executable_product(tr)
JLLWrappers.@declare_executable_product(true_bin)
JLLWrappers.@declare_executable_product(truncate_bin)
JLLWrappers.@declare_executable_product(tsort)
JLLWrappers.@declare_executable_product(tty)
JLLWrappers.@declare_executable_product(uname)
JLLWrappers.@declare_executable_product(unexpand)
JLLWrappers.@declare_executable_product(uniq)
JLLWrappers.@declare_executable_product(unlink_bin)
JLLWrappers.@declare_executable_product(uptime)
JLLWrappers.@declare_executable_product(users)
JLLWrappers.@declare_executable_product(vdir)
JLLWrappers.@declare_executable_product(wc)
JLLWrappers.@declare_executable_product(who)
JLLWrappers.@declare_executable_product(whoami)
JLLWrappers.@declare_executable_product(yes)
function __init__()
    JLLWrappers.@generate_init_header()
    JLLWrappers.@init_library_product(
        libstdbuf,
        "libexec/coreutils/libstdbuf.so",
        RTLD_LAZY | RTLD_DEEPBIND,
    )

    JLLWrappers.@init_executable_product(
        b2sum,
        "bin/b2sum",
    )

    JLLWrappers.@init_executable_product(
        base32,
        "bin/base32",
    )

    JLLWrappers.@init_executable_product(
        base64,
        "bin/base64",
    )

    JLLWrappers.@init_executable_product(
        basename_bin,
        "bin/basename",
    )

    JLLWrappers.@init_executable_product(
        basenc,
        "bin/basenc",
    )

    JLLWrappers.@init_executable_product(
        cat_bin,
        "bin/cat",
    )

    JLLWrappers.@init_executable_product(
        chcon,
        "bin/chcon",
    )

    JLLWrappers.@init_executable_product(
        chgrp,
        "bin/chgrp",
    )

    JLLWrappers.@init_executable_product(
        chmod_bin,
        "bin/chmod",
    )

    JLLWrappers.@init_executable_product(
        chown_bin,
        "bin/chown",
    )

    JLLWrappers.@init_executable_product(
        chroot,
        "bin/chroot",
    )

    JLLWrappers.@init_executable_product(
        cksum,
        "bin/cksum",
    )

    JLLWrappers.@init_executable_product(
        comm,
        "bin/comm",
    )

    JLLWrappers.@init_executable_product(
        cp_bin,
        "bin/cp",
    )

    JLLWrappers.@init_executable_product(
        csplit,
        "bin/csplit",
    )

    JLLWrappers.@init_executable_product(
        cut,
        "bin/cut",
    )

    JLLWrappers.@init_executable_product(
        date,
        "bin/date",
    )

    JLLWrappers.@init_executable_product(
        dd,
        "bin/dd",
    )

    JLLWrappers.@init_executable_product(
        df,
        "bin/df",
    )

    JLLWrappers.@init_executable_product(
        dir,
        "bin/dir",
    )

    JLLWrappers.@init_executable_product(
        dircolors,
        "bin/dircolors",
    )

    JLLWrappers.@init_executable_product(
        dirname_bin,
        "bin/dirname",
    )

    JLLWrappers.@init_executable_product(
        du,
        "bin/du",
    )

    JLLWrappers.@init_executable_product(
        echo,
        "bin/echo",
    )

    JLLWrappers.@init_executable_product(
        env,
        "bin/env",
    )

    JLLWrappers.@init_executable_product(
        expand,
        "bin/expand",
    )

    JLLWrappers.@init_executable_product(
        expr,
        "bin/expr",
    )

    JLLWrappers.@init_executable_product(
        factor,
        "bin/factor",
    )

    JLLWrappers.@init_executable_product(
        false_bin,
        "bin/false",
    )

    JLLWrappers.@init_executable_product(
        fmt,
        "bin/fmt",
    )

    JLLWrappers.@init_executable_product(
        fold,
        "bin/fold",
    )

    JLLWrappers.@init_executable_product(
        groups,
        "bin/groups",
    )

    JLLWrappers.@init_executable_product(
        head,
        "bin/head",
    )

    JLLWrappers.@init_executable_product(
        hostid,
        "bin/hostid",
    )

    JLLWrappers.@init_executable_product(
        id,
        "bin/id",
    )

    JLLWrappers.@init_executable_product(
        install,
        "bin/install",
    )

    JLLWrappers.@init_executable_product(
        join_bin,
        "bin/join",
    )

    JLLWrappers.@init_executable_product(
        kill_bin,
        "bin/kill",
    )

    JLLWrappers.@init_executable_product(
        left_bracket,
        "bin/[",
    )

    JLLWrappers.@init_executable_product(
        link,
        "bin/link",
    )

    JLLWrappers.@init_executable_product(
        ln,
        "bin/ln",
    )

    JLLWrappers.@init_executable_product(
        logname,
        "bin/logname",
    )

    JLLWrappers.@init_executable_product(
        ls,
        "bin/ls",
    )

    JLLWrappers.@init_executable_product(
        md5sum,
        "bin/md5sum",
    )

    JLLWrappers.@init_executable_product(
        mkdir_bin,
        "bin/mkdir",
    )

    JLLWrappers.@init_executable_product(
        mkfifo,
        "bin/mkfifo",
    )

    JLLWrappers.@init_executable_product(
        mknod,
        "bin/mknod",
    )

    JLLWrappers.@init_executable_product(
        mktemp_bin,
        "bin/mktemp",
    )

    JLLWrappers.@init_executable_product(
        mv_bin,
        "bin/mv",
    )

    JLLWrappers.@init_executable_product(
        nice,
        "bin/nice",
    )

    JLLWrappers.@init_executable_product(
        nl,
        "bin/nl",
    )

    JLLWrappers.@init_executable_product(
        nohup,
        "bin/nohup",
    )

    JLLWrappers.@init_executable_product(
        nproc,
        "bin/nproc",
    )

    JLLWrappers.@init_executable_product(
        numfmt,
        "bin/numfmt",
    )

    JLLWrappers.@init_executable_product(
        od,
        "bin/od",
    )

    JLLWrappers.@init_executable_product(
        paste,
        "bin/paste",
    )

    JLLWrappers.@init_executable_product(
        pathchk,
        "bin/pathchk",
    )

    JLLWrappers.@init_executable_product(
        pinky,
        "bin/pinky",
    )

    JLLWrappers.@init_executable_product(
        pr,
        "bin/pr",
    )

    JLLWrappers.@init_executable_product(
        printenv,
        "bin/printenv",
    )

    JLLWrappers.@init_executable_product(
        printf,
        "bin/printf",
    )

    JLLWrappers.@init_executable_product(
        ptx,
        "bin/ptx",
    )

    JLLWrappers.@init_executable_product(
        pwd_bin,
        "bin/pwd",
    )

    JLLWrappers.@init_executable_product(
        readlink_bin,
        "bin/readlink",
    )

    JLLWrappers.@init_executable_product(
        realpath_bin,
        "bin/realpath",
    )

    JLLWrappers.@init_executable_product(
        rm_bin,
        "bin/rm",
    )

    JLLWrappers.@init_executable_product(
        rmdir,
        "bin/rmdir",
    )

    JLLWrappers.@init_executable_product(
        runcon,
        "bin/runcon",
    )

    JLLWrappers.@init_executable_product(
        seq,
        "bin/seq",
    )

    JLLWrappers.@init_executable_product(
        sha1sum,
        "bin/sha1sum",
    )

    JLLWrappers.@init_executable_product(
        sha224sum,
        "bin/sha224sum",
    )

    JLLWrappers.@init_executable_product(
        sha256sum,
        "bin/sha256sum",
    )

    JLLWrappers.@init_executable_product(
        sha384sum,
        "bin/sha384sum",
    )

    JLLWrappers.@init_executable_product(
        sha512sum,
        "bin/sha512sum",
    )

    JLLWrappers.@init_executable_product(
        shred,
        "bin/shred",
    )

    JLLWrappers.@init_executable_product(
        shuf,
        "bin/shuf",
    )

    JLLWrappers.@init_executable_product(
        sleep_bin,
        "bin/sleep",
    )

    JLLWrappers.@init_executable_product(
        sort_bin,
        "bin/sort",
    )

    JLLWrappers.@init_executable_product(
        split_bin,
        "bin/split",
    )

    JLLWrappers.@init_executable_product(
        stat_bin,
        "bin/stat",
    )

    JLLWrappers.@init_executable_product(
        stdbuf,
        "bin/stdbuf",
    )

    JLLWrappers.@init_executable_product(
        stty,
        "bin/stty",
    )

    JLLWrappers.@init_executable_product(
        sum_bin,
        "bin/sum",
    )

    JLLWrappers.@init_executable_product(
        sync,
        "bin/sync",
    )

    JLLWrappers.@init_executable_product(
        tac,
        "bin/tac",
    )

    JLLWrappers.@init_executable_product(
        tail_bin,
        "bin/tail",
    )

    JLLWrappers.@init_executable_product(
        tee,
        "bin/tee",
    )

    JLLWrappers.@init_executable_product(
        test,
        "bin/test",
    )

    JLLWrappers.@init_executable_product(
        timeout,
        "bin/timeout",
    )

    JLLWrappers.@init_executable_product(
        touch_bin,
        "bin/touch",
    )

    JLLWrappers.@init_executable_product(
        tr,
        "bin/tr",
    )

    JLLWrappers.@init_executable_product(
        true_bin,
        "bin/true",
    )

    JLLWrappers.@init_executable_product(
        truncate_bin,
        "bin/truncate",
    )

    JLLWrappers.@init_executable_product(
        tsort,
        "bin/tsort",
    )

    JLLWrappers.@init_executable_product(
        tty,
        "bin/tty",
    )

    JLLWrappers.@init_executable_product(
        uname,
        "bin/uname",
    )

    JLLWrappers.@init_executable_product(
        unexpand,
        "bin/unexpand",
    )

    JLLWrappers.@init_executable_product(
        uniq,
        "bin/uniq",
    )

    JLLWrappers.@init_executable_product(
        unlink_bin,
        "bin/unlink",
    )

    JLLWrappers.@init_executable_product(
        uptime,
        "bin/uptime",
    )

    JLLWrappers.@init_executable_product(
        users,
        "bin/users",
    )

    JLLWrappers.@init_executable_product(
        vdir,
        "bin/vdir",
    )

    JLLWrappers.@init_executable_product(
        wc,
        "bin/wc",
    )

    JLLWrappers.@init_executable_product(
        who,
        "bin/who",
    )

    JLLWrappers.@init_executable_product(
        whoami,
        "bin/whoami",
    )

    JLLWrappers.@init_executable_product(
        yes,
        "bin/yes",
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()


text:
Base.@kwdef struct ParallelQuickSortAlg{Alg,SmallSize,BaseSize} <: ParallelSortAlgorithm
    smallsort::Alg = Base.Sort.DEFAULT_UNSTABLE
    smallsize::SmallSize = nothing  # lazily determined
    basesize::BaseSize = 10_000
end
# `basesize` is tuned using `Float64`.  Make it `eltype`-aware?
#
# Note: Small `smallsize` (e.g., 32) is beneficial on random data. However, it
# is slower on sorted data with small number of worker threads.
# TODO: Implement some heuristics for `smallsize`.

function Base.sort!(
    v::AbstractVector,
    lo::Integer,
    hi::Integer,
    a::ParallelQuickSortAlg,
    o::Ordering,
)
    if a.basesize === nothing
        a = @set a.basesize = default_basesize(hi - lo + 1)
    end
    if a.smallsize === nothing
        a = @set a.smallsize = a.basesize
    end
    ys = view(v, lo:hi)
    _quicksort!(
        similar(ys),
        ys,
        a,
        o,
        Vector{Int8}(undef, length(ys)),
        false,  # ys_is_result
        true,   # mutable_xs
    )
    return v
end

function _quicksort!(
    ys,
    xs,
    alg,
    order,
    cs = Vector{Int8}(undef, length(ys)),
    ys_is_result = true,
    mutable_xs = false,
)
    @check length(ys) == length(xs)
    if length(ys) <= max(8, alg.basesize)
        return _quicksort_serial!(ys, xs, alg, order, cs, ys_is_result, mutable_xs)
    end
    pivot = _median(
        order,
        (
            xs[1],
            xs[end÷8],
            xs[end÷4],
            xs[3*(end÷8)],
            xs[end÷2],
            xs[5*(end÷8)],
            xs[3*(end÷4)],
            xs[7*(end÷8)],
            xs[end],
        ),
    )
    chunksize = alg.basesize

    # TODO: Calculate extrema during the first pass if it's possible
    # to use counting sort.
    # TODO: When recursing, fuse copying _from_ `ys` to `xs` with the
    # first pass.

    # Compute sizes of each partition for each chunks.
    xs_chunk_list = _partition(xs, chunksize)
    cs_chunk_list = _partition(cs, chunksize)
    nchunks = cld(length(xs), chunksize)
    nbelows = Vector{Int}(undef, nchunks)
    nequals = Vector{Int}(undef, nchunks)
    naboves = Vector{Int}(undef, nchunks)
    @DBG begin
        VERSION >= v"1.4" &&
            @check length(xs_chunk_list) == length(cs_chunk_list) == nchunks
        fill!(nbelows, -1)
        fill!(nequals, -1)
        fill!(naboves, -1)
    end
    @sync for (nb, ne, na, xs_chunk, cs_chunk) in zip(
        referenceable(nbelows),
        referenceable(nequals),
        referenceable(naboves),
        xs_chunk_list,
        cs_chunk_list,
    )
        @spawn partition_sizes!(nb, ne, na, xs_chunk, cs_chunk, pivot, order)
    end
    @DBG begin
        @check all(>=(0), nbelows)
        @check all(>=(0), nequals)
        @check all(>=(0), naboves)
    end

    below_offsets = nbelows
    equal_offsets = nequals
    above_offsets = naboves
    acc = exclusive_cumsum!(below_offsets)
    acc = exclusive_cumsum!(equal_offsets, acc)
    acc = exclusive_cumsum!(above_offsets, acc)
    @check acc == length(xs)

    @inline function singleton_chunkid(i)
        nb = @inbounds get(below_offsets, i + 1, equal_offsets[1]) - below_offsets[i]
        ne = @inbounds get(equal_offsets, i + 1, above_offsets[1]) - equal_offsets[i]
        na = @inbounds get(above_offsets, i + 1, length(ys)) - above_offsets[i]
        if (nb > 0) + (ne > 0) + (na > 0) == 1
            return 1 * (nb > 0) + 2 * (ne > 0) + 3 * (na > 0)
        else
            return 0
        end
    end

    @sync begin
        for (i, (xs_chunk, cs_chunk)) in enumerate(zip(xs_chunk_list, cs_chunk_list))
            singleton_chunkid(i) > 0 && continue
            @spawn unsafe_quicksort_scatter!(
                ys,
                xs_chunk,
                cs_chunk,
                below_offsets[i],
                equal_offsets[i],
                above_offsets[i],
            )
        end
        for (i, xs_chunk) in enumerate(xs_chunk_list)
            sid = singleton_chunkid(i)
            sid > 0 || continue
            idx = (
                below_offsets[i]+1:get(below_offsets, i + 1, equal_offsets[1]),
                equal_offsets[i]+1:get(equal_offsets, i + 1, above_offsets[1]),
                above_offsets[i]+1:get(above_offsets, i + 1, length(ys)),
            )[sid]
            # There is only one partition. Short-circuit scattering.
            ys_chunk = view(ys, idx)
            copyto!(ys_chunk, xs_chunk)
            # Is it better to multi-thread this?
        end
    end

    partitions = (1:equal_offsets[1], above_offsets[1]+1:length(xs))
    @sync begin
        for idx in partitions
            length(idx) <= alg.smallsize && continue
            ys_new = view(ys, idx)
            xs_new = view(xs, idx)
            cs_new = view(cs, idx)
            @spawn let zs
                if mutable_xs
                    zs = xs_new
                else
                    zs = similar(ys_new)
                end
                _quicksort!(zs, ys_new, alg, order, cs_new, !ys_is_result, true)
            end
        end
        for idx in partitions
            length(idx) <= alg.smallsize || continue
            if ys_is_result
                ys_new = view(ys, idx)
            else
                ys_new = copyto!(view(xs, idx), view(ys, idx))
            end
            sort!(ys_new, alg.smallsort, order)
        end
        if !ys_is_result
            let idx = equal_offsets[1]+1:above_offsets[1]
                copyto!(view(xs, idx), view(ys, idx))
            end
        end
    end

    return ys_is_result ? ys : xs
end

function _quicksort_serial!(
    ys,
    xs,
    alg,
    order,
    cs = Vector{Int8}(undef, length(ys)),
    ys_is_result = true,
    mutable_xs = false,
)
    # @check length(ys) == length(xs)
    if length(ys) <= max(8, alg.smallsize)
        if ys_is_result
            zs = copyto!(ys, xs)
        else
            zs = xs
        end
        return sort!(zs, alg.smallsort, order)
    end
    pivot = _median(
        order,
        (
            xs[1],
            xs[end÷8],
            xs[end÷4],
            xs[3*(end÷8)],
            xs[end÷2],
            xs[5*(end÷8)],
            xs[3*(end÷4)],
            xs[7*(end÷8)],
            xs[end],
        ),
    )

    (nbelows, nequals) = partition_sizes!(xs, cs, pivot, order)
    if nequals == length(xs)
        if ys_is_result
            copyto!(ys, xs)
            return ys
        else
            return xs
        end
    end
    @assert nequals > 0
    below_offset = 0
    equal_offset = nbelows
    above_offset = nbelows + nequals
    unsafe_quicksort_scatter!(ys, xs, cs, below_offset, equal_offset, above_offset)

    below = 1:equal_offset
    above = above_offset+1:length(xs)
    ya = view(ys, above)
    yb = view(ys, below)
    ca = view(cs, above)
    cb = view(cs, below)
    if mutable_xs
        _quicksort_serial!(view(xs, above), ya, alg, order, ca, !ys_is_result, true)
        _quicksort_serial!(view(xs, below), yb, alg, order, cb, !ys_is_result, true)
    else
        let zs = similar(ys)
            _quicksort_serial!(view(zs, above), ya, alg, order, ca, !ys_is_result, true)
            _quicksort_serial!(view(zs, below), yb, alg, order, cb, !ys_is_result, true)
        end
    end
    if !ys_is_result
        let idx = equal_offset+1:above_offset
            copyto!(view(xs, idx), view(ys, idx))
        end
    end

    return ys_is_result ? ys : xs
end

function partition_sizes!(nbelows, nequals, naboves, xs, cs, pivot, order)
    (nb, ne) = partition_sizes!(xs, cs, pivot, order)
    nbelows[] = nb
    nequals[] = ne
    naboves[] = length(xs) - (nb + ne)
    return
end

function partition_sizes!(xs, cs, pivot, order)
    nbelows = 0
    nequals = 0
    @inbounds for i in eachindex(xs, cs)
        x = xs[i]
        b = Base.lt(order, x, pivot)
        a = Base.lt(order, pivot, x)
        cs[i] = ifelse(b, -Int8(1), ifelse(a, Int8(1), Int8(0)))
        nbelows += Int(b)
        nequals += Int(!(a | b))
    end
    return (nbelows, nequals)
end

function unsafe_quicksort_scatter!(
    ys,
    xs_chunk,
    cs_chunk,
    below_offset,
    equal_offset,
    above_offset,
)
    b = below_offset
    e = equal_offset
    a = above_offset
    _foldl((b, a, e), Unroll{4}(eachindex(xs_chunk, cs_chunk))) do (b, a, e), i
        @inbounds x = xs_chunk[i]
        @inbounds c = cs_chunk[i]
        is_equal = c == 0
        is_above = c > 0
        is_below = c < 0
        e += Int(is_equal)
        a += Int(is_above)
        b += Int(is_below)
        @inbounds ys[ifelse(is_equal, e, ifelse(is_above, a, b))] = x
        (b, a, e)
    end
    return
end


text:
include("$(@__DIR__)/../src/FiniteGroups.jl")
using .FiniteGroups, LinearAlgebra, Test

g = pointgroup("Oh")

function get_coeff(g, χ)
    n = order(g)
    ct = Matrix{Int64}(undef, 2n, 2n)
    ct[1:2n,1:n] .= 0
    for i = 1:n
        p = χ[inclass(g, i)]
        ct[i, n+1:2n] .= (1-p)÷2
        ct[n+i, n+1:2n] .= (1-p)÷2
    end
    ct
end

function bigger_group(g)
    n = order(g)
    mt = Matrix{Int64}(undef, 2n, 2n)
    for i=1:n, j=1:n
        k = g[i,j]
        mt[i, j] = k
        mt[n+i, j] = n+k
        mt[i, n+j] = n+k
        mt[n+i, n+j] = k
    end
    FiniteGroup(mt)
end

function check_group(g)
    n = order(g)
    g.multab[1, :] = 1:n
    g.multab[:, 1] = 1:n
    for i=1:n, j=1:n, k=1:n
        g[g[i,j],k] == g[i,g[j,k]] || error("$i, $j, $k not associate")
    end
    for i = 1:n
        sort(g.multab[i, :]) == 1:n || error("Not cycle")
    end
    true
end

ctab = charactertable(g); display(ctab)
χ = ctab[4,:];
coeff = get_coeff(g, χ);
reps = proj_reps(biggergroup(g), coeff, 2)
"""
cg = cover_group(bigger_group(g), coeff, 2);
reps = irreps(cg);
check = [check_proj_coeff(bigger_group(g), rep, coeff, 2, tol=1e-7) for rep in reps]

ct = charactertable(cg)
reg = regular_rep(cg)
p = proj_operator(ct[])
"""


text:
export cnp_1d

"""
    cnp_1d(;
        dim_embedding::Integer,
        num_encoder_layers::Integer,
        num_decoder_layers::Integer
    )

# Arguments
- `dim_embedding::Integer`: Dimensionality of the embedding.
- `num_encoder_layers::Integer`: Number of layers in the encoder.
- `num_decoder_layers::Integer`: Number of layers in the decoder.

# Returns
- `Model`: Corresponding model.
"""
function cnp_1d(;
    dim_embedding::Integer,
    num_encoder_layers::Integer,
    num_decoder_layers::Integer
)
    dim_x = 1
    dim_y = 1
    return Model(
        Parallel(
            Chain(
                InputsCoder(),
                DeterministicLikelihood()
            ),
            Chain(
                MLPCoder(
                    batched_mlp(
                        dim_in    =dim_x + dim_y,
                        dim_hidden=dim_embedding,
                        dim_out   =dim_embedding,
                        num_layers=num_encoder_layers
                    ),
                    batched_mlp(
                        dim_in    =dim_embedding,
                        dim_hidden=dim_embedding,
                        dim_out   =dim_embedding,
                        num_layers=num_encoder_layers
                    )
                ),
                DeterministicLikelihood()
            )
        ),
        Chain(
            Materialise(),
            batched_mlp(
                dim_in    =dim_x + dim_embedding,
                dim_hidden=dim_embedding,
                dim_out   =2dim_y,
                num_layers=num_decoder_layers
            ),
            HeterogeneousGaussianLikelihood()
        )
    )
end


text:
using BigFloatX

N = 100
X = randn(N,N)
Y = randn(N,N)

Xj = map(BigFloat,X);
Yj = map(BigFloat,Y);
@time Xj*Yj;
@time Xj*Yj;
@time Xj*Yj;

Xa = map(BigFloatA,X);
Ya = map(BigFloatA,Y);
@time Xa*Ya;
@time Xa*Ya;
@time Xa*Ya;

Xb = map(BigFloatB,X);
Yb = map(BigFloatB,Y);
@time Xb*Yb;
@time Xb*Yb;
@time Xb*Yb;


text:
export project

"""
    project(ϕ, f::GP, g)

The projection operation `f′(x) = Σₘ ϕₘ(x) f(m) + g(x)`. Assumes `f` is finite.
"""
function project(ϕ, f::GP, g)
    @assert isfinite(length(f))
    return GP(project, ϕ, f, g)
end

function project(ϕ, f::GP)
    if size(ϕ, 2) < Inf
        return project(ϕ, f, FiniteZeroMean(eachindex(ϕ, 2)))
    else
        return project(ϕ, f, ZeroMean{Float64}())
    end
end

function μ_p′(::typeof(project), ϕ, f, g)
    return iszero(mean(f)) ? g : DeltaSumMean(ϕ, mean(f), g)
end
function k_p′(::typeof(project), ϕ, f, g)
    if iszero(kernel(f))
        return get_zero(size(ϕ, 2), size(ϕ, 2))
    else
        return DeltaSumKernel(ϕ, kernel(f))
    end
end
function k_p′p(::typeof(project), ϕ, f, g, fp)
    if iszero(kernel(f, fp))
        return get_zero(size(ϕ, 2), length(fp))
    else
        return LhsDeltaSumCrossKernel(ϕ, kernel(f, fp))
    end
end
function k_pp′(fp, ::typeof(project), ϕ, f, g)
    if iszero(kernel(f, fp))
        return get_zero(length(fp), size(ϕ, 2))
    else
        return RhsDeltaSumCrossKernel(kernel(f, fp), ϕ)
    end
end


text:
cleanit(x::T) where T <: AbstractArray{String} = filter(!isempty, strip.(vec(x)))
function assert_metadata(folder::String)
    @assert isdir(folder) "$folder does not exist"
    @assert isdir(joinpath(folder, "metadata")) "No metadata folder found"
    @assert isfile(joinpath(folder, "metadata", "poi.csv")) "No `poi.csv` file found"
    x = cleanit(readcsv(joinpath(folder, "metadata", "poi.csv"), String))
    @assert !all(isempty.(x)) "There are no enteries in the `poi.csv` file"
    @assert allunique(x) "`poi.csv` contains duplicate items"

    @assert isfile(joinpath(folder, "metadata", "run.csv")) "No `run.csv` file found"
    x = readcsv(joinpath(folder, "metadata", "run.csv"), String)
    for i in 1:size(x,1)
        xi = cleanit(x[i,:])
        @assert !all(isempty.(xi)) "There are no enteries in the `run.csv` file at line $i"
        @assert allunique(xi) "`run.csv` contains duplicate items in line $i"
    end
end

function assert_log(folder::String)
    @assert isdir(folder) "$folder does not exist"
    isdir(joinpath(folder, "log")) || return nothing
    for i in ["associations.csv", "files.csv", "pois.csv", "runs.csv"]
        @assert isfile(joinpath(folder, "log", i)) "No `$i` file found"
    x = cleanit(readcsv(joinpath(folder, "log", "poi.csv"), String))
    @assert !all(isempty.(x)) "There are no enteries in the `poi.csv` file"
    @assert allunique(x) "`poi.csv` contains duplicate items"

    @assert isfile(joinpath(folder, "log", "run.csv")) "No `run.csv` file found"
    x = readcsv(joinpath(folder, "log", "run.csv"), String)
    for i in 1:size(x,1)
        xi = cleanit(x[i,:])
        @assert !all(isempty.(xi)) "There are no enteries in the `run.csv` file at line $i"
        @assert allunique(xi) "`run.csv` contains duplicate items in line $i"
    end
end

function second2hms(x::Second)::Dict{DataType, Int}
    ps = Dates.canonicalize(Dates.CompoundPeriod(x))
    a = Dict{DataType, Int}(k => 0 for k in [Hour, Minute, Second])
    ts = [Day, Week, Month, Year]
    for p in ps.periods
        if typeof(p) in ts
            a[Hour] += Hour(p).value
        else
            a[typeof(p)] += p.value
        end
    end
    return a
end
shorten(s::String, k::Int) = length(s) > 2k + 1 ? s[1:k]*"…"*s[(end-k + 1):end] : s
function shorten(vfs::OrderedSet{String}, m)
    nmax = maximum(map(length, vfs))
    n = min(m, nmax) - 1
    while n < nmax
        n += 1
        if allunique(shorten(vf, n) for vf in vfs)
            break
        end
    end
    return OrderedDict(shorten(vf, n) => vf for vf in vfs)
end

function openit(f::String)
    isfile(f) || systemerror("$f not found", true)
    cmd = if is_windows()
        `explorer $f`
    elseif is_linux()
        `xdg-open $f`
    elseif is_apple()
        try
            `open -a "quicktime player" $f`
        catch
            `open $f`
        end
    else
        error("Unknown OS")
    end
    return run(cmd)
    # return run(ignorestatus(cmd))
        #stream, proc = open(cmd)
        # return proc
        # try to see if you can kill the spawned process (closing the movie player). this will be useful for testing this, and for managing shit once the user is done (not sure if all the players automatically close when the user quits julia)
end
function findshortfile(v::String, d::Dict{String, String})::String
    for k in keys(d)
        d[k] == v && return k
    end
    error("Couldn't find $v in $d")
end
function validargs(_, y, m, d, rest...)
    0 < m < 13 || return Nullable{ArgumentError}(ArgumentError("Month: $m out of range (1:12)"))
    0 < d < daysinmonth(y,m)+1 || return Nullable{ArgumentError}(ArgumentError("Day: $d out of range (1:$(daysinmonth(y,m)))"))
    return Nullable{ArgumentError}()
end


text:
using CUDA: CUDA, @cuda
using DataFrames: DataFrame, Tables
using LoopVectorization: @avx

export satisfies, satisfies_all

#####################
# Circuit logical evaluation
#####################

# evaluate a circuit as a function
function (root::LogicCircuit)(data...)
    satisfies(root, data...)
end

"Evaluate satisfaction of the circuit bottom-up for a given input"
satisfies(root::LogicCircuit, data::Real...) =
    satisfies(root, collect(data))

satisfies(root::LogicCircuit, data::AbstractVector{Bool}) =
    satisfies(root, DataFrame(reshape(BitVector(data), 1, :), :auto))[1]

satisfies(root::LogicCircuit, data::AbstractVector{<:AbstractFloat}) =
    satisfies(root, DataFrame(reshape(data, 1, :), :auto))[1]

satisfies(circuit::LogicCircuit, data::DataFrame) =
    satisfies(same_device(BitCircuit(circuit, data), data) , data)

function satisfies(circuit::BitCircuit, data::DataFrame)::AbstractVector
    output = satisfies_all(circuit,data)[:,end]
    if isbinarydata(data)
        return AbstractBitVector(output, num_examples(data))
    else
        @assert isfpdata(data) "Only floating point and binary flows are supported"
        return output
    end
end

#####################
# Circuit evaluation of *all* nodes in circuit
#####################

"Evaluate the circuit bottom-up for a given input and return the value of all nodes"
function satisfies_all(circuit::BitCircuit, data, reuse=nothing)
    @assert num_features(data) == num_features(circuit)
    @assert iscomplete(data)
    values = init_satisfies(data, reuse, num_nodes(circuit))
    satisfies_layers(circuit, values)
    return values
end

"Initialize values from the data (data frames)"
function init_satisfies(data, reuse, num_nodes)
    if isbinarydata(data)
        flowtype = isgpu(data) ? CuMatrix{UInt} : Matrix{UInt}
        values = similar!(reuse, flowtype, num_chunks(data), num_nodes)

        @views values[:,TRUE_BITS] .= typemax(UInt)
        @views values[:,FALSE_BITS] .= typemin(UInt)
        nfeatures = num_features(data)
        for i=1:nfeatures
            @views values[:,2+i] .= chunks(data,i)
            @views values[:,2+nfeatures+i] .= .~ chunks(data,i)
        end
        # warning: there are now some 1 bits beyond the BitVector mask.
        # therefore, reset those to 0
        @views values[end,1:2+2*num_features(data)] .&= _msk_end(data)

    else
        @assert isfpdata(data) "Only floating point and binary flows are supported"
        flowtype = isgpu(data) ? CuMatrix{eltype(data)} : Matrix{eltype(data)}
        values = similar!(reuse, flowtype, num_examples(data), num_nodes)
        @views values[:,TRUE_BITS] .= one(eltype(values))
        @views values[:,FALSE_BITS] .= zero(eltype(values))
        nfeatures = num_features(data)
        for i=1:nfeatures
            @views values[:,2+i] .= feature_values(data,i)
            @views values[:,2+nfeatures+i] .= one(eltype(values)) .- feature_values(data,i)
        end

    end
    return values
end

# upward pass helpers on CPU

"Evaluate the layers of a bit circuit on the CPU (SIMD & multi-threaded)"
function satisfies_layers(circuit::BitCircuit, values::Matrix)
    els = circuit.elements
    for layer in circuit.layers[2:end]
        Threads.@threads for dec_id in layer
            j = @inbounds circuit.nodes[1,dec_id]
            els_end = @inbounds circuit.nodes[2,dec_id]
            if j == els_end
                assign_value(values, dec_id, els[2,j], els[3,j])
                j += 1
            else
                assign_value(values, dec_id, els[2,j], els[3,j], els[2,j+1], els[3,j+1])
                j += 2
            end
            while j <= els_end
                if j == els_end
                    accum_value(values, dec_id, els[2,j], els[3,j])
                    j += 1
                else
                    accum_value(values, dec_id, els[2,j], els[3,j], els[2,j+1], els[3,j+1])
                    j += 2
                end
            end
        end
    end
end

assign_value(v::Matrix{<:AbstractFloat}, i, e1p, e1s) =
    @views @. @avx v[:,i] = v[:,e1p] * v[:,e1s]

accum_value(v::Matrix{<:AbstractFloat}, i, e1p, e1s) =
    @views @. v[:,i] += v[:,e1p] * v[:,e1s] # adding @avx crashes macro

assign_value(v::Matrix{<:AbstractFloat}, i, e1p, e1s, e2p, e2s) =
    @views @. @avx v[:,i] = v[:,e1p] * v[:,e1s] + v[:,e2p] * v[:,e2s]

accum_value(v::Matrix{<:AbstractFloat}, i, e1p, e1s, e2p, e2s) =
    @views @. v[:,i] += v[:,e1p] * v[:,e1s] + v[:,e2p] * v[:,e2s] # adding @avx crashes macro

assign_value(v::Matrix{<:Unsigned}, i, e1p, e1s) =
    @views @. @avx v[:,i] = v[:,e1p] & v[:,e1s]

accum_value(v::Matrix{<:Unsigned}, i, e1p, e1s) =
    @views @. v[:,i] |= v[:,e1p] & v[:,e1s] # adding @avx crashes macro

assign_value(v::Matrix{<:Unsigned}, i, e1p, e1s, e2p, e2s) =
    @views @. @avx v[:,i] = v[:,e1p] & v[:,e1s] | v[:,e2p] & v[:,e2s]

accum_value(v::Matrix{<:Unsigned}, i, e1p, e1s, e2p, e2s) =
    @views @. v[:,i] |= v[:,e1p] & v[:,e1s] | v[:,e2p] & v[:,e2s] # adding @avx crashes macro

# upward pass helpers on GPU

"Evaluate the layers of a bit circuit on the GPU"
function satisfies_layers(circuit::BitCircuit, values::CuMatrix)
    CUDA.@sync for layer in circuit.layers[2:end]
        num_examples = size(values, 1)
        num_decision_sets = length(layer)

        kernel = @cuda name="satisfies_layers_cuda" launch=false satisfies_layers_cuda(layer, circuit.nodes, circuit.elements, values)
        config = launch_configuration(kernel.fun)
        threads, blocks = balance_threads_2d(num_examples, num_decision_sets, config.threads)
        kernel(layer, circuit.nodes, circuit.elements, values; threads, blocks)
    end
end

"assign threads to examples and decision nodes so that everything is a power of 2"
function balance_threads(num_examples, num_decisions, total_log2)
    ratio = num_examples / num_decisions
    k = ceil(Int, (log2(ratio) + total_log2)/2)
    k = min(max(0, k), total_log2)
    l = total_log2-k
    (2^k, 2^l)
end

"""

"""
function balance_threads_2d(num_examples, num_decisions, total_threads_per_block)
    lsb(n) = n ⊻ ( n& (n-1))
    ratio_diff(a, b) =  ceil(Int, num_examples/a) - ceil(Int, num_decisions/b)

    n_lsb = lsb(total_threads_per_block)
    options_d1 = [Int32(2^i) for i = 0 : log2(n_lsb)]
    append!(options_d1, [total_threads_per_block / n_lsb * Int32(2^i) for i = 0 : log2(n_lsb)])
    options_d2 = [Int32(total_threads_per_block / d1) for d1 in options_d1]

    best_d1 = options_d1[1]
    best_d2 = options_d2[1]
    best_ratio = ratio_diff(best_d1, best_d2)

    for (d1, d2) in zip(options_d1, options_d2)
        cur_ratio = ratio_diff(d1, d2)
        if abs(best_ratio) > abs(cur_ratio)
            best_d1 = d1
            best_d2 = d2
            best_ratio = cur_ratio
        end
    end
    threads = (best_d1, best_d2)
    blocks = (ceil(Int, num_examples / threads[1]),
                ceil(Int, num_decisions / threads[2]))

    threads, blocks
end

"CUDA kernel for circuit evaluation"
function satisfies_layers_cuda(layer, nodes, elements, values)
    index_x = (blockIdx().x - 1) * blockDim().x + threadIdx().x
    index_y = (blockIdx().y - 1) * blockDim().y + threadIdx().y
    stride_x = blockDim().x * gridDim().x
    stride_y = blockDim().y * gridDim().y
    for j = index_x:stride_x:size(values,1)
        for i = index_y:stride_y:length(layer)
            decision_id = @inbounds layer[i]
            k = @inbounds nodes[1,decision_id]
            els_end = @inbounds nodes[2,decision_id]
            @inbounds values[j, decision_id] =
                el_value(values[j, elements[2,k]], values[j, elements[3,k]])
            while k < els_end
                k += 1
                @inbounds accum_el_value(values, j, decision_id, values[j,
                            elements[2,k]], values[j, elements[3,k]])
            end # would loop unrolling help here as on CPU? probably not?
        end
    end
    return nothing
end

el_value(p::AbstractFloat, s) = p * s
el_value(p::Unsigned, s) = p & s

accum_el_value(values, j, decision_id, p::AbstractFloat, s) =
    @inbounds values[j, decision_id] += el_value(p, s)
accum_el_value(values, j, decision_id, p::Unsigned, s) =
    @inbounds values[j, decision_id] |= el_value(p, s)

max_issues_repo_name: Juice-jl/LogicCircuits

max_forks_repo_name: Juice-jl/LogicCircuits


text:
#!/usr/bin/env julia

#==============================================================================
 = Example of utilization of the library 'ros_interface'
 =
 = Maintainer: Sidney Carvalho - sydney.rdc@gmail.com
 = Last Change: 2021 Mar 21 23:28:41
 = Info: Send and receive information from a node in the ROS environment.
 =============================================================================#

# add ros_interface wrapper to ros_interface.hpp
include("../julia/ros_interface.jl")

# to call ros_interface functions
using ROS_INTERFACE

# start ROS interface
ros_interface("exe3")

# insert the robots (if do you want to use another type of robot, change
# the second parameter of 'add_node' with the following values - for now, julia
# does not support macros from C, so use its equivalence in integer):
# T_REAL for rosaria robots (use 1)
# T_STAGE for stageros robots (use 2)
# T_TURTLE for turtlesim robots (use 3)
add_bot(1, 2, "r1", [0.0; 0.0; 0.0])

# control saturations
const MAX_VX = 1                  # maximum linear velocity
const MAX_VA = 1                  # maximum angular velocity

# dead zones
const MAX_DX = 10                 # maximum linear variation
const MIN_DX = 0.1                # minimum linear variation
const MAX_DA = pi/4               # maximum angular variation
const MIN_DA = 0.01               # minimum angular variation

# variables of proportionality
kx = 0                            # to linear control
ka = 0                            # to angular control

# reference (x, y, θ)
r = [3.0, 4.5, 0.0]

# control signals (vx, vθ)
u = zeros(2)

# pose
x = zeros(3)

# main loop
while ros_ok()
    # set the data capture frequency (s)
    clock(0.01)

    # read positions from ROS
    x = get_pose(1)

    # sides of triangle formed by current point and reference
    op = r[2] - x[2]                            # opposite side
    ad = r[1] - x[1]                            # adjacent side

    # get the reference angle (the angle that the line initiating in the current
    # point and finishing in the reference gets with the frame x axis)
    # check the quadrants and hold the reference angle between (-π, π)
    if ad >= 0 && op >= 0                       # 1st quadrant
        r[3] = atan(op/ad)
    elseif ad <= 0 && op >= 0                   # 2nd quadrant
        r[3] = pi - atan(abs(op/ad))
    elseif ad <= 0 && op <= 0                   # 3rd quadrant
        r[3] = atan(abs(op/ad)) - pi
    elseif ad >= 0 && op <= 0                   # 4th quadrant
        r[3] = atan(op/ad)
    end

    # calculate the angular variation
    da = r[3] - x[3]

    # get the shorter arc and with opposite signal (to perform the best movement)
    da > pi ? da = da - 2*pi : da < -pi ? da = 2*pi - abs(da) : 0

    # calculate the angular control gain
    if da > -MIN_DA && da < MIN_DA
        ka = 0
    elseif da > MAX_DA
        ka = 1
    elseif da < -MAX_DA
        ka = -1
    elseif abs(da) >= MIN_DA && abs(da) <= MAX_DA
        ka = da/MAX_DA
    end

    # calculate the angular control signal
    u[2] = MAX_VA*ka

    # linear distance between the reference and the current point
    dx = norm(r[1:2] - x[1:2])

    # calculate the linear control gain
    if dx < MIN_DX
        kx = 0
    elseif dx > MAX_DX
        kx = 1
    elseif dx >= MIN_DX && dx <= MAX_DX
        kx = dx/MAX_DX
    end

    # calculate the linear control signal
    u[1] = MAX_VX*kx*(1 - abs(ka))

    # send velocities to ROS
    set_vel(1, [u[1]; 0.0; u[2]])

    # print robot's position and control actions
    @printf("x_r:%.2f y_r:%.2f x:%.2f y:%.2f k_x:%.2f u_x:%.2f\n", r[1], r[2], x[1], x[2], kx, u[1])
    @printf("θ_r:%.2f θ:%.2f k_θ:%.2f u_θ:%.2f\n", r[3]*180/pi, x[3]*180/pi, ka, u[2]*180/pi)
    println("time->$(get_time())")
end