slatec-bessel-cpp/f2c/d9b0mp.cpp

/* d9b0mp.f -- translated by f2c (version 20100827).
   This file no longer depends on f2c.
*/

#include "slatec-internal.hpp"

/* Table of constant values */

static integer const c__3 = 3;
static integer const c__37 = 37;
static integer const c__39 = 39;
static integer const c__40 = 40;
static integer const c__44 = 44;
static integer const c__4 = 4;
static integer const c__1 = 1;
static integer const c__2 = 2;

/* Initialized data */

static double const bm0cs[37] = { .09211656246827742712573767730182,
                                      -.001050590997271905102480716371755,
                                      1.470159840768759754056392850952e-5,
                                      -5.058557606038554223347929327702e-7,
                                      2.787254538632444176630356137881e-8,
                                      -2.062363611780914802618841018973e-9,
                                      1.870214313138879675138172596261e-10,
                                      -1.969330971135636200241730777825e-11,
                                      2.325973793999275444012508818052e-12,
                                      -3.009520344938250272851224734482e-13,
                                      4.194521333850669181471206768646e-14,
                                      -6.219449312188445825973267429564e-15,
                                      9.718260411336068469601765885269e-16,
                                      -1.588478585701075207366635966937e-16,
                                      2.700072193671308890086217324458e-17,
                                      -4.750092365234008992477504786773e-18,
                                      8.61512816260437087319170374656e-19,
                                      -1.605608686956144815745602703359e-19,
                                      3.066513987314482975188539801599e-20,
                                      -5.987764223193956430696505617066e-21,
                                      1.192971253748248306489069841066e-21,
                                      -2.420969142044805489484682581333e-22,
                                      4.996751760510616453371002879999e-23,
                                      -1.047493639351158510095040511999e-23,
                                      2.227786843797468101048183466666e-24,
                                      -4.801813239398162862370542933333e-25,
                                      1.047962723470959956476996266666e-25,
                                      -2.3138581656786153251012608e-26,
                                      5.164823088462674211635199999999e-27,
                                      -1.164691191850065389525401599999e-27,
                                      2.651788486043319282958336e-28,
                                      -6.092559503825728497691306666666e-29,
                                      1.411804686144259308038826666666e-29,
                                      -3.298094961231737245750613333333e-30,
                                      7.763931143074065031714133333333e-31,
                                      -1.841031343661458478421333333333e-31,
                                      4.395880138594310737100799999999e-32 };
static double const bth0cs[44] = { -.24901780862128936717709793789967,
                                       4.8550299609623749241048615535485e-4,
                                       -5.4511837345017204950656273563505e-6,
                                       1.3558673059405964054377445929903e-7,
                                       -5.569139890222762622758321841492e-9,
                                       3.2609031824994335304004205719468e-10,
                                       -2.4918807862461341125237903877993e-11,
                                       2.3449377420882520554352413564891e-12,
                                       -2.6096534444310387762177574766136e-13,
                                       3.3353140420097395105869955014923e-14,
                                       -4.7890000440572684646750770557409e-15,
                                       7.5956178436192215972642568545248e-16,
                                       -1.3131556016891440382773397487633e-16,
                                       2.4483618345240857495426820738355e-17,
                                       -4.8805729810618777683256761918331e-18,
                                       1.0327285029786316149223756361204e-18,
                                       -2.3057633815057217157004744527025e-19,
                                       5.4044443001892693993017108483765e-20,
                                       -1.3240695194366572724155032882385e-20,
                                       3.3780795621371970203424792124722e-21,
                                       -8.9457629157111779003026926292299e-22,
                                       2.4519906889219317090899908651405e-22,
                                       -6.9388422876866318680139933157657e-23,
                                       2.0228278714890138392946303337791e-23,
                                       -6.0628500002335483105794195371764e-24,
                                       1.864974896403763538182378839627e-24,
                                       -5.8783732384849894560245036530867e-25,
                                       1.8958591447999563485531179503513e-25,
                                       -6.2481979372258858959291620728565e-26,
                                       2.1017901684551024686638633529074e-26,
                                       -7.2084300935209253690813933992446e-27,
                                       2.5181363892474240867156405976746e-27,
                                       -8.9518042258785778806143945953643e-28,
                                       3.2357237479762298533256235868587e-28,
                                       -1.1883010519855353657047144113796e-28,
                                       4.4306286907358104820579231941731e-29,
                                       -1.6761009648834829495792010135681e-29,
                                       6.4292946921207466972532393966088e-30,
                                       -2.4992261166978652421207213682763e-30,
                                       9.8399794299521955672828260355318e-31,
                                       -3.9220375242408016397989131626158e-31,
                                       1.5818107030056522138590618845692e-31,
                                       -6.4525506144890715944344098365426e-32,
                                       2.6611111369199356137177018346367e-32 };
