Index: explore/precision.py
===================================================================
--- explore/precision.py	(revision 2a7e20e9c4ba72da5c6088e743a64f3345ee554d)
+++ explore/precision.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -357,10 +357,12 @@
 )
 add_function(
-    name="debye",
+    name="gauss_coil",
     mp_function=lambda x: 2*(mp.exp(-x**2) + x**2 - 1)/x**4,
     np_function=lambda x: 2*(np.expm1(-x**2) + x**2)/x**4,
     ocl_function=make_ocl("""
     const double qsq = q*q;
-    if (qsq < 1.0) { // Pade approximation
+    // For double: use O(5) Pade with 0.5 cutoff (10 mad + 1 divide)
+    // For single: use O(7) Taylor with 0.8 cutoff (7 mad)
+    if (qsq < 0.0) {
         const double x = qsq;
         if (0) { // 0.36 single
@@ -372,5 +374,5 @@
             const double B1=3./8., B2=3./56., B3=1./336.;
             return (((A3*x + A2)*x + A1)*x + 1.)/(((B3*x + B2)*x + B1)*x + 1.);
-        } else if (1) { // 1.0 for single, 0.25 for double
+        } else if (0) { // 1.0 for single, 0.25 for double
             // PadeApproximant[2*Exp[-x^2] + x^2-1)/x^4, {x, 0, 8}]
             const double A1=1./15., A2=1./60, A3=0., A4=1./75600.;
@@ -385,5 +387,5 @@
                   /(((((B5*x + B4)*x + B3)*x + B2)*x + B1)*x + 1.);
         }
-    } else if (qsq < 1.) { // Taylor series; 0.9 for single, 0.25 for double
+    } else if (qsq < 0.8) {
         const double x = qsq;
         const double C0 = +1.;
Index: sasmodels/core.py
===================================================================
--- sasmodels/core.py	(revision 53998094412e4e529314e9fc45914cc59e148317)
+++ sasmodels/core.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -364,5 +364,6 @@
     actual = [p.name for p in model.info.parameters.kernel_parameters]
     target = ("sld sld_solvent radius length theta phi"
-              " radius_effective volfraction structure_factor_mode"
+              " radius_effective volfraction "
+              " structure_factor_mode radius_effective_mode"
               " A_sld A_sld_solvent A_radius").split()
     assert target == actual, "%s != %s"%(target, actual)
Index: sasmodels/kernel.py
===================================================================
--- sasmodels/kernel.py	(revision 53998094412e4e529314e9fc45914cc59e148317)
+++ sasmodels/kernel.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -92,3 +92,3 @@
     def _call_kernel(self, call_details, values, cutoff, magnetic, effective_radius_type):
         # type: (CallDetails, np.ndarray, np.ndarray, float, bool, int) -> np.ndarray
-        raise NotImpmentedError()
+        raise NotImplementedError()
Index: sasmodels/models/barbell.c
===================================================================
--- sasmodels/models/barbell.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/barbell.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -79,11 +79,12 @@
 effective_radius(int mode, double radius_bell, double radius, double length)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(radius_bell, radius , length);
-    } else if (mode == 2) {
+    case 2: // radius
         return radius;
-    } else if (mode == 3) {
+    case 3: // half length
         return 0.5*length;
-    } else {
+    case 4: // half total length
         return radius_from_totallength(radius_bell,radius,length);
     }
Index: sasmodels/models/barbell.py
===================================================================
--- sasmodels/models/barbell.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/barbell.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -117,5 +117,7 @@
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "barbell.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","radius","half length","half total length"]
+effective_radius_type = [
+    "equivalent sphere", "radius", "half length", "half total length",
+    ]
 
 def random():
Index: sasmodels/models/capped_cylinder.c
===================================================================
--- sasmodels/models/capped_cylinder.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/capped_cylinder.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -101,11 +101,12 @@
 effective_radius(int mode, double radius, double radius_cap, double length)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(radius, radius_cap, length);
-    } else if (mode == 2) {
+    case 2: // radius
         return radius;
-    } else if (mode == 3) {
+    case 3: // half length
         return 0.5*length;
-    } else {
+    case 4: // half total length
         return radius_from_totallength(radius, radius_cap,length);
     }
Index: sasmodels/models/capped_cylinder.py
===================================================================
--- sasmodels/models/capped_cylinder.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/capped_cylinder.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -137,5 +137,7 @@
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "capped_cylinder.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","radius","half length","half total length"]
+effective_radius_type = [
+    "equivalent sphere", "radius", "half length", "half total length",
+    ]
 
 def random():
Index: sasmodels/models/core_multi_shell.c
===================================================================
--- sasmodels/models/core_multi_shell.c	(revision a94046f61e19c19ba17cb07580d048db76ec896c)
+++ sasmodels/models/core_multi_shell.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -6,15 +6,4 @@
   const double vol = M_4PI_3 * cube(r);
   return sld * vol * bes;
-}
-
-static double
-form_volume(double core_radius, double fp_n, double thickness[])
-{
-  double r = core_radius;
-  int n = (int)(fp_n+0.5);
-  for (int i=0; i < n; i++) {
-    r += thickness[i];
-  }
-  return M_4PI_3 * cube(r);
 }
 
@@ -31,22 +20,18 @@
 
 static double
+form_volume(double core_radius, double fp_n, double thickness[])
+{
+  return M_4PI_3 * cube(outer_radius(core_radius, fp_n, thickness));
+}
+
+static double
 effective_radius(int mode, double core_radius, double fp_n, double thickness[])
-// this seems regardless to always give the result for outer radius for n=1 shells; why??
-// printf shows fp_n is always 1, not 0,1,2
 {
-//        printf("fp_n =%g \n",fp_n);
-        if (mode == 1) {
-        double r = core_radius;
-        int n = (int)(fp_n+0.5);
-        if ( n > 0) {
-            for (int i=0; i < n; i++) {
-                r += thickness[i];
-            }
-        }
-        return r;
-        //return outer_radius(core_radius,fp_n,thickness);
-    } else {
-        return core_radius;
-    }
+  switch (mode) {
+  case 1: // outer radius
+    return outer_radius(core_radius, fp_n, thickness);
+  case 2: // core radius
+    return core_radius;
+  }
 }
 
Index: sasmodels/models/core_multi_shell.py
===================================================================
--- sasmodels/models/core_multi_shell.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_multi_shell.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -145,9 +145,4 @@
     return np.asarray(z), np.asarray(rho)
 
-#def ER(radius, n, thickness):
-#    """Effective radius"""
-#    n = int(n[0]+0.5)  # n is a control parameter and is not polydisperse
-#    return np.sum(thickness[:n], axis=0) + radius
-
 demo = dict(sld_core=6.4,
             radius=60,
Index: sasmodels/models/core_shell_bicelle.c
===================================================================
--- sasmodels/models/core_shell_bicelle.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_shell_bicelle.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -54,11 +54,12 @@
 effective_radius(int mode, double radius, double thick_rim, double thick_face, double length)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(radius, thick_rim, thick_face, length);
-    } else if (mode == 2) {
+    case 2: // outer rim radius
         return radius + thick_rim;
-    } else if (mode == 3) {
+    case 3: // half outer thickness
         return 0.5*length + thick_face;
-    } else {
+    case 4: // half diagonal
         return radius_from_diagonal(radius,thick_rim,thick_face,length);
     }
Index: sasmodels/models/core_shell_bicelle.py
===================================================================
--- sasmodels/models/core_shell_bicelle.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_shell_bicelle.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -155,5 +155,8 @@
           "core_shell_bicelle.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","outer rim radius", "half outer thickness","half diagonal"]
+effective_radius_type = [
+    "equivalent sphere", "outer rim radius",
+    "half outer thickness", "half diagonal",
+    ]
 
 def random():
Index: sasmodels/models/core_shell_bicelle_elliptical.c
===================================================================
--- sasmodels/models/core_shell_bicelle_elliptical.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_shell_bicelle_elliptical.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -29,15 +29,16 @@
 effective_radius(int mode, double r_minor, double x_core, double thick_rim, double thick_face, double length)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(r_minor, x_core, thick_rim, thick_face, length);
-    } else if (mode == 2) {
+    case 2: // outer rim average radius
         return 0.5*r_minor*(1.0 + x_core) + thick_rim;
-    } else if (mode == 3) {
+    case 3: // outer rim min radius
         return (x_core < 1.0 ? x_core*r_minor+thick_rim : r_minor+thick_rim);
-    } else if (mode == 4) {
+    case 4: // outer max radius
         return (x_core > 1.0 ? x_core*r_minor+thick_rim : r_minor+thick_rim);
-    } else if (mode ==5) {
+    case 5: // half outer thickness
         return 0.5*length + thick_face;
-    } else {
+    case 6: // half diagonal
         return radius_from_diagonal(r_minor,x_core,thick_rim,thick_face,length);
     }
Index: sasmodels/models/core_shell_bicelle_elliptical.py
===================================================================
--- sasmodels/models/core_shell_bicelle_elliptical.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_shell_bicelle_elliptical.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -147,6 +147,8 @@
           "core_shell_bicelle_elliptical.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","outer rim average radius","outer rim min radius",
-                         "outer max radius", "half outer thickness","half diagonal"]
+effective_radius_type = [
+    "equivalent sphere", "outer rim average radius", "outer rim min radius",
+    "outer max radius", "half outer thickness", "half diagonal",
+    ]
 
 def random():
Index: sasmodels/models/core_shell_bicelle_elliptical_belt_rough.c
===================================================================
--- sasmodels/models/core_shell_bicelle_elliptical_belt_rough.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_shell_bicelle_elliptical_belt_rough.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -30,15 +30,16 @@
 effective_radius(int mode, double r_minor, double x_core, double thick_rim, double thick_face, double length)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(r_minor, x_core, thick_rim, thick_face, length);
-    } else if (mode == 2) {
+    case 2: // outer rim average radius
         return 0.5*r_minor*(1.0 + x_core) + thick_rim;
-    } else if (mode == 3) {
+    case 3: // outer rim min radius
         return (x_core < 1.0 ? x_core*r_minor+thick_rim : r_minor+thick_rim);
-    } else if (mode == 4) {
+    case 4: // outer max radius
         return (x_core > 1.0 ? x_core*r_minor+thick_rim : r_minor+thick_rim);
-    } else if (mode ==5) {
+    case 5: // half outer thickness
         return 0.5*length + thick_face;
-    } else {
+    case 6: // half diagonal
         return radius_from_diagonal(r_minor,x_core,thick_rim,thick_face,length);
     }
Index: sasmodels/models/core_shell_bicelle_elliptical_belt_rough.py
===================================================================
--- sasmodels/models/core_shell_bicelle_elliptical_belt_rough.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_shell_bicelle_elliptical_belt_rough.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -160,6 +160,8 @@
           "core_shell_bicelle_elliptical_belt_rough.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","outer rim average radius","outer rim min radius",
-                         "outer max radius", "half outer thickness","half diagonal"]
+effective_radius_type = [
+    "equivalent sphere", "outer rim average radius", "outer rim min radius",
+    "outer max radius", "half outer thickness", "half diagonal",
+    ]
 
 demo = dict(scale=1, background=0,
Index: sasmodels/models/core_shell_cylinder.c
===================================================================
--- sasmodels/models/core_shell_cylinder.c	(revision a94046f61e19c19ba17cb07580d048db76ec896c)
+++ sasmodels/models/core_shell_cylinder.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -30,18 +30,17 @@
 static double
 effective_radius(int mode, double radius, double thickness, double length)
-//effective_radius_type = ["equivalent sphere","outer radius","half outer length","half min outer dimension",
-//                         "half max outer dimension","half outer diagonal"]
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(radius, thickness, length);
-    } else if (mode == 2) {
+    case 2: // outer radius
         return radius + thickness;
-    } else if (mode == 3) {
+    case 3: // half outer length
         return 0.5*length + thickness;
-    } else if (mode == 4) {
+    case 4: // half min outer length
         return (radius < 0.5*length ? radius + thickness : 0.5*length + thickness);
-    } else if (mode == 5) {
+    case 5: // half max outer length
         return (radius > 0.5*length ? radius + thickness : 0.5*length + thickness);
-    } else {
+    case 6: // half outer diagonal
         return radius_from_diagonal(radius,thickness,length);
     }
Index: sasmodels/models/core_shell_cylinder.py
===================================================================
--- sasmodels/models/core_shell_cylinder.py	(revision e31b19a46729ff0f8f1a5eb9847f69524a6ed1d3)
+++ sasmodels/models/core_shell_cylinder.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -132,15 +132,9 @@
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "core_shell_cylinder.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","outer radius","half outer length","half min outer dimension",
-                         "half max outer dimension","half outer diagonal"]
-
-#def ER(radius, thickness, length):
-#    """
-#    Returns the effective radius used in the S*P calculation
-#    """
-#    radius = radius + thickness
-#    length = length + 2 * thickness
-#    ddd = 0.75 * radius * (2 * radius * length + (length + radius) * (length + pi * radius))
-#    return 0.5 * (ddd) ** (1. / 3.)
+effective_radius_type = [
+    "equivalent sphere", "outer radius", "half outer length",
+    "half min outer dimension", "half max outer dimension",
+    "half outer diagonal",
+    ]
 
 def random():
Index: sasmodels/models/core_shell_ellipsoid.c
===================================================================
--- sasmodels/models/core_shell_ellipsoid.c	(revision 3c60146f44ec3507ebbe0e1e7c740546e003f5b2)
+++ sasmodels/models/core_shell_ellipsoid.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -73,11 +73,12 @@
     const double radius_polar_tot = radius_equat_core*x_core + thick_shell*x_polar_shell;
 
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(radius_equat_core, x_core, thick_shell, x_polar_shell);
-    } else if (mode == 2) {
+    case 2: // average outer curvature
         return radius_from_curvature(radius_equat_core, x_core, thick_shell, x_polar_shell);
-    } else if (mode == 3) {
+    case 3: // min outer radius
         return (radius_polar_tot < radius_equat_tot ? radius_polar_tot : radius_equat_tot);
-    } else {
+    case 4: // max outer radius
         return (radius_polar_tot > radius_equat_tot ? radius_polar_tot : radius_equat_tot);
     }
Index: sasmodels/models/core_shell_ellipsoid.py
===================================================================
--- sasmodels/models/core_shell_ellipsoid.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_shell_ellipsoid.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -146,14 +146,8 @@
 source = ["lib/sas_3j1x_x.c", "lib/gauss76.c", "core_shell_ellipsoid.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","average outer curvature", "min outer radius", "max outer radius"]
-
-#def ER(radius_equat_core, x_core, thick_shell, x_polar_shell):
-#    """
-#        Returns the effective radius used in the S*P calculation
-#    """
-#    from .ellipsoid import ER as ellipsoid_ER
-#    polar_outer = radius_equat_core*x_core + thick_shell*x_polar_shell
-#    equat_outer = radius_equat_core + thick_shell
-#    return ellipsoid_ER(polar_outer, equat_outer)
+effective_radius_type = [
+    "equivalent sphere", "average outer curvature",
+    "min outer radius", "max outer radius",
+    ]
 
 def random():
Index: sasmodels/models/core_shell_parallelepiped.c
===================================================================
--- sasmodels/models/core_shell_parallelepiped.c	(revision a94046f61e19c19ba17cb07580d048db76ec896c)
+++ sasmodels/models/core_shell_parallelepiped.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -31,6 +31,6 @@
                    double thick_rim_a, double thick_rim_b, double thick_rim_c)
 {
-    const double volume_paral = form_volume(length_a, length_b, length_c, thick_rim_a, thick_rim_b, thick_rim_c);
-    return cbrt(0.75*volume_paral/M_PI);
+    const double volume = form_volume(length_a, length_b, length_c, thick_rim_a, thick_rim_b, thick_rim_c);
+    return cbrt(volume/M_4PI_3);
 }
 
