source: sans/Analysis/trunk/Put in User Procedures/SANS_Models_v3.00/HardSphereStruct.ipf @ 56

Last change on this file since 56 was 42, checked in by srkline, 16 years ago

initial checkin of Analysis v.3.00 files

File size: 3.0 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2
3////////////////////////////////////////////////
4// GaussUtils.ipf is not required (no smearing calculation done) PlotUtils.ipf should be included
5////////////////////////////////////////////////
6// this function is for the hard sphere structure factor
7//
8//
9// 06 NOV 98 SRK
10////////////////////////////////////////////////
11
12Proc PlotHardSphereStruct(num,qmin,qmax)
13        Variable num=128,qmin=0.001,qmax=0.3
14        Prompt num "Enter number of data points for model: "
15        Prompt qmin "Enter minimum q-value (^-1) for model: "
16        Prompt qmax "Enter maximum q-value (^-1) for model: "
17       
18        Make/O/D/n=(num) xwave_hss,ywave_hss
19        //xwave_hss = (x+1)*((qmax-qmin)/num)
20        xwave_hss = alog(log(qmin) + x*((log(qmax)-log(qmin))/num))     
21        Make/O/D coef_hss = {50,0.2}
22        make/o/t parameters_hss = {"Radius (A)","vol fraction"}
23        Edit parameters_hss,coef_hss
24        ywave_hss := HardSphereStruct(coef_hss,xwave_hss)
25        Display ywave_hss vs xwave_hss
26        ModifyGraph marker=29,msize=2,mode=4
27        Label bottom "q (\\S-1\\M)"
28        Label left "Structure Factor"
29        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
30End
31
32///////////////////////////////////////////////////////////
33// no smeared version is calculated - it is simply not appropriate
34
35///////////////////////////////////////////////////////////////
36// unsmeared model calculation
37///////////////////////////
38Function HardSphereStruct(w,x) : FitFunc
39        Wave w
40        Variable x
41                       
42//     SUBROUTINE HSSTRCT: CALCULATES THE STRUCTURE FACTOR FOR A
43//                         DISPERSION OF MONODISPERSE HARD SPHERES
44//                         IN THE PERCUS-YEVICK APPROXIMATION
45//
46//     REFS:  PERCUS,YEVICK PHYS. REV. 110 1 (1958)
47//            THIELE J. CHEM PHYS. 39 474 (1968)
48//            WERTHEIM  PHYS. REV. LETT. 47 1462 (1981)
49//
50// Input variables are:
51        //[0] radius
52        //[1] volume fraction
53        //Variable timer=StartMSTimer
54       
55        Variable r,phi,struc
56        r = w[0]
57        phi = w[1]
58
59// Local variables
60        Variable denom,dnum,alpha,beta,gamm,q,a,asq,ath,afor,rca,rsa
61        Variable calp,cbeta,cgam,prefac,c,vstruc
62//  COMPUTE CONSTANTS
63//
64      DENOM = (1.0-PHI)^4
65      DNUM = (1.0 + 2.0*PHI)^2
66      ALPHA = DNUM/DENOM
67      BETA = -6.0*PHI*((1.0 + PHI/2.0)^2)/DENOM
68      GAMM = 0.50*PHI*DNUM/DENOM
69//
70//  CALCULATE THE STRUCTURE FACTOR
71//     
72// loop over q-values used to be here
73//      DO 10 I=1,NPTSM
74        Q = x           // q-value for the calculation is passed in as variable x
75        A = 2.0*Q*R
76//        IF(A.LT.1.0D-10)  A = 1.0D-10
77        ASQ = A*A
78        ATH = ASQ*A
79        AFOR = ATH*A
80        RCA = COS(A)
81        RSA = SIN(A)
82        CALP = ALPHA*(RSA/ASQ - RCA/A)
83        CBETA = BETA*(2.0*RSA/ASQ - (ASQ - 2.0)*RCA/ATH - 2.0/ATH)
84        CGAM = GAMM*(-RCA/A + (4.0/A)*((3.0*ASQ - 6.0)*RCA/AFOR + (ASQ - 6.0)*RSA/ATH + 6.0/AFOR))
85        PREFAC = -24.0*PHI/A
86        C = PREFAC*(CALP + CBETA + CGAM)
87        VSTRUC = 1.0/(1.0-C)
88        STRUC = VSTRUC
89//   10 CONTINUE
90        //Variable elapse=StopMSTimer(timer)
91      //Print "HS struct eval time (s) = ",elapse
92      RETURN Struc
93End
94// End of HardSphereStruct
95
Note: See TracBrowser for help on using the repository browser.