Index: sasmodels/models/bcc_paracrystal.py
===================================================================
--- sasmodels/models/bcc_paracrystal.py	(revision b0c42719ac9096c35137bf644db59e1c5915b06b)
+++ sasmodels/models/bcc_paracrystal.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -128,5 +128,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "lib/gauss150.c", "lib/sphere_form.c", "bcc_paracrystal.c"]
+source = ["lib/sas_3j1x_x.c", "lib/gauss150.c", "lib/sphere_form.c", "bcc_paracrystal.c"]
 
 # parameters for demo
Index: sasmodels/models/binary_hard_sphere.c
===================================================================
--- sasmodels/models/binary_hard_sphere.c	(revision 4f79d94ef5835e630e5531fb9a69ea5d02515e60)
+++ sasmodels/models/binary_hard_sphere.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -58,6 +58,6 @@
     qr2 = r2*q;
 
-    sc1 = sph_j1c(qr1);
-    sc2 = sph_j1c(qr2);
+    sc1 = sas_3j1x_x(qr1);
+    sc2 = sas_3j1x_x(qr2);
     b1 = r1*r1*r1*(rho1-rhos)*M_4PI_3*sc1;
     b2 = r2*r2*r2*(rho2-rhos)*M_4PI_3*sc2;
Index: sasmodels/models/binary_hard_sphere.py
===================================================================
--- sasmodels/models/binary_hard_sphere.py	(revision b0c42719ac9096c35137bf644db59e1c5915b06b)
+++ sasmodels/models/binary_hard_sphere.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -108,5 +108,5 @@
              ]
 
-source = ["lib/sph_j1c.c", "binary_hard_sphere.c"]
+source = ["lib/sas_3j1x_x.c", "binary_hard_sphere.c"]
 
 # parameters for demo and documentation
Index: sasmodels/models/core_multi_shell.c
===================================================================
--- sasmodels/models/core_multi_shell.c	(revision c5ac2b2355ae3705fd04e87490d99bfe4140d545)
+++ sasmodels/models/core_multi_shell.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -3,5 +3,5 @@
 f_constant(double q, double r, double sld)
 {
-  const double bes = sph_j1c(q * r);
+  const double bes = sas_3j1x_x(q * r);
   const double vol = M_4PI_3 * cube(r);
   return sld * vol * bes;
@@ -33,9 +33,9 @@
   f = 0.;
   for (int i=0; i<n; i++) {
-    f += M_4PI_3 * cube(r) * (sld[i] - last_sld) * sph_j1c(q*r);
+    f += M_4PI_3 * cube(r) * (sld[i] - last_sld) * sas_3j1x_x(q*r);
     last_sld = sld[i];
     r += thickness[i];
   }
-  f += M_4PI_3 * cube(r) * (solvent_sld - last_sld) * sph_j1c(q*r);
+  f += M_4PI_3 * cube(r) * (solvent_sld - last_sld) * sas_3j1x_x(q*r);
   return f * f * 1.0e-4;
 }
Index: sasmodels/models/core_multi_shell.py
===================================================================
--- sasmodels/models/core_multi_shell.py	(revision 2d73a53de4a4fb58419bdaa7cf5c5e8cddaa8c5c)
+++ sasmodels/models/core_multi_shell.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -101,5 +101,5 @@
              ]
 
-source = ["lib/sph_j1c.c", "core_multi_shell.c"]
+source = ["lib/sas_3j1x_x.c", "core_multi_shell.c"]
 
 def profile(sld_core, radius, sld_solvent, n, sld, thickness):
Index: sasmodels/models/core_shell_ellipsoid.py
===================================================================
--- sasmodels/models/core_shell_ellipsoid.py	(revision fcb33e4672814b4a8fe1ee865d430090b5fc03e6)
+++ sasmodels/models/core_shell_ellipsoid.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -137,5 +137,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "lib/gfn.c", "lib/gauss76.c",
+source = ["lib/sas_3j1x_x.c", "lib/gfn.c", "lib/gauss76.c",
           "core_shell_ellipsoid.c"]
 
