Changeset 6196487 in sasmodels
- Timestamp:
- Mar 8, 2017 6:24:32 AM (8 years ago)
- Branches:
- master, core_shell_microgels, costrafo411, magnetic_model, ticket-1257-vesicle-product, ticket_1156, ticket_1265_superball, ticket_822_more_unit_tests
- Children:
- ceb8a8c
- Parents:
- d24e390 (diff), 9ae85f0 (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. - Location:
- sasmodels
- Files:
-
- 28 edited
Legend:
- Unmodified
- Added
- Removed
-
sasmodels/compare.py
rfe25eda rd504bcd 340 340 if pars['radius'] < pars['thick_string']: 341 341 pars['radius'], pars['thick_string'] = pars['thick_string'], pars['radius'] 342 pars['num_pearls'] = math.ceil(pars['num_pearls'])343 342 pass 344 343 … … 353 352 for c in '1234': 354 353 pars['Phi'+c] /= total 355 356 elif name == 'stacked_disks':357 pars['n_stacking'] = math.ceil(pars['n_stacking'])358 354 359 355 def parlist(model_info, pars, is2d): -
sasmodels/models/core_multi_shell.c
r925ad6e rc3ccaec 8 8 } 9 9 10 double 11 form_volume(double core_radius, double n, double thickness[]); 12 double 13 form_volume(double core_radius, double n, double thickness[]) 10 static double 11 form_volume(double core_radius, double fp_n, double thickness[]) 14 12 { 15 13 double r = core_radius; 14 int n = (int)(fp_n+0.5); 16 15 for (int i=0; i < n; i++) { 17 16 r += thickness[i]; … … 20 19 } 21 20 22 double21 static double 23 22 Iq(double q, double core_sld, double core_radius, 24 double solvent_sld, double num_shells, double sld[], double thickness[]); 25 double 26 Iq(double q, double core_sld, double core_radius, 27 double solvent_sld, double num_shells, double sld[], double thickness[]) 23 double solvent_sld, double fp_n, double sld[], double thickness[]) 28 24 { 29 const int n = (int) ceil(num_shells);25 const int n = (int)(fp_n+0.5); 30 26 double f, r, last_sld; 31 27 r = core_radius; -
sasmodels/models/core_multi_shell.py
r925ad6e r5a0b3d7 107 107 Returns the SLD profile *r* (Ang), and *rho* (1e-6/Ang^2). 108 108 """ 109 n = int(n+0.5) 109 110 z = [] 110 111 rho = [] … … 133 134 def ER(radius, n, thickness): 134 135 """Effective radius""" 135 n = int(n[0] ) # n cannot bepolydisperse136 n = int(n[0]+0.5) # n is a control parameter and is not polydisperse 136 137 return np.sum(thickness[:n], axis=0) + radius 137 138 -
sasmodels/models/flexible_cylinder.c
r592343f rc3ccaec 1 1 static double 2 form_volume( length, kuhn_length,radius)2 form_volume(double length, double kuhn_length, double radius) 3 3 { 4 4 return 1.0; -
sasmodels/models/lamellar_hg_stack_caille.c
ra807206 r1c6286d 3 3 */ 4 4 5 double Iq(double qval, 6 double length_tail, 7 double length_head, 8 double Nlayers, 9 double dd, 10 double Cp, 11 double tail_sld, 12 double head_sld, 13 double solvent_sld); 14 15 double Iq(double qval, 16 double length_tail, 17 double length_head, 18 double Nlayers, 19 double dd, 20 double Cp, 21 double tail_sld, 22 double head_sld, 23 double solvent_sld) 5 static double 6 Iq(double qval, 7 double length_tail, 8 double length_head, 9 double fp_Nlayers, 10 double dd, 11 double Cp, 12 double tail_sld, 13 double head_sld, 14 double solvent_sld) 24 15 { 25 double NN; //local variables of coefficient wave16 int Nlayers = (int)(fp_Nlayers+0.5); //cast to an integer for the loop 26 17 double inten,Pq,Sq,alpha,temp,t2; 27 18 //double dQ, dQDefault, t1, t3; 28 int ii,NNint;29 19 // from wikipedia 0.577215664901532860606512090082402431042159335 30 20 const double Euler = 0.577215664901533; // Euler's constant, increased sig figs for new models Feb 2015 … … 32 22 //dQ = dQDefault; // REMOVED UNUSED dQ calculations for new models Feb 2015 33 23 34 NN = trunc(Nlayers); //be sure that NN is an integer 35 36 Pq = (head_sld-solvent_sld)*(sin(qval*(length_head+length_tail))-sin(qval*length_tail)) + 37 (tail_sld-solvent_sld)*sin(qval*length_tail); 24 Pq = (head_sld-solvent_sld)*(sin(qval*(length_head+length_tail))-sin(qval*length_tail)) 25 + (tail_sld-solvent_sld)*sin(qval*length_tail); 38 26 Pq *= Pq; 39 27 Pq *= 4.0/(qval*qval); 40 28 41 NNint = (int)NN; //cast to an integer for the loop42 ii=0;43 29 Sq = 0.0; 44 for(ii=1;ii<=(NNint-1);ii+=1) { 45 46 //fii = (double)ii; //do I really need to do this? - unused variable, removed 18Feb2015 47 30 for(int ii=1; ii < Nlayers; ii++) { 48 31 temp = 0.0; 49 32 alpha = Cp/4.0/M_PI/M_PI*(log(M_PI*ii) + Euler); … … 52 35 //t3 = dQ*dQ*dd*dd*ii*ii; 53 36 54 temp = 1.0- ii/NN;37 temp = 1.0-(double)ii/(double)Nlayers; 55 38 //temp *= cos(dd*qval*ii/(1.0+t1)); 56 39 temp *= cos(dd*qval*ii); 57 40 //if (temp < 0) printf("q=%g: ii=%d, cos(dd*q*ii)=cos(%g) < 0\n",qval,ii,dd*qval*ii); 58 41 //temp *= exp(-1.0*(t2 + t3)/(2.0*(1.0+t1)) ); 59 temp *= exp(-t2/2.0 42 temp *= exp(-t2/2.0); 60 43 //temp /= sqrt(1.0+t1); 61 44 … … 71 54 72 55 inten *= 1.0e-04; // 1/A to 1/cm 73 return (inten);56 return inten; 74 57 } 75 58 -
sasmodels/models/lamellar_hg_stack_caille.py
r7c57861 ra57b31d 98 98 ["length_head", "Ang", 2, [0, inf], "volume", 99 99 "head thickness"], 100 ["Nlayers", "", 30, [ 0, inf], "",100 ["Nlayers", "", 30, [1, inf], "", 101 101 "Number of layers"], 102 102 ["d_spacing", "Ang", 40., [0.0, inf], "volume", -
sasmodels/models/lamellar_stack_caille.c
r0bef47b r1c6286d 3 3 */ 4 4 5 double Iq(double qval, 6 double del, 7 double Nlayers, 8 double dd, 9 double Cp, 10 double sld, 11 double solvent_sld); 12 13 double Iq(double qval, 14 double del, 15 double Nlayers, 16 double dd, 17 double Cp, 18 double sld, 19 double solvent_sld) 5 static double 6 Iq(double qval, 7 double del, 8 double fp_Nlayers, 9 double dd, 10 double Cp, 11 double sld, 12 double solvent_sld) 20 13 { 21 double contr,NN; //local variables of coefficient wave 14 int Nlayers = (int)(fp_Nlayers+0.5); //cast to an integer for the loop 15 double contr; //local variables of coefficient wave 22 16 double inten,Pq,Sq,alpha,temp,t2; 23 17 //double dQ, dQDefault, t1, t3; 24 int ii,NNint;25 18 // from wikipedia 0.577215664901532860606512090082402431042159335 26 19 const double Euler = 0.577215664901533; // Euler's constant, increased sig figs for new models Feb 2015 … … 28 21 //dQ = dQDefault; // REMOVED UNUSED dQ calculations for new models Feb 2015 29 22 30 NN = trunc(Nlayers); //be sure that NN is an integer31 32 23 contr = sld - solvent_sld; 33 24 34 25 Pq = 2.0*contr*contr/qval/qval*(1.0-cos(qval*del)); 35 26 36 NNint = (int)NN; //cast to an integer for the loop37 ii=0;38 27 Sq = 0.0; 39 // the vital "=" in ii<= added March 2015 40 for(ii=1;ii<=(NNint-1);ii+=1) { 41 42 //fii = (double)ii; //do I really need to do this? - unused variable, removed 18Feb2015 43 28 for (int ii=1; ii < Nlayers; ii++) { 44 29 temp = 0.0; 45 30 alpha = Cp/4.0/M_PI/M_PI*(log(M_PI*ii) + Euler); … … 48 33 //t3 = dQ*dQ*dd*dd*ii*ii; 49 34 50 temp = 1.0 -ii/NN;35 temp = 1.0 - (double)ii / (double)Nlayers; 51 36 //temp *= cos(dd*qval*ii/(1.0+t1)); 52 37 temp *= cos(dd*qval*ii); 53 38 //temp *= exp(-1.0*(t2 + t3)/(2.0*(1.0+t1)) ); 54 temp *= exp(-t2/2.0 39 temp *= exp(-t2/2.0); 55 40 //temp /= sqrt(1.0+t1); 56 41 -
sasmodels/models/lamellar_stack_caille.py
r7c57861 ra57b31d 88 88 parameters = [ 89 89 ["thickness", "Ang", 30.0, [0, inf], "volume", "sheet thickness"], 90 ["Nlayers", "", 20, [ 0, inf], "", "Number of layers"],90 ["Nlayers", "", 20, [1, inf], "", "Number of layers"], 91 91 ["d_spacing", "Ang", 400., [0.0,inf], "volume", "lamellar d-spacing of Caille S(Q)"], 92 92 ["Caille_parameter", "1/Ang^2", 0.1, [0.0,0.8], "", "Caille parameter"], -
sasmodels/models/lamellar_stack_paracrystal.c
r4962519 r5467cd8 2 2 3 3 */ 4 double Iq(double qval, 5 double th, 6 double Nlayers, 7 double davg, 8 double pd, 9 double sld, 10 double solvent_sld); 11 double paraCryst_sn(double ww, double qval, double davg, long Nlayers, double an); 12 double paraCryst_an(double ww, double qval, double davg, long Nlayers); 4 double paraCryst_sn(double ww, double qval, double davg, int Nlayers, double an); 5 double paraCryst_an(double ww, double qval, double davg, int Nlayers); 13 6 14 double Iq(double qval, 15 double th, 16 double Nlayers, 17 double davg, 18 double pd, 19 double sld, 20 double solvent_sld) 7 static double 8 Iq(double qval, 9 double th, 10 double fp_Nlayers, 11 double davg, 12 double pd, 13 double sld, 14 double solvent_sld) 21 15 { 22 23 double inten,contr,xn;24 double xi,ww,Pbil,Znq,Snq,an;25 long n1,n2;16 //get the fractional part of Nlayers, to determine the "mixing" of N's 17 int n1 = (int)(fp_Nlayers); //truncate towards zero 18 int n2 = n1 + 1; 19 const double xn = (double)n2 - fp_Nlayers; //fractional contribution of n1 26 20 27 contr = sld - solvent_sld; 28 //get the fractional part of Nlayers, to determine the "mixing" of N's 29 30 n1 = (long)trunc(Nlayers); //rounds towards zero 31 n2 = n1 + 1; 32 xn = (double)n2 - Nlayers; //fractional contribution of n1 33 34 ww = exp(-qval*qval*pd*pd*davg*davg/2.0); 21 const double ww = exp(-0.5*square(qval*pd*davg)); 35 22 36 23 //calculate the n1 contribution 24 double Znq,Snq,an; 37 25 an = paraCryst_an(ww,qval,davg,n1); 38 26 Snq = paraCryst_sn(ww,qval,davg,n1,an); 39 27 40 28 Znq = xn*Snq; 41 29 … … 52 40 // Zq = (1-ww^2)/(1+ww^2-2*ww*cos(qval*davg)) 53 41 54 xi = th/2.0; //use 1/2 the bilayer thickness55 Pbil = (sin(qval*xi)/(qval*xi))*(sin(qval*xi)/(qval*xi));42 const double xi = th/2.0; //use 1/2 the bilayer thickness 43 const double Pbil = square(sas_sinx_x(qval*xi)); 56 44 57 inten = 2.0*M_PI*contr*contr*Pbil*Znq/(qval*qval);58 inten *= 1.0e-04;45 const double contr = sld - solvent_sld; 46 const double inten = 2.0*M_PI*contr*contr*Pbil*Znq/(qval*qval); 59 47 //printf("q=%.7e wwm1=%g ww=%.5e an=% 12.5e Snq=% 12.5e Znq=% 12.5e Pbil=% 12.5e\n",qval,wwm1,ww,an,Snq,Znq,Pbil); 60 return (inten);48 return 1.0e-4*inten; 61 49 } 62 50 63 51 // functions for the lamellar paracrystal model 64 52 double 65 paraCryst_sn(double ww, double qval, double davg, longNlayers, double an) {53 paraCryst_sn(double ww, double qval, double davg, int Nlayers, double an) { 66 54 67 55 double Snq; … … 69 57 Snq = an/( (double)Nlayers*square(1.0+ww*ww-2.0*ww*cos(qval*davg)) ); 70 58 71 return (Snq);59 return Snq; 72 60 } 73 61 74 62 double 75 paraCryst_an(double ww, double qval, double davg, long Nlayers) { 76 63 paraCryst_an(double ww, double qval, double davg, int Nlayers) { 77 64 double an; 78 65 … … 82 69 an += 2.0*pow(ww,(Nlayers+1))*cos((double)(Nlayers+1)*qval*davg); 83 70 84 return (an);71 return an; 85 72 } 86 73 -
sasmodels/models/lamellar_stack_paracrystal.py
r7c57861 ra0168e8 113 113 parameters = [["thickness", "Ang", 33.0, [0, inf], "volume", 114 114 "sheet thickness"], 115 ["Nlayers", "", 20, [ 0, inf], "",115 ["Nlayers", "", 20, [1, inf], "", 116 116 "Number of layers"], 117 117 ["d_spacing", "Ang", 250., [0.0, inf], "", -
sasmodels/models/linear_pearls.c
r925ad6e rc3ccaec 4 4 double radius, 5 5 double edge_sep, 6 double num_pearls,6 double fp_num_pearls, 7 7 double pearl_sld, 8 8 double solvent_sld); … … 11 11 double radius, 12 12 double edge_sep, 13 doublenum_pearls,13 int num_pearls, 14 14 double pearl_sld, 15 15 double solvent_sld); 16 16 17 17 18 double form_volume(double radius, double num_pearls)18 double form_volume(double radius, double fp_num_pearls) 19 19 { 20 int num_pearls = (int)(fp_num_pearls + 0.5); 20 21 // Pearl volume 21 22 double pearl_vol = M_4PI_3 * cube(radius); … … 27 28 double radius, 28 29 double edge_sep, 29 doublenum_pearls,30 int num_pearls, 30 31 double pearl_sld, 31 32 double solvent_sld) 32 33 { 33 double n_contrib;34 34 //relative sld 35 35 double contrast_pearl = pearl_sld - solvent_sld; … … 46 46 double psi = sas_3j1x_x(q * radius); 47 47 48 // N pearls contribution 49 int n_max = num_pearls - 1; 50 n_contrib = num_pearls; 51 for(int num=1; num<=n_max; num++) { 52 n_contrib += (2.0*(num_pearls-num)*sas_sinx_x(q*separation*num)); 48 // N pearls interaction terms 49 double structure_factor = (double)num_pearls; 50 for(int num=1; num<num_pearls; num++) { 51 structure_factor += 2.0*(num_pearls-num)*sas_sinx_x(q*separation*num); 53 52 } 54 53 // form factor for num_pearls 55 double form_factor = 1.0e-4 * n_contrib* square(m_s*psi) / tot_vol;54 double form_factor = 1.0e-4 * structure_factor * square(m_s*psi) / tot_vol; 56 55 57 56 return form_factor; … … 61 60 double radius, 62 61 double edge_sep, 63 double num_pearls,62 double fp_num_pearls, 64 63 double pearl_sld, 65 64 double solvent_sld) 66 65 { 67 66 67 int num_pearls = (int)(fp_num_pearls + 0.5); 68 68 double result = linear_pearls_kernel(q, 69 69 radius, -
sasmodels/models/linear_pearls.py
r925ad6e rc3ccaec 16 16 .. math:: 17 17 18 P(Q) = \frac{ scale}{V}\left[ m_{p}^219 \left(N+2\sum_{n-1}^{N-1}(N-n)\frac{ sin(qnl)}{qnl}\right)20 \left( 3\frac{ sin(qR)-qRcos(qR)}{(qr)^3}\right)^2\right]18 P(Q) = \frac{\text{scale}}{V}\left[ m_{p}^2 19 \left(N+2\sum_{n-1}^{N-1}(N-n)\frac{\sin(qnl)}{qnl}\right) 20 \left( 3\frac{\sin(qR)-qR\cos(qR)}{(qr)^3}\right)^2\right] 21 21 22 22 where the mass $m_p$ is $(SLD_{pearl}-SLD_{solvent})*(volume\ of\ N\ pearls)$. … … 56 56 ["radius", "Ang", 80.0, [0, inf], "", "Radius of the pearls"], 57 57 ["edge_sep", "Ang", 350.0, [0, inf], "", "Length of the string segment - surface to surface"], 58 ["num_pearls", "", 3.0, [ 0, inf], "", "Number of the pearls"],58 ["num_pearls", "", 3.0, [1, inf], "", "Number of the pearls"], 59 59 ["sld", "1e-6/Ang^2", 1.0, [-inf, inf], "sld", "SLD of the pearl spheres"], 60 60 ["sld_solvent", "1e-6/Ang^2", 6.3, [-inf, inf], "sld", "SLD of the solvent"], -
sasmodels/models/multilayer_vesicle.c
r925ad6e rc3ccaec 7 7 double sld_solvent, 8 8 double sld, 9 doublen_pairs)9 int n_pairs) 10 10 { 11 11 //calculate with a loop, two shells at a time … … 47 47 double sld_solvent, 48 48 double sld, 49 double n_pairs)49 double fp_n_pairs) 50 50 { 51 int n_pairs = (int)(fp_n_pairs + 0.5); 51 52 return multilayer_vesicle_kernel(q, 52 53 volfraction, -
sasmodels/models/multilayer_vesicle.py
r925ad6e rc3ccaec 19 19 20 20 .. math:: 21 P(q) = \text{scale} \cdot \frac{V_f}{V_t} F^2(q) + \text{background} 21 22 22 P(q) = \frac{\text{scale.volfraction}}{V_t} F^2(q) + \text{background} 23 24 where 23 for 25 24 26 25 .. math:: 26 F(q) = (\rho_\text{shell}-\rho_\text{solv}) \sum_{i=1}^{n_\text{pairs}} 27 \left[ 28 3V(R_i)\frac{\sin(qR_i)-qR_i\cos(qR_i)}{(qR_i)^3} \\ 29 - 3V(R_i+t_s)\frac{\sin(q(R_i+t_s))-q(R_i+t_s)\cos(q(R_i+t_s))}{(q(R_i+t_s))^3} 30 \right] 27 31 28 F(q) = (\rho_{shell}-\rho_{solv}) \sum_{i=1}^{n\_pairs} \left[ 29 3V(R_i)\frac{\sin(qR_i)-qR_i\cos(qR_i)}{(qR_i)^3} \\ 30 - 3V(R_i+t_s)\frac{\sin(q(R_i+t_s))-q(R_i+t_s)\cos(q(R_i+t_s))}{(q(R_i+t_s))^3} 31 \right] 32 and 32 33 34 .. math:: 35 R_i = r_c + (i-1)(t_s + t_w) 33 36 34 where $R_i = r_c + (i-1)(t_s + t_w)$ 35 36 where $V_t$ is the volume of the whole particle, $V(R)$ is the volume of a sphere 37 of radius $R$, $r_c$ is the radius of the core, $\rho_{shell}$ is the scattering length 38 density of a shell, $\rho_{solv}$ is the scattering length density of the solvent. 37 where $V_f$ is the volume fraction of particles, $V_t$ is the volume of the 38 whole particle, $V(r)$ is the volume of a sphere of radius $r$, $r_c$ is the 39 radius of the core, $\rho_\text{shell}$ is the scattering length density of a 40 shell, $\rho_\text{solv}$ is the scattering length density of the solvent. 39 41 42 The outer most radius, $r_o = R_n + t_s$, is used for both the volume fraction 43 normalization and for the effective radius for *S(Q)* when $P(Q) * S(Q)$ 44 is applied. 40 45 41 46 The 2D scattering intensity is the same as 1D, regardless of the orientation … … 45 50 46 51 q = \sqrt{q_x^2 + q_y^2} 47 48 49 The outer most radius50 51 $radius + n\_pairs * thick\_shell + (n\_pairs- 1) * thick\_solvent$52 53 is used for both the volume fraction normalization and for the54 effective radius for *S(Q)* when $P(Q) * S(Q)$ is applied.55 52 56 53 For information about polarised and magnetic scattering, see … … 70 67 **Author:** NIST IGOR/DANSE **on:** pre 2010 71 68 72 **Last Modified by:** Piotr Rozyczko **on:** Feb 24, 201669 **Last Modified by:** Piotr Rozyczko **on:** Feb 24, 2016 73 70 74 71 **Last Reviewed by:** Paul Butler **on:** March 20, 2016 … … 109 106 source = ["lib/sas_3j1x_x.c", "multilayer_vesicle.c"] 110 107 108 # TODO: the following line does nothing 111 109 polydispersity = ["radius", "n_pairs"] 112 110 -
sasmodels/models/onion.c
r925ad6e rc3ccaec 30 30 31 31 static double 32 form_volume(double radius_core, double n , double thickness[])32 form_volume(double radius_core, double n_shells, double thickness[]) 33 33 { 34 int i;34 int n = (int)(n_shells+0.5); 35 35 double r = radius_core; 36 for (i =0; i < n; i++) {36 for (int i=0; i < n; i++) { 37 37 r += thickness[i]; 38 38 } -
sasmodels/models/onion.py
r925ad6e rc3ccaec 323 323 Returns shape profile with x=radius, y=SLD. 324 324 """ 325 325 n_shells = int(n_shells+0.5) 326 326 total_radius = 1.25*(sum(thickness[:n_shells]) + radius_core + 1) 327 327 dz = total_radius/400 # 400 points for a smooth plot … … 366 366 return np.asarray(z), np.asarray(rho) 367 367 368 def ER(radius_core, n , thickness):368 def ER(radius_core, n_shells, thickness): 369 369 """Effective radius""" 370 return np.sum(thickness[:int(n[0])], axis=0) + radius_core 370 n = int(n_shells[0]+0.5) 371 return np.sum(thickness[:n], axis=0) + radius_core 371 372 372 373 demo = { -
sasmodels/models/pearl_necklace.c
r4b541ac r3f853beb 1 1 double form_volume(double radius, double edge_sep, 2 double thick_string, double num_pearls);2 double thick_string, double fp_num_pearls); 3 3 double Iq(double q, double radius, double edge_sep, 4 double thick_string, double num_pearls, double sld,4 double thick_string, double fp_num_pearls, double sld, 5 5 double string_sld, double solvent_sld); 6 6 … … 9 9 // From Igor library 10 10 static double 11 _pearl_necklace_kernel(double q, double radius, double edge_sep, double thick_string,12 doublenum_pearls, double sld_pearl, double sld_string, double sld_solv)11 pearl_necklace_kernel(double q, double radius, double edge_sep, double thick_string, 12 int num_pearls, double sld_pearl, double sld_string, double sld_solv) 13 13 { 14 14 // number of string segments 15 num_pearls = floor(num_pearls + 0.5); //Force integer number of pearls 16 const double num_strings = num_pearls - 1.0; 15 const int num_strings = num_pearls - 1; 17 16 18 17 //each masses: contrast * volume … … 69 68 70 69 double form_volume(double radius, double edge_sep, 71 double thick_string, double num_pearls)70 double thick_string, double fp_num_pearls) 72 71 { 73 num_pearls = floor(num_pearls + 0.5); //Force integer number of pearls 74 75 const double num_strings = num_pearls - 1.0; 72 const int num_pearls = (int)(fp_num_pearls + 0.5); //Force integer number of pearls 73 const int num_strings = num_pearls - 1; 76 74 const double string_vol = edge_sep * M_PI_4 * thick_string * thick_string; 77 75 const double pearl_vol = M_4PI_3 * radius * radius * radius; … … 82 80 83 81 double Iq(double q, double radius, double edge_sep, 84 double thick_string, double num_pearls, double sld,82 double thick_string, double fp_num_pearls, double sld, 85 83 double string_sld, double solvent_sld) 86 84 { 87 const double form = _pearl_necklace_kernel(q, radius, edge_sep, 85 const int num_pearls = (int)(fp_num_pearls + 0.5); //Force integer number of pearls 86 const double form = pearl_necklace_kernel(q, radius, edge_sep, 88 87 thick_string, num_pearls, sld, string_sld, solvent_sld); 89 88 -
sasmodels/models/pearl_necklace.py
r4b541ac r1bd1ea2 82 82 ["thick_string", "Ang", 2.5, [0, inf], "volume", 83 83 "Thickness of the chain linkage"], 84 ["num_pearls", "none", 3, [ 0, inf], "volume",84 ["num_pearls", "none", 3, [1, inf], "volume", 85 85 "Number of pearls in the necklace (must be integer)"], 86 86 ["sld", "1e-6/Ang^2", 1.0, [-inf, inf], "sld", … … 100 100 Redundant with form_volume. 101 101 """ 102 num_pearls = int(num_pearls + 0.5) 102 103 number_of_strings = num_pearls - 1.0 103 104 string_vol = edge_sep * pi * pow((thick_string / 2.0), 2.0) … … 111 112 Calculation for effective radius. 112 113 """ 114 num_pearls = int(num_pearls + 0.5) 113 115 tot_vol = volume(radius, edge_sep, thick_string, num_pearls) 114 116 rad_out = (tot_vol/(4.0/3.0*pi)) ** (1./3.) -
sasmodels/models/raspberry.py
r8e68ea0 rc3ccaec 10 10 Schematic of the raspberry model 11 11 12 In order to calculate the form factor of the entire complex, the self-13 correlation of the large droplet, the self-correlation of the particles, the 14 correlation terms between different particles and the cross terms between large 15 droplet and small particles all need to be calculated.12 In order to calculate the form factor of the entire complex, the 13 self-correlation of the large droplet, the self-correlation of the particles, 14 the correlation terms between different particles and the cross terms between 15 large droplet and small particles all need to be calculated. 16 16 17 Consider two infinitely thin shells of radii R1 and R2 separated by distance r.18 The general structure of the equation is then the form factor of the two shells 19 multiplied by the phase factor that accounts for the separation of their 20 centers.17 Consider two infinitely thin shells of radii $R_1$ and $R_2$ separated by 18 distance $r$. The general structure of the equation is then the form factor 19 of the two shells multiplied by the phase factor that accounts for the 20 separation of their centers. 21 21 22 22 .. math:: -
sasmodels/models/rpa.c
racfb094 reb63414 1 double Iq(double q, double case_num,1 double Iq(double q, double fp_case_num, 2 2 double N[], double Phi[], double v[], double L[], double b[], 3 3 double Kab, double Kac, double Kad, … … 5 5 ); 6 6 7 double Iq(double q, double case_num,7 double Iq(double q, double fp_case_num, 8 8 double N[], // DEGREE OF POLYMERIZATION 9 9 double Phi[], // VOL FRACTION … … 15 15 ) 16 16 { 17 int icase = (int) case_num;17 int icase = (int)(fp_case_num+0.5); 18 18 19 19 double Nab,Nac,Nad,Nbc,Nbd,Ncd; … … 309 309 310 310 //calculate contrast where L[i] is the scattering length of i and D is the matrix 311 //need to verify why the sqrt of Nav rather than just Nav (assuming v is molar volume) 311 //Note that should multiply by Nav to get units of SLD which will become 312 // Nav*2 in the next line (SLD^2) but then normalization in that line would 313 //need to divide by Nav leaving only Nav or sqrt(Nav) before squaring. 312 314 Nav=6.022045e+23; 313 315 Lad=(L[0]/v[0]-L[3]/v[3])*sqrt(Nav); … … 317 319 Intg=Lad*Lad*S11+Lbd*Lbd*S22+Lcd*Lcd*S33+2.0*Lad*Lbd*S12+2.0*Lbd*Lcd*S23+2.0*Lad*Lcd*S13; 318 320 319 //rescale for units of Lij^2 ( in 10e-12 m^2 to m^2 ?)320 Intg *= 1.0e-2 4;321 //rescale for units of Lij^2 (fm^2 to cm^2) 322 Intg *= 1.0e-26; 321 323 322 324 return Intg; -
sasmodels/models/rpa.py
rbb73096 reb63414 1 1 r""" 2 Definition 3 ---------- 4 2 5 Calculates the macroscopic scattering intensity for a multi-component 3 6 homogeneous mixture of polymers using the Random Phase Approximation. … … 24 27 Case 9: A-B-C-D tetra-block copolymer 25 28 26 **NB: these case numbers are different from those in the NIST SANS package!** 29 .. note:: 30 These case numbers are different from those in the NIST SANS package! 27 31 28 Only one case can be used at any one time. 32 USAGE NOTES: 29 33 30 The RPA (mean field) formalism only applies only when the multicomponent 31 polymer mixture is in the homogeneous mixed-phase region. 32 33 **Component D is assumed to be the "background" component (ie, all contrasts 34 are calculated with respect to component D).** So the scattering contrast 35 for a C/D blend = [SLD(component C) - SLD(component D)]\ :sup:`2`. 36 37 Depending on which case is being used, the number of fitting parameters - the 38 segment lengths (ba, bb, etc) and $\chi$ parameters (Kab, Kac, etc) - vary. 39 The *scale* parameter should be held equal to unity. 40 41 The input parameters are the degrees of polymerization, the volume fractions, 42 the specific volumes, and the neutron scattering length densities for each 43 component. 34 * Only one case can be used at any one time. 35 * The RPA (mean field) formalism only applies only when the multicomponent 36 polymer mixture is in the homogeneous mixed-phase region. 37 * **Component D is assumed to be the "background" component (ie, all contrasts 38 are calculated with respect to component D).** So the scattering contrast 39 for a C/D blend = [SLD(component C) - SLD(component D)]\ :sup:`2`. 40 * Depending on which case is being used, the number of fitting parameters can 41 vary. Note that in general the degrees of polymerization, the volume 42 fractions, the molar volumes, and the neutron scattering lengths for each 43 component are obtained from other methods and held fixed while the segment 44 lengths (b\ :sub:`a`, b\ :sub:`b`, etc) and $\chi$ parameters (K\ :sub:`ab`, 45 K\ :sub:`ac`, etc). The *scale* parameter should be held equal to unity. 44 46 45 47 … … 47 49 ---------- 48 50 49 A Z Akcasu, R Klein and B Hammouda, *Macromolecules*, 26 (1993) 413651 .. [#] A Z Akcasu, R Klein and B Hammouda, *Macromolecules*, 26 (1993) 4136 50 52 """ 51 53 … … 53 55 54 56 name = "rpa" 55 title = "Random Phase Approximation - unfinished work in progress"57 title = "Random Phase Approximation" 56 58 description = """ 57 59 This formalism applies to multicomponent polymer mixtures in the … … 90 92 ["N[4]", "", 1000.0, [1, inf], "", "Degree of polymerization"], 91 93 ["Phi[4]", "", 0.25, [0, 1], "", "volume fraction"], 92 ["v[4]", "mL/mol", 100.0, [0, inf], "", " specificvolume"],94 ["v[4]", "mL/mol", 100.0, [0, inf], "", "molar volume"], 93 95 ["L[4]", "fm", 10.0, [-inf, inf], "", "scattering length"], 94 96 ["b[4]", "Ang", 5.0, [0, inf], "", "segment length"], … … 114 116 Return a list of parameters to hide depending on the multiplicity parameter. 115 117 """ 118 case_num = int(case_num+0.5) 116 119 if case_num < 2: 117 120 return HIDE_AB -
sasmodels/models/spherical_sld.c
r925ad6e rc3ccaec 1 1 static double form_volume( 2 intn_shells,2 double fp_n_shells, 3 3 double thickness[], 4 4 double interface[]) 5 5 { 6 int n_shells= (int)(fp_n_shells + 0.5); 6 7 double r = 0.0; 7 8 for (int i=0; i < n_shells; i++) { … … 20 21 return pow(z, nu); 21 22 } else if (shape==2) { 22 return 1.0 - pow(1. - z, nu);23 return 1.0 - pow(1.0 - z, nu); 23 24 } else if (shape==3) { 24 25 return expm1(-nu*z)/expm1(-nu); … … 44 45 static double Iq( 45 46 double q, 46 intn_shells,47 double fp_n_shells, 47 48 double sld_solvent, 48 49 double sld[], … … 51 52 double shape[], 52 53 double nu[], 53 intn_steps)54 double fp_n_steps) 54 55 { 55 56 // iteration for # of shells + core + solvent 57 int n_shells = (int)(fp_n_shells + 0.5); 58 int n_steps = (int)(fp_n_steps + 0.5); 56 59 double f=0.0; 57 60 double r=0.0; -
sasmodels/models/spherical_sld.py
r925ad6e rc3ccaec 233 233 """ 234 234 235 n_shells = int(n_shells + 0.5) 236 n_steps = int(n_steps + 0.5) 235 237 z = [] 236 238 rho = [] … … 240 242 rho.append(sld[0]) 241 243 242 for i in range(0, int(n_shells)):244 for i in range(0, n_shells): 243 245 z_next += thickness[i] 244 246 z.append(z_next) … … 261 263 def ER(n_shells, thickness, interface): 262 264 """Effective radius""" 263 n_shells = int(n_shells )265 n_shells = int(n_shells + 0.5) 264 266 total = (np.sum(thickness[:n_shells], axis=1) 265 267 + np.sum(interface[:n_shells], axis=1)) -
sasmodels/models/stacked_disks.c
r6c3e266 r19f996b 1 double form_volume(double thick_core, 2 double thick_layer, 3 double radius, 4 double n_stacking); 5 6 double Iq(double q, 7 double thick_core, 8 double thick_layer, 9 double radius, 10 double n_stacking, 11 double sigma_dnn, 12 double core_sld, 13 double layer_sld, 14 double solvent_sld); 15 16 double Iqxy(double qx, double qy, 17 double thick_core, 18 double thick_layer, 19 double radius, 20 double n_stacking, 21 double sigma_dnn, 22 double core_sld, 23 double layer_sld, 24 double solvent_sld, 25 double theta, 26 double phi); 27 28 static 29 double _kernel(double q, 30 double radius, 31 double core_sld, 32 double layer_sld, 33 double solvent_sld, 34 double halfheight, 35 double thick_layer, 36 double sin_alpha, 37 double cos_alpha, 38 double sigma_dnn, 39 double d, 40 double n_stacking) 1 static double stacked_disks_kernel( 2 double q, 3 double halfheight, 4 double thick_layer, 5 double radius, 6 int n_stacking, 7 double sigma_dnn, 8 double core_sld, 9 double layer_sld, 10 double solvent_sld, 11 double sin_alpha, 12 double cos_alpha, 13 double d) 41 14 42 15 { … … 88 61 89 62 90 static 91 double stacked_disks_kernel(double q,92 93 94 95 doublen_stacking,96 97 98 99 63 static double stacked_disks_1d( 64 double q, 65 double thick_core, 66 double thick_layer, 67 double radius, 68 int n_stacking, 69 double sigma_dnn, 70 double core_sld, 71 double layer_sld, 72 double solvent_sld) 100 73 { 101 74 /* StackedDiscsX : calculates the form factor of a stacked "tactoid" of core shell disks … … 111 84 double sin_alpha, cos_alpha; // slots to hold sincos function output 112 85 SINCOS(zi, sin_alpha, cos_alpha); 113 double yyy = _kernel(q, 86 double yyy = stacked_disks_kernel(q, 87 halfheight, 88 thick_layer, 114 89 radius, 90 n_stacking, 91 sigma_dnn, 115 92 core_sld, 116 93 layer_sld, 117 94 solvent_sld, 118 halfheight,119 thick_layer,120 95 sin_alpha, 121 96 cos_alpha, 122 sigma_dnn, 123 d, 124 n_stacking); 97 d); 125 98 summ += Gauss76Wt[i] * yyy * sin_alpha; 126 99 } … … 132 105 } 133 106 134 double form_volume(double thick_core, 135 double thick_layer, 136 double radius, 137 double n_stacking){ 107 static double form_volume( 108 double thick_core, 109 double thick_layer, 110 double radius, 111 double fp_n_stacking) 112 { 113 int n_stacking = (int)(fp_n_stacking + 0.5); 138 114 double d = 2.0 * thick_layer + thick_core; 139 115 return M_PI * radius * radius * d * n_stacking; 140 116 } 141 117 142 double Iq(double q, 143 double thick_core, 144 double thick_layer, 145 double radius, 146 double n_stacking, 147 double sigma_dnn, 148 double core_sld, 149 double layer_sld, 150 double solvent_sld) 118 static double Iq( 119 double q, 120 double thick_core, 121 double thick_layer, 122 double radius, 123 double fp_n_stacking, 124 double sigma_dnn, 125 double core_sld, 126 double layer_sld, 127 double solvent_sld) 151 128 { 152 return stacked_disks_kernel(q, 129 int n_stacking = (int)(fp_n_stacking + 0.5); 130 return stacked_disks_1d(q, 153 131 thick_core, 154 132 thick_layer, … … 162 140 163 141 164 double 165 Iqxy(double qx, double qy, 166 double thick_core, 167 double thick_layer, 168 double radius, 169 double n_stacking, 170 double sigma_dnn, 171 double core_sld, 172 double layer_sld, 173 double solvent_sld, 174 double theta, 175 double phi) 142 static double Iqxy(double qx, double qy, 143 double thick_core, 144 double thick_layer, 145 double radius, 146 double fp_n_stacking, 147 double sigma_dnn, 148 double core_sld, 149 double layer_sld, 150 double solvent_sld, 151 double theta, 152 double phi) 176 153 { 154 int n_stacking = (int)(fp_n_stacking + 0.5); 177 155 double q, sin_alpha, cos_alpha; 178 156 ORIENT_SYMMETRIC(qx, qy, theta, phi, q, sin_alpha, cos_alpha); … … 180 158 double d = 2.0 * thick_layer + thick_core; 181 159 double halfheight = 0.5*thick_core; 182 double answer = _kernel(q, 160 double answer = stacked_disks_kernel(q, 161 halfheight, 162 thick_layer, 183 163 radius, 164 n_stacking, 165 sigma_dnn, 184 166 core_sld, 185 167 layer_sld, 186 168 solvent_sld, 187 halfheight,188 thick_layer,189 169 sin_alpha, 190 170 cos_alpha, 191 sigma_dnn, 192 d, 193 n_stacking); 171 d); 194 172 195 173 //convert to [cm-1] -
sasmodels/models/stacked_disks.py
rb7e8b94 rc3ccaec 126 126 ["thick_layer", "Ang", 10.0, [0, inf], "volume", "Thickness of layer each side of core"], 127 127 ["radius", "Ang", 15.0, [0, inf], "volume", "Radius of the stacked disk"], 128 ["n_stacking", "", 1.0, [ 0, inf], "volume", "Number of stacked layer/core/layer disks"],128 ["n_stacking", "", 1.0, [1, inf], "volume", "Number of stacked layer/core/layer disks"], 129 129 ["sigma_d", "Ang", 0, [0, inf], "", "Sigma of nearest neighbor spacing"], 130 130 ["sld_core", "1e-6/Ang^2", 4, [-inf, inf], "sld", "Core scattering length density"], -
sasmodels/models/star_polymer.c
r3a48772 r2586093f 3 3 double Iq(double q, double radius2, double arms); 4 4 5 static double _mass_fractal_kernel(double q, double radius2, double arms)5 static double star_polymer_kernel(double q, double radius2, double arms) 6 6 { 7 7 … … 23 23 double Iq(double q, double radius2, double arms) 24 24 { 25 return _mass_fractal_kernel(q, radius2, arms);25 return star_polymer_kernel(q, radius2, arms); 26 26 } -
sasmodels/models/unified_power_Rg.py
r66ca2a6 rc3ccaec 97 97 98 98 def Iq(q, level, rg, power, B, G): 99 ilevel = int(level)100 if ilevel == 0:99 level = int(level + 0.5) 100 if level == 0: 101 101 with errstate(divide='ignore'): 102 102 return 1./q … … 104 104 with errstate(divide='ignore', invalid='ignore'): 105 105 result = np.zeros(q.shape, 'd') 106 for i in range( ilevel):106 for i in range(level): 107 107 exp_now = exp(-(q*rg[i])**2/3.) 108 108 pow_now = (erf(q*rg[i]/sqrt(6.))**3/q)**power[i] 109 if i < ilevel-1:109 if i < level-1: 110 110 exp_next = exp(-(q*rg[i+1])**2/3.) 111 111 else: … … 113 113 result += G[i]*exp_now + B[i]*exp_next*pow_now 114 114 115 result[q == 0] = np.sum(G[: ilevel])115 result[q == 0] = np.sum(G[:level]) 116 116 return result 117 117 -
sasmodels/resolution.py
rb397165 r7b7fcf0 84 84 self.weight_matrix = pinhole_resolution(self.q_calc, self.q, 85 85 np.maximum(q_width, MINIMUM_RESOLUTION)) 86 self.q_calc = abs(self.q_calc) 86 87 87 88 def apply(self, theory): … … 125 126 self.weight_matrix = \ 126 127 slit_resolution(self.q_calc, self.q, qx_width, qy_width) 128 self.q_calc = abs(self.q_calc) 127 129 128 130 def apply(self, theory): … … 153 155 # neither trapezoid nor Simpson's rule improved the accuracy. 154 156 edges = bin_edges(q_calc) 155 edges[edges < 0.0] = 0.0 # clip edges below zero157 #edges[edges < 0.0] = 0.0 # clip edges below zero 156 158 cdf = erf((edges[:, None] - q[None, :]) / (sqrt(2.0)*q_width)[None, :]) 157 159 weights = cdf[1:] - cdf[:-1] … … 286 288 # The current algorithm is a midpoint rectangle rule. 287 289 q_edges = bin_edges(q_calc) # Note: requires q > 0 288 q_edges[q_edges < 0.0] = 0.0 # clip edges below zero290 #q_edges[q_edges < 0.0] = 0.0 # clip edges below zero 289 291 weights = np.zeros((len(q), len(q_calc)), 'd') 290 292
Note: See TracChangeset
for help on using the changeset viewer.