@@ -45,22 +45,19 @@
 effective_radius(int mode, double length_a, double length_b, double length_c,
                  double thick_rim_a, double thick_rim_b, double thick_rim_c)
-//effective_radius_type = ["equivalent sphere","half outer length_a", "half outer length_b", "half outer length_c",
-//                         "equivalent circular cross-section","half outer ab diagonal","half outer diagonal"]
-// note the core box is A*B*C with slabs ta, tb & tc on each face but missing the corners, though that fact is ignored here
-// in the equvalent sphere option
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(length_a, length_b, length_c, thick_rim_a, thick_rim_b, thick_rim_c);
-    } else if (mode == 2) {
+    case 2: // half outer length a
         return 0.5 * length_a + thick_rim_a;
-    } else if (mode == 3) {
+    case 3: // half outer length b
         return 0.5 * length_b + thick_rim_b;
-    } else if (mode == 4) {
+    case 4: // half outer length c
         return 0.5 * length_c + thick_rim_c;
-    } else if (mode == 5) {
+    case 5: // equivalent circular cross-section
         return radius_from_crosssection(length_a, length_b, thick_rim_a, thick_rim_b);
-    } else if (mode == 6) {
+    case 6: // half outer ab diagonal
         return 0.5*sqrt(square(length_a+ 2.0*thick_rim_a) + square(length_b+ 2.0*thick_rim_b));
-    } else {
+    case 7: // half outer diagonal
         return 0.5*sqrt(square(length_a+ 2.0*thick_rim_a) + square(length_b+ 2.0*thick_rim_b) + square(length_c+ 2.0*thick_rim_c));
     }
