# Changeset db472a5 in sasmodels

Ignore:
Timestamp:
Sep 18, 2018 9:19:29 AM (17 months ago)
Branches:
master, core_shell_microgels, magnetic_model, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
Children:
fb7c176, c5f7aa9
Parents:
0159b5e (diff), d299327 (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.
Message:

Merge branch 'beta_approx_new_R_eff' into beta_approx

Files:
68 edited

Unmodified
Removed

• ## explore/precision.py

 r2a7e20e ) 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 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.; /(((((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.;
• ## sasmodels/core.py

 r71b751d model = load_model("cylinder@hardsphere*sphere") actual = [p.name for p in model.info.parameters.kernel_parameters] target = ("sld sld_solvent radius length theta phi volfraction" " beta_mode A_sld A_sld_solvent A_radius").split() target = ("sld sld_solvent radius length theta phi" " radius_effective volfraction " " structure_factor_mode radius_effective_mode" " A_sld A_sld_solvent A_radius").split() assert target == actual, "%s != %s"%(target, actual)
• ## sasmodels/generate.py

 rc88f983 for p in partable.kernel_parameters)) # Define the function calls call_effective_radius = "#define CALL_EFFECTIVE_RADIUS(mode, v) 0.0" if partable.form_volume_parameters: refs = _call_pars("_v.", partable.form_volume_parameters) call_volume = "#define CALL_VOLUME(_v) form_volume(%s)"%(",".join(refs)) if model_info.effective_radius_type: call_effective_radius = "#define CALL_EFFECTIVE_RADIUS(mode, _v) effective_radius(mode, %s)"%(",".join(refs)) else: # Model doesn't have volume.  We could make the kernel run a little call_volume = "#define CALL_VOLUME(v) 1.0" source.append(call_volume) source.append(call_effective_radius) model_refs = _call_pars("_v.", partable.iq_parameters)
• ## sasmodels/kernel.py

 r108e70e #define NEED_EXPM1 #define NEED_TGAMMA #define NEED_CBRT // expf missing from windows? #define expf exp #  define pown(a,b) pow(a,b) #endif // !USE_OPENCL #if defined(NEED_CBRT) #define cbrt(_x) pow(_x, 0.33333333333333333333333) #endif #if defined(NEED_EXPM1)

• ## sasmodels/kernelcl.py

 rc036ddb # at this point, so instead using 32, which is good on the set of # architectures tested so far. extra_q = 3  # total weight, weighted volume and weighted radius if self.is_2d: # Note: 16 rather than 15 because result is 1 longer than input. width = ((self.nq+16)//16)*16 width = ((self.nq+15+extra_q)//16)*16 self.q = np.empty((width, 2), dtype=dtype) self.q[:self.nq, 0] = q_vectors[0] self.q[:self.nq, 1] = q_vectors[1] else: # Note: 32 rather than 31 because result is 1 longer than input. width = ((self.nq+32)//32)*32 width = ((self.nq+31+extra_q)//32)*32 self.q = np.empty(width, dtype=dtype) self.q[:self.nq] = q_vectors[0] self.dim = '2d' if q_input.is_2d else '1d' # leave room for f1/f2 results in case we need to compute beta for 1d models num_returns = 1 if self.dim == '2d' else 2  # # plus 1 for the normalization value self.result = np.empty((q_input.nq+1)*num_returns, dtype) nout = 2 if self.info.have_Fq and self.dim == '1d' else 1 # plus 3 weight, volume, radius self.result = np.empty(q_input.nq*nout + 3, self.dtype) # Inputs and outputs for each kernel call self.result_b = cl.Buffer(self.queue.context, mf.READ_WRITE, q_input.global_size[0] * num_returns * dtype.itemsize) q_input.global_size[0] * nout * dtype.itemsize) self.q_input = q_input # allocated by GpuInput above else np.float32)  # will never get here, so use np.float32 def Iq(self, call_details, values, cutoff, magnetic): # type: (CallDetails, np.ndarray, np.ndarray, float, bool) -> np.ndarray self._call_kernel(call_details, values, cutoff, magnetic) #print("returned",self.q_input.q, self.result) pd_norm = self.result[self.q_input.nq] scale = values[0]/(pd_norm if pd_norm != 0.0 else 1.0) background = values[1] #print("scale",scale,background) return scale*self.result[:self.q_input.nq] + background __call__ = Iq def beta(self, call_details, values, cutoff, magnetic): # type: (CallDetails, np.ndarray, np.ndarray, float, bool) -> np.ndarray if self.dim == '2d': raise NotImplementedError("beta not yet supported for 2D") self._call_kernel(call_details, values, cutoff, magnetic) w_norm = self.result[2*self.q_input.nq + 1] pd_norm = self.result[self.q_input.nq] if w_norm == 0.: w_norm = 1. F2 = self.result[:self.q_input.nq]/w_norm F1 = self.result[self.q_input.nq+1:2*self.q_input.nq+1]/w_norm volume_avg = pd_norm/w_norm return F1, F2, volume_avg def _call_kernel(self, call_details, values, cutoff, magnetic): def _call_kernel(self, call_details, values, cutoff, magnetic, effective_radius_type): # type: (CallDetails, np.ndarray, np.ndarray, float, bool) -> np.ndarray context = self.queue.context details_b, values_b, self.q_input.q_b, self.result_b, self.real(cutoff), np.uint32(effective_radius_type), ] #print("Calling OpenCL")
• ## sasmodels/kerneldll.py

 r2f8cbb9 # int, int, int, int*, double*, double*, double*, double*, double argtypes = [ct.c_int32]*3 + [ct.c_void_p]*4 + [float_type] argtypes = [ct.c_int32]*3 + [ct.c_void_p]*4 + [float_type, ct.c_int32] names = [generate.kernel_name(self.info, variant) for variant in ("Iq", "Iqxy", "Imagnetic")] self.dim = '2d' if q_input.is_2d else '1d' # leave room for f1/f2 results in case we need to compute beta for 1d models num_returns = 1 if self.dim == '2d' else 2  # # plus 1 for the normalization value self.result = np.empty((q_input.nq+1)*num_returns, self.dtype) nout = 2 if self.info.have_Fq else 1 # plus 3 weight, volume, radius self.result = np.empty(q_input.nq*nout + 3, self.dtype) self.real = (np.float32 if self.q_input.dtype == generate.F32 else np.float64 if self.q_input.dtype == generate.F64 else np.float128) def Iq(self, call_details, values, cutoff, magnetic): # type: (CallDetails, np.ndarray, np.ndarray, float, bool) -> np.ndarray self._call_kernel(call_details, values, cutoff, magnetic) #print("returned",self.q_input.q, self.result) pd_norm = self.result[self.q_input.nq] scale = values[0]/(pd_norm if pd_norm != 0.0 else 1.0) background = values[1] #print("scale",scale,background) return scale*self.result[:self.q_input.nq] + background __call__ = Iq def beta(self, call_details, values, cutoff, magnetic): # type: (CallDetails, np.ndarray, np.ndarray, float, bool) -> np.ndarray if self.dim == '2d': raise NotImplementedError("beta not yet supported for 2D") self._call_kernel(call_details, values, cutoff, magnetic) w_norm = self.result[2*self.q_input.nq + 1] pd_norm = self.result[self.q_input.nq] if w_norm == 0.: w_norm = 1. F2 = self.result[:self.q_input.nq]/w_norm F1 = self.result[self.q_input.nq+1:2*self.q_input.nq+1]/w_norm volume_avg = pd_norm/w_norm return F1, F2, volume_avg def _call_kernel(self, call_details, values, cutoff, magnetic): def _call_kernel(self, call_details, values, cutoff, magnetic, effective_radius_type): # type: (CallDetails, np.ndarray, np.ndarray, float, bool, int) -> np.ndarray kernel = self.kernel[1 if magnetic else 0] args = [ self.result.ctypes.data,   # results self.real(cutoff), # cutoff effective_radius_type, # cutoff ] #print("Calling DLL")

• ## sasmodels/modelinfo.py

 rc88f983 info.structure_factor = getattr(kernel_module, 'structure_factor', False) # TODO: find Fq by inspection info.effective_radius_type = getattr(kernel_module, 'effective_radius_type', None) info.have_Fq = getattr(kernel_module, 'have_Fq', False) info.profile_axes = getattr(kernel_module, 'profile_axes', ['x', 'y']) info.source = getattr(kernel_module, 'source', []) info.c_code = getattr(kernel_module, 'c_code', None) info.effective_radius = getattr(kernel_module, 'effective_radius', None) info.ER = None  # CRUFT info.VR = None  # CRUFT # TODO: check the structure of the tests info.tests = getattr(kernel_module, 'tests', []) info.ER = getattr(kernel_module, 'ER', None) # type: ignore info.VR = getattr(kernel_module, 'VR', None) # type: ignore info.form_volume = getattr(kernel_module, 'form_volume', None) # type: ignore info.Iq = getattr(kernel_module, 'Iq', None) # type: ignore #: use those functions.  Form factors are indicated by providing #: an :attr:ER function. effective_radius_type = None   # type: List[str] #: Returns the occupied volume and the total volume for each parameter set. #: See :attr:ER for details on the parameters. source = None           # type: List[str] #: The set of tests that must pass.  The format of the tests is described #: each parameter set.  Multiplicity parameters will be received as #: arrays, with one row per polydispersity condition. ER = None               # type: Optional[Callable[[np.ndarray], np.ndarray]] #: Returns the occupied volume and the total volume for each parameter set. #: See :attr:ER for details on the parameters. VR = None               # type: Optional[Callable[[np.ndarray], Tuple[np.ndarray, np.ndarray]]] #: Arbitrary C code containing supporting functions, etc., to be inserted #: after everything in source.  This can include Iq and Iqxy functions with #: the full function signature, including all parameters. c_code = None #: Returns the form volume for python-based models.  Form volume is needed

• ## sasmodels/models/barbell.py

 r71b751d 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", ] def random():

• ## sasmodels/models/capped_cylinder.py

 r71b751d 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", ] def random():

• ## sasmodels/models/core_multi_shell.py

 r71b751d source = ["lib/sas_3j1x_x.c", "core_multi_shell.c"] have_Fq = True effective_radius_type = ["outer radius", "core radius"] def random(): 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,

• ## sasmodels/models/core_shell_bicelle.py

 r71b751d "core_shell_bicelle.c"] have_Fq = True effective_radius_type = [ "equivalent sphere", "outer rim radius", "half outer thickness", "half diagonal", ] def random():

• ## sasmodels/models/core_shell_bicelle_elliptical.py

 r71b751d "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", ] def random():

• ## sasmodels/models/core_shell_bicelle_elliptical_belt_rough.py

 r71b751d "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", ] demo = dict(scale=1, background=0,

• ## sasmodels/models/core_shell_cylinder.py

 r0159b5e source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "core_shell_cylinder.c"] have_Fq = True 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():

