[754c454] | 1 | double form_volume(double radius); |
---|
| 2 | double Iq(double q,double dnn,double d_factor, double radius,double sld, double solvent_sld); |
---|
| 3 | double Iqxy(double qx, double qy, double dnn, |
---|
| 4 | double d_factor, double radius,double sld, double solvent_sld, |
---|
| 5 | double theta, double phi, double psi); |
---|
| 6 | |
---|
| 7 | double _BCC_Integrand(double q, double dnn, double d_factor, double theta, double phi); |
---|
| 8 | double _BCCeval(double Theta, double Phi, double temp1, double temp3); |
---|
| 9 | double _sphereform(double q, double radius, double sld, double solvent_sld); |
---|
| 10 | |
---|
| 11 | |
---|
| 12 | double _BCC_Integrand(double q, double dnn, double d_factor, double theta, double phi) { |
---|
| 13 | |
---|
| 14 | const double Da = d_factor*dnn; |
---|
| 15 | const double temp1 = q*q*Da*Da; |
---|
| 16 | const double temp3 = q*dnn; |
---|
| 17 | |
---|
| 18 | double retVal = _BCCeval(theta,phi,temp1,temp3)/(4.0*M_PI); |
---|
| 19 | return(retVal); |
---|
| 20 | } |
---|
| 21 | |
---|
| 22 | double _BCCeval(double Theta, double Phi, double temp1, double temp3) { |
---|
| 23 | |
---|
| 24 | double result; |
---|
[82d239a] | 25 | double sin_theta,cos_theta,sin_phi,cos_phi; |
---|
| 26 | SINCOS(Theta, sin_theta, cos_theta); |
---|
| 27 | SINCOS(Phi, sin_phi, cos_phi); |
---|
| 28 | |
---|
| 29 | const double temp6 = sin_theta; |
---|
| 30 | const double temp7 = sin_theta*cos_phi + sin_theta*sin_phi + cos_theta; |
---|
| 31 | const double temp8 = -sin_theta*cos_phi - sin_theta*sin_phi + cos_theta; |
---|
| 32 | const double temp9 = -sin_theta*cos_phi + sin_theta*sin_phi - cos_theta; |
---|
[11ca2ab] | 33 | |
---|
[82d239a] | 34 | const double temp10 = exp((-1.0/8.0)*temp1*(temp7*temp7 + temp8*temp8 + temp9*temp9)); |
---|
[4962519] | 35 | result = cube(1.0 - (temp10*temp10))*temp6 |
---|
[82d239a] | 36 | / ( (1.0 - 2.0*temp10*cos(0.5*temp3*temp7) + temp10*temp10) |
---|
| 37 | * (1.0 - 2.0*temp10*cos(0.5*temp3*temp8) + temp10*temp10) |
---|
| 38 | * (1.0 - 2.0*temp10*cos(0.5*temp3*temp9) + temp10*temp10)); |
---|
[754c454] | 39 | |
---|
| 40 | return (result); |
---|
| 41 | } |
---|
| 42 | |
---|
| 43 | double form_volume(double radius){ |
---|
[ad90df9] | 44 | return sphere_volume(radius); |
---|
[754c454] | 45 | } |
---|
| 46 | |
---|
| 47 | |
---|
| 48 | double Iq(double q, double dnn, |
---|
| 49 | double d_factor, double radius, |
---|
| 50 | double sld, double solvent_sld){ |
---|
| 51 | |
---|
| 52 | //Volume fraction calculated from lattice symmetry and sphere radius |
---|
| 53 | const double s1 = dnn/sqrt(0.75); |
---|
[ad90df9] | 54 | const double latticescale = 2.0*sphere_volume(radius/s1); |
---|
[754c454] | 55 | |
---|
| 56 | const double va = 0.0; |
---|
| 57 | const double vb = 2.0*M_PI; |
---|
| 58 | const double vaj = 0.0; |
---|
| 59 | const double vbj = M_PI; |
---|
| 60 | |
---|
| 61 | double summ = 0.0; |
---|
| 62 | double answer = 0.0; |
---|
| 63 | for(int i=0; i<150; i++) { |
---|
| 64 | //setup inner integral over the ellipsoidal cross-section |
---|
| 65 | double summj=0.0; |
---|
[c95dc908] | 66 | const double zphi = ( Gauss150Z[i]*(vb-va) + va + vb )/2.0; //the outer dummy is phi |
---|
[754c454] | 67 | for(int j=0;j<150;j++) { |
---|
| 68 | //20 gauss points for the inner integral |
---|
[c95dc908] | 69 | double ztheta = ( Gauss150Z[j]*(vbj-vaj) + vaj + vbj )/2.0; //the inner dummy is theta |
---|
| 70 | double yyy = Gauss150Wt[j] * _BCC_Integrand(q,dnn,d_factor,ztheta,zphi); |
---|
[754c454] | 71 | summj += yyy; |
---|
| 72 | } |
---|
| 73 | //now calculate the value of the inner integral |
---|
| 74 | double answer = (vbj-vaj)/2.0*summj; |
---|
| 75 | |
---|
| 76 | //now calculate outer integral |
---|
| 77 | summ = summ+(Gauss150Wt[i] * answer); |
---|
| 78 | } //final scaling is done at the end of the function, after the NT_FP64 case |
---|
| 79 | |
---|
| 80 | answer = (vb-va)/2.0*summ; |
---|
[9aac25d] | 81 | answer = answer*sphere_form(q,radius,sld,solvent_sld)*latticescale; |
---|
[754c454] | 82 | |
---|
| 83 | return answer; |
---|
| 84 | } |
---|
| 85 | |
---|
| 86 | |
---|
[11ca2ab] | 87 | double Iqxy(double qx, double qy, |
---|
| 88 | double dnn, double d_factor, double radius, |
---|
| 89 | double sld, double solvent_sld, |
---|
| 90 | double theta, double phi, double psi) |
---|
| 91 | { |
---|
| 92 | double q, cos_a1, cos_a2, cos_a3; |
---|
| 93 | ORIENT_ASYMMETRIC(qx, qy, theta, phi, psi, q, cos_a3, cos_a2, cos_a1); |
---|
| 94 | |
---|
| 95 | const double a1 = +cos_a3 - cos_a1 + cos_a2; |
---|
| 96 | const double a2 = +cos_a3 + cos_a1 - cos_a2; |
---|
| 97 | const double a3 = -cos_a3 + cos_a1 + cos_a2; |
---|
| 98 | |
---|
| 99 | const double qd = 0.5*q*dnn; |
---|
[0b717c5] | 100 | const double arg = 0.5*square(qd*d_factor)*(a1*a1 + a2*a2 + a3*a3); |
---|
| 101 | const double tanh_qd = tanh(arg); |
---|
| 102 | const double cosh_qd = cosh(arg); |
---|
| 103 | const double Zq = tanh_qd/(1. - cos(qd*a1)/cosh_qd) |
---|
| 104 | * tanh_qd/(1. - cos(qd*a2)/cosh_qd) |
---|
| 105 | * tanh_qd/(1. - cos(qd*a3)/cosh_qd); |
---|
[11ca2ab] | 106 | |
---|
| 107 | const double Fq = sphere_form(q,radius,sld,solvent_sld)*Zq; |
---|
| 108 | //the occupied volume of the lattice |
---|
| 109 | const double lattice_scale = 2.0*sphere_volume(sqrt(0.75)*radius/dnn); |
---|
| 110 | return lattice_scale * Fq; |
---|
| 111 | } |
---|