summaryrefslogtreecommitdiff
path: root/gsl-1.9/specfunc/airy_der.c
diff options
context:
space:
mode:
Diffstat (limited to 'gsl-1.9/specfunc/airy_der.c')
-rw-r--r--gsl-1.9/specfunc/airy_der.c888
1 files changed, 888 insertions, 0 deletions
diff --git a/gsl-1.9/specfunc/airy_der.c b/gsl-1.9/specfunc/airy_der.c
new file mode 100644
index 0000000..3830875
--- /dev/null
+++ b/gsl-1.9/specfunc/airy_der.c
@@ -0,0 +1,888 @@
+/* specfunc/airy_der.c
+ *
+ * Copyright (C) 1996, 1997, 1998, 1999, 2000 Gerard Jungman
+ *
+ * This program is free software; you can redistribute it and/or modify
+ * it under the terms of the GNU General Public License as published by
+ * the Free Software Foundation; either version 2 of the License, or (at
+ * your option) any later version.
+ *
+ * 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. See the GNU
+ * General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software
+ * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ */
+
+/* Author: G. Jungman */
+
+#include <config.h>
+#include <gsl/gsl_math.h>
+#include <gsl/gsl_errno.h>
+#include <gsl/gsl_sf_exp.h>
+#include <gsl/gsl_sf_airy.h>
+
+#include "error.h"
+
+#include "chebyshev.h"
+#include "cheb_eval_mode.c"
+
+/*-*-*-*-*-*-*-*-*-*-*-* Private Section *-*-*-*-*-*-*-*-*-*-*-*/
+
+
+/* based on SLATEC aide.f, bide.f, aid.f, bid.f, r9admp.f */
+
+/*
+ series for aif on the interval -1.00000e+00 to 1.00000e+00
+ with weighted error 5.22e-18
+ log weighted error 17.28
+ significant figures required 16.01
+ decimal places required 17.73
+*/
+static double aif_data[8] = {
+ 0.10527461226531408809,
+ 0.01183613628152997844,
+ 0.00012328104173225664,
+ 0.00000062261225638140,
+ 0.00000000185298887844,
+ 0.00000000000363328873,
+ 0.00000000000000504622,
+ 0.00000000000000000522
+};
+static cheb_series aif_cs = {
+ aif_data,
+ 7,
+ -1, 1,
+ 7
+};
+
+/*
+ series for aig on the interval -1.00000e+00 to 1.00000e+00
+ with weighted error 3.14e-19
+ log weighted error 18.50
+ significant figures required 17.44
+ decimal places required 18.98
+*/
+static double aig_data[9] = {
+ 0.021233878150918666852,
+ 0.086315930335214406752,
+ 0.001797594720383231358,
+ 0.000014265499875550693,
+ 0.000000059437995283683,
+ 0.000000000152403366479,
+ 0.000000000000264587660,
+ 0.000000000000000331562,
+ 0.000000000000000000314
+};
+static cheb_series aig_cs = {
+ aig_data,
+ 8,
+ -1, 1,
+ 8
+};
+
+/*
+ series for aip2 on the interval 0.00000e+00 to 1.25000e-01
+ with weighted error 2.15e-17
+ log weighted error 16.67
+ significant figures required 14.27
+ decimal places required 17.26
+*/
+static double aip2_data[15] = {
+ 0.0065457691989713757,
+ 0.0023833724120774592,
+ -0.0000430700770220586,
+ 0.0000015629125858629,
+ -0.0000000815417186163,
+ 0.0000000054103738057,
+ -0.0000000004284130883,
+ 0.0000000000389497963,
+ -0.0000000000039623161,
+ 0.0000000000004428184,
+ -0.0000000000000536297,
+ 0.0000000000000069650,
+ -0.0000000000000009620,
+ 0.0000000000000001403,
+ -0.0000000000000000215
+};
+static cheb_series aip2_cs = {
+ aip2_data,
+ 14,
+ -1, 1,
+ 9
+};
+
+/*
+ series for aip1 on the interval 1.25000e-01 to 1.00000e+00
+ with weighted error 2.60e-17
+ log weighted error 16.58
+ significant figures required 14.91
+ decimal places required 17.28
+*/
+static double aip1_data[25] = {
+ 0.0358865097808301538,
+ 0.0114668575627764899,
+ -0.0007592073583861400,
+ 0.0000869517610893841,
+ -0.0000128237294298592,
+ 0.0000022062695681038,
+ -0.0000004222295185921,
+ 0.0000000874686415726,
+ -0.0000000192773588418,
+ 0.0000000044668460054,
+ -0.0000000010790108052,
+ 0.0000000002700029447,
+ -0.0000000000696480108,
+ 0.0000000000184489907,
+ -0.0000000000050027817,
+ 0.0000000000013852243,
+ -0.0000000000003908218,
+ 0.0000000000001121536,
+ -0.0000000000000326862,
+ 0.0000000000000096619,
+ -0.0000000000000028935,
+ 0.0000000000000008770,
+ -0.0000000000000002688,
+ 0.0000000000000000832,
+ -0.0000000000000000260
+};
+static cheb_series aip1_cs = {
+ aip1_data,
+ 24,
+ -1, 1,
+ 14
+};
+
+
+/*
+ series for bif on the interval -1.00000e+00 to 1.00000e+00
+ with weighted error 9.05e-18
+ log weighted error 17.04
+ significant figures required 15.83
+ decimal places required 17.49
+*/
+static double bif_data[8] = {
+ 0.1153536790828570243,
+ 0.0205007894049192875,
+ 0.0002135290278902876,
+ 0.0000010783960614677,
+ 0.0000000032094708833,
+ 0.0000000000062930407,
+ 0.0000000000000087403,
+ 0.0000000000000000090
+};
+static cheb_series bif_cs = {
+ bif_data,
+ 7,
+ -1, 1,
+ 7
+};
+
+/*
+ series for big on the interval -1.00000e+00 to 1.00000e+00
+ with weighted error 5.44e-19
+ log weighted error 18.26
+ significant figures required 17.46
+ decimal places required 18.74
+*/
+static double big_data[9] = {
+ -0.097196440416443537390,
+ 0.149503576843167066571,
+ 0.003113525387121326042,
+ 0.000024708570579821297,
+ 0.000000102949627731379,
+ 0.000000000263970373987,
+ 0.000000000000458279271,
+ 0.000000000000000574283,
+ 0.000000000000000000544
+};
+static cheb_series big_cs = {
+ big_data,
+ 8,
+ -1, 1,
+ 8
+};
+
+/*
+ series for bif2 on the interval 1.00000e+00 to 8.00000e+00
+ with weighted error 3.82e-19
+ log weighted error 18.42
+ significant figures required 17.68
+ decimal places required 18.92
+*/
+static double bif2_data[10] = {
+ 0.323493987603522033521,
+ 0.086297871535563559139,
+ 0.002994025552655397426,
+ 0.000051430528364661637,
+ 0.000000525840250036811,
+ 0.000000003561751373958,
+ 0.000000000017146864007,
+ 0.000000000000061663520,
+ 0.000000000000000171911,
+ 0.000000000000000000382
+};
+static cheb_series bif2_cs = {
+ bif2_data,
+ 9,
+ -1, 1,
+ 9
+};
+
+/*
+ series for big2 on the interval 1.00000e+00 to 8.00000e+00
+ with weighted error 3.35e-17
+ log weighted error 16.48
+ significant figures required 16.52
+ decimal places required 16.98
+*/
+static double big2_data[10] = {
+ 1.6062999463621294578,
+ 0.7449088819876088652,
+ 0.0470138738610277380,
+ 0.0012284422062548239,
+ 0.0000173222412256624,
+ 0.0000001521901652368,
+ 0.0000000009113560249,
+ 0.0000000000039547918,
+ 0.0000000000000130017,
+ 0.0000000000000000335
+};
+static cheb_series big2_cs = {
+ big2_data,
+ 9,
+ -1, 1,
+ 9
+};
+
+/*
+ series for bip2 on the interval 0.00000e+00 to 1.25000e-01
+ with weighted error 2.07e-18
+ log weighted error 17.69
+ significant figures required 16.51
+ decimal places required 18.42
+*/
+static double bip2_data[29] = {
+ -0.13269705443526630495,
+ -0.00568443626045977481,
+ -0.00015643601119611610,
+ -0.00001136737203679562,
+ -0.00000143464350991284,
+ -0.00000018098531185164,
+ 0.00000000926177343611,
+ 0.00000001710005490721,
+ 0.00000000476698163504,
+ -0.00000000035195022023,
+ -0.00000000058890614316,
+ -0.00000000006678499608,
+ 0.00000000006395565102,
+ 0.00000000001554529427,
+ -0.00000000000792397000,
+ -0.00000000000258326243,
+ 0.00000000000121655048,
+ 0.00000000000038707207,
+ -0.00000000000022487045,
+ -0.00000000000004953477,
+ 0.00000000000004563782,
+ 0.00000000000000332998,
+ -0.00000000000000921750,
+ 0.00000000000000094157,
+ 0.00000000000000167154,
+ -0.00000000000000055134,
+ -0.00000000000000022369,
+ 0.00000000000000017487,
+ 0.00000000000000000207
+};
+static cheb_series bip2_cs = {
+ bip2_data,
+ 28,
+ -1, 1,
+ 14
+};
+
+/*
+ series for bip1 on the interval 1.25000e-01 to 3.53553e-01
+ with weighted error 1.86e-17
+ log weighted error 16.73
+ significant figures required 15.67
+ decimal places required 17.42
+*/
+static double bip1_data[24] = {
+ -0.1729187351079553719,
+ -0.0149358492984694364,
+ -0.0005471104951678566,
+ 0.0001537966292958408,
+ 0.0000154353476192179,
+ -0.0000065434113851906,
+ 0.0000003728082407879,
+ 0.0000002072078388189,
+ -0.0000000658173336470,
+ 0.0000000074926746354,
+ 0.0000000011101336884,
+ -0.0000000007265140553,
+ 0.0000000001782723560,
+ -0.0000000000217346352,
+ -0.0000000000020302035,
+ 0.0000000000019311827,
+ -0.0000000000006044953,
+ 0.0000000000001209450,
+ -0.0000000000000125109,
+ -0.0000000000000019917,
+ 0.0000000000000015154,
+ -0.0000000000000004977,
+ 0.0000000000000001155,
+ -0.0000000000000000186
+};
+static cheb_series bip1_cs = {
+ bip1_data,
+ 23,
+ -1, 1,
+ 13
+};
+
+/*
+ series for an22 on the interval -1.00000e+00 to -1.25000e-01
+ with weighted error 3.30e-17
+ log weighted error 16.48
+ significant figures required 14.95
+ decimal places required 17.24
+*/
+static double an22_data[33] = {
+ 0.0537418629629794329,
+ -0.0126661435859883193,
+ -0.0011924334106593007,
+ -0.0002032327627275655,
+ -0.0000446468963075164,
+ -0.0000113359036053123,
+ -0.0000031641352378546,
+ -0.0000009446708886149,
+ -0.0000002966562236472,
+ -0.0000000969118892024,
+ -0.0000000326822538653,
+ -0.0000000113144618964,
+ -0.0000000040042691002,
+ -0.0000000014440333684,
+ -0.0000000005292853746,
+ -0.0000000001967763374,
+ -0.0000000000740800096,
+ -0.0000000000282016314,
+ -0.0000000000108440066,
+ -0.0000000000042074801,
+ -0.0000000000016459150,
+ -0.0000000000006486827,
+ -0.0000000000002574095,
+ -0.0000000000001027889,
+ -0.0000000000000412846,
+ -0.0000000000000166711,
+ -0.0000000000000067657,
+ -0.0000000000000027585,
+ -0.0000000000000011296,
+ -0.0000000000000004645,
+ -0.0000000000000001917,
+ -0.0000000000000000794,
+ -0.0000000000000000330
+};
+static cheb_series an22_cs = {
+ an22_data,
+ 32,
+ -1, 1,
+ 18
+};
+
+/*
+ series for an21 on the interval -1.25000e-01 to -1.56250e-02
+ with weighted error 3.43e-17
+ log weighted error 16.47
+ significant figures required 14.48
+ decimal places required 17.16
+*/
+static double an21_data[24] = {
+ 0.0198313155263169394,
+ -0.0029376249067087533,
+ -0.0001136260695958196,
+ -0.0000100554451087156,
+ -0.0000013048787116563,
+ -0.0000002123881993151,
+ -0.0000000402270833384,
+ -0.0000000084996745953,
+ -0.0000000019514839426,
+ -0.0000000004783865344,
+ -0.0000000001236733992,
+ -0.0000000000334137486,
+ -0.0000000000093702824,
+ -0.0000000000027130128,
+ -0.0000000000008075954,
+ -0.0000000000002463214,
+ -0.0000000000000767656,
+ -0.0000000000000243883,
+ -0.0000000000000078831,
+ -0.0000000000000025882,
+ -0.0000000000000008619,
+ -0.0000000000000002908,
+ -0.0000000000000000993,
+ -0.0000000000000000343
+};
+static cheb_series an21_cs = {
+ an21_data,
+ 23,
+ -1, 1,
+ 12
+};
+
+/*
+ series for an20 on the interval -1.56250e-02 to 0.00000e+00
+ with weighted error 4.41e-17
+ log weighted error 16.36
+ significant figures required 14.16
+ decimal places required 16.96
+*/
+static double an20_data[16] = {
+ 0.0126732217145738027,
+ -0.0005212847072615621,
+ -0.0000052672111140370,
+ -0.0000001628202185026,
+ -0.0000000090991442687,
+ -0.0000000007438647126,
+ -0.0000000000795494752,
+ -0.0000000000104050944,
+ -0.0000000000015932426,
+ -0.0000000000002770648,
+ -0.0000000000000535343,
+ -0.0000000000000113062,
+ -0.0000000000000025772,
+ -0.0000000000000006278,
+ -0.0000000000000001621,
+ -0.0000000000000000441
+};
+static cheb_series an20_cs = {
+ an20_data,
+ 15,
+ -1, 1,
+ 8
+};
+
+/*
+ series for aph2 on the interval -1.00000e+00 to -1.25000e-01
+ with weighted error 2.94e-17
+ log weighted error 16.53
+ significant figures required 15.58
+ decimal places required 17.28
+*/
+static double aph2_data[32] = {
+ -0.2057088719781465107,
+ 0.0422196961357771922,
+ 0.0020482560511207275,
+ 0.0002607800735165006,
+ 0.0000474824268004729,
+ 0.0000105102756431612,
+ 0.0000026353534014668,
+ 0.0000007208824863499,
+ 0.0000002103236664473,
+ 0.0000000644975634555,
+ 0.0000000205802377264,
+ 0.0000000067836273921,
+ 0.0000000022974015284,
+ 0.0000000007961306765,
+ 0.0000000002813860610,
+ 0.0000000001011749057,
+ 0.0000000000369306738,
+ 0.0000000000136615066,
+ 0.0000000000051142751,
+ 0.0000000000019351689,
+ 0.0000000000007393607,
+ 0.0000000000002849792,
+ 0.0000000000001107281,
+ 0.0000000000000433412,
+ 0.0000000000000170801,
+ 0.0000000000000067733,
+ 0.0000000000000027017,
+ 0.0000000000000010835,
+ 0.0000000000000004367,
+ 0.0000000000000001769,
+ 0.0000000000000000719,
+ 0.0000000000000000294
+};
+static cheb_series aph2_cs = {
+ aph2_data,
+ 31,
+ -1, 1,
+ 16
+};
+
+/*
+ series for aph1 on the interval -1.25000e-01 to -1.56250e-02
+ with weighted error 6.38e-17
+ log weighted error 16.20
+ significant figures required 14.91
+ decimal places required 16.87
+*/
+static double aph1_data[22] = {
+ -0.1024172908077571694,
+ 0.0071697275146591248,
+ 0.0001209959363122329,
+ 0.0000073361512841220,
+ 0.0000007535382954272,
+ 0.0000001041478171741,
+ 0.0000000174358728519,
+ 0.0000000033399795033,
+ 0.0000000007073075174,
+ 0.0000000001619187515,
+ 0.0000000000394539982,
+ 0.0000000000101192282,
+ 0.0000000000027092778,
+ 0.0000000000007523806,
+ 0.0000000000002156369,
+ 0.0000000000000635283,
+ 0.0000000000000191757,
+ 0.0000000000000059143,
+ 0.0000000000000018597,
+ 0.0000000000000005950,
+ 0.0000000000000001934,
+ 0.0000000000000000638
+};
+static cheb_series aph1_cs = {
+ aph1_data,
+ 21,
+ -1, 1,
+ 10
+};
+
+/*
+ series for aph0 on the interval -1.56250e-02 to 0.00000e+00
+ with weighted error 2.29e-17
+ log weighted error 16.64
+ significant figures required 15.27
+ decimal places required 17.23
+*/
+static double aph0_data[15] = {
+ -0.0855849241130933257,
+ 0.0011214378867065261,
+ 0.0000042721029353664,
+ 0.0000000817607381483,
+ 0.0000000033907645000,
+ 0.0000000002253264423,
+ 0.0000000000206284209,
+ 0.0000000000023858763,
+ 0.0000000000003301618,
+ 0.0000000000000527010,
+ 0.0000000000000094555,
+ 0.0000000000000018709,
+ 0.0000000000000004024,
+ 0.0000000000000000930,
+ 0.0000000000000000229
+};
+static cheb_series aph0_cs = {
+ aph0_data,
+ 14,
+ -1, 1,
+ 7
+};
+
+
+static
+int
+airy_deriv_mod_phase(const double x, gsl_mode_t mode,
+ gsl_sf_result * ampl, gsl_sf_result * phi)
+{
+ const double pi34 = 2.356194490192344928847;
+ gsl_sf_result result_a;
+ gsl_sf_result result_p;
+ double a, p;
+ double sqx;
+ double x32;
+
+ if(x <= -4.0) {
+ double z = 128.0/(x*x*x) + 1.0;
+ cheb_eval_mode_e(&an20_cs, z, mode, &result_a);
+ cheb_eval_mode_e(&aph0_cs, z, mode, &result_p);
+ }
+ else if(x <= -2.0) {
+ double z = (128.0/(x*x*x) + 9.0) / 7.0;
+ cheb_eval_mode_e(&an21_cs, z, mode, &result_a);
+ cheb_eval_mode_e(&aph1_cs, z, mode, &result_p);
+ }
+ else if(x <= -1.0) {
+ double z = (16.0/(x*x*x) + 9.0) / 7.0;
+ cheb_eval_mode_e(&an22_cs, z, mode, &result_a);
+ cheb_eval_mode_e(&aph2_cs, z, mode, &result_p);
+ }
+ else {
+ ampl->val = 0.0;
+ ampl->err = 0.0;
+ phi->val = 0.0;
+ phi->err = 0.0;
+ GSL_ERROR ("x is greater than 1.0", GSL_EDOM);
+ }
+
+ a = 0.3125 + result_a.val;
+ p = -0.625 + result_p.val;
+
+ sqx = sqrt(-x);
+ x32 = x*sqx;
+
+ ampl->val = sqrt(a * sqx);
+ ampl->err = fabs(ampl->val) * (GSL_DBL_EPSILON + fabs(result_a.err/result_a.val));
+ phi->val = pi34 - x * sqx * p;
+ phi->err = fabs(phi->val) * (GSL_DBL_EPSILON + fabs(result_p.err/result_p.val));
+
+ return GSL_SUCCESS;
+}
+
+
+/*-*-*-*-*-*-*-*-*-*-*-* Functions with Error Codes *-*-*-*-*-*-*-*-*-*-*-*/
+
+int
+gsl_sf_airy_Ai_deriv_scaled_e(const double x, gsl_mode_t mode, gsl_sf_result * result)
+{
+ /* CHECK_POINTER(result) */
+
+ if(x < -1.0) {
+ gsl_sf_result a;
+ gsl_sf_result p;
+ int status_ap = airy_deriv_mod_phase(x, mode, &a, &p);
+ double c = cos(p.val);
+ result->val = a.val * c;
+ result->err = fabs(result->val * p.err) + fabs(c * a.err);
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return status_ap;
+ }
+ else if(x <= 1.0) {
+ const double x3 = x*x*x;
+ const double x2 = x*x;
+ gsl_sf_result result_c0;
+ gsl_sf_result result_c1;
+ cheb_eval_mode_e(&aif_cs, x3, mode, &result_c0);
+ cheb_eval_mode_e(&aig_cs, x3, mode, &result_c1);
+
+ result->val = (x2*(0.125 + result_c0.val) - result_c1.val) - 0.25;
+ result->err = fabs(x2*result_c0.val) + result_c1.err;
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+
+ if(x > GSL_ROOT3_DBL_EPSILON*GSL_ROOT3_DBL_EPSILON) {
+ /* scale only if x is positive */
+ double s = exp(2.0*x*sqrt(x)/3.0);
+ result->val *= s;
+ result->err *= s;
+ }
+
+ return GSL_SUCCESS;
+ }
+ else if(x <= 4.0) {
+ const double sqrtx = sqrt(x);
+ const double z = (16.0/(x*sqrtx) - 9.0)/7.0;
+ const double s = sqrt(sqrtx);
+ gsl_sf_result result_c0;
+ cheb_eval_mode_e(&aip1_cs, z, mode, &result_c0);
+ result->val = -(0.28125 + result_c0.val) * s;
+ result->err = result_c0.err * s;
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return GSL_SUCCESS;
+ }
+ else {
+ const double sqrtx = sqrt(x);
+ const double z = 16.0/(x*sqrtx) - 1.0;
+ const double s = sqrt(sqrtx);
+ gsl_sf_result result_c0;
+ cheb_eval_mode_e(&aip2_cs, z, mode, &result_c0);
+ result->val = -(0.28125 + result_c0.val) * s;
+ result->err = result_c0.err * s;
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return GSL_SUCCESS;
+ }
+}
+
+
+int
+gsl_sf_airy_Ai_deriv_e(const double x, gsl_mode_t mode, gsl_sf_result * result)
+{
+ /* CHECK_POINTER(result) */
+
+ if(x < -1.0) {
+ gsl_sf_result a;
+ gsl_sf_result p;
+ int status_ap = airy_deriv_mod_phase(x, mode, &a, &p);
+ double c = cos(p.val);
+ result->val = a.val * c;
+ result->err = fabs(result->val * p.err) + fabs(c * a.err);
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return status_ap;
+ }
+ else if(x < 1.0) {
+ const double x3 = x*x*x;
+ gsl_sf_result result_c1;
+ gsl_sf_result result_c2;
+ cheb_eval_mode_e(&aif_cs, x3, mode, &result_c1);
+ cheb_eval_mode_e(&aig_cs, x3, mode, &result_c2);
+ result->val = (x*x*(0.125 + result_c1.val) - result_c2.val) - 0.25;
+ result->err = fabs(x*x*result_c1.err) + result_c2.err;
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return GSL_SUCCESS;
+ }
+ else if(x*x*x < 9.0/4.0 * GSL_LOG_DBL_MIN*GSL_LOG_DBL_MIN) {
+ gsl_sf_result result_aps;
+ const double arg = -2.0*x*sqrt(x)/3.0;
+ const int stat_a = gsl_sf_airy_Ai_deriv_scaled_e(x, mode, &result_aps);
+ const int stat_e = gsl_sf_exp_mult_err_e(arg, 1.5*fabs(arg*GSL_DBL_EPSILON),
+ result_aps.val, result_aps.err,
+ result);
+ return GSL_ERROR_SELECT_2(stat_e, stat_a);
+ }
+ else {
+ UNDERFLOW_ERROR(result);
+ }
+}
+
+
+int
+gsl_sf_airy_Bi_deriv_scaled_e(const double x, gsl_mode_t mode, gsl_sf_result * result)
+{
+ const double atr = 8.7506905708484345; /* 16./(sqrt(8)-1) */
+ const double btr = -2.0938363213560543; /* -(sqrt(8)+1)/(sqrt(8)-1) */
+
+ /* CHECK_POINTER(result) */
+
+ if(x < -1.0) {
+ gsl_sf_result a;
+ gsl_sf_result p;
+ int status_ap = airy_deriv_mod_phase(x, mode, &a, &p);
+ double s = sin(p.val);
+ result->val = a.val * s;
+ result->err = fabs(result->val * p.err) + fabs(s * a.err);
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return status_ap;
+ }
+ else if(x < 1.0) {
+ const double x3 = x*x*x;
+ const double x2 = x*x;
+ gsl_sf_result result_c1;
+ gsl_sf_result result_c2;
+ cheb_eval_mode_e(&bif_cs, x3, mode, &result_c1);
+ cheb_eval_mode_e(&big_cs, x3, mode, &result_c2);
+ result->val = x2 * (result_c1.val + 0.25) + result_c2.val + 0.5;
+ result->err = x2 * result_c1.err + result_c2.err;
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+
+ if(x > GSL_ROOT3_DBL_EPSILON*GSL_ROOT3_DBL_EPSILON) {
+ /* scale only if x is positive */
+ const double s = exp(-2.0*x*sqrt(x)/3.0);
+ result->val *= s;
+ result->err *= s;
+ }
+
+ return GSL_SUCCESS;
+ }
+ else if(x < 2.0) {
+ const double z = (2.0*x*x*x - 9.0) / 7.0;
+ const double s = exp(-2.0*x*sqrt(x)/3.0);
+ gsl_sf_result result_c0;
+ gsl_sf_result result_c1;
+ cheb_eval_mode_e(&bif2_cs, z, mode, &result_c0);
+ cheb_eval_mode_e(&big2_cs, z, mode, &result_c1);
+ result->val = s * (x*x * (0.25 + result_c0.val) + 0.5 + result_c1.val);
+ result->err = s * (x*x * result_c0.err + result_c1.err);
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return GSL_SUCCESS;
+ }
+ else if(x < 4.0) {
+ const double sqrtx = sqrt(x);
+ const double z = atr/(x*sqrtx) + btr;
+ const double s = sqrt(sqrtx);
+ gsl_sf_result result_c0;
+ cheb_eval_mode_e(&bip1_cs, z, mode, &result_c0);
+ result->val = s * (0.625 + result_c0.val);
+ result->err = s * result_c0.err;
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return GSL_SUCCESS;
+ }
+ else {
+ const double sqrtx = sqrt(x);
+ const double z = 16.0/(x*sqrtx) - 1.0;
+ const double s = sqrt(sqrtx);
+ gsl_sf_result result_c0;
+ cheb_eval_mode_e(&bip2_cs, z, mode, &result_c0);
+ result->val = s * (0.625 + result_c0.val);
+ result->err = s * result_c0.err;
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return GSL_SUCCESS;
+ }
+}
+
+
+int
+gsl_sf_airy_Bi_deriv_e(const double x, gsl_mode_t mode, gsl_sf_result * result)
+{
+ /* CHECK_POINTER(result) */
+
+ if(x < -1.0) {
+ gsl_sf_result a;
+ gsl_sf_result p;
+ int status_ap = airy_deriv_mod_phase(x, mode, &a, &p);
+ double s = sin(p.val);
+ result->val = a.val * s;
+ result->err = fabs(result->val * p.err) + fabs(s * a.err);
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return status_ap;
+ }
+ else if(x < 1.0) {
+ const double x3 = x*x*x;
+ const double x2 = x*x;
+ gsl_sf_result result_c1;
+ gsl_sf_result result_c2;
+ cheb_eval_mode_e(&bif_cs, x3, mode, &result_c1);
+ cheb_eval_mode_e(&big_cs, x3, mode, &result_c2);
+ result->val = x2 * (result_c1.val + 0.25) + result_c2.val + 0.5;
+ result->err = x2 * result_c1.err + result_c2.err;
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return GSL_SUCCESS;
+ }
+ else if(x < 2.0) {
+ const double z = (2.0*x*x*x - 9.0) / 7.0;
+ gsl_sf_result result_c1;
+ gsl_sf_result result_c2;
+ cheb_eval_mode_e(&bif2_cs, z, mode, &result_c1);
+ cheb_eval_mode_e(&big2_cs, z, mode, &result_c2);
+ result->val = x*x * (result_c1.val + 0.25) + result_c2.val + 0.5;
+ result->err = x*x * result_c1.err + result_c2.err;
+ result->err += GSL_DBL_EPSILON * fabs(result->val);
+ return GSL_SUCCESS;
+ }
+ else if(x < GSL_ROOT3_DBL_MAX*GSL_ROOT3_DBL_MAX) {
+ gsl_sf_result result_bps;
+ const double arg = 2.0*(x*sqrt(x)/3.0);
+ int stat_b = gsl_sf_airy_Bi_deriv_scaled_e(x, mode, &result_bps);
+ int stat_e = gsl_sf_exp_mult_err_e(arg, 1.5*fabs(arg*GSL_DBL_EPSILON),
+ result_bps.val, result_bps.err,
+ result);
+ return GSL_ERROR_SELECT_2(stat_e, stat_b);
+ }
+ else {
+ OVERFLOW_ERROR(result);
+ }
+}
+
+/*-*-*-*-*-*-*-*-*-* Functions w/ Natural Prototypes *-*-*-*-*-*-*-*-*-*-*/
+
+#include "eval.h"
+
+double gsl_sf_airy_Ai_deriv_scaled(const double x, gsl_mode_t mode)
+{
+ EVAL_RESULT(gsl_sf_airy_Ai_deriv_scaled_e(x, mode, &result));
+}
+
+double gsl_sf_airy_Ai_deriv(const double x, gsl_mode_t mode)
+{
+ EVAL_RESULT(gsl_sf_airy_Ai_deriv_e(x, mode, &result));
+}
+
+double gsl_sf_airy_Bi_deriv_scaled(const double x, gsl_mode_t mode)
+{
+ EVAL_RESULT(gsl_sf_airy_Bi_deriv_scaled_e(x, mode, &result));
+}
+
+double gsl_sf_airy_Bi_deriv(const double x, gsl_mode_t mode)
+{
+ EVAL_RESULT(gsl_sf_airy_Bi_deriv_e(x, mode, &result));
+}
generated by cgit v1.2.3 (git 2.25.1) at 2024-05-01 10:28:31 +0000