Changeset b3f4831 in sasmodels for explore


Ignore:
Timestamp:
Oct 30, 2018 11:07:41 AM (6 years ago)
Author:
Paul Kienzle <pkienzle@…>
Branches:
ticket_1156
Children:
cc8b183
Parents:
5778279 (diff), c6084f1 (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.
Message:

Merge branch 'master' into ticket_1156

Location:
explore
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • explore/precision.py

    r2a7e20e rfba9ca0  
    9595            neg:    [-100,100] 
    9696 
     97        For arbitrary range use "start:stop:steps:scale" where scale is 
     98        one of log, lin, or linear. 
     99 
    97100        *diff* is "relative", "absolute" or "none" 
    98101 
    … …  
    102105        linear = not xrange.startswith("log") 
    103106        if xrange == "zoom": 
    104             lin_min, lin_max, lin_steps = 1000, 1010, 2000 
     107            start, stop, steps = 1000, 1010, 2000 
    105108        elif xrange == "neg": 
    106             lin_min, lin_max, lin_steps = -100.1, 100.1, 2000 
     109            start, stop, steps = -100.1, 100.1, 2000 
    107110        elif xrange == "linear": 
    108             lin_min, lin_max, lin_steps = 1, 1000, 2000 
    109             lin_min, lin_max, lin_steps = 0.001, 2, 2000 
     111            start, stop, steps = 1, 1000, 2000 
     112            start, stop, steps = 0.001, 2, 2000 
    110113        elif xrange == "log": 
    111             log_min, log_max, log_steps = -3, 5, 400 
     114            start, stop, steps = -3, 5, 400 
    112115        elif xrange == "logq": 
    113             log_min, log_max, log_steps = -4, 1, 400 
     116            start, stop, steps = -4, 1, 400 
     117        elif ':' in xrange: 
     118            parts = xrange.split(':') 
     119            linear = parts[3] != "log" if len(parts) == 4 else True 
     120            steps = int(parts[2]) if len(parts) > 2 else 400 
     121            start = float(parts[0]) 
     122            stop = float(parts[1]) 
     123 
    114124        else: 
    115125            raise ValueError("unknown range "+xrange) 
    … …  
    121131            # value to x in the given precision. 
    122132            if linear: 
    123                 lin_min = max(lin_min, self.limits[0]) 
    124                 lin_max = min(lin_max, self.limits[1]) 
    125                 qrf = np.linspace(lin_min, lin_max, lin_steps, dtype='single') 
    126                 #qrf = np.linspace(lin_min, lin_max, lin_steps, dtype='double') 
     133                start = max(start, self.limits[0]) 
     134                stop = min(stop, self.limits[1]) 
     135                qrf = np.linspace(start, stop, steps, dtype='single') 
     136                #qrf = np.linspace(start, stop, steps, dtype='double') 
    127137                qr = [mp.mpf(float(v)) for v in qrf] 
    128                 #qr = mp.linspace(lin_min, lin_max, lin_steps) 
     138                #qr = mp.linspace(start, stop, steps) 
    129139            else: 
    130                 log_min = np.log10(max(10**log_min, self.limits[0])) 
    131                 log_max = np.log10(min(10**log_max, self.limits[1])) 
    132                 qrf = np.logspace(log_min, log_max, log_steps, dtype='single') 
    133                 #qrf = np.logspace(log_min, log_max, log_steps, dtype='double') 
     140                start = np.log10(max(10**start, self.limits[0])) 
     141                stop = np.log10(min(10**stop, self.limits[1])) 
     142                qrf = np.logspace(start, stop, steps, dtype='single') 
     143                #qrf = np.logspace(start, stop, steps, dtype='double') 
    134144                qr = [mp.mpf(float(v)) for v in qrf] 
    135                 #qr = [10**v for v in mp.linspace(log_min, log_max, log_steps)] 
     145                #qr = [10**v for v in mp.linspace(start, stop, steps)] 
    136146 
    137147        target = self.call_mpmath(qr, bits=500) 
    … …  
    176186    """ 
    177187    if diff == "relative": 
    178         err = np.array([abs((t-a)/t) for t, a in zip(target, actual)], 'd') 
     188        err = np.array([(abs((t-a)/t) if t != 0 else a) for t, a in zip(target, actual)], 'd') 
    179189        #err = np.clip(err, 0, 1) 
    180190        pylab.loglog(x, err, '-', label=label) 
    … …  
    197207    return model_info 
    198208 
     209# Hack to allow second parameter A in two parameter functions 
     210A = 1 
     211def parse_extra_pars(): 
     212    global A 
     213 
     214    A_str = str(A) 
     215    pop = [] 
     216    for k, v in enumerate(sys.argv[1:]): 
     217        if v.startswith("A="): 
     218            A_str = v[2:] 
     219            pop.append(k+1) 
     220    if pop: 
     221        sys.argv = [v for k, v in enumerate(sys.argv) if k not in pop] 
     222        A = float(A_str) 
     223 
     224parse_extra_pars() 
     225 
    199226 
    200227# =============== FUNCTION DEFINITIONS ================ 
    … …  
    297324    ocl_function=make_ocl("return sas_gamma(q);", "sas_gamma", ["lib/sas_gamma.c"]), 
    298325    limits=(-3.1, 10), 
     326) 
     327add_function( 
     328    name="gammaln(x)", 
     329    mp_function=mp.loggamma, 
     330    np_function=scipy.special.gammaln, 
     331    ocl_function=make_ocl("return sas_gammaln(q);", "sas_gammaln", ["lib/sas_gammainc.c"]), 
     332    #ocl_function=make_ocl("return lgamma(q);", "sas_gammaln"), 
     333) 
     334add_function( 
     335    name="gammainc(x)", 
     336    mp_function=lambda x, a=A: mp.gammainc(a, a=0, b=x)/mp.gamma(a), 
     337    np_function=lambda x, a=A: scipy.special.gammainc(a, x), 
     338    ocl_function=make_ocl("return sas_gammainc(%.15g,q);"%A, "sas_gammainc", ["lib/sas_gammainc.c"]), 
     339) 
     340add_function( 
     341    name="gammaincc(x)", 
     342    mp_function=lambda x, a=A: mp.gammainc(a, a=x, b=mp.inf)/mp.gamma(a), 
     343    np_function=lambda x, a=A: scipy.special.gammaincc(a, x), 
     344    ocl_function=make_ocl("return sas_gammaincc(%.15g,q);"%A, "sas_gammaincc", ["lib/sas_gammainc.c"]), 
    299345) 
    300346add_function( 
    … …  
    463509lanczos_gamma = """\ 
    464510    const double coeff[] = { 
    465             76.18009172947146,     -86.50532032941677, 
    466             24.01409824083091,     -1.231739572450155, 
     511            76.18009172947146, -86.50532032941677, 
     512            24.01409824083091, -1.231739572450155, 
    467513            0.1208650973866179e-2,-0.5395239384953e-5 
    468514            }; 
    … …  
    475521""" 
    476522add_function( 
    477     name="log gamma(x)", 
     523    name="loggamma(x)", 
    478524    mp_function=mp.loggamma, 
    479525    np_function=scipy.special.gammaln, 
    … …  
    599645 
    600646ALL_FUNCTIONS = set(FUNCTIONS.keys()) 
    601 ALL_FUNCTIONS.discard("loggamma")  # OCL version not ready yet 
     647ALL_FUNCTIONS.discard("loggamma")  # use cephes-based gammaln instead 
    602648ALL_FUNCTIONS.discard("3j1/x:taylor") 
    603649ALL_FUNCTIONS.discard("3j1/x:trig") 
    … …  
    615661    -r indicates that the relative error should be plotted (default), 
    616662    -x<range> indicates the steps in x, where <range> is one of the following 
    617       log indicates log stepping in [10^-3, 10^5] (default) 
    618       logq indicates log stepping in [10^-4, 10^1] 
    619       linear indicates linear stepping in [1, 1000] 
    620       zoom indicates linear stepping in [1000, 1010] 
    621       neg indicates linear stepping in [-100.1, 100.1] 
    622 and name is "all" or one of: 
     663        log indicates log stepping in [10^-3, 10^5] (default) 
     664        logq indicates log stepping in [10^-4, 10^1] 
     665        linear indicates linear stepping in [1, 1000] 
     666        zoom indicates linear stepping in [1000, 1010] 
     667        neg indicates linear stepping in [-100.1, 100.1] 
     668        start:stop:n[:stepping] indicates an n-step plot in [start, stop] 
     669            or [10^start, 10^stop] if stepping is "log" (default n=400) 
     670Some functions (notably gammainc/gammaincc) have an additional parameter A 
     671which can be set from the command line as A=value.  Default is A=1. 
     672 
     673Name is one of: 
    623674    """+names) 
    624675    sys.exit(1) 

  • explore/realspace.py

    r362a66f r5778279  
    99import numpy as np 
    1010from numpy import pi, radians, sin, cos, sqrt 
    11 from numpy.random import poisson, uniform, randn, rand 
     11from numpy.random import poisson, uniform, randn, rand, randint 
    1212from numpy.polynomial.legendre import leggauss 
    1313from scipy.integrate import simps 
    … …  
    7878 
    7979 
     80I3 = np.matrix([[1., 0, 0], [0, 1, 0], [0, 0, 1]]) 
     81 
    8082class Shape: 
    81     rotation = np.matrix([[1., 0, 0], [0, 1, 0], [0, 0, 1]]) 
     83    rotation = I3 
    8284    center = np.array([0., 0., 0.])[:, None] 
    8385    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 = "" 
    8491 
    8592    def volume(self): 
    … …  
    96103 
    97104    def rotate(self, theta, phi, psi): 
    98         self.rotation = rotation(theta, phi, psi) * self.rotation 
     105        if theta != 0. or phi != 0. or psi != 0.: 
     106            self.rotation = rotation(theta, phi, psi) * self.rotation 
    99107        return self 
    100108 
    … …  
    103111        return self 
    104112 
     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 
    105121    def _adjust(self, points): 
    106         points = np.asarray(self.rotation * np.matrix(points.T)) + self.center 
     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) 
    107128        return points.T 
    108129 
    109     def r_bins(self, q, over_sampling=1, r_step=0.): 
    110         r_max = min(2 * pi / q[0], self.r_max) 
     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) 
    111136        if r_step == 0.: 
    112137            r_step = 2 * pi / q[-1] / over_sampling 
    113138        #r_step = 0.01 
    114139        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 
    115188 
    116189class Composite(Shape): 
    … …  
    669742    Iq = 100 * np.ones_like(qx) 
    670743    data = Data2D(x=qx, y=qy, z=Iq, dx=None, dy=None, dz=np.sqrt(Iq)) 
    671     data.x_bins = qx[0,:] 
    672     data.y_bins = qy[:,0] 
     744    data.x_bins = qx[0, :] 
     745    data.y_bins = qy[:, 0] 
    673746    data.filename = "fake data" 
    674747 
    … …  
    695768    return shape, fn, fn_xy 
    696769 
    697 def build_sphere(radius=125, rho=2): 
     770DEFAULT_SPHERE_RADIUS = 125 
     771DEFAULT_SPHERE_CONTRAST = 2 
     772def build_sphere(radius=DEFAULT_SPHERE_RADIUS, rho=DEFAULT_SPHERE_CONTRAST): 
    698773    shape = TriaxialEllipsoid(radius, radius, radius, rho) 
    699774    fn = lambda q: sphere_Iq(q, radius)*rho**2 
    … …  
    751826    return shape, fn, fn_xy 
    752827 
    753 def build_cubic_lattice(shape, nx=1, ny=1, nz=1, dx=2, dy=2, dz=2, 
    754                   shuffle=0, rotate=0): 
     828def build_sc_lattice(shape, nx=1, ny=1, nz=1, dx=2, dy=2, dz=2, 
     829                        shuffle=0, rotate=0): 
    755830    a, b, c = shape.dims 
    756     shapes = [copy(shape) 
     831    corners= [copy(shape) 
    757832              .shift((ix+(randn() if shuffle < 0.3 else rand())*shuffle)*dx*a, 
    758833                     (iy+(randn() if shuffle < 0.3 else rand())*shuffle)*dy*b, 
    … …  
    762837              for iy in range(ny) 
    763838              for iz in range(nz)] 
    764     lattice = Composite(shapes) 
     839    lattice = Composite(corners) 
    765840    return lattice 
    766841 
     842def 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 
     864def 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 
    767901 
    768902SHAPE_FUNCTIONS = OrderedDict([ 
    … …  
    775909]) 
    776910SHAPES = list(SHAPE_FUNCTIONS.keys()) 
     911LATTICE_FUNCTIONS = OrderedDict([ 
     912    ("sc", build_sc_lattice), 
     913    ("bcc", build_bcc_lattice), 
     914    ("fcc", build_fcc_lattice), 
     915]) 
     916LATTICE_TYPES = list(LATTICE_FUNCTIONS.keys()) 
    777917 
    778918def check_shape(title, shape, fn=None, show_points=False, 
    … …  
    783923    r = shape.r_bins(q, r_step=r_step) 
    784924    sampling_density = samples / shape.volume 
     925    print("sampling points") 
    785926    rho, points = shape.sample(sampling_density) 
     927    print("calculating Pr") 
    786928    t0 = time.time() 
    787929    Pr = calc_Pr(r, rho-rho_solvent, points) 
    … …  
    792934    import pylab 
    793935    if show_points: 
    794          plot_points(rho, points); pylab.figure() 
     936        plot_points(rho, points); pylab.figure() 
    795937    plot_calc(r, Pr, q, Iq, theory=theory, title=title) 
    796938    pylab.gcf().canvas.set_window_title(title) 
    … …  
    806948    Qx, Qy = np.meshgrid(qx, qy) 
    807949    sampling_density = samples / shape.volume 
     950    print("sampling points") 
    808951    t0 = time.time() 
    809952    rho, points = shape.sample(sampling_density) 
    … …  
    844987                        help='lattice size') 
    845988    parser.add_argument('-z', '--spacing', type=str, default='2,2,2', 
    846                         help='lattice spacing') 
     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') 
    847993    parser.add_argument('-r', '--rotate', type=float, default=0., 
    848994                        help="rotation relative to lattice, gaussian < 30 degrees, uniform otherwise") 
    … …  
    8581004    nx, ny, nz = [int(v) for v in opts.lattice.split(',')] 
    8591005    dx, dy, dz = [float(v) for v in opts.spacing.split(',')] 
    860     shuffle, rotate = opts.shuffle, opts.rotate 
     1006    distortion, rotation = opts.shuffle, opts.rotate 
    8611007    shape, fn, fn_xy = SHAPE_FUNCTIONS[opts.shape](**pars) 
    862     if nx > 1 or ny > 1 or nz > 1: 
    863         shape = build_cubic_lattice(shape, nx, ny, nz, dx, dy, dz, shuffle, rotate) 
     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 
    8641038    title = "%s(%s)" % (opts.shape, " ".join(opts.pars)) 
    8651039    if opts.dim == 1: 
    … …  
    8671041                    mesh=opts.mesh, qmax=opts.qmax, samples=opts.samples) 
    8681042    else: 
    869         view = tuple(float(v) for v in opts.view.split(',')) 
    8701043        check_shape_2d(title, shape, fn_xy, view=view, show_points=opts.plot, 
    8711044                       mesh=opts.mesh, qmax=opts.qmax, samples=opts.samples) 
Note: See TracChangeset for help on using the changeset viewer.