Changeset 961


Ignore:
Timestamp:
Jun 20, 2013 11:57:22 AM (8 years ago)
Author:
vondreele
Message:

fix rb constraints - useCount
fix Uij2Ueqv
mods to anneal algorithms
allow save of MC/SA results from one set of runs to another
use flat shading for polyhedra
a callafter in UpdatePDFGrid - may need more
fixes for when RBs defined but not used

Location:
trunk
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIIconstrGUI.py

    r934 r961  
    987987                del data['Vector'][rbId]
    988988                data['RBIds']['Vector'].remove(rbId)
     989                rbData['useCount'] -= 1
    989990                wx.CallAfter(UpdateVectorRB)
    990991               
  • trunk/GSASIIlattice.py

    r959 r961  
    284284    U = np.inner(Amat,np.inner(U,Amat).T)
    285285    E,R = nl.eigh(U)
    286     return (E[0]+E[2]+E[5])/3.      #lower triangle?
     286    return np.sum(E)/3.
    287287       
    288288def CosSinAngle(U,V,G):
  • trunk/GSASIImath.py

    r960 r961  
    20872087
    20882088
    2089 #  A schedule due to Lester Ingber
     2089#  A schedule due to Lester Ingber modified to use bounds - OK
    20902090class fast_sa(base_schedule):
    20912091    def init(self, **options):
     
    21012101        u = squeeze(random.uniform(0.0, 1.0, size=self.dims))
    21022102        T = self.T
    2103         y = sign(u-0.5)*T*((1+1.0/T)**abs(2*u-1)-1.0)+1.0
    2104         xc = y*(self.upper - self.lower)/2.0+self.lower
    2105         return xc
     2103        xc = (sign(u-0.5)*T*((1+1.0/T)**abs(2*u-1)-1.0)+1.0)/2.0
     2104        xnew = xc*(self.upper - self.lower)+self.lower
     2105        return xnew
    21062106#        y = sign(u-0.5)*T*((1+1.0/T)**abs(2*u-1)-1.0)
    21072107#        xc = y*(self.upper - self.lower)
     
    21142114        return
    21152115
    2116 class cauchy_sa(base_schedule):
    2117 #    def update_guess(self, x0):
    2118 #        x0 = asarray(x0)
     2116class cauchy_sa(base_schedule):     #modified to use bounds - not good
     2117    def update_guess(self, x0):
     2118        x0 = asarray(x0)
     2119        numbers = squeeze(random.uniform(-pi/4, pi/4, size=self.dims))
     2120        xc = (1.+(self.learn_rate * self.T * tan(numbers))%1.)
     2121        xnew = xc*(self.upper - self.lower)+self.lower
     2122        return xnew
    21192123#        numbers = squeeze(random.uniform(-pi/2, pi/2, size=self.dims))
    21202124#        xc = self.learn_rate * self.T * tan(numbers)
     
    21412145        return
    21422146
    2143 class log_sa(base_schedule):
     2147class log_sa(base_schedule):        #OK
    21442148
    21452149    def init(self,**options):
     
    23612365            if abs(af[-1]-best_state.cost) > feps*10:
    23622366                retval = 5
    2363                 print "Warning: Cooled to %f at %s but this is not" \
    2364                       % (squeeze(last_state.cost), str(squeeze(last_state.x))) \
    2365                       + " the smallest point found."
     2367#                print "Warning: Cooled to %f at %s but this is not" \
     2368#                      % (squeeze(last_state.cost), str(squeeze(last_state.x))) \
     2369#                      + " the smallest point found."
    23662370            break
    23672371        if (Tf is not None) and (schedule.T < Tf):
     
    23952399    '''
    23962400    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.)
     2401   
     2402    twopi = 2.0*np.pi
     2403    global tsum
     2404    tsum = 0.
    23972405   
    23982406    def getMDparms(item,pfx,parmDict,varyList):
     
    24632471    def GetAtomM(Xdata,SGData):
    24642472        Mdata = []
    2465         for xyz in Xdata.T:
     2473        for xyz in Xdata:
    24662474            Mdata.append(float(len(G2spc.GenAtom(xyz,SGData))))
    24672475        return np.array(Mdata)
     
    24852493            if parmDict[pfx+'Type'] in ['Vector','Residue']:
    24862494                if parmDict[pfx+'Type'] == 'Vector':
    2487                     RBId = parmDict[pfx+'RBId']
    2488                     RBRes = RBdata['Vector'][RBId]
     2495                    RBRes = RBdata['Vector'][parmDict[pfx+'RBId']]
    24892496                    aTypes = RBRes['rbTypes']
    24902497                    vecs = RBRes['rbVect']
     
    24942501                        Cart += vec*mag
    24952502                elif parmDict[pfx+'Type'] == 'Residue':
    2496                     RBId = parmDict[pfx+'RBId']
    2497                     RBRes = RBdata['Residue'][RBId]
     2503                    RBRes = RBdata['Residue'][parmDict[pfx+'RBId']]
    24982504                    aTypes = RBRes['rbTypes']
    24992505                    Cart = np.array(RBRes['rbXYZ'])
     
    25192525                    if parm in parmDict:
    25202526                        keys[key][atNo] = parmDict[parm]
     2527                iatm += 1
    25212528            else:
    25222529                continue        #skips March Dollase
    2523         return Tdata,Xdata
     2530        return Tdata,Xdata.T
    25242531   
    25252532    def calcMDcorr(MDval,MDaxis,Uniq,G):
     
    25372544#            ParmDict:
    25382545#        puts result F^2 in each ref[8] in refList
    2539 #        '''        
    2540         twopi = 2.0*np.pi
     2546#        '''       
     2547        global tsum
    25412548        parmDict.update(dict(zip(varyList,values)))
    25422549        Tdata,Xdata = GetAtomTX(RBdata,parmDict)
     
    25482555        Srefs = np.zeros(len(refList))
    25492556        sumFcsq = 0.
     2557        t0 = time.time()
    25502558        for refl in refList:
    2551             fbs = 0.
    2552             H = refl[:3]
    2553             for i,El in enumerate(Tdata):
     2559            for i,El in enumerate(Tdata):       #NB: generator here is slower!
    25542560                FF[i] = refl[7][El]
    2555             Uniq = refl[8]
    2556             phi = refl[9]
    2557             phase = twopi*(np.inner(Uniq,(Xdata.T))+phi[:,np.newaxis])
    2558             sinp = np.sin(phase)
    2559             cosp = np.cos(phase)
    2560             occ = Mdata/len(Uniq)
    2561             fa = np.array(FF*occ*cosp)
    2562             fas = np.sum(fa)
    2563             if not ifInv:
    2564                 fb = np.array(FF*occ*sinp)
    2565                 fbs = np.sum(fb)
    2566             fcsq = (fas**2+fbs**2)*refl[3]      #*calcMDcorr(MDval,MDaxis,Uniq,Gmat)
    2567             sumFcsq += fcsq
    2568             refl[5] = fcsq
     2561            FF *= Mdata/len(refl[8])            #FF*occ
     2562            phase = np.inner(refl[8],Xdata)     #hx+ky+lz
     2563            phase += refl[9][:,np.newaxis]      #hx+ky+lz+phi
     2564            cosp = np.cos(twopi*phase)
     2565            fas = np.sum(FF*cosp)
     2566            if ifInv:
     2567                fbs = 0.
     2568            else:
     2569                sinp = np.sin(twopi*phase)
     2570                fbs = np.sum(FF*sinp)
     2571            refl[5] = (fas**2+fbs**2)*refl[3]  #*calcMDcorr(MDval,MDaxis,Uniq,Gmat)
     2572            sumFcsq += refl[5]
     2573        tsum += (time.time()-t0)
    25692574        scale = (parmDict['sumFosq']/sumFcsq)
    25702575        for iref,refl in enumerate(refList):
     
    27002705        quench=MCSA['fast parms'][0], m=MCSA['fast parms'][1], n=MCSA['fast parms'][2],
    27012706        lower=lower, upper=upper, slope=MCSA['log slope'],dlg=pgbar)
    2702     return [False,results[1],results[2],results[0],varyList]
     2707    Result = [False,False,results[1],results[2],]+list(results[0])
     2708    Result.append(varyList)
     2709    return Result,tsum
    27032710
    27042711       
  • trunk/GSASIIphsGUI.py

    r953 r961  
    42084208                    resultsGrid.ForceRefresh()
    42094209                    result = Results[r]
    4210                     for key,val in zip(result[4],result[3]):
     4210                    for key,val in zip(result[-1],result[4:-1]):
    42114211                        vals = key.split(':')
    42124212                        nObj,name = int(vals[0]),vals[1]
     
    42344234                    wx.CallAfter(UpdateMCSA)
    42354235                    G2plt.PlotStructure(G2frame,data)
     4236                elif c == 1:
     4237                    if Results[r][1]:
     4238                        Results[r][1] = False
     4239                    else:
     4240                        Results[r][1] = True
     4241                    resultsTable.SetValue(r,c,Results[r][1])
     4242                    resultsGrid.ForceRefresh()
     4243                       
    42364244               
    42374245            resultsSizer = wx.BoxSizer(wx.VERTICAL)
     
    42394247            resultVals = []
    42404248            for result in Results:
    4241                 maxVary = max(maxVary,len(result[3]))
    4242                 resultVals.append(result[:3]+list(result[3]))
     4249                maxVary = max(maxVary,len(result[-1]))
     4250                resultVals.append(result[:-1])
    42434251            rowLabels = []
    42444252            for i in range(len(Results)): rowLabels.append(str(i))
    4245             colLabels = ['Select','Residual','Tmin',]
    4246             for item in result[4]: colLabels.append(item)
    4247 #            for i in range(maxVary): colLabels.append('variable:'+str(i))
    4248             Types = [wg.GRID_VALUE_BOOL,wg.GRID_VALUE_FLOAT+':10,4',
     4253            colLabels = ['Select','Keep','Residual','Tmin',]
     4254            for item in result[-1]: colLabels.append(item)   #from last result from for loop above
     4255            Types = [wg.GRID_VALUE_BOOL,wg.GRID_VALUE_BOOL,wg.GRID_VALUE_FLOAT+':10,4',
    42494256                wg.GRID_VALUE_FLOAT+':10,4',]+maxVary*[wg.GRID_VALUE_FLOAT+':10,5',]
    42504257            resultsTable = G2gd.Table(resultVals,rowLabels=rowLabels,colLabels=colLabels,types=Types)
     
    42554262            for r in range(resultsGrid.GetNumberRows()):
    42564263                for c in range(resultsGrid.GetNumberCols()):
    4257                     if c == 0:
     4264                    if c in [0,1]:
    42584265                        resultsGrid.SetReadOnly(r,c,isReadOnly=False)
    42594266                    else:
     
    43074314        if G2frame.dataFrame.PhaseUserSize is None:
    43084315            mainSizer.FitInside(G2frame.dataFrame)
    4309             Size = mainSizer.GetMinSize()
     4316            Size = mainSizer.Fit()
    43104317            Size[0] += 40
    4311             Size[1] = max(Size[1],290) + 35
     4318            Size[1] = max(Size[1],350) + 35
    43124319            MCSA.SetSize(Size)
    43134320            MCSA.SetScrollbars(10,10,Size[0]/10-4,Size[1]/10-1)
     
    43264333        phaseName = generalData['Name']
    43274334        MCSAdata = data['MCSA']
     4335        saveResult = []
     4336        for result in MCSAdata['Results']:
     4337            if result[1]:       #keep?
     4338                saveResult.append(result)
     4339        MCSAdata['Results'] = saveResult           
    43284340        covData = {}
    43294341        if 'PWDR' in reflName:
     
    43564368        pgbar.SetPosition(wx.Point(screenSize[2]-Size[0]-305,screenSize[1]+5))
    43574369        pgbar.SetSize(Size)
     4370        time1 = time.time()
    43584371        try:
     4372            tsf = 0.
    43594373            for i in range(mcsaControls['Cycles']):
    4360                 MCSAdata['Results'].append(G2mth.mcsaSearch(data,RBdata,reflType,reflData,covData,pgbar))
    4361                 print ' MC/SA runs completed: ',i
     4374                Result,tsum = G2mth.mcsaSearch(data,RBdata,reflType,reflData,covData,pgbar)
     4375                MCSAdata['Results'].append(Result)
     4376                print ' MC/SA runs completed: %d residual: %.3f%%'%(i,100*Result[2])
     4377                tsf += tsum
     4378            print ' MC/SA run time: %.2f'%(time.time()-time1)
     4379            print ' Structure factor time: %.2f'%(tsf)
    43624380        finally:
    43634381            pgbar.Destroy()
     4382        MCSAdata['Results'] = G2mth.sortArray(MCSAdata['Results'],2,reverse=False)
    43644383        UpdateMCSA()
    43654384        G2plt.PlotStructure(G2frame,data)
  • trunk/GSASIIplot.py

    r953 r961  
    32473247       
    32483248    def RenderPolyhedra(x,y,z,Faces,color):
     3249        glShadeModel(GL_FLAT)
    32493250        glPushMatrix()
    32503251        glTranslate(x,y,z)
     
    32613262            glEnd()
    32623263        glPopMatrix()
     3264        glShadeModel(GL_SMOOTH)
    32633265
    32643266    def RenderMapPeak(x,y,z,color,den):
  • trunk/GSASIIpwdGUI.py

    r945 r961  
    21162116            data['Form Vol'] = max(10.0,SumElementVolumes())
    21172117            formVol.SetValue('%.2f'%(data['Form Vol']))
    2118             UpdatePDFGrid(G2frame,data)
     2118            wx.CallAfter(UpdatePDFGrid,G2frame,data)
    21192119            auxPlot = ComputePDF(data)
    21202120            G2plt.PlotISFG(G2frame,newPlot=True)       
  • trunk/GSASIIstrMath.py

    r954 r961  
    5858        VRBData = RBData['Vector']
    5959        for i,rbId in enumerate(VRBIds):
    60             for j in range(len(VRBData[rbId]['VectMag'])):
    61                 name = '::RBV;'+str(j)+':'+str(i)
    62                 VRBData[rbId]['VectMag'][j] = parmDict[name]
     60            if VRBData[rbId]['useCount']:
     61                for j in range(len(VRBData[rbId]['VectMag'])):
     62                    name = '::RBV;'+str(j)+':'+str(i)
     63                    VRBData[rbId]['VectMag'][j] = parmDict[name]
    6364    for phase in Phases:
    6465        Phase = Phases[phase]
Note: See TracChangeset for help on using the changeset viewer.