source: sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/SquareWellStruct_v40.ipf @ 510

Last change on this file since 510 was 510, checked in by srkline, 14 years ago

Simple change in all of the model function files to include the name of the parameter wave in the Keyword=list that is generated when a model is plotted. This is becoming an issue where the proper parameter wave can't be deduced from just the suffix, then there is nothing to put in the table.

I should have added this when I initially wrote the wrapper...

File size: 3.8 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion = 6.0
3
4////////////////////////////////////////////////
5// GaussUtils.ipf is not necessary (no smearing calculations done) PlotUtils.ipf is recommended
6////////////////////////////////////////////////
7//
8// this function calculates the interparticle structure factor for spherical particles interacting
9// through a square well potential.
10//
11// 06 NOV 98 SRK
12////////////////////////////////////////////////
13
14Proc PlotSquareWellStruct(num,qmin,qmax)
15        Variable num=128,qmin=0.001,qmax=0.3
16        Prompt num "Enter number of data points for model: "
17        Prompt qmin "Enter minimum q-value (A^-1) for model: "
18        Prompt qmax "Enter maximum q-value (A^-1) for model: "
19       
20        Make/O/D/n=(num) xwave_sws,ywave_sws
21        xwave_sws = alog(log(qmin) + x*((log(qmax)-log(qmin))/num))             
22        Make/O/D coef_sws = {50,0.04,1.5,1.2}
23        make/o/t parameters_sws = {"Radius (A)","vol fraction","well depth (kT)","well width (diameters)"}
24        Edit parameters_sws,coef_sws
25        Variable/G root:g_sws
26        g_sws := SquareWellStruct(coef_sws,ywave_sws,xwave_sws)
27//      ywave_sws := SquareWellStruct(coef_sws,xwave_sws)
28        Display ywave_sws vs xwave_sws
29        ModifyGraph marker=29,msize=2,mode=4
30        Label bottom "q (A\\S-1\\M)"
31        Label left "Structure Factor"
32        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
33       
34        AddModelToStrings("SquareWellStruct","coef_sws","parameters_sws","sws")
35End
36
37//AAO version
38Function SquareWellStruct(cw,yw,xw) : FitFunc
39        Wave cw,yw,xw
40       
41#if exists("SquareWellStructX")
42        yw = SquareWellStructX(cw,xw)
43#else
44        yw = fSquareWellStruct(cw,xw)
45#endif
46        return(0)
47End
48///////////////////////////////////////////////////////////
49// smearing a structure factor is not appropriate, so it's not done
50///////////////////////////////////////////////////////////////
51// unsmeared model calculation
52///////////////////////////
53Function fSquareWellStruct(w,x) : FitFunc
54        Wave w
55        Variable x
56
57//     SUBROUTINE SQWELL: CALCULATES THE STRUCTURE FACTOR FOR A
58//                        DISPERSION OF MONODISPERSE HARD SPHERES
59//     IN THE Mean Spherical APPROXIMATION ASSUMING THE SPHERES
60//     INTERACT THROUGH A SQUARE WELL POTENTIAL.
61//** not the best choice of closure ** see note below
62//     REFS:  SHARMA,SHARMA, PHYSICA 89A,(1977),212
63//
64//     
65 
66// NOTE - depths >1.5kT and volume fractions > 0.08 give UNPHYSICAL RESULTS
67// when compared to Monte Carlo simulations
68
69// Input variables are:
70        //[0] radius
71        //[1] volume fraction
72        //[2] well depth e/kT, dimensionless, +ve depths are attractive
73        //[3] well width, multiples of diameter
74       
75        Variable req,phis,edibkb,lambda,struc
76        req = w[0]
77        phis = w[1]
78        edibkb = w[2]
79        lambda = w[3]
80       
81//  COMPUTE CONSTANTS
82//  local variables
83        Variable sigma,eta,eta2,eta3,eta4,etam1,etam14,alpha,beta,gamma
84        Variable qvs,sk,sk2,sk3,sk4,t1,t2,t3,t4,ck
85         
86      SIGMA = req*2.
87      ETA = phis
88      ETA2 = ETA*ETA
89      ETA3 = ETA*ETA2
90      ETA4 = ETA*ETA3       
91      ETAM1 = 1. - ETA
92      ETAM14 = ETAM1*ETAM1*ETAM1*ETAM1
93      ALPHA = (  ( (1. + 2.*ETA)^2 ) + ETA3*( ETA-4.0 )  )/ETAM14
94      BETA = -(ETA/3.0) * ( 18. + 20.*ETA - 12.*ETA2 + ETA4 )/ETAM14
95      GAMMA = 0.5*ETA*( (1. + 2.*ETA)^2 + ETA3*(ETA-4.) )/ETAM14
96//
97//  CALCULATE THE STRUCTURE FACTOR
98//
99// the loop over q-values used to be here     
100//      DO 20 I=1,NPTSM
101        QVS = x
102        SK = x*SIGMA
103        SK2 = SK*SK
104        SK3 = SK*SK2
105        SK4 = SK3*SK
106        T1 = ALPHA * SK3 * ( SIN(SK) - SK * COS(SK) )
107        T2 = BETA * SK2 * ( 2.*SK*SIN(SK) - (SK2-2.)*COS(SK) - 2.0 )
108        T3 =   ( 4.0*SK3 - 24.*SK ) * SIN(SK) 
109        T3 = T3 - ( SK4 - 12.0*SK2 + 24.0 )*COS(SK) + 24.0   
110        T3 = GAMMA*T3
111        T4 = -EDIBKB*SK3*(SIN(LAMBDA*SK) - LAMBDA*SK*COS(LAMBDA*SK)+ SK*COS(SK) - SIN(SK) )
112        CK =  -24.0*ETA*( T1 + T2 + T3 + T4 )/SK3/SK3
113        STRUC  = 1./(1.-CK)
114//   20 CONTINUE
115      Return struc
116End
Note: See TracBrowser for help on using the repository browser.