• ## sasmodels/models/core_shell_ellipsoid.py

 r71b751d source = ["lib/sas_3j1x_x.c", "lib/gauss76.c", "core_shell_ellipsoid.c"] have_Fq = True 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():
• ## sasmodels/models/core_shell_parallelepiped.c

 r71b751d 2.0 * length_a * length_b * thick_rim_c; #endif } static double radius_from_volume(double length_a, double length_b, double length_c, double thick_rim_a, double thick_rim_b, double thick_rim_c) { 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); } static double radius_from_crosssection(double length_a, double length_b, double thick_rim_a, double thick_rim_b) { const double area_xsec_paral = length_a*length_b + 2.0*thick_rim_a*length_b + 2.0*thick_rim_b*length_a; return sqrt(area_xsec_paral/M_PI); } static double 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) { 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); case 2: // half outer length a return 0.5 * length_a + thick_rim_a; case 3: // half outer length b return 0.5 * length_b + thick_rim_b; case 4: // half outer length c return 0.5 * length_c + thick_rim_c; case 5: // equivalent circular cross-section return radius_from_crosssection(length_a, length_b, thick_rim_a, thick_rim_b); 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)); 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)); } }
• ## sasmodels/models/core_shell_parallelepiped.py

 r71b751d is the scattering length of the solvent. .. note:: .. note:: the code actually implements two substitutions: $d(cos\alpha)$ is source = ["lib/gauss76.c", "core_shell_parallelepiped.c"] have_Fq = True 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():
