add some docs to explore/check_adaptive.py

explore/check_adaptive.py

@@ -1,3 +1,75 @@
+"""
+Check the speed and accuracy of the 1-D model integration routines,
+comparing them to a fixed 76-point gaussian.
+
+You can look at the overall performance by fitting the latex_sphere
+sasdata example with a triaxial ellipsoid. This can be done in sasview
+for a real world check, but it includes the overhead of the sasview shim.
+You can skip over the shim by running a fit directly in bumps:
+
+    $ bumps example/simul_fit.py --fit=dream --pop=-100
+
+Modify simul_fit.py to choose SANS, USANS or both for
+triaxial_ellipsoid, ellipsoid or sphere.
+
+Sasmodels compare help:
+
+    $ python -m sasmodels.compare -?
+
+Environment variables are listed in the sasmodels.compare help. For example,
+to control use of the GPU in the execution engine:
+
+    SAS_OPENCL=vendor:device|cuda:device|none sets the target GPU device
+
+Check speed and accuracy via explore/check_adaptive.py (this file).
+
+There are multiple control parameters for running this program, but you need
+to change the code to set them. Search for "runtime control flags" below.
+
+speed_target=2 warns if adaptive is 2x slower than gauss-76.
+
+speed_only only compares the calculation speed, not the accuracy
+
+speed_check is what to compare adaptive against: "gauss-76" fixed grid
+or adaptive running on the "gpu", or None for no speed check
+
+In the code UNNESTED are 1D integrals and NESTED are 2d integrals.
+
+big_n is the grid size for the target value. Time scales as n² for NESTED models,
+so only a few q values are tested.
+
+I've checked a couple of models using big_n=15000 by listing them on the command line:
+
+    $ python explore/check_adaptive.py capped_cylinder barbell
+
+Target sizes are 20 μm diameter for big (USANS) and 20 nm diameter for small (SANS/SAXS).
+Small models are tested over the SANS q range, 1e-3/Ang to 1/Ang. USANS only tests to 0.1/Ang
+
+For each model the generic a,b,c dimensions are translated to model specific parameters.
+
+There are other model parameters such as "sld_core" for core-shell models, which is set to
+sld_solvent=0, and the thickness/core_size ratio "shell" which defaults to 0.1. Maybe a
+future version could take these parameters on the command line. The latter is relevant
+to USANS since large objects with thin shells have a 1/q rather than 1/q^4 falloff, and
+so slit resolution values will be higher.
+
+Three different aspect ratios are tested: rod-like, disk-like and cube-like
+
+You can check random parameter sets for a model using:
+
+    $ python -m sasmodels.compare background=0 -ngauss=0,1000 -sets=5 -double triaxial_ellipsoid
+
+This compares adaptive (ngauss=0) to a 1000x1000 fixed grid.
+
+If you find a dataset that is anomalous then look for the -random=n output on the console
+and add it to the command line. He will reproduce a model with the same parameters.
+
+You may want something like "-random=n -ngauss=0,12000 -nq=5 q=0.001:0.1" to look at
+five q values over two decades with 12000x12000 grid
+
+Really small models are faster on the CPU. Large models are faster on the GPU.
+"""
+
 from math import sqrt
 
 from sasmodels import core, data, generate
@@ -125,7 +197,7 @@ def capped_cylinder(a, b, c):
     return pars
 
 @register(2)
-def core_shell_bicelle_elliptical_belt_rough(a, b, c, shell=0.1, sld_core=0, roughness=0.1):
+def core_shell_bicelle_elliptical_belt_rough(a, b, c, shell=0.1, sld_core=0, roughness=0.02):
     thick_rim = thick_face = a*shell/2 if shell < 1 else shell
     length = c - 2*thick_face
     radius = a/2 - thick_rim
@@ -322,6 +394,10 @@ def print_label():
 
         k_accurate = get_kernel(model_name, n_gauss=test_n_gauss, dtype="double", platform="dll")
         Iq_test_accurate = DirectModel(test_data, k_accurate)(**pars)
+        if not np.isfinite(Iq_test_accurate).any() or (Iq_test_accurate == 0).any():
+            print(f">>> Bad target values for\n    {model_name} {par_str}\n    q   : {test_q}\n    I(q): {Iq_test_accurate}")
+        if not np.isfinite(Iq_test_adaptive).any() or (Iq_test_adaptive == 0).any():
+            print(f">>> Bad target values for\n    {model_name} {par_str}\n    q   : {test_q}\n    I(q): {Iq_test_adaptive}")
         relerr = abs(Iq_test_adaptive - Iq_test_accurate)/Iq_test_accurate
         index = relerr > tol
         if index.any():
@@ -334,6 +410,8 @@ def print_label():
         # Compare against 76-point gaussian
         k_76 = get_kernel(model_name, n_gauss=76, dtype="double", platform="dll")
         Iq_test_76 = DirectModel(test_data, k_76)(**pars)
+        if not np.isfinite(Iq_test_76).any() or (Iq_test_76 == 0).any():
+            print(f">>> Bad target values for\n    {model_name} {par_str}\n    q   : {test_q}\n    I(q): {Iq_test_76}")
         relerr_76 = abs(Iq_test_76 - Iq_test_accurate)/Iq_test_accurate
         index = (relerr_76 > tol) #& (relerr > tol) & (relerr > 10*relerr_76)
         if index.any():
@@ -344,11 +422,13 @@ def print_label():
             print("  relerr", relerr[index])
             print("  rel-76", relerr_76[index])
 
+
+    # === runtime control flags ===
     speed_target = 2
-    #speed_only = True
-    #speed_check = None
     speed_only = False
     speed_check = "gauss-76"  # None | "" | "gpu" | "gauss-76"
+    #speed_only = True
+    #speed_check = None
     big_n = 5000
     #small_n = 1000
     small_n = big_n