Changeset f9245d4 in sasmodels
- Timestamp:
- Mar 20, 2016 9:43:34 AM (9 years ago)
- Branches:
- master, core_shell_microgels, costrafo411, magnetic_model, release_v0.94, release_v0.95, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
- Children:
- f78a2a1
- Parents:
- 0cbcec4
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
sasmodels/kernel_iq.c
r0cbcec4 rf9245d4 13 13 14 14 /* 15 The environment needs to provide the following #defines 15 The environment needs to provide the following #defines: 16 16 17 USE_OPENCL is defined if running in opencl 17 18 KERNEL declares a function to be available externally 18 19 KERNEL_NAME is the name of the function being declared 19 20 NPARS is the number of parameters in the kernel 20 PARAMETER_ TABLEis the declaration of the parameters to the kernel.21 PARAMETER_DECL is the declaration of the parameters to the kernel. 21 22 22 23 Cylinder: 23 24 24 #define PARAMETER_ TABLE\25 #define PARAMETER_DECL \ 25 26 double length; \ 26 27 double radius; \ 27 28 double sld; \ 28 29 double sld_solvent 30 31 Note: scale and background are not included 29 32 30 33 Multi-shell cylinder (10 shell max): … … 65 68 #define INVALID(var) (var.bell_radius > var.radius) 66 69 70 Model with complicated constraints: 71 72 inline bool constrained(p1, p2, p3) { return expression; } 73 #define INVALID(var) constrained(var.p1, var.p2, var.p3) 74 75 IQ_FUNC could be Iq or Iqxy 76 IQ_PARS could be q[i] or q[2*i],q[2*i+1] 77 67 78 Our design supports a limited number of polydispersity loops, wherein 68 79 we need to cycle through the values of the polydispersity, calculate … … 72 83 calculator treats the parameter set as one long vector. 73 84 74 Let's assume we have 6 polydisperse parameters 75 0 : sld_solvent {s1 = constant} 76 1: sld_material {s2 = constant} 77 2: radius {r = vector of 10pts} 78 3: lenghth {l = vector of 30pts} 79 4: background {b = constant} 80 5: scale {s = constant} 85 Let's assume we have 6 parameters in the model, with two polydisperse:: 86 87 0: scale {scl = constant} 88 1: background {bkg = constant} 89 5: length {l = vector of 30pts} 90 4: radius {r = vector of 10pts} 91 3: sld {s = constant/(radius**2*length)} 92 2: sld_solvent {s2 = constant} 93 94 This generates the following call to the kernel (where x stands for an 95 arbitrary value that is not used by the kernel evaluator): 96 97 NPARS = 4 // scale and background are in all models 98 problem { 99 pd_par = {5, 4, x, x} // parameters *radius* and *length* vary 100 pd_length = {30,10,0,0} // *length* has more, so it is first 101 pd_offset = {10,0,x,x} // *length* starts at index 10 in weights 102 pd_stride = {1,30,300,300} // cumulative product of pd length 103 par_offset = {2, 3, 303, 313} // parameter offsets 104 par_coord = {0, 3, 2, 1} // bitmap of parameter dependencies 105 fast_coord_count = 2 // two parameters vary with *length* distribution 106 fast_coord_index = {5, 3, x, x, x, x} 107 } 108 109 weight = { l0, .., l29, r0, .., r9} 110 pars = { scl, bkg, l0, ..., l29, r0, r1, ..., r9, 111 s[l0,r0], ... s[l0,r9], s[l1,r0], ... s[l29,r9] , s2} 112 113 nq = 130 114 q = { q0, q1, ..., q130, x, x } # pad to 8 element boundary 115 result = {r1, ..., r130, norm, vol, vol_norm, x, x, x, x, x, x, x} 116 81 117 82 118 The polydisperse parameters are stored in as an array of parameter 83 119 indices, one for each polydisperse parameter, stored in pd_par[n]. 84 For the above mentioned expample that will give pd_par = {3, 2, x, x},85 where x stands for abitrary number and 3 corresponds to longest parameter (lenght).86 120 Non-polydisperse parameters do not appear in this array. Each polydisperse 87 121 parameter has a weight vector whose length is stored in pd_length[n], 88 In our case pd_length = {30,10,0,0}. The weights are stored in 89 a contiguous vector of weights for all parameters, with the starting position 90 for the each parameter stored in pd_offset[n] (pd_offset = {10,0,x,x}. 91 The values corresponding to the weights are stored together in a 92 separate weights[] vector, with offset stored in par_offset[pd_par[n]]. 93 In the above mentioned example weight = {r0, .., r9, l0, .., l29}. 94 Polydisperse parameters should be stored in decreasing order of length 95 for highest efficiency. 122 The weights are stored in a contiguous vector of weights for all 123 parameters, with the starting position for the each parameter stored 124 in pd_offset[n]. The values corresponding to the weights are stored 125 together in a separate weights[] vector, with offset stored in 126 par_offset[pd_par[n]]. Polydisperse parameters should be stored in 127 decreasing order of length for highest efficiency. 96 128 97 129 We limit the number of polydisperse dimensions to MAX_PD (currently 4). … … 136 168 137 169 If there is no polydispersity we pretend that it is polydisperisty with one 138 parameter. 170 parameter, pd_start=0 and pd_stop=1. We may or may not short circuit the 171 calculation in this case, depending on how much time it saves. 139 172 140 173 The problem details structure can be allocated and sent in as an integer … … 184 217 } ParameterBlock; 185 218 186 219 #define FULL_KERNEL_NAME KERNEL_NAME ## _ ## IQ_FUNC 187 220 KERNEL 188 void KERNEL_NAME(221 void FULL_KERNEL_NAME( 189 222 int nq, // number of q values 190 223 global const ProblemDetails *problem, … … 193 226 global const double *q, // nq q values, with padding to boundary 194 227 global double *result, // nq return values, again with padding 195 global int *offset, // nq+padding space for current parameter index196 228 const double cutoff, // cutoff in the polydispersity weight product 197 229 const int pd_start, // where we are in the polydispersity loop … … 202 234 // Storage for the current parameter values. These will be updated as we 203 235 // walk the polydispersity cube. 204 local ParameterBlock local_pars; 236 local ParameterBlock local_pars; // current parameter values 205 237 const double *parvec = &local_pars; // Alias named parameters with a vector 238 239 local int offset[NPARS-2]; 206 240 207 241 #if defined(USE_SHORTCUT_OPTIMIZATION) … … 209 243 // Shouldn't need to copy!! 210 244 for (int k=0; k < NPARS; k++) { 211 parvec[k] = pars[k ];245 parvec[k] = pars[k+2]; // skip scale and background 212 246 } 213 247 … … 218 252 { 219 253 const double scattering = IQ_FUNC(IQ_PARS, IQ_PARAMETERS); 220 result[i] += local_pars.scale*scattering + local_pars.background;254 result[i] += pars[0]*scattering + pars[1]; 221 255 } 222 256 return; … … 254 288 local int pd_index[PD_MAX]; 255 289 256 // Set the initial index greater than its vector length in order to257 // trigger the reset behaviour that happens at the end the fast loop.290 // Trigger the reset behaviour that happens at the end the fast loop 291 // by setting the initial index >= weight vector length. 258 292 pd_index[0] = pd_length[0]; 259 293 … … 329 363 result[i] *= norm_vol/vol; 330 364 } 331 result[i] = local_pars.scale*result[i]/norm+local_pars.background;365 result[i] = pars[0]*result[i]/norm + pars[1]; 332 366 } 333 367 }
Note: See TracChangeset
for help on using the changeset viewer.