• ## sasmodels/models/core_shell_sphere.c

 r71b751d { return M_4PI_3 * cube(radius + thickness); } static double effective_radius(int mode, double radius, double thickness) { switch (mode) { case 1: // outer radius return radius + thickness; case 2: // core radius return radius; } }

• ## sasmodels/models/cylinder.py

 r71b751d source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "cylinder.c"] have_Fq = True 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():

• ## sasmodels/models/fractal_core_shell.py

 reb3eb38 return radius + thickness tests = [[{'radius': 20.0, 'thickness': 10.0}, 'ER', 30.0], #tests = [[{'radius': 20.0, 'thickness': 10.0}, 'ER', 30.0], tests = [ #         # At some point the SasView 3.x test result was deemed incorrect. The #following tests were verified against NIST IGOR macros ver 7.850.

• ## sasmodels/models/hollow_rectangular_prism.c

 r71b751d // 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; 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); } static double effective_radius(int mode, double length_a, double b2a_ratio, double c2a_ratio, double thickness) //effective_radius_type = ["equivalent sphere","half length_a", "half length_b", "half length_c", //                         "equivalent outer circular cross-section","half ab diagonal","half diagonal"] // NOTE length_a is external dimension! { switch (mode) { case 1: // equivalent sphere return cbrt(0.75*cube(length_a)*b2a_ratio*c2a_ratio/M_PI); case 2: // half length_a return 0.5 * length_a; case 3: // half length_b return 0.5 * length_a*b2a_ratio; case 4: // half length_c return 0.5 * length_a*c2a_ratio; case 5: // equivalent outer circular cross-section return length_a*sqrt(b2a_ratio/M_PI); case 6: // half ab diagonal return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio))); case 7: // half diagonal return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio) + square(c2a_ratio))); } }
• ## sasmodels/models/hollow_rectangular_prism.py

 r0159b5e "Solvent scattering length density"], ["length_a", "Ang", 35, [0, inf], "volume", "Shorter side of the parallelepiped"], "Shortest, external, size of the parallelepiped"], ["b2a_ratio", "Ang", 1, [0, inf], "volume", "Ratio sides b/a"], source = ["lib/gauss76.c", "hollow_rectangular_prism.c"] have_Fq = True 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():