Index: sasmodels/models/core_shell_parallelepiped.py
===================================================================
--- sasmodels/models/core_shell_parallelepiped.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_shell_parallelepiped.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -98,5 +98,5 @@
 is the scattering length of the solvent.
 
-.. note:: 
+.. note::
 
    the code actually implements two substitutions: $d(cos\alpha)$ is
@@ -227,19 +227,10 @@
 source = ["lib/gauss76.c", "core_shell_parallelepiped.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","half outer length_a", "half outer length_b", "half outer length_c",
-                         "equivalent circular cross-section","half outer ab diagonal","half outer diagonal"]
-
-#def ER(length_a, length_b, length_c, thick_rim_a, thick_rim_b, thick_rim_c):
-#    """
-#        Return equivalent radius (ER)
-#    """
-#    from .parallelepiped import ER as ER_p
-#
-#    a = length_a + 2*thick_rim_a
-#    b = length_b + 2*thick_rim_b
-#    c = length_c + 2*thick_rim_c
-#    return ER_p(a, b, c)
-
-# VR defaults to 1.0
+effective_radius_type = [
+    "equivalent sphere",
+    "half outer length_a", "half outer length_b", "half outer length_c",
+    "equivalent circular cross-section",
+    "half outer ab diagonal", "half outer diagonal",
+    ]
 
 def random():
Index: sasmodels/models/core_shell_sphere.c
===================================================================
--- sasmodels/models/core_shell_sphere.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/core_shell_sphere.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -8,7 +8,8 @@
 effective_radius(int mode, double radius, double thickness)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // outer radius
         return radius + thickness;
-    } else {
+    case 2: // core radius
         return radius;
     }
Index: sasmodels/models/core_shell_sphere.py
===================================================================
--- sasmodels/models/core_shell_sphere.py	(revision 2cc8aa2a798267dbefb6309f9f93a254fa180b82)
+++ sasmodels/models/core_shell_sphere.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -83,12 +83,4 @@
             sld_core=1.0, sld_shell=2.0, sld_solvent=0.0)
 
-#def ER(radius, thickness):
-#    """
-#        Equivalent radius
-#        @param radius: core radius
-#        @param thickness: shell thickness
-#    """
-#    return radius + thickness
-
 def random():
     outer_radius = 10**np.random.uniform(1.3, 4.3)
Index: sasmodels/models/cylinder.c
===================================================================
--- sasmodels/models/cylinder.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/cylinder.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -28,15 +28,16 @@
 effective_radius(int mode, double radius, double length)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1:
         return radius_from_volume(radius, length);
-    } else if (mode == 2) {
+    case 2:
         return radius;
-    } else if (mode == 3) {
+    case 3:
         return 0.5*length;
-    } else if (mode == 4) {
+    case 4:
         return (radius < 0.5*length ? radius : 0.5*length);
-    } else if (mode == 5) {
+    case 5:
         return (radius > 0.5*length ? radius : 0.5*length);
-    } else {
+    case 6:
         return radius_from_diagonal(radius,length);
     }
Index: sasmodels/models/cylinder.py
===================================================================
--- sasmodels/models/cylinder.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/cylinder.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -139,12 +139,8 @@
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "cylinder.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","radius","half length","half min dimension","half max dimension","half diagonal"]
-
-#def ER(radius, length):
-#    """
-#        Return equivalent radius (ER)
-#    """
-#    ddd = 0.75 * radius * (2 * radius * length + (length + radius) * (length + pi * radius))
-#    return 0.5 * (ddd) ** (1. / 3.)
+effective_radius_type = [
+    "equivalent sphere", "radius",
+    "half length", "half min dimension", "half max dimension", "half diagonal",
+    ]
 
 def random():
Index: sasmodels/models/ellipsoid.c
===================================================================
--- sasmodels/models/ellipsoid.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/ellipsoid.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -34,11 +34,12 @@
 effective_radius(int mode, double radius_polar, double radius_equatorial)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(radius_polar, radius_equatorial);
