source: sasmodels/doc/guide/plugin.rst @ 93fbc34

core_shell_microgelsmagnetic_modelticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since 93fbc34 was c654160, checked in by GitHub <noreply@…>, 6 years ago

Corrected typo in plugin.rst

  • Property mode set to 100644
File size: 45.9 KB
RevLine 
[990d8df]1.. _Writing_a_Plugin:
2
3Writing a Plugin Model
4======================
5
6Overview
7^^^^^^^^
8
9In addition to the models provided with the sasmodels package, you are free to
10create your own models.
11
12Models can be of three types:
13
14- A pure python model : Example -
15  `broadpeak.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/broad_peak.py>`_
16
17- A python model with embedded C : Example -
18  `sphere.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/sphere.py>`_
19
20- A python wrapper with separate C code : Example -
21  `cylinder.py <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/cylinder.py>`_,
22  `cylinder.c <https://github.com/SasView/sasmodels/blob/master/sasmodels/models/cylinder.c>`_
23
24When using SasView, plugin models should be saved to the SasView
25*plugin_models* folder *C:\\Users\\{username}\\.sasview\\plugin_models*
26(on Windows) or */Users/{username}/.sasview\\plugin_models* (on Mac).
27The next time SasView is started it will compile the plugin and add
28it to the list of *Plugin Models* in a FitPage.  Scripts can load
29the models from anywhere.
30
31The built-in modules are available in the *models* subdirectory
32of the sasmodels package.  For SasView on Windows, these will
33be found in *C:\\Program Files (x86)\\SasView\\sasmodels-data\\models*.
34On Mac OSX, these will be within the application bundle as
35*/Applications/SasView 4.0.app/Contents/Resources/sasmodels-data/models*.
36
37Other models are available for download from the
38`Model Marketplace <http://marketplace.sasview.org/>`_. You can contribute your
39own models to the Marketplace as well.
40
41Create New Model Files
42^^^^^^^^^^^^^^^^^^^^^^
43
44Copy the appropriate files to your plugin models directory (we recommend
45using the examples above as templates) as mymodel.py (and mymodel.c, etc)
46as required, where "mymodel" is the name for the model you are creating.
47
48*Please follow these naming rules:*
49
50- No capitalization and thus no CamelCase
51- If necessary use underscore to separate words (i.e. barbell not BarBell or
52  broad_peak not BroadPeak)
53- Do not include "model" in the name (i.e. barbell not BarBellModel)
54
55
56Edit New Model Files
57^^^^^^^^^^^^^^^^^^^^
58
59Model Contents
60..............
61
62The model interface definition is in the .py file.  This file contains:
63
64- a **model name**:
65   - this is the **name** string in the *.py* file
66   - titles should be:
67
68    - all in *lower* case
69    - without spaces (use underscores to separate words instead)
70    - without any capitalization or CamelCase
71    - without incorporating the word "model"
72    - examples: *barbell* **not** *BarBell*; *broad_peak* **not** *BroadPeak*;
73      *barbell* **not** *BarBellModel*
74
75- a **model title**:
76   - this is the **title** string in the *.py* file
77   - this is a one or two line description of the model, which will appear
78     at the start of the model documentation and as a tooltip in the SasView GUI
79
[3048ec6]80- a **short description**:
[990d8df]81   - this is the **description** string in the *.py* file
82   - this is a medium length description which appears when you click
83     *Description* on the model FitPage
84
85- a **parameter table**:
86   - this will be auto-generated from the *parameters* in the *.py* file
87
88- a **long description**:
89   - this is ReStructuredText enclosed between the r""" and """ delimiters
90     at the top of the *.py* file
91   - what you write here is abstracted into the SasView help documentation
92   - this is what other users will refer to when they want to know what
93     your model does; so please be helpful!
94
95- a **definition** of the model:
96   - as part of the **long description**
97
98- a **formula** defining the function the model calculates:
99   - as part of the **long description**
100
101- an **explanation of the parameters**:
102   - as part of the **long description**
103   - explaining how the symbols in the formula map to the model parameters
104
105- a **plot** of the function, with a **figure caption**:
106   - this is automatically generated from your default parameters
107
108- at least one **reference**:
109   - as part of the **long description**
110   - specifying where the reader can obtain more information about the model
111
112- the **name of the author**
113   - as part of the **long description**
114   - the *.py* file should also contain a comment identifying *who*
115     converted/created the model file
116
117Models that do not conform to these requirements will *never* be incorporated
118into the built-in library.
119
120
121Model Documentation
122...................
123
124The *.py* file starts with an r (for raw) and three sets of quotes
125to start the doc string and ends with a second set of three quotes.
126For example::
127
128    r"""
129    Definition
130    ----------
131
132    The 1D scattering intensity of the sphere is calculated in the following
133    way (Guinier, 1955)
134
135    .. math::
136
137        I(q) = \frac{\text{scale}}{V} \cdot \left[
138            3V(\Delta\rho) \cdot \frac{\sin(qr) - qr\cos(qr))}{(qr)^3}
139            \right]^2 + \text{background}
140
141    where *scale* is a volume fraction, $V$ is the volume of the scatterer,
142    $r$ is the radius of the sphere and *background* is the background level.
143    *sld* and *sld_solvent* are the scattering length densities (SLDs) of the
144    scatterer and the solvent respectively, whose difference is $\Delta\rho$.
145
146    You can included figures in your documentation, as in the following
147    figure for the cylinder model.
148
149    .. figure:: img/cylinder_angle_definition.jpg
150
151        Definition of the angles for oriented cylinders.
152
153    References
154    ----------
155
156    A Guinier, G Fournet, *Small-Angle Scattering of X-Rays*,
157    John Wiley and Sons, New York, (1955)
158    """
159
160This is where the FULL documentation for the model goes (to be picked up by
161the automatic documentation system).  Although it feels odd, you
162should start the documentation immediately with the **definition**---the model
163name, a brief description and the parameter table are automatically inserted
164above the definition, and the a plot of the model is automatically inserted
165before the **reference**.
166
167Figures can be included using the *figure* command, with the name
168of the *.png* file containing the figure and a caption to appear below the
169figure.  Figure numbers will be added automatically.
170
171See this `Sphinx cheat sheet <http://matplotlib.org/sampledoc/cheatsheet.html>`_
172for a quick guide to the documentation layout commands, or the
173`Sphinx Documentation <http://www.sphinx-doc.org/en/stable/>`_ for
174complete details.
175
176The model should include a **formula** written using LaTeX markup.
177The example above uses the *math* command to make a displayed equation.  You
178can also use *\$formula\$* for an inline formula. This is handy for defining
179the relationship between the model parameters and formula variables, such
180as the phrase "\$r\$ is the radius" used above.  The live demo MathJax
181page `<http://www.mathjax.org/>`_ is handy for checking that the equations
182will look like you intend.
183
184Math layout uses the `amsmath <http://www.ams.org/publications/authors/tex/amslatex>`_
185package for aligning equations (see amsldoc.pdf on that page for complete
186documentation). You will automatically be in an aligned environment, with
187blank lines separating the lines of the equation.  Place an ampersand before
188the operator on which to align.  For example::
189
190    .. math::
191
192      x + y &= 1 \\
193      y &= x - 1
194
195produces
196
197.. math::
198
199      x + y &= 1 \\
200      y &= x - 1
201
202If you need more control, use::
203
204    .. math::
205        :nowrap:
206
207
208Model Definition
209................
210
211Following the documentation string, there are a series of definitions::
212
213    name = "sphere"  # optional: defaults to the filename without .py
214
215    title = "Spheres with uniform scattering length density"
216
217    description = """\
218    P(q)=(scale/V)*[3V(sld-sld_solvent)*(sin(qr)-qr cos(qr))
219                    /(qr)^3]^2 + background
220        r: radius of sphere
221        V: The volume of the scatter
222        sld: the SLD of the sphere
223        sld_solvent: the SLD of the solvent
224    """
225
226    category = "shape:sphere"
227
228    single = True   # optional: defaults to True
229
230    opencl = False  # optional: defaults to False
231
232    structure_factor = False  # optional: defaults to False
233
234**name = "mymodel"** defines the name of the model that is shown to the user.
[3048ec6]235If it is not provided it will use the name of the model file. The name must
236be a valid variable name, starting with a letter and contains only letters,
237numbers or underscore.  Spaces, dashes, and other symbols are not permitted.
[990d8df]238
239**title = "short description"** is short description of the model which
240is included after the model name in the automatically generated documentation.
241The title can also be used for a tooltip.
242
243**description = """doc string"""** is a longer description of the model. It
244shows up when you press the "Description" button of the SasView FitPage.
245It should give a brief description of the equation and the parameters
246without the need to read the entire model documentation. The triple quotes
247allow you to write the description over multiple lines. Keep the lines
248short since the GUI will wrap each one separately if they are too long.
249**Make sure the parameter names in the description match the model definition!**
250
251**category = "shape:sphere"** defines where the model will appear in the
252model documentation.  In this example, the model will appear alphabetically
253in the list of spheroid models in the *Shape* category.
254
255**single = True** indicates that the model can be run using single
256precision floating point values.  Set it to False if the numerical
257calculation for the model is unstable, which is the case for about 20 of
258the built in models.  It is worthwhile modifying the calculation to support
259single precision, allowing models to run up to 10 times faster.  The
260section `Test_Your_New_Model`_  describes how to compare model values for
261single vs. double precision so you can decide if you need to set
262single to False.
263
264**opencl = False** indicates that the model should not be run using OpenCL.
265This may be because the model definition includes code that cannot be
266compiled for the GPU (for example, goto statements).  It can also be used
267for large models which can't run on most GPUs.  This flag has not been
268used on any of the built in models; models which were failing were
269streamlined so this flag was not necessary.
270
271**structure_factor = True** indicates that the model can be used as a
272structure factor to account for interactions between particles.  See
273`Form_Factors`_ for more details.
274
275Model Parameters
276................
277
278Next comes the parameter table.  For example::
279
280    # pylint: disable=bad-whitespace, line-too-long
281    #   ["name",        "units", default, [min, max], "type",    "description"],
282    parameters = [
283        ["sld",         "1e-6/Ang^2",  1, [-inf, inf], "sld",    "Layer scattering length density"],
284        ["sld_solvent", "1e-6/Ang^2",  6, [-inf, inf], "sld",    "Solvent scattering length density"],
285        ["radius",      "Ang",        50, [0, inf],    "volume", "Sphere radius"],
286    ]
287    # pylint: enable=bad-whitespace, line-too-long
288
289**parameters = [["name", "units", default, [min,max], "type", "tooltip"],...]**
290defines the parameters that form the model.
291
292**Note: The order of the parameters in the definition will be the order of the
293parameters in the user interface and the order of the parameters in Iq(),
[108e70e]294Iqac(), Iqabc() and form_volume(). And** *scale* **and** *background*
295**parameters are implicit to all models, so they do not need to be included
296in the parameter table.**
[990d8df]297
298- **"name"** is the name of the parameter shown on the FitPage.
299
[3048ec6]300  - the name must be a valid variable name, starting with a letter and
301    containing only letters, numbers and underscore.
302
[990d8df]303  - parameter names should follow the mathematical convention; e.g.,
304    *radius_core* not *core_radius*, or *sld_solvent* not *solvent_sld*.
305
306  - model parameter names should be consistent between different models,
307    so *sld_solvent*, for example, should have exactly the same name
308    in every model.
309
310  - to see all the parameter names currently in use, type the following in the
311    python shell/editor under the Tools menu::
312
313       import sasmodels.list_pars
314       sasmodels.list_pars.list_pars()
315
316    *re-use* as many as possible!!!
317
318  - use "name[n]" for multiplicity parameters, where *n* is the name of
319    the parameter defining the number of shells/layers/segments, etc.
320
321- **"units"** are displayed along with the parameter name
322
323  - every parameter should have units; use "None" if there are no units.
324
325  - **sld's should be given in units of 1e-6/Ang^2, and not simply
326    1/Ang^2 to be consistent with the builtin models.  Adjust your formulas
327    appropriately.**
328
329  - fancy units markup is available for some units, including::
330
331        Ang, 1/Ang, 1/Ang^2, 1e-6/Ang^2, degrees, 1/cm, Ang/cm, g/cm^3, mg/m^2
332
333  - the list of units is defined in the variable *RST_UNITS* within
334    `sasmodels/generate.py <https://github.com/SasView/sasmodels/tree/master/sasmodels/generate.py>`_
335
336    - new units can be added using the macros defined in *doc/rst_prolog*
337      in the sasmodels source.
338    - units should be properly formatted using sub-/super-scripts
339      and using negative exponents instead of the / operator, though
340      the unit name should use the / operator for consistency.
341    - please post a message to the SasView developers mailing list with your changes.
342
343- **default** is the initial value for the parameter.
344
345  - **the parameter default values are used to auto-generate a plot of
346    the model function in the documentation.**
347
348- **[min, max]** are the lower and upper limits on the parameter.
349
350  - lower and upper limits can be any number, or *-inf* or *inf*.
351
352  - the limits will show up as the default limits for the fit making it easy,
353    for example, to force the radius to always be greater than zero.
354
355  - these are hard limits defining the valid range of parameter values;
356    polydisperity distributions will be truncated at the limits.
357
358- **"type"** can be one of: "", "sld", "volume", or "orientation".
359
360  - "sld" parameters can have magnetic moments when fitting magnetic models;
361    depending on the spin polarization of the beam and the $q$ value being
362    examined, the effective sld for that material will be used to compute the
363    scattered intensity.
364
[108e70e]365  - "volume" parameters are passed to Iq(), Iqac(), Iqabc() and form_volume(),
366    and have polydispersity loops generated automatically.
[990d8df]367
[108e70e]368  - "orientation" parameters are not passed, but instead are combined with
369    orientation dispersity to translate *qx* and *qy* to *qa*, *qb* and *qc*.
370    These parameters should appear at the end of the table with the specific
371    names *theta*, *phi* and for asymmetric shapes *psi*, in that order.
[990d8df]372
[9844c3a]373Some models will have integer parameters, such as number of pearls in the
374pearl necklace model, or number of shells in the multi-layer vesicle model.
375The optimizers in BUMPS treat all parameters as floating point numbers which
376can take arbitrary values, even for integer parameters, so your model should
377round the incoming parameter value to the nearest integer inside your model
378you should round to the nearest integer.  In C code, you can do this using::
379
380    static double
381    Iq(double q, ..., double fp_n, ...)
382    {
383        int n = (int)(fp_n + 0.5);
384        ...
385    }
386
387in python::
388
389    def Iq(q, ..., n, ...):
390        n = int(n+0.5)
391        ...
392
[3048ec6]393Derivative based optimizers such as Levenberg-Marquardt will not work
[9844c3a]394for integer parameters since the partial derivative is always zero, but
395the remaining optimizers (DREAM, differential evolution, Nelder-Mead simplex)
396will still function.
397
[990d8df]398Model Computation
399.................
400
401Models can be defined as pure python models, or they can be a mixture of
402python and C models.  C models are run on the GPU if it is available,
403otherwise they are compiled and run on the CPU.
404
405Models are defined by the scattering kernel, which takes a set of parameter
406values defining the shape, orientation and material, and returns the
407expected scattering. Polydispersity and angular dispersion are defined
408by the computational infrastructure.  Any parameters defined as "volume"
409parameters are polydisperse, with polydispersity defined in proportion
410to their value.  "orientation" parameters use angular dispersion defined
411in degrees, and are not relative to the current angle.
412
413Based on a weighting function $G(x)$ and a number of points $n$, the
414computed value is
415
416.. math::
417
418     \hat I(q)
419     = \frac{\int G(x) I(q, x)\,dx}{\int G(x) V(x)\,dx}
420     \approx \frac{\sum_{i=1}^n G(x_i) I(q,x_i)}{\sum_{i=1}^n G(x_i) V(x_i)}
421
[3048ec6]422That is, the individual models do not need to include polydispersity
[990d8df]423calculations, but instead rely on numerical integration to compute the
[108e70e]424appropriately smeared pattern.
[990d8df]425
426Python Models
427.............
428
429For pure python models, define the *Iq* function::
430
431      import numpy as np
432      from numpy import cos, sin, ...
433
434      def Iq(q, par1, par2, ...):
435          return I(q, par1, par2, ...)
436      Iq.vectorized = True
437
438The parameters *par1, par2, ...* are the list of non-orientation parameters
439to the model in the order that they appear in the parameter table.
[3048ec6]440**Note that the auto-generated model file uses** *x* **rather than** *q*.
[990d8df]441
442The *.py* file should import trigonometric and exponential functions from
443numpy rather than from math.  This lets us evaluate the model for the whole
444range of $q$ values at once rather than looping over each $q$ separately in
445python.  With $q$ as a vector, you cannot use if statements, but must instead
446do tricks like
447
448::
449
450     a = x*q*(q>0) + y*q*(q<=0)
451
452or
453
454::
455
456     a = np.empty_like(q)
457     index = q>0
458     a[index] = x*q[index]
459     a[~index] = y*q[~index]
460
461which sets $a$ to $q \cdot x$ if $q$ is positive or $q \cdot y$ if $q$
462is zero or negative. If you have not converted your function to use $q$
463vectors, you can set the following and it will only receive one $q$
464value at a time::
465
466    Iq.vectorized = False
467
468Return np.NaN if the parameters are not valid (e.g., cap_radius < radius in
469barbell).  If I(q; pars) is NaN for any $q$, then those parameters will be
470ignored, and not included in the calculation of the weighted polydispersity.
471
472Models should define *form_volume(par1, par2, ...)* where the parameter
473list includes the *volume* parameters in order.  This is used for a weighted
474volume normalization so that scattering is on an absolute scale.  If
475*form_volume* is not defined, then the default *form_volume = 1.0* will be
476used.
477
478Embedded C Models
479.................
480
481Like pure python models, inline C models need to define an *Iq* function::
482
483    Iq = """
484        return I(q, par1, par2, ...);
485    """
486
487This expands into the equivalent C code::
488
489    #include <math.h>
490    double Iq(double q, double par1, double par2, ...);
491    double Iq(double q, double par1, double par2, ...)
492    {
493        return I(q, par1, par2, ...);
494    }
495
496*form_volume* defines the volume of the shape. As in python models, it
497includes only the volume parameters.
498
499**source=['fn.c', ...]** includes the listed C source files in the
[108e70e]500program before *Iq* and *form_volume* are defined. This allows you to
[ef85a09]501extend the library of C functions available to your model.
502
503*c_code* includes arbitrary C code into your kernel, which can be
504handy for defining helper functions for *Iq* and *form_volume*. Note that
[108e70e]505you can put the full function definition for *Iq* and *form_volume*
[ef85a09]506(include function declaration) into *c_code* as well, or put them into an
507external C file and add that file to the list of sources.
[990d8df]508
509Models are defined using double precision declarations for the
510parameters and return values.  When a model is run using single
511precision or long double precision, each variable is converted
512to the target type, depending on the precision requested.
513
514**Floating point constants must include the decimal point.**  This allows us
515to convert values such as 1.0 (double precision) to 1.0f (single precision)
516so that expressions that use these values are not promoted to double precision
517expressions.  Some graphics card drivers are confused when functions
518that expect floating point values are passed integers, such as 4*atan(1); it
519is safest to not use integers in floating point expressions.  Even better,
520use the builtin constant M_PI rather than 4*atan(1); it is faster and smaller!
521
522The C model operates on a single $q$ value at a time.  The code will be
523run in parallel across different $q$ values, either on the graphics card
524or the processor.
525
526Rather than returning NAN from Iq, you must define the *INVALID(v)*.  The
527*v* parameter lets you access all the parameters in the model using
528*v.par1*, *v.par2*, etc. For example::
529
530    #define INVALID(v) (v.bell_radius < v.radius)
531
[ef85a09]532The INVALID define can go into *Iq*, or *c_code*, or an external C file
533listed in *source*.
534
[108e70e]535Oriented Shapes
536...............
537
538If the scattering is dependent on the orientation of the shape, then you
539will need to include *orientation* parameters *theta*, *phi* and *psi*
540at the end of the parameter table.  Shape orientation uses *a*, *b* and *c*
541axes, corresponding to the *x*, *y* and *z* axes in the laboratory coordinate
542system, with *z* along the beam and *x*-*y* in the detector plane, with *x*
543horizontal and *y* vertical.  The *psi* parameter rotates the shape
544about its *c* axis, the *theta* parameter then rotates the *c* axis toward
545the *x* axis of the detector, then *phi* rotates the shape in the detector
546plane.  (Prior to these rotations, orientation dispersity will be applied
547as roll-pitch-yaw, rotating *c*, then *b* then *a* in the shape coordinate
548system.)  A particular *qx*, *qy* point on the detector, then corresponds
549to *qa*, *qb*, *qc* with respect to the shape.
550
551The oriented C model is called as *Iqabc(qa, qb, qc, par1, par2, ...)* where
552*par1*, etc. are the parameters to the model.  If the shape is rotationally
553symmetric about *c* then *psi* is not needed, and the model is called
554as *Iqac(qab, qc, par1, par2, ...)*.  In either case, the orientation
555parameters are not included in the function call.
556
557For 1D oriented shapes, an integral over all angles is usually needed for
[b85227d]558the *Iq* function. Given symmetry and the substitution $u = \cos(\alpha)$,
[108e70e]559$du = -\sin(\alpha)\,d\alpha$ this becomes
560
561.. math::
562
[b85227d]563    I(q) &= \frac{1}{4\pi} \int_{-\pi/2}^{pi/2} \int_{-pi}^{pi}
564            F(q_a, q_b, q_c)^2 \sin(\alpha)\,d\beta\,d\alpha \\
565        &= \frac{8}{4\pi} \int_{0}^{pi/2} \int_{0}^{\pi/2}
566            F^2 \sin(\alpha)\,d\beta\,d\alpha \\
567        &= \frac{8}{4\pi} \int_1^0 \int_{0}^{\pi/2} - F^2 \,d\beta\,du \\
568        &= \frac{8}{4\pi} \int_0^1 \int_{0}^{\pi/2} F^2 \,d\beta\,du
569
570for
571
572.. math::
573
574    q_a &= q \sin(\alpha)\sin(\beta) = q \sqrt{1-u^2} \sin(\beta) \\
575    q_b &= q \sin(\alpha)\cos(\beta) = q \sqrt{1-u^2} \cos(\beta) \\
576    q_c &= q \cos(\alpha) = q u
[108e70e]577
578Using the $z, w$ values for Gauss-Legendre integration in "lib/gauss76.c", the
579numerical integration is then::
580
581    double outer_sum = 0.0;
582    for (int i = 0; i < GAUSS_N; i++) {
583        const double cos_alpha = 0.5*GAUSS_Z[i] + 0.5;
584        const double sin_alpha = sqrt(1.0 - cos_alpha*cos_alpha);
585        const double qc = cos_alpha * q;
586        double inner_sum = 0.0;
587        for (int j = 0; j < GAUSS_N; j++) {
588            const double beta = M_PI_4 * GAUSS_Z[j] + M_PI_4;
589            double sin_beta, cos_beta;
590            SINCOS(beta, sin_beta, cos_beta);
591            const double qa = sin_alpha * sin_beta * q;
[b85227d]592            const double qb = sin_alpha * cos_beta * q;
593            const double form = Fq(qa, qb, qc, ...);
594            inner_sum += GAUSS_W[j] * form * form;
[108e70e]595        }
596        outer_sum += GAUSS_W[i] * inner_sum;
597    }
598    outer_sum *= 0.25; // = 8/(4 pi) * outer_sum * (pi/2) / 4
599
600The *z* values for the Gauss-Legendre integration extends from -1 to 1, so
601the double sum of *w[i]w[j]* explains the factor of 4.  Correcting for the
602average *dz[i]dz[j]* gives $(1-0) \cdot (\pi/2-0) = \pi/2$.  The $8/(4 \pi)$
603factor comes from the integral over the quadrant.  With less symmetry (eg.,
604in the bcc and fcc paracrystal models), then an integral over the entire
605sphere may be necessary.
606
607For simpler models which are rotationally symmetric a single integral
608suffices:
609
610.. math::
611
[b85227d]612    I(q) &= \frac{1}{\pi}\int_{-\pi/2}^{\pi/2}
613            F(q_{ab}, q_c)^2 \sin(\alpha)\,d\alpha/\pi \\
614        &= \frac{2}{\pi} \int_0^1 F^2\,du
615
616for
617
618.. math::
619
620    q_{ab} &= q \sin(\alpha) = q \sqrt{1 - u^2} \\
621    q_c &= q \cos(\alpha) = q u
622
[108e70e]623
624with integration loop::
625
626    double sum = 0.0;
627    for (int i = 0; i < GAUSS_N; i++) {
628        const double cos_alpha = 0.5*GAUSS_Z[i] + 0.5;
629        const double sin_alpha = sqrt(1.0 - cos_alpha*cos_alpha);
630        const double qab = sin_alpha * q;
[b85227d]631        const double qc = cos_alpha * q;
632        const double form = Fq(qab, qc, ...);
633        sum += GAUSS_W[j] * form * form;
[108e70e]634    }
635    sum *= 0.5; // = 2/pi * sum * (pi/2) / 2
636
637Magnetism
638.........
639
640Magnetism is supported automatically for all shapes by modifying the
641effective SLD of particle according to the Halpern-Johnson vector
[c654160]642describing the interaction between neutron spin and magnetic field.  All
[108e70e]643parameters marked as type *sld* in the parameter table are treated as
644possibly magnetic particles with magnitude *M0* and direction
645*mtheta* and *mphi*.  Polarization parameters are also provided
646automatically for magnetic models to set the spin state of the measurement.
647
648For more complicated systems where magnetism is not uniform throughout
649the individual particles, you will need to write your own models.
650You should not mark the nuclear sld as type *sld*, but instead leave
651them unmarked and provide your own magnetism and polarization parameters.
652For 2D measurements you will need $(q_x, q_y)$ values for the measurement
653to compute the proper magnetism and orientation, which you can implement
654using *Iqxy(qx, qy, par1, par2, ...)*.
655
[990d8df]656Special Functions
657.................
658
659The C code follows the C99 standard, with the usual math functions,
660as defined in
661`OpenCL <https://www.khronos.org/registry/cl/sdk/1.1/docs/man/xhtml/mathFunctions.html>`_.
662This includes the following:
663
664    M_PI, M_PI_2, M_PI_4, M_SQRT1_2, M_E:
665        $\pi$, $\pi/2$, $\pi/4$, $1/\sqrt{2}$ and Euler's constant $e$
[d0dc9a3]666    exp, log, pow(x,y), expm1, log1p, sqrt, cbrt:
667        Power functions $e^x$, $\ln x$, $x^y$, $e^x - 1$, $\ln 1 + x$,
668        $\sqrt{x}$, $\sqrt[3]{x}$. The functions expm1(x) and log1p(x)
669        are accurate across all $x$, including $x$ very close to zero.
[990d8df]670    sin, cos, tan, asin, acos, atan:
671        Trigonometry functions and inverses, operating on radians.
672    sinh, cosh, tanh, asinh, acosh, atanh:
673        Hyperbolic trigonometry functions.
674    atan2(y,x):
675        Angle from the $x$\ -axis to the point $(x,y)$, which is equal to
676        $\tan^{-1}(y/x)$ corrected for quadrant.  That is, if $x$ and $y$ are
677        both negative, then atan2(y,x) returns a value in quadrant III where
678        atan(y/x) would return a value in quadrant I. Similarly for
679        quadrants II and IV when $x$ and $y$ have opposite sign.
[d0dc9a3]680    fabs(x), fmin(x,y), fmax(x,y), trunc, rint:
[990d8df]681        Floating point functions.  rint(x) returns the nearest integer.
682    NAN:
683        NaN, Not a Number, $0/0$.  Use isnan(x) to test for NaN.  Note that
684        you cannot use :code:`x == NAN` to test for NaN values since that
[d0dc9a3]685        will always return false.  NAN does not equal NAN!  The alternative,
686        :code:`x != x` may fail if the compiler optimizes the test away.
[990d8df]687    INFINITY:
688        $\infty, 1/0$.  Use isinf(x) to test for infinity, or isfinite(x)
689        to test for finite and not NaN.
690    erf, erfc, tgamma, lgamma:  **do not use**
691        Special functions that should be part of the standard, but are missing
692        or inaccurate on some platforms. Use sas_erf, sas_erfc and sas_gamma
693        instead (see below). Note: lgamma(x) has not yet been tested.
694
695Some non-standard constants and functions are also provided:
696
697    M_PI_180, M_4PI_3:
698        $\frac{\pi}{180}$, $\frac{4\pi}{3}$
699    SINCOS(x, s, c):
700        Macro which sets s=sin(x) and c=cos(x). The variables *c* and *s*
701        must be declared first.
702    square(x):
703        $x^2$
704    cube(x):
705        $x^3$
706    sas_sinx_x(x):
707        $\sin(x)/x$, with limit $\sin(0)/0 = 1$.
708    powr(x, y):
709        $x^y$ for $x \ge 0$; this is faster than general $x^y$ on some GPUs.
710    pown(x, n):
711        $x^n$ for $n$ integer; this is faster than general $x^n$ on some GPUs.
712    FLOAT_SIZE:
713        The number of bytes in a floating point value.  Even though all
714        variables are declared double, they may be converted to single
715        precision float before running. If your algorithm depends on
716        precision (which is not uncommon for numerical algorithms), use
717        the following::
718
719            #if FLOAT_SIZE>4
720            ... code for double precision ...
721            #else
722            ... code for single precision ...
723            #endif
724    SAS_DOUBLE:
725        A replacement for :code:`double` so that the declared variable will
726        stay double precision; this should generally not be used since some
727        graphics cards do not support double precision.  There is no provision
728        for forcing a constant to stay double precision.
729
730The following special functions and scattering calculations are defined in
731`sasmodels/models/lib <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib>`_.
732These functions have been tuned to be fast and numerically stable down
733to $q=0$ even in single precision.  In some cases they work around bugs
734which appear on some platforms but not others, so use them where needed.
735Add the files listed in :code:`source = ["lib/file.c", ...]` to your *model.py*
736file in the order given, otherwise these functions will not be available.
737
738    polevl(x, c, n):
739        Polynomial evaluation $p(x) = \sum_{i=0}^n c_i x^i$ using Horner's
740        method so it is faster and more accurate.
741
742        $c = \{c_n, c_{n-1}, \ldots, c_0 \}$ is the table of coefficients,
743        sorted from highest to lowest.
744
745        :code:`source = ["lib/polevl.c", ...]` (`link to code <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/polevl.c>`_)
746
747    p1evl(x, c, n):
748        Evaluation of normalized polynomial $p(x) = x^n + \sum_{i=0}^{n-1} c_i x^i$
749        using Horner's method so it is faster and more accurate.
750
751        $c = \{c_{n-1}, c_{n-2} \ldots, c_0 \}$ is the table of coefficients,
752        sorted from highest to lowest.
753
754        :code:`source = ["lib/polevl.c", ...]`
[870a2f4]755        (`polevl.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/polevl.c>`_)
[990d8df]756
757    sas_gamma(x):
[30b60d2]758        Gamma function sas_gamma\ $(x) = \Gamma(x)$.
[990d8df]759
760        The standard math function, tgamma(x) is unstable for $x < 1$
761        on some platforms.
762
[870a2f4]763        :code:`source = ["lib/sas_gamma.c", ...]`
764        (`sas_gamma.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/sas_gamma.c>`_)
[990d8df]765
766    sas_erf(x), sas_erfc(x):
767        Error function
[30b60d2]768        sas_erf\ $(x) = \frac{2}{\sqrt\pi}\int_0^x e^{-t^2}\,dt$
[990d8df]769        and complementary error function
[30b60d2]770        sas_erfc\ $(x) = \frac{2}{\sqrt\pi}\int_x^{\infty} e^{-t^2}\,dt$.
[990d8df]771
772        The standard math functions erf(x) and erfc(x) are slower and broken
773        on some platforms.
774
775        :code:`source = ["lib/polevl.c", "lib/sas_erf.c", ...]`
[870a2f4]776        (`sas_erf.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/sas_erf.c>`_)
[990d8df]777
778    sas_J0(x):
[30b60d2]779        Bessel function of the first kind sas_J0\ $(x)=J_0(x)$ where
[990d8df]780        $J_0(x) = \frac{1}{\pi}\int_0^\pi \cos(x\sin(\tau))\,d\tau$.
781
782        The standard math function j0(x) is not available on all platforms.
783
784        :code:`source = ["lib/polevl.c", "lib/sas_J0.c", ...]`
[870a2f4]785        (`sas_J0.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/sas_J0.c>`_)
[990d8df]786
787    sas_J1(x):
[30b60d2]788        Bessel function of the first kind  sas_J1\ $(x)=J_1(x)$ where
[990d8df]789        $J_1(x) = \frac{1}{\pi}\int_0^\pi \cos(\tau - x\sin(\tau))\,d\tau$.
790
791        The standard math function j1(x) is not available on all platforms.
792
793        :code:`source = ["lib/polevl.c", "lib/sas_J1.c", ...]`
[870a2f4]794        (`sas_J1.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/sas_J1.c>`_)
[990d8df]795
796    sas_JN(n, x):
[30b60d2]797        Bessel function of the first kind and integer order $n$,
798        sas_JN\ $(n, x) =J_n(x)$ where
[990d8df]799        $J_n(x) = \frac{1}{\pi}\int_0^\pi \cos(n\tau - x\sin(\tau))\,d\tau$.
[30b60d2]800        If $n$ = 0 or 1, it uses sas_J0($x$) or sas_J1($x$), respectively.
[990d8df]801
802        The standard math function jn(n, x) is not available on all platforms.
803
804        :code:`source = ["lib/polevl.c", "lib/sas_J0.c", "lib/sas_J1.c", "lib/sas_JN.c", ...]`
[870a2f4]805        (`sas_JN.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/sas_JN.c>`_)
[990d8df]806
807    sas_Si(x):
[30b60d2]808        Sine integral Si\ $(x) = \int_0^x \tfrac{\sin t}{t}\,dt$.
[990d8df]809
810        This function uses Taylor series for small and large arguments:
811
812        For large arguments,
813
814        .. math::
815
816             \text{Si}(x) \sim \frac{\pi}{2}
817             - \frac{\cos(x)}{x}\left(1 - \frac{2!}{x^2} + \frac{4!}{x^4} - \frac{6!}{x^6} \right)
818             - \frac{\sin(x)}{x}\left(\frac{1}{x} - \frac{3!}{x^3} + \frac{5!}{x^5} - \frac{7!}{x^7}\right)
819
820        For small arguments,
821
822        .. math::
823
824           \text{Si}(x) \sim x
825           - \frac{x^3}{3\times 3!} + \frac{x^5}{5 \times 5!} - \frac{x^7}{7 \times 7!}
826           + \frac{x^9}{9\times 9!} - \frac{x^{11}}{11\times 11!}
827
828        :code:`source = ["lib/Si.c", ...]`
[870a2f4]829        (`Si.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/Si.c>`_)
[990d8df]830
831    sas_3j1x_x(x):
832        Spherical Bessel form
[30b60d2]833        sph_j1c\ $(x) = 3 j_1(x)/x = 3 (\sin(x) - x \cos(x))/x^3$,
[990d8df]834        with a limiting value of 1 at $x=0$, where $j_1(x)$ is the spherical
835        Bessel function of the first kind and first order.
836
837        This function uses a Taylor series for small $x$ for numerical accuracy.
838
839        :code:`source = ["lib/sas_3j1x_x.c", ...]`
[870a2f4]840        (`sas_3j1x_x.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/sas_3j1x_x.c>`_)
[990d8df]841
842
843    sas_2J1x_x(x):
[30b60d2]844        Bessel form sas_J1c\ $(x) = 2 J_1(x)/x$, with a limiting value
[990d8df]845        of 1 at $x=0$, where $J_1(x)$ is the Bessel function of first kind
846        and first order.
847
848        :code:`source = ["lib/polevl.c", "lib/sas_J1.c", ...]`
[870a2f4]849        (`sas_J1.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/sas_J1.c>`_)
[990d8df]850
851
852    Gauss76Z[i], Gauss76Wt[i]:
853        Points $z_i$ and weights $w_i$ for 76-point Gaussian quadrature, respectively,
854        computing $\int_{-1}^1 f(z)\,dz \approx \sum_{i=1}^{76} w_i\,f(z_i)$.
855
856        Similar arrays are available in :code:`gauss20.c` for 20-point
857        quadrature and in :code:`gauss150.c` for 150-point quadrature.
[d0dc9a3]858        The macros :code:`GAUSS_N`, :code:`GAUSS_Z` and :code:`GAUSS_W` are
859        defined so that you can change the order of the integration by
860        selecting an different source without touching the C code.
[990d8df]861
862        :code:`source = ["lib/gauss76.c", ...]`
[870a2f4]863        (`gauss76.c <https://github.com/SasView/sasmodels/tree/master/sasmodels/models/lib/gauss76.c>`_)
[990d8df]864
865
866
867Problems with C models
868......................
869
870The graphics processor (GPU) in your computer is a specialized computer tuned
871for certain kinds of problems.  This leads to strange restrictions that you
872need to be aware of.  Your code may work fine on some platforms or for some
873models, but then return bad values on other platforms.  Some examples of
874particular problems:
875
876  **(1) Code is too complex, or uses too much memory.** GPU devices only
877  have a limited amount of memory available for each processor. If you run
878  programs which take too much memory, then rather than running multiple
879  values in parallel as it usually does, the GPU may only run a single
880  version of the code at a time, making it slower than running on the CPU.
881  It may fail to run on some platforms, or worse, cause the screen to go
882  blank or the system to reboot.
883
884  **(2) Code takes too long.** Because GPU devices are used for the computer
885  display, the OpenCL drivers are very careful about the amount of time they
886  will allow any code to run. For example, on OS X, the model will stop
887  running after 5 seconds regardless of whether the computation is complete.
888  You may end up with only some of your 2D array defined, with the rest
889  containing random data. Or it may cause the screen to go blank or the
890  system to reboot.
891
892  **(3) Memory is not aligned**. The GPU hardware is specialized to operate
893  on multiple values simultaneously. To keep the GPU simple the values in
894  memory must be aligned with the different GPU compute engines. Not
895  following these rules can lead to unexpected values being loaded into
896  memory, and wrong answers computed. The conclusion from a very long and
897  strange debugging session was that any arrays that you declare in your
898  model should be a multiple of four. For example::
899
900      double Iq(q, p1, p2, ...)
901      {
902          double vector[8];  // Only going to use seven slots, but declare 8
903          ...
904      }
905
906The first step when your model is behaving strangely is to set
907**single=False**. This automatically restricts the model to only run on the
908CPU, or on high-end GPU cards. There can still be problems even on high-end
909cards, so you can force the model off the GPU by setting **opencl=False**.
910This runs the model as a normal C program without any GPU restrictions so
911you know that strange results are probably from your code rather than the
912environment. Once the code is debugged, you can compare your output to the
913output on the GPU.
914
915Although it can be difficult to get your model to work on the GPU, the reward
916can be a model that runs 1000x faster on a good card.  Even your laptop may
917show a 50x improvement or more over the equivalent pure python model.
918
919
920.. _Form_Factors:
921
922Form Factors
923............
924
925Away from the dilute limit you can estimate scattering including
926particle-particle interactions using $I(q) = P(q)*S(q)$ where $P(q)$
927is the form factor and $S(q)$ is the structure factor.  The simplest
928structure factor is the *hardsphere* interaction, which
929uses the effective radius of the form factor as an input to the structure
930factor model.  The effective radius is the average radius of the
931form averaged over all the polydispersity values.
932
933::
934
935    def ER(radius, thickness):
936        """Effective radius of a core-shell sphere."""
937        return radius + thickness
938
939Now consider the *core_shell_sphere*, which has a simple effective radius
940equal to the radius of the core plus the thickness of the shell, as
941shown above. Given polydispersity over *(r1, r2, ..., rm)* in radius and
942*(t1, t2, ..., tn)* in thickness, *ER* is called with a mesh
943grid covering all possible combinations of radius and thickness.
944That is, *radius* is *(r1, r2, ..., rm, r1, r2, ..., rm, ...)*
945and *thickness* is *(t1, t1, ... t1, t2, t2, ..., t2, ...)*.
946The *ER* function returns one effective radius for each combination.
947The effective radius calculator weights each of these according to
948the polydispersity distributions and calls the structure factor
949with the average *ER*.
950
951::
952
953    def VR(radius, thickness):
954        """Sphere and shell volumes for a core-shell sphere."""
955        whole = 4.0/3.0 * pi * (radius + thickness)**3
956        core = 4.0/3.0 * pi * radius**3
957        return whole, whole - core
958
959Core-shell type models have an additional volume ratio which scales
960the structure factor.  The *VR* function returns the volume of
961the whole sphere and the volume of the shell. Like *ER*, there is
962one return value for each point in the mesh grid.
963
964*NOTE: we may be removing or modifying this feature soon. As of the
965time of writing, core-shell sphere returns (1., 1.) for VR, giving a volume
966ratio of 1.0.*
967
968Unit Tests
969..........
970
971THESE ARE VERY IMPORTANT. Include at least one test for each model and
972PLEASE make sure that the answer value is correct (i.e. not a random number).
973
974::
975
976    tests = [
977        [{}, 0.2, 0.726362],
978        [{"scale": 1., "background": 0., "sld": 6., "sld_solvent": 1.,
979          "radius": 120., "radius_pd": 0.2, "radius_pd_n":45},
980         0.2, 0.228843],
981        [{"radius": 120., "radius_pd": 0.2, "radius_pd_n":45}, "ER", 120.],
982        [{"radius": 120., "radius_pd": 0.2, "radius_pd_n":45}, "VR", 1.],
983    ]
984
985
986**tests=[[{parameters}, q, result], ...]** is a list of lists.
987Each list is one test and contains, in order:
988
989- a dictionary of parameter values. This can be *{}* using the default
990  parameters, or filled with some parameters that will be different from the
991  default, such as *{"radius":10.0, "sld":4}*. Unlisted parameters will
992  be given the default values.
993- the input $q$ value or tuple of $(q_x, q_y)$ values.
994- the output $I(q)$ or $I(q_x,q_y)$ expected of the model for the parameters
995  and input value given.
996- input and output values can themselves be lists if you have several
997  $q$ values to test for the same model parameters.
998- for testing *ER* and *VR*, give the inputs as "ER" and "VR" respectively;
999  the output for *VR* should be the sphere/shell ratio, not the individual
1000  sphere and shell values.
1001
1002.. _Test_Your_New_Model:
1003
1004Test Your New Model
1005^^^^^^^^^^^^^^^^^^^
1006
1007Minimal Testing
1008...............
1009
1010From SasView either open the Python shell (*Tools* > *Python Shell/Editor*)
1011or the plugin editor (*Fitting* > *Plugin Model Operations* > *Advanced
1012Plugin Editor*), load your model, and then select *Run > Check Model* from
1013the menu bar. An *Info* box will appear with the results of the compilation
1014and a check that the model runs.
1015
1016If you are not using sasmodels from SasView, skip this step.
1017
1018Recommended Testing
1019...................
1020
1021If the model compiles and runs, you can next run the unit tests that
1022you have added using the **test =** values.
1023
1024From SasView, switch to the *Shell* tab and type the following::
1025
1026    from sasmodels.model_test import run_one
1027    run_one("~/.sasview/plugin_models/model.py")
1028
1029This should print::
1030
1031    test_model_python (sasmodels.model_test.ModelTestCase) ... ok
1032
1033To check whether single precision is good enough, type the following::
1034
1035    from sasmodels.compare import main as compare
1036    compare("~/.sasview/plugin_models/model.py")
1037
1038This will pop up a plot showing the difference between single precision
1039and double precision on a range of $q$ values.
1040
1041::
1042
1043  demo = dict(scale=1, background=0,
1044              sld=6, sld_solvent=1,
1045              radius=120,
1046              radius_pd=.2, radius_pd_n=45)
1047
1048**demo={'par': value, ...}** in the model file sets the default values for
1049the comparison. You can include polydispersity parameters such as
1050*radius_pd=0.2, radius_pd_n=45* which would otherwise be zero.
1051
1052These commands can also be run directly in the python interpreter:
1053
1054    $ python -m sasmodels.model_test -v ~/.sasview/plugin_models/model.py
1055    $ python -m sasmodels.compare ~/.sasview/plugin_models/model.py
1056
1057The options to compare are quite extensive; type the following for help::
1058
1059    compare()
1060
1061Options will need to be passed as separate strings.
1062For example to run your model with a random set of parameters::
1063
1064    compare("-random", "-pars", "~/.sasview/plugin_models/model.py")
1065
1066For the random models,
1067
1068- *sld* will be in the range (-0.5,10.5),
1069- angles (*theta, phi, psi*) will be in the range (-180,180),
1070- angular dispersion will be in the range (0,45),
1071- polydispersity will be in the range (0,1)
1072- other values will be in the range (0, 2\ *v*), where *v* is the value
1073  of the parameter in demo.
1074
1075Dispersion parameters *n*\, *sigma* and *type* will be unchanged from
1076demo so that run times are more predictable (polydispersity calculated
1077across multiple parameters can be very slow).
1078
[3048ec6]1079If your model has 2D orientation calculation, then you should also
[990d8df]1080test with::
1081
1082    compare("-2d", "~/.sasview/plugin_models/model.py")
1083
1084Check The Docs
1085^^^^^^^^^^^^^^
1086
1087You can get a rough idea of how the documentation will look using the
1088following::
1089
1090    compare("-help", "~/.sasview/plugin_models/model.py")
1091
1092This does not use the same styling as the rest of the docs, but it will
1093allow you to check that your ReStructuredText and LaTeX formatting.
1094Here are some tools to help with the inevitable syntax errors:
1095
1096- `Sphinx cheat sheet <http://matplotlib.org/sampledoc/cheatsheet.html>`_
1097- `Sphinx Documentation <http://www.sphinx-doc.org/en/stable/>`_
1098- `MathJax <http://www.mathjax.org/>`_
1099- `amsmath <http://www.ams.org/publications/authors/tex/amslatex>`_
1100
1101There is also a neat online WYSIWYG ReStructuredText editor at
1102http://rst.ninjs.org\ .
1103
1104
1105Clean Lint - (Developer Version Only)
1106^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
1107
1108**NB: For now we are not providing pylint with the installer version
1109of SasView; so unless you have a SasView build environment available,
1110you can ignore this section!**
1111
1112Run the lint check with::
1113
1114    python -m pylint --rcfile=extra/pylint.rc ~/.sasview/plugin_models/model.py
1115
1116We are not aiming for zero lint just yet, only keeping it to a minimum.
1117For now, don't worry too much about *invalid-name*. If you really want a
1118variable name *Rg* for example because $R_g$ is the right name for the model
1119parameter then ignore the lint errors.  Also, ignore *missing-docstring*
[108e70e]1120for standard model functions *Iq*, *Iqac*, etc.
[990d8df]1121
1122We will have delinting sessions at the SasView Code Camps, where we can
1123decide on standards for model files, parameter names, etc.
1124
1125For now, you can tell pylint to ignore things.  For example, to align your
1126parameters in blocks::
1127
1128    # pylint: disable=bad-whitespace,line-too-long
1129    #   ["name",                  "units", default, [lower, upper], "type", "description"],
1130    parameters = [
1131        ["contrast_factor",       "barns",    10.0,  [-inf, inf], "", "Contrast factor of the polymer"],
1132        ["bjerrum_length",        "Ang",       7.1,  [0, inf],    "", "Bjerrum length"],
1133        ["virial_param",          "1/Ang^2",  12.0,  [-inf, inf], "", "Virial parameter"],
1134        ["monomer_length",        "Ang",      10.0,  [0, inf],    "", "Monomer length"],
1135        ["salt_concentration",    "mol/L",     0.0,  [-inf, inf], "", "Concentration of monovalent salt"],
1136        ["ionization_degree",     "",          0.05, [0, inf],    "", "Degree of ionization"],
1137        ["polymer_concentration", "mol/L",     0.7,  [0, inf],    "", "Polymer molar concentration"],
1138        ]
1139    # pylint: enable=bad-whitespace,line-too-long
1140
1141Don't put in too many pylint statements, though, since they make the code ugly.
1142
1143Share Your Model!
1144^^^^^^^^^^^^^^^^^
1145
1146Once compare and the unit test(s) pass properly and everything is done,
1147consider adding your model to the
1148`Model Marketplace <http://marketplace.sasview.org/>`_ so that others may use it!
1149
1150.. ZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZZ
1151
1152*Document History*
1153
1154| 2016-10-25 Steve King
[c654160]1155| 2017-05-07 Paul Kienzle - Moved from sasview to sasmodels docs
Note: See TracBrowser for help on using the repository browser.