Changeset 5031ca3 in sasmodels
 Timestamp:
 Oct 1, 2016 10:47:44 AM (4 years ago)
 Branches:
 master, core_shell_microgels, costrafo411, magnetic_model, release_v0.94, release_v0.95, ticket1257vesicleproduct, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
 Children:
 e0de72f
 Parents:
 a807206
 Location:
 sasmodels/models
 Files:

 2 deleted
 2 edited
Legend:
 Unmodified
 Added
 Removed

sasmodels/models/core_shell_ellipsoid.c
r2222134 r5031ca3 1 1 double form_volume(double radius_equat_core, 2 double radius_polar_core,3 double radius_equat_shell,4 double radius_polar_shell);2 double polar_core, 3 double equat_shell, 4 double polar_shell); 5 5 double Iq(double q, 6 6 double radius_equat_core, 7 double radius_polar_core,8 double radius_equat_shell,9 double radius_polar_shell,10 double sld_core,11 double s ld_shell,12 double s ld_solvent);7 double x_core, 8 double thick_shell, 9 double x_polar_shell, 10 double core_sld, 11 double shell_sld, 12 double solvent_sld); 13 13 14 14 15 15 double Iqxy(double qx, double qy, 16 16 double radius_equat_core, 17 double radius_polar_core,18 double radius_equat_shell,19 double radius_polar_shell,20 double sld_core,21 double s ld_shell,22 double s ld_solvent,17 double x_core, 18 double thick_shell, 19 double x_polar_shell, 20 double core_sld, 21 double shell_sld, 22 double solvent_sld, 23 23 double theta, 24 24 double phi); … … 26 26 27 27 double form_volume(double radius_equat_core, 28 double radius_polar_core,29 double radius_equat_shell,30 double radius_polar_shell)28 double x_core, 29 double thick_shell, 30 double x_polar_shell) 31 31 { 32 double vol = 4.0*M_PI/3.0*radius_equat_shell*radius_equat_shell*radius_polar_shell; 32 const double equat_shell = radius_equat_core + thick_shell; 33 const double polar_shell = radius_equat_core*x_core + thick_shell*x_polar_shell; 34 double vol = 4.0*M_PI/3.0*equat_shell*equat_shell*polar_shell; 33 35 return vol; 34 36 } 35 37 36 38 static double 37 core_shell_ellipsoid_ kernel(double q,39 core_shell_ellipsoid_xt_kernel(double q, 38 40 double radius_equat_core, 39 double radius_polar_core,40 double radius_equat_shell,41 double radius_polar_shell,42 double sld_core,43 double s ld_shell,44 double s ld_solvent)41 double x_core, 42 double thick_shell, 43 double x_polar_shell, 44 double core_sld, 45 double shell_sld, 46 double solvent_sld) 45 47 { 46 47 //upper and lower integration limits48 48 const double lolim = 0.0; 49 49 const double uplim = 1.0; … … 51 51 double summ = 0.0; //initialize intergral 52 52 53 const double delpc = sld_core  sld_shell; //core  shell 54 const double delps = sld_shell  sld_solvent; //shell  solvent 53 const double delpc = core_sld  shell_sld; //core  shell 54 const double delps = shell_sld  solvent_sld; //shell  solvent 55 56 57 const double polar_core = radius_equat_core*x_core; 58 const double equat_shell = radius_equat_core + thick_shell; 59 const double polar_shell = radius_equat_core*x_core + thick_shell*x_polar_shell; 55 60 56 61 for(int i=0;i<N_POINTS_76;i++) { … … 58 63 double yyy = Gauss76Wt[i] * gfn4(zi, 59 64 radius_equat_core, 60 radius_polar_core,61 radius_equat_shell,62 radius_polar_shell,65 polar_core, 66 equat_shell, 67 polar_shell, 63 68 delpc, 64 69 delps, … … 68 73 69 74 double answer = (uplimlolim)/2.0*summ; 70 71 75 //convert to [cm1] 72 76 answer *= 1.0e4; … … 76 80 77 81 static double 78 core_shell_ellipsoid_ kernel_2d(double q, double q_x, double q_y,82 core_shell_ellipsoid_xt_kernel_2d(double q, double q_x, double q_y, 79 83 double radius_equat_core, 80 double radius_polar_core,81 double radius_equat_shell,82 double radius_polar_shell,83 double sld_core,84 double s ld_shell,85 double s ld_solvent,84 double x_core, 85 double thick_shell, 86 double x_polar_shell, 87 double core_sld, 88 double shell_sld, 89 double solvent_sld, 86 90 double theta, 87 91 double phi) … … 91 95 phi = phi * M_PI_180; 92 96 93 94 97 // ellipsoid orientation, the axis of the rotation is consistent with the ploar axis. 95 98 const double cyl_x = cos(theta) * cos(phi); 96 99 const double cyl_y = sin(theta); 97 100 98 const double sldcs = sld_core  sld_shell;99 const double sldss = s ld_shell sld_solvent;101 const double sldcs = core_sld  shell_sld; 102 const double sldss = shell_sld solvent_sld; 100 103 101 104 // Compute the angle btw vector q and the … … 103 106 const double cos_val = cyl_x*q_x + cyl_y*q_y; 104 107 105 // Call the IGOR library function to get the kernel: MUST use gfn4 not gf2 because of the def of params. 108 const double polar_core = radius_equat_core*x_core; 109 const double equat_shell = radius_equat_core + thick_shell; 110 const double polar_shell = radius_equat_core*x_core + thick_shell*x_polar_shell; 111 112 // Call the IGOR library function to get the kernel: 113 // MUST use gfn4 not gf2 because of the def of params. 106 114 double answer = gfn4(cos_val, 107 115 radius_equat_core, 108 radius_polar_core,109 radius_equat_shell,110 radius_polar_shell,116 polar_core, 117 equat_shell, 118 polar_shell, 111 119 sldcs, 112 120 sldss, … … 121 129 double Iq(double q, 122 130 double radius_equat_core, 123 double radius_polar_core,124 double radius_equat_shell,125 double radius_polar_shell,126 double sld_core,127 double s ld_shell,128 double s ld_solvent)131 double x_core, 132 double thick_shell, 133 double x_polar_shell, 134 double core_sld, 135 double shell_sld, 136 double solvent_sld) 129 137 { 130 double intensity = core_shell_ellipsoid_ kernel(q,138 double intensity = core_shell_ellipsoid_xt_kernel(q, 131 139 radius_equat_core, 132 radius_polar_core,133 radius_equat_shell,134 radius_polar_shell,135 sld_core,136 s ld_shell,137 s ld_solvent);140 x_core, 141 thick_shell, 142 x_polar_shell, 143 core_sld, 144 shell_sld, 145 solvent_sld); 138 146 139 147 return intensity; … … 143 151 double Iqxy(double qx, double qy, 144 152 double radius_equat_core, 145 double radius_polar_core,146 double radius_equat_shell,147 double radius_polar_shell,148 double sld_core,149 double s ld_shell,150 double s ld_solvent,153 double x_core, 154 double thick_shell, 155 double x_polar_shell, 156 double core_sld, 157 double shell_sld, 158 double solvent_sld, 151 159 double theta, 152 160 double phi) … … 154 162 double q; 155 163 q = sqrt(qx*qx+qy*qy); 156 double intensity = core_shell_ellipsoid_ kernel_2d(q, qx/q, qy/q,164 double intensity = core_shell_ellipsoid_xt_kernel_2d(q, qx/q, qy/q, 157 165 radius_equat_core, 158 radius_polar_core,159 radius_equat_shell,160 radius_polar_shell,161 sld_core,162 s ld_shell,163 s ld_solvent,166 x_core, 167 thick_shell, 168 x_polar_shell, 169 core_sld, 170 shell_sld, 171 solvent_sld, 164 172 theta, 165 173 phi); 
sasmodels/models/core_shell_ellipsoid.py
r2222134 r5031ca3 1 1 r""" 2 This model provides the form factor, $P(q)$, for a core shell ellipsoid (below) 3 where the form factor is normalized by the volume of the outer [CHECK]. 2 An alternative version of $P(q)$ for the core_shell_ellipsoid 3 having as parameters the core axial ratio X and a shell thickness, 4 which are more often what we would like to determine. 4 5 5 .. math:: 6 7 P(q) = \text{scale} * \left<f^2\right>/V + \text{background} 8 9 where the volume $V = (4/3)\pi(r_\text{major outer} r_\text{minor outer}^2)$ 10 and the averaging $< >$ is applied over all orientations for 1D. 11 12 .. figure:: img/core_shell_ellipsoid_geometry.png 13 14 The returned value is in units of cm^1, on absolute scale. 6 This model is also better behaved when polydispersity is applied than the four 7 independent radii in core_shell_ellipsoid model. 15 8 16 9 Definition 17 10  18 11 19 The form factor calculated is 12 .. figure:: img/core_shell_ellipsoid_geometry.png 20 13 21 .. math:: 14 The geometric parameters of this model are 22 15 23 P(q) &= \frac{\text{scale}}{V}\int_0^1 24 \leftF(q,r_\text{minor core},r_\text{major core},\alpha) 25 + F(q,r_\text{minor outer},r_\text{major outer},\alpha)\right^2 26 d\alpha 27 + \text{background} 16 *radius_equat_core =* equatorial core radius *= R_minor_core* 28 17 29 \leftF(q,r_\text{minor},r_\text{major},\alpha)\right 30 &=(4\pi/3)r_\text{major}r_\text{minor}^2 \Delta \rho \cdot (3j_1(u)/u) 18 *X_core = polar_core / radius_equat_core = Rmajor_core / Rminor_core* 31 19 32 u &= q\left[ r_\text{major}^2\alpha ^2 33 + r_\text{minor}^2(1\alpha ^2)\right]^{1/2} 20 *Thick_shell = equat_outer  radius_equat_core = Rminor_outer  Rminor_core* 34 21 35 where 22 *XpolarShell = Tpolar_shell / Thick_shell = (Rmajor_outer  Rmajor_core)/ 23 (Rminor_outer  Rminor_core)* 36 24 37 .. math:: 25 In terms of the original radii 38 26 39 j_1(u)=(\sin x  x \cos x)/x^2 27 *polar_core = radius_equat_core * X_core* 40 28 41 To provide easy access to the orientation of the coreshell ellipsoid, 42 we define the axis of the solid ellipsoid using two angles $\theta$ and $\phi$. 43 These angles are defined as for 44 :ref:`cylinder orientation <cylinderangledefinition>`. 45 The contrast is defined as SLD(core)  SLD(shell) and SLD(shell)  SLD(solvent). 29 *equat_shell = radius_equat_core + Thick_shell* 46 30 47 Note: It is the users' responsibility to ensure that shell radii are larger than 48 the core radii, especially if both are polydisperse, in which case the 49 core_shell_ellipsoid_xt model may be much better. 31 *polar_shell = radius_equat_core * X_core + Thick_shell * XpolarShell* 50 32 33 (where we note that "shell" perhaps confusingly, relates to the outer radius) 34 When *X_core < 1* the core is oblate; when *X_core > 1* it is prolate. 35 *X_core = 1* is a spherical core. 51 36 52 .. note:: 53 The 2nd virial coefficient of the solid ellipsoid is calculated based on 54 the *radius_a* (= *radius_polar_shell)* and *radius_b (= radius_equat_shell)* values, 55 and used as the effective radius for *S(Q)* when $P(Q) * S(Q)$ is applied. 37 For a fixed shell thickness *XpolarShell = 1*, to scale the shell thickness 38 prorata with the radius *XpolarShell = X_core*. 56 39 57 .. figure:: img/core_shell_ellipsoid_angle_projection.jpg 40 When including an $S(q)$, the radius in $S(q)$ is calculated to be that of 41 a sphere with the same 2nd virial coefficient of the outer surface of the 42 ellipsoid. This may have some undesirable effects if the aspect ratio of the 43 ellipsoid is large (ie, if $X << 1$ or $X >> 1$ ), when the $S(q)$ 44  which assumes spheres  will not in any case be valid. 58 45 59 The angles for oriented core_shell_ellipsoid. 60 61 Our model uses the form factor calculations implemented in a clibrary provided 62 by the NIST Center for Neutron Research (Kline, 2006).46 If SAS data are in absolute units, and the SLDs are correct, then scale should 47 be the total volume fraction of the "outer particle". When $S(q)$ is introduced 48 this moves to the $S(q)$ volume fraction, and scale should then be 1.0, 49 or contain some other units conversion factor (for example, if you have SAXS data). 63 50 64 51 References 65 52  66 53 67 M Kotlarchyk, S H Chen, *J. Chem. Phys.*, 79 (1983) 2461 54 R K Heenan, 2015, reparametrised the core_shell_ellipsoid model 68 55 69 S J Berr, *Phys. Chem.*, 91 (1987) 476070 56 """ 71 57 72 58 from numpy import inf, sin, cos, pi 73 59 74 name = "core_shell_ellipsoid "60 name = "core_shell_ellipsoid_xt" 75 61 title = "Form factor for an spheroid ellipsoid particle with a core shell structure." 76 62 description = """ 77 [SpheroidCoreShellModel] Calculates the form factor for an spheroid78 ellipsoid particle with a core_shell structure.79 The form factor is averaged over all possible80 orientations of the ellipsoid such that P(q)81 = scale*<f^2>/Vol + bkg, where f is the82 single particle scattering amplitude.83 [Parameters]:84 radius_equat_core = equatorial radius of core, Rminor_core,85 radius_polar_core = polar radius of core, Rmajor_core,86 radius_equat_shell = equatorial radius of shell, Rminor_outer,87 radius_polar_shell = polar radius of shell, Rmajor_outer,88 sld_core = scattering length density of core,89 sld_shell = scattering length density of shell,90 sld_solvent = scattering length density of solvent,91 background = Incoherent bkg92 scale =scale93 Note:It is the users' responsibility to ensure94 that shell radii are larger than core radii,95 especially if both are polydisperse.96 oblate: polar radius < equatorial radius97 prolate : polar radius > equatorial radius63 [core_shell_ellipsoid_xt] Calculates the form factor for an spheroid 64 ellipsoid particle with a core_shell structure. 65 The form factor is averaged over all possible 66 orientations of the ellipsoid such that P(q) 67 = scale*<f^2>/Vol + bkg, where f is the 68 single particle scattering amplitude. 69 [Parameters]: 70 radius_equat_core = equatorial radius of core, 71 x_core = ratio of core polar/equatorial radii, 72 thick_shell = equatorial radius of outer surface, 73 x_polar_shell = ratio of polar shell thickness to equatorial shell thickness, 74 sld_core = SLD_core 75 sld_shell = SLD_shell 76 sld_solvent = SLD_solvent 77 background = Incoherent bkg 78 scale =scale 79 Note:It is the users' responsibility to ensure 80 that shell radii are larger than core radii. 81 oblate: polar radius < equatorial radius 82 prolate : polar radius > equatorial radius  this new model will make this easier 83 and polydispersity integrals more logical (as previously the shell could disappear). 98 84 """ 99 85 category = "shape:ellipsoid" 100 86 101 87 # pylint: disable=badwhitespace, linetoolong 102 # ["name", "units", default, [lower, upper], "type", "description"],88 # ["name", "units", default, [lower, upper], "type", "description"], 103 89 parameters = [ 104 ["radius_equat_core", "Ang", 200, [0, inf], "volume", "Equatorial radius of core, r minorcore"],105 [" radius_polar_core", "Ang", 10, [0, inf], "volume", "Polar radius of core, r major core"],106 [" radius_equat_shell", "Ang", 250, [0, inf], "volume", "Equatorial radius of shell, r minor outer"],107 [" radius_polar_shell", "Ang", 30, [0, inf], "volume", "Polar radius of shell, r major outer"],108 ["sld_core", "1e6/Ang^2", 2, [inf, inf], "sld", "Core scattering length density"],109 ["sld_shell", "1e6/Ang^2", 1, [inf, inf], "sld", "Shell scattering length density"],110 ["sld_solvent", "1e6/Ang^2", 6.3, [inf, inf], "sld", "Solvent scattering length density"],111 ["theta", "degrees", 0, [inf, inf], "orientation", "Oblate orientation wrt incoming beam"],112 ["phi", "degrees", 0, [inf, inf], "orientation", "Oblate orientation in the plane of the detector"],90 ["radius_equat_core","Ang", 20, [0, inf], "volume", "Equatorial radius of core"], 91 ["x_core", "None", 3, [0, inf], "volume", "axial ratio of core, X = r_polar/r_equatorial"], 92 ["thick_shell", "Ang", 30, [0, inf], "volume", "thickness of shell at equator"], 93 ["x_polar_shell", "", 1, [0, inf], "volume", "ratio of thickness of shell at pole to that at equator"], 94 ["sld_core", "1e6/Ang^2", 2, [inf, inf], "sld", "Core scattering length density"], 95 ["sld_shell", "1e6/Ang^2", 1, [inf, inf], "sld", "Shell scattering length density"], 96 ["sld_solvent", "1e6/Ang^2", 6.3, [inf, inf], "sld", "Solvent scattering length density"], 97 ["theta", "degrees", 0, [inf, inf], "orientation", "Oblate orientation wrt incoming beam"], 98 ["phi", "degrees", 0, [inf, inf], "orientation", "Oblate orientation in the plane of the detector"], 113 99 ] 114 100 # pylint: enable=badwhitespace, linetoolong 115 101 116 source = ["lib/sph_j1c.c", "lib/gfn.c", "lib/gauss76.c", "core_shell_ellipsoid.c"] 102 source = ["lib/sph_j1c.c", "lib/gfn.c", "lib/gauss76.c", 103 "core_shell_ellipsoid_xt.c"] 117 104 118 def ER(radius_equat_core, radius_polar_core, radius_equat_shell, radius_polar_shell):105 def ER(radius_equat_core, x_core, thick_shell, x_polar_shell): 119 106 """ 120 107 Returns the effective radius used in the S*P calculation 121 108 """ 122 109 from .ellipsoid import ER as ellipsoid_ER 123 return ellipsoid_ER(radius_polar_shell, radius_equat_shell) 110 polar_outer = radius_equat_core*x_core + thick_shell*x_polar_shell 111 equat_outer = radius_equat_core + thick_shell 112 return ellipsoid_ER(polar_outer, equat_outer) 124 113 125 114 126 demo = dict(scale= 1, background=0.001,127 radius_equat_core=20 0.0,128 radius_polar_core=10.0,129 radius_equat_shell=250.0,130 radius_polar_shell=30.0,115 demo = dict(scale=0.05, background=0.001, 116 radius_equat_core=20.0, 117 x_core=3.0, 118 thick_shell=30.0, 119 x_polar_shell=1.0, 131 120 sld_core=2.0, 132 121 sld_shell=1.0, … … 141 130 142 131 tests = [ 143 # Accuracy tests based on content in test/utest_ other_models.py132 # Accuracy tests based on content in test/utest_coreshellellipsoidXTmodel.py 144 133 [{'radius_equat_core': 200.0, 145 ' radius_polar_core': 20.0,146 ' radius_equat_shell': 250.0,147 ' radius_polar_shell': 30.0,134 'x_core': 0.1, 135 'thick_shell': 50.0, 136 'x_polar_shell': 0.2, 148 137 'sld_core': 2.0, 149 138 'sld_shell': 1.0, … … 154 143 155 144 # Additional tests with larger range of parameters 156 [{'background': 0.01}, 0.1, 8.86741],145 [{'background': 0.01}, 0.1, 11.6915], 157 146 158 147 [{'radius_equat_core': 20.0, 159 ' radius_polar_core': 200.0,160 ' radius_equat_shell': 54.0,161 ' radius_polar_shell': 3.0,148 'x_core': 200.0, 149 'thick_shell': 54.0, 150 'x_polar_shell': 3.0, 162 151 'sld_core': 20.0, 163 152 'sld_shell': 10.0, … … 165 154 'background': 0.0, 166 155 'scale': 1.0, 167 }, 0.01, 26150.4],156 }, 0.01, 8688.53], 168 157 169 [{'background': 0.001}, (0.4, 0.5), 0.00 170471],158 [{'background': 0.001}, (0.4, 0.5), 0.00690673], 170 159 171 160 [{'radius_equat_core': 20.0, 172 ' radius_polar_core': 200.0,173 ' radius_equat_shell': 54.0,174 ' radius_polar_shell': 3.0,161 'x_core': 200.0, 162 'thick_shell': 54.0, 163 'x_polar_shell': 3.0, 175 164 'sld_core': 20.0, 176 165 'sld_shell': 10.0, … … 178 167 'background': 0.01, 179 168 'scale': 0.01, 180 }, (qx, qy), 0. 105764],169 }, (qx, qy), 0.0100002], 181 170 ]
Note: See TracChangeset
for help on using the changeset viewer.