-    } else if (mode == 2) {
+    case 2: // average curvature
         return radius_from_curvature(radius_polar, radius_equatorial);
-    } else if (mode == 3) {
+    case 3: // min radius
         return (radius_polar < radius_equatorial ? radius_polar : radius_equatorial);
-    } else {
+    case 4: // max radius
         return (radius_polar > radius_equatorial ? radius_polar : radius_equatorial);
     }
Index: sasmodels/models/ellipsoid.py
===================================================================
--- sasmodels/models/ellipsoid.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/ellipsoid.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -169,5 +169,7 @@
 source = ["lib/sas_3j1x_x.c", "lib/gauss76.c", "ellipsoid.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","average curvature", "min radius", "max radius"]
+effective_radius_type = [
+    "equivalent sphere", "average curvature", "min radius", "max radius",
+    ]
 
 def random():
Index: sasmodels/models/elliptical_cylinder.c
===================================================================
--- sasmodels/models/elliptical_cylinder.c	(revision fbaef045776fd0f829e54df2bac54fda45ab5ea3)
+++ sasmodels/models/elliptical_cylinder.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -35,24 +35,23 @@
 static double
 effective_radius(int mode, double radius_minor, double r_ratio, double length)
-//effective_radius_type = ["equivalent sphere","average radius","min radius","max radius",
-//                         "equivalent circular cross-section","half length","half min dimension","half max dimension","half diagonal"]
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(radius_minor, r_ratio, length);
-    } else if (mode == 2) {
+    case 2: // average radius
         return 0.5*radius_minor*(1.0 + r_ratio);
-    } else if (mode == 3) {
+    case 3: // min radius
         return (r_ratio > 1.0 ? radius_minor : r_ratio*radius_minor);
-    } else if (mode == 4) {
+    case 4: // max radius
         return (r_ratio < 1.0 ? radius_minor : r_ratio*radius_minor);
-    } else if (mode == 5) {
+    case 5: // equivalent circular cross-section
         return sqrt(radius_minor*radius_minor*r_ratio);
-    } else if (mode == 6) {
+    case 6: // half length
         return 0.5*length;
-    } else if (mode == 7) {
+    case 7: // half min dimension
         return radius_from_min_dimension(radius_minor,r_ratio,0.5*length);
-    } else if (mode == 8) {
+    case 8: // half max dimension
         return radius_from_max_dimension(radius_minor,r_ratio,0.5*length);
-    } else {
+    case 9: // half diagonal
         return radius_from_diagonal(radius_minor,r_ratio,length);
     }
Index: sasmodels/models/elliptical_cylinder.py
===================================================================
--- sasmodels/models/elliptical_cylinder.py	(revision 2cc8aa2a798267dbefb6309f9f93a254fa180b82)
+++ sasmodels/models/elliptical_cylinder.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -123,22 +123,13 @@
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "elliptical_cylinder.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","average radius","min radius","max radius",
-                         "equivalent circular cross-section","half length","half min dimension","half max dimension","half diagonal"]
+effective_radius_type = [
+    "equivalent sphere", "average radius", "min radius", "max radius",
+    "equivalent circular cross-section",
+    "half length", "half min dimension", "half max dimension", "half diagonal",
+    ]
 
 demo = dict(scale=1, background=0, radius_minor=100, axis_ratio=1.5, length=400.0,
             sld=4.0, sld_solvent=1.0, theta=10.0, phi=20, psi=30,
             theta_pd=10, phi_pd=2, psi_pd=3)
-
-#def ER(radius_minor, axis_ratio, length):
-#    """
-#        Equivalent radius
-#        @param radius_minor: Ellipse minor radius
-#        @param axis_ratio: Ratio of major radius over minor radius
-#        @param length: Length of the cylinder
-#    """
-#    radius = sqrt(radius_minor * radius_minor * axis_ratio)
-#    ddd = 0.75 * radius * (2 * radius * length
-#                           + (length + radius) * (length + pi * radius))
-#    return 0.5 * (ddd) ** (1. / 3.)
 
 def random():
Index: sasmodels/models/fuzzy_sphere.c
===================================================================
--- sasmodels/models/fuzzy_sphere.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/fuzzy_sphere.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -7,7 +7,8 @@
 effective_radius(int mode, double radius, double fuzziness)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // radius
         return radius;
-    } else {
+    case 2: // radius + fuzziness
         return radius + fuzziness;
     }
Index: sasmodels/models/fuzzy_sphere.py
===================================================================
--- sasmodels/models/fuzzy_sphere.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/fuzzy_sphere.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -84,13 +84,5 @@
 source = ["lib/sas_3j1x_x.c","fuzzy_sphere.c"]
 have_Fq = True
-effective_radius_type = ["radius","radius + fuzziness"]
-
-#def ER(radius):
-#    """
-#    Return radius
-#    """
-#    return radius
-
-# VR defaults to 1.0
+effective_radius_type = ["radius", "radius + fuzziness"]
 
 def random():
Index: sasmodels/models/hollow_cylinder.c
===================================================================
--- sasmodels/models/hollow_cylinder.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/hollow_cylinder.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -15,9 +15,17 @@
 }
 
+// TODO: interface to form_volume/shell_volume not yet settled
+static double
+shell_volume(double *total, double radius, double thickness, double length)
+{
+    *total = M_PI*length*square(radius+thickness);
+    return *total - M_PI*length*radius*radius;
+}
+
 static double
 form_volume(double radius, double thickness, double length)
 {
-    double v_shell = M_PI*length*(square(radius+thickness) - radius*radius);
-    return v_shell;
+    double total;
+    return shell_volume(&total, radius, thickness, length);
 }
 
@@ -38,15 +46,16 @@
 effective_radius(int mode, double radius, double thickness, double length)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(radius, thickness, length);
-    } else if (mode == 2) {
+    case 2: // outer radius
         return radius + thickness;
-    } else if (mode == 3) {
+    case 3: // half length
         return 0.5*length;
-    } else if (mode == 4) {
+    case 4: // half outer min dimension
         return (radius + thickness < 0.5*length ? radius + thickness : 0.5*length);
-    } else if (mode == 5) {
+    case 5: // half outer max dimension
         return (radius + thickness > 0.5*length ? radius + thickness : 0.5*length);
-    } else {
+    case 6: // half outer diagonal
         return radius_from_diagonal(radius,thickness,length);
     }
