Changeset a08c47f in sasmodels for sasmodels


Ignore:
Timestamp:
Mar 11, 2019 12:43:15 PM (5 years ago)
Author:
Adam Washington <adam.washington@…>
Branches:
master, core_shell_microgels, magnetic_model, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
Children:
c11d09f
Parents:
b9c7379 (diff), e589e9a (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 remote-tracking branch 'upstream/beta_approx' into test_args

Location:
sasmodels
Files:
1 added
53 edited

Legend:

Unmodified
Added
Removed

  • sasmodels/compare.py

    r07646b6 rc1799d3  
    11521152        'rel_err'   : True, 
    11531153        'explore'   : False, 
    1154         'use_demo'  : True, 
     1154        'use_demo'  : False, 
    11551155        'zero'      : False, 
    11561156        'html'      : False, 

  • sasmodels/direct_model.py

    r304c775 rc1799d3  
    3232from .details import make_kernel_args, dispersion_mesh 
    3333from .modelinfo import DEFAULT_BACKGROUND 
     34from .product import RADIUS_MODE_ID 
    3435 
    3536# pylint: disable=unused-import 
    … …  
    6768    # type: (Kernel, ParameterSet, float, bool) -> np.ndarray 
    6869    """ 
    69     Like :func:`call_kernel`, but returning F, F^2, R_eff, V, V_form/V_shell. 
    70     """ 
    71     R_eff_type = int(pars.pop('radius_effective_type', 1.0)) 
     70    Like :func:`call_kernel`, but returning F, F^2, R_eff, V_shell, V_form/V_shell. 
     71 
     72    For solid objects V_shell is equal to V_form and the volume ratio is 1. 
     73 
     74    Use parameter *radius_effective_mode* to select the effective radius 
     75    calculation. 
     76    """ 
     77    R_eff_type = int(pars.pop(RADIUS_MODE_ID, 1.0)) 
    7278    mesh = get_mesh(calculator.info, pars, dim=calculator.dim, mono=mono) 
    7379    #print("pars", list(zip(*mesh))[0]) 
    … …  
    7682    return calculator.Fq(call_details, values, cutoff, is_magnetic, R_eff_type) 
    7783 
    78 def call_profile(model_info, **pars): 
    79     # type: (ModelInfo, ...) -> Tuple[np.ndarray, np.ndarray, Tuple[str, str]] 
     84def call_profile(model_info, pars=None): 
     85    # type: (ModelInfo, ParameterSet) -> Tuple[np.ndarray, np.ndarray, Tuple[str, str]] 
    8086    """ 
    8187    Returns the profile *x, y, (xlabel, ylabel)* representing the model. 
    8288    """ 
     89    if pars is None: 
     90        pars = {} 
    8391    args = {} 
    8492    for p in model_info.parameters.kernel_parameters: 
    … …  
    324332 
    325333        # Need to pull background out of resolution for multiple scattering 
    326         background = pars.get('background', DEFAULT_BACKGROUND) 
     334        default_background = self._model.info.parameters.common_parameters[1].default 
     335        background = pars.get('background', default_background) 
    327336        pars = pars.copy() 
    328337        pars['background'] = 0. 
    … …  
    377386        Generate a plottable profile. 
    378387        """ 
    379         return call_profile(self.model.info, **pars) 
     388        return call_profile(self.model.info, pars) 
    380389 
    381390def main(): 

  • sasmodels/generate.py

    r39a06c9 ra8a1d48  
    703703    """ 
    704704    for code in source: 
    705         m = _FQ_PATTERN.search(code) 
    706         if m is not None: 
     705        if _FQ_PATTERN.search(code) is not None: 
    707706            return True 
    708707    return False 
    … …  
    712711    # type: (List[str]) -> bool 
    713712    """ 
    714     Return True if C source defines "void Fq(". 
     713    Return True if C source defines "double shell_volume(". 
    715714    """ 
    716715    for code in source: 
    717         m = _SHELL_VOLUME_PATTERN.search(code) 
    718         if m is not None: 
     716        if _SHELL_VOLUME_PATTERN.search(code) is not None: 
    719717            return True 
    720718    return False 
    … …  
    10081006        pars = model_info.parameters.kernel_parameters 
    10091007    else: 
    1010         pars = model_info.parameters.COMMON + model_info.parameters.kernel_parameters 
     1008        pars = (model_info.parameters.common_parameters 
     1009                + model_info.parameters.kernel_parameters) 
    10111010    partable = make_partable(pars) 
    10121011    subst = dict(id=model_info.id.replace('_', '-'), 

  • sasmodels/jitter.py

    r1198f90 r7d97437  
    1515    pass 
    1616 
     17import matplotlib as mpl 
    1718import matplotlib.pyplot as plt 
    1819from matplotlib.widgets import Slider 
    … …  
    746747        pass 
    747748 
    748     axcolor = 'lightgoldenrodyellow' 
     749    # CRUFT: use axisbg instead of facecolor for matplotlib<2 
     750    facecolor_prop = 'facecolor' if mpl.__version__ > '2' else 'axisbg' 
     751    props = {facecolor_prop: 'lightgoldenrodyellow'} 
    749752 
    750753    ## add control widgets to plot 
    751     axes_theta = plt.axes([0.1, 0.15, 0.45, 0.04], axisbg=axcolor) 
    752     axes_phi = plt.axes([0.1, 0.1, 0.45, 0.04], axisbg=axcolor) 
    753     axes_psi = plt.axes([0.1, 0.05, 0.45, 0.04], axisbg=axcolor) 
     754    axes_theta = plt.axes([0.1, 0.15, 0.45, 0.04], **props) 
     755    axes_phi = plt.axes([0.1, 0.1, 0.45, 0.04], **props) 
     756    axes_psi = plt.axes([0.1, 0.05, 0.45, 0.04], **props) 
    754757    stheta = Slider(axes_theta, 'Theta', -90, 90, valinit=theta) 
    755758    sphi = Slider(axes_phi, 'Phi', -180, 180, valinit=phi) 
    756759    spsi = Slider(axes_psi, 'Psi', -180, 180, valinit=psi) 
    757760 
    758     axes_dtheta = plt.axes([0.75, 0.15, 0.15, 0.04], axisbg=axcolor) 
    759     axes_dphi = plt.axes([0.75, 0.1, 0.15, 0.04], axisbg=axcolor) 
    760     axes_dpsi = plt.axes([0.75, 0.05, 0.15, 0.04], axisbg=axcolor) 
     761    axes_dtheta = plt.axes([0.75, 0.15, 0.15, 0.04], **props) 
     762    axes_dphi = plt.axes([0.75, 0.1, 0.15, 0.04], **props) 
     763    axes_dpsi = plt.axes([0.75, 0.05, 0.15, 0.04], **props) 
    761764    # Note: using ridiculous definition of rectangle distribution, whose width 
    762765    # in sasmodels is sqrt(3) times the given width.  Divide by sqrt(3) to keep 

  • sasmodels/kernel.py

    re44432d rcd28947  
    133133        nout = 2 if self.info.have_Fq and self.dim == '1d' else 1 
    134134        total_weight = self.result[nout*self.q_input.nq + 0] 
     135        # Note: total_weight = sum(weight > cutoff), with cutoff >= 0, so it 
     136        # is okay to test directly against zero.  If weight is zero then I(q), 
     137        # etc. must also be zero. 
    135138        if total_weight == 0.: 
    136139            total_weight = 1. 

  • sasmodels/kernelcl.py

    rf872fd1 r8064d5e  
    5858import time 
    5959 
     60try: 
     61    from time import perf_counter as clock 
     62except ImportError: # CRUFT: python < 3.3 
     63    import sys 
     64    if sys.platform.count("darwin") > 0: 
     65        from time import time as clock 
     66    else: 
     67        from time import clock 
     68 
    6069import numpy as np  # type: ignore 
    61  
    6270 
    6371# Attempt to setup opencl. This may fail if the pyopencl package is not 
    … …  
    611619        #Call kernel and retrieve results 
    612620        wait_for = None 
    613         last_nap = time.clock() 
     621        last_nap = clock() 
    614622        step = 1000000//self.q_input.nq + 1 
    615623        for start in range(0, call_details.num_eval, step): 
    … …  
    622630                # Allow other processes to run 
    623631                wait_for[0].wait() 
    624                 current_time = time.clock() 
     632                current_time = clock() 
    625633                if current_time - last_nap > 0.5: 
    626634                    time.sleep(0.001) 

  • sasmodels/kernelpy.py

    re44432d raa8c6e0  
    4242        self.info = model_info 
    4343        self.dtype = np.dtype('d') 
     44        logger.info("make python model " + self.info.name) 
    4445 
    4546    def make_kernel(self, q_vectors): 

  • sasmodels/model_test.py

    rb9c7379 ra08c47f  
    4848import sys 
    4949import unittest 
     50import traceback 
    5051 
    5152try: 
    … …  
    7778# pylint: enable=unused-import 
    7879 
    79  
    8080def make_suite(loaders, models): 
    8181    # type: (List[str], List[str]) -> unittest.TestSuite 
    … …  
    8888    *models* is the list of models to test, or *["all"]* to test all models. 
    8989    """ 
    90     ModelTestCase = _hide_model_case_from_nose() 
    9190    suite = unittest.TestSuite() 
    9291 
    … …  
    9796        skip = [] 
    9897    for model_name in models: 
    99         if model_name in skip: 
    100             continue 
    101         model_info = load_model_info(model_name) 
    102  
    103         #print('------') 
    104         #print('found tests in', model_name) 
    105         #print('------') 
    106  
    107         # if ispy then use the dll loader to call pykernel 
    108         # don't try to call cl kernel since it will not be 
    109         # available in some environmentes. 
    110         is_py = callable(model_info.Iq) 
    111  
    112         # Some OpenCL drivers seem to be flaky, and are not producing the 
    113         # expected result.  Since we don't have known test values yet for 
    114         # all of our models, we are instead going to compare the results 
    115         # for the 'smoke test' (that is, evaluation at q=0.1 for the default 
    116         # parameters just to see that the model runs to completion) between 
    117         # the OpenCL and the DLL.  To do this, we define a 'stash' which is 
    118         # shared between OpenCL and DLL tests.  This is just a list.  If the 
    119         # list is empty (which it will be when DLL runs, if the DLL runs 
    120         # first), then the results are appended to the list.  If the list 
    121         # is not empty (which it will be when OpenCL runs second), the results 
    122         # are compared to the results stored in the first element of the list. 
    123         # This is a horrible stateful hack which only makes sense because the 
    124         # test suite is thrown away after being run once. 
    125         stash = [] 
    126  
    127         if is_py:  # kernel implemented in python 
    128             test_name = "%s-python"%model_name 
    129             test_method_name = "test_%s_python" % model_info.id 
     98        if model_name not in skip: 
     99            model_info = load_model_info(model_name) 
     100            _add_model_to_suite(loaders, suite, model_info) 
     101 
     102    return suite 
     103 
     104def _add_model_to_suite(loaders, suite, model_info): 
     105    ModelTestCase = _hide_model_case_from_nose() 
     106 
     107    #print('------') 
     108    #print('found tests in', model_name) 
     109    #print('------') 
     110 
     111    # if ispy then use the dll loader to call pykernel 
     112    # don't try to call cl kernel since it will not be 
     113    # available in some environmentes. 
     114    is_py = callable(model_info.Iq) 
     115 
     116    # Some OpenCL drivers seem to be flaky, and are not producing the 
     117    # expected result.  Since we don't have known test values yet for 
     118    # all of our models, we are instead going to compare the results 
     119    # for the 'smoke test' (that is, evaluation at q=0.1 for the default 
     120    # parameters just to see that the model runs to completion) between 
     121    # the OpenCL and the DLL.  To do this, we define a 'stash' which is 
     122    # shared between OpenCL and DLL tests.  This is just a list.  If the 
     123    # list is empty (which it will be when DLL runs, if the DLL runs 
     124    # first), then the results are appended to the list.  If the list 
     125    # is not empty (which it will be when OpenCL runs second), the results 
     126    # are compared to the results stored in the first element of the list. 
     127    # This is a horrible stateful hack which only makes sense because the 
     128    # test suite is thrown away after being run once. 
     129    stash = [] 
     130 
     131    if is_py:  # kernel implemented in python 
     132        test_name = "%s-python"%model_info.name 
     133        test_method_name = "test_%s_python" % model_info.id 
     134        test = ModelTestCase(test_name, model_info, 
     135                                test_method_name, 
     136                                platform="dll",  # so that 
     137                                dtype="double", 
     138                                stash=stash) 
     139        suite.addTest(test) 
     140    else:   # kernel implemented in C 
     141 
     142        # test using dll if desired 
     143        if 'dll' in loaders or not use_opencl(): 
     144            test_name = "%s-dll"%model_info.name 
     145            test_method_name = "test_%s_dll" % model_info.id 
    130146            test = ModelTestCase(test_name, model_info, 
    131                                  test_method_name, 
    132                                  platform="dll",  # so that 
    133                                  dtype="double", 
    134                                  stash=stash) 
     147                                    test_method_name, 
     148                                    platform="dll", 
     149                                    dtype="double", 
     150                                    stash=stash) 
    135151            suite.addTest(test) 
    136         else:   # kernel implemented in C 
    137  
    138             # test using dll if desired 
    139             if 'dll' in loaders: 
    140                 test_name = "%s-dll"%model_name 
    141                 test_method_name = "test_%s_dll" % model_info.id 
    142                 test = ModelTestCase(test_name, model_info, 
    143                                      test_method_name, 
    144                                      platform="dll", 
    145                                      dtype="double", 
    146                                      stash=stash) 
    147                 suite.addTest(test) 
    148  
    149             # test using opencl if desired and available 
    150             if 'opencl' in loaders and use_opencl(): 
    151                 test_name = "%s-opencl"%model_name 
    152                 test_method_name = "test_%s_opencl" % model_info.id 
    153                 # Using dtype=None so that the models that are only 
    154                 # correct for double precision are not tested using 
    155                 # single precision.  The choice is determined by the 
    156                 # presence of *single=False* in the model file. 
    157                 test = ModelTestCase(test_name, model_info, 
    158                                      test_method_name, 
    159                                      platform="ocl", dtype=None, 
    160                                      stash=stash) 
    161                 #print("defining", test_name) 
    162                 suite.addTest(test) 
    163  
    164             # test using cuda if desired and available 
    165             if 'cuda' in loaders and use_cuda(): 
    166                 test_name = "%s-cuda"%model_name 
    167                 test_method_name = "test_%s_cuda" % model_info.id 
    168                 # Using dtype=None so that the models that are only 
    169                 # correct for double precision are not tested using 
    170                 # single precision.  The choice is determined by the 
    171                 # presence of *single=False* in the model file. 
    172                 test = ModelTestCase(test_name, model_info, 
    173                                      test_method_name, 
    174                                      platform="cuda", dtype=None, 
    175                                      stash=stash) 
    176                 #print("defining", test_name) 
    177                 suite.addTest(test) 
    178  
    179     return suite 
     152 
     153        # test using opencl if desired and available 
     154        if 'opencl' in loaders and use_opencl(): 
     155            test_name = "%s-opencl"%model_info.name 
     156            test_method_name = "test_%s_opencl" % model_info.id 
     157            # Using dtype=None so that the models that are only 
     158            # correct for double precision are not tested using 
     159            # single precision.  The choice is determined by the 
     160            # presence of *single=False* in the model file. 
     161            test = ModelTestCase(test_name, model_info, 
     162                                    test_method_name, 
     163                                    platform="ocl", dtype=None, 
     164                                    stash=stash) 
     165            #print("defining", test_name) 
     166            suite.addTest(test) 
     167 
     168        # test using cuda if desired and available 
     169        if 'cuda' in loaders and use_cuda(): 
     170            test_name = "%s-cuda"%model_name 
     171            test_method_name = "test_%s_cuda" % model_info.id 
     172            # Using dtype=None so that the models that are only 
     173            # correct for double precision are not tested using 
     174            # single precision.  The choice is determined by the 
     175            # presence of *single=False* in the model file. 
     176            test = ModelTestCase(test_name, model_info, 
     177                                    test_method_name, 
     178                                    platform="cuda", dtype=None, 
     179                                    stash=stash) 
     180            #print("defining", test_name) 
     181            suite.addTest(test) 
     182 
    180183 
    181184def _hide_model_case_from_nose(): 
    … …  
    384387    for par in sorted(pars.keys()): 
    385388        # special handling of R_eff mode, which is not a usual parameter 
    386         if par == 'radius_effective_type': 
     389        if par == product.RADIUS_MODE_ID: 
    387390            continue 
    388391        parts = par.split('_pd') 
    … …  
    404407    return abs(target-actual)/shift < 1.5*10**-digits 
    405408 
    406 def run_one(model): 
    407     # type: (str) -> str 
    408     """ 
    409     Run the tests for a single model, printing the results to stdout. 
    410  
    411     *model* can by a python file, which is handy for checking user defined 
    412     plugin models. 
     409# CRUFT: old interface; should be deprecated and removed 
     410def run_one(model_name): 
     411    # msg = "use check_model(model_info) rather than run_one(model_name)" 
     412    # warnings.warn(msg, category=DeprecationWarning, stacklevel=2) 
     413    try: 
     414        model_info = load_model_info(model_name) 
     415    except Exception: 
     416        output = traceback.format_exc() 
     417        return output 
     418 
     419    success, output = check_model(model_info) 
     420    return output 
     421 
     422def check_model(model_info): 
     423    # type: (ModelInfo) -> str 
     424    """ 
     425    Run the tests for a single model, capturing the output. 
     426 
     427    Returns success status and the output string. 
    413428    """ 
    414429    # Note that running main() directly did not work from within the 
    … …  
    424439 
    425440    # Build a test suite containing just the model 
    426     loader = 'opencl' if use_opencl() else 'cuda' if use_cuda() else 'dll' 
    427     models = [model] 
    428     try: 
    429         suite = make_suite([loader], models) 
    430     except Exception: 
    431         import traceback 
    432         stream.writeln(traceback.format_exc()) 
    433         return 
     441    loaders = ['opencl' if use_opencl() else 'cuda' if use_cuda() else 'dll'] 
     442    suite = unittest.TestSuite() 
     443    _add_model_to_suite(loaders, suite, model_info) 
    434444 
    435445    # Warn if there are no user defined tests. 
    … …  
    446456    for test in suite: 
    447457        if not test.info.tests: 
    448             stream.writeln("Note: %s has no user defined tests."%model) 
     458            stream.writeln("Note: %s has no user defined tests."%model_info.name) 
    449459        break 
    450460    else: 
    … …  
    462472    output = stream.getvalue() 
    463473    stream.close() 
    464     return output 
     474    return result.wasSuccessful(), output 
    465475 
    466476 

  • sasmodels/modelinfo.py

    r39a06c9 rc1799d3  
    404404      parameters counted as n individual parameters p1, p2, ... 
    405405 
     406    * *common_parameters* is the list of common parameters, with a unique 
     407      copy for each model so that structure factors can have a default 
     408      background of 0.0. 
     409 
    406410    * *call_parameters* is the complete list of parameters to the kernel, 
    407411      including scale and background, with vector parameters recorded as 
    … …  
    416420    parameters don't use vector notation, and instead use p1, p2, ... 
    417421    """ 
    418     # scale and background are implicit parameters 
    419     COMMON = [Parameter(*p) for p in COMMON_PARAMETERS] 
    420  
    421422    def __init__(self, parameters): 
    422423        # type: (List[Parameter]) -> None 
     424 
     425        # scale and background are implicit parameters 
     426        # Need them to be unique to each model in case they have different 
     427        # properties, such as default=0.0 for structure factor backgrounds. 
     428        self.common_parameters = [Parameter(*p) for p in COMMON_PARAMETERS] 
     429 
    423430        self.kernel_parameters = parameters 
    424431        self._set_vector_lengths() 
    … …  
    468475                         if p.polydisperse and p.type not in ('orientation', 'magnetic')) 
    469476        self.pd_2d = set(p.name for p in self.call_parameters if p.polydisperse) 
     477 
     478    def set_zero_background(self): 
     479        """ 
     480        Set the default background to zero for this model.  This is done for 
     481        structure factor models. 
     482        """ 
     483        # type: () -> None 
     484        # Make sure background is the second common parameter. 
     485        assert self.common_parameters[1].id == "background" 
     486        self.common_parameters[1].default = 0.0 
     487        self.defaults = self._get_defaults() 
    470488 
    471489    def check_angles(self): 
    … …  
    567585    def _get_call_parameters(self): 
    568586        # type: () -> List[Parameter] 
    569         full_list = self.COMMON[:] 
     587        full_list = self.common_parameters[:] 
    570588        for p in self.kernel_parameters: 
    571589            if p.length == 1: 
    … …  
    670688 
    671689        # Gather the user parameters in order 
    672         result = control + self.COMMON 
     690        result = control + self.common_parameters 
    673691        for p in self.kernel_parameters: 
    674692            if not is2d and p.type in ('orientation', 'magnetic'): 
    … …  
    770788 
    771789    info = ModelInfo() 
     790 
     791    # Build the parameter table 
    772792    #print("make parameter table", kernel_module.parameters) 
    773793    parameters = make_parameter_table(getattr(kernel_module, 'parameters', [])) 
     794 
     795    # background defaults to zero for structure factor models 
     796    structure_factor = getattr(kernel_module, 'structure_factor', False) 
     797    if structure_factor: 
     798        parameters.set_zero_background() 
     799 
     800    # TODO: remove demo parameters 
     801    # The plots in the docs are generated from the model default values. 
     802    # Sascomp set parameters from the command line, and so doesn't need 
     803    # demo values for testing. 
    774804    demo = expand_pars(parameters, getattr(kernel_module, 'demo', None)) 
     805 
    775806    filename = abspath(kernel_module.__file__).replace('.pyc', '.py') 
    776807    kernel_id = splitext(basename(filename))[0] 

  • sasmodels/models/barbell.c

    rd42dd4a r99658f6  
    6363 
    6464static double 
     65radius_from_excluded_volume(double radius_bell, double radius, double length) 
     66{ 
     67    const double hdist = sqrt(square(radius_bell) - square(radius)); 
     68    const double length_tot = length + 2.0*(hdist+ radius); 
     69    return 0.5*cbrt(0.75*radius_bell*(2.0*radius_bell*length_tot + (radius_bell + length_tot)*(M_PI*radius_bell + length_tot))); 
     70} 
     71 
     72static double 
    6573radius_from_volume(double radius_bell, double radius, double length) 
    6674{ 
    … …  
    8189    switch (mode) { 
    8290    default: 
    83     case 1: // equivalent sphere 
     91    case 1: // equivalent cylinder excluded volume 
     92        return radius_from_excluded_volume(radius_bell, radius , length); 
     93    case 2: // equivalent volume sphere 
    8494        return radius_from_volume(radius_bell, radius , length); 
    85     case 2: // radius 
     95    case 3: // radius 
    8696        return radius; 
    87     case 3: // half length 
     97    case 4: // half length 
    8898        return 0.5*length; 
    89     case 4: // half total length 
     99    case 5: // half total length 
    90100        return radius_from_totallength(radius_bell,radius,length); 
    91101    } 

  • sasmodels/models/barbell.py

    ree60aa7 r99658f6  
    7979.. [#] H Kaya and N R deSouza, *J. Appl. Cryst.*, 37 (2004) 508-509 (addenda 
    8080   and errata) 
     81L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    8182 
    8283Authorship and Verification 
    … …  
    116117 
    117118source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "barbell.c"] 
    118 have_Fq = True 
     119have_Fq = True  
    119120effective_radius_type = [ 
    120     "equivalent sphere", "radius", "half length", "half total length", 
     121    "equivalent cylinder excluded volume","equivalent volume sphere", "radius", "half length", "half total length", 
    121122    ] 
    122123 

  • sasmodels/models/capped_cylinder.c

    rd42dd4a r99658f6  
    8585 
    8686static double 
     87radius_from_excluded_volume(double radius, double radius_cap, double length) 
     88{ 
     89    const double hc = radius_cap - sqrt(radius_cap*radius_cap - radius*radius); 
     90    const double length_tot = length + 2.0*hc; 
     91    return 0.5*cbrt(0.75*radius*(2.0*radius*length_tot + (radius + length_tot)*(M_PI*radius + length_tot))); 
     92} 
     93 
     94static double 
    8795radius_from_volume(double radius, double radius_cap, double length) 
    8896{ 
    … …  
    103111    switch (mode) { 
    104112    default: 
    105     case 1: // equivalent sphere 
     113    case 1: // equivalent cylinder excluded volume 
     114        return radius_from_excluded_volume(radius, radius_cap, length); 
     115    case 2: // equivalent volume sphere 
    106116        return radius_from_volume(radius, radius_cap, length); 
    107     case 2: // radius 
     117    case 3: // radius 
    108118        return radius; 
    109     case 3: // half length 
     119    case 4: // half length 
    110120        return 0.5*length; 
    111     case 4: // half total length 
     121    case 5: // half total length 
    112122        return radius_from_totallength(radius, radius_cap,length); 
    113123    } 

  • sasmodels/models/capped_cylinder.py

    ree60aa7 r99658f6  
    8282.. [#] H Kaya and N-R deSouza, *J. Appl. Cryst.*, 37 (2004) 508-509 (addenda 
    8383   and errata) 
     84L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    8485 
    8586Authorship and Verification 
    … …  
    138139have_Fq = True 
    139140effective_radius_type = [ 
    140     "equivalent sphere", "radius", "half length", "half total length", 
     141    "equivalent cylinder excluded volume", "equivalent volume sphere", "radius", "half length", "half total length", 
    141142    ] 
    142143 

  • sasmodels/models/core_shell_bicelle.c

    rd42dd4a r99658f6  
    3737 
    3838static double 
     39radius_from_excluded_volume(double radius, double thick_rim, double thick_face, double length) 
     40{ 
     41    const double radius_tot = radius + thick_rim; 
     42    const double length_tot = length + 2.0*thick_face; 
     43    return 0.5*cbrt(0.75*radius_tot*(2.0*radius_tot*length_tot + (radius_tot + length_tot)*(M_PI*radius_tot + length_tot))); 
     44} 
     45 
     46static double 
    3947radius_from_volume(double radius, double thick_rim, double thick_face, double length) 
    4048{ 
    … …  
    5664    switch (mode) { 
    5765    default: 
    58     case 1: // equivalent sphere 
     66    case 1: // equivalent cylinder excluded volume 
     67        return radius_from_excluded_volume(radius, thick_rim, thick_face, length); 
     68    case 2: // equivalent sphere 
    5969        return radius_from_volume(radius, thick_rim, thick_face, length); 
    60     case 2: // outer rim radius 
     70    case 3: // outer rim radius 
    6171        return radius + thick_rim; 
    62     case 3: // half outer thickness 
     72    case 4: // half outer thickness 
    6373        return 0.5*length + thick_face; 
    64     case 4: // half diagonal 
     74    case 5: // half diagonal 
    6575        return radius_from_diagonal(radius,thick_rim,thick_face,length); 
    6676    } 

  • sasmodels/models/core_shell_bicelle.py

    ree60aa7 r99658f6  
    8989   from Proquest <http://search.proquest.com/docview/304915826?accountid 
    9090   =26379>`_ 
     91    
     92   L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    9193 
    9294Authorship and Verification 
    … …  
    156158have_Fq = True 
    157159effective_radius_type = [ 
    158     "equivalent sphere", "outer rim radius", 
     160    "excluded volume","equivalent volume sphere", "outer rim radius", 
    159161    "half outer thickness", "half diagonal", 
    160162    ] 

  • sasmodels/models/core_shell_bicelle_elliptical.c

    rd42dd4a r99658f6  
    88{ 
    99    return M_PI*(r_minor+thick_rim)*(r_minor*x_core+thick_rim)*(length+2.0*thick_face); 
     10} 
     11 
     12static double 
     13radius_from_excluded_volume(double r_minor, double x_core, double thick_rim, double thick_face, double length) 
     14{ 
     15    const double r_equiv     = sqrt((r_minor + thick_rim)*(r_minor*x_core + thick_rim)); 
     16    const double length_tot  = length + 2.0*thick_face; 
     17    return 0.5*cbrt(0.75*r_equiv*(2.0*r_equiv*length_tot + (r_equiv + length_tot)*(M_PI*r_equiv + length_tot))); 
    1018} 
    1119 
    … …  
    3139    switch (mode) { 
    3240    default: 
    33     case 1: // equivalent sphere 
     41    case 1: // equivalent cylinder excluded volume 
     42        return radius_from_excluded_volume(r_minor, x_core, thick_rim, thick_face, length); 
     43    case 2: // equivalent volume sphere 
    3444        return radius_from_volume(r_minor, x_core, thick_rim, thick_face, length); 
    35     case 2: // outer rim average radius 
     45    case 3: // outer rim average radius 
    3646        return 0.5*r_minor*(1.0 + x_core) + thick_rim; 
    37     case 3: // outer rim min radius 
     47    case 4: // outer rim min radius 
    3848        return (x_core < 1.0 ? x_core*r_minor+thick_rim : r_minor+thick_rim); 
    39     case 4: // outer max radius 
     49    case 5: // outer max radius 
    4050        return (x_core > 1.0 ? x_core*r_minor+thick_rim : r_minor+thick_rim); 
    41     case 5: // half outer thickness 
     51    case 6: // half outer thickness 
    4252        return 0.5*length + thick_face; 
    43     case 6: // half diagonal 
     53    case 7: // half diagonal 
    4454        return radius_from_diagonal(r_minor,x_core,thick_rim,thick_face,length); 
    4555    } 

  • sasmodels/models/core_shell_bicelle_elliptical.py

    r304c775 r99658f6  
    100100 
    101101.. [#] 
     102L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    102103 
    103104Authorship and Verification 
    … …  
    148149have_Fq = True 
    149150effective_radius_type = [ 
    150     "equivalent sphere", "outer rim average radius", "outer rim min radius", 
     151    "equivalent cylinder excluded volume", "equivalent volume sphere", "outer rim average radius", "outer rim min radius", 
    151152    "outer max radius", "half outer thickness", "half diagonal", 
    152153    ] 

  • sasmodels/models/core_shell_bicelle_elliptical_belt_rough.c

    rd42dd4a r99658f6  
    99    return M_PI*(  (r_minor + thick_rim)*(r_minor*x_core + thick_rim)* length + 
    1010                 square(r_minor)*x_core*2.0*thick_face  ); 
     11} 
     12 
     13static double 
     14radius_from_excluded_volume(double r_minor, double x_core, double thick_rim, double thick_face, double length) 
     15{ 
     16    const double r_equiv     = sqrt((r_minor + thick_rim)*(r_minor*x_core + thick_rim)); 
     17    const double length_tot  = length + 2.0*thick_face; 
     18    return 0.5*cbrt(0.75*r_equiv*(2.0*r_equiv*length_tot + (r_equiv + length_tot)*(M_PI*r_equiv + length_tot))); 
    1119} 
    1220 
    … …  
    3240    switch (mode) { 
    3341    default: 
    34     case 1: // equivalent sphere 
     42    case 1: // equivalent cylinder excluded volume 
     43        return radius_from_excluded_volume(r_minor, x_core, thick_rim, thick_face, length); 
     44    case 2: // equivalent sphere 
    3545        return radius_from_volume(r_minor, x_core, thick_rim, thick_face, length); 
    36     case 2: // outer rim average radius 
     46    case 3: // outer rim average radius 
    3747        return 0.5*r_minor*(1.0 + x_core) + thick_rim; 
    38     case 3: // outer rim min radius 
     48    case 4: // outer rim min radius 
    3949        return (x_core < 1.0 ? x_core*r_minor+thick_rim : r_minor+thick_rim); 
    40     case 4: // outer max radius 
     50    case 5: // outer max radius 
    4151        return (x_core > 1.0 ? x_core*r_minor+thick_rim : r_minor+thick_rim); 
    42     case 5: // half outer thickness 
     52    case 6: // half outer thickness 
    4353        return 0.5*length + thick_face; 
    44     case 6: // half diagonal (this ignores the missing "corners", so may give unexpected answer if thick_face 
     54    case 7: // half diagonal (this ignores the missing "corners", so may give unexpected answer if thick_face 
    4555            //  or thick_rim is extremely large) 
    4656        return radius_from_diagonal(r_minor,x_core,thick_rim,thick_face,length); 

  • sasmodels/models/core_shell_bicelle_elliptical_belt_rough.py

    r304c775 r99658f6  
    112112 
    113113.. [#] 
     114L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    114115 
    115116Authorship and Verification 
    … …  
    161162have_Fq = True 
    162163effective_radius_type = [ 
    163     "equivalent sphere", "outer rim average radius", "outer rim min radius", 
     164    "equivalent cylinder excluded volume", "equivalent volume sphere", "outer rim average radius", "outer rim min radius", 
    164165    "outer max radius", "half outer thickness", "half diagonal", 
    165166    ] 

  • sasmodels/models/core_shell_cylinder.c

    rd42dd4a r99658f6  
    1111{ 
    1212    return M_PI*square(radius+thickness)*(length+2.0*thickness); 
     13} 
     14 
     15static double 
     16radius_from_excluded_volume(double radius, double thickness, double length) 
     17{ 
     18    const double radius_tot = radius + thickness; 
     19    const double length_tot = length + 2.0*thickness; 
     20    return 0.5*cbrt(0.75*radius_tot*(2.0*radius_tot*length_tot + (radius_tot + length_tot)*(M_PI*radius_tot + length_tot))); 
    1321} 
    1422 
    … …  
    3341    switch (mode) { 
    3442    default: 
    35     case 1: // equivalent sphere 
     43    case 1: //cylinder excluded volume 
     44        return radius_from_excluded_volume(radius, thickness, length); 
     45    case 2: // equivalent volume sphere 
    3646        return radius_from_volume(radius, thickness, length); 
    37     case 2: // outer radius 
     47    case 3: // outer radius 
    3848        return radius + thickness; 
    39     case 3: // half outer length 
     49    case 4: // half outer length 
    4050        return 0.5*length + thickness; 
    41     case 4: // half min outer length 
     51    case 5: // half min outer length 
    4252        return (radius < 0.5*length ? radius + thickness : 0.5*length + thickness); 
    43     case 5: // half max outer length 
     53    case 6: // half max outer length 
    4454        return (radius > 0.5*length ? radius + thickness : 0.5*length + thickness); 
    45     case 6: // half outer diagonal 
     55    case 7: // half outer diagonal 
    4656        return radius_from_diagonal(radius,thickness,length); 
    4757    } 

  • sasmodels/models/core_shell_cylinder.py

    ree60aa7 r99658f6  
    7070   1445-1452 
    7171.. [#kline] S R Kline, *J Appl. Cryst.*, 39 (2006) 895 
     72L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    7273 
    7374Authorship and Verification 
    … …  
    133134have_Fq = True 
    134135effective_radius_type = [ 
    135     "equivalent sphere", "outer radius", "half outer length", 
    136     "half min outer dimension", "half max outer dimension", 
    137     "half outer diagonal", 
     136    "excluded volume", "equivalent volume sphere", "outer radius", "half outer length", 
     137    "half min outer dimension", "half max outer dimension", "half outer diagonal", 
    138138    ] 
    139139 

  • sasmodels/models/core_shell_ellipsoid.c

    rd42dd4a r99658f6  
    7575    switch (mode) { 
    7676    default: 
    77     case 1: // equivalent sphere 
     77    case 1: // average outer curvature 
     78        return radius_from_curvature(radius_equat_core, x_core, thick_shell, x_polar_shell); 
     79    case 2: // equivalent volume sphere 
    7880        return radius_from_volume(radius_equat_core, x_core, thick_shell, x_polar_shell); 
    79     case 2: // average outer curvature 
    80         return radius_from_curvature(radius_equat_core, x_core, thick_shell, x_polar_shell); 
    8181    case 3: // min outer radius 
    8282        return (radius_polar_tot < radius_equat_tot ? radius_polar_tot : radius_equat_tot); 

  • sasmodels/models/core_shell_ellipsoid.py

    ree60aa7 r99658f6  
    147147have_Fq = True 
    148148effective_radius_type = [ 
    149     "equivalent sphere", "average outer curvature", 
     149    "average outer curvature", "equivalent volume sphere",      
    150150    "min outer radius", "max outer radius", 
    151151    ] 

  • sasmodels/models/core_shell_parallelepiped.c

    rd42dd4a r99658f6  
    2828 
    2929static double 
     30radius_from_excluded_volume(double length_a, double length_b, double length_c, 
     31                   double thick_rim_a, double thick_rim_b, double thick_rim_c) 
     32{ 
     33    double r_equiv, length; 
     34    double lengths[3] = {length_a+thick_rim_a, length_b+thick_rim_b, length_c+thick_rim_c}; 
     35    double lengthmax = fmax(lengths[0],fmax(lengths[1],lengths[2])); 
     36    double length_1 = lengthmax; 
     37    double length_2 = lengthmax; 
     38    double length_3 = lengthmax; 
     39 
     40    for(int ilen=0; ilen<3; ilen++) { 
     41        if (lengths[ilen] < length_1) { 
     42            length_2 = length_1; 
     43            length_1 = lengths[ilen]; 
     44            } else { 
     45                if (lengths[ilen] < length_2) { 
     46                        length_2 = lengths[ilen]; 
     47                } 
     48            } 
     49    } 
     50    if(length_2-length_1 > length_3-length_2) { 
     51        r_equiv = sqrt(length_2*length_3/M_PI); 
     52        length  = length_1; 
     53    } else  { 
     54        r_equiv = sqrt(length_1*length_2/M_PI); 
     55        length  = length_3; 
     56    } 
     57 
     58    return 0.5*cbrt(0.75*r_equiv*(2.0*r_equiv*length + (r_equiv + length)*(M_PI*r_equiv + length))); 
     59} 
     60 
     61static double 
    3062radius_from_volume(double length_a, double length_b, double length_c, 
    3163                   double thick_rim_a, double thick_rim_b, double thick_rim_c) 
    … …  
    4880    switch (mode) { 
    4981    default: 
    50     case 1: // equivalent sphere 
     82    case 1: // equivalent cylinder excluded volume 
     83        return radius_from_excluded_volume(length_a, length_b, length_c, thick_rim_a, thick_rim_b, thick_rim_c); 
     84    case 2: // equivalent volume sphere 
    5185        return radius_from_volume(length_a, length_b, length_c, thick_rim_a, thick_rim_b, thick_rim_c); 
    52     case 2: // half outer length a 
     86    case 3: // half outer length a 
    5387        return 0.5 * length_a + thick_rim_a; 
    54     case 3: // half outer length b 
     88    case 4: // half outer length b 
    5589        return 0.5 * length_b + thick_rim_b; 
    56     case 4: // half outer length c 
     90    case 5: // half outer length c 
    5791        return 0.5 * length_c + thick_rim_c; 
    58     case 5: // equivalent circular cross-section 
     92    case 6: // equivalent circular cross-section 
    5993        return radius_from_crosssection(length_a, length_b, thick_rim_a, thick_rim_b); 
    60     case 6: // half outer ab diagonal 
     94    case 7: // half outer ab diagonal 
    6195        return 0.5*sqrt(square(length_a+ 2.0*thick_rim_a) + square(length_b+ 2.0*thick_rim_b)); 
    62     case 7: // half outer diagonal 
     96    case 8: // half outer diagonal 
    6397        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)); 
    6498    } 

  • sasmodels/models/core_shell_parallelepiped.py

    ree60aa7 r99658f6  
    173173   from Proquest <http://search.proquest.com/docview/304915826?accountid 
    174174   =26379>`_ 
     175L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    175176 
    176177Authorship and Verification 
    … …  
    228229have_Fq = True 
    229230effective_radius_type = [ 
    230     "equivalent sphere", 
     231    "equivalent cylinder excluded volume",  
     232    "equivalent volume sphere", 
    231233    "half outer length_a", "half outer length_b", "half outer length_c", 
    232234    "equivalent circular cross-section", 

  • sasmodels/models/cylinder.c

    rd42dd4a r99658f6  
    1111{ 
    1212    return sas_2J1x_x(qab*radius) * sas_sinx_x(qc*0.5*length); 
     13} 
     14 
     15static double 
     16radius_from_excluded_volume(double radius, double length) 
     17{ 
     18    return 0.5*cbrt(0.75*radius*(2.0*radius*length + (radius + length)*(M_PI*radius + length))); 
    1319} 
    1420 
    … …  
    3137    default: 
    3238    case 1: 
     39        return radius_from_excluded_volume(radius, length); 
     40    case 2: 
    3341        return radius_from_volume(radius, length); 
    34     case 2: 
     42    case 3: 
    3543        return radius; 
    36     case 3: 
     44    case 4: 
    3745        return 0.5*length; 
    38     case 4: 
     46    case 5: 
    3947        return (radius < 0.5*length ? radius : 0.5*length); 
    40     case 5: 
     48    case 6: 
    4149        return (radius > 0.5*length ? radius : 0.5*length); 
    42     case 6: 
     50    case 7: 
    4351        return radius_from_diagonal(radius,length); 
    4452    } 

  • sasmodels/models/cylinder.py

    r304c775 r99658f6  
    9898J. S. Pedersen, Adv. Colloid Interface Sci. 70, 171-210 (1997). 
    9999G. Fournet, Bull. Soc. Fr. Mineral. Cristallogr. 74, 39-113 (1951). 
     100L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    100101""" 
    101102 
    … …  
    140141have_Fq = True 
    141142effective_radius_type = [ 
    142     "equivalent sphere", "radius", 
     143    "excluded volume", "equivalent volume sphere", "radius", 
    143144    "half length", "half min dimension", "half max dimension", "half diagonal", 
    144145    ] 
    … …  
    185186radius, length = parameters[2][2], parameters[3][2] 
    186187tests.extend([ 
    187     ({'radius_effective_type': 0}, 0.1, None, None, 0., pi*radius**2*length, 1.0), 
    188     ({'radius_effective_type': 1}, 0.1, None, None, (0.75*radius**2*length)**(1./3.), None, None), 
    189     ({'radius_effective_type': 2}, 0.1, None, None, radius, None, None), 
    190     ({'radius_effective_type': 3}, 0.1, None, None, length/2., None, None), 
    191     ({'radius_effective_type': 4}, 0.1, None, None, min(radius, length/2.), None, None), 
    192     ({'radius_effective_type': 5}, 0.1, None, None, max(radius, length/2.), None, None), 
    193     ({'radius_effective_type': 6}, 0.1, None, None, np.sqrt(4*radius**2 + length**2)/2., None, None), 
     188    ({'radius_effective_mode': 0}, 0.1, None, None, 0., pi*radius**2*length, 1.0),    
     189    ({'radius_effective_mode': 1}, 0.1, None, None, 0.5*(0.75*radius*(2.0*radius*length + (radius + length)*(pi*radius + length)))**(1./3.), None, None),     
     190    ({'radius_effective_mode': 2}, 0.1, None, None, (0.75*radius**2*length)**(1./3.), None, None), 
     191    ({'radius_effective_mode': 3}, 0.1, None, None, radius, None, None), 
     192    ({'radius_effective_mode': 4}, 0.1, None, None, length/2., None, None), 
     193    ({'radius_effective_mode': 5}, 0.1, None, None, min(radius, length/2.), None, None), 
     194    ({'radius_effective_mode': 6}, 0.1, None, None, max(radius, length/2.), None, None), 
     195    ({'radius_effective_mode': 7}, 0.1, None, None, np.sqrt(4*radius**2 + length**2)/2., None, None), 
    194196]) 
    195197del radius, length 

  • sasmodels/models/ellipsoid.c

    rd42dd4a r99658f6  
    3636    switch (mode) { 
    3737    default: 
    38     case 1: // equivalent sphere 
     38    case 1: // average curvature 
     39        return radius_from_curvature(radius_polar, radius_equatorial); 
     40    case 2: // equivalent volume sphere 
    3941        return radius_from_volume(radius_polar, radius_equatorial); 
    40     case 2: // average curvature 
    41         return radius_from_curvature(radius_polar, radius_equatorial); 
    4242    case 3: // min radius 
    4343        return (radius_polar < radius_equatorial ? radius_polar : radius_equatorial); 

  • sasmodels/models/ellipsoid.py

    ree60aa7 r99658f6  
    170170have_Fq = True 
    171171effective_radius_type = [ 
    172     "equivalent sphere", "average curvature", "min radius", "max radius", 
     172    "average curvature", "equivalent volume sphere", "min radius", "max radius", 
    173173    ] 
    174174 

  • sasmodels/models/elliptical_cylinder.c

    rd42dd4a r99658f6  
    33{ 
    44    return M_PI * radius_minor * radius_minor * r_ratio * length; 
     5} 
     6 
     7static double 
     8radius_from_excluded_volume(double radius_minor, double r_ratio, double length) 
     9{ 
     10    const double r_equiv = sqrt(radius_minor*radius_minor*r_ratio); 
     11    return 0.5*cbrt(0.75*r_equiv*(2.0*r_equiv*length + (r_equiv + length)*(M_PI*r_equiv + length))); 
    512} 
    613 
    … …  
    3845    switch (mode) { 
    3946    default: 
    40     case 1: // equivalent sphere 
     47    case 1: // equivalent cylinder excluded volume 
     48        return radius_from_excluded_volume(radius_minor, r_ratio, length); 
     49    case 2: // equivalent volume sphere 
    4150        return radius_from_volume(radius_minor, r_ratio, length); 
    42     case 2: // average radius 
     51    case 3: // average radius 
    4352        return 0.5*radius_minor*(1.0 + r_ratio); 
    44     case 3: // min radius 
     53    case 4: // min radius 
    4554        return (r_ratio > 1.0 ? radius_minor : r_ratio*radius_minor); 
    46     case 4: // max radius 
     55    case 5: // max radius 
    4756        return (r_ratio < 1.0 ? radius_minor : r_ratio*radius_minor); 
    48     case 5: // equivalent circular cross-section 
     57    case 6: // equivalent circular cross-section 
    4958        return sqrt(radius_minor*radius_minor*r_ratio); 
    50     case 6: // half length 
     59    case 7: // half length 
    5160        return 0.5*length; 
    52     case 7: // half min dimension 
     61    case 8: // half min dimension 
    5362        return radius_from_min_dimension(radius_minor,r_ratio,0.5*length); 
    54     case 8: // half max dimension 
     63    case 9: // half max dimension 
    5564        return radius_from_max_dimension(radius_minor,r_ratio,0.5*length); 
    56     case 9: // half diagonal 
     65    case 10: // half diagonal 
    5766        return radius_from_diagonal(radius_minor,r_ratio,length); 
    5867    } 

  • sasmodels/models/elliptical_cylinder.py

    ree60aa7 r99658f6  
    8888L A Feigin and D I Svergun, *Structure Analysis by Small-Angle X-Ray and 
    8989Neutron Scattering*, Plenum, New York, (1987) [see table 3.4] 
     90L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    9091 
    9192Authorship and Verification 
    … …  
    124125have_Fq = True 
    125126effective_radius_type = [ 
    126     "equivalent sphere", "average radius", "min radius", "max radius", 
     127    "equivalent cylinder excluded volume", "equivalent volume sphere", "average radius", "min radius", "max radius", 
    127128    "equivalent circular cross-section", 
    128129    "half length", "half min dimension", "half max dimension", "half diagonal", 

  • sasmodels/models/hardsphere.py

    r304c775 rc1799d3  
    162162    return pars 
    163163 
    164 demo = dict(radius_effective=200, volfraction=0.2, 
    165             radius_effective_pd=0.1, radius_effective_pd_n=40) 
    166164# Q=0.001 is in the Taylor series, low Q part, so add Q=0.1, 
    167165# assuming double precision sasview is correct 

  • sasmodels/models/hollow_cylinder.c

    rd42dd4a r99658f6  
    1111{ 
    1212    return M_PI*length*square(radius+thickness); 
     13} 
     14 
     15static double 
     16radius_from_excluded_volume(double radius, double thickness, double length) 
     17{ 
     18    const double radius_tot = radius + thickness; 
     19    return 0.5*cbrt(0.75*radius_tot*(2.0*radius_tot*length + (radius_tot + length)*(M_PI*radius_tot + length))); 
    1320} 
    1421 
    … …  
    3138    switch (mode) { 
    3239    default: 
    33     case 1: // equivalent sphere 
     40    case 1: // excluded volume 
     41        return radius_from_excluded_volume(radius, thickness, length); 
     42    case 2: // equivalent volume sphere 
    3443        return radius_from_volume(radius, thickness, length); 
    35     case 2: // outer radius 
     44    case 3: // outer radius 
    3645        return radius + thickness; 
    37     case 3: // half length 
     46    case 4: // half length 
    3847        return 0.5*length; 
    39     case 4: // half outer min dimension 
     48    case 5: // half outer min dimension 
    4049        return (radius + thickness < 0.5*length ? radius + thickness : 0.5*length); 
    41     case 5: // half outer max dimension 
     50    case 6: // half outer max dimension 
    4251        return (radius + thickness > 0.5*length ? radius + thickness : 0.5*length); 
    43     case 6: // half outer diagonal 
     52    case 7: // half outer diagonal 
    4453        return radius_from_diagonal(radius,thickness,length); 
    4554    } 

  • sasmodels/models/hollow_cylinder.py

    r304c775 r99658f6  
    6060.. [#] L A Feigin and D I Svergun, *Structure Analysis by Small-Angle X-Ray and 
    6161   Neutron Scattering*, Plenum Press, New York, (1987) 
     62L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    6263 
    6364Authorship and Verification 
    … …  
    102103have_Fq = True 
    103104effective_radius_type = [ 
    104     "equivalent sphere", "outer radius", "half length", 
     105    "excluded volume", "equivalent outer volume sphere", "outer radius", "half length", 
    105106    "half outer min dimension", "half outer max dimension", 
    106107    "half outer diagonal", 
    … …  
    140141    [{}, 0.00005, 1764.926], 
    141142    [{}, 0.1, None, None, 
    142      (3./4*(radius+thickness)**2*length)**(1./3),  # R_eff from volume 
     143     0.5*(0.75*(radius+thickness)*(2.0*(radius+thickness)*length + ((radius+thickness) + length)*(pi*(radius+thickness) + length)))**(1./3.),  # R_eff from excluded volume 
    143144     pi*((radius+thickness)**2-radius**2)*length,  # shell volume 
    144145     (radius+thickness)**2/((radius+thickness)**2 - radius**2), # form:shell ratio 

  • sasmodels/models/hollow_rectangular_prism.c

    rd42dd4a r99658f6  
    2222 
    2323static double 
     24radius_from_excluded_volume(double length_a, double b2a_ratio, double c2a_ratio) 
     25{ 
     26    const double r_equiv = sqrt(length_a*length_a*b2a_ratio/M_PI); 
     27    const double length_c = length_a*c2a_ratio; 
     28    return 0.5*cbrt(0.75*r_equiv*(2.0*r_equiv*length_c + (r_equiv + length_c)*(M_PI*r_equiv + length_c))); 
     29} 
     30 
     31static double 
    2432effective_radius(int mode, double length_a, double b2a_ratio, double c2a_ratio, double thickness) 
    2533// NOTE length_a is external dimension! 
    … …  
    2735    switch (mode) { 
    2836    default: 
    29     case 1: // equivalent sphere 
     37    case 1: // equivalent cylinder excluded volume 
     38        return radius_from_excluded_volume(length_a, b2a_ratio, c2a_ratio); 
     39    case 2: // equivalent outer volume sphere 
    3040        return cbrt(cube(length_a)*b2a_ratio*c2a_ratio/M_4PI_3); 
    31     case 2: // half length_a 
     41    case 3: // half length_a 
    3242        return 0.5 * length_a; 
    33     case 3: // half length_b 
     43    case 4: // half length_b 
    3444        return 0.5 * length_a*b2a_ratio; 
    35     case 4: // half length_c 
     45    case 5: // half length_c 
    3646        return 0.5 * length_a*c2a_ratio; 
    37     case 5: // equivalent outer circular cross-section 
     47    case 6: // equivalent outer circular cross-section 
    3848        return length_a*sqrt(b2a_ratio/M_PI); 
    39     case 6: // half ab diagonal 
     49    case 7: // half ab diagonal 
    4050        return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio))); 
    41     case 7: // half diagonal 
     51    case 8: // half diagonal 
    4252        return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio) + square(c2a_ratio))); 
    4353    } 

  • sasmodels/models/hollow_rectangular_prism.py

    ree60aa7 r99658f6  
    9898 
    9999.. [#Nayuk2012] R Nayuk and K Huber, *Z. Phys. Chem.*, 226 (2012) 837-854 
     100L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    100101 
    101102 
    … …  
    150151have_Fq = True 
    151152effective_radius_type = [ 
    152     "equivalent sphere", "half length_a", "half length_b", "half length_c", 
     153    "equivalent cylinder excluded volume", "equivalent outer volume sphere",  
     154    "half length_a", "half length_b", "half length_c", 
    153155    "equivalent outer circular cross-section", 
    154156    "half ab diagonal", "half diagonal", 

  • sasmodels/models/hollow_rectangular_prism_thin_walls.c

    rd42dd4a r99658f6  
    1818 
    1919static double 
     20radius_from_excluded_volume(double length_a, double b2a_ratio, double c2a_ratio) 
     21{ 
     22    const double r_equiv = sqrt(length_a*length_a*b2a_ratio/M_PI); 
     23    const double length_c = length_a*c2a_ratio; 
     24    return 0.5*cbrt(0.75*r_equiv*(2.0*r_equiv*length_c + (r_equiv + length_c)*(M_PI*r_equiv + length_c))); 
     25} 
     26 
     27static double 
    2028effective_radius(int mode, double length_a, double b2a_ratio, double c2a_ratio) 
    2129{ 
    2230    switch (mode) { 
    2331    default: 
    24     case 1: // equivalent sphere 
     32    case 1: // equivalent cylinder excluded volume 
     33        return radius_from_excluded_volume(length_a, b2a_ratio, c2a_ratio); 
     34    case 2: // equivalent outer volume sphere 
    2535        return cbrt(cube(length_a)*b2a_ratio*c2a_ratio/M_4PI_3); 
    26     case 2: // half length_a 
     36    case 3: // half length_a 
    2737        return 0.5 * length_a; 
    28     case 3: // half length_b 
     38    case 4: // half length_b 
    2939        return 0.5 * length_a*b2a_ratio; 
    30     case 4: // half length_c 
     40    case 5: // half length_c 
    3141        return 0.5 * length_a*c2a_ratio; 
    32     case 5: // equivalent outer circular cross-section 
     42    case 6: // equivalent outer circular cross-section 
    3343        return length_a*sqrt(b2a_ratio/M_PI); 
    34     case 6: // half ab diagonal 
     44    case 7: // half ab diagonal 
    3545        return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio))); 
    36     case 7: // half diagonal 
     46    case 8: // half diagonal 
    3747        return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio) + square(c2a_ratio))); 
    3848    } 

  • sasmodels/models/hollow_rectangular_prism_thin_walls.py

    ree60aa7 r99658f6  
    7272 
    7373.. [#Nayuk2012] R Nayuk and K Huber, *Z. Phys. Chem.*, 226 (2012) 837-854 
     74L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    7475 
    7576 
    … …  
    110111have_Fq = True 
    111112effective_radius_type = [ 
    112     "equivalent sphere", "half length_a", "half length_b", "half length_c", 
     113    "equivalent cylinder excluded volume", "equivalent outer volume sphere",  
     114    "half length_a", "half length_b", "half length_c", 
    113115    "equivalent outer circular cross-section", 
    114116    "half ab diagonal", "half diagonal", 

  • sasmodels/models/parallelepiped.c

    rd42dd4a r99658f6  
    66 
    77static double 
     8radius_from_excluded_volume(double length_a, double length_b, double length_c) 
     9{ 
     10    double r_equiv, length; 
     11    double lengths[3] = {length_a, length_b, length_c}; 
     12    double lengthmax = fmax(lengths[0],fmax(lengths[1],lengths[2])); 
     13    double length_1 = lengthmax; 
     14    double length_2 = lengthmax; 
     15    double length_3 = lengthmax; 
     16 
     17    for(int ilen=0; ilen<3; ilen++) { 
     18        if (lengths[ilen] < length_1) { 
     19            length_2 = length_1; 
     20            length_1 = lengths[ilen]; 
     21            } else { 
     22                if (lengths[ilen] < length_2) { 
     23                        length_2 = lengths[ilen]; 
     24                } 
     25            } 
     26    } 
     27    if(length_2-length_1 > length_3-length_2) { 
     28        r_equiv = sqrt(length_2*length_3/M_PI); 
     29        length  = length_1; 
     30    } else  { 
     31        r_equiv = sqrt(length_1*length_2/M_PI); 
     32        length  = length_3; 
     33    } 
     34 
     35    return 0.5*cbrt(0.75*r_equiv*(2.0*r_equiv*length + (r_equiv + length)*(M_PI*r_equiv + length))); 
     36} 
     37 
     38static double 
    839effective_radius(int mode, double length_a, double length_b, double length_c) 
    940{ 
    1041    switch (mode) { 
    1142    default: 
    12     case 1: // equivalent sphere 
     43    case 1: // equivalent cylinder excluded volume 
     44        return radius_from_excluded_volume(length_a,length_b,length_c); 
     45    case 2: // equivalent volume sphere 
    1346        return cbrt(length_a*length_b*length_c/M_4PI_3); 
    14     case 2: // half length_a 
     47    case 3: // half length_a 
    1548        return 0.5 * length_a; 
    16     case 3: // half length_b 
     49    case 4: // half length_b 
    1750        return 0.5 * length_b; 
    18     case 4: // half length_c 
     51    case 5: // half length_c 
    1952        return 0.5 * length_c; 
    20     case 5: // equivalent circular cross-section 
     53    case 6: // equivalent circular cross-section 
    2154        return sqrt(length_a*length_b/M_PI); 
    22     case 6: // half ab diagonal 
     55    case 7: // half ab diagonal 
    2356        return 0.5*sqrt(length_a*length_a + length_b*length_b); 
    24     case 7: // half diagonal 
     57    case 8: // half diagonal 
    2558        return 0.5*sqrt(length_a*length_a + length_b*length_b + length_c*length_c); 
    2659    } 

  • sasmodels/models/parallelepiped.py

    ree60aa7 r99658f6  
    180180   14 (1961) 185-211 
    181181.. [#] R Nayuk and K Huber, *Z. Phys. Chem.*, 226 (2012) 837-854 
     182L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    182183 
    183184Authorship and Verification 
    … …  
    232233have_Fq = True 
    233234effective_radius_type = [ 
    234     "equivalent sphere", "half length_a", "half length_b", "half length_c", 
     235    "equivalent cylinder excluded volume", "equivalent volume sphere",  
     236    "half length_a", "half length_b", "half length_c", 
    235237    "equivalent circular cross-section", "half ab diagonal", "half diagonal", 
    236238    ] 

  • sasmodels/models/pearl_necklace.c

    r3f853beb r9b5fd42  
    4040    const double si = sas_sinx_x(q*A_s); 
    4141    const double omsi = 1.0 - si; 
    42     const double pow_si = pow(si, num_pearls); 
     42    const double pow_si = pown(si, num_pearls); 
    4343 
    4444    // form factor for num_pearls 
    … …  
    6767} 
    6868 
    69 double form_volume(double radius, double edge_sep, 
    70     double thick_string, double fp_num_pearls) 
     69double form_volume(double radius, double edge_sep, double thick_string, double fp_num_pearls) 
    7170{ 
    7271    const int num_pearls = (int)(fp_num_pearls + 0.5); //Force integer number of pearls 
    … …  
    7776 
    7877    return volume; 
     78} 
     79 
     80static double 
     81radius_from_volume(double radius, double edge_sep, double thick_string, double fp_num_pearls) 
     82{ 
     83    const double vol_tot = form_volume(radius, edge_sep, thick_string, fp_num_pearls); 
     84    return cbrt(vol_tot/M_4PI_3); 
     85} 
     86 
     87static double 
     88effective_radius(int mode, double radius, double edge_sep, double thick_string, double fp_num_pearls) 
     89{ 
     90    switch (mode) { 
     91    default: 
     92    case 1: 
     93        return radius_from_volume(radius, edge_sep, thick_string, fp_num_pearls); 
     94    } 
    7995} 
    8096 

  • sasmodels/models/pearl_necklace.py

    r2cc8aa2 rcf3d0ce  
    5353R Schweins and K Huber, *Particle Scattering Factor of Pearl Necklace Chains*, 
    5454*Macromol. Symp.* 211 (2004) 25-42 2004 
     55L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
    5556""" 
    5657 
    … …  
    9596source = ["lib/sas_Si.c", "lib/sas_3j1x_x.c", "pearl_necklace.c"] 
    9697single = False  # use double precision unless told otherwise 
    97  
    98 def volume(radius, edge_sep, thick_string, num_pearls): 
    99     """ 
    100     Calculates the total particle volume of the necklace. 
    101     Redundant with form_volume. 
    102     """ 
    103     num_pearls = int(num_pearls + 0.5) 
    104     number_of_strings = num_pearls - 1.0 
    105     string_vol = edge_sep * pi * pow((thick_string / 2.0), 2.0) 
    106     pearl_vol = 4.0 /3.0 * pi * pow(radius, 3.0) 
    107     total_vol = number_of_strings * string_vol 
    108     total_vol += num_pearls * pearl_vol 
    109     return total_vol 
    110  
    111 def ER(radius, edge_sep, thick_string, num_pearls): 
    112     """ 
    113     Calculation for effective radius. 
    114     """ 
    115     num_pearls = int(num_pearls + 0.5) 
    116     tot_vol = volume(radius, edge_sep, thick_string, num_pearls) 
    117     rad_out = (tot_vol/(4.0/3.0*pi)) ** (1./3.) 
    118     return rad_out 
    119  
     98effective_radius_type = [ 
     99    "equivalent volume sphere",  
     100    ] 
     101     
    120102def random(): 
    121103    radius = 10**np.random.uniform(1, 3) # 1 - 1000 

  • sasmodels/models/pringle.c

    rd42dd4a r99658f6  
    105105 
    106106static double 
     107radius_from_excluded_volume(double radius, double thickness) 
     108{ 
     109    return 0.5*cbrt(0.75*radius*(2.0*radius*thickness + (radius + thickness)*(M_PI*radius + thickness))); 
     110} 
     111 
     112static double 
    107113effective_radius(int mode, double radius, double thickness, double alpha, double beta) 
    108114{ 
    109115    switch (mode) { 
    110116    default: 
    111     case 1: // equivalent sphere 
     117    case 1: // equivalent cylinder excluded volume 
     118        return radius_from_excluded_volume(radius, thickness); 
     119    case 2: // equivalent volume sphere 
    112120        return cbrt(M_PI*radius*radius*thickness/M_4PI_3); 
    113     case 2: // radius 
     121    case 3: // radius 
    114122        return radius; 
    115123    } 

  • sasmodels/models/pringle.py

    ree60aa7 r99658f6  
    4242Karen Edler, Universtiy of Bath, Private Communication. 2012. 
    4343Derivation by Stefan Alexandru Rautu. 
     44L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949). 
    4445 
    4546* **Author:** Andrew Jackson **Date:** 2008 
    … …  
    7475source = ["lib/polevl.c", "lib/sas_J0.c", "lib/sas_J1.c", 
    7576          "lib/sas_JN.c", "lib/gauss76.c", "pringle.c"] 
    76 effective_radius_type = ["equivalent sphere", "radius"] 
     77effective_radius_type = ["equivalent cylinder excluded volume", "equivalent volume sphere", "radius"] 
    7778 
    7879def random(): 

  • sasmodels/models/rectangular_prism.c

    rd42dd4a r99658f6  
    66 
    77static double 
     8radius_from_excluded_volume(double length_a, double b2a_ratio, double c2a_ratio) 
     9{ 
     10    double const r_equiv   = sqrt(length_a*length_a*b2a_ratio/M_PI); 
     11    double const length_c  = c2a_ratio*length_a; 
     12    return 0.5*cbrt(0.75*r_equiv*(2.0*r_equiv*length_c + (r_equiv + length_c)*(M_PI*r_equiv + length_c))); 
     13} 
     14 
     15static double 
    816effective_radius(int mode, double length_a, double b2a_ratio, double c2a_ratio) 
    917{ 
    1018    switch (mode) { 
    1119    default: 
    12     case 1: // equivalent sphere 
     20    case 1: // equivalent cylinder excluded volume 
     21        return radius_from_excluded_volume(length_a,b2a_ratio,c2a_ratio); 
     22    case 2: // equivalent volume sphere 
    1323        return cbrt(cube(length_a)*b2a_ratio*c2a_ratio/M_4PI_3); 
    14     case 2: // half length_a 
     24    case 3: // half length_a 
    1525        return 0.5 * length_a; 
    16     case 3: // half length_b 
     26    case 4: // half length_b 
    1727        return 0.5 * length_a*b2a_ratio; 
    18     case 4: // half length_c 
     28    case 5: // half length_c 
    1929        return 0.5 * length_a*c2a_ratio; 
    20     case 5: // equivalent circular cross-section 
     30    case 6: // equivalent circular cross-section 
    2131        return length_a*sqrt(b2a_ratio/M_PI); 
    22     case 6: // half ab diagonal 
     32    case 7: // half ab diagonal 
    2333        return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio))); 
    24     case 7: // half diagonal 
     34    case 8: // half diagonal 
    2535        return 0.5*sqrt(square(length_a) * (1.0 + square(b2a_ratio) + square(c2a_ratio))); 
    2636    } 

  • sasmodels/models/rectangular_prism.py

    ree60aa7 r99658f6  
    9999 
    100100R Nayuk and K Huber, *Z. Phys. Chem.*, 226 (2012) 837-854 
     101 
     102L. Onsager, Ann. New York Acad. Sci. 51, 627-659 (1949).  
     103 
    101104""" 
    102105 
    … …  
    137140have_Fq = True 
    138141effective_radius_type = [ 
    139     "equivalent sphere", "half length_a", "half length_b", "half length_c", 
     142    "equivalent cylinder excluded volume", "equivalent volume sphere",  
     143    "half length_a", "half length_b", "half length_c", 
    140144    "equivalent circular cross-section", "half ab diagonal", "half diagonal", 
    141145    ] 

  • sasmodels/models/triaxial_ellipsoid.c

    rd42dd4a r99658f6  
    55{ 
    66    return M_4PI_3*radius_equat_minor*radius_equat_major*radius_polar; 
     7} 
     8 
     9static double 
     10radius_from_curvature(double radius_equat_minor, double radius_equat_major, double radius_polar) 
     11{ 
     12    // Trivial cases 
     13    if (radius_equat_minor == radius_equat_major == radius_polar) return radius_polar; 
     14    if (radius_equat_minor * radius_equat_major * radius_polar == 0.)  return 0.; 
     15 
     16 
     17    double r_equat_equiv, r_polar_equiv; 
     18    double radii[3] = {radius_equat_minor, radius_equat_major, radius_polar}; 
     19    double radmax = fmax(radii[0],fmax(radii[1],radii[2])); 
     20 
     21    double radius_1 = radmax; 
     22    double radius_2 = radmax; 
     23    double radius_3 = radmax; 
     24 
     25    for(int irad=0; irad<3; irad++) { 
     26        if (radii[irad] < radius_1) { 
     27            radius_3 = radius_2; 
     28            radius_2 = radius_1; 
     29            radius_1 = radii[irad]; 
     30            } else { 
     31                if (radii[irad] < radius_2) { 
     32                        radius_2 = radii[irad]; 
     33                } 
     34            } 
     35    } 
     36    if(radius_2-radius_1 > radius_3-radius_2) { 
     37        r_equat_equiv = sqrt(radius_2*radius_3); 
     38        r_polar_equiv = radius_1; 
     39    } else  { 
     40        r_equat_equiv = sqrt(radius_1*radius_2); 
     41        r_polar_equiv = radius_3; 
     42    } 
     43 
     44    // see equation (26) in A.Isihara, J.Chem.Phys. 18(1950)1446-1449 
     45    const double ratio = (r_polar_equiv < r_equat_equiv 
     46                          ? r_polar_equiv / r_equat_equiv 
     47                          : r_equat_equiv / r_polar_equiv); 
     48    const double e1 = sqrt(1.0 - ratio*ratio); 
     49    const double b1 = 1.0 + asin(e1) / (e1 * ratio); 
     50    const double bL = (1.0 + e1) / (1.0 - e1); 
     51    const double b2 = 1.0 + 0.5 * ratio * ratio / e1 * log(bL); 
     52    const double delta = 0.75 * b1 * b2; 
     53    const double ddd = 2.0 * (delta + 1.0) * r_polar_equiv * r_equat_equiv * r_equat_equiv; 
     54    return 0.5 * cbrt(ddd); 
    755} 
    856 
    … …  
    3280    switch (mode) { 
    3381    default: 
    34     case 1: // equivalent sphere 
     82    case 1: // equivalent biaxial ellipsoid average curvature 
     83        return radius_from_curvature(radius_equat_minor,radius_equat_major, radius_polar); 
     84    case 2: // equivalent volume sphere 
    3585        return radius_from_volume(radius_equat_minor,radius_equat_major, radius_polar); 
    36     case 2: // min radius 
     86    case 3: // min radius 
    3787        return radius_from_min_dimension(radius_equat_minor,radius_equat_major, radius_polar); 
    38     case 3: // max radius 
     88    case 4: // max radius 
    3989        return radius_from_max_dimension(radius_equat_minor,radius_equat_major, radius_polar); 
    4090    } 

  • sasmodels/models/triaxial_ellipsoid.py

    ree60aa7 r99658f6  
    158158source = ["lib/sas_3j1x_x.c", "lib/gauss76.c", "triaxial_ellipsoid.c"] 
    159159have_Fq = True 
    160 effective_radius_type = ["equivalent sphere", "min radius", "max radius"] 
     160effective_radius_type = ["equivalent biaxial ellipsoid average curvature", "equivalent volume sphere", "min radius", "max radius"] 
    161161 
    162162def random(): 

  • sasmodels/product.py

    r39a06c9 r99658f6  
    3737#] 
    3838 
     39STRUCTURE_MODE_ID = "structure_factor_mode" 
     40RADIUS_MODE_ID = "radius_effective_mode" 
    3941RADIUS_ID = "radius_effective" 
    4042VOLFRAC_ID = "volfraction" 
    … …  
    4345    if p_info.have_Fq: 
    4446        par = parse_parameter( 
    45                 "structure_factor_mode", 
     47                STRUCTURE_MODE_ID, 
    4648                "", 
    4749                0, 
    … …  
    5254    if p_info.effective_radius_type is not None: 
    5355        par = parse_parameter( 
    54                 "radius_effective_mode", 
     56                RADIUS_MODE_ID, 
    5557                "", 
    56                 0, 
     58                1, 
    5759                [["unconstrained"] + p_info.effective_radius_type], 
    5860                "", 

  • sasmodels/resolution.py

    r9e7837a re2592f0  
    445445    q = np.sort(q) 
    446446    if q_min + 2*MINIMUM_RESOLUTION < q[0]: 
    447         n_low = np.ceil((q[0]-q_min) / (q[1]-q[0])) if q[1] > q[0] else 15 
     447        n_low = int(np.ceil((q[0]-q_min) / (q[1]-q[0]))) if q[1] > q[0] else 15 
    448448        q_low = np.linspace(q_min, q[0], n_low+1)[:-1] 
    449449    else: 
    450450        q_low = [] 
    451451    if q_max - 2*MINIMUM_RESOLUTION > q[-1]: 
    452         n_high = np.ceil((q_max-q[-1]) / (q[-1]-q[-2])) if q[-1] > q[-2] else 15 
     452        n_high = int(np.ceil((q_max-q[-1]) / (q[-1]-q[-2]))) if q[-1] > q[-2] else 15 
    453453        q_high = np.linspace(q[-1], q_max, n_high+1)[1:] 
    454454    else: 
    … …  
    499499            q_min = q[0]*MINIMUM_ABSOLUTE_Q 
    500500        n_low = log_delta_q * (log(q[0])-log(q_min)) 
    501         q_low = np.logspace(log10(q_min), log10(q[0]), np.ceil(n_low)+1)[:-1] 
     501        q_low = np.logspace(log10(q_min), log10(q[0]), int(np.ceil(n_low))+1)[:-1] 
    502502    else: 
    503503        q_low = [] 
    504504    if q_max > q[-1]: 
    505505        n_high = log_delta_q * (log(q_max)-log(q[-1])) 
    506         q_high = np.logspace(log10(q[-1]), log10(q_max), np.ceil(n_high)+1)[1:] 
     506        q_high = np.logspace(log10(q[-1]), log10(q_max), int(np.ceil(n_high))+1)[1:] 
    507507    else: 
    508508        q_high = [] 

  • sasmodels/sasview_model.py

    r39a06c9 ra8a1d48  
    382382            hidden.add('scale') 
    383383            hidden.add('background') 
    384             self._model_info.parameters.defaults['background'] = 0. 
    385384 
    386385        # Update the parameter lists to exclude any hidden parameters 
    … …  
    695694            return self._calculate_Iq(qx, qy) 
    696695 
    697     def _calculate_Iq(self, qx, qy=None, Fq=False, effective_radius_type=1): 
     696    def _calculate_Iq(self, qx, qy=None): 
    698697        if self._model is None: 
    699698            self._model = core.build_model(self._model_info) 
    … …  
    715714        #print("values", values) 
    716715        #print("is_mag", is_magnetic) 
    717         if Fq: 
    718             result = calculator.Fq(call_details, values, cutoff=self.cutoff, 
    719                                    magnetic=is_magnetic, 
    720                                    effective_radius_type=effective_radius_type) 
    721716        result = calculator(call_details, values, cutoff=self.cutoff, 
    722717                            magnetic=is_magnetic) 
    … …  
    736731        Calculate the effective radius for P(q)*S(q) 
    737732 
     733        *mode* is the R_eff type, which defaults to 1 to match the ER 
     734        calculation for sasview models from version 3.x. 
     735 
    738736        :return: the value of the effective radius 
    739737        """ 
    740         Fq = self._calculate_Iq([0.1], True, mode) 
    741         return Fq[2] 
     738        # ER and VR are only needed for old multiplication models, based on 
     739        # sas.sascalc.fit.MultiplicationModel.  Fail for now.  If we want to 
     740        # continue supporting them then add some test cases so that the code 
     741        # is exercised.  We can access ER/VR using the kernel Fq function by 
     742        # extending _calculate_Iq so that it calls: 
     743        #    if er_mode > 0: 
     744        #        res = calculator.Fq(call_details, values, cutoff=self.cutoff, 
     745        #                            magnetic=False, effective_radius_type=mode) 
     746        #        R_eff, form_shell_ratio = res[2], res[4] 
     747        #        return R_eff, form_shell_ratio 
     748        # Then use the following in calculate_ER: 
     749        #    ER, VR = self._calculate_Iq(q=[0.1], er_mode=mode) 
     750        #    return ER 
     751        # Similarly, for calculate_VR: 
     752        #    ER, VR = self._calculate_Iq(q=[0.1], er_mode=1) 
     753        #    return VR 
     754        # Obviously a combined calculate_ER_VR method would be better, but 
     755        # we only need them to support very old models, so ignore the 2x 
     756        # performance hit. 
     757        raise NotImplementedError("ER function is no longer available.") 
    742758 
    743759    def calculate_VR(self): 
    … …  
    748764        :return: the value of the form:shell volume ratio 
    749765        """ 
    750         Fq = self._calculate_Iq([0.1], True, mode) 
    751         return Fq[4] 
     766        # See comments in calculate_ER. 
     767        raise NotImplementedError("VR function is no longer available.") 
    752768 
    753769    def set_dispersion(self, parameter, dispersion): 
    … …  
    819835            return value, [value], [1.0] 
    820836 
     837    @classmethod 
     838    def runTests(cls): 
     839        """ 
     840        Run any tests built into the model and captures the test output. 
     841 
     842        Returns success flag and output 
     843        """ 
     844        from .model_test import check_model 
     845        return check_model(cls._model_info) 
     846 
    821847def test_cylinder(): 
    822848    # type: () -> float 
    … …  
    849875    P = _make_standard_model('cylinder')() 
    850876    model = MultiplicationModel(P, S) 
    851     model.setParam('radius_effective_mode', 1.0) 
     877    model.setParam(product.RADIUS_MODE_ID, 1.0) 
    852878    value = model.evalDistribution([0.1, 0.1]) 
    853879    if np.isnan(value): 
    … …  
    904930    CylinderModel().evalDistribution([0.1, 0.1]) 
    905931 
     932def test_structure_factor_background(): 
     933    # type: () -> None 
     934    """ 
     935    Check that sasview model and direct model match, with background=0. 
     936    """ 
     937    from .data import empty_data1D 
     938    from .core import load_model_info, build_model 
     939    from .direct_model import DirectModel 
     940 
     941    model_name = "hardsphere" 
     942    q = [0.0] 
     943 
     944    sasview_model = _make_standard_model(model_name)() 
     945    sasview_value = sasview_model.evalDistribution(np.array(q))[0] 
     946 
     947    data = empty_data1D(q) 
     948    model_info = load_model_info(model_name) 
     949    model = build_model(model_info) 
     950    direct_model = DirectModel(data, model) 
     951    direct_value_zero_background = direct_model(background=0.0) 
     952 
     953    assert sasview_value == direct_value_zero_background 
     954 
     955    # Additionally check that direct value background defaults to zero 
     956    direct_value_default = direct_model() 
     957    assert sasview_value == direct_value_default 
     958 
     959 
    906960def magnetic_demo(): 
    907961    Model = _make_standard_model('sphere') 
    … …  
    924978    #print("rpa:", test_rpa()) 
    925979    #test_empty_distribution() 
     980    #test_structure_factor_background() 

  • sasmodels/special.py

    rdf69efa rfba9ca0  
    113113        The standard math function, tgamma(x) is unstable for $x < 1$ 
    114114        on some platforms. 
     115 
     116    sas_gammaln(x): 
     117        log gamma function sas_gammaln\ $(x) = \log \Gamma(|x|)$. 
     118 
     119        The standard math function, lgamma(x), is incorrect for single 
     120        precision on some platforms. 
     121 
     122    sas_gammainc(a, x), sas_gammaincc(a, x): 
     123        Incomplete gamma function 
     124        sas_gammainc\ $(a, x) = \int_0^x t^{a-1}e^{-t}\,dt / \Gamma(a)$ 
     125        and complementary incomplete gamma function 
     126        sas_gammaincc\ $(a, x) = \int_x^\infty t^{a-1}e^{-t}\,dt / \Gamma(a)$ 
    115127 
    116128    sas_erf(x), sas_erfc(x): 
    … …  
    207219from numpy import pi, nan, inf 
    208220from scipy.special import gamma as sas_gamma 
     221from scipy.special import gammaln as sas_gammaln 
     222from scipy.special import gammainc as sas_gammainc 
     223from scipy.special import gammaincc as sas_gammaincc 
    209224from scipy.special import erf as sas_erf 
    210225from scipy.special import erfc as sas_erfc 

  • sasmodels/weights.py

    r3d58247 re2592f0  
    2323    default = dict(npts=35, width=0, nsigmas=3) 
    2424    def __init__(self, npts=None, width=None, nsigmas=None): 
    25         self.npts = self.default['npts'] if npts is None else npts 
     25        self.npts = self.default['npts'] if npts is None else int(npts) 
    2626        self.width = self.default['width'] if width is None else width 
    2727        self.nsigmas = self.default['nsigmas'] if nsigmas is None else nsigmas 
Note: See TracChangeset for help on using the changeset viewer.