diff --git a/.gitattributes b/.gitattributes index d784d52568..ddb6a52573 100644 --- a/.gitattributes +++ b/.gitattributes @@ -4379,7 +4379,6 @@ components/tachart/numlib_fix/direct.inc svneol=native#text/pascal components/tachart/numlib_fix/ipf.pas svneol=native#text/pascal components/tachart/numlib_fix/mdt.pas svneol=native#text/pascal components/tachart/numlib_fix/sle.pas svneol=native#text/pascal -components/tachart/numlib_fix/spe.pas svneol=native#text/pascal components/tachart/taanimatedsource.pas svneol=native#text/pascal components/tachart/taaxissource.pas svneol=native#text/pascal components/tachart/tabgrautils.pas svneol=native#text/pascal diff --git a/components/tachart/numlib_fix/spe.pas b/components/tachart/numlib_fix/spe.pas deleted file mode 100644 index 3a58d9bfb8..0000000000 --- a/components/tachart/numlib_fix/spe.pas +++ /dev/null @@ -1,1295 +0,0 @@ -{ - This file is part of the Numlib package. - Copyright (c) 1986-2000 by - Kees van Ginneken, Wil Kortsmit and Loek van Reij of the - Computational centre of the Eindhoven University of Technology - - FPC port Code by Marco van de Voort (marco@freepascal.org) - documentation by Michael van Canneyt (Michael@freepascal.org) - - Special functions. (Currently, most of these functions only work for - ArbFloat=REAL/DOUBLE, not for Extended(at least not at full - precision, implement the tables in the diverse functions - for extended) - - See the file COPYING.FPC, included in this distribution, - for details about the license. - - This program is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. - - **********************************************************************} -unit spe; -{$I DIRECT.INC} - -interface - -uses typ; - -{ Calculate modified Besselfunction "of the first kind" I0(x) } -function spebi0(x: ArbFloat): ArbFloat; - -{ Calculate modified Besselfunction "of the first kind" I1(x) } -function spebi1(x: ArbFloat): ArbFloat; - -{ Calculate Besselfunction "of the first kind" J0(x) } -function spebj0(x: ArbFloat): ArbFloat; - -{ Calculate Besselfunction "of the first kind" J1(x) } -function spebj1(x: ArbFloat): ArbFloat; - -{ Calculate modified Besselfunction "of the second kind" K0(x) } -function spebk0(x: ArbFloat): ArbFloat; - -{ Calculate modified Besselfunction "of the second kind" K1(x) } -function spebk1(x: ArbFloat): ArbFloat; - -{ Calculate Besselfunction "of the second kind" Y0(x) } -function speby0(x: ArbFloat): ArbFloat; - -{ Calculate Besselfunction "of the second kind" Y1(x) } -function speby1(x: ArbFloat): ArbFloat; - -{ Entier function, calculates first integer greater or equal than X} -function speent(x: ArbFloat): longint; - -{ Errorfunction ( 2/sqrt(pi)* Int(t,0,pi,exp(sqr(t)) )} -function speerf(x: ArbFloat): ArbFloat; - -{ Errorfunction's complement ( 2/sqrt(pi)* Int(t,pi,inf,exp(sqr(t)) )} -function speefc(x: ArbFloat): ArbFloat; - -{ Function to calculate the Gamma function ( int(t,0,inf,t^(x-1)*exp(-t)) } -function spegam(x: ArbFloat): ArbFloat; - -{ Function to calculate the natural logaritm of the Gamma function} -function spelga(x: ArbFloat): ArbFloat; - -{ "Calculates" the maximum of two ArbFloat values } -function spemax(a, b: ArbFloat): ArbFloat; - -{ Calculates the functionvalue of a polynomalfunction with n coefficients in a -for variable X } -function spepol(x: ArbFloat; var a: ArbFloat; n: ArbInt): ArbFloat; - -{ Calc a^b with a and b real numbers} -function spepow(a, b: ArbFloat): ArbFloat; - -{ Returns sign of x (-1 for x<0, 0 for x=0 and 1 for x>0) } -function spesgn(x: ArbFloat): ArbInt; - -{ ArcSin(x) } -function spears(x: ArbFloat): ArbFloat; - -{ ArcCos(x) } -function spearc(x: ArbFloat): ArbFloat; - -{ Sinh(x) } -function spesih(x: ArbFloat): ArbFloat; - -{ Cosh(x) } -function specoh(x: ArbFloat): ArbFloat; - -{ Tanh(x) } -function spetah(x: ArbFloat): ArbFloat; - -{ ArcSinH(x) } -function speash(x: ArbFloat): ArbFloat; - -{ ArcCosh(x) } -function speach(x: ArbFloat): ArbFloat; - -{ ArcTanH(x) } -function speath(x: ArbFloat): ArbFloat; - -implementation - -function spebi0(x: ArbFloat): ArbFloat; - -const - - xvsmall = 3.2e-9; - a1 : array[0..23] of ArbFloat = - ( 3.08508322553671039e-1, -1.86478066609466760e-1, - 1.57686843969995904e-1, -1.28895621330524993e-1, - 9.41616340200868389e-2, -6.04316795007737183e-2, - 3.41505388391452157e-2, -1.71317947935716536e-2, - 7.70061052263382555e-3, -3.12923286656374358e-3, - 1.15888319775791686e-3, -3.93934532072526720e-4, - 1.23682594989692688e-4, -3.60645571444886286e-5, - 9.81395862769787105e-6, -2.50298975966588680e-6, - 6.00566861079330132e-7, -1.36042013507151017e-7, - 2.92096163521178835e-8, -5.94856273204259507e-9, - 1.13415934215369209e-9, -2.10071360134551962e-10, - 4.44484446637868974e-11,-7.48150165756234957e-12); - a2 : array[0..26] of ArbFloat = - ( 1.43431781856850311e-1, -3.71571542566085323e-2, - 1.44861237337359455e-2, -6.30121694459896307e-3, - 2.89362046530968701e-3, -1.37638906941232170e-3, - 6.72508592273773611e-4, -3.35833513200679384e-4, - 1.70524543267970595e-4, -8.74354291104467762e-5, - 4.48739019580173804e-5, -2.28278155280668483e-5, - 1.14032404021741277e-5, -5.54917762110482949e-6, - 2.61457634142262604e-6, -1.18752840689765504e-6, - 5.18632519069546106e-7, -2.17653548816447667e-7, - 8.75291839187305722e-8, -3.34900221934314738e-8, - 1.24131668344616429e-8, -4.66215489983794905e-9, - 1.58599776268172290e-9, -3.80370174256271589e-10, - 1.23188158175419302e-10,-8.46900307934754898e-11, - 2.45185252963941089e-11); - a3: array[0..19] of ArbFloat = - ( 4.01071065066847416e-1, 2.18216817211694382e-3, - 5.59848253337377763e-5, 2.79770701849785597e-6, - 2.17160501061222148e-7, 2.36884434055843528e-8, - 3.44345025431425567e-9, 6.47994117793472057e-10, - 1.56147127476528831e-10, 4.82726630988879388e-11, - 1.89599322920800794e-11, 1.05863621425699789e-11, - 8.27719401266046976e-12, 2.82807056475555021e-12, - -4.34624739357691085e-12,-4.29417106720584499e-12, - 4.30812577328136192e-13, 1.44572313799118029e-12, - 4.73229306831831040e-14,-1.95679809047625728e-13); - - -var t : ArbFloat; - -begin - t:=abs(x); - if t <=xvsmall - then - spebi0:=1 - else - if t <= 4 - then - spebi0 := exp(t)*spepol(t/2-1, a1[0], SizeOf(a1) div SizeOf(ArbFloat) -1) - else - if t <= 12 - then - spebi0:=exp(t)*spepol(t/4-2, a2[0], SizeOf(a2) div SizeOf(ArbFloat) -1) - else { t > 12} - spebi0:=(exp(t)/sqrt(t))* - spepol(24/t-1, a3[0], SizeOf(a3) div SizeOf(ArbFloat) -1) -end; {spebi0} - -function spebi1(x: ArbFloat): ArbFloat; - - -const xvsmall = 3.2e-9; - a1: array[0..11] of ArbFloat = - ( 1.19741654963670236e+0, 9.28758890114609554e-1, - 2.68657659522092832e-1, 4.09286371827770484e-2, - 3.84763940423809498e-3, 2.45224314039278904e-4, - 1.12849795779951847e-5, 3.92368710996392755e-7, - 1.06662712314503955e-8, 2.32856921884663846e-10, - 4.17372709788222413e-12,6.24387910353848320e-14); - - a2: array[0..26] of ArbFloat = - ( 1.34142493292698178e-1, -2.99140923897405570e-2, - 9.76021102528646704e-3, -3.40759647928956354e-3, - 1.17313412855965374e-3, -3.67626180992174570e-4, - 8.47999438119288094e-5, 5.21557319070236939e-6, - -2.62051678511418163e-5, 2.47493270133518925e-5, - -1.79026222757948636e-5, 1.13818992442463952e-5, - -6.63144162982509821e-6, 3.60186151617732531e-6, - -1.83910206626348772e-6, 8.86951515545183908e-7, - -4.05456611578551130e-7, 1.76305222240064495e-7, - -7.28978293484163628e-8, 2.84961041291017650e-8, - -1.07563514207617768e-8, 4.11321223904934809e-9, - -1.41575617446629553e-9, 3.38883570696523350e-10, - -1.10970391104678003e-10, 7.79929176497056645e-11, - -2.27061376122617856e-11); - - a3: array[0..19] of ArbFloat = - ( 3.92624494204116555e-1, -6.40545360348237412e-3, - -9.12475535508497109e-5, -3.82795135453556215e-6, - -2.72684545741400871e-7, -2.82537120880041703e-8, - -3.96757162863209348e-9, -7.28107961041827952e-10, - -1.72060490748583241e-10,-5.23524129533553498e-11, - -2.02947854602758139e-11,-1.11795516742222899e-11, - -8.69631766630563635e-12,-3.05957293450420224e-12, - 4.42966462319664333e-12, 4.47735589657057690e-12, - -3.95353303949377536e-13,-1.48765082315961139e-12, - -5.77176811730370560e-14, 1.99448557598015488e-13); - -var t : ArbFloat; - -begin - t:=abs(x); - if t <= xvsmall - then - spebi1:=x/2 - else - if t <= 4 - then - spebi1:=x*spepol(sqr(t)/8-1, a1[0], sizeof(a1) div sizeof(ArbFloat)-1) - else - if t <= 12 - then - spebi1:= - exp(t)*spepol(t/4-2, a2[0], sizeof(a2) div sizeof(ArbFloat)-1)*spesgn(x) - else { t > 12} - spebi1:= - (exp(t)/sqrt(t))* - spepol(24/t-1, a3[0], sizeof(a3) div sizeof(ArbFloat)-1)*spesgn(x) -end; {spebi1} - -function spebj0(x: ArbFloat): ArbFloat; -const - - xvsmall = 3.2e-9; - tbpi = 6.36619772367581343e-1; - a1 : array[0..5] of ArbFloat = - ( 1.22200000000000000e-17, -1.94383469000000000e-12, - 7.60816359241900000e-8, -4.60626166206275050e-4, - 1.58067102332097261e-1, -8.72344235285222129e-3); - - b1 : array[0..5] of ArbFloat = - ( - 7.58850000000000000e-16, 7.84869631400000000e-11, - - 1.76194690776215000e-6, 4.81918006946760450e-3, - - 3.70094993872649779e-1, 1.57727971474890120e-1); - - c1 : array[0..4] of ArbFloat = - ( 4.12532100000000000e-14, - 2.67925353056000000e-9, - 3.24603288210050800e-5, - 3.48937694114088852e-2, - 2.65178613203336810e-1); - - d1 : array[0..13] of ArbFloat = - ( 9.99457275788251954e-1, -5.36367319213004570e-4, - 6.13741608010926000e-6, -2.05274481565160000e-7, - 1.28037614434400000e-8, -1.21211819632000000e-9, - 1.55005642880000000e-10,-2.48827276800000000e-11, - 4.78702080000000000e-12,-1.06365696000000000e-12, - 2.45294080000000000e-13,-6.41843200000000000e-14, - 3.34028800000000000e-14,-1.17964800000000000e-14); - - d2 : array[0..16] of ArbFloat = - ( -1.55551138795135187e-2, 6.83314909934390000e-5, - -1.47713883264594000e-6, 7.10621485930000000e-8, - -5.66871613024000000e-9, 6.43278173280000000e-10, - -9.47034774400000000e-11, 1.70330918400000000e-11, - -3.59094272000000000e-12, 8.59855360000000000e-13, - -2.28807680000000000e-13, 6.95193600000000000e-14, - -2.27942400000000000e-14, 4.75136000000000000e-15, - -1.14688000000000000e-15, 2.12992000000000000e-15, - -9.83040000000000000e-16); - -var t, g, y, t2, a, b, c, cx, sx : ArbFloat; - i, bov : ArbInt; - -begin - t:=abs(x); - if t<=xvsmall - then - spebj0:=1 - else - if t<=8 - then - begin - t:=0.03125*sqr(t)-1; t2:=2*t; - b:=0; c:=0; - bov:=sizeof(a1) div sizeof(ArbFloat) - 1; - for i:=0 to bov do - begin - a:=t2*c-b+a1[i]; - if i<5 - then - b:=t2*a-c+b1[i] - else - spebj0:=t*a-c+b1[i]; - if i8} - -end {spebj0}; - -function spebj1(x: ArbFloat): ArbFloat; -const - - xvsmall = 3.2e-9; - tbpi = 6.36619772367581343e-1; - a1 : array[0..5] of ArbFloat = - ( 2.95000000000000000e-18, -5.77740420000000000e-13, - 2.94970700727800000e-8, -2.60444389348580680e-4, - 1.77709117239728283e-1, -1.19180116054121687e+0); - - b1 : array[0..5] of ArbFloat = - ( -1.95540000000000000e-16, 2.52812366400000000e-11, - -7.61758780540030000e-7, 3.24027018268385747e-3, - -6.61443934134543253e-1, 6.48358770605264921e-1); - - c1 : array[0..4] of ArbFloat = - ( 1.13857200000000000e-14, -9.42421298160000000e-10, - 1.58870192399321300e-5, -2.91755248061542077e-2, - 1.28799409885767762e+0); - - d1 : array[0..13] of ArbFloat = - ( 1.00090702627808217e+0, 8.98804941670557880e-4, - -7.95969469843846000e-6, 2.45367662227560000e-7, - -1.47085129889600000e-8, 1.36030580128000000e-9, - -1.71310758400000000e-10, 2.72040729600000000e-11, - -5.19113984000000000e-12, 1.14622464000000000e-12, - -2.63372800000000000e-13, 6.86387200000000000e-14, - -3.54508800000000000e-14, 1.24928000000000000e-14); - - d2 : array[0..15] of ArbFloat = - ( 4.67768740274489776e-2, -9.62145882205441600e-5, - 1.82120185123076000e-6, -8.29196070929200000e-8, - 6.42013250344000000e-9, -7.15110504800000000e-10, - 1.03950931840000000e-10, -1.85248000000000000e-11, - 3.87554432000000000e-12, -9.23228160000000000e-13, - 2.50224640000000000e-13, -7.43936000000000000e-14, - 1.75718400000000000e-14, -4.83328000000000000e-15, - 5.32480000000000000e-15, -2.29376000000000000e-15); - -var t, g, y, t2, a, b, c, cx, sx : ArbFloat; - i, bov : ArbInt; - -begin - t:=abs(x); - if t8} -end {spebj1}; - -function spebk0(x: ArbFloat): ArbFloat; - -const - - egam = 0.57721566490153286; - xvsmall = 3.2e-9; - highexp = 745; - - a0: array[0..7] of ArbFloat = - ( 1.12896092945412762e+0, 1.32976966478338191e-1, - 4.07157485171389048e-3, 5.59702338227915383e-5, - 4.34562671546158210e-7, 2.16382411824721532e-9, - 7.49110736894134794e-12, 1.90674197514561280e-14); - - a1: array[0..8] of ArbFloat = - ( 2.61841879258687055e-1, 1.52436921799395196e-1, - 6.63513979313943827e-3, 1.09534292632401542e-4, - 9.57878493265929443e-7, 5.19906865800665633e-9, - 1.92405264219706684e-11, 5.16867886946332160e-14, - 1.05407718191360000e-16); - - a2: array[0..22] of ArbFloat = - ( 9.58210053294896496e-1, -1.42477910128828254e-1, - 3.23582010649653009e-2, -8.27780350351692662e-3, - 2.24709729617770471e-3, -6.32678357460594866e-4, - 1.82652460089342789e-4, -5.37101208898441760e-5, - 1.60185974149720562e-5, -4.83134250336922161e-6, - 1.47055796078231691e-6, -4.51017292375200017e-7, - 1.39217270224614153e-7, -4.32185089841834127e-8, - 1.34790467361340101e-8, -4.20597329258249948e-9, - 1.32069362385968867e-9, -4.33326665618780914e-10, - 1.37999268074442719e-10, -3.19241059198852137e-11, - 9.74410152270679245e-12, -7.83738609108569293e-12, - 2.57466288575820595e-12); - - a3: array[0..22] of ArbFloat = - ( 6.97761598043851776e-1, -1.08801882084935132e-1, - 2.56253646031960321e-2, -6.74459607940169198e-3, - 1.87292939725962385e-3, -5.37145622971910027e-4, - 1.57451516235860573e-4, -4.68936653814896712e-5, - 1.41376509343622727e-5, -4.30373871727268511e-6, - 1.32052261058932425e-6, -4.07851207862189007e-7, - 1.26672629417567360e-7, -3.95403255713518420e-8, - 1.23923137898346852e-8, -3.88349705250555658e-9, - 1.22424982779432970e-9, -4.03424607871960089e-10, - 1.28905587479980147e-10,-2.97787564633235128e-11, - 9.11109430833001267e-12,-7.39672783987933184e-12, - 2.43538242247537459e-12); - a4: array[0..16] of ArbFloat = - ( 1.23688664769425422e+0, -1.72683652385321641e-2, - -9.25551464765637133e-4, -9.02553345187404564e-5, - -6.31692398333746470e-6, -7.69177622529272933e-7, - -4.16044811174114579e-8, -9.41555321137176073e-9, - 1.75359321273580603e-10, -2.22829582288833265e-10, - 3.49564293256545992e-11, -1.11391758572647639e-11, - 2.85481235167705907e-12, -7.31344482663931904e-13, - 2.06328892562554880e-13, -1.28108310826991616e-13, - 4.43741979886551040e-14); - - -var t: ArbFloat; - -begin - if x<=0 - then - RunError(401); - if x<=xvsmall - then - spebk0:=-(ln(x/2)+egam) - else - if x<=1 - then - begin - t:=2*sqr(x)-1; - spebk0:=-ln(x)*spepol(t, a0[0], sizeof(a0) div sizeof(ArbFloat) - 1) - + spepol(t, a1[0], sizeof(a1) div sizeof(ArbFloat) - 1) - end - else - if x<=2 - then - spebk0:=exp(-x)*spepol(2*x-3, a2[0], sizeof(a2) div sizeof(ArbFloat) - 1) - else - if x<=4 - then - spebk0:=exp(-x)*spepol(x-3, a3[0], sizeof(a3) div sizeof(ArbFloat) - 1) - else - if x <= highexp - then - spebk0:=exp(-x)* - spepol(10/(1+x)-1, a4[0], sizeof(a4) div sizeof(ArbFloat) - 1)/sqrt(x) - else - spebk0:=0 -end; {spebk0} - -function spebk1(x: ArbFloat): ArbFloat; - -const - - xsmall = 7.9e-10; - highexp = 745; - a0: array[0..7] of ArbFloat = - ( 5.31907865913352762e-1, 3.25725988137110495e-2, - 6.71642805873498653e-4, 6.95300274548206237e-6, - 4.32764823642997753e-8, 1.79784792380155752e-10, - 5.33888268665658944e-13, 1.18964962439910400e-15); - - a1: array[0..7] of ArbFloat = - ( 3.51825828289325536e-1, 4.50490442966943726e-2, - 1.20333585658219028e-3, 1.44612432533006139e-5, - 9.96686689273781531e-8, 4.46828628435618679e-10, - 1.40917103024514301e-12, 3.29881058019865600e-15); - - a2: array[0..23] of ArbFloat = - ( 1.24316587355255299e+0, -2.71910714388689413e-1, - 8.20250220860693888e-2, -2.62545818729427417e-2, - 8.57388087067410089e-3, -2.82450787841655951e-3, - 9.34594154387642940e-4, -3.10007681013626626e-4, - 1.02982746700060730e-4, -3.42424912211942134e-5, - 1.13930169202553526e-5, -3.79227698821142908e-6, - 1.26265578331941923e-6, -4.20507152338934956e-7, - 1.40138351985185509e-7, -4.66928912168020101e-8, - 1.54456653909012693e-8, -5.13783508140332214e-9, - 1.82808381381205361e-9, -6.15211416898895086e-10, - 1.28044023949946257e-10, -4.02591066627023831e-11, - 4.27404330568767242e-11, -1.46639291782948454e-11); - - a3: array[0..23] of ArbFloat = - ( 8.06563480128786903e-1, -1.60052611291327173e-1, - 4.58591528414023064e-2, -1.42363136684423646e-2, - 4.55865751206724687e-3, -1.48185472032688523e-3, - 4.85707174778663652e-4, -1.59994873621599146e-4, - 5.28712919123131781e-5, -1.75089594354079944e-5, - 5.80692311842296724e-6, -1.92794586996432593e-6, - 6.40581814037398274e-7, -2.12969229346310343e-7, - 7.08723366696569880e-8, -2.35855618461025265e-8, - 7.79421651144832709e-9, -2.59039399308009059e-9, - 9.20781685906110546e-10, -3.09667392343245062e-10, - 6.44913423545894175e-11, -2.02680401514735862e-11, - 2.14736751065133220e-11, -7.36478297050421658e-12); - - a4: array[0..16] of ArbFloat = - ( 1.30387573604230402e+0, 5.44845254318931612e-2, - 4.31639434283445364e-3, 4.29973970898766831e-4, - 4.04720631528495020e-5, 4.32776409784235211e-6, - 4.07563856931843484e-7, 4.86651420008153956e-8, - 3.82717692121438315e-9, 6.77688943857588882e-10, - 6.97075379117731379e-12, 1.72026097285930936e-11, - -2.60774502020271104e-12, 8.58211523713560576e-13, - -2.19287104441802752e-13, 1.39321122940600320e-13, - -4.77850238111580160e-14); - -var t: ArbFloat; - -begin - if x<=0 - then - RunError(402); - if x<=xsmall - then - spebk1:=1/x - else - if x<=1 - then - begin - t:=2*sqr(x)-1; - spebk1:=( ln(x)*spepol(t, a0[0], sizeof(a0) div sizeof(ArbFloat) - 1) - -spepol(t, a1[0], sizeof(a1) div sizeof(ArbFloat) -1) )*x + 1/x - end - else - if x<=2 - then - spebk1:=exp(-x)*spepol(2*x-3, a2[0], sizeof(a2) div sizeof(ArbFloat) - 1) - else - if x<=4 - then - spebk1:=exp(-x)*spepol(x-3, a3[0], sizeof(a3) div sizeof(ArbFloat) - 1) - else - if x <= highexp - then - spebk1:=exp(-x)*spepol(10/(1+x)-1, a4[0], - sizeof(a4) div sizeof(ArbFloat) - 1)/sqrt(x) - else - spebk1:=0 -end; {spebk1} - -function speby0(x: ArbFloat): ArbFloat; - -const - - tbpi = 6.36619772367581343e-1; - egam = 5.77215664901532861e-1; - xvsmall = 3.2e-9; - a1 : array[0..5] of ArbFloat = - ( 3.90000000000000000e-19, -8.74734100000000000e-14, - 5.24879478733000000e-9, -5.63207914105698700e-5, - 4.71966895957633869e-2, 1.79034314077182663e-1); - - b1 : array[0..5] of ArbFloat = - ( -2.69800000000000000e-17, 4.02633082000000000e-12, - -1.44072332740190000e-7, 7.53113593257774230e-4, - -1.77302012781143582e-1, -2.74474305529745265e-1); - - c1 : array[0..5] of ArbFloat = - ( 1.64349000000000000e-15, -1.58375525420000000e-10, - 3.20653253765480000e-6, -7.28796247955207918e-3, - 2.61567346255046637e-1, -3.31461132032849417e-2); - - d1 : array[0..13] of ArbFloat = - ( 9.99457275788251954e-1, -5.36367319213004570e-4, - 6.13741608010926000e-6, -2.05274481565160000e-7, - 1.28037614434400000e-8, -1.21211819632000000e-9, - 1.55005642880000000e-10,-2.48827276800000000e-11, - 4.78702080000000000e-12,-1.06365696000000000e-12, - 2.45294080000000000e-13,-6.41843200000000000e-14, - 3.34028800000000000e-14,-1.17964800000000000e-14); - - d2 : array[0..16] of ArbFloat = - (-1.55551138795135187e-2, 6.83314909934390000e-5, - -1.47713883264594000e-6, 7.10621485930000000e-8, - -5.66871613024000000e-9, 6.43278173280000000e-10, - -9.47034774400000000e-11, 1.70330918400000000e-11, - -3.59094272000000000e-12, 8.59855360000000000e-13, - -2.28807680000000000e-13, 6.95193600000000000e-14, - -2.27942400000000000e-14, 4.75136000000000000e-15, - -1.14688000000000000e-15, 2.12992000000000000e-15, - -9.83040000000000000e-16); - -var t, g, y, t2, a, b, c, cx, sx : ArbFloat; - i, bov : ArbInt; - -begin - if x<=0 - then - RunError(403); - if x<=xvsmall - then - speby0:=(ln(x/2)+egam)*tbpi - else - if x<=8 - then - begin - t:=0.03125*sqr(x)-1; t2:=2*t; - b:=0; c:=0; - bov:=sizeof(a1) div sizeof(ArbFloat) - 1; - for i:=0 to bov do - begin - a:=t2*c-b+a1[i]; - b:=t2*a-c+b1[i]; - if i8} -end {speby0}; - -function speby1(x: ArbFloat): ArbFloat; - -const - tbpi = 6.36619772367581343e-1; - xsmall = 7.9e-10; - a1 : array[0..5] of ArbFloat = - (-6.58000000000000000e-18, 1.21143321000000000e-12, - -5.68844003991900000e-8, 4.40478629867099510e-4, - -2.26624991556754924e-1, -1.28697384381350001e-1); - - b1 : array[0..5] of ArbFloat = - ( 4.27730000000000000e-16,-5.17212147300000000e-11, - 1.41662436449235000e-6, -5.13164116106108479e-3, - 6.75615780772187667e-1, 2.03041058859342538e-2); - - c1 : array[0..4] of ArbFloat = - (-2.44094900000000000e-14, 1.87547032473000000e-9, - -2.83046401495148000e-5, 4.23191803533369041e-2, - -7.67296362886645940e-1); - - d1 : array[0..13] of ArbFloat = - ( 1.00090702627808217e+0, 8.98804941670557880e-4, - -7.95969469843846000e-6, 2.45367662227560000e-7, - -1.47085129889600000e-8, 1.36030580128000000e-9, - -1.71310758400000000e-10, 2.72040729600000000e-11, - -5.19113984000000000e-12, 1.14622464000000000e-12, - -2.63372800000000000e-13, 6.86387200000000000e-14, - -3.54508800000000000e-14, 1.24928000000000000e-14); - - d2 : array[0..15] of ArbFloat = - ( 4.67768740274489776e-2, -9.62145882205441600e-5, - 1.82120185123076000e-6, -8.29196070929200000e-8, - 6.42013250344000000e-9, -7.15110504800000000e-10, - 1.03950931840000000e-10,-1.85248000000000000e-11, - 3.87554432000000000e-12,-9.23228160000000000e-13, - 2.50224640000000000e-13,-7.43936000000000000e-14, - 1.75718400000000000e-14,-4.83328000000000000e-15, - 5.32480000000000000e-15,-2.29376000000000000e-15); - -var t, g, y, t2, a, b, c, cx, sx : ArbFloat; - i, bov : ArbInt; - -begin - if x<=0 - then - RunError(404); - if x<=xsmall - then - speby1:=-tbpi/x - else - if x<=8 - then - begin - t:=0.03125*sqr(x)-1; t2:=2*t; - b:=0; c:=0; - bov:=sizeof(a1) div sizeof(ArbFloat) - 1; - for i:=0 to bov do - begin - a:=t2*c-b+a1[i]; - if i8} -end {speby1}; - -function speent(x : ArbFloat): longint; - -var tx : longint; - -begin - tx:=trunc(x); - if x>=0 - then - speent:=tx - else - if x=tx - then - speent:=tx - else - speent:=tx-1 -end; {speent} - -function speerf(x : ArbFloat) : ArbFloat; -const - - xup = 6.25; - sqrtpi = 1.7724538509055160; - c : array[1..18] of ArbFloat = - ( 1.9449071068178803e0, 4.20186582324414e-2, -1.86866103976769e-2, - 5.1281061839107e-3, -1.0683107461726e-3, 1.744737872522e-4, - -2.15642065714e-5, 1.7282657974e-6, -2.00479241e-8, - -1.64782105e-8, 2.0008475e-9, 2.57716e-11, - -3.06343e-11, 1.9158e-12, 3.703e-13, - -5.43e-14, -4.0e-15, 1.2e-15); - - d : array[1..17] of ArbFloat = - ( 1.4831105640848036e0, -3.010710733865950e-1, 6.89948306898316e-2, - -1.39162712647222e-2, 2.4207995224335e-3, -3.658639685849e-4, - 4.86209844323e-5, -5.7492565580e-6, 6.113243578e-7, - -5.89910153e-8, 5.2070091e-9, -4.232976e-10, - 3.18811e-11, -2.2361e-12, 1.467e-13, - -9.0e-15, 5.0e-16); - - var t, s, s1, s2, x2: ArbFloat; - bovc, bovd, j: ArbInt; -begin - bovc:=sizeof(c) div sizeof(ArbFloat); - bovd:=sizeof(d) div sizeof(ArbFloat); - t:=abs(x); - if t <= 2 - then - begin - x2:=sqr(x)-2; - s1:=d[bovd]; s2:=0; j:=bovd-1; - s:=x2*s1-s2+d[j]; - while j > 1 do - begin - s2:=s1; s1:=s; j:=j-1; - s:=x2*s1-s2+d[j] - end; - speerf:=(s-s2)*x/2 - end - else - if t < xup - then - begin - x2:=2-20/(t+3); - s1:=c[bovc]; s2:=0; j:=bovc-1; - s:=x2*s1-s2+c[j]; - while j > 1 do - begin - s2:=s1; s1:=s; j:=j-1; - s:=x2*s1-s2+c[j] - end; - x2:=((s-s2)/(2*t))*exp(-sqr(x))/sqrtpi; - speerf:=(1-x2)*spesgn(x) - end - else - speerf:=spesgn(x) -end; {speerf} - -function spemax(a, b: ArbFloat): ArbFloat; -begin - if a>b - then - spemax:=a - else - spemax:=b -end; {spemax} - -function speefc(x : ArbFloat) : ArbFloat; -const - - xlow = -6.25; - xhigh = 27.28; - c : array[0..22] of ArbFloat = - ( 1.455897212750385e-1, -2.734219314954260e-1, - 2.260080669166197e-1, -1.635718955239687e-1, - 1.026043120322792e-1, -5.480232669380236e-2, - 2.414322397093253e-2, -8.220621168415435e-3, - 1.802962431316418e-3, -2.553523453642242e-5, - -1.524627476123466e-4, 4.789838226695987e-5, - 3.584014089915968e-6, -6.182369348098529e-6, - 7.478317101785790e-7, 6.575825478226343e-7, - -1.822565715362025e-7, -7.043998994397452e-8, - 3.026547320064576e-8, 7.532536116142436e-9, - -4.066088879757269e-9, -5.718639670776992e-10, - 3.328130055126039e-10); - - var t, s : ArbFloat; -begin - if x <= xlow - then - speefc:=2 - else - if x >= xhigh - then - speefc:=0 - else - begin - t:=1-7.5/(abs(x)+3.75); - s:=exp(-x*x)*spepol(t, c[0], sizeof(c) div sizeof(ArbFloat) - 1); - if x < 0 - then - speefc:=2-s - else - speefc:=s - end -end {speefc}; - -function spegam(x: ArbFloat): ArbFloat; -const - - tmax = 170; - a: array[0..23] of ArbFloat = - ( 8.86226925452758013e-1, 1.61691987244425092e-2, - 1.03703363422075456e-1, -1.34118505705967765e-2, - 9.04033494028101968e-3, -2.42259538436268176e-3, - 9.15785997288933120e-4, -2.96890121633200000e-4, - 1.00928148823365120e-4, -3.36375833240268800e-5, - 1.12524642975590400e-5, -3.75499034136576000e-6, - 1.25281466396672000e-6, -4.17808776355840000e-7, - 1.39383522590720000e-7, -4.64774927155200000e-8, - 1.53835215257600000e-8, -5.11961333760000000e-9, - 1.82243164160000000e-9, -6.13513953280000000e-10, - 1.27679856640000000e-10,-4.01499750400000000e-11, - 4.26560716800000000e-11,-1.46381209600000000e-11); - -var tvsmall, t, g: ArbFloat; - m, i: ArbInt; -begin - if sizeof(ArbFloat) = 6 - then - {%H-}tvsmall:=2*midget - else - tvsmall:=midget; - t:=abs(x); - if t > tmax - then - RunError(407); - if t < macheps - then - begin - if t < tvsmall - then - RunError(407); - spegam:=1/x - end - else { abs(x) >= macheps } - begin - m:=trunc(x); - if x > 0 - then - begin - t:=x-m; m:=m-1; g:=1; - if m<0 - then - g:=g/x - else - if m>0 - then - for i:=1 to m do - g:=(x-i)*g - end - else { x < 0 } - begin - t:=x-m+1; - if t=1 - then - RunError(407); - m:=1-m; - g:=x; - for i:=1 to m do - g:=(i+x)*g; - g:=1/g - end; - spegam:=spepol(2*t-1, a[0], sizeof(a) div sizeof(ArbFloat) - 1)*g - end { abs(x) >= macheps } -end; {spegam} - -function spelga(x: ArbFloat): ArbFloat; - -const - - xbig = 7.7e7; - xmax = 2.559e305; - lnr2pi = 9.18938533204672742e-1; - a: array[0..23] of ArbFloat = - ( 8.86226925452758013e-1, 1.61691987244425092e-2, - 1.03703363422075456e-1, -1.34118505705967765e-2, - 9.04033494028101968e-3, -2.42259538436268176e-3, - 9.15785997288933120e-4, -2.96890121633200000e-4, - 1.00928148823365120e-4, -3.36375833240268800e-5, - 1.12524642975590400e-5, -3.75499034136576000e-6, - 1.25281466396672000e-6, -4.17808776355840000e-7, - 1.39383522590720000e-7, -4.64774927155200000e-8, - 1.53835215257600000e-8, -5.11961333760000000e-9, - 1.82243164160000000e-9, -6.13513953280000000e-10, - 1.27679856640000000e-10,-4.01499750400000000e-11, - 4.26560716800000000e-11,-1.46381209600000000e-11); - b: array[0..5] of ArbFloat = - ( 8.33271644065786580e-2, -6.16502049453716986e-6, - 3.89978899876484712e-9, -6.45101975779653651e-12, - 2.00201927337982364e-14, -9.94561064728159347e-17); - - -var t, g : ArbFloat; - m, i : ArbInt; - -begin - if x <= 0 - then - RunError(408); - if x <= macheps - then - spelga:=-ln(x) - else - if x <= 15 - then - begin - m:=trunc(x); t:=x-m; m:=m-1; g:=1; - if m < 0 - then - g:=g/x - else - if m > 0 - then - for i:=1 to m do - g:=(x-i)*g; - spelga:=ln(g*spepol(2*t-1, a[0], sizeof(a) div sizeof(ArbFloat) - 1)) - end - else { x > 15 } - if x <= xbig - then - spelga:=(x-0.5)*ln(x) - x + lnr2pi - + spepol(450/sqr(x)-1, b[0], sizeof(b) div sizeof(ArbFloat) - 1)/x - else { x > xbig } - if x <= xmax - then - spelga:=(x-0.5)*ln(x) - x + lnr2pi - else { x > xmax => x*ln(x) > giant } - RunError(408) -end; {spelga} - -function spepol(x: ArbFloat; var a: ArbFloat; n: ArbInt): ArbFloat; -var pa : ^arfloat0; - i : ArbInt; - polx : ArbFloat; -begin - pa:=@a; - polx:=0; - for i:=n downto 0 do - polx:=polx*x+pa^[i]; - spepol:=polx -end {spepol}; - -function spepow(a, b: ArbFloat): ArbFloat; - - function PowInt(a: double; n: longint): double; - var a1 : double; - begin - if n=0 then PowInt := 1 else - begin - a1 := 1; - if n<0 then begin a := 1/a; n := -n end; - while n>0 - do if Odd(n) - then begin Dec(n); a1 := a1*a end - else begin n := n div 2; a := sqr(a) end; - PowInt := a1 - end - end; - -var tb : longint; - fb : double; -begin - - { (a < 0, b niet geheel) of (a = 0, b <= 0), dan afbreken} - if (a=0) then if (b<=0) then RunError(400) else begin SpePow :=0; exit end; - tb := Trunc(b); fb := b-tb; - if (a<0) and (fb<>0) then RunError(400); - - if a>0 - then if fb=0 then SpePow := PowInt(a, tb) - else SpePow := PowInt(a, tb)*exp(fb*ln(a)) - else if odd(tb) then Spepow := -PowInt(-a, tb) - else SpePow := PowInt(-a, tb) - -end; {spepow} - -function spesgn(x: ArbFloat): ArbInt; - -begin - if x<0 - then - spesgn:=-1 - else - if x=0 - then - spesgn:=0 - else - spesgn:=1 -end; {spesgn} - -function spears(x: ArbFloat): ArbFloat; -const - - pi2 = 1.570796326794897; - a : array[0..17] of ArbFloat = - ( 1.047197551196598e+0, 5.375149359132719e-2, 7.798902238957732e-3, - 1.519668539582420e-3, 3.408637238430600e-4, 8.302317819598986e-5, - 2.134554822576075e-5, 5.701781046148566e-6, 1.566985123962741e-6, - 4.402076871418002e-7, 1.257811162594110e-7, 3.646577948300129e-8, - 1.081021746966715e-8, 3.212744286269388e-9, 8.515014302985799e-10, - 2.513296398553196e-10, 1.342121568282535e-10, 4.210346761190271e-11); - -var y, u, t, s : ArbFloat; - uprang : boolean; -begin - if abs(x) > 1 - then - RunError(401); - u:=sqr(x); uprang:= u > 0.5; - if uprang - then - u:=1-u; - t:=4*u-1; y:=spepol(t, a[0], sizeof(a) div sizeof(ArbFloat) - 1); - if uprang - then - begin - s:=pi2-sqrt(u)*y; - if x < 0 - then - s:=-s; - spears:=s - end - else - spears:=x*y -end; {spears} - -function spearc(x: ArbFloat): ArbFloat; -const - - pi2 = 1.570796326794897; - a : array[0..17] of ArbFloat = - ( 1.047197551196598e+0, 5.375149359132719e-2, 7.798902238957732e-3, - 1.519668539582420e-3, 3.408637238430600e-4, 8.302317819598986e-5, - 2.134554822576075e-5, 5.701781046148566e-6, 1.566985123962741e-6, - 4.402076871418002e-7, 1.257811162594110e-7, 3.646577948300129e-8, - 1.081021746966715e-8, 3.212744286269388e-9, 8.515014302985799e-10, - 2.513296398553196e-10, 1.342121568282535e-10, 4.210346761190271e-11); - -var u, t, y, s : ArbFloat; - uprang : boolean; -begin - if abs(x)>1.0 - then - RunError(402); - u:=sqr(x); uprang:=u>0.5; - if uprang - then - u:=1-u; - t:=4*u-1; y:=spepol(t, a[0], sizeof(a) div sizeof(ArbFloat) - 1); - if uprang - then - begin - s:=sqrt(u)*y; - if x<0 - then - s:=2*pi2-s; - spearc:=s - end - else - spearc:=pi2-x*y -end; {spearc} - -function spesih(x: ArbFloat): ArbFloat; -const - - a : array[0..6] of ArbFloat = - ( 1.085441641272607e+0, 8.757509762437522e-2, 2.158779361257021e-3, - 2.549839945498292e-5, 1.761854853281383e-7, 7.980704288665359e-10, - 2.551377137317034e-12); - -var u : ArbFloat; -begin - if abs(x)<=1.0 - then - begin - u:=2*sqr(x)-1; - spesih:=x*spepol(u, a[0], sizeof(a) div sizeof(ArbFloat) - 1) - end - else - begin - u:=exp(x); spesih:=(u-1/u)/2 - end -end; {spesih} - -function specoh(x: ArbFloat): ArbFloat; -var u: ArbFloat; -begin - u:=exp(x); specoh:=(u+1/u)/2 -end; {specoh} - -function spetah(x: ArbFloat): ArbFloat; -const - xhi = 18.50; - a : array[0..15] of ArbFloat = - ( 8.610571715805476e-1, -1.158834489728470e-1, 1.918072383973950e-2, - -3.225255180728459e-3, 5.433071386922689e-4, -9.154289983175165e-5, - 1.542469328074432e-5, -2.599022539340038e-6, 4.379282308765732e-7, - -7.378980192173815e-8, 1.243517352745986e-8, -2.095373768837420e-9, - 3.509758916273561e-10,-5.908745181531817e-11, 1.124199312776748e-11, - -1.907888434471600e-12); - -var t, y: ArbFloat; - -begin - t:=abs(x); - if t <= 1 - then - begin - y:=2*sqr(x)-1; - spetah:=x*spepol(y, a[0], sizeof(a) div sizeof(ArbFloat) - 1) - end - else - if t < xhi - then - begin - y:=exp(2*x); spetah:=(y-1)/(y+1) - end - else - spetah:=spesgn(x) -end; {spetah} - -function speash(x: ArbFloat): ArbFloat; -const - - xhi = 1e9; - c : array[0..18] of ArbFloat = - ( 9.312298594527122e-1, -5.736663926249348e-2, - 9.004288574881897e-3, -1.833458667045431e-3, - 4.230023450529706e-4, -1.050715136470630e-4, - 2.740790473603819e-5, -7.402952157663977e-6, - 2.052474396638805e-6, -5.807433412373489e-7, - 1.670117348345774e-7, -4.863477336087045e-8, - 1.432753532351304e-8, -4.319978113584910e-9, - 1.299779213740398e-9, -3.394726871170490e-10, - 1.008344962167889e-10, -5.731943029121004e-11, - 1.810792296549804e-11); - - -var t : ArbFloat; - -begin - t:=abs(x); - if t <= 1 then - speash:=x*spepol(2*sqr(x)-1, c[0], sizeof(c) div sizeof(ArbFloat) - 1) - else - if t < xhi then - speash:=ln(sqrt(sqr(x)+1)+t)*spesgn(x) - else - speash:=ln(2*t)*spesgn(x) -end; {speash} - -function speach(x: ArbFloat): ArbFloat; -const - - xhi = 1e9; - -begin - if x<1 then - RunError(405); - if x=1 then - speach:=0 - else - if x<=xhi then - speach:=ln(x+sqrt(sqr(x)-1)) - else - speach:=ln(2*x) -end; {speach} - -function speath(x: ArbFloat): ArbFloat; -const - - c : array[0..19] of ArbFloat = - ( 1.098612288668110e+0, 1.173605223326117e-1, 2.309071936165689e-2, - 5.449091889986991e-3, 1.404884102286929e-3, 3.816948426588058e-4, - 1.073604335435426e-4, 3.095027782918129e-5, 9.088050814470148e-6, - 2.706881064641104e-6, 8.155200644023077e-7, 2.479830612463254e-7, - 7.588067811453948e-8, 2.339295963220429e-8, 7.408486568719348e-9, - 2.319454882064018e-9, 5.960921368486746e-10, 1.820410351379402e-10, - 1.184977617320312e-10, 3.856235316559190e-11); - -var t, u: ArbFloat; -begin - t:=abs(x); - if t >= 1 then - RunError(406); - u:=sqr(x); - if u < 0.5 then - speath:=x*spepol(4*u-1, c[0], sizeof(c) div sizeof(ArbFloat) - 1) - else { 0.5 < x*x < 1 } - speath:=ln((1+t)/(1-t))/2*spesgn(x) -end; {speath} - -var exitsave : pointer; - -procedure MyExit; -{ -const ErrorS : array[400..408,1..6] of char = - ('spepow', - 'spebk0', - 'spebk1', - 'speby0', - 'speby1', - 'speach', - 'speath', - 'spegam', - 'spelga'); -} -begin - ExitProc := ExitSave; - if (ExitCode>=400) AND (ExitCode<=408) then - begin - //write(ErrFil, 'critical error in ', ErrorS[ExitCode]); - ExitCode := 201 - end; -end; - -begin - ExitSave := ExitProc; - ExitProc := @MyExit; -end.