Changeset a5af4e1 in sasmodels for sasmodels/models/lib/j1_cephes.c


Ignore:
Timestamp:
Mar 18, 2016 8:39:55 AM (8 years ago)
Author:
wojciech
Branches:
master, core_shell_microgels, costrafo411, magnetic_model, release_v0.94, release_v0.95, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
Children:
7904790
Parents:
c094758
Message:

Polevl fixed

File:
1 edited

Legend:

Unmodified
Added
Removed

  • sasmodels/models/lib/j1_cephes.c

    re2af2a9 ra5af4e1  
    3939Copyright 1984, 1987, 1989, 2000 by Stephen L. Moshier 
    4040*/ 
    41 double j1(double ); 
     41double J1(double ); 
    4242 
    43 double j1(double x) { 
     43double J1(double x) { 
    4444 
    4545//Cephes double pression function 
    … …  
    4848    double w, z, p, q, xn; 
    4949 
    50     const double DR1 = 5.78318596294678452118E0; 
    51     const double DR2 = 3.04712623436620863991E1; 
    5250    const double Z1 = 1.46819706421238932572E1; 
    5351    const double Z2 = 4.92184563216946036703E1; 
    5452    const double THPIO4 =  2.35619449019234492885; 
    5553    const double SQ2OPI = 0.79788456080286535588; 
    56  
    57     double RP[8] = { 
    58     -8.99971225705559398224E8, 
    59     4.52228297998194034323E11, 
    60     -7.27494245221818276015E13, 
    61     3.68295732863852883286E15, 
    62     0.0, 
    63     0.0, 
    64     0.0, 
    65     0.0 
    66     }; 
    67  
    68     double RQ[8] = { 
    69     /* 1.00000000000000000000E0,*/ 
    70     6.20836478118054335476E2, 
    71     2.56987256757748830383E5, 
    72     8.35146791431949253037E7, 
    73     2.21511595479792499675E10, 
    74     4.74914122079991414898E12, 
    75     7.84369607876235854894E14, 
    76     8.95222336184627338078E16, 
    77     5.32278620332680085395E18, 
    78     }; 
    79  
    80     double PP[8] = { 
    81     7.62125616208173112003E-4, 
    82     7.31397056940917570436E-2, 
    83     1.12719608129684925192E0, 
    84     5.11207951146807644818E0, 
    85     8.42404590141772420927E0, 
    86     5.21451598682361504063E0, 
    87     1.00000000000000000254E0, 
    88     0.0, 
    89     }; 
    90     double PQ[8] = { 
    91     5.71323128072548699714E-4, 
    92     6.88455908754495404082E-2, 
    93     1.10514232634061696926E0, 
    94     5.07386386128601488557E0, 
    95     8.39985554327604159757E0, 
    96     5.20982848682361821619E0, 
    97     9.99999999999999997461E-1, 
    98     0.0, 
    99     }; 
    100  
    101     double QP[8] = { 
    102     5.10862594750176621635E-2, 
    103     4.98213872951233449420E0, 
    104     7.58238284132545283818E1, 
    105     3.66779609360150777800E2, 
    106     7.10856304998926107277E2, 
    107     5.97489612400613639965E2, 
    108     2.11688757100572135698E2, 
    109     2.52070205858023719784E1, 
    110     }; 
    111  
    112     double QQ[8] = { 
    113     /* 1.00000000000000000000E0,*/ 
    114     7.42373277035675149943E1, 
    115     1.05644886038262816351E3, 
    116     4.98641058337653607651E3, 
    117     9.56231892404756170795E3, 
    118     7.99704160447350683650E3, 
    119     2.82619278517639096600E3, 
    120     3.36093607810698293419E2, 
    121     0.0, 
    122     }; 
    12354 
    12455    w = x; 
    … …  
    12960        { 
    13061            z = x * x; 
    131             w = polevl( z, RP, 3 ) / p1evl( z, RQ, 8 ); 
     62            w = polevlRP( z, 3 ) / p1evlRQ( z, 8 ); 
    13263            w = w * x * (z - Z1) * (z - Z2); 
    13364            return( w ); 
    … …  
    13667    w = 5.0/x; 
    13768    z = w * w; 
    138     p = polevl( z, PP, 6)/polevl( z, PQ, 6 ); 
    139     q = polevl( z, QP, 7)/p1evl( z, QQ, 7 ); 
     69 
     70    p = polevlPP( z, 6)/polevlPQ( z, 6 ); 
     71    q = polevlQP( z, 7)/p1evlQQ( z, 7 ); 
     72 
    14073    xn = x - THPIO4; 
    14174 
    … …  
    15588 
    15689 
    157     double JP[8] = { 
    158         -4.878788132172128E-009, 
    159         6.009061827883699E-007, 
    160         -4.541343896997497E-005, 
    161         1.937383947804541E-003, 
    162         -3.405537384615824E-002, 
    163         0.0, 
    164         0.0, 
    165         0.0 
    166     }; 
    167  
    168     double MO1[8] = { 
    169         6.913942741265801E-002, 
    170         -2.284801500053359E-001, 
    171         3.138238455499697E-001, 
    172         -2.102302420403875E-001, 
    173         5.435364690523026E-003, 
    174         1.493389585089498E-001, 
    175         4.976029650847191E-006, 
    176         7.978845453073848E-001 
    177     }; 
    178  
    179     double PH1[8] = { 
    180         -4.497014141919556E+001, 
    181         5.073465654089319E+001, 
    182         -2.485774108720340E+001, 
    183         7.222973196770240E+000, 
    184         -1.544842782180211E+000, 
    185         3.503787691653334E-001, 
    186         -1.637986776941202E-001, 
    187         3.749989509080821E-001 
    188     }; 
    189  
    19090    xx = x; 
    19191    if( xx < 0 ) 
    … …  
    19595        { 
    19696            z = xx * xx; 
    197             p = (z-Z1) * xx * polevl( z, JP, 4 ); 
     97            p = (z-Z1) * xx * polevlJP( z, 4 ); 
    19898            return( p ); 
    19999        } 
    … …  
    202102    w = sqrt(q); 
    203103 
    204     p = w * polevl( q, MO1, 7); 
     104    p = w * polevlMO1( q, 7); 
    205105    w = q*q; 
    206     xn = q * polevl( w, PH1, 7) - THPIO4F; 
     106    xn = q * polevlPH1( w, 7) - THPIO4F; 
    207107    p = p * cos(xn + xx); 
    208108 
    … …  
    210110#endif 
    211111} 
    212  
Note: See TracChangeset for help on using the changeset viewer.