static double const bm02cs[40] = { .0950041514522838136933086133556,
                                       -3.801864682365670991748081566851e-4,
                                       2.258339301031481192951829927224e-6,
                                       -3.895725802372228764730621412605e-8,
                                       1.246886416512081697930990529725e-9,
                                       -6.065949022102503779803835058387e-11,
                                       4.008461651421746991015275971045e-12,
                                       -3.350998183398094218467298794574e-13,
                                       3.377119716517417367063264341996e-14,
                                       -3.964585901635012700569356295823e-15,
                                       5.286111503883857217387939744735e-16,
                                       -7.852519083450852313654640243493e-17,
                                       1.280300573386682201011634073449e-17,
                                       -2.263996296391429776287099244884e-18,
                                       4.300496929656790388646410290477e-19,
                                       -8.705749805132587079747535451455e-20,
                                       1.86586271396209514118144277205e-20,
                                       -4.210482486093065457345086972301e-21,
                                       9.956676964228400991581627417842e-22,
                                       -2.457357442805313359605921478547e-22,
                                       6.307692160762031568087353707059e-23,
                                       -1.678773691440740142693331172388e-23,
                                       4.620259064673904433770878136087e-24,
                                       -1.311782266860308732237693402496e-24,
                                       3.834087564116302827747922440276e-25,
                                       -1.151459324077741271072613293576e-25,
                                       3.547210007523338523076971345213e-26,
                                       -1.119218385815004646264355942176e-26,
                                       3.611879427629837831698404994257e-27,
                                       -1.190687765913333150092641762463e-27,
                                       4.005094059403968131802476449536e-28,
                                       -1.373169422452212390595193916017e-28,
                                       4.794199088742531585996491526437e-29,
                                       -1.702965627624109584006994476452e-29,
                                       6.149512428936330071503575161324e-30,
                                       -2.255766896581828349944300237242e-30,
                                       8.3997075092942994860616583532e-31,
                                       -3.172997595562602355567423936152e-31,
                                       1.215205298881298554583333026514e-31,
                                       -4.715852749754438693013210568045e-32 };
static double const bt02cs[39] = { -.24548295213424597462050467249324,
                                       .0012544121039084615780785331778299,
                                       -3.1253950414871522854973446709571e-5,
                                       1.4709778249940831164453426969314e-6,
                                       -9.9543488937950033643468850351158e-8,
                                       8.5493166733203041247578711397751e-9,
                                       -8.6989759526554334557985512179192e-10,
                                       1.0052099533559791084540101082153e-10,
                                       -1.2828230601708892903483623685544e-11,
                                       1.7731700781805131705655750451023e-12,
                                       -2.6174574569485577488636284180925e-13,
                                       4.0828351389972059621966481221103e-14,
                                       -6.6751668239742720054606749554261e-15,
                                       1.1365761393071629448392469549951e-15,
                                       -2.0051189620647160250559266412117e-16,
                                       3.6497978794766269635720591464106e-17,
                                       -6.83096375645823031693558437888e-18,
                                       1.3107583145670756620057104267946e-18,
                                       -2.5723363101850607778757130649599e-19,
                                       5.1521657441863959925267780949333e-20,
                                       -1.0513017563758802637940741461333e-20,
                                       2.1820381991194813847301084501333e-21,
                                       -4.6004701210362160577225905493333e-22,
                                       9.8407006925466818520953651199999e-23,
                                       -2.1334038035728375844735986346666e-23,
                                       4.6831036423973365296066286933333e-24,
                                       -1.0400213691985747236513382399999e-24,
                                       2.33491056773015100517777408e-25,
                                       -5.2956825323318615788049749333333e-26,
                                       1.2126341952959756829196287999999e-26,
                                       -2.8018897082289428760275626666666e-27,
                                       6.5292678987012873342593706666666e-28,
                                       -1.5337980061873346427835733333333e-28,
                                       3.6305884306364536682359466666666e-29,
                                       -8.6560755713629122479172266666666e-30,
                                       2.0779909972536284571238399999999e-30,
                                       -5.0211170221417221674325333333333e-31,
                                       1.2208360279441714184191999999999e-31,
                                       -2.9860056267039913454250666666666e-32 };
