Ignore:
Timestamp:
Feb 4, 2010 10:21:53 PM (13 years ago)
Author:
srkline
Message:

1 - added function in HPMSA to calculate U(r). See the function for details of usage.

2 - fixed bug in NSORT to trim the _res wave at the point when it is loaded.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/HPMSA_v40.ipf

    r570 r629  
    793793      return SofQ 
    794794end 
     795 
     796 
     797////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 
     798// 
     799// calculate U(r) given the HPMSA parameters 
     800// 
     801// pass in the x wave as scaled r/diam values 
     802// 
     803Function fHPMSA_Ur(w,x) 
     804        wave w 
     805        variable x 
    795806       
     807//      variable timer=StartMSTimer 
     808       
     809        variable Elcharge=1.602189e-19          // electron charge in Coulombs (C) 
     810        variable kB=1.380662e-23                                // Boltzman constant in J/K 
     811        variable FrSpPerm=8.85418782E-12        //Permittivity of free space in C^2/(N m^2) 
     812 
     813        variable SofQ, QQ, Qdiam, gammaek, Vp, csalt, ss 
     814        variable VolFrac, SIdiam, diam, Kappa, cs, IonSt 
     815        variable dialec, Perm, beta1, Temp, zz, charge, ierr 
     816        Variable ur,gam,kapSig 
     817 
     818        diam=w[0]               //in A  (not SI .. should force people to think in nm!!!) 
     819        zz = w[1]               //# of charges 
     820        VolFrac=w[2]     
     821        QQ=x                    //in A^-1 (not SI .. should force people to think in nm^-1!!!) 
     822        Temp=w[3]               //in degrees Kelvin 
     823        csalt=w[4]              //in molarity 
     824        dialec=w[5]             // unitless 
     825         
     826 
     827//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 
     828//////////////////////////// convert to USEFUL inputs in SI units                                                // 
     829////////////////////////////    NOTE: easiest to do EVERYTHING in SI units                               // 
     830//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// 
     831        beta1=1/(kB*Temp)               // in Joules^-1 
     832        Perm=dialec*FrSpPerm    //in C^2/(N  m^2) 
     833        charge=zz*Elcharge              //in Coulomb (C) 
     834        SIdiam = diam*1E-10             //in m 
     835        Vp=4*pi/3*(SIdiam/2)^3  //in m^3 
     836        cs=csalt*6.022E23*1E3   //# salt molecules/m^3 
     837 
     838//         Compute the derived values of : 
     839//                       Ionic strength IonSt (in C^2/m^3)   
     840//                      Kappa (Debye-Huckel screening length in m) 
     841//      and             gamma Exp(-k) 
     842        IonSt=0.5 * Elcharge^2*(zz*VolFrac/Vp+2*cs) 
     843        Kappa=sqrt(2*beta1*IonSt/Perm)     //Kappa calc from Ionic strength 
     844//      Kappa=2/SIdiam                                  // Use to compare with HP paper 
     845 
     846        kapSig = kappa*SIDiam 
     847         
     848        gam=beta1*charge^2/(pi*Perm*SIdiam*(2+Kappa*SIdiam)^2)*exp(kapSig) 
     849        ur = gam*exp(-kapSig*x)/x 
     850 
     851 
     852        Return(ur) 
     853End 
Note: See TracChangeset for help on using the changeset viewer.