Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed
  • sasmodels/models/lib/j1_cephes.c

    ra5af4e1 rfad5dc1  
    3939Copyright 1984, 1987, 1989, 2000 by Stephen L. Moshier 
    4040*/ 
    41 double J1(double ); 
     41 
     42constant double RP[8] = { 
     43    -8.99971225705559398224E8, 
     44    4.52228297998194034323E11, 
     45    -7.27494245221818276015E13, 
     46    3.68295732863852883286E15, 
     47    0.0, 
     48    0.0, 
     49    0.0, 
     50    0.0 }; 
     51 
     52constant double RQ[8] = { 
     53        6.20836478118054335476E2, 
     54        2.56987256757748830383E5, 
     55        8.35146791431949253037E7, 
     56        2.21511595479792499675E10, 
     57        4.74914122079991414898E12, 
     58        7.84369607876235854894E14, 
     59        8.95222336184627338078E16, 
     60        5.32278620332680085395E18 
     61    }; 
     62 
     63constant double PP[8] = { 
     64    7.62125616208173112003E-4, 
     65    7.31397056940917570436E-2, 
     66    1.12719608129684925192E0, 
     67    5.11207951146807644818E0, 
     68    8.42404590141772420927E0, 
     69    5.21451598682361504063E0, 
     70    1.00000000000000000254E0, 
     71    0.0} ; 
     72 
     73 
     74constant double PQ[8] = { 
     75    5.71323128072548699714E-4, 
     76    6.88455908754495404082E-2, 
     77    1.10514232634061696926E0, 
     78    5.07386386128601488557E0, 
     79    8.39985554327604159757E0, 
     80    5.20982848682361821619E0, 
     81    9.99999999999999997461E-1, 
     82    0.0 }; 
     83 
     84constant double QP[8] = { 
     85    5.10862594750176621635E-2, 
     86    4.98213872951233449420E0, 
     87    7.58238284132545283818E1, 
     88    3.66779609360150777800E2, 
     89    7.10856304998926107277E2, 
     90    5.97489612400613639965E2, 
     91    2.11688757100572135698E2, 
     92    2.52070205858023719784E1 }; 
     93 
     94constant double QQ[8] = { 
     95    7.42373277035675149943E1, 
     96    1.05644886038262816351E3, 
     97    4.98641058337653607651E3, 
     98    9.56231892404756170795E3, 
     99    7.99704160447350683650E3, 
     100    2.82619278517639096600E3, 
     101    3.36093607810698293419E2, 
     102    0.0 }; 
     103 
     104constant double JP[8] = { 
     105        -4.878788132172128E-009, 
     106        6.009061827883699E-007, 
     107        -4.541343896997497E-005, 
     108        1.937383947804541E-003, 
     109        -3.405537384615824E-002, 
     110        0.0, 
     111        0.0, 
     112        0.0 
     113    }; 
     114 
     115constant double MO1[8] = { 
     116        6.913942741265801E-002, 
     117        -2.284801500053359E-001, 
     118        3.138238455499697E-001, 
     119        -2.102302420403875E-001, 
     120        5.435364690523026E-003, 
     121        1.493389585089498E-001, 
     122        4.976029650847191E-006, 
     123        7.978845453073848E-001 
     124    }; 
     125 
     126constant double PH1[8] = { 
     127        -4.497014141919556E+001, 
     128        5.073465654089319E+001, 
     129        -2.485774108720340E+001, 
     130        7.222973196770240E+000, 
     131        -1.544842782180211E+000, 
     132        3.503787691653334E-001, 
     133        -1.637986776941202E-001, 
     134        3.749989509080821E-001 
     135    }; 
    42136 
    43137double J1(double x) { 
     
    60154        { 
    61155            z = x * x; 
    62             w = polevlRP( z, 3 ) / p1evlRQ( z, 8 ); 
     156            w = polevl( z, RP, 3 ) / p1evl( z, RQ, 8 ); 
    63157            w = w * x * (z - Z1) * (z - Z2); 
    64158            return( w ); 
     
    68162    z = w * w; 
    69163 
    70     p = polevlPP( z, 6)/polevlPQ( z, 6 ); 
    71     q = polevlQP( z, 7)/p1evlQQ( z, 7 ); 
     164    p = polevl( z, PP, 6)/polevl( z, PQ, 6 ); 
     165    q = polevl( z, QP, 7)/p1evl( z, QQ, 7 ); 
    72166 
    73167    xn = x - THPIO4; 
     
    95189        { 
    96190            z = xx * xx; 
    97             p = (z-Z1) * xx * polevlJP( z, 4 ); 
     191            p = (z-Z1) * xx * polevl( z, JP, 4 ); 
    98192            return( p ); 
    99193        } 
     
    102196    w = sqrt(q); 
    103197 
    104     p = w * polevlMO1( q, 7); 
     198    p = w * polevl( q, MO1, 7); 
    105199    w = q*q; 
    106     xn = q * polevlPH1( w, 7) - THPIO4F; 
     200    xn = q * polevl( w, PH1, 7) - THPIO4F; 
    107201    p = p * cos(xn + xx); 
    108202 
     
    110204#endif 
    111205} 
     206 
     207 
     208/*double J1c(double x) { 
     209    return 
     210}*/ 
Note: See TracChangeset for help on using the changeset viewer.