• ## sasmodels/models/hollow_rectangular_prism_thin_walls.c

 r71b751d // 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) { switch (mode) { case 1: // equivalent sphere return cbrt(0.75*cube(length_a)*b2a_ratio*c2a_ratio/M_PI); case 2: // half length_a return 0.5 * length_a; case 3: // half length_b return 0.5 * length_a*b2a_ratio; case 4: // half length_c return 0.5 * length_a*c2a_ratio; case 5: // equivalent outer circular cross-section return length_a*sqrt(b2a_ratio/M_PI); case 6: // half ab diagonal return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio))); case 7: // half diagonal return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio) + square(c2a_ratio))); } }
• ## sasmodels/models/hollow_rectangular_prism_thin_walls.py

 r0159b5e source = ["lib/gauss76.c", "hollow_rectangular_prism_thin_walls.c"] have_Fq = True 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", ]
• ## sasmodels/models/mono_gauss_coil.py

 r2d81cfe import numpy as np from numpy import inf, exp, errstate from numpy import inf name = "mono_gauss_coil" parameters = [ ["i_zero", "1/cm", 70.0, [0.0, inf], "", "Intensity at q=0"], ["rg", "Ang", 75.0, [0.0, inf], "", "Radius of gyration"], ["rg", "Ang", 75.0, [0.0, inf], "volume", "Radius of gyration"], ] # 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 source = ["mono_gauss_coil.c"] have_Fq = False effective_radius_type = ["R_g", "2R_g", "3R_g", "(5/3)^0.5*R_g"] 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():
