Changeset 903


Ignore:
Timestamp:
May 13, 2013 3:18:21 PM (10 years ago)
Author:
toby
Message:

Fix problem with large constraints; Change binding of data item menus so they can be used from any window (Mac only); Start on AUI notebook (not in use); fix periodic table on Mac; capture error if argument file is not found; start on sphinx documentation formatting

Location:
trunk
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • trunk/GSASII.py

    r891 r903  
    99# $Id$
    1010########### SVN repository information ###################
     11'''
     12*GSAS-II Main Module*
     13=====================
     14
     15Main routines for the GSAS-II program
     16'''
    1117
    1218import os
    … …  
    101107
    102108class GSASII(wx.Frame):
    103    
     109    '''Define the main GSAS-II frame and its associated menu items
     110    '''
    104111    def _Add_FileMenuItems(self, parent):
    105112        item = parent.Append(
    … …  
    11011108            self.dirname = os.path.dirname(arg[1])
    11021109            if self.dirname: os.chdir(self.dirname)
    1103             G2IO.ProjFileOpen(self)
    1104             self.PatternTree.Expand(self.root)
    1105             for item in self.Refine: item.Enable(True)
    1106             for item in self.SeqRefine: item.Enable(True)
     1110            try:
     1111                G2IO.ProjFileOpen(self)
     1112                self.PatternTree.Expand(self.root)
     1113                for item in self.Refine: item.Enable(True)
     1114                for item in self.SeqRefine: item.Enable(True)
     1115            except:
     1116                print 'Error opening file',arg[1]
    11071117
    11081118    def OnSize(self,event):
    … …  
    23272337
    23282338class GSASIImain(wx.App):
     2339    '''Defines a wxApp for GSAS-II
     2340
     2341    Creates a wx frame (self.main) which contains the display of the
     2342    data tree.
     2343    '''
    23292344    def OnInit(self):
     2345        '''Called automatically when the app is created.'''
    23302346        self.main = GSASII(None)
    23312347        self.main.Show()
    … …  
    23342350
    23352351def main():
     2352    '''Start up the GSAS-II application'''
    23362353    application = GSASIImain(0)
    23372354    if wxInspector: wxeye.InspectionTool().Show()
    … …  
    23412358   
    23422359if __name__ == '__main__':
    2343     main()
     2360    main() # start the GUI

  • trunk/GSASIIElemGUI.py

    r762 r903  
    2121              style=wx.DEFAULT_DIALOG_STYLE, title='Pick Element')
    2222        import ElementTable as ET
    23         self.butWid = 55
     23        self.butWid = 60
    2424        if 'nt' in os.name:
    2525            self.butWid = 40

  • trunk/GSASIIgrid.py

    r899 r903  
    1111import wx.grid as wg
    1212import wx.wizard as wz
     13import wx.aui
    1314import time
    1415import cPickle
    … …  
    766767    '''Create the dataframe window and its menus
    767768    '''
     769    def Bind(self,*args,**kwargs):
     770        '''Override the Bind() function on the Mac the binding is to
     771        the main window, so that menus operate with any window on top.
     772        For other platforms, call the default class Bind()
     773        '''
     774        if sys.platform == "darwin": # mac
     775            self.G2frame.Bind(*args,**kwargs)
     776        else:
     777            wx.Frame.Bind(self,*args,**kwargs)     
     778       
    768779    def PrefillDataMenu(self,menu,helpType,helpLbl=None,empty=False):
    769780        '''Create the "standard" part of data frame menus. Note that on Linux and
    … …  
    13371348            if self.GetPageText(page) == name:
    13381349                return page
     1350class GSauiNoteBook(wx.aui.AuiNotebook):
     1351    '''Notebook implemented with wx.aui extension'''
     1352    def __init__(self, parent, name='',size = None):
     1353        wx.aui.AuiNotebook.__init__(self, parent, -1)
     1354        #if size: self.SetSize(size)
     1355                                                     
     1356    def Clear(self):       
     1357        GSNoteBook.DeleteAllPages(self)
     1358       
     1359    def FindPage(self,name):
     1360        numPage = self.GetPageCount()
     1361        for page in range(numPage):
     1362            if self.GetPageText(page) == name:
     1363                return page
     1364
    13391365       
    13401366################################################################################

  • trunk/GSASIImapvars.py

    r762 r903  
    77# $Id$
    88########### SVN repository information ###################
    9 """Module that implements algebraic contraints, parameter redefinition
     9"""
     10*GSASIImapvars: Parameter constraints*
     11======================================
     12
     13Module to implements algebraic contraints, parameter redefinition
    1014and parameter simplification contraints.
    1115
    1216Parameter redefinition (new vars) is done by creating one or more relationships
    1317between a set of parameters
     18
     19::
     20
    1421   Mx1 * Px + My1 * Py +...
    1522   Mx2 * Px + Mz2 * Pz + ...
     23
    1624where Pj is a parameter name and Mjk is a constant.
    1725
    1826Constant constraint Relations can also be supplied in the form of an equation:
     27
     28::
     29
    1930  nx1 * Px + ny1 * Py +... = C1
     31
    2032where Cn is a constant. These equations define an algebraic
    2133constrant.
    … …  
    2335Parameters can also be "fixed" (held), which prevents them from being refined.
    2436
    25 All of the above three cases is supplied as input using routines
     37All of the above three cases are input using routines
    2638GroupConstraints and GenerateConstraints. The input consists of a list of
    2739relationship dictionaries:
     40
     41.. code-block:: python
     42
    2843    constrDict = [
    2944         {'0:12:Scale': 2.0, '0:14:Scale': 4.0, '0:13:Scale': 3.0, '0:0:Scale': 0.5},
    3045         {'2::C(10,6,1)': 1.0, '1::C(10,6,1)': 1.0},
    3146         {'0::A0': 0.0}]
    32 
    3347    fixedList = ['5.0', None, '0']
    3448
    … …  
    4155(independent) parameter is equated to several (now dependent) parameters. In
    4256algebraic form this is:
     57
     58::
     59
    4360   P0 = M1 * P1 = M2 * P2 = ...
     61
    4462Thus parameters P0, P1 and P2,... are linearly equivalent. Routine StoreEquivalence is
    4563used to specify these equivalences.
    … …  
    6886varied (appear in varyList) or none can be varied. This is checked in GenerateConstraints
    6987(as well as for generated relationships in SetVaryFlags).
     88
    7089* If all parameters in a parameter redefinition (new var) relationship are varied, the
    7190  parameter assigned to this expression (::constr:n, see paramPrefix) newly generated
    7291  parameter is varied. Note that any generated "missing" relations are not varied. Only
    7392  the input relations are varied.
     93 
    7494* If all parameters in a fixed constraint equation are varied, the generated "missing"
    7595  relations in the group are all varied. This provides the N-C degrees of freedom.
    7696
    77 External Routines:
    78    To define a set of constrained and unconstrained relations, one
    79      calls InitVars, GroupConstraints and GenerateConstraints. It
    80      is best to supply all relations/equations in one call to these
    81      routines so that grouping is done correctly.
    82 
    83    To implement parameter redefinition, one calls
    84      StoreEquivalence. This should be called for every set of
    85      equivalence relationships. There is no harm in using
    86      StoreEquivalence with the same independent variable:
     97*External Routines*
     98-------------------
     99
     100To define a set of constrained and unconstrained relations, one
     101defines a list of dictionary defining constraint parameters and their
     102values, a list of fixed values for each constraint and a list of
     103parameters to be varied. In addition, one uses
     104:func:`StoreEquivalence` to define parameters that are equivalent. One
     105can then use :func:`CheckConstraints` to check that the input is
     106internally consistent and finally :func:`GroupConstraints` and
     107:func:`GenerateConstraints` to generate the internally used
     108tables. Routines :func:`Map2Dict` is used to initialize the parameter
     109dictionary and :func:`Dict2Map`, :func:`Dict2Deriv`, and
     110:func:`ComputeDepESD` are used to apply constraints. Routine
     111:func:`VarRemapShow` is used to print out the constraint information,
     112as stored by :func:`GenerateConstraints`.
     113
     114:func:`InitVars`
     115  This is optionally used to clear out all defined previously defined constraint information
     116 
     117:func:`StoreEquivalence`
     118  To implement parameter redefinition, one calls StoreEquivalence. This should be called for every set of
     119  equivalence relationships. There is no harm in using StoreEquivalence with the same independent variable:
     120
     121  .. code-block:: python
     122
    87123       StoreEquivalence('x',('y',))
    88124       StoreEquivalence('x',('z',))
    89      (though StoreEquivalence('x',('y','z')) will run more
    90      efficiently) but mixing independent and dependent variables is
    91      problematic. This is not allowed:
     125
     126  or equivalently
     127
     128  .. code-block:: python
     129
     130       StoreEquivalence('x',('y','z'))
     131
     132  The latter will run more efficiently. Note that mixing independent and dependent variables is
     133  problematic. This is not allowed:
     134
     135  .. code-block:: python
     136
    92137        StoreEquivalence('x',('y',))
    93138        StoreEquivalence('y',('z',))
    94    Use StoreEquivalence before calling GenerateConstraints or
    95       CheckConstraints
    96 
     139       
     140  Use StoreEquivalence before calling GenerateConstraints or CheckConstraints
     141
     142:func:`CheckConstraints`
    97143   To check that input in internally consistent, use CheckConstraints
    98144
    99    To show a summary of the parameter remapping, one calls
    100       VarRemapShow
    101 
     145:func:`Map2Dict`
    102146   To determine values for the parameters created in this module, one
    103       calls Map2Dict. This will not apply contraints.
    104 
     147   calls Map2Dict. This will not apply contraints.
     148
     149:func:`Dict2Map`
    105150   To take values from the new independent parameters and constraints,
    106       one calls Dict2Map. This will apply contraints.
    107 
     151   one calls Dict2Map. This will apply contraints.
     152
     153:func:`Dict2Deriv`
    108154   Use Dict2Deriv to determine derivatives on independent parameters
    109       from those on dependent ones
    110      
     155   from those on dependent ones
     156
     157:func:`ComputeDepESD`     
    111158   Use ComputeDepESD to compute uncertainties on dependent variables
    112159
    113 Global Variables:
    114    dependentParmList: contains a list by group of lists of
    115      parameters used in the group. Note that parameters listed in
    116      dependentParmList should not be refined as they will not affect
    117      the model
    118    indParmList: a list of groups of Independent parameters defined in
     160:func:`VarRemapShow`
     161   To show a summary of the parameter remapping, one calls VarRemapShow
     162
     163*Global Variables*
     164------------------
     165
     166dependentParmList:
     167   contains a list by group of lists of
     168   parameters used in the group. Note that parameters listed in
     169   dependentParmList should not be refined as they will not affect
     170   the model
     171
     172indParmList:
     173     a list of groups of Independent parameters defined in
    119174     each group. This contains both parameters used in parameter
    120175     redefinitions as well as names of generated new parameters.
    121176
    122    fixedVarList: a list of variables that have been 'fixed'
     177fixedVarList:
     178     a list of variables that have been 'fixed'
    123179     by defining them as equal to a constant (::var: = 0). Note that
    124180     the constant value is ignored at present. These variables are
    125181     later removed from varyList which prevents them from being refined.
    126182     Unlikely to be used externally.
    127    arrayList: a list by group of relationship matrices to relate
     183
     184arrayList:
     185     a list by group of relationship matrices to relate
    128186     parameters in dependentParmList to those in indParmList. Unlikely
    129187     to be used externally.
    130    invarrayList: a list by group of relationship matrices to relate
     188
     189invarrayList:
     190     a list by group of relationship matrices to relate
    131191     parameters in indParmList to those in dependentParmList. Unlikely
    132192     to be used externally.
    133    fixedDict: a dictionary containing the fixed values corresponding
     193
     194fixedDict:
     195     a dictionary containing the fixed values corresponding
    134196     to parameter equations.  The dict key is an ascii string, but the
    135197     dict value is a float.  Unlikely to be used externally.
     198
     199*Routines*
     200----------
     201
     202Note that parameter names in GSAS-II are strings of form ``<ph>:<hst>:<nam>``
     203
    136204"""
    137205
    … …  
    167235def GroupConstraints(constrDict):
    168236    """divide the constraints into groups that share no parameters.
    169     Input
    170        constrDict: a list of dicts defining relationships/constraints
     237
     238    :param dict constrDict: a list of dicts defining relationships/constraints
     239
     240    ::
     241   
    171242       constrDict = [{<constr1>}, {<constr2>}, ...]
    172        where {<constr1>} is {'param1': mult1, 'param2': mult2,...}
    173     Returns two lists of lists:
    174       a list of relationship list indices for each group pointing to
    175         elements in constrDict and
    176       a list containing the parameters used in each group
     243
     244    where {<constr1>} is {'param1': mult1, 'param2': mult2,...}
     245
     246    :returns: two lists of lists:
     247   
     248      * a list of grouped contraints where each constraint grouped containts a list of indices for constraint constrDict entries
     249      * a list containing lists of parameter names contained in each group
     250     
    177251      """
    178252    assignedlist = [] # relationships that have been used
    … …  
    205279    constraints and checks for inconsistencies such as conflicts in
    206280    parameter/variable definitions and or inconsistently varied parameters.
    207     Input: see GenerateConstraints
    208     Output: returns two strings
    209       the first lists conflicts internal to the specified constraints
    210       the second lists conflicts where the varyList specifies some
     281
     282    :param list varyList: a list of parameters names that will be varied
     283
     284    :param dict constrDict: a list of dicts defining relationships/constraints (as defined in :func:`GroupConstraints`)
     285
     286    :param list fixedList: a list of values specifying a fixed value for each dict in constrDict. Values are
     287      either strings that can be converted to floats or None if the constraint defines a new parameter rather
     288      than a constant.
     289
     290    :returns: two strings:
     291
     292      * the first lists conflicts internal to the specified constraints
     293      * the second lists conflicts where the varyList specifies some
    211294        parameters in a constraint, but not all
     295       
    212296      If there are no errors, both strings will be empty
    213297    '''
    … …  
    267351                errmsg += str(mv) + " => " + s + '\n'
    268352
    269     #print 'indepVarList',indepVarList
    270     #print 'depVarList',depVarList
    271353    # check for errors:
    272354    inboth = set(indepVarList).intersection(set(depVarList))
    … …  
    313395                if var in fixVlist:
    314396                    errmsg += '\nParameter '+var+" is Fixed and used in a constraint:\n\t"
    315                     errmsg += FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
    316 #                if var in equivVarList:
    317 #                    errmsg += '\nParameter '+var+" is Equivalenced and used in a constraint:\n\t"
    318 #                    errmsg += FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
     397                    errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
    319398            if varied > 0 and varied != len(constrDict[rel]):
    320399                warnmsg += "\nNot all variables refined in constraint:\n\t"
    321                 warnmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     400                warnmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
    322401                warnmsg += '\nNot refined: ' + notvaried + '\n'
    323     if errmsg or warnmsg: return errmsg,warnmsg
     402    if errmsg or warnmsg:
     403        return errmsg,warnmsg
    324404
    325405    # now look for process each group and create the relations that are needed to form
    326406    # non-singular square matrix
    327407    for group,varlist in zip(groups,parmlist):
    328         if len(varlist) == 1: continue
    329         try:
    330             arr = MakeArray(constrDict, group, varlist)
     408        if len(varlist) == 1: continue # a constraint group with a single variable can be ignored
     409        if len(varlist) < len(group): # too many relationships -- no can do
     410            errmsg += "\nOver-constrained input. "
     411            errmsg += "There are more constraints " + str(len(group))
     412            errmsg += "\n\tthan variables " + str(len(varlist)) + "\n"
     413            for rel in group:
     414                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
     415                errmsg += "\n"
     416                continue
     417        try:
     418            multarr = _FillArray(group,constrDict,varlist)
     419            _RowEchelon(len(group),multarr,varlist)
    331420        except:
    332             errmsg += "\nOver-constrained input. "
    333             errmsg += "There are more constraints" + str(len(group))
    334             errmsg += "\n\tthan variables" + str(len(varlist)) + "\n"
     421            errmsg += "\nSingular input. "
     422            errmsg += "There are internal inconsistencies in these constraints\n"
    335423            for rel in group:
    336                 errmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     424                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
    337425                errmsg += "\n"
     426            continue
    338427        try:
    339             constrArr = FillInMissingRelations(arr,len(group))
    340         except Exception,msg:
    341             if msg.message.startswith('VectorProduct'):
    342                 errmsg += "\nSingular input. "
    343                 errmsg += "This set of constraints is not linearly independent:\n\n"
    344             else:
    345                 errmsg += "\nInconsistent input. "
    346                 errmsg += "Cound not generate a set of non-singular constraints"
    347                 errmsg += "\n\tfor this set of input:\n\n"
     428            multarr = _FillArray(group,constrDict,varlist,FillDiagonals=True)
     429            GramSchmidtOrtho(multarr,len(group))
     430        except:
     431            errmsg += "\nUnexpected singularity with constraints (in Gram-Schmidt)\n"
    348432            for rel in group:
    349                 errmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     433                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
    350434                errmsg += "\n"
    351             #import traceback
    352             #print traceback.format_exc()
    353435            continue
    354        
    355436        mapvar = []
    356437        group = group[:]
    … …  
    386467                warnmsg += '\nPlease report this unexpected error\n'
    387468                for rel in group:
    388                     warnmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     469                    warnmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
    389470                    warnmsg += "\n"
    390471                warnmsg += '\n\tNot refined: ' + notvaried + '\n'
    391472        try:
    392             np.linalg.inv(constrArr)
     473            np.linalg.inv(multarr)           
    393474        except:
    394475            errmsg += "\nSingular input. "
    … …  
    398479            errmsg += 'Please report this unexpected error\n'
    399480            for rel in group:
    400                 errmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     481                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
    401482                errmsg += "\n"
    402483    return errmsg,warnmsg
    … …  
    407488    and stores them in global variables Also checks for internal
    408489    conflicts or inconsistencies in parameter/variable definitions.
    409     Input:
    410        groups,parmlist: see GroupConstraints
    411        varyList: a list of parameters that will be varied
    412        constrDict: a list of dicts defining relationships/constraints
    413          constrDict = [{<constr1>}, {<constr2>}, ...]
    414          where {<constr1>} is {'param1': mult1, 'param2': mult2,...}
    415        fixedList: a list of values for each constraint/variable in
    416           constrDict, value is either a float (contraint) or None (for
    417           a new variable).
     490
     491    :param list groups: a list of grouped contraints where each constraint grouped containts a list of
     492      indices for constraint constrDict entries, created in :func:`GroupConstraints` (returned as 1st value)
     493
     494    :param list parmlist: a list containing lists of parameter names contained in each group,
     495      created in :func:`GroupConstraints` (returned as 1st value)
     496
     497    :param list varyList: a list of parameters names (strings of form ``<ph>:<hst>:<nam>``) that will be varied
     498   
     499    :param dict constrDict: a list of dicts defining relationships/constraints (as defined in :func:`GroupConstraints`)
     500
     501    :param list fixedList: a list of values specifying a fixed value for each dict in constrDict. Values are
     502      either strings that can be converted to floats, float values or None if the constraint defines a new
     503      parameter
     504     
     505    :param dict constrDict: a list of dicts defining relationships/constraints
     506
    418507    '''
    419508    global dependentParmList,arrayList,invarrayList,indParmList,consNum
    … …  
    424513        if len(cdict) == 1:
    425514            fixedVarList.append(cdict.keys()[0])
    426     #print 'fixedVarList',fixedVarList
    427515   
    428516    # process equivalences: make a list of dependent and independent vars
    … …  
    525613                if var in fixedVarList:
    526614                    msg += '\nError: parameter '+var+" is Fixed and used in a constraint:\n\t"
    527                     msg += FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
     615                    msg += _FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
    528616#                if var in equivVarList:
    529617#                    msg += '\nError: parameter '+var+" is Equivalenced and used in a constraint:\n\t"
    530 #                    msg += FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
     618#                    msg += _FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
    531619            if varied > 0 and varied != len(constrDict[rel]):
    532620                msg += "\nNot all variables refined in constraint:\n\t"
    533                 msg += FormatConstraint(constrDict[rel],fixedList[rel])
     621                msg += _FormatConstraint(constrDict[rel],fixedList[rel])
    534622                msg += '\nNot refined: ' + notvaried + '\n'
    535623            if fixedList[rel] is not None and varied > 0:
    … …  
    545633    for group,varlist in zip(groups,parmlist):
    546634        if len(varlist) == 1: continue
    547         arr = MakeArray(constrDict, group, varlist)
    548         constrArr = FillInMissingRelations(arr,len(group))
     635        if len(varlist) < len(group): # too many relationships -- no can do
     636            msg = 'too many relationships'
     637        try:
     638            arr = _FillArray(group,constrDict,varlist)
     639            _RowEchelon(len(group),arr,varlist)
     640            constrArr = _FillArray(group,constrDict,varlist,FillDiagonals=True)
     641            GramSchmidtOrtho(constrArr,len(group))
     642        except:
     643            msg = 'Singular relationships'
     644
    549645        mapvar = []
    550646        group = group[:]
    … …  
    613709    '''Takes a list of dependent parameter(s) and stores their
    614710    relationship to a single independent parameter (independentVar)
    615     input:
    616        independentVar -- name of parameter that will set others (may
    617          be varied)
    618        dependentList -- list of parameter that will set from
    619          independentVar (may not be varied). Each item can be a parameter
     711
     712    :param str independentVar: name of master parameter that will be used to determine the value
     713      to set the dependent variables
     714
     715    :param list dependentList: a list of parameters that will set from
     716         independentVar. Each item in the list can be a string with the parameter
    620717         name or a tuple containing a name and multiplier:
    621          ('parm1',('parm2',.5),)
     718         ``['parm1',('parm2',.5),]``
     719
    622720    '''
    623721   
    … …  
    643741def GetDependentVars():
    644742    '''Return a list of dependent variables: e.g. variables that are
    645     constrained in terms of other variables'''
     743    constrained in terms of other variables
     744
     745    :returns: a list of variable names
     746
     747    '''
    646748    dependentVars = []
    647749    global dependentParmList
    … …  
    652754def GetIndependentVars():
    653755    '''Return a list of independent variables: e.g. variables that are
    654     created by constraints of other variables'''
     756    created by constraints of other variables
     757
     758    :returns: a list of variable names
     759
     760    '''
    655761    independentVars = []
    656762    global indParmList,fixedDict
    … …  
    717823    return sigmaDict
    718824
    719 def FormatConstraint(RelDict,RelVal):
     825def _FormatConstraint(RelDict,RelVal):
    720826    '''Formats a Constraint or Function for use in a convenient way'''
    721827    linelen = 45
    … …  
    742848
    743849def VarRemapShow(varyList):
    744     '''List out the saved relationships.
    745     Returns a string containing the details.
     850    '''List out the saved relationships. This should be done after the constraints have been
     851    defined using :func:`StoreEquivalence`, :func:`GroupConstraints` and :func:`GenerateConstraints`.
     852
     853    :returns: a string containing the details of the contraint relationships
    746854    '''
    747855    s = ''
    … …  
    792900    '''Compute derivatives for Independent Parameters from the
    793901    derivatives for the original parameters
     902
     903    :param list varyList: a list of parameters names that will be varied
     904
     905    :param dict derivDict: a dict containing derivatives for parameter values keyed by the
     906      parameter names.
     907
     908    :param dict dMdv: a dict containing derivatives for dependent parameter computed from
     909      derivDict
     910
    794911    '''
    795912    global dependentParmList,arrayList,invarrayList,indParmList,invarrayList
    … …  
    816933    Removes dependent variables from the varyList
    817934
    818     This should be done once, before any variable refinement is done
    819     to complete the parameter dictionary with the Independent Parameters
     935    This should be done once, after the constraints have been
     936    defined using :func:`StoreEquivalence`,
     937    :func:`GroupConstraints` and :func:`GenerateConstraints` and
     938    before any variable refinement is done
     939    to complete the parameter dictionary by defining independent
     940    parameters and satisfying the constraint equations.
     941
     942    :param dict parmDict: a dict containing parameter values keyed by the
     943      parameter names.
     944      This will contain updated values for both dependent and independent
     945      parameters after Dict2Map is called. It will also contain some
     946      unexpected entries of every constant value {'0':0.0} & {'1.0':1.0},
     947      which do not cause any problems.
     948
     949    :param list varyList: a list of parameters names that will be varied
     950   
     951
    820952    '''
    821953    # process the Independent vars: remove dependent ones from varylist
    … …  
    844976
    845977def Dict2Map(parmDict,varyList):
    846     '''Convert the remapped values back to the original parameters
     978    '''Applies the constraints defined using :func:`StoreEquivalence`,
     979    :func:`GroupConstraints` and :func:`GenerateConstraints` by changing
     980    values in a dict containing the parameters. This should be
     981    done before the parameters are used for any computations
     982
     983    :param dict parmDict: a dict containing parameter values keyed by the
     984      parameter names.
     985      This will contain updated values for both dependent and independent
     986      parameters after Dict2Map is called. It will also contain some
     987      unexpected entries of every constant value {'0':0.0} & {'1.0':1.0},
     988      which do not cause any problems.
     989
     990    :param list varyList: a list of parameters names that will be varied
    847991   
    848     This should be done to apply constraints to parameter values (use
    849     Map2Dict once first). It should be done as part of the Model function
    850     evaluation, before any computation is done
    851     '''
    852     #I think this needs fixing to update parmDict with new values
    853     #   from the constraints based on what is in varyList - RVD. Don't think so -- BHT
     992    '''
    854993    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict
    855994    # reset fixed values (should not be needed, but very quick)
    … …  
    8671006# internal routines follow (these routines are unlikely to be called
    8681007# from outside the module)
    869 def FillInMissingRelations(arrayin,nvars):
    870     '''Fill in unspecified rows in array with non-colinear unit vectors
    871     arrayin is a square array, where only the first nvars rows are defined.
    872    
    873     Returns a new array where all rows are defined
    874 
    875     Throws an exception if there is no non-signular result
    876     (likely because two or more input rows are co-linear)
    877     '''
    878     a = arrayin.copy()
    879     n = nvars
    880     # fill in the matrix with basis vectors not colinear with any of the starting set
    881     for i in range(n,len(a)):
    882         try:
    883             a[n] = VectorProduct(a[:n])
    884         except:
    885             raise Exception,"VectorProduct failed, singular input?"
    886         n += 1
    887     # use the G-S algorithm to compute a complete set of unit vectors orthogonal
    888     # to the first in array
    889     gs = GramSchmidtOrtho(a)
    890     n = nvars
    891     # now replace the created vectors with the ones least correlated to the
    892     # initial set
    893     vlist = [v for v in gs[1:]] # drop the first row
    894     for j in range(n,len(a)):
    895         minavg = None
    896         droplist = []
    897         for k in range(len(vlist)):
    898             v = vlist[k]
    899             avgcor = 0
    900             for i in range(n):
    901                 cor = np.dot(a[i],vlist[k])/(np.linalg.norm(a[i]) * np.linalg.norm(vlist[k]) )
    902                 if np.allclose(cor, 1.0):
    903                     droplist.append(k) # found a vector co-linear to one we already have
    904                                        # don't need it again
    905                     #vlist.pop(k)
    906                     break
    907                 avgcor += cor
    908             else:
    909                 if minavg == None:
    910                     minavg = abs(avgcor/n)
    911                     minnum = k
    912                 elif abs(avgcor/n) < minavg:
    913                     minavg = abs(avgcor/n)
    914                     minnum = k
    915         if minavg == None:
    916             raise Exception,"Failed to find a non-colinear G-S vector for row %d. Should not happen!" % n
    917         a[j] = vlist[minnum]
    918         droplist.append(minnum) # don't need to use this vector again
    919         for i in sorted(droplist,reverse=True):
    920             vlist.pop(i)
    921         n += 1
    922     return a
    923 
    924 
    925 def MakeArray(constrDict, group, varlist):
    926     """Convert the multipliers in a constraint group to an array of
    927     coefficients and place in a square numpy array, adding empty rows as
    928     needed.
    929 
    930     Throws an exception if some sanity checks fail.
    931     """
    932     # do some error checks
    933     if len(varlist) < len(group): # too many relationships -- no can do
    934         raise Exception,'The number of relationships (%d) exceeds the number of parameters (%d):\n\t%s\n\t%s'% (
    935             len(group),len(varlist),group,varlist)
    936     # put the coefficients into an array
    937     multarr = np.zeros(2*[len(varlist),])
    938     i = 0
    939     for cnum in group:
    940         cdict = constrDict[cnum]
    941         j = 0
    942         for var in varlist:
    943             m = cdict.get(var)
    944             if m != None:
    945                 multarr[i,j] = m
    946             j += 1
    947         i += 1
    948     return multarr
    949 
    950 def GramSchmidtOrtho(arrayin):
     1008
     1009def GramSchmidtOrtho(a,nkeep=0):
    9511010    '''Use the Gram-Schmidt process (http://en.wikipedia.org/wiki/Gram-Schmidt) to
    9521011    find orthonormal unit vectors relative to first row.
     1012
     1013    If nkeep is non-zero, the first nkeep rows in the array are not changed
     1014   
    9531015    input:
    954        arrayin: a 2-D non-signular square array
     1016       arrayin: a 2-D non-singular square array
    9551017    returns:
    9561018       a orthonormal set of unit vectors as a square array
    … …  
    9591021        'Projection operator'
    9601022        return a*(np.dot(a,b)/np.dot(a,a))
    961     a = arrayin.copy()
    962     for j in range (len(a)):
     1023    for j in range(nkeep,len(a)):
    9631024        for i in range(j):
    964             a[j] = a[j] - proj(a[i],a[j])
     1025            a[j] -= proj(a[i],a[j])
    9651026        if np.allclose(np.linalg.norm(a[j]),0.0):
    9661027            raise Exception,"Singular input to GramSchmidtOrtho"
    967         a[j] = a[j]/np.linalg.norm(a[j])
     1028        a[j] /= np.linalg.norm(a[j])
    9681029    return a
    9691030
    970 def VectorProduct(vl):
    971     '''Evaluate the n-dimensional "cross" product of the list of vectors in vl
    972     vl can be any length. The length of each vector is arbitrary, but
    973     all must be the same length. See http://en.wikipedia.org/wiki/Levi-Civita_symbol
    974 
    975     This appears to return a vector perpendicular to the supplied set.
    976 
    977     Throws an exception if a vector can not be obtained because the input set of
    978     vectors is already complete or contains any redundancies.
    979    
    980     Uses LeviCitvitaMult recursively to obtain all permutations of vector elements
    981     '''
    982     i = 0
    983     newvec = []
    984     for e in vl[0]:
    985         i += 1
    986         tl = [([i,],1),]
    987         seq = LeviCitvitaMult(tl ,vl)
    988         tsum = 0
    989         for terms,prod in seq:
    990             tsum += EvalLCterm(terms) * prod
    991         newvec.append(tsum)
    992     if np.allclose(newvec,0.0):
    993         raise Exception,"VectorProduct failed, singular or already complete input"
    994     return newvec
    995 
    996 def LeviCitvitaMult(tin,vl):
    997     '''Recursion formula to compute all permutations of vector element products
    998     The first term in each list is the indices of the Levi-Civita symbol, note
    999     that if an index is repeated, the value is zero, so the element is not included.
    1000     The second term is the product of the vector elements
    1001     '''
    1002     v = vl[0]
    1003     vl1 = vl[1:]       
    1004 
    1005     j = 0
    1006     tl = []
    1007     for e in vl[0]:
    1008         j += 1
    1009         for ind,vals in tin:
    1010             if j in ind: continue
    1011             tl.append((ind+[j,],vals*e))
    1012     if len(vl1):
    1013         return LeviCitvitaMult(tl,vl1)
     1031def _FillArray(sel,dict,collist,FillDiagonals=False):
     1032    '''Construct a n by n matrix (n = len(collist)
     1033    filling in the rows using the relationships defined
     1034    in the dictionaries selected by sel
     1035
     1036    If FillDiagonals is True, diagonal elements in the
     1037    array are set to 1.0
     1038    '''
     1039    n = len(collist)
     1040    if FillDiagonals:
     1041        arr = np.eye(n)
    10141042    else:
    1015         return tl
    1016 
    1017 def EvalLCterm(term):
    1018     '''Evaluate the Levi-Civita symbol Epsilon(i,j,k,...)'''
    1019     p = 1
    1020     for i in range(len(term)):
    1021         for j in range(i+1,len(term)):
    1022             p *= (term[j] - term[i])
    1023             if p == 0: return 0
    1024     return p/abs(p)
    1025 
     1043        arr = np.zeros(2*[n,])
     1044    # fill the top rows
     1045    for i,cnum in enumerate(sel):
     1046        for j,var in enumerate(collist):
     1047            arr[i,j] = dict[cnum].get(var,0)
     1048    return arr
     1049
     1050def _SwapColumns(i,m,v):
     1051    '''Swap columns in matrix m as well as the labels in v
     1052    so that element (i,i) is replaced by the first non-zero element in row i after that element
     1053
     1054    Throws an exception if there are no non-zero elements in that row
     1055    '''
     1056    for j in range(i+1,len(v)):
     1057        if not np.allclose(m[i,j],0):
     1058            m[:,(i,j)] = m[:,(j,i)]
     1059            v[i],v[j] = v[j],v[i]
     1060            return
     1061    else:
     1062        raise Exception,'Singular input'
     1063
     1064def _RowEchelon(m,arr,collist):
     1065    '''Convert the first m rows in Matrix arr to row-echelon form
     1066    exchanging columns in the matrix and collist as needed.
     1067
     1068    throws an exception if the matrix is singular because
     1069    the first m rows are not linearly independent
     1070    '''
     1071    n = len(collist)
     1072    for i in range(m):
     1073        if np.allclose(arr[i,i],0):
     1074            _SwapColumns(i,arr,collist)
     1075        arr[i,:] /= arr[i,i] # normalize row
     1076        # subtract current row from subsequent rows to set values to left of diagonal to 0
     1077        for j in range(i+1,m):
     1078            arr[j,:] -= arr[i,:] * arr[j,i]
    10261079
    10271080if __name__ == "__main__":
    1028     import sys
    10291081    parmdict = {}
    10301082    constrDict = [
    … …  
    10471099                '2::atomy:3', '2::atomz:3',
    10481100                '0:12:Scale', '0:11:Scale', '0:14:Scale', '0:13:Scale', '0:0:Scale']
    1049     varyList = ['0::A0', '0::AUiso:0', '0::Afrac:1', '0::Afrac:2', '0::Afrac:3', '0::Afrac:4', '0::dAx:5', '0::dAy:5', '0::dAz:5', '0::AUiso:5', ':0:Back:0', ':0:Back:1', ':0:Back:2', ':0:Back:3', ':0:Back:4', ':0:Back:5', ':0:Back:6', ':0:Back:7', ':0:Back:8', ':0:Back:9', ':0:Back:10', ':0:Back:11', ':0:U', ':0:V', ':0:W', ':0:X', ':0:Y', ':0:Scale', ':0:DisplaceX', ':0:DisplaceY']
    1050     constrDict = [
    1051         {'0::Afrac:4': 24.0, '0::Afrac:1': 16.0, '0::Afrac:3': 24.0, '0::Afrac:2': 16.0},
    1052         {'0::Afrac:1': 1.0, '0::Afrac:2': 1.0},
    1053         {'0::Afrac:4': 1.0, '0::Afrac:3': 1.0}]
    1054     fixedList = ['40.0', '1.0', '1.0']
     1101#    e,w = CheckConstraints([,
     1102#                     [{'2:0:Scale': 1.0, '5:0:Scale': 1.0, '10:0:Scale': 1.0, '6:0:Scale': 1.0, '9:0:Scale': 1.0, '8:0:Scale': 1.0,# '3:0:Scale': 1.0, '4:0:Scale': 1.0, '7:0:Scale': 1.0, '1:0:Scale': 1.0, '0:0:Scale': 1.0}],
     1103#                     ['1.0'])
     1104#    if e: print 'error=',e
     1105#    if w: print 'error=',w
     1106#    varyList = ['0::A0', '0::AUiso:0', '0::Afrac:1', '0::Afrac:2', '0::Afrac:3', '0::Afrac:4', '0::dAx:5', '0::dAy:5', '0::dAz:5', '0::AUiso:5', ':0:Back:0', ':0:Back:1', ':0:Back:2', ':0:Back:3', ':0:Back:4', ':0:Back:5', ':0:Back:6', ':0:Back:7', ':0:Back:8', ':0:Back:9', ':0:Back:10', ':0:Back:11', ':0:U', ':0:V', ':0:W', ':0:X', ':0:Y', ':0:Scale', ':0:DisplaceX', ':0:DisplaceY']
     1107#    constrDict = [
     1108#        {'0::Afrac:4': 24.0, '0::Afrac:1': 16.0, '0::Afrac:3': 24.0, '0::Afrac:2': 16.0},
     1109#        {'0::Afrac:1': 1.0, '0::Afrac:2': 1.0},
     1110#        {'0::Afrac:4': 1.0, '0::Afrac:3': 1.0}]
     1111#    fixedList = ['40.0', '1.0', '1.0']
    10551112
    10561113    errmsg, warnmsg = CheckConstraints(varylist,constrDict,fixedList)
    … …  
    10901147    for key in sorted(parmdict.keys()):
    10911148        print '  '+key,'\t',before.get(key),'\t',parmdict[key]
    1092     dMdv = len(varylist)*[0]
    1093     deriv = {}
    1094     for i,v in enumerate(parmdict.keys()): deriv[v]=i
    1095     Dict2Deriv(varylist,deriv,dMdv)
     1149#    dMdv = len(varylist)*[0]
     1150#    deriv = {}
     1151#    for i,v in enumerate(parmdict.keys()): deriv[v]=i
     1152#    Dict2Deriv(varylist,deriv,dMdv)
Note: See TracChangeset for help on using the changeset viewer.