Changeset 950


Ignore:
Timestamp:
Jun 13, 2013 9:20:50 AM (8 years ago)
Author:
vondreele
Message:

a few MC/SA changes
some work on GSASIImath docs
dummy string replaced with a fleshed out default block
Some filled in.

Location:
trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIImath.py

    r949 r950  
    6565        until other limits are met (ftol, xtol)
    6666
    67     :Returns: (x,cov_x,infodict) where
     67    :returns: (x,cov_x,infodict) where
    6868
    6969      * x : ndarray
     
    147147   
    148148def getVCov(varyNames,varyList,covMatrix):
    149     'Needs a doc string'
     149    '''obtain variance-covariance terms for a set of variables. NB: the varyList
     150    and covMatrix were saved by the last least squares refineemnt so they must match
     151   
     152    :param list varyNames: variable names to find v-cov matric for
     153    :param list varyList: full list of all variables in v-cov matrix
     154    :param nparray covMatrix: full variance-covariance matrix from the last
     155     least squares refinement
     156   
     157    :returns: nparray vcov: variance-covariance matrix for the variables given
     158     in varyNames
     159   
     160    '''
    150161    vcov = np.zeros((len(varyNames),len(varyNames)))
    151162    for i1,name1 in enumerate(varyNames):
     
    158169
    159170def FindAtomIndexByIDs(atomData,IDs,Draw=True):
    160     'Needs a doc string'
     171    '''finds the set of atom array indices for a list of atom IDs. Will search
     172    either the Atom table or the drawAtom table.
     173   
     174    :param list atomData: Atom or drawAtom table containting coordinates, etc.
     175    :param list IDs: atom IDs to be found
     176    :param bool Draw: True if drawAtom table to be searched; False if Atom table
     177    is searched
     178   
     179    :returns: list indx: atom (or drawAtom) indices
     180   
     181    '''
    161182    indx = []
    162183    for i,atom in enumerate(atomData):
     
    168189
    169190def FillAtomLookUp(atomData):
    170     'Needs a doc string'
     191    '''create a dictionary of atom indexes with atom IDs as keys
     192   
     193    :param list atomData: Atom table to be used
     194   
     195    :returns: dict atomLookUp: dictionary of atom indexes with atom IDs as keys
     196   
     197    '''
    171198    atomLookUp = {}
    172199    for iatm,atom in enumerate(atomData):
     
    175202
    176203def GetAtomsById(atomData,atomLookUp,IdList):
    177     'Needs a doc string'
     204    '''gets a list of atoms from Atom table that match a set of atom IDs
     205   
     206    :param list atomData: Atom table to be used
     207    :param dict atomLookUp: dictionary of atom indexes with atom IDs as keys
     208    :param list IdList: atom IDs to be found
     209   
     210    :returns: list atoms: list of atoms found
     211   
     212    '''
    178213    atoms = []
    179214    for id in IdList:
     
    182217   
    183218def GetAtomItemsById(atomData,atomLookUp,IdList,itemLoc,numItems=1):
    184     'Needs a doc string'
     219    '''gets atom parameters for atoms using atom IDs
     220   
     221    :param list atomData: Atom table to be used
     222    :param dict atomLookUp: dictionary of atom indexes with atom IDs as keys
     223    :param list IdList: atom IDs to be found
     224    :param int itemLoc: pointer to desired 1st item in an atom table entry
     225    :param int numItems: number of items to be retrieved
     226   
     227    :returns: type name: description
     228   
     229    '''
    185230    Items = []
    186231    if not isinstance(IdList,list):
     
    194239   
    195240def GetAtomCoordsByID(pId,parmDict,AtLookup,indx):
    196     'Needs a doc string'
     241    '''default doc string
     242   
     243    :param type name: description
     244   
     245    :returns: type name: description
     246   
     247    '''
    197248    pfx = [str(pId)+'::A'+i+':' for i in ['x','y','z']]
    198249    dpfx = [str(pId)+'::dA'+i+':' for i in ['x','y','z']]
     
    205256
    206257def AtomUij2TLS(atomData,atPtrs,Amat,Bmat,rbObj):   #unfinished & not used
    207     'Needs a doc string; unfinished & not used'
     258    '''default doc string
     259   
     260    :param type name: description
     261   
     262    :returns: type name: description
     263   
     264    '''
    208265    for atom in atomData:
    209266        XYZ = np.inner(Amat,atom[cx:cx+3])
     
    212269               
    213270def TLS2Uij(xyz,g,Amat,rbObj):    #not used anywhere, but could be?
    214     'Needs a doc string; not used anywhere'
     271    '''default doc string
     272   
     273    :param type name: description
     274   
     275    :returns: type name: description
     276   
     277    '''
    215278    TLStype,TLS = rbObj['ThermalMotion'][:2]
    216279    Tmat = np.zeros((3,3))
     
    232295       
    233296def AtomTLS2UIJ(atomData,atPtrs,Amat,rbObj):    #not used anywhere, but could be?
    234     'Needs a doc string; not used anywhere'
     297    '''default doc string
     298   
     299    :param type name: description
     300   
     301    :returns: type name: description
     302   
     303    '''
    235304    cx,ct,cs,cia = atPtrs
    236305    TLStype,TLS = rbObj['ThermalMotion'][:2]
     
    259328
    260329def GetXYZDist(xyz,XYZ,Amat):
    261     ''' gets distance from position xyz to all XYZ, xyz & XYZ are np.array
     330    '''gets distance from position xyz to all XYZ, xyz & XYZ are np.array
    262331        and are in crystal coordinates; Amat is crystal to Cart matrix
     332   
     333    :param type name: description
     334   
     335    :returns: type name: description
     336   
    263337    '''
    264338    return np.sqrt(np.sum(np.inner(Amat,XYZ-xyz)**2,axis=0))
    265339
    266340def getAtomXYZ(atoms,cx):
    267     'Needs a doc string; not used anywhere'
     341    '''default doc string
     342   
     343    :param type name: description
     344   
     345    :returns: type name: description
     346   
     347    '''
    268348    XYZ = []
    269349    for atom in atoms:
     
    272352
    273353def UpdateRBXYZ(Bmat,RBObj,RBData,RBType):
    274     ''' Returns crystal coordinates for atoms described by RBObj
     354    '''default doc string
     355   
     356    :param type name: description
     357   
     358    :returns: type name: description
     359   
     360    '''
     361    ''' returns crystal coordinates for atoms described by RBObj
    275362    '''
    276363    RBRes = RBData[RBType][RBObj['RBId']]
     
    294381
    295382def UpdateMCSAxyz(Bmat,MCSA):
    296     'Needs a doc string'
     383    '''default doc string
     384   
     385    :param type name: description
     386   
     387    :returns: type name: description
     388   
     389    '''
    297390    xyz = []
    298391    atTypes = []
     
    329422   
    330423def SetMolCent(model,RBData):
    331     'Needs a doc string'
     424    '''default doc string
     425   
     426    :param type name: description
     427   
     428    :returns: type name: description
     429   
     430    '''
    332431    rideList = []
    333432    RBRes = RBData[model['Type']][model['RBId']]
     
    349448   
    350449def UpdateRBUIJ(Bmat,Cart,RBObj):
    351     ''' Returns atom I/A, Uiso or UIJ for atoms at XYZ as described by RBObj
     450    '''default doc string
     451   
     452    :param type name: description
     453   
     454    :returns: type name: description
     455   
     456    '''
     457    ''' returns atom I/A, Uiso or UIJ for atoms at XYZ as described by RBObj
    352458    '''
    353459    TLStype,TLS = RBObj['ThermalMotion'][:2]
     
    389495   
    390496def GetSHCoeff(pId,parmDict,SHkeys):
    391     'Needs a doc string'
     497    '''default doc string
     498   
     499    :param type name: description
     500   
     501    :returns: type name: description
     502   
     503    '''
    392504    SHCoeff = {}
    393505    for shkey in SHkeys:
     
    397509       
    398510def getMass(generalData):
     511    '''default doc string
     512   
     513    :param type name: description
     514   
     515    :returns: type name: description
     516   
     517    '''
    399518    'Computes mass of unit cell contents'
    400519    mass = 0.
     
    404523
    405524def getDensity(generalData):
    406     'Computes density of unit cell contents'   
     525    '''default doc string
     526   
     527    :param type name: description
     528   
     529    :returns: type name: description
     530   
     531    '''
    407532    mass = getMass(generalData)
    408533    Volume = generalData['Cell'][7]
     
    411536   
    412537def getWave(Parms):
     538    '''default doc string
     539   
     540    :param type name: description
     541   
     542    :returns: type name: description
     543   
     544    '''
    413545    try:
    414546        return Parms['Lam'][1]
     
    421553
    422554def getSyXYZ(XYZ,ops,SGData):
    423     'Needs a doc string'
     555    '''default doc string
     556   
     557    :param type name: description
     558   
     559    :returns: type name: description
     560   
     561    '''
    424562    XYZout = np.zeros_like(XYZ)
    425563    for i,[xyz,op] in enumerate(zip(XYZ,ops)):
     
    442580   
    443581def getRestDist(XYZ,Amat):
    444     'Needs a doc string'
     582    '''default doc string
     583   
     584    :param type name: description
     585   
     586    :returns: type name: description
     587   
     588    '''
    445589    return np.sqrt(np.sum(np.inner(Amat,(XYZ[1]-XYZ[0]))**2))
    446590   
    447591def getRestDeriv(Func,XYZ,Amat,ops,SGData):
    448     'Needs a doc string'
     592    '''default doc string
     593   
     594    :param type name: description
     595   
     596    :returns: type name: description
     597   
     598    '''
    449599    deriv = np.zeros((len(XYZ),3))
    450600    dx = 0.00001
     
    460610
    461611def getRestAngle(XYZ,Amat):
    462     'Needs a doc string'
     612    '''default doc string
     613   
     614    :param type name: description
     615   
     616    :returns: type name: description
     617   
     618    '''
    463619   
    464620    def calcVec(Ox,Tx,Amat):
     
    476632   
    477633def getRestPlane(XYZ,Amat):
    478     'Needs a doc string'
     634    '''default doc string
     635   
     636    :param type name: description
     637   
     638    :returns: type name: description
     639   
     640    '''
    479641    sumXYZ = np.zeros(3)
    480642    for xyz in XYZ:
     
    492654   
    493655def getRestChiral(XYZ,Amat):   
    494     'Needs a doc string'
     656    '''default doc string
     657   
     658    :param type name: description
     659   
     660    :returns: type name: description
     661   
     662    '''
    495663    VecA = np.empty((3,3))   
    496664    VecA[0] = np.inner(XYZ[1]-XYZ[0],Amat)
     
    500668   
    501669def getRestTorsion(XYZ,Amat):
    502     'Needs a doc string'
     670    '''default doc string
     671   
     672    :param type name: description
     673   
     674    :returns: type name: description
     675   
     676    '''
    503677    VecA = np.empty((3,3))
    504678    VecA[0] = np.inner(XYZ[1]-XYZ[0],Amat)
     
    517691   
    518692def calcTorsionEnergy(TOR,Coeff=[]):
    519     'Needs a doc string'
     693    '''default doc string
     694   
     695    :param type name: description
     696   
     697    :returns: type name: description
     698   
     699    '''
    520700    sum = 0.
    521701    if len(Coeff):
     
    532712
    533713def getTorsionDeriv(XYZ,Amat,Coeff):
    534     'Needs a doc string'
     714    '''default doc string
     715   
     716    :param type name: description
     717   
     718    :returns: type name: description
     719   
     720    '''
    535721    deriv = np.zeros((len(XYZ),3))
    536722    dx = 0.00001
     
    548734
    549735def getRestRama(XYZ,Amat):
    550     'Needs a doc string'
     736    '''Computes a pair of torsion angles in a 5 atom string
     737   
     738    :param nparray XYZ: crystallographic coordinates of 5 atoms
     739    :param nparray Amat: crystal to cartesian transformation matrix
     740   
     741    :returns: list (phi,psi) two torsion angles in degrees
     742   
     743    '''
    551744    phi = getRestTorsion(XYZ[:5],Amat)
    552745    psi = getRestTorsion(XYZ[1:],Amat)
     
    554747   
    555748def calcRamaEnergy(phi,psi,Coeff=[]):
    556     'Needs a doc string'
     749    '''Computes pseudo potential energy from a pair of torsion angles and a
     750     numerical description of the potential energy surface. Used to create
     751     penalty function in LS refinement.
     752     
     753     Eval(phi,psi) = C[0]*exp(-V/1000) where
     754     
     755     V = -C[3]*(phi-C[1])**2-C[4]*(psi-C[2])**2-2*(phi-C[1])*(psi-C[2])
     756   
     757    :param float phi: first torsion angle
     758    :param float psi: second torsion angle
     759    :param list Coeff: pseudo potential coefficients
     760   
     761    :returns: list (sum,Eval): pseudo-potential difference from minimum
     762     & value; sum is used for penalty function.
     763   
     764    '''
    557765    sum = 0.
     766    Eval = 0.
    558767    if len(Coeff):
    559768        cof = Coeff.T
     
    572781
    573782def getRamaDeriv(XYZ,Amat,Coeff):
    574     'Needs a doc string'
     783    '''Computes numerical derivatives of torsion angle pair pseudo potential
     784    with respect of crystallographic atom coordinates of the 5 atom sequence
     785   
     786    :param nparray XYZ: crystallographic coordinates of 5 atoms
     787    :param nparray Amat: crystal to cartesian transformation matrix
     788    :param list Coeff: pseudo potential coefficients
     789   
     790    :returns: list (deriv) derivatives of pseudopotential with respect to 5 atom
     791     crystallographic xyz coordinates.
     792   
     793    '''
    575794    deriv = np.zeros((len(XYZ),3))
    576795    dx = 0.00001
     
    588807
    589808def getRestPolefig(ODFln,SamSym,Grid):
    590     'Needs a doc string'
     809    '''default doc string
     810   
     811    :param type name: description
     812   
     813    :returns: type name: description
     814   
     815    '''
    591816    X,Y = np.meshgrid(np.linspace(1.,-1.,Grid),np.linspace(-1.,1.,Grid))
    592817    R,P = np.sqrt(X**2+Y**2).flatten(),atan2d(Y,X).flatten()
     
    598823
    599824def getRestPolefigDerv(HKL,Grid,SHCoeff):
    600     'Needs a doc string'
     825    '''default doc string
     826   
     827    :param type name: description
     828   
     829    :returns: type name: description
     830   
     831    '''
    601832    pass
    602833       
    603834def getDistDerv(Oxyz,Txyz,Amat,Tunit,Top,SGData):
    604     'Needs a doc string'
     835    '''default doc string
     836   
     837    :param type name: description
     838   
     839    :returns: type name: description
     840   
     841    '''
    605842    def calcDist(Ox,Tx,U,inv,C,M,T,Amat):
    606843        TxT = inv*(np.inner(M,Tx)+T)+C+U
     
    628865   
    629866def getAngSig(VA,VB,Amat,SGData,covData={}):
    630     'Needs a doc string'   
     867    '''default doc string
     868   
     869    :param type name: description
     870   
     871    :returns: type name: description
     872   
     873    '''
    631874    def calcVec(Ox,Tx,U,inv,C,M,T,Amat):
    632875        TxT = inv*(np.inner(M,Tx)+T)+C
     
    692935
    693936def GetDistSig(Oatoms,Atoms,Amat,SGData,covData={}):
    694     'Needs a doc string'
     937    '''default doc string
     938   
     939    :param type name: description
     940   
     941    :returns: type name: description
     942   
     943    '''
    695944    def calcDist(Atoms,SyOps,Amat):
    696945        XYZ = []
     
    737986
    738987def GetAngleSig(Oatoms,Atoms,Amat,SGData,covData={}):
    739     'Needs a doc string'       
     988    '''default doc string
     989   
     990    :param type name: description
     991   
     992    :returns: type name: description
     993   
     994    '''
     995
    740996    def calcAngle(Atoms,SyOps,Amat):
    741997        XYZ = []
     
    7871043
    7881044def GetTorsionSig(Oatoms,Atoms,Amat,SGData,covData={}):
    789     'Needs a doc string'
     1045    '''default doc string
     1046   
     1047    :param type name: description
     1048   
     1049    :returns: type name: description
     1050   
     1051    '''
     1052
    7901053    def calcTorsion(Atoms,SyOps,Amat):
    7911054       
     
    8461109       
    8471110def GetDATSig(Oatoms,Atoms,Amat,SGData,covData={}):
    848     'Needs a doc string'   
     1111    '''default doc string
     1112   
     1113    :param type name: description
     1114   
     1115    :returns: type name: description
     1116   
     1117    '''
     1118
    8491119    def calcDist(Atoms,SyOps,Amat):
    8501120        XYZ = []
     
    10101280
    10111281def adjHKLmax(SGData,Hmax):
    1012     'Needs a doc string'       
     1282    '''default doc string
     1283   
     1284    :param type name: description
     1285   
     1286    :returns: type name: description
     1287   
     1288    '''
    10131289    if SGData['SGLaue'] in ['3','3m1','31m','6/m','6/mmm']:
    10141290        Hmax[0] = ((Hmax[0]+3)/6)*6
     
    10211297
    10221298def OmitMap(data,reflData):
    1023     'Needs a doc string'       
     1299    '''default doc string
     1300   
     1301    :param type name: description
     1302   
     1303    :returns: type name: description
     1304   
     1305    '''
    10241306    generalData = data['General']
    10251307    if not generalData['Map']['MapType']:
     
    10591341
    10601342def FourierMap(data,reflData):
    1061     'Needs a doc string'           
     1343    '''default doc string
     1344   
     1345    :param type name: description
     1346   
     1347    :returns: type name: description
     1348   
     1349    '''
    10621350    generalData = data['General']
    10631351    if not generalData['Map']['MapType']:
     
    11211409# map printing for testing purposes
    11221410def printRho(SGLaue,rho,rhoMax):                         
    1123     'Needs a doc string'       
     1411    '''default doc string
     1412   
     1413    :param type name: description
     1414   
     1415    :returns: type name: description
     1416   
     1417    '''
    11241418    dim = len(rho.shape)
    11251419    if dim == 2:
     
    11481442               
    11491443def findOffset(SGData,A,Fhkl):   
    1150     'Needs a doc string'       
     1444    '''default doc string
     1445   
     1446    :param type name: description
     1447   
     1448    :returns: type name: description
     1449   
     1450    '''
    11511451    if SGData['SpGrp'] == 'P 1':
    11521452        return [0,0,0]   
     
    12051505   
    12061506def ChargeFlip(data,reflData,pgbar):
    1207     'Needs a doc string'       
     1507    '''default doc string
     1508   
     1509    :param type name: description
     1510   
     1511    :returns: type name: description
     1512   
     1513    '''
    12081514    generalData = data['General']
    12091515    mapData = generalData['Map']
     
    12871593   
    12881594def SearchMap(data):
    1289     'Needs a doc string'       
     1595    '''Does a search of a density map for peaks meeting the criterion of peak
     1596    height is greater than mapData['cutOff']/100 of mapData['rhoMax'] where
     1597    mapData is data['General']['mapData']; the map is also in mapData.
     1598
     1599    :param data: the phase data structure
     1600    :returns: (peaks,mags,dzeros) where
     1601        * peaks : ndarray
     1602            x,y,z positions of the peaks found in the map
     1603        * mags : ndarray
     1604            the magnitudes of the peaks
     1605        * dzeros : ndarray
     1606            the distance of the peaks from  the unit cell origin
     1607
     1608    '''       
    12901609    rollMap = lambda rho,roll: np.roll(np.roll(np.roll(rho,roll[0],axis=0),roll[1],axis=1),roll[2],axis=2)
    12911610   
     
    14051724
    14061725def PeaksEquiv(data,Ind):
    1407     'Needs a doc string'       
     1726    '''Find the equivalent map peaks for those selected. Works on the
     1727    contents of data['Map Peaks'].
     1728
     1729    :param data: the phase data structure
     1730    :param list Ind: list of selected peak indices
     1731    :returns: augmented list of peaks including those related by symmetry to the
     1732      ones in Ind
     1733
     1734    '''       
    14081735    def Duplicate(xyz,peaks,Amat):
    14091736        if True in [np.allclose(np.inner(Amat,xyz),np.inner(Amat,peak),atol=0.5) for peak in peaks]:
     
    14321759               
    14331760def PeaksUnique(data,Ind):
    1434     'Needs a doc string'       
     1761    '''Finds the symmetry unique set of peaks from those selected. Works on the
     1762    contents of data['Map Peaks'].
     1763
     1764    :param data: the phase data structure
     1765    :param list Ind: list of selected peak indices
     1766    :returns: the list of symmetry unique peaks from among those given in Ind
     1767
     1768    '''       
    14351769#    XYZE = np.array([[equiv[0] for equiv in G2spc.GenAtom(xyz[1:4],SGData,Move=True)] for xyz in mapPeaks]) #keep this!!
    14361770
     
    14701804
    14711805def getCWsig(ins,pos):
     1806    '''default doc string
     1807   
     1808    :param type name: description
     1809   
     1810    :returns: type name: description
     1811   
     1812    '''
    14721813    tp = tand(pos/2.0)
    14731814    return ins['U']*tp**2+ins['V']*tp+ins['W']
    14741815   
    14751816def getCWsigDeriv(pos):
     1817    '''default doc string
     1818   
     1819    :param type name: description
     1820   
     1821    :returns: type name: description
     1822   
     1823    '''
    14761824    tp = tand(pos/2.0)
    14771825    return tp**2,tp,1.0
    14781826   
    14791827def getCWgam(ins,pos):
     1828    '''default doc string
     1829   
     1830    :param type name: description
     1831   
     1832    :returns: type name: description
     1833   
     1834    '''
    14801835    return ins['X']/cosd(pos/2.0)+ins['Y']*tand(pos/2.0)
    14811836   
    14821837def getCWgamDeriv(pos):
     1838    '''default doc string
     1839   
     1840    :param type name: description
     1841   
     1842    :returns: type name: description
     1843   
     1844    '''
    14831845    return 1./cosd(pos/2.0),tand(pos/2.0)
    14841846   
    14851847def getTOFsig(ins,dsp):
     1848    '''default doc string
     1849   
     1850    :param type name: description
     1851   
     1852    :returns: type name: description
     1853   
     1854    '''
    14861855    return ins['sig-0']+ins['sig-1']*dsp**2+ins['sig-q']*dsp
    14871856   
    14881857def getTOFsigDeriv(dsp):
     1858    '''default doc string
     1859   
     1860    :param type name: description
     1861   
     1862    :returns: type name: description
     1863   
     1864    '''
    14891865    return 1.0,dsp**2,dsp
    14901866   
    14911867def getTOFgamma(ins,dsp):
     1868    '''default doc string
     1869   
     1870    :param type name: description
     1871   
     1872    :returns: type name: description
     1873   
     1874    '''
    14921875    return ins['X']*dsp+ins['Y']*dsp**2
    14931876   
    14941877def getTOFgammaDeriv(dsp):
     1878    '''default doc string
     1879   
     1880    :param type name: description
     1881   
     1882    :returns: type name: description
     1883   
     1884    '''
    14951885    return dsp,dsp**2
    14961886   
    14971887def getTOFbeta(ins,dsp):
     1888    '''default doc string
     1889   
     1890    :param type name: description
     1891   
     1892    :returns: type name: description
     1893   
     1894    '''
    14981895    return ins['beta-0']+ins['beta-1']/dsp**4+ins['beta-q']/dsp
    14991896   
    15001897def getTOFbetaDeriv(dsp):
     1898    '''default doc string
     1899   
     1900    :param type name: description
     1901   
     1902    :returns: type name: description
     1903   
     1904    '''
    15011905    return 1.0,1./dsp**4,1./dsp
    15021906   
    15031907def getTOFalpha(ins,dsp):
     1908    '''default doc string
     1909   
     1910    :param type name: description
     1911   
     1912    :returns: type name: description
     1913   
     1914    '''
    15041915    return ins['alpha']/dsp
    15051916   
    15061917def getTOFalphaDeriv(dsp):
     1918    '''default doc string
     1919   
     1920    :param type name: description
     1921   
     1922    :returns: type name: description
     1923   
     1924    '''
    15071925    return 1./dsp
    15081926   
    15091927def setPeakparms(Parms,Parms2,pos,mag,ifQ=False,useFit=False):
    1510     'Needs a doc string'       
     1928    '''default doc string
     1929   
     1930    :param type name: description
     1931   
     1932    :returns: type name: description
     1933   
     1934    '''
    15111935    ind = 0
    15121936    if useFit:
     
    15942018            A _state object to store the function value and x0 found.
    15952019
    1596         Returns
     2020        returns
    15972021        -------
    15982022        x0 : array
     
    16482072        self.c = self.m * exp(-self.n * self.quench)
    16492073
    1650 #    def update_guess(self, x0):
    1651 #        x0 = asarray(x0)
    1652 #        u = squeeze(random.uniform(0.0, 1.0, size=self.dims))
    1653 #        T = self.T
     2074    def update_guess(self, x0):
     2075        x0 = asarray(x0)
     2076        u = squeeze(random.uniform(0.0, 1.0, size=self.dims))
     2077        T = self.T
     2078        y = sign(u-0.5)*T*((1+1.0/T)**abs(2*u-1)-1.0)+1.0
     2079        xc = y*(self.upper - self.lower)/2.0+self.lower
     2080        return xc
    16542081#        y = sign(u-0.5)*T*((1+1.0/T)**abs(2*u-1)-1.0)
    16552082#        xc = y*(self.upper - self.lower)
     
    19352362
    19362363def mcsaSearch(data,RBdata,reflType,reflData,covData,pgbar):
     2364    '''default doc string
     2365   
     2366    :param type name: description
     2367   
     2368    :returns: type name: description
     2369   
     2370    '''
    19372371    gamFW = lambda s,g: math.exp(math.log(s**5+2.69269*s**4*g+2.42843*s**3*g**2+4.47163*s**2*g**3+0.07842*s*g**4+g**5)/5.)
    19382372   
     
    21412575    parmDict['nfixAt'] = len(fixAtoms)       
    21422576    MCSA = generalData['MCSA controls']
     2577    Results = MCSA.get('Results',[])
    21432578    reflName = MCSA['Data source']
    21442579    phaseName = generalData['Name']
     
    22382673    x0 = [parmDict[val] for val in varyList]
    22392674    ifInv = SGData['SGInv']
    2240     results = anneal(mcsaCalc,x0,args=(refs,rcov,ifInv,RBdata,varyList,parmDict),
    2241         schedule=MCSA['Algorithm'], full_output=True,maxiter=MCSA['nRuns'],
    2242         T0=MCSA['Annealing'][0], Tf=MCSA['Annealing'][1],dwell=MCSA['Annealing'][2],
    2243         boltzmann=MCSA['boltzmann'], learn_rate=0.5, feps=MCSA['Annealing'][3],
    2244         quench=MCSA['fast parms'][0], m=MCSA['fast parms'][1], n=MCSA['fast parms'][2],
    2245         lower=lower, upper=upper, slope=MCSA['log slope'],dlg=pgbar)
    2246     print results
    2247            
    2248 #    parmDict.update(zip(varylist,results[0]))           
    2249                        
    2250     return {}
     2675    for i in range(MCSA['Cycles']):     
     2676        results = anneal(mcsaCalc,x0,args=(refs,rcov,ifInv,RBdata,varyList,parmDict),
     2677            schedule=MCSA['Algorithm'], full_output=True,maxiter=MCSA['nRuns'],
     2678            T0=MCSA['Annealing'][0], Tf=MCSA['Annealing'][1],dwell=MCSA['Annealing'][2],
     2679            boltzmann=MCSA['boltzmann'], learn_rate=0.5, feps=MCSA['Annealing'][3],
     2680            quench=MCSA['fast parms'][0], m=MCSA['fast parms'][1], n=MCSA['fast parms'][2],
     2681            lower=lower, upper=upper, slope=MCSA['log slope'],dlg=pgbar)
     2682        Results.append([results[1],results[2],results[0],varyList])
    22512683
    22522684       
     
    23822814
    23832815        :param np.array A,B,C: Cartesian 3-vectors
    2384         :Returns: quaternion & rotation angle in radians q=r+ai+bj+ck
     2816        :returns: quaternion & rotation angle in radians q=r+ai+bj+ck
    23852817    '''
    23862818
  • trunk/GSASIIphsGUI.py

    r949 r950  
    128128            generalData['MCSA controls'] = {'Data source':'','Annealing':[50.,0.001,50,1.e-6],
    129129            'dmin':2.0,'Algorithm':'fast','Jump coeff':[0.95,0.5],'nRuns':50,'boltzmann':1.0,
    130             'fast parms':[1.0,1.0,1.0],'log slope':0.9}
     130            'fast parms':[1.0,1.0,1.0],'log slope':0.9,'Cycles':1,'Results':[]}
    131131# end of patches
    132132        generalData['NoAtoms'] = {}
     
    657657                    pass
    658658                dmin.SetValue("%.3f"%(MCSA['dmin']))          #reset in case of error
    659                
     659
     660            def OnCycles(event):
     661                MCSA['Cycles'] = int(noRuns.GetValue())
     662                               
    660663            def OnNoRuns(event):
    661664                MCSA['nRuns'] = int(noRuns.GetValue())
     
    723726            dmin.Bind(wx.EVT_KILL_FOCUS,OnDmin)
    724727            lineSizer.Add(dmin,0,wx.ALIGN_CENTER_VERTICAL)
     728            lineSizer.Add(wx.StaticText(General,label=' Cycles: '),0,wx.ALIGN_CENTER_VERTICAL)
     729            Cchoice = ['1','2','3','5','10','15','20','30']
     730            cycles = wx.ComboBox(General,-1,value=str(MCSA.get('Cycles',1)),choices=Cchoice,
     731                style=wx.CB_READONLY|wx.CB_DROPDOWN)
     732            cycles.Bind(wx.EVT_COMBOBOX,OnCycles)       
     733            lineSizer.Add(cycles,0,wx.ALIGN_CENTER_VERTICAL)
    725734            mcsaSizer.Add(lineSizer)
    726735            mcsaSizer.Add((5,5),)
     
    42854294        pgbar.SetSize(Size)
    42864295        try:
    4287             MCSAdata['Results'] = G2mth.mcsaSearch(data,RBdata,reflType,reflData,covData,pgbar)
     4296            G2mth.mcsaSearch(data,RBdata,reflType,reflData,covData,pgbar)
    42884297        finally:
    42894298            pgbar.Destroy()
    42904299        if not data['Drawing']:                 #if new drawing - no drawing data!
    42914300            SetupDrawingData()
    4292         print ' MC/SA run finished: best model residual: %f.3'%(0.01)
    42934301        UpdateMCSA()
    42944302
Note: See TracChangeset for help on using the changeset viewer.