• ## sasmodels/models/multilayer_vesicle.c

 r71b751d } static double 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; } static void Fq(double q,
• ## sasmodels/models/multilayer_vesicle.py

 r71b751d source = ["lib/sas_3j1x_x.c", "multilayer_vesicle.c"] have_Fq = True 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 effective_radius_type = ["outer radius"] def random():

• ## sasmodels/models/onion.py

 r71b751d single = False have_Fq = True effective_radius_type = ["outer radius"] profile_axes = ['Radius (A)', 'SLD (1e-6/A^2)'] 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 = {
• ## sasmodels/models/parallelepiped.c

 r71b751d } static double effective_radius(int mode, double length_a, double length_b, double length_c) { switch (mode) { case 1: // equivalent sphere return cbrt(0.75*length_a*length_b*length_c/M_PI); case 2: // half length_a return 0.5 * length_a; case 3: // half length_b return 0.5 * length_b; case 4: // half length_c return 0.5 * length_c; case 5: // equivalent circular cross-section return sqrt(length_a*length_b/M_PI); case 6: // half ab diagonal return 0.5*sqrt(length_a*length_a + length_b*length_b); case 7: // half diagonal return 0.5*sqrt(length_a*length_a + length_b*length_b + length_c*length_c); } } static void
• ## sasmodels/models/parallelepiped.py

 r71b751d 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: source = ["lib/gauss76.c", "parallelepiped.c"] have_Fq = True 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():
• ## sasmodels/models/pearl_necklace.py

 r2d81cfe thick_string_pd=0.2, thick_string_pd_n=5, ) tests = [[{}, 0.001, 17380.245], [{}, 'ER', 115.39502]] # ER function is not being used here, not that it is likely very sensible to # include an S(Q) with this model, the default in sasview 5.0 would be to the # "unconstrained" radius_effective. #tests = [[{}, 0.001, 17380.245], [{}, 'ER', 115.39502]] tests = [[{}, 0.001, 17380.245]]
• ## sasmodels/models/pringle.c

 r74768cb } static double effective_radius(int mode, double radius, double thickness, double alpha, double beta) { switch (mode) { case 1: // equivalent sphere return cbrt(0.75*radius*radius*thickness); case 2: // radius return radius; } } double Iq( double q,
• ## sasmodels/models/pringle.py

 ref07e95 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"] 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():

• ## sasmodels/models/raspberry.py

 ref07e95 ["radius_lg", "Ang", 5000, [0, inf], "volume", "radius of large spheres"], ["radius_sm", "Ang", 100, [0, inf], "", ["radius_sm", "Ang", 100, [0, inf], "volume", "radius of small spheres"], ["penetration", "Ang", 0, [-1, 1], "", ["penetration", "Ang", 0, [-1, 1], "volume", "fractional penetration depth of small spheres into large sphere"], ] source = ["lib/sas_3j1x_x.c", "raspberry.c"] effective_radius_type = ["radius_large", "radius_outer"] def random():
• ## sasmodels/models/rectangular_prism.c

 r71b751d { return length_a * (length_a*b2a_ratio) * (length_a*c2a_ratio); } static double effective_radius(int mode, double length_a, double b2a_ratio, double c2a_ratio) { switch (mode) { case 1: // equivalent sphere return cbrt(0.75*cube(length_a)*b2a_ratio*c2a_ratio/M_PI); case 2: // half length_a return 0.5 * length_a; case 3: // half length_b return 0.5 * length_a*b2a_ratio; case 4: // half length_c return 0.5 * length_a*c2a_ratio; case 5: // equivalent circular cross-section return length_a*sqrt(b2a_ratio/M_PI); case 6: // half ab diagonal return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio))); case 7: // half diagonal return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio) + square(c2a_ratio))); } }
• ## sasmodels/models/rectangular_prism.py

 r71b751d source = ["lib/gauss76.c", "rectangular_prism.c"] have_Fq = True 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():

• ## sasmodels/models/spherical_sld.c

 r71b751d static double form_volume( double fp_n_shells, double thickness[], double interface[]) static double outer_radius(double fp_n_shells, double thickness[], double interface[]) { int n_shells= (int)(fp_n_shells + 0.5); r += thickness[i] + interface[i]; } return M_4PI_3*cube(r); 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); }
• ## sasmodels/models/spherical_sld.py

 r0159b5e 0: erf($\nu z$) 1: Rpow($z^\nu$) 2: Lpow($z^\nu$) 3: Rexp($-\nu z$) 4: Lexp($-\nu z$) ["interface[n_shells]",  "Ang",        50.0,   [0, inf],       "volume", "thickness of the interface"], ["shape[n_shells]",      "",           0,      [SHAPES],       "", "interface shape"], ["nu[n_shells]",         "",           2.5,    [0, inf],       "", "interface shape exponent"], ["nu[n_shells]",         "",           2.5,    [1, inf],       "", "interface shape exponent"], ["n_steps",              "",           35,     [0, inf],       "", "number of steps in each interface (must be an odd integer)"], ] single = False  # TODO: fix low q behaviour have_Fq = True effective_radius_type = ["outer radius"] profile_axes = ['Radius (A)', 'SLD (1e-6/A^2)'] 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,