static double const pi4 = .785398163397448309615660845819876;

static float const eta = (float) d1mach_(3) * (float).1;
static integer const nbm0 = initds_(bm0cs, &c__37, &eta);
static integer const nbm02 = initds_(bm02cs, &c__40, &eta);
static integer const nbt02 = initds_(bt02cs, &c__39, &eta);
static double const xmax = 1. / d1mach_(4);
static integer const nbth0 = initds_(bth0cs, &c__44, &eta);

int d9b0mp_(double const *x, double *ampl, double *theta)
{
    /* Local variables */
    double z__;

/* ***BEGIN PROLOGUE  D9B0MP */
/* ***SUBSIDIARY */
/* ***PURPOSE  Evaluate the modulus and phase for the J0 and Y0 Bessel */
/*            functions. */
/* ***LIBRARY   SLATEC (FNLIB) */
/* ***CATEGORY  C10A1 */
/* ***TYPE      DOUBLE PRECISION (D9B0MP-D) */
/* ***KEYWORDS  BESSEL FUNCTION, FNLIB, MODULUS, PHASE, SPECIAL FUNCTIONS */
/* ***AUTHOR  Fullerton, W., (LANL) */
/* ***DESCRIPTION */

/* Evaluate the modulus and phase for the Bessel J0 and Y0 functions. */

/* Series for BM0        on the interval  1.56250E-02 to  6.25000E-02 */
/*                                        with weighted error   4.40E-32 */
/*                                         log weighted error  31.36 */
/*                               significant figures required  30.02 */
/*                                    decimal places required  32.14 */

/* Series for BTH0       on the interval  0.          to  1.56250E-02 */
/*                                        with weighted error   2.66E-32 */
/*                                         log weighted error  31.57 */
/*                               significant figures required  30.67 */
/*                                    decimal places required  32.40 */

/* Series for BM02       on the interval  0.          to  1.56250E-02 */
/*                                        with weighted error   4.72E-32 */
/*                                         log weighted error  31.33 */
/*                               significant figures required  30.00 */
/*                                    decimal places required  32.13 */

/* Series for BT02       on the interval  1.56250E-02 to  6.25000E-02 */
/*                                        with weighted error   2.99E-32 */
/*                                         log weighted error  31.52 */
/*                               significant figures required  30.61 */
/*                                    decimal places required  32.32 */

/* ***REFERENCES  (NONE) */
/* ***ROUTINES CALLED  D1MACH, DCSEVL, INITDS, XERMSG */
/* ***REVISION HISTORY  (YYMMDD) */
/*   770701  DATE WRITTEN */
/*   890531  Changed all specific intrinsics to generic.  (WRB) */
/*   890531  REVISION DATE from Version 3.2 */
/*   891214  Prologue converted to Version 4.0 format.  (BAB) */
/*   900315  CALLs to XERROR changed to CALLs to XERMSG.  (THJ) */
/*   900720  Routine changed from user-callable to subsidiary.  (WRB) */
/*   920618  Removed space from variable names.  (RWC, WRB) */
/* ***END PROLOGUE  D9B0MP */
/* ***FIRST EXECUTABLE STATEMENT  D9B0MP */
    if (*x < 4.) {
	xermsg_("SLATEC", "D9B0MP", "X MUST BE GE 4", &c__1, &c__2, (ftnlen)6,
		 (ftnlen)6, (ftnlen)14);
    }

    if (*x > 8.) {
	goto L20;
    }
    z__ = (128. / (*x * *x) - 5.) / 3.;
    *ampl = (dcsevl_(&z__, bm0cs, &nbm0) + .75) / sqrt(*x);
    *theta = *x - pi4 + dcsevl_(&z__, bt02cs, &nbt02) / *x;
    return 0;

L20:
    if (*x > xmax) {
	xermsg_("SLATEC", "D9B0MP", "NO PRECISION BECAUSE X IS BIG", &c__2, &
		c__2, (ftnlen)6, (ftnlen)6, (ftnlen)29);
    }

    z__ = 128. / (*x * *x) - 1.;
    *ampl = (dcsevl_(&z__, bm02cs, &nbm02) + .75) / sqrt(*x);
    *theta = *x - pi4 + dcsevl_(&z__, bth0cs, &nbth0) / *x;
    return 0;

} /* d9b0mp_ */