Index: sasmodels/models/hollow_cylinder.py
===================================================================
--- sasmodels/models/hollow_cylinder.py	(revision 2cc8aa2a798267dbefb6309f9f93a254fa180b82)
+++ sasmodels/models/hollow_cylinder.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -100,38 +100,9 @@
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "hollow_cylinder.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","outer radius","half length",
-                         "half outer min dimension","half outer max dimension","half outer diagonal"]
-
-# pylint: disable=W0613
-#def ER(radius, thickness, length):
-#    """
-#    :param radius:      Cylinder core radius
-#    :param thickness:   Cylinder wall thickness
-#    :param length:      Cylinder length
-#    :return:            Effective radius
-#    """
-#    router = radius + thickness
-#    if router == 0 or length == 0:
-#        return 0.0
-#    len1 = router
-#    len2 = length/2.0
-#    term1 = len1*len1*2.0*len2/2.0
-#    term2 = 1.0 + (len2/len1)*(1.0 + 1/len2/2.0)*(1.0 + pi*len1/len2/2.0)
-#    ddd = 3.0*term1*term2
-#    diam = pow(ddd, (1.0/3.0))
-#    return diam
-
-def VR(radius, thickness, length):
-    """
-    :param radius:      Cylinder radius
-    :param thickness:   Cylinder wall thickness
-    :param length:      Cylinder length
-    :return:            Volf ratio for P(q)*S(q)
-    """
-    router = radius + thickness
-    vol_core = pi*radius*radius*length
-    vol_total = pi*router*router*length
-    vol_shell = vol_total - vol_core
-    return vol_total, vol_shell
+effective_radius_type = [
+    "equivalent sphere", "outer radius", "half length",
+    "half outer min dimension", "half outer max dimension",
+    "half outer diagonal",
+    ]
 
 def random():
Index: sasmodels/models/hollow_rectangular_prism.c
===================================================================
--- sasmodels/models/hollow_rectangular_prism.c	(revision a94046f61e19c19ba17cb07580d048db76ec896c)
+++ sasmodels/models/hollow_rectangular_prism.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -1,4 +1,5 @@
+// TODO: interface to form_volume/shell_volume not yet settled
 static double
-form_volume(double length_a, double b2a_ratio, double c2a_ratio, double thickness)
+shell_volume(double *total, double length_a, double b2a_ratio, double c2a_ratio, double thickness)
 {
     double length_b = length_a * b2a_ratio;
@@ -8,7 +9,13 @@
     double c_core = length_c - 2.0*thickness;
     double vol_core = a_core * b_core * c_core;
-    double vol_total = length_a * length_b * length_c;
-    double vol_shell = vol_total - vol_core;
-    return vol_shell;
+    *total = length_a * length_b * length_c;
+    return *total - vol_core;
+}
+
+static double
+form_volume(double length_a, double b2a_ratio, double c2a_ratio, double thickness)
+{
+    double total;
+    return shell_volume(&total, length_a, b2a_ratio, c2a_ratio, thickness);
 }
 
@@ -19,17 +26,18 @@
 // NOTE length_a is external dimension!
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return cbrt(0.75*cube(length_a)*b2a_ratio*c2a_ratio/M_PI);
-    } else if (mode == 2) {
+    case 2: // half length_a
         return 0.5 * length_a;
-    } else if (mode == 3) {
+    case 3: // half length_b
         return 0.5 * length_a*b2a_ratio;
-    } else if (mode == 4) {
+    case 4: // half length_c
         return 0.5 * length_a*c2a_ratio;
-    } else if (mode == 5) {
+    case 5: // equivalent outer circular cross-section
         return length_a*sqrt(b2a_ratio/M_PI);
-    } else if (mode == 6) {
+    case 6: // half ab diagonal
         return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio)));
-    } else {
+    case 7: // half diagonal
         return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio) + square(c2a_ratio)));
     }
Index: sasmodels/models/hollow_rectangular_prism.py
===================================================================
--- sasmodels/models/hollow_rectangular_prism.py	(revision a94046f61e19c19ba17cb07580d048db76ec896c)
+++ sasmodels/models/hollow_rectangular_prism.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -149,35 +149,9 @@
 source = ["lib/gauss76.c", "hollow_rectangular_prism.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","half length_a", "half length_b", "half length_c",
-                         "equivalent outer circular cross-section","half ab diagonal","half diagonal"]
-
-#def ER(length_a, b2a_ratio, c2a_ratio, thickness):
-#    """
-#    Return equivalent radius (ER)
-#    thickness parameter not used
-#    """
-#    b_side = length_a * b2a_ratio
-#    c_side = length_a * c2a_ratio
-#
-#    # surface average radius (rough approximation)
-#    surf_rad = sqrt(length_a * b_side / pi)
-#
-#    ddd = 0.75 * surf_rad * (2 * surf_rad * c_side + (c_side + surf_rad) * (c_side + pi * surf_rad))
-#    return 0.5 * (ddd) ** (1. / 3.)
-
-def VR(length_a, b2a_ratio, c2a_ratio, thickness):
-    """
-    Return shell volume and total volume
-    """
-    b_side = length_a * b2a_ratio
-    c_side = length_a * c2a_ratio
-    a_core = length_a - 2.0*thickness
-    b_core = b_side - 2.0*thickness
-    c_core = c_side - 2.0*thickness
-    vol_core = a_core * b_core * c_core
-    vol_total = length_a * b_side * c_side
-    vol_shell = vol_total - vol_core
-    return vol_total, vol_shell
-
+effective_radius_type = [
+    "equivalent sphere", "half length_a", "half length_b", "half length_c",
+    "equivalent outer circular cross-section",
+    "half ab diagonal", "half diagonal",
+    ]
 
 def random():
Index: sasmodels/models/hollow_rectangular_prism_thin_walls.c
===================================================================
--- sasmodels/models/hollow_rectangular_prism_thin_walls.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/hollow_rectangular_prism_thin_walls.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -1,27 +1,38 @@
+// TODO: interface to form_volume/shell_volume not yet settled
 static double
-form_volume(double length_a, double b2a_ratio, double c2a_ratio)
+shell_volume(double *total, double length_a, double b2a_ratio, double c2a_ratio)
 {
     double length_b = length_a * b2a_ratio;
     double length_c = length_a * c2a_ratio;
     double vol_shell = 2.0 * (length_a*length_b + length_a*length_c + length_b*length_c);
+    *total = length_a * length_b * length_c;
     return vol_shell;
 }
 
 static double
+form_volume(double length_a, double b2a_ratio, double c2a_ratio)
+{
+    double total;
+    return shell_volume(&total, length_a, b2a_ratio, c2a_ratio);
+}
+
+
+static double
 effective_radius(int mode, double length_a, double b2a_ratio, double c2a_ratio)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return cbrt(0.75*cube(length_a)*b2a_ratio*c2a_ratio/M_PI);
-    } else if (mode == 2) {
+    case 2: // half length_a
         return 0.5 * length_a;
-    } else if (mode == 3) {
+    case 3: // half length_b
         return 0.5 * length_a*b2a_ratio;
-    } else if (mode == 4) {
+    case 4: // half length_c
         return 0.5 * length_a*c2a_ratio;
-    } else if (mode == 5) {
+    case 5: // equivalent outer circular cross-section
         return length_a*sqrt(b2a_ratio/M_PI);
-    } else if (mode == 6) {
+    case 6: // half ab diagonal
         return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio)));
-    } else {
+    case 7: // half diagonal
         return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio) + square(c2a_ratio)));
     }