Index: sasmodels/models/core_shell_sphere.py
===================================================================
--- sasmodels/models/core_shell_sphere.py	(revision 4962519cc1910422a467e4f643dba223623e0e09)
+++ sasmodels/models/core_shell_sphere.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -75,5 +75,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "lib/core_shell.c", "core_shell_sphere.c"]
+source = ["lib/sas_3j1x_x.c", "lib/core_shell.c", "core_shell_sphere.c"]
 
 demo = dict(scale=1, background=0, radius=60, thickness=10,
Index: sasmodels/models/ellipsoid.c
===================================================================
--- sasmodels/models/ellipsoid.c	(revision 73e08aee242fc777e9b6a88ee4ab9fea7732e685)
+++ sasmodels/models/ellipsoid.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -18,5 +18,5 @@
     const double r = radius_equatorial
                      * sqrt(1.0 + sin_alpha*sin_alpha*(ratio*ratio - 1.0));
-    const double f = sph_j1c(q*r);
+    const double f = sas_3j1x_x(q*r);
 
     return f*f;
Index: sasmodels/models/ellipsoid.py
===================================================================
--- sasmodels/models/ellipsoid.py	(revision 0d6e86566d51a1fba22933bdbac8c40d36f9bd76)
+++ sasmodels/models/ellipsoid.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -135,5 +135,5 @@
              ]
 
-source = ["lib/sph_j1c.c", "lib/gauss76.c", "ellipsoid.c"]
+source = ["lib/sas_3j1x_x.c", "lib/gauss76.c", "ellipsoid.c"]
 
 def ER(radius_polar, radius_equatorial):
Index: sasmodels/models/fcc_paracrystal.py
===================================================================
--- sasmodels/models/fcc_paracrystal.py	(revision 0bef47b278e819df7ea6b572375dd7bc9607227e)
+++ sasmodels/models/fcc_paracrystal.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -116,5 +116,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "lib/gauss150.c", "lib/sphere_form.c", "fcc_paracrystal.c"]
+source = ["lib/sas_3j1x_x.c", "lib/gauss150.c", "lib/sphere_form.c", "fcc_paracrystal.c"]
 
 # parameters for demo
Index: sasmodels/models/fractal.c
===================================================================
--- sasmodels/models/fractal.c	(revision 217590bce5ffcf2c4de704a4770a9cab003dd287)
+++ sasmodels/models/fractal.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -14,5 +14,5 @@
     //calculate P(q) for the spherical subunits
     const double V = M_4PI_3*cube(radius);
-    const double pq = V * square((sld_block-sld_solvent)*sph_j1c(q*radius));
+    const double pq = V * square((sld_block-sld_solvent)*sas_3j1x_x(q*radius));
 
     // scale to units cm-1 sr-1 (assuming data on absolute scale)
