Changeset 5bddd89 in sasmodels for sasmodels/models/core_shell_bicelle.c
- Timestamp:
- Oct 14, 2016 2:05:39 PM (7 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:
- 9068f4c
- Parents:
- 0b717c5
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
sasmodels/models/core_shell_bicelle.c
r0d6e865 r5bddd89 26 26 double form_volume(double radius, double thick_rim, double thick_face, double length) 27 27 { 28 return M_PI*(radius+thick_rim)*(radius+thick_rim)*(length+2 *thick_face);28 return M_PI*(radius+thick_rim)*(radius+thick_rim)*(length+2.0*thick_face); 29 29 } 30 30 … … 39 39 double rhor, 40 40 double rhosolv, 41 double dum) 41 double sin_alpha, 42 double cos_alpha) 42 43 { 43 44 double si1,si2,be1,be2; … … 49 50 const double vol2 = M_PI*(rad+radthick)*(rad+radthick)*2.0*(length+facthick); 50 51 const double vol3 = M_PI*rad*rad*2.0*(length+facthick); 51 double sn,cn; 52 SINCOS(dum, sn, cn); 53 double besarg1 = qq*rad*sn; 54 double besarg2 = qq*(rad+radthick)*sn; 55 double sinarg1 = qq*length*cn; 56 double sinarg2 = qq*(length+facthick)*cn; 52 double besarg1 = qq*rad*sin_alpha; 53 double besarg2 = qq*(rad+radthick)*sin_alpha; 54 double sinarg1 = qq*length*cos_alpha; 55 double sinarg2 = qq*(length+facthick)*cos_alpha; 57 56 58 57 be1 = sas_J1c(besarg1); … … 65 64 vol3*dr3*si2*be1; 66 65 67 const double retval = t*t*s n;66 const double retval = t*t*sin_alpha; 68 67 69 return (retval);68 return retval; 70 69 71 70 } … … 83 82 { 84 83 // set up the integration end points 85 const double uplim = M_PI /4;86 const double halfheight = length/2.0;84 const double uplim = M_PI_4; 85 const double halfheight = 0.5*length; 87 86 88 87 double summ = 0.0; 89 88 for(int i=0;i<N_POINTS_76;i++) { 90 double zi = (Gauss76Z[i] + 1.0)*uplim; 89 double alpha = (Gauss76Z[i] + 1.0)*uplim; 90 double sin_alpha, cos_alpha; // slots to hold sincos function output 91 SINCOS(alpha, sin_alpha, cos_alpha); 91 92 double yyy = Gauss76Wt[i] * bicelle_kernel(qq, rad, radthick, facthick, 92 halfheight, rhoc, rhoh, rhor,rhosolv, zi); 93 halfheight, rhoc, rhoh, rhor, rhosolv, 94 sin_alpha, cos_alpha); 93 95 summ += yyy; 94 96 } … … 96 98 // calculate value of integral to return 97 99 double answer = uplim*summ; 98 return (answer);100 return answer; 99 101 } 100 102 101 103 static double 102 bicelle_kernel_2d(double q , double q_x, double q_y,104 bicelle_kernel_2d(double qx, double qy, 103 105 double radius, 104 106 double thick_rim, … … 112 114 double phi) 113 115 { 114 //convert angle degree to radian 115 theta *= M_PI_180; 116 phi *= M_PI_180; 116 double q, sin_alpha, cos_alpha; 117 ORIENT_SYMMETRIC(qx, qy, theta, phi, q, sin_alpha, cos_alpha); 117 118 118 // Cylinder orientation119 const double cyl_x = sin(theta) * cos(phi);120 const double cyl_y = sin(theta) * sin(phi);121 122 // Compute the angle btw vector q and the axis of the cylinder123 const double cos_val = cyl_x*q_x + cyl_y*q_y;124 const double alpha = acos( cos_val );125 126 // Get the kernel127 119 double answer = bicelle_kernel(q, radius, thick_rim, thick_face, 128 length/2.0, core_sld, face_sld, rim_sld,129 solvent_sld, alpha) / fabs(sin(alpha));120 0.5*length, core_sld, face_sld, rim_sld, 121 solvent_sld, sin_alpha, cos_alpha) / fabs(sin_alpha); 130 122 131 123 answer *= 1.0e-4; … … 162 154 double phi) 163 155 { 164 double q; 165 q = sqrt(qx*qx+qy*qy); 166 double intensity = bicelle_kernel_2d(q, qx/q, qy/q, 156 double intensity = bicelle_kernel_2d(qx, qy, 167 157 radius, 168 158 thick_rim,
Note: See TracChangeset
for help on using the changeset viewer.