Index: sasmodels/models/hollow_rectangular_prism_thin_walls.py
===================================================================
--- sasmodels/models/hollow_rectangular_prism_thin_walls.py	(revision 6e7d7b66fc1ec444c4e30ac47db7bbcdc554ecc3)
+++ sasmodels/models/hollow_rectangular_prism_thin_walls.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -109,29 +109,9 @@
 source = ["lib/gauss76.c", "hollow_rectangular_prism_thin_walls.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","half length_a", "half length_b", "half length_c",
-                         "equivalent outer circular cross-section","half ab diagonal","half diagonal"]
-
-#def ER(length_a, b2a_ratio, c2a_ratio):
-#    """
-#        Return equivalent radius (ER)
-#    """
-#    b_side = length_a * b2a_ratio
-#    c_side = length_a * c2a_ratio
-#
-#    # surface average radius (rough approximation)
-#    surf_rad = sqrt(length_a * b_side / pi)
-#
-#    ddd = 0.75 * surf_rad * (2 * surf_rad * c_side + (c_side + surf_rad) * (c_side + pi * surf_rad))
-#    return 0.5 * (ddd) ** (1. / 3.)
-
-def VR(length_a, b2a_ratio, c2a_ratio):
-    """
-        Return shell volume and total volume
-    """
-    b_side = length_a * b2a_ratio
-    c_side = length_a * c2a_ratio
-    vol_total = length_a * b_side * c_side
-    vol_shell = 2.0 * (length_a*b_side + length_a*c_side + b_side*c_side)
-    return vol_shell, vol_total
+effective_radius_type = [
+    "equivalent sphere", "half length_a", "half length_b", "half length_c",
+    "equivalent outer circular cross-section",
+    "half ab diagonal", "half diagonal",
+    ]
 
 
Index: sasmodels/models/mono_gauss_coil.c
===================================================================
--- sasmodels/models/mono_gauss_coil.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/mono_gauss_coil.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -7,25 +7,54 @@
 effective_radius(int mode, double rg)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // R_g
         return rg;
-    } else if (mode == 2) {
+    case 2: // 2R_g
         return 2.0*rg;
-    } else if (mode == 3) {
+    case 3: // 3R_g
         return 3.0*rg;
-    } else {
+    case 4: // (5/3)^0.5*R_g
         return sqrt(5.0/3.0)*rg;
     }
 }
 
+static double
+gauss_coil(double qr)
+{
+    const double x = qr*qr;
+
+    // Use series expansion at low q for higher accuracy. We could use
+    // smaller polynomials if we sacrifice some digits of precision or
+    // introduce an additional series expansion around x == 1.
+    // See explore/precision.py, gauss_coil function.
+#if FLOAT_SIZE>4 // DOUBLE_PRECISION
+    // For double precision: use O(5) Pade with 0.5 cutoff (10 mad + 1 divide)
+    if (x < 0.5) {
+        // PadeApproximant[2*Exp[-x^2] + x^2-1)/x^4, {x, 0, 8}]
+        const double A1=1./12., A2=2./99., A3=1./2640., A4=1./23760., A5=-1./1995840.;
+        const double B1=5./12., B2=5./66., B3=1./132., B4=1./2376., B5=1./95040.;
+        return (((((A5*x + A4)*x + A3)*x + A2)*x + A1)*x + 1.)
+                / (((((B5*x + B4)*x + B3)*x + B2)*x + B1)*x + 1.);
+    }
+#else
+    // For single precision: use O(7) Taylor with 0.8 cutoff (7 mad)
+    if (x < 0.8) {
+        const double C0 = +1.;
+        const double C1 = -1./3.;
+        const double C2 = +1./12.;
+        const double C3 = -1./60.;
+        const double C4 = +1./360.;
+        const double C5 = -1./2520.;
+        const double C6 = +1./20160.;
+        const double C7 = -1./181440.;
+        return ((((((C7*x + C6)*x + C5)*x + C4)*x + C3)*x + C2)*x + C1)*x + C0;
+    }
+#endif
+
+    return 2.0 * (expm1(-x) + x)/(x*x);
+}
+
 double Iq(double q, double i_zero, double rg)
 {
-    const double uarg = square(q*rg);
-    const double inten;
-    if (q == 0) {
-        inten = i_zero;
-    } else {
-        inten = 2.0*i_zero * (exp(-uarg) + uarg - 1.0)/square(uarg);
-    }
-
-    return inten;
+    return i_zero * gauss_coil(q*rg);
 }
Index: sasmodels/models/mono_gauss_coil.py
===================================================================
--- sasmodels/models/mono_gauss_coil.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/mono_gauss_coil.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -71,22 +71,10 @@
     ["rg", "Ang", 75.0, [0.0, inf], "volume", "Radius of gyration"],
     ]
+# pylint: enable=bad-whitespace, line-too-long
 
 source = ["mono_gauss_coil.c"]
 have_Fq = False
-effective_radius_type = ["R_g","2R_g","3R_g","(5/3)^0.5*R_g"]
+effective_radius_type = ["R_g", "2R_g", "3R_g", "(5/3)^0.5*R_g"]
 
-
-# pylint: enable=bad-whitespace, line-too-long
-
-## NB: Scale and Background are implicit parameters on every model
-#def Iq(q, i_zero, rg):
-#    # pylint: disable = missing-docstring
-#    z = (q * rg)**2
-#
-#    with errstate(invalid='ignore'):
-#        inten = (i_zero * 2.0) * (exp(-z) + z - 1.0)/z**2
-#        inten[q == 0] = i_zero
-#    return inten
-#Iq.vectorized = True # Iq accepts an array of q values
 
 def random():
Index: sasmodels/models/multilayer_vesicle.c
===================================================================
--- sasmodels/models/multilayer_vesicle.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/multilayer_vesicle.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -50,4 +50,5 @@
 effective_radius(int mode, double radius, double thick_shell, double thick_solvent, double fp_n_shells)
 {
+    // case 1: outer radius
     return radius + fp_n_shells*thick_shell + (fp_n_shells - 1.0)*thick_solvent;
 }
Index: sasmodels/models/multilayer_vesicle.py
===================================================================
--- sasmodels/models/multilayer_vesicle.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/multilayer_vesicle.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -148,8 +148,4 @@
 effective_radius_type = ["outer radius"]
 
-#def ER(radius, thick_shell, thick_solvent, n_shells):
-#    n_shells = int(n_shells+0.5)
-#    return radius + n_shells * (thick_shell + thick_solvent) - thick_solvent
-
 def random():
     volfraction = 10**np.random.uniform(-3, -0.5)  # scale from 0.1% to 30%
Index: sasmodels/models/onion.c
===================================================================
--- sasmodels/models/onion.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/onion.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -30,16 +30,5 @@
 
 static double
-form_volume(double radius_core, double n_shells, double thickness[])
-{
-  int n = (int)(n_shells+0.5);
-  double r = radius_core;
-  for (int i=0; i < n; i++) {
-    r += thickness[i];
-  }
-  return M_4PI_3*cube(r);
-}
-
-static double
-effective_radius(int mode, double radius_core, double n_shells, double thickness[])
+outer_radius(double radius_core, double n_shells, double thickness[])
 {
   int n = (int)(n_shells+0.5);
@@ -49,4 +38,17 @@
   }
   return r;
+}
+
+static double
+form_volume(double radius_core, double n_shells, double thickness[])
+{
+  return M_4PI_3*cube(outer_radius(radius_core, n_shells, thickness));
+}
+
+static double
+effective_radius(int mode, double radius_core, double n_shells, double thickness[])
+{
+  // case 1: outer radius
+  return outer_radius(radius_core, n_shells, thickness);
 }
 