Index: sasmodels/models/fractal.py
===================================================================
--- sasmodels/models/fractal.py	(revision d1cfa86b03d71f574114b242dc91eefade357ba2)
+++ sasmodels/models/fractal.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -97,5 +97,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "lib/sas_gamma.c", "lib/fractal_sq.c", "fractal.c"]
+source = ["lib/sas_3j1x_x.c", "lib/sas_gamma.c", "lib/fractal_sq.c", "fractal.c"]
 
 demo = dict(volfraction=0.05,
Index: sasmodels/models/fractal_core_shell.py
===================================================================
--- sasmodels/models/fractal_core_shell.py	(revision d6f60c35e7286f24ca0fdd9d5f740e642d007b8a)
+++ sasmodels/models/fractal_core_shell.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -95,5 +95,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "lib/sas_gamma.c", "lib/core_shell.c",
+source = ["lib/sas_3j1x_x.c", "lib/sas_gamma.c", "lib/core_shell.c",
           "lib/fractal_sq.c", "fractal_core_shell.c"]
 
Index: sasmodels/models/fuzzy_sphere.py
===================================================================
--- sasmodels/models/fuzzy_sphere.py	(revision 3a48772f136eb2795903c63707474e7527b4a325)
+++ sasmodels/models/fuzzy_sphere.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -81,5 +81,5 @@
 # pylint: enable=bad-whitespace,line-too-long
 
-source = ["lib/sph_j1c.c"]
+source = ["lib/sas_3j1x_x.c"]
 
 # No volume normalization despite having a volume parameter
@@ -91,5 +91,5 @@
 Iq = """
     const double qr = q*radius;
-    const double bes = sph_j1c(qr);
+    const double bes = sas_3j1x_x(qr);
     const double qf = q*fuzziness;
     const double fq = bes * (sld - sld_solvent) * form_volume(radius) * exp(-0.5*qf*qf);
Index: sasmodels/models/lib/core_shell.c
===================================================================
--- sasmodels/models/lib/core_shell.c	(revision 3a48772f136eb2795903c63707474e7527b4a325)
+++ sasmodels/models/lib/core_shell.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -16,5 +16,5 @@
     const double core_qr = q * radius;
     const double core_contrast = core_sld - shell_sld;
-    const double core_bes = sph_j1c(core_qr);
+    const double core_bes = sas_3j1x_x(core_qr);
     const double core_volume = M_4PI_3 * cube(radius);
     double f = core_volume * core_bes * core_contrast;
@@ -23,5 +23,5 @@
     const double shell_qr = q * (radius + thickness);
     const double shell_contrast = shell_sld - solvent_sld;
-    const double shell_bes = sph_j1c(shell_qr);
+    const double shell_bes = sas_3j1x_x(shell_qr);
     const double shell_volume = M_4PI_3 * cube(radius + thickness);
     f += shell_volume * shell_bes * shell_contrast;
Index: sasmodels/models/lib/gfn.c
===================================================================
--- sasmodels/models/lib/gfn.c	(revision 3a48772f136eb2795903c63707474e7527b4a325)
+++ sasmodels/models/lib/gfn.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -18,5 +18,5 @@
     // changing to more accurate sph_j1c since the following inexplicably fails on Radeon Nano.
     //const double siq = (uq == 0.0 ? 1.0 : 3.0*(sin(uq)/uq/uq - cos(uq)/uq)/uq);
-    const double siq = sph_j1c(uq);
+    const double siq = sas_3j1x_x(uq);
     const double vc = M_4PI_3*aa*aa*bb;
     const double gfnc = siq*vc*delpc;
@@ -26,5 +26,5 @@
     const double vt = M_4PI_3*trmaj*trmaj*trmin;
     //const double sit = (ut == 0.0 ? 1.0 : 3.0*(sin(ut)/ut/ut - cos(ut)/ut)/ut);
-    const double sit = sph_j1c(ut);
+    const double sit = sas_3j1x_x(ut);
     const double gfnt = sit*vt*delps;
 
Index: sasmodels/models/lib/sphere_form.c
===================================================================
--- sasmodels/models/lib/sphere_form.c	(revision ba32cddebeff04b66d828ca58b1843244d4425fb)
+++ sasmodels/models/lib/sphere_form.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -9,5 +9,5 @@
 double sphere_form(double q, double radius, double sld, double solvent_sld)
 {
-    const double fq = sphere_volume(radius) * sph_j1c(q*radius);
+    const double fq = sphere_volume(radius) * sas_3j1x_x(q*radius);
     const double contrast = (sld - solvent_sld);
     return 1.0e-4*square(contrast * fq);
Index: sasmodels/models/linear_pearls.c
===================================================================
--- sasmodels/models/linear_pearls.c	(revision 1e7b0db0b7b4315c62b59e1d659332bcbc962b8f)
+++ sasmodels/models/linear_pearls.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -44,5 +44,5 @@
 
     //sine functions of a pearl
-    double psi = sph_j1c(q * radius);
+    double psi = sas_3j1x_x(q * radius);
 
     // N pearls contribution
Index: sasmodels/models/linear_pearls.py
===================================================================
--- sasmodels/models/linear_pearls.py	(revision 4962519cc1910422a467e4f643dba223623e0e09)
+++ sasmodels/models/linear_pearls.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -63,5 +63,5 @@
 single = False
 
-source = ["lib/sph_j1c.c", "linear_pearls.c"]
+source = ["lib/sas_3j1x_x.c", "linear_pearls.c"]
 
 demo = dict(scale=1.0, background=0.0,
Index: sasmodels/models/mass_fractal.c
===================================================================
--- sasmodels/models/mass_fractal.c	(revision 6d96b66e1e58166fbc135a8df1ddf8d37f64b54a)
+++ sasmodels/models/mass_fractal.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -5,5 +5,5 @@
 {
     //calculate P(q)
-    const double pq = square(sph_j1c(q*radius));
+    const double pq = square(sas_3j1x_x(q*radius));
 
     //calculate S(q)
Index: sasmodels/models/mass_fractal.py
===================================================================
--- sasmodels/models/mass_fractal.py	(revision 6d96b66e1e58166fbc135a8df1ddf8d37f64b54a)
+++ sasmodels/models/mass_fractal.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -86,5 +86,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "lib/sas_gamma.c", "mass_fractal.c"]
+source = ["lib/sas_3j1x_x.c", "lib/sas_gamma.c", "mass_fractal.c"]
 
 demo = dict(scale=1, background=0,
Index: sasmodels/models/multilayer_vesicle.c
===================================================================
--- sasmodels/models/multilayer_vesicle.c	(revision 3a48772f136eb2795903c63707474e7527b4a325)
+++ sasmodels/models/multilayer_vesicle.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -20,5 +20,5 @@
         // layer 1
         voli = M_4PI_3*ri*ri*ri;
-        fval += voli*sldi*sph_j1c(ri*q);
+        fval += voli*sldi*sas_3j1x_x(ri*q);
 
         ri += thick_shell;
@@ -26,5 +26,5 @@
         // layer 2
         voli = M_4PI_3*ri*ri*ri;
-        fval -= voli*sldi*sph_j1c(ri*q);
+        fval -= voli*sldi*sas_3j1x_x(ri*q);
 
         //do 2 layers at a time
Index: sasmodels/models/multilayer_vesicle.py
===================================================================
--- sasmodels/models/multilayer_vesicle.py	(revision 1a6cd57bd7b12db89d69b335e984c55413295845)
+++ sasmodels/models/multilayer_vesicle.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -107,5 +107,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "multilayer_vesicle.c"]
+source = ["lib/sas_3j1x_x.c", "multilayer_vesicle.c"]
 
 polydispersity = ["radius", "n_pairs"]
Index: sasmodels/models/onion.c
===================================================================
--- sasmodels/models/onion.c	(revision 97623415b6465256171bc225ecfa9d65aac160ce)
+++ sasmodels/models/onion.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -6,5 +6,5 @@
   const double vol = M_4PI_3 * cube(r);
   const double qr = q * r;
-  const double bes = sph_j1c(qr);
+  const double bes = sas_3j1x_x(qr);
   const double alpha = A * r/thickness;
   double result;
Index: sasmodels/models/onion.py
===================================================================
--- sasmodels/models/onion.py	(revision 97623415b6465256171bc225ecfa9d65aac160ce)
+++ sasmodels/models/onion.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -314,5 +314,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "onion.c"]
+source = ["lib/sas_3j1x_x.c", "onion.c"]
 single = False
 
Index: sasmodels/models/pearl_necklace.c
===================================================================
--- sasmodels/models/pearl_necklace.c	(revision 1e7b0db0b7b4315c62b59e1d659332bcbc962b8f)
+++ sasmodels/models/pearl_necklace.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -36,5 +36,5 @@
     const double beta = (Si(q*(A_s-radius)) - Si(q*radius)) / q_edge;
     const double gamma = Si(q_edge) / q_edge;
-    const double psi = sph_j1c(q*radius);
+    const double psi = sas_3j1x_x(q*radius);
 
     // Precomputed sinc terms
Index: sasmodels/models/pearl_necklace.py
===================================================================
--- sasmodels/models/pearl_necklace.py	(revision 21261311ccad3df50fbfc3f21eb04942110f1c9d)
+++ sasmodels/models/pearl_necklace.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -92,5 +92,5 @@
              ]
 
-source = ["lib/Si.c", "lib/sph_j1c.c", "pearl_necklace.c"]
+source = ["lib/Si.c", "lib/sas_3j1x_x.c", "pearl_necklace.c"]
 single = False  # use double precision unless told otherwise
 
Index: sasmodels/models/polymer_micelle.c
===================================================================
--- sasmodels/models/polymer_micelle.c	(revision 1e7b0db0b7b4315c62b59e1d659332bcbc962b8f)
+++ sasmodels/models/polymer_micelle.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -31,5 +31,5 @@
 
     // Self-correlation term of the core
-    const double bes_core = sph_j1c(q*radius_core);
+    const double bes_core = sas_3j1x_x(q*radius_core);
     const double term1 = square(n_aggreg*beta_core*bes_core);
 
Index: sasmodels/models/polymer_micelle.py
===================================================================
--- sasmodels/models/polymer_micelle.py	(revision bba9361c45987d76f3673d7d0a4eb16d28a8b701)
+++ sasmodels/models/polymer_micelle.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -104,5 +104,5 @@
 single = False
 
-source = ["lib/sph_j1c.c", "polymer_micelle.c"]
+source = ["lib/sas_3j1x_x.c", "polymer_micelle.c"]
 
 demo = dict(scale=1, background=0,
Index: sasmodels/models/raspberry.c
===================================================================
--- sasmodels/models/raspberry.c	(revision 1e7b0db0b7b4315c62b59e1d659332bcbc962b8f)
+++ sasmodels/models/raspberry.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -51,6 +51,6 @@
 
     //Form factors for each particle
-    psiL = sph_j1c(q*rL);
-    psiS = sph_j1c(q*rS);
+    psiL = sas_3j1x_x(q*rL);
+    psiS = sas_3j1x_x(q*rS);
 
     //Cross term between large and small particles
Index: sasmodels/models/raspberry.py
===================================================================
--- sasmodels/models/raspberry.py	(revision 40a87fa5096610e24c3ead10299f609e99f093f2)
+++ sasmodels/models/raspberry.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -152,5 +152,5 @@
              ]
 
-source = ["lib/sph_j1c.c", "raspberry.c"]
+source = ["lib/sas_3j1x_x.c", "raspberry.c"]
 
 # parameters for demo
Index: sasmodels/models/sc_paracrystal.py
===================================================================
--- sasmodels/models/sc_paracrystal.py	(revision 0bef47b278e819df7ea6b572375dd7bc9607227e)
+++ sasmodels/models/sc_paracrystal.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -133,5 +133,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "lib/sphere_form.c", "lib/gauss150.c", "sc_paracrystal.c"]
+source = ["lib/sas_3j1x_x.c", "lib/sphere_form.c", "lib/gauss150.c", "sc_paracrystal.c"]
 
 demo = dict(scale=1, background=0,
Index: sasmodels/models/sphere.py
===================================================================
--- sasmodels/models/sphere.py	(revision 7e6bea817616b3ff05f1c96c2a11a94853dc6276)
+++ sasmodels/models/sphere.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -66,5 +66,5 @@
              ]
 
-source = ["lib/sph_j1c.c", "lib/sphere_form.c"]
+source = ["lib/sas_3j1x_x.c", "lib/sphere_form.c"]
 
 # No volume normalization despite having a volume parameter
Index: sasmodels/models/spherical_sld.c
===================================================================
--- sasmodels/models/spherical_sld.c	(revision 54bcd4abeaba4faf03e24319c050286fc10969da)
+++ sasmodels/models/spherical_sld.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -34,5 +34,5 @@
     const double qr = q * r;
     const double qrsq = qr * qr;
-    const double bes = sph_j1c(qr);
+    const double bes = sas_3j1x_x(qr);
     double sinqr, cosqr;
     SINCOS(qr, sinqr, cosqr);
@@ -60,7 +60,7 @@
 
         // uniform shell; r=0 => r^3=0 => f=0, so works for core as well.
-        f -= M_4PI_3 * cube(r) * sld_l * sph_j1c(q*r);
+        f -= M_4PI_3 * cube(r) * sld_l * sas_3j1x_x(q*r);
         r += thickness[shell];
-        f += M_4PI_3 * cube(r) * sld_l * sph_j1c(q*r);
+        f += M_4PI_3 * cube(r) * sld_l * sas_3j1x_x(q*r);
 
         // iterate over sub_shells in the interface
@@ -92,5 +92,5 @@
     }
     // add in solvent effect
-    f -= M_4PI_3 * cube(r) * sld_solvent * sph_j1c(q*r);
+    f -= M_4PI_3 * cube(r) * sld_solvent * sas_3j1x_x(q*r);
 
     const double f2 = f * f * 1.0e-4;
Index: sasmodels/models/spherical_sld.py
===================================================================
--- sasmodels/models/spherical_sld.py	(revision 2d65d51719530a31ea7870f34610b834f36a1319)
+++ sasmodels/models/spherical_sld.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -214,5 +214,5 @@
              ]
 # pylint: enable=bad-whitespace, line-too-long
-source = ["lib/polevl.c", "lib/sas_erf.c", "lib/sph_j1c.c", "spherical_sld.c"]
+source = ["lib/polevl.c", "lib/sas_erf.c", "lib/sas_3j1x_x.c", "spherical_sld.c"]
 single = False  # TODO: fix low q behaviour
 
Index: sasmodels/models/surface_fractal.c
===================================================================
--- sasmodels/models/surface_fractal.c	(revision b716cc63882c24230b143debfbf45e545ba6bac2)
+++ sasmodels/models/surface_fractal.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -15,5 +15,5 @@
 {
     // calculate P(q)
-    const double pq = square(sph_j1c(q*radius));
+    const double pq = square(sas_3j1x_x(q*radius));
 
     // calculate S(q)
Index: sasmodels/models/surface_fractal.py
===================================================================
--- sasmodels/models/surface_fractal.py	(revision 5c94f414f50d8c0d9e3d33e06b3e085b7730049a)
+++ sasmodels/models/surface_fractal.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -72,5 +72,5 @@
 # pylint: enable=bad-whitespace, line-too-long
 
-source = ["lib/sph_j1c.c", "lib/sas_gamma.c", "surface_fractal.c"]
+source = ["lib/sas_3j1x_x.c", "lib/sas_gamma.c", "surface_fractal.c"]
 
 demo = dict(scale=1, background=1e-5,
Index: sasmodels/models/triaxial_ellipsoid.c
===================================================================
--- sasmodels/models/triaxial_ellipsoid.c	(revision 3a48772f136eb2795903c63707474e7527b4a325)
+++ sasmodels/models/triaxial_ellipsoid.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -37,5 +37,5 @@
             const double ysq = square(Gauss76Z[j]*zm + zb);
             const double t = q*sqrt(acosx2 + bsinx2*(1.0-ysq) + c2*ysq);
-            const double fq = sph_j1c(t);
+            const double fq = sas_3j1x_x(t);
             inner += Gauss76Wt[j] * fq * fq ;
         }
@@ -64,5 +64,5 @@
                           + radius_equat_major*radius_equat_major*cmu*cmu
                           + radius_polar*radius_polar*calpha*calpha);
-    const double fq = sph_j1c(t);
+    const double fq = sas_3j1x_x(t);
     const double s = (sld - sld_solvent) * form_volume(radius_equat_minor, radius_equat_major, radius_polar);
 
Index: sasmodels/models/triaxial_ellipsoid.py
===================================================================
--- sasmodels/models/triaxial_ellipsoid.py	(revision 416f5c7f5af38d5e5d926a5b6b140a94e7e60502)
+++ sasmodels/models/triaxial_ellipsoid.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -104,5 +104,5 @@
              ]
 
-source = ["lib/sph_j1c.c", "lib/gauss76.c", "triaxial_ellipsoid.c"]
+source = ["lib/sas_3j1x_x.c", "lib/gauss76.c", "triaxial_ellipsoid.c"]
 
 def ER(radius_equat_minor, radius_equat_major, radius_polar):
Index: sasmodels/models/vesicle.c
===================================================================
--- sasmodels/models/vesicle.c	(revision 3a48772f136eb2795903c63707474e7527b4a325)
+++ sasmodels/models/vesicle.c	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -30,10 +30,10 @@
     contrast = sld_solvent-sld;
     vol = M_4PI_3*cube(radius);
-    f = vol * sph_j1c(q*radius) * contrast;
+    f = vol * sas_3j1x_x(q*radius) * contrast;
  
     //now the shell. No volume normalization as this is done by the caller
     contrast = sld-sld_solvent;
     vol = M_4PI_3*cube(radius+thickness);
-    f += vol * sph_j1c(q*(radius+thickness)) * contrast;
+    f += vol * sas_3j1x_x(q*(radius+thickness)) * contrast;
 
     //rescale to [cm-1]. 
Index: sasmodels/models/vesicle.py
===================================================================
--- sasmodels/models/vesicle.py	(revision 3a48772f136eb2795903c63707474e7527b4a325)
+++ sasmodels/models/vesicle.py	(revision 925ad6edca2969db8fb763eb02beee2c5a5f2872)
@@ -94,5 +94,5 @@
              ]
 
-source = ["lib/sph_j1c.c", "vesicle.c"]
+source = ["lib/sas_3j1x_x.c", "vesicle.c"]
 
 def ER(radius, thickness):
