Changes in / [a3412a6:9150036] in sasmodels


Ignore:
Files:
8 edited

Legend:

Unmodified
Added
Removed

  • explore/realspace.py

    r5778279 r362a66f  
    99import numpy as np 
    1010from numpy import pi, radians, sin, cos, sqrt 
    11 from numpy.random import poisson, uniform, randn, rand, randint 
     11from numpy.random import poisson, uniform, randn, rand 
    1212from numpy.polynomial.legendre import leggauss 
    1313from scipy.integrate import simps 
    … …  
    7878 
    7979 
    80 I3 = np.matrix([[1., 0, 0], [0, 1, 0], [0, 0, 1]]) 
    81  
    8280class Shape: 
    83     rotation = I3 
     81    rotation = np.matrix([[1., 0, 0], [0, 1, 0], [0, 0, 1]]) 
    8482    center = np.array([0., 0., 0.])[:, None] 
    8583    r_max = None 
    86     lattice_size = np.array((1, 1, 1)) 
    87     lattice_spacing = np.array((1., 1., 1.)) 
    88     lattice_distortion = 0.0 
    89     lattice_rotation = 0.0 
    90     lattice_type = "" 
    9184 
    9285    def volume(self): 
    … …  
    10396 
    10497    def rotate(self, theta, phi, psi): 
    105         if theta != 0. or phi != 0. or psi != 0.: 
    106             self.rotation = rotation(theta, phi, psi) * self.rotation 
     98        self.rotation = rotation(theta, phi, psi) * self.rotation 
    10799        return self 
    108100 
    … …  
    111103        return self 
    112104 
    113     def lattice(self, size=(1, 1, 1), spacing=(2, 2, 2), type="sc", 
    114                 distortion=0.0, rotation=0.0): 
    115         self.lattice_size = np.asarray(size, 'i') 
    116         self.lattice_spacing = np.asarray(spacing, 'd') 
    117         self.lattice_type = type 
    118         self.lattice_distortion = distortion 
    119         self.lattice_rotation = rotation 
    120  
    121105    def _adjust(self, points): 
    122         if self.rotation is I3: 
    123             points = points.T + self.center 
    124         else: 
    125             points = np.asarray(self.rotation * np.matrix(points.T)) + self.center 
    126         if self.lattice_type: 
    127             points = self._apply_lattice(points) 
     106        points = np.asarray(self.rotation * np.matrix(points.T)) + self.center 
    128107        return points.T 
    129108 
    130     def r_bins(self, q, over_sampling=10, r_step=0.): 
    131         if self.lattice_type: 
    132             r_max = np.sqrt(np.sum(self.lattice_size*self.lattice_spacing*self.dims)**2)/2 
    133         else: 
    134             r_max = self.r_max 
    135         #r_max = min(2 * pi / q[0], r_max) 
     109    def r_bins(self, q, over_sampling=1, r_step=0.): 
     110        r_max = min(2 * pi / q[0], self.r_max) 
    136111        if r_step == 0.: 
    137112            r_step = 2 * pi / q[-1] / over_sampling 
    138113        #r_step = 0.01 
    139114        return np.arange(r_step, r_max, r_step) 
    140  
    141     def _apply_lattice(self, points): 
    142         """Spread points to different lattice positions""" 
    143         size = self.lattice_size 
    144         spacing = self.lattice_spacing 
    145         shuffle = self.lattice_distortion 
    146         rotate = self.lattice_rotation 
    147         lattice = self.lattice_type 
    148  
    149         if rotate != 0: 
    150             # To vectorize the rotations we will need to unwrap the matrix multiply 
    151             raise NotImplementedError("don't handle rotations yet") 
    152  
    153         # Determine the number of lattice points in the lattice 
    154         shapes_per_cell = 2 if lattice == "bcc" else 4 if lattice == "fcc" else 1 
    155         number_of_lattice_points = np.prod(size) * shapes_per_cell 
    156  
    157         # For each point in the original shape, figure out which lattice point 
    158         # to translate it to.  This is both cell index (i*ny*nz + j*nz  + k) as 
    159         # well as the point in the cell (corner, body center or face center). 
    160         nsamples = points.shape[1] 
    161         lattice_point = randint(number_of_lattice_points, size=nsamples) 
    162  
    163         # Translate the cell index into the i,j,k coordinates of the senter 
    164         cell_index = lattice_point // shapes_per_cell 
    165         center = np.vstack((cell_index//(size[1]*size[2]), 
    166                             (cell_index%(size[1]*size[2]))//size[2], 
    167                             cell_index%size[2])) 
    168         center = np.asarray(center, dtype='d') 
    169         if lattice == "bcc": 
    170             center[:, lattice_point % shapes_per_cell == 1] += [[0.5], [0.5], [0.5]] 
    171         elif lattice == "fcc": 
    172             center[:, lattice_point % shapes_per_cell == 1] += [[0.0], [0.5], [0.5]] 
    173             center[:, lattice_point % shapes_per_cell == 2] += [[0.5], [0.0], [0.5]] 
    174             center[:, lattice_point % shapes_per_cell == 3] += [[0.5], [0.5], [0.0]] 
    175  
    176         # Each lattice point has its own displacement from the ideal position. 
    177         # Not checking that shapes do not overlap if displacement is too large. 
    178         offset = shuffle*(randn(3, number_of_lattice_points) if shuffle < 0.3 
    179                           else rand(3, number_of_lattice_points)) 
    180         center += offset[:, cell_index] 
    181  
    182         # Each lattice point has its own rotation.  Rotate the point prior to 
    183         # applying any displacement. 
    184         # rotation = rotate*(randn(size=(shapes, 3)) if shuffle < 30 else rand(size=(nsamples, 3))) 
    185         # for k in shapes: points[k] = rotation[k]*points[k] 
    186         points += center*(np.array([spacing])*np.array(self.dims)).T 
    187         return points 
    188115 
    189116class Composite(Shape): 
    … …  
    742669    Iq = 100 * np.ones_like(qx) 
    743670    data = Data2D(x=qx, y=qy, z=Iq, dx=None, dy=None, dz=np.sqrt(Iq)) 
    744     data.x_bins = qx[0, :] 
    745     data.y_bins = qy[:, 0] 
     671    data.x_bins = qx[0,:] 
     672    data.y_bins = qy[:,0] 
    746673    data.filename = "fake data" 
    747674 
    … …  
    768695    return shape, fn, fn_xy 
    769696 
    770 DEFAULT_SPHERE_RADIUS = 125 
    771 DEFAULT_SPHERE_CONTRAST = 2 
    772 def build_sphere(radius=DEFAULT_SPHERE_RADIUS, rho=DEFAULT_SPHERE_CONTRAST): 
     697def build_sphere(radius=125, rho=2): 
    773698    shape = TriaxialEllipsoid(radius, radius, radius, rho) 
    774699    fn = lambda q: sphere_Iq(q, radius)*rho**2 
    … …  
    826751    return shape, fn, fn_xy 
    827752 
    828 def build_sc_lattice(shape, nx=1, ny=1, nz=1, dx=2, dy=2, dz=2, 
    829                         shuffle=0, rotate=0): 
     753def build_cubic_lattice(shape, nx=1, ny=1, nz=1, dx=2, dy=2, dz=2, 
     754                  shuffle=0, rotate=0): 
    830755    a, b, c = shape.dims 
    831     corners= [copy(shape) 
     756    shapes = [copy(shape) 
    832757              .shift((ix+(randn() if shuffle < 0.3 else rand())*shuffle)*dx*a, 
    833758                     (iy+(randn() if shuffle < 0.3 else rand())*shuffle)*dy*b, 
    … …  
    837762              for iy in range(ny) 
    838763              for iz in range(nz)] 
    839     lattice = Composite(corners) 
     764    lattice = Composite(shapes) 
    840765    return lattice 
    841766 
    842 def build_bcc_lattice(shape, nx=1, ny=1, nz=1, dx=2, dy=2, dz=2, 
    843                       shuffle=0, rotate=0): 
    844     a, b, c = shape.dims 
    845     corners = [copy(shape) 
    846                .shift((ix+(randn() if shuffle < 0.3 else rand())*shuffle)*dx*a, 
    847                       (iy+(randn() if shuffle < 0.3 else rand())*shuffle)*dy*b, 
    848                       (iz+(randn() if shuffle < 0.3 else rand())*shuffle)*dz*c) 
    849                .rotate(*((randn(3) if rotate < 30 else rand(3))*rotate)) 
    850                for ix in range(nx) 
    851                for iy in range(ny) 
    852                for iz in range(nz)] 
    853     centers = [copy(shape) 
    854                .shift((ix+0.5+(randn() if shuffle < 0.3 else rand())*shuffle)*dx*a, 
    855                       (iy+0.5+(randn() if shuffle < 0.3 else rand())*shuffle)*dy*b, 
    856                       (iz+0.5+(randn() if shuffle < 0.3 else rand())*shuffle)*dz*c) 
    857                .rotate(*((randn(3) if rotate < 30 else rand(3))*rotate)) 
    858                for ix in range(nx) 
    859                for iy in range(ny) 
    860                for iz in range(nz)] 
    861     lattice = Composite(corners + centers) 
    862     return lattice 
    863  
    864 def build_fcc_lattice(shape, nx=1, ny=1, nz=1, dx=2, dy=2, dz=2, 
    865                       shuffle=0, rotate=0): 
    866     a, b, c = shape.dims 
    867     corners = [copy(shape) 
    868                .shift((ix+(randn() if shuffle < 0.3 else rand())*shuffle)*dx*a, 
    869                       (iy+(randn() if shuffle < 0.3 else rand())*shuffle)*dy*b, 
    870                       (iz+(randn() if shuffle < 0.3 else rand())*shuffle)*dz*c) 
    871                .rotate(*((randn(3) if rotate < 30 else rand(3))*rotate)) 
    872                for ix in range(nx) 
    873                for iy in range(ny) 
    874                for iz in range(nz)] 
    875     faces_a = [copy(shape) 
    876                .shift((ix+0.0+(randn() if shuffle < 0.3 else rand())*shuffle)*dx*a, 
    877                       (iy+0.5+(randn() if shuffle < 0.3 else rand())*shuffle)*dy*b, 
    878                       (iz+0.5+(randn() if shuffle < 0.3 else rand())*shuffle)*dz*c) 
    879                .rotate(*((randn(3) if rotate < 30 else rand(3))*rotate)) 
    880                for ix in range(nx) 
    881                for iy in range(ny) 
    882                for iz in range(nz)] 
    883     faces_b = [copy(shape) 
    884                .shift((ix+0.5+(randn() if shuffle < 0.3 else rand())*shuffle)*dx*a, 
    885                       (iy+0.0+(randn() if shuffle < 0.3 else rand())*shuffle)*dy*b, 
    886                       (iz+0.5+(randn() if shuffle < 0.3 else rand())*shuffle)*dz*c) 
    887                .rotate(*((randn(3) if rotate < 30 else rand(3))*rotate)) 
    888                for ix in range(nx) 
    889                for iy in range(ny) 
    890                for iz in range(nz)] 
    891     faces_c = [copy(shape) 
    892                .shift((ix+0.5+(randn() if shuffle < 0.3 else rand())*shuffle)*dx*a, 
    893                       (iy+0.5+(randn() if shuffle < 0.3 else rand())*shuffle)*dy*b, 
    894                       (iz+0.0+(randn() if shuffle < 0.3 else rand())*shuffle)*dz*c) 
    895                .rotate(*((randn(3) if rotate < 30 else rand(3))*rotate)) 
    896                for ix in range(nx) 
    897                for iy in range(ny) 
    898                for iz in range(nz)] 
    899     lattice = Composite(corners + faces_a + faces_b + faces_c) 
    900     return lattice 
    901767 
    902768SHAPE_FUNCTIONS = OrderedDict([ 
    … …  
    909775]) 
    910776SHAPES = list(SHAPE_FUNCTIONS.keys()) 
    911 LATTICE_FUNCTIONS = OrderedDict([ 
    912     ("sc", build_sc_lattice), 
    913     ("bcc", build_bcc_lattice), 
    914     ("fcc", build_fcc_lattice), 
    915 ]) 
    916 LATTICE_TYPES = list(LATTICE_FUNCTIONS.keys()) 
    917777 
    918778def check_shape(title, shape, fn=None, show_points=False, 
    … …  
    923783    r = shape.r_bins(q, r_step=r_step) 
    924784    sampling_density = samples / shape.volume 
    925     print("sampling points") 
    926785    rho, points = shape.sample(sampling_density) 
    927     print("calculating Pr") 
    928786    t0 = time.time() 
    929787    Pr = calc_Pr(r, rho-rho_solvent, points) 
    … …  
    934792    import pylab 
    935793    if show_points: 
    936         plot_points(rho, points); pylab.figure() 
     794         plot_points(rho, points); pylab.figure() 
    937795    plot_calc(r, Pr, q, Iq, theory=theory, title=title) 
    938796    pylab.gcf().canvas.set_window_title(title) 
    … …  
    948806    Qx, Qy = np.meshgrid(qx, qy) 
    949807    sampling_density = samples / shape.volume 
    950     print("sampling points") 
    951808    t0 = time.time() 
    952809    rho, points = shape.sample(sampling_density) 
    … …  
    987844                        help='lattice size') 
    988845    parser.add_argument('-z', '--spacing', type=str, default='2,2,2', 
    989                         help='lattice spacing (relative to shape)') 
    990     parser.add_argument('-t', '--type', choices=LATTICE_TYPES, 
    991                         default=LATTICE_TYPES[0], 
    992                         help='lattice type') 
     846                        help='lattice spacing') 
    993847    parser.add_argument('-r', '--rotate', type=float, default=0., 
    994848                        help="rotation relative to lattice, gaussian < 30 degrees, uniform otherwise") 
    … …  
    1004858    nx, ny, nz = [int(v) for v in opts.lattice.split(',')] 
    1005859    dx, dy, dz = [float(v) for v in opts.spacing.split(',')] 
    1006     distortion, rotation = opts.shuffle, opts.rotate 
     860    shuffle, rotate = opts.shuffle, opts.rotate 
    1007861    shape, fn, fn_xy = SHAPE_FUNCTIONS[opts.shape](**pars) 
    1008     view = tuple(float(v) for v in opts.view.split(',')) 
    1009     # If comparing a sphere in a cubic lattice, compare against the 
    1010     # corresponding paracrystalline model. 
    1011     if opts.shape == "sphere" and dx == dy == dz and nx*ny*nz > 1: 
    1012         radius = pars.get('radius', DEFAULT_SPHERE_RADIUS) 
    1013         model_name = opts.type + "_paracrystal" 
    1014         model_pars = { 
    1015             "scale": 1., 
    1016             "background": 0., 
    1017             "lattice_spacing": 2*radius*dx, 
    1018             "lattice_distortion": distortion, 
    1019             "radius": radius, 
    1020             "sld": pars.get('rho', DEFAULT_SPHERE_CONTRAST), 
    1021             "sld_solvent": 0., 
    1022             "theta": view[0], 
    1023             "phi": view[1], 
    1024             "psi": view[2], 
    1025         } 
    1026         fn, fn_xy = wrap_sasmodel(model_name, **model_pars) 
    1027     if nx*ny*nz > 1: 
    1028         if rotation != 0: 
    1029             print("building %s lattice"%opts.type) 
    1030             build_lattice = LATTICE_FUNCTIONS[opts.type] 
    1031             shape = build_lattice(shape, nx, ny, nz, dx, dy, dz, 
    1032                                   distortion, rotation) 
    1033         else: 
    1034             shape.lattice(size=(nx, ny, nz), spacing=(dx, dy, dz), 
    1035                           type=opts.type, 
    1036                           rotation=rotation, distortion=distortion) 
    1037  
     862    if nx > 1 or ny > 1 or nz > 1: 
     863        shape = build_cubic_lattice(shape, nx, ny, nz, dx, dy, dz, shuffle, rotate) 
    1038864    title = "%s(%s)" % (opts.shape, " ".join(opts.pars)) 
    1039865    if opts.dim == 1: 
    … …  
    1041867                    mesh=opts.mesh, qmax=opts.qmax, samples=opts.samples) 
    1042868    else: 
     869        view = tuple(float(v) for v in opts.view.split(',')) 
    1043870        check_shape_2d(title, shape, fn_xy, view=view, show_points=opts.plot, 
    1044871                       mesh=opts.mesh, qmax=opts.qmax, samples=opts.samples) 

  • sasmodels/convert.py

    r610ef23 r610ef23  
    166166        oldpars = _hand_convert_3_1_2_to_4_1(name, oldpars) 
    167167    if version < (4, 2, 0): 
    168         oldpars = _hand_convert_4_1_to_4_2(name, oldpars) 
    169168        oldpars = _rename_magnetic_pars(oldpars) 
    170     return oldpars 
    171  
    172 def _hand_convert_4_1_to_4_2(name, oldpars): 
    173     if name in ('bcc_paracrystal', 'fcc_paracrystal', 'sc_paracrystal'): 
    174         oldpars['lattice_spacing'] = oldpars.pop('dnn') 
    175         oldpars['lattice_distortion'] = oldpars.pop('d_factor') 
    176169    return oldpars 
    177170 

  • sasmodels/models/bcc_paracrystal.c

    r642046e r108e70e  
    11static double 
    2 bcc_Zq(double qa, double qb, double qc, double lattice_spacing, double lattice_distortion) 
     2bcc_Zq(double qa, double qb, double qc, double dnn, double d_factor) 
    33{ 
    44    // Equations from Matsuoka 26-27-28, multiplied by |q| 
    … …  
    1717    //         => exp(a)^2 - 2 cos(d ak) exp(a) + 1) 
    1818    //         => (exp(a) - 2 cos(d ak)) * exp(a) + 1 
    19     const double arg = -0.5*square(lattice_spacing*lattice_distortion)*(a1*a1 + a2*a2 + a3*a3); 
     19    const double arg = -0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
    2020    const double exp_arg = exp(arg); 
    2121    const double Zq = -cube(expm1(2.0*arg)) 
    22         / ( ((exp_arg - 2.0*cos(lattice_spacing*a1))*exp_arg + 1.0) 
    23           * ((exp_arg - 2.0*cos(lattice_spacing*a2))*exp_arg + 1.0) 
    24           * ((exp_arg - 2.0*cos(lattice_spacing*a3))*exp_arg + 1.0)); 
     22        / ( ((exp_arg - 2.0*cos(dnn*a1))*exp_arg + 1.0) 
     23          * ((exp_arg - 2.0*cos(dnn*a2))*exp_arg + 1.0) 
     24          * ((exp_arg - 2.0*cos(dnn*a3))*exp_arg + 1.0)); 
    2525 
    2626#elif 0 
    … …  
    3636    //            = tanh(-a) / [1 - cos(d a_k)/cosh(-a)] 
    3737    // 
    38     const double arg = 0.5*square(lattice_spacing*lattice_distortion)*(a1*a1 + a2*a2 + a3*a3); 
     38    const double arg = 0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
    3939    const double sinh_qd = sinh(arg); 
    4040    const double cosh_qd = cosh(arg); 
    41     const double Zq = sinh_qd/(cosh_qd - cos(lattice_spacing*a1)) 
    42                     * sinh_qd/(cosh_qd - cos(lattice_spacing*a2)) 
    43                     * sinh_qd/(cosh_qd - cos(lattice_spacing*a3)); 
     41    const double Zq = sinh_qd/(cosh_qd - cos(dnn*a1)) 
     42                    * sinh_qd/(cosh_qd - cos(dnn*a2)) 
     43                    * sinh_qd/(cosh_qd - cos(dnn*a3)); 
    4444#else 
    45     const double arg = 0.5*square(lattice_spacing*lattice_distortion)*(a1*a1 + a2*a2 + a3*a3); 
     45    const double arg = 0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
    4646    const double tanh_qd = tanh(arg); 
    4747    const double cosh_qd = cosh(arg); 
    48     const double Zq = tanh_qd/(1.0 - cos(lattice_spacing*a1)/cosh_qd) 
    49                     * tanh_qd/(1.0 - cos(lattice_spacing*a2)/cosh_qd) 
    50                     * tanh_qd/(1.0 - cos(lattice_spacing*a3)/cosh_qd); 
     48    const double Zq = tanh_qd/(1.0 - cos(dnn*a1)/cosh_qd) 
     49                    * tanh_qd/(1.0 - cos(dnn*a2)/cosh_qd) 
     50                    * tanh_qd/(1.0 - cos(dnn*a3)/cosh_qd); 
    5151#endif 
    5252 
    … …  
    5757// occupied volume fraction calculated from lattice symmetry and sphere radius 
    5858static double 
    59 bcc_volume_fraction(double radius, double lattice_spacing) 
     59bcc_volume_fraction(double radius, double dnn) 
    6060{ 
    61     return 2.0*sphere_volume(radius/lattice_spacing); 
     61    return 2.0*sphere_volume(sqrt(0.75)*radius/dnn); 
    6262} 
    6363 
    … …  
    6969 
    7070 
    71 static double Iq(double q, double lattice_spacing, 
    72     double lattice_distortion, double radius, 
     71static double Iq(double q, double dnn, 
     72    double d_factor, double radius, 
    7373    double sld, double solvent_sld) 
    7474{ 
    … …  
    9494            const double qa = qab*cos_phi; 
    9595            const double qb = qab*sin_phi; 
    96             const double form = bcc_Zq(qa, qb, qc, lattice_spacing, lattice_distortion); 
     96            const double form = bcc_Zq(qa, qb, qc, dnn, d_factor); 
    9797            inner_sum += GAUSS_W[j] * form; 
    9898        } 
    … …  
    103103    const double Zq = outer_sum/(4.0*M_PI); 
    104104    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    105     return bcc_volume_fraction(radius, lattice_spacing) * Pq * Zq; 
     105    return bcc_volume_fraction(radius, dnn) * Pq * Zq; 
    106106} 
    107107 
    108108 
    109109static double Iqabc(double qa, double qb, double qc, 
    110     double lattice_spacing, double lattice_distortion, double radius, 
     110    double dnn, double d_factor, double radius, 
    111111    double sld, double solvent_sld) 
    112112{ 
    113113    const double q = sqrt(qa*qa + qb*qb + qc*qc); 
    114     const double Zq = bcc_Zq(qa, qb, qc, lattice_spacing, lattice_distortion); 
     114    const double Zq = bcc_Zq(qa, qb, qc, dnn, d_factor); 
    115115    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    116     return bcc_volume_fraction(radius, lattice_spacing) * Pq * Zq; 
     116    return bcc_volume_fraction(radius, dnn) * Pq * Zq; 
    117117} 

  • sasmodels/models/bcc_paracrystal.py

    rda7b26b rda7b26b  
    3434.. math:: 
    3535 
    36     V_\text{lattice} = \frac{8\pi}{3} \frac{R^3}{\left(2D/\sqrt{3}\right)^3} 
     36    V_\text{lattice} = \frac{16\pi}{3} \frac{R^3}{\left(D\sqrt{2}\right)^3} 
    3737 
    3838 
    … …  
    104104 
    105105* **Author:** NIST IGOR/DANSE **Date:** pre 2010 
    106 * **Last Modified by:** Paul Butler **Date:** September 16, 2018 
    107 * **Last Reviewed by:** Paul Butler **Date:** September 16, 2018 
     106* **Last Modified by:** Paul Butler **Date:** September 29, 2016 
     107* **Last Reviewed by:** Richard Heenan **Date:** March 21, 2016 
    108108""" 
    109109 
    … …  
    127127# pylint: disable=bad-whitespace, line-too-long 
    128128#             ["name", "units", default, [lower, upper], "type","description" ], 
    129 parameters = [["lattice_spacing",         "Ang",       220,    [-inf, inf], "",            "Lattice spacing"], 
    130               ["lattice_distortion",    "",            0.06, [-inf, inf], "",            "Paracrystal distortion factor"], 
     129parameters = [["dnn",         "Ang",       220,    [-inf, inf], "",            "Nearest neighbour distance"], 
     130              ["d_factor",    "",            0.06, [-inf, inf], "",            "Paracrystal distortion factor"], 
    131131              ["radius",      "Ang",        40,    [0, inf],    "volume",      "Particle radius"], 
    132132              ["sld",         "1e-6/Ang^2",  4,    [-inf, inf], "sld",         "Particle scattering length density"], 
    … …  
    149149    # in this range 'cuz its easy. 
    150150    radius = 10**np.random.uniform(1.3, 4) 
    151     lattice_distortion = 10**np.random.uniform(-2, -0.7)  # sigma_d in 0.01-0.7 
    152     lattice_spacing_fraction = np.random.beta(a=10, b=1) 
    153     lattice_spacing = radius*4/np.sqrt(3)/lattice_spacing_fraction 
     151    d_factor = 10**np.random.uniform(-2, -0.7)  # sigma_d in 0.01-0.7 
     152    dnn_fraction = np.random.beta(a=10, b=1) 
     153    dnn = radius*4/np.sqrt(3)/dnn_fraction 
    154154    pars = dict( 
    155155        #sld=1, sld_solvent=0, scale=1, background=1e-32, 
    156         lattice_spacing=lattice_spacing, 
    157         lattice_distortion=lattice_distortion, 
     156        dnn=dnn, 
     157        d_factor=d_factor, 
    158158        radius=radius, 
    159159    ) 

  • sasmodels/models/fcc_paracrystal.c

    r642046e r71b751d  
    11static double 
    2 fcc_Zq(double qa, double qb, double qc, double lattice_spacing, double lattice_distortion) 
     2fcc_Zq(double qa, double qb, double qc, double dnn, double d_factor) 
    33{ 
    44    // Equations from Matsuoka 17-18-19, multiplied by |q| 
    … …  
    1616    //         => exp(a)^2 - 2 cos(d ak) exp(a) + 1) 
    1717    //         => (exp(a) - 2 cos(d ak)) * exp(a) + 1 
    18     const double arg = -0.5*square(lattice_spacing*lattice_distortion)*(a1*a1 + a2*a2 + a3*a3); 
     18    const double arg = -0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
    1919    const double exp_arg = exp(arg); 
    2020    const double Zq = -cube(expm1(2.0*arg)) 
    21         / ( ((exp_arg - 2.0*cos(lattice_spacing*a1))*exp_arg + 1.0) 
    22           * ((exp_arg - 2.0*cos(lattice_spacing*a2))*exp_arg + 1.0) 
    23           * ((exp_arg - 2.0*cos(lattice_spacing*a3))*exp_arg + 1.0)); 
     21        / ( ((exp_arg - 2.0*cos(dnn*a1))*exp_arg + 1.0) 
     22          * ((exp_arg - 2.0*cos(dnn*a2))*exp_arg + 1.0) 
     23          * ((exp_arg - 2.0*cos(dnn*a3))*exp_arg + 1.0)); 
    2424 
    2525    return Zq; 
    … …  
    2929// occupied volume fraction calculated from lattice symmetry and sphere radius 
    3030static double 
    31 fcc_volume_fraction(double radius, double lattice_spacing) 
     31fcc_volume_fraction(double radius, double dnn) 
    3232{ 
    33     return 4.0*sphere_volume(radius/lattice_spacing); 
     33    return 4.0*sphere_volume(M_SQRT1_2*radius/dnn); 
    3434} 
    3535 
    … …  
    4141 
    4242 
    43 static double Iq(double q, double lattice_spacing, 
    44   double lattice_distortion, double radius, 
     43static double Iq(double q, double dnn, 
     44  double d_factor, double radius, 
    4545  double sld, double solvent_sld) 
    4646{ 
    … …  
    6666            const double qa = qab*cos_phi; 
    6767            const double qb = qab*sin_phi; 
    68             const double form = fcc_Zq(qa, qb, qc, lattice_spacing, lattice_distortion); 
     68            const double form = fcc_Zq(qa, qb, qc, dnn, d_factor); 
    6969            inner_sum += GAUSS_W[j] * form; 
    7070        } 
    … …  
    7676    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    7777 
    78     return fcc_volume_fraction(radius, lattice_spacing) * Pq * Zq; 
     78    return fcc_volume_fraction(radius, dnn) * Pq * Zq; 
    7979} 
    8080 
    8181static double Iqabc(double qa, double qb, double qc, 
    82     double lattice_spacing, double lattice_distortion, double radius, 
     82    double dnn, double d_factor, double radius, 
    8383    double sld, double solvent_sld) 
    8484{ 
    8585    const double q = sqrt(qa*qa + qb*qb + qc*qc); 
    8686    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    87     const double Zq = fcc_Zq(qa, qb, qc, lattice_spacing, lattice_distortion); 
    88     return fcc_volume_fraction(radius, lattice_spacing) * Pq * Zq; 
     87    const double Zq = fcc_Zq(qa, qb, qc, dnn, d_factor); 
     88    return fcc_volume_fraction(radius, dnn) * Pq * Zq; 
    8989} 

  • sasmodels/models/fcc_paracrystal.py

    rda7b26b rda7b26b  
    33#note - calculation requires double precision 
    44r""" 
    5 .. warning:: This model and this model description are under review following 
    6              concerns raised by SasView users. If you need to use this model, 
    7              please email help@sasview.org for the latest situation. *The 
     5.. warning:: This model and this model description are under review following  
     6             concerns raised by SasView users. If you need to use this model,  
     7             please email help@sasview.org for the latest situation. *The  
    88             SasView Developers. September 2018.* 
    99 
    … …  
    100100 
    101101Authorship and Verification 
    102 ---------------------------- 
     102--------------------------- 
    103103 
    104104* **Author:** NIST IGOR/DANSE **Date:** pre 2010 
    105 * **Last Modified by:** Paul Butler **Date:** September 16, 2018 
    106 * **Last Reviewed by:** Paul Butler **Date:** September 16, 2018 
     105* **Last Modified by:** Paul Butler **Date:** September 29, 2016 
     106* **Last Reviewed by:** Richard Heenan **Date:** March 21, 2016 
    107107""" 
    108108 
    … …  
    123123# pylint: disable=bad-whitespace, line-too-long 
    124124#             ["name", "units", default, [lower, upper], "type","description"], 
    125 parameters = [["lattice_spacing", "Ang", 220, [-inf, inf], "", "Lattice spacing"], 
    126               ["lattice_distortion", "", 0.06, [-inf, inf], "", "Paracrystal distortion factor"], 
     125parameters = [["dnn", "Ang", 220, [-inf, inf], "", "Nearest neighbour distance"], 
     126              ["d_factor", "", 0.06, [-inf, inf], "", "Paracrystal distortion factor"], 
    127127              ["radius", "Ang", 40, [0, inf], "volume", "Particle radius"], 
    128128              ["sld", "1e-6/Ang^2", 4, [-inf, inf], "sld", "Particle scattering length density"], 
    … …  
    139139    # copied from bcc_paracrystal 
    140140    radius = 10**np.random.uniform(1.3, 4) 
    141     lattice_distortion = 10**np.random.uniform(-2, -0.7)  # sigma_d in 0.01-0.7 
    142     lattice_spacing_fraction = np.random.beta(a=10, b=1) 
    143     lattice_spacing = radius*4/np.sqrt(2)/lattice_spacing_fraction 
     141    d_factor = 10**np.random.uniform(-2, -0.7)  # sigma_d in 0.01-0.7 
     142    dnn_fraction = np.random.beta(a=10, b=1) 
     143    dnn = radius*4/np.sqrt(2)/dnn_fraction 
    144144    pars = dict( 
    145145        #sld=1, sld_solvent=0, scale=1, background=1e-32, 
    146         latice_spacing=lattice_spacing, 
    147         lattice_distortion=d_factor, 
     146        dnn=dnn, 
     147        d_factor=d_factor, 
    148148        radius=radius, 
    149149    ) 

  • sasmodels/models/sc_paracrystal.c

    r6530963 r71b751d  
    11static double 
    2 sc_Zq(double qa, double qb, double qc, double lattice_spacing, double lattice_distortion) 
     2sc_Zq(double qa, double qb, double qc, double dnn, double d_factor) 
    33{ 
    44    // Equations from Matsuoka 9-10-11, multiplied by |q| 
    … …  
    1616    //         => exp(a)^2 - 2 cos(d ak) exp(a) + 1) 
    1717    //         => (exp(a) - 2 cos(d ak)) * exp(a) + 1 
    18     const double arg = -0.5*square(lattice_spacing*lattice_distortion)*(a1*a1 + a2*a2 + a3*a3); 
     18    const double arg = -0.5*square(dnn*d_factor)*(a1*a1 + a2*a2 + a3*a3); 
    1919    const double exp_arg = exp(arg); 
    2020    const double Zq = -cube(expm1(2.0*arg)) 
    21         / ( ((exp_arg - 2.0*cos(lattice_spacing*a1))*exp_arg + 1.0) 
    22           * ((exp_arg - 2.0*cos(lattice_spacing*a2))*exp_arg + 1.0) 
    23           * ((exp_arg - 2.0*cos(lattice_spacing*a3))*exp_arg + 1.0)); 
     21        / ( ((exp_arg - 2.0*cos(dnn*a1))*exp_arg + 1.0) 
     22          * ((exp_arg - 2.0*cos(dnn*a2))*exp_arg + 1.0) 
     23          * ((exp_arg - 2.0*cos(dnn*a3))*exp_arg + 1.0)); 
    2424 
    2525    return Zq; 
    … …  
    2828// occupied volume fraction calculated from lattice symmetry and sphere radius 
    2929static double 
    30 sc_volume_fraction(double radius, double lattice_spacing) 
     30sc_volume_fraction(double radius, double dnn) 
    3131{ 
    32     return sphere_volume(radius/lattice_spacing); 
     32    return sphere_volume(radius/dnn); 
    3333} 
    3434 
    … …  
    4141 
    4242static double 
    43 Iq(double q, double lattice_spacing, 
    44     double lattice_distortion, double radius, 
     43Iq(double q, double dnn, 
     44    double d_factor, double radius, 
    4545    double sld, double solvent_sld) 
    4646{ 
    … …  
    6767            const double qa = qab*cos_phi; 
    6868            const double qb = qab*sin_phi; 
    69             const double form = sc_Zq(qa, qb, qc, lattice_spacing, lattice_distortion); 
     69            const double form = sc_Zq(qa, qb, qc, dnn, d_factor); 
    7070            inner_sum += GAUSS_W[j] * form; 
    7171        } 
    … …  
    7777    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    7878 
    79     return sc_volume_fraction(radius, lattice_spacing) * Pq * Zq; 
     79    return sc_volume_fraction(radius, dnn) * Pq * Zq; 
    8080} 
    8181 
    8282static double 
    8383Iqabc(double qa, double qb, double qc, 
    84     double lattice_spacing, double lattice_distortion, double radius, 
     84    double dnn, double d_factor, double radius, 
    8585    double sld, double solvent_sld) 
    8686{ 
    8787    const double q = sqrt(qa*qa + qb*qb + qc*qc); 
    8888    const double Pq = sphere_form(q, radius, sld, solvent_sld); 
    89     const double Zq = sc_Zq(qa, qb, qc, lattice_spacing, lattice_distortion); 
    90     return sc_volume_fraction(radius, lattice_spacing) * Pq * Zq; 
     89    const double Zq = sc_Zq(qa, qb, qc, dnn, d_factor); 
     90    return sc_volume_fraction(radius, dnn) * Pq * Zq; 
    9191} 

  • sasmodels/models/sc_paracrystal.py

    rda7b26b rda7b26b  
    11r""" 
    2 .. warning:: This model and this model description are under review following 
    3              concerns raised by SasView users. If you need to use this model, 
    4              please email help@sasview.org for the latest situation. *The 
     2.. warning:: This model and this model description are under review following  
     3             concerns raised by SasView users. If you need to use this model,  
     4             please email help@sasview.org for the latest situation. *The  
    55             SasView Developers. September 2018.* 
    6  
     6              
    77Definition 
    88---------- 
    … …  
    104104 
    105105Authorship and Verification 
    106 ---------------------------- 
     106--------------------------- 
    107107 
    108108* **Author:** NIST IGOR/DANSE **Date:** pre 2010 
    109 * **Last Modified by:** Paul Butler **Date:** September 16, 2018 
    110 * **Last Reviewed by:** Paul Butler **Date:** September 16, 2018 
     109* **Last Modified by:** Paul Butler **Date:** September 29, 2016 
     110* **Last Reviewed by:** Richard Heenan **Date:** March 21, 2016 
    111111""" 
    112112 
    … …  
    129129        scale: volume fraction of spheres 
    130130        bkg:background, R: radius of sphere 
    131         lattice_spacing: Nearest neighbor distance 
    132         lattice_distortion: Paracrystal distortion factor 
     131        dnn: Nearest neighbor distance 
     132        d_factor: Paracrystal distortion factor 
    133133        radius: radius of the spheres 
    134134        sldSph: SLD of the sphere 
    … …  
    139139# pylint: disable=bad-whitespace, line-too-long 
    140140#             ["name", "units", default, [lower, upper], "type","description"], 
    141 parameters = [["lattice_spacing",         "Ang",       220.0, [0.0, inf],  "",  "Lattice spacing"], 
    142               ["lattice_distortion",    "",           0.06, [-inf, inf], "",   "Paracrystal distortion factor"], 
     141parameters = [["dnn",         "Ang",       220.0, [0.0, inf],  "",            "Nearest neighbor distance"], 
     142              ["d_factor",    "",           0.06, [-inf, inf], "",            "Paracrystal distortion factor"], 
    143143              ["radius",      "Ang",        40.0, [0.0, inf],  "volume",      "Radius of sphere"], 
    144144              ["sld",  "1e-6/Ang^2",         3.0, [0.0, inf],  "sld",         "Sphere scattering length density"], 
    … …  
    155155    # copied from bcc_paracrystal 
    156156    radius = 10**np.random.uniform(1.3, 4) 
    157     lattice_distortion = 10**np.random.uniform(-2, -0.7)  # sigma_d in 0.01-0.7 
    158     lattice_spacing_fraction = np.random.beta(a=10, b=1) 
    159     lattice_spacing = radius*4/np.sqrt(4)/lattice_spacing_fraction 
     157    d_factor = 10**np.random.uniform(-2, -0.7)  # sigma_d in 0.01-0.7 
     158    dnn_fraction = np.random.beta(a=10, b=1) 
     159    dnn = radius*4/np.sqrt(4)/dnn_fraction 
    160160    pars = dict( 
    161161        #sld=1, sld_solvent=0, scale=1, background=1e-32, 
    162         lattice_spacing=lattice_spacing, 
    163         lattice_distortion=lattice_distortion, 
     162        dnn=dnn, 
     163        d_factor=d_factor, 
    164164        radius=radius, 
    165165    ) 
Note: See TracChangeset for help on using the changeset viewer.