Index: sasmodels/models/onion.py
===================================================================
--- sasmodels/models/onion.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/onion.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -367,9 +367,4 @@
     return np.asarray(z), np.asarray(rho)
 
-#def ER(radius_core, n_shells, thickness):
-#    """Effective radius"""
-#    n = int(n_shells[0]+0.5)
-#    return np.sum(thickness[:n], axis=0) + radius_core
-
 demo = {
     "sld_solvent": 2.2,
Index: sasmodels/models/parallelepiped.c
===================================================================
--- sasmodels/models/parallelepiped.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/parallelepiped.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -8,21 +8,21 @@
 effective_radius(int mode, double length_a, double length_b, double length_c)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return cbrt(0.75*length_a*length_b*length_c/M_PI);
-    } else if (mode == 2) {
+    case 2: // half length_a
         return 0.5 * length_a;
-    } else if (mode == 3) {
+    case 3: // half length_b
         return 0.5 * length_b;
-    } else if (mode == 4) {
+    case 4: // half length_c
         return 0.5 * length_c;
-    } else if (mode == 5) {
+    case 5: // equivalent circular cross-section
         return sqrt(length_a*length_b/M_PI);
-    } else if (mode == 6) {
+    case 6: // half ab diagonal
         return 0.5*sqrt(length_a*length_a + length_b*length_b);
-    } else {
+    case 7: // half diagonal
         return 0.5*sqrt(length_a*length_a + length_b*length_b + length_c*length_c);
     }
 }
-
 
 static void
Index: sasmodels/models/parallelepiped.py
===================================================================
--- sasmodels/models/parallelepiped.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/parallelepiped.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -140,5 +140,5 @@
    ensuring that the inequality $A < B < C$ is not violated,  The calculation
    will not report an error, but the results may be not correct.
-   
+
 .. _parallelepiped-orientation:
 
@@ -231,25 +231,8 @@
 source = ["lib/gauss76.c", "parallelepiped.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","half length_a", "half length_b", "half length_c",
-                         "equivalent circular cross-section","half ab diagonal","half diagonal"]
-
-#def ER(length_a, length_b, length_c):
-#    """
-#    Return effective radius (ER) for P(q)*S(q)
-#    """
-#    # now that axes can be in any size order, need to sort a,b,c
-#    # where a~b and c is either much smaller or much larger
-#    abc = np.vstack((length_a, length_b, length_c))
-#    abc = np.sort(abc, axis=0)
-#    selector = (abc[1] - abc[0]) > (abc[2] - abc[1])
-#    length = np.where(selector, abc[0], abc[2])
-#    # surface average radius (rough approximation)
-#    radius = sqrt(np.where(~selector, abc[0]*abc[1], abc[1]*abc[2]) / pi)
-#
-#    ddd = 0.75 * radius * (2*radius*length + (length + radius)*(length + pi*radius))
-#    return 0.5 * (ddd) ** (1. / 3.)
-
-# VR defaults to 1.0
-
+effective_radius_type = [
+    "equivalent sphere", "half length_a", "half length_b", "half length_c",
+    "equivalent circular cross-section", "half ab diagonal", "half diagonal",
+    ]
 
 def random():
Index: sasmodels/models/pringle.c
===================================================================
--- sasmodels/models/pringle.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/pringle.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -107,7 +107,8 @@
 effective_radius(int mode, double radius, double thickness, double alpha, double beta)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return cbrt(0.75*radius*radius*thickness);
-    } else {
+    case 2: // radius
         return radius;
     }
Index: sasmodels/models/pringle.py
===================================================================
--- sasmodels/models/pringle.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/pringle.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -72,15 +72,7 @@
 
 
-source = ["lib/polevl.c", "lib/sas_J0.c", "lib/sas_J1.c", \
+source = ["lib/polevl.c", "lib/sas_J0.c", "lib/sas_J1.c",
           "lib/sas_JN.c", "lib/gauss76.c", "pringle.c"]
-effective_radius_type = ["equivalent sphere","radius"]
-
-#def ER(radius, thickness, alpha, beta):
-#    """
-#    Return equivalent radius (ER)
-#    """
-#    ddd = 0.75 * radius * (2 * radius * thickness + (thickness + radius) \
-#                           * (thickness + pi * radius))
-#    return 0.5 * (ddd) ** (1. / 3.)
+effective_radius_type = ["equivalent sphere", "radius"]
 
 def random():
Index: sasmodels/models/raspberry.c
===================================================================
--- sasmodels/models/raspberry.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/raspberry.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -17,7 +17,8 @@
 effective_radius(int mode, double radius_lg, double radius_sm, double penetration)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // radius_large
         return radius_lg;
-    } else {
+    case 2: // radius_outer
         return radius_lg + 2.0*radius_sm - penetration;
     }
Index: sasmodels/models/raspberry.py
===================================================================
--- sasmodels/models/raspberry.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/raspberry.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -152,5 +152,5 @@
 
 source = ["lib/sas_3j1x_x.c", "raspberry.c"]
-effective_radius_type = ["radius_large","radius_outer"]
+effective_radius_type = ["radius_large", "radius_outer"]
 
 def random():
Index: sasmodels/models/rectangular_prism.c
===================================================================
--- sasmodels/models/rectangular_prism.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/rectangular_prism.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -8,17 +8,18 @@
 effective_radius(int mode, double length_a, double b2a_ratio, double c2a_ratio)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return cbrt(0.75*cube(length_a)*b2a_ratio*c2a_ratio/M_PI);
-    } else if (mode == 2) {
+    case 2: // half length_a
         return 0.5 * length_a;
-    } else if (mode == 3) {
+    case 3: // half length_b
         return 0.5 * length_a*b2a_ratio;
-    } else if (mode == 4) {
+    case 4: // half length_c
         return 0.5 * length_a*c2a_ratio;
-    } else if (mode == 5) {
+    case 5: // equivalent circular cross-section
         return length_a*sqrt(b2a_ratio/M_PI);
-    } else if (mode == 6) {
+    case 6: // half ab diagonal
         return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio)));
-    } else {
+    case 7: // half diagonal
         return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio) + square(c2a_ratio)));
     }
Index: sasmodels/models/rectangular_prism.py
===================================================================
--- sasmodels/models/rectangular_prism.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/rectangular_prism.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -136,19 +136,8 @@
 source = ["lib/gauss76.c", "rectangular_prism.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","half length_a", "half length_b", "half length_c",
-                         "equivalent circular cross-section","half ab diagonal","half diagonal"]
-
-#def ER(length_a, b2a_ratio, c2a_ratio):
-#    """
-#        Return equivalent radius (ER)
-#    """
-#    b_side = length_a * b2a_ratio
-#    c_side = length_a * c2a_ratio
-#
-#    # surface average radius (rough approximation)
-#    surf_rad = sqrt(length_a * b_side / pi)
-#
-#    ddd = 0.75 * surf_rad * (2 * surf_rad * c_side + (c_side + surf_rad) * (c_side + pi * surf_rad))
-#    return 0.5 * (ddd) ** (1. / 3.)
+effective_radius_type = [
+    "equivalent sphere", "half length_a", "half length_b", "half length_c",
+    "equivalent circular cross-section", "half ab diagonal", "half diagonal",
+    ]
 
 def random():
Index: sasmodels/models/sphere.c
===================================================================
--- sasmodels/models/sphere.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/sphere.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -7,4 +7,5 @@
 effective_radius(int mode, double radius)
 {
+    // case 1: radius
     return radius;
 }
Index: sasmodels/models/sphere.py
===================================================================
--- sasmodels/models/sphere.py	(revision 2cc8aa2a798267dbefb6309f9f93a254fa180b82)
+++ sasmodels/models/sphere.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -71,12 +71,4 @@
 effective_radius_type = ["radius"]
 
-#def ER(radius):
-#    """
-#    Return equivalent radius (ER)
-#    """
-#    return radius
-
-# VR defaults to 1.0
-
 def random():
     radius = 10**np.random.uniform(1.3, 4)
Index: sasmodels/models/spherical_sld.c
===================================================================
--- sasmodels/models/spherical_sld.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/spherical_sld.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -1,17 +1,4 @@
-static double form_volume(
-    double fp_n_shells,
-    double thickness[],
-    double interface[])
-{
-    int n_shells= (int)(fp_n_shells + 0.5);
-    double r = 0.0;
-    for (int i=0; i < n_shells; i++) {
-        r += thickness[i] + interface[i];
-    }
-    return M_4PI_3*cube(r);
-}
-
 static double
-effective_radius(int mode, double fp_n_shells, double thickness[], double interface[])
+outer_radius(double fp_n_shells, double thickness[], double interface[])
 {
     int n_shells= (int)(fp_n_shells + 0.5);
@@ -21,4 +8,19 @@
     }
     return r;
+}
+
+static double form_volume(
+    double fp_n_shells,
+    double thickness[],
+    double interface[])
+{
+    return M_4PI_3*cube(outer_radius(fp_n_shells, thickness, interface));
+}
+
+static double
+effective_radius(int mode, double fp_n_shells, double thickness[], double interface[])
+{
+    // case 1: outer radius
+    return outer_radius(fp_n_shells, thickness, interface);
 }
 
Index: sasmodels/models/spherical_sld.py
===================================================================
--- sasmodels/models/spherical_sld.py	(revision e8eff7b6b415fc238f4ec9a17539b740a8df0411)
+++ sasmodels/models/spherical_sld.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -22,11 +22,11 @@
 
     0: erf($\nu z$)
-    
+
     1: Rpow($z^\nu$)
-    
+
     2: Lpow($z^\nu$)
-    
+
     3: Rexp($-\nu z$)
-    
+
     4: Lexp($-\nu z$)
 
@@ -270,12 +270,4 @@
 
 
-#def ER(n_shells, thickness, interface):
-#    """Effective radius"""
-#    n_shells = int(n_shells + 0.5)
-#    total = (np.sum(thickness[:n_shells], axis=1)
-#             + np.sum(interface[:n_shells], axis=1))
-#    return total
-
-
 demo = {
     "n_shells": 5,
Index: sasmodels/models/triaxial_ellipsoid.c
===================================================================
--- sasmodels/models/triaxial_ellipsoid.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/triaxial_ellipsoid.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -30,9 +30,10 @@
 effective_radius(int mode, double radius_equat_minor, double radius_equat_major, double radius_polar)
 {
-    if (mode == 1) {
+    switch (mode) {
+    case 1: // equivalent sphere
         return radius_from_volume(radius_equat_minor,radius_equat_major, radius_polar);
-    } else if (mode == 2) {
+    case 2: // min radius
         return radius_from_min_dimension(radius_equat_minor,radius_equat_major, radius_polar);
-    } else {
+    case 3: // max radius
         return radius_from_max_dimension(radius_equat_minor,radius_equat_major, radius_polar);
     }
Index: sasmodels/models/triaxial_ellipsoid.py
===================================================================
--- sasmodels/models/triaxial_ellipsoid.py	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/triaxial_ellipsoid.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -158,20 +158,5 @@
 source = ["lib/sas_3j1x_x.c", "lib/gauss76.c", "triaxial_ellipsoid.c"]
 have_Fq = True
-effective_radius_type = ["equivalent sphere","min radius", "max radius"]
-
-def ER(radius_equat_minor, radius_equat_major, radius_polar):
-    """
-    Returns the effective radius used in the S*P calculation
-    """
-    from .ellipsoid import ER as ellipsoid_ER
-
-    # now that radii can be in any size order, radii need sorting a,b,c
-    # where a~b and c is either much smaller or much larger
-    radii = np.vstack((radius_equat_major, radius_equat_minor, radius_polar))
-    radii = np.sort(radii, axis=0)
-    selector = (radii[1] - radii[0]) > (radii[2] - radii[1])
-    polar = np.where(selector, radii[0], radii[2])
-    equatorial = np.sqrt(np.where(~selector, radii[0]*radii[1], radii[1]*radii[2]))
-    return ellipsoid_ER(polar, equatorial)
+effective_radius_type = ["equivalent sphere", "min radius", "max radius"]
 
 def random():
Index: sasmodels/models/vesicle.c
===================================================================
--- sasmodels/models/vesicle.c	(revision d277229f02cb7770fbb4ea2ad506d96bed705d75)
+++ sasmodels/models/vesicle.c	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -1,7 +1,17 @@
+// TODO: interface to form_volume/shell_volume not yet settled
+static double
+shell_volume(double *total, double radius, double thickness)
+{
+    //note that for the vesicle model, the volume is ONLY the shell volume
+    *total = M_4PI_3 * cube(radius+thickness);
+    return *total - M_4PI_3*cube(radius);
+}
+
 static double
 form_volume(double radius, double thickness)
 {
     //note that for the vesicle model, the volume is ONLY the shell volume
-    return M_4PI_3*(cube(radius+thickness) - cube(radius));
+    double total;
+    return shell_volume(&total, radius, thickness);
 }
 
@@ -9,4 +19,5 @@
 effective_radius(int mode, double radius, double thickness)
 {
+    // case 1: outer radius
     return radius + thickness;
 }
Index: sasmodels/models/vesicle.py
===================================================================
--- sasmodels/models/vesicle.py	(revision 2cc8aa2a798267dbefb6309f9f93a254fa180b82)
+++ sasmodels/models/vesicle.py	(revision ee60aa768b29c58e20b19e94578ec3a94a0889d7)
@@ -102,28 +102,4 @@
 effective_radius_type = ["outer radius"]
 
-#def ER(radius, thickness):
-#    '''
-#    returns the effective radius used in the S*P calculation
-#
-#    :param radius: core radius
-#    :param thickness: shell thickness
-#    '''
-#    return radius + thickness
-
-def VR(radius, thickness):
-    '''
-    returns the volumes of the shell and of the whole sphere including the
-    core plus shell - is used to normalize when including polydispersity.
-
-    :param radius: core radius
-    :param thickness: shell thickness
-    :return whole: volume of core and shell
-    :return whole-core: volume of the shell
-    '''
-
-    whole = 4./3. * pi * (radius + thickness)**3
-    core = 4./3. * pi * radius**3
-    return whole, whole - core
-
 def random():
     total_radius = 10**np.random.uniform(1.3, 5)
