source: trunk/GSASIImapvars.py @ 3056

Last change on this file since 3056 was 3056, checked in by toby, 4 years ago

expand error reporting with constraints (no fixes yet to prevent them)

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 52.8 KB
Line 
1# -*- coding: utf-8 -*-
2########### SVN repository information ###################
3# $Date: 2017-09-08 19:17:29 +0000 (Fri, 08 Sep 2017) $
4# $Author: toby $
5# $Revision: 3056 $
6# $URL: trunk/GSASIImapvars.py $
7# $Id: GSASIImapvars.py 3056 2017-09-08 19:17:29Z toby $
8########### SVN repository information ###################
9"""
10*GSASIImapvars: Parameter constraints*
11======================================
12
13Module to implements algebraic contraints, parameter redefinition
14and parameter simplification contraints.
15
16Parameter redefinition (new vars) is done by creating one or more relationships
17between a set of parameters
18
19::
20
21   Mx1 * Px + My1 * Py +...
22   Mx2 * Px + Mz2 * Pz + ...
23
24where Pj is a parameter name and Mjk is a constant.
25
26Constant constraint Relations can also be supplied in the form of an equation:
27
28::
29
30  nx1 * Px + ny1 * Py +... = C1
31
32where Cn is a constant. These equations define an algebraic
33constrant.
34
35Parameters can also be "fixed" (held), which prevents them from being refined.
36
37All of the above three cases are input using routines
38GroupConstraints and GenerateConstraints. The input consists of a list of
39relationship dictionaries:
40
41.. code-block:: python
42
43    constrDict = [
44         {'0:12:Scale': 2.0, '0:14:Scale': 4.0, '0:13:Scale': 3.0, '0:0:Scale': 0.5},
45         {'2::C(10,6,1)': 1.0, '1::C(10,6,1)': 1.0},
46         {'0::A0': 0.0}]
47    fixedList = ['5.0', None, '0']
48
49Where the dictionary defines the first part of an expression and the corresponding fixedList
50item is either None (for parameter redefinition) or the constant value, for a constant
51constraint equation. A dictionary that contains a single term defines a variable that
52will be fixed (held). The multiplier and the fixedList value in this case are ignored.
53
54Parameters can also be equivalenced or "slaved" to another parameter, such that one
55(independent) parameter is equated to several (now dependent) parameters. In
56algebraic form this is:
57
58::
59
60   P0 = M1 * P1 = M2 * P2 = ...
61
62Thus parameters P0, P1 and P2,... are linearly equivalent. Routine StoreEquivalence is
63used to specify these equivalences.
64
65Parameter redefinition (new vars) describes a new, independent, parameter, which
66is defined in terms of dependent parameters that are defined in the
67Model, while fixed constrained relations effectively reduce the complexity
68of the Model by removing a degree of freedom. It is possible for a parameter to appear
69in both a parameter redefinition expression and a fixed constraint equation, but a
70parameter cannot be used a parameter equivalance cannot be used elsewhere (not fixed,
71constrained or redefined). Likewise a fixed parameter cannot be used elsewhere (not
72equivalanced, constrained or redefined).
73
74Relationships are grouped so that a set of dependent parameters appear
75in only one group (done in routine GroupConstraints.) Note that if a
76group contains relationships/equations that involve N dependent
77parameters, there must exist N-C new parameters, where C is the number
78of contraint equations in the group. Routine GenerateConstraints takes
79the output from GroupConstraints and generates the
80"missing" relationships and saves that information in the module's
81global variables. Each generated parameter is named sequentially using paramPrefix.
82
83A list of parameters that will be varied is specified as input to GenerateConstraints
84(varyList). A fixed parameter will simply be removed from this list preventing that
85parameter from being varied. Note that all parameters in a constraint relationship
86must specified as varied (appear in varyList) or none can be varied. This is
87checked in GenerateConstraints. Likewise, if all parameters in a constraint are
88not referenced in a refinement, the constraint is ignored, but if some parameters
89in a constraint group are not referenced in a refinement, but others are this
90constitutes and error.
91
92* When a new variable is created, the variable is assigned the name associated
93  in the constraint definition or it is assigned a default name of form
94  ``::constr<n>`` (see paramPrefix). The vary setting for variables used in the
95  constraint are ignored.
96  Note that any generated "missing" relations are not varied. Only
97  the input relations can be are varied.
98 
99* If all parameters in a fixed constraint equation are varied, the generated "missing"
100  relations in the group are all varied. This provides the N-C degrees of freedom.
101
102*External Routines*
103-------------------
104
105To define a set of constrained and unconstrained relations, one
106defines a list of dictionary defining constraint parameters and their
107values, a list of fixed values for each constraint and a list of
108parameters to be varied. In addition, one uses
109:func:`StoreEquivalence` to define parameters that are equivalent. One
110can then use :func:`CheckConstraints` to check that the input is
111internally consistent and finally :func:`GroupConstraints` and
112:func:`GenerateConstraints` to generate the internally used
113tables. Routines :func:`Map2Dict` is used to initialize the parameter
114dictionary and :func:`Dict2Map`, :func:`Dict2Deriv`, and
115:func:`ComputeDepESD` are used to apply constraints. Routine
116:func:`VarRemapShow` is used to print out the constraint information,
117as stored by :func:`GenerateConstraints`.
118
119:func:`InitVars`
120  This is optionally used to clear out all defined previously defined constraint information
121 
122:func:`StoreEquivalence`
123  To implement parameter redefinition, one calls StoreEquivalence. This should be called for every set of
124  equivalence relationships. There is no harm in using StoreEquivalence with the same independent variable:
125
126  .. code-block:: python
127
128       StoreEquivalence('x',('y',))
129       StoreEquivalence('x',('z',))
130
131  or equivalently
132
133  .. code-block:: python
134
135       StoreEquivalence('x',('y','z'))
136
137  The latter will run more efficiently. Note that mixing independent and dependent variables would
138  require assignments to be done in a particular order and thus is not is not allowed. This will cause an
139  error:
140
141  .. code-block:: python
142
143        StoreEquivalence('x',('y',))
144        StoreEquivalence('y',('z',))
145       
146  Use StoreEquivalence before calling GenerateConstraints or CheckConstraints
147
148:func:`CheckConstraints`
149   To check that input in internally consistent, use CheckConstraints
150
151:func:`Map2Dict`
152   To determine values for the parameters created in this module, one
153   calls Map2Dict. This will not apply contraints.
154
155:func:`Dict2Map`
156   To take values from the new independent parameters and constraints,
157   one calls Dict2Map. This will apply contraints.
158
159:func:`Dict2Deriv`
160   Use Dict2Deriv to determine derivatives on independent parameters
161   from those on dependent ones
162
163:func:`ComputeDepESD`     
164   Use ComputeDepESD to compute uncertainties on dependent variables
165
166:func:`VarRemapShow`
167   To show a summary of the parameter remapping, one calls VarRemapShow
168
169*Global Variables*
170------------------
171
172dependentParmList:
173   contains a list by group of lists of
174   parameters used in the group. Note that parameters listed in
175   dependentParmList should not be refined as they will not affect
176   the model
177
178indParmList:
179     a list of groups of Independent parameters defined in
180     each group. This contains both parameters used in parameter
181     redefinitions as well as names of generated new parameters.
182
183fixedVarList:
184     a list of variables that have been 'fixed'
185     by defining them as equal to a constant (::var: = 0). Note that
186     the constant value is ignored at present. These variables are
187     later removed from varyList which prevents them from being refined.
188     Unlikely to be used externally.
189
190arrayList:
191     a list by group of relationship matrices to relate
192     parameters in dependentParmList to those in indParmList. Unlikely
193     to be used externally.
194
195invarrayList:
196     a list by group of relationship matrices to relate
197     parameters in indParmList to those in dependentParmList. Unlikely
198     to be used externally.
199
200fixedDict:
201     a dictionary containing the fixed values corresponding
202     to parameter equations.  The dict key is an ascii string, but the
203     dict value is a float.  Unlikely to be used externally.
204
205symGenList:
206     a list of boolean values that will be True to indicate that a constraint
207     (only equivalences) is generated by symmetry (or Pawley overlap)
208     
209problemVars:
210     a list containing variables that show up in constraints producing errors
211
212*Routines*
213----------
214
215Note that parameter names in GSAS-II are strings of form ``<ph>:<hst>:<nam>``
216
217"""
218
219import numpy as np
220import sys
221import GSASIIpath
222GSASIIpath.SetVersionNumber("$Revision: 3056 $")
223# data used for constraints;
224debug = False # turns on printing as constraint input is processed
225# note that constraints are stored listed by contraint groups, where each constraint
226# group contains those parameters that must be handled together
227dependentParmList = [] # contains a list of parameters in each group
228# note that parameters listed in dependentParmList should not be refined
229arrayList = [] # a list of of relationship matrices
230invarrayList = [] # a list of inverse relationship matrices
231indParmList = [] # a list of names for the new parameters
232fixedDict = {} # a dictionary containing the fixed values corresponding to defined parameter equations
233               # key is original ascii string, value is float
234fixedVarList = [] # List of variables that should not be refined
235symGenList = [] # Flag if constraint is generated by symmetry
236problemVars = [] # variables causing errors
237
238# prefix for parameter names
239paramPrefix = "::constr"
240consNum = 0 # number of the next constraint to be created
241
242def InitVars():
243    '''Initializes all constraint information'''
244    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict,consNum,symGenList
245    dependentParmList = [] # contains a list of parameters in each group
246    arrayList = [] # a list of of relationship matrices
247    invarrayList = [] # a list of inverse relationship matrices
248    indParmList = [] # a list of names for the new parameters
249    fixedDict = {} # a dictionary containing the fixed values corresponding to defined parameter equations
250    symGenList = [] # Flag if constraint is generated by symmetry
251    consNum = 0 # number of the next constraint to be created
252
253def VarKeys(constr):
254    """Finds the keys in a constraint that represent variables
255    e.g. eliminates any that start with '_'
256
257    :param dict constr: a single constraint entry of form::
258
259        {'var1': mult1, 'var2': mult2,... '_notVar': val,...}
260
261        (see :func:`GroupConstraints`)
262    :returns: a list of keys where any keys beginning with '_' are
263      removed.
264    """
265    return [i for i in constr.keys() if not i.startswith('_')]
266
267
268def GroupConstraints(constrDict):
269    """divide the constraints into groups that share no parameters.
270
271    :param dict constrDict: a list of dicts defining relationships/constraints
272
273    ::
274   
275       constrDict = [{<constr1>}, {<constr2>}, ...]
276
277    where {<constr1>} is {'var1': mult1, 'var2': mult2,... }
278
279    :returns: two lists of lists:
280   
281      * a list of grouped contraints where each constraint grouped containts a list
282        of indices for constraint constrDict entries
283      * a list containing lists of parameter names contained in each group
284     
285      """
286    assignedlist = [] # relationships that have been used
287    groups = [] # contains a list of grouplists
288    ParmList = []
289    for i,consi in enumerate(constrDict):
290        if i in assignedlist: continue # already in a group, skip
291        # starting a new group
292        grouplist = [i,]
293        assignedlist.append(i)
294        groupset = set(VarKeys(consi))
295        changes = True # always loop at least once
296        while(changes): # loop until we can't find anything to add to the current group
297            changes = False # but don't loop again unless we find something
298            for j,consj in enumerate(constrDict):
299                if j in assignedlist: continue # already in a group, skip
300                if len(set(VarKeys(consj)) & groupset) > 0: # true if this needs to be added
301                    changes = True
302                    grouplist.append(j)
303                    assignedlist.append(j)
304                    groupset = groupset | set(VarKeys(consj))
305        group = sorted(grouplist)
306        varlist = sorted(list(groupset))
307        groups.append(group)
308        ParmList.append(varlist)
309    return groups,ParmList
310
311def CheckConstraints(varyList,constrDict,fixedList):
312    '''Takes a list of relationship entries comprising a group of
313    constraints and checks for inconsistencies such as conflicts in
314    parameter/variable definitions and or inconsistently varied parameters.
315
316    :param list varyList: a list of parameters names that will be varied
317
318    :param dict constrDict: a list of dicts defining relationships/constraints
319      (as created in :func:`GSASIIstrIO.ProcessConstraints` and
320      documented in :func:`GroupConstraints`)
321
322    :param list fixedList: a list of values specifying a fixed value for each
323      dict in constrDict. Values are either strings that can be converted to
324      floats or ``None`` if the constraint defines a new parameter rather
325      than a constant.
326
327    :returns: two strings:
328
329      * the first lists conflicts internal to the specified constraints
330      * the second lists conflicts where the varyList specifies some
331        parameters in a constraint, but not all
332       
333      If there are no errors, both strings will be empty
334    '''
335    import re
336    global dependentParmList,arrayList,invarrayList,indParmList,consNum
337    errmsg = ''
338    warnmsg = ''
339    global problemVars
340    problemVars = []
341    fixVlist = []
342    # process fixed variables (holds)
343    for cdict in constrDict:
344        # N.B. No "_" names in holds
345        if len(cdict) == 1:
346            fixVlist.append(cdict.keys()[0])
347   
348    # process equivalences: make a list of dependent and independent vars
349    #    and check for repeated uses (repetition of a parameter as an
350    #    independent var is OK)
351    indepVarList = []
352    depVarList = []
353    multdepVarList = []
354    for varlist,mapvars,multarr,invmultarr in zip(
355        dependentParmList,indParmList,arrayList,invarrayList):
356        if multarr is None: # an equivalence
357            zeromult = False
358            for mv in mapvars:
359                varied = 0
360                notvaried = ''
361                if mv in varyList:
362                    varied += 1
363                else:
364                    if notvaried: notvaried += ', '
365                    notvaried += mv
366                if mv not in indepVarList: indepVarList.append(mv)
367                for v,m in zip(varlist,invmultarr):
368                    if v in indepVarList:
369                        errmsg += '\nVariable '+v+' is used to set values in a constraint before its value is set in another constraint\n'
370                        if v not in problemVars: problemVars.append(v)
371                    if m == 0: zeromult = True
372                    if v in varyList:
373                        varied += 1
374                    else:
375                        if notvaried: notvaried += ', '
376                        notvaried += v
377                    if v in depVarList:
378                        multdepVarList.append(v)
379                    else:
380                        depVarList.append(v)
381            if varied > 0 and varied != len(varlist)+1:
382                warnmsg += "\nNot all variables refined in equivalence:\n\t"
383                s = ""
384                for v in varlist:
385                    if s != "": s+= " & "
386                    s += str(v)           
387                warnmsg += str(mv) + " => " + s
388                warnmsg += '\nNot refined: ' + notvaried + '\n'
389            if zeromult:
390                errmsg += "\nZero multiplier is invalid in equivalence:\n\t"
391                s = ""
392                for v in varlist:
393                    if s != "": s+= " & "
394                    s += str(v)           
395                errmsg += str(mv) + " => " + s + '\n'
396
397    # check for errors:
398    if len(multdepVarList) > 0:
399        errmsg += "\nThe following parameters(s) are used in conflicting Equivalence relations as dependent variables:\n"
400        s = ''
401        for var in sorted(set(multdepVarList)):
402            if v not in problemVars: problemVars.append(v)
403            if s != "": s+= ", "
404            s += str(var)           
405        errmsg += '\t'+ s + '\n'
406    equivVarList = list(set(indepVarList).union(set(depVarList)))
407    if debug: print 'equivVarList',equivVarList
408    inboth = set(fixVlist).intersection(set(equivVarList))
409    if len(inboth) > 0:
410        errmsg += "\nThe following parameter(s) are used in both Equivalence and Fixed constraints:\n"
411        s = ''
412        for var in sorted(inboth):
413            if var not in problemVars: problemVars.append(var)
414            if s != "": s+= ", "
415            s += str(var)
416        errmsg += '\t'+ s + '\n'
417
418    groups,parmlist = GroupConstraints(constrDict)
419    # scan through parameters in each relationship. Are all varied? If only some are
420    # varied, create a warning message.
421    for group,varlist in zip(groups,parmlist):
422        if len(varlist) == 1: continue
423        for rel in group:
424            varied = 0
425            notvaried = ''
426            for var in constrDict[rel]:
427                if var.startswith('_'): continue
428                if not re.match('[0-9]*:[0-9\*]*:',var):
429                    warnmsg += "\nVariable "+str(var)+" does not begin with a ':'"
430                if var in varyList:
431                    varied += 1
432                else:
433                    if notvaried: notvaried += ', '
434                    notvaried += var
435                if var in fixVlist:
436                    errmsg += '\nParameter '+var+" is Fixed and used in a constraint:\n\t"
437                    if var not in problemVars: problemVars.append(var)
438                    errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
439            if varied > 0 and varied != len(VarKeys(constrDict[rel])):
440                warnmsg += "\nNot all variables refined in constraint:\n\t"
441                warnmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
442                warnmsg += '\nNot refined: ' + notvaried + '\n'
443    if errmsg or warnmsg:
444        return errmsg,warnmsg
445
446    # now look for process each group and create the relations that are needed to form
447    # non-singular square matrix
448    for group,varlist in zip(groups,parmlist):
449        if len(varlist) == 1: continue # a constraint group with a single variable can be ignored
450        if len(varlist) < len(group): # too many relationships -- no can do
451            errmsg += "\nOver-constrained input. "
452            errmsg += "There are more constraints " + str(len(group))
453            errmsg += "\n\tthan variables " + str(len(varlist)) + "\n"
454            for rel in group:
455                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
456                errmsg += "\n"
457                continue
458        try:
459            multarr = _FillArray(group,constrDict,varlist)
460            _RowEchelon(len(group),multarr,varlist)
461        except:
462            errmsg += "\nSingular input. "
463            errmsg += "There are internal inconsistencies in these constraints\n"
464            for rel in group:
465                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
466                errmsg += "\n"
467            continue
468        try:
469            multarr = _FillArray(group,constrDict,varlist,FillDiagonals=True)
470            GramSchmidtOrtho(multarr,len(group))
471        except:
472            errmsg += "\nUnexpected singularity with constraints (in Gram-Schmidt)\n"
473            for rel in group:
474                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
475                errmsg += "\n"
476            continue
477        mapvar = []
478        group = group[:]
479        # scan through all generated and input variables
480        # Check again for inconsistent variable use
481        # for new variables -- where varied and unvaried parameters get grouped
482        # together. I don't think this can happen when not flagged before, but
483        # it does not hurt to check again.
484        for i in range(len(varlist)):
485            varied = 0
486            notvaried = ''
487            if len(group) > 0:
488                rel = group.pop(0)
489                fixedval = fixedList[rel]
490                for var in VarKeys(constrDict[rel]):
491                    if var in varyList:
492                        varied += 1
493                    else:
494                        if notvaried: notvaried += ', '
495                        notvaried += var
496            else:
497                fixedval = None
498            if fixedval is None:
499                varname = paramPrefix + str(consNum) # assign a name to a variable
500                mapvar.append(varname)
501                consNum += 1
502            else:
503                mapvar.append(fixedval)
504            if varied > 0 and notvaried != '':
505                warnmsg += "\nNot all variables refined in generated constraint"
506                warnmsg += '\nPlease report this unexpected error\n'
507                for rel in group:
508                    warnmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
509                    warnmsg += "\n"
510                warnmsg += '\n\tNot refined: ' + notvaried + '\n'
511        try:
512            np.linalg.inv(multarr)           
513        except:
514            errmsg += "\nSingular input. "
515            errmsg += "The following constraints are not "
516            errmsg += "linearly independent\n\tor do not "
517            errmsg += "allow for generation of a non-singular set\n"
518            errmsg += 'Please report this unexpected error\n'
519            for rel in group:
520                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
521                errmsg += "\n"
522    return errmsg,warnmsg
523
524def GenerateConstraints(groups,parmlist,varyList,constrDict,fixedList,parmDict=None,SeqHist=None):
525    '''Takes a list of relationship entries comprising a group of
526    constraints and builds the relationship lists and their inverse
527    and stores them in global variables Also checks for internal
528    conflicts or inconsistencies in parameter/variable definitions.
529
530    :param list groups: a list of grouped contraints where each constraint
531      grouped containts a list of indices for constraint constrDict entries,
532      created in :func:`GroupConstraints` (returned as 1st value)
533
534    :param list parmlist: a list containing lists of parameter names
535      contained in each group, created in :func:`GroupConstraints`
536      (returned as 2nd value)
537
538    :param list varyList: a list of parameters names (strings of form
539      ``<ph>:<hst>:<nam>``) that will be varied. Note that this is changed here.
540   
541    :param dict constrDict: a list of dicts defining relationships/constraints
542      (as defined in :func:`GroupConstraints`)
543
544    :param list fixedList: a list of values specifying a fixed value for each
545      dict in constrDict. Values are either strings that can be converted to
546      floats, float values or None if the constraint defines a new parameter.
547     
548    :param dict parmDict: a dict containing all parameters defined in current
549      refinement.
550
551    :param int SeqHist: number of current histogram, when used in a sequential
552      refinement. None (default) otherwise. Wildcard variable names are
553      set to the current histogram, when found if not None.
554    '''
555    global dependentParmList,arrayList,invarrayList,indParmList,consNum
556    msg = ''
557    shortmsg = ''
558
559    # process fixed (held) variables
560    for cdict in constrDict:
561        if len(cdict) == 1:
562            fixedVarList.append(cdict.keys()[0])
563   
564    # process equivalences: make a list of dependent and independent vars
565    #    and check for repeated uses (repetition of a parameter as an
566    #    independent var is OK [A=B; A=C], but chaining: [A=B; B=C] is not good)
567    dropVarList = []
568    translateTable = {} # lookup table for wildcard referenced variables
569    for varlist,mapvars,multarr,invmultarr in zip(       # process equivalences
570        dependentParmList,indParmList,arrayList,invarrayList):
571        if multarr is None: # true only if an equivalence
572            zeromult = False
573            for mv in mapvars:
574                #s = ''
575                varied = 0
576                notvaried = ''
577                if mv in varyList:
578                    varied += 1
579                else:
580                    if notvaried: notvaried += ', '
581                    notvaried += mv
582                if parmDict is not None and mv not in parmDict:
583                    print "Dropping equivalence for variable "+str(mv)+". Not defined in this refinement"
584                    if mv not in dropVarList: dropVarList.append(mv)
585                    #msg += "\nCannot equivalence to variable "+str(mv)+". Not defined in this refinement"
586                    #continue
587            for v,m in zip(varlist,invmultarr):
588                if parmDict is not None and v not in parmDict:
589                    print "Dropping equivalence for dep. variable "+str(v)+". Not defined in this refinement"
590                    if v not in dropVarList: dropVarList.append(v)
591                    continue
592                if m == 0: zeromult = True
593                if v in varyList:
594                    varied += 1
595                else:
596                    if notvaried: notvaried += ', '
597                    notvaried += v
598            if varied > 0 and varied != len(varlist)+1:
599                msg += "\nNot all variables refined in equivalence:\n\t"
600                s = ""
601                for v in varlist:
602                    if s != "": s+= " & "
603                    s += str(v)           
604                msg += str(mv) + " => " + s
605                msg += '\nNot refined: ' + notvaried + '\n'
606            if zeromult:
607                msg += "\nZero multiplier is invalid in equivalence:\n\t"
608                s = ""
609                for v in varlist:
610                    if s != "": s+= " & "
611                    s += str(v)           
612                msg += str(mv) + " => " + s + '\n'
613
614    # scan through parameters in each relationship. Are all varied? If only some are
615    # varied, create an error message.
616    for group,varlist in zip(groups,parmlist):
617        if len(varlist) == 1: continue
618        for rel in group:
619            varied = 0
620            notvaried = ''
621            unused = 0
622            notused = ''
623            for var in constrDict[rel]:
624                if var.startswith('_'): continue
625                if var.split(':')[1] == '*' and SeqHist is not None:
626                    # convert wildcard var to reference current histogram; save translation in table
627                    sv = var.split(':')
628                    sv[1] = str(SeqHist)
629                    translateTable[var] = ':'.join(sv)
630                    var = translateTable[var]
631                if parmDict is not None and var not in parmDict:
632                    unused += 1
633                    if notvaried: notused += ', '
634                    notused += var
635                if var in varyList:
636                    varied += 1
637                else:
638                    if notvaried: notvaried += ', '
639                    notvaried += var
640                if var in fixedVarList:
641                    msg += '\nError: parameter '+var+" is Fixed and used in a constraint:\n\t"
642                    msg += _FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
643            #if unused > 0:# and unused != len(VarKeys(constrDict[rel])):
644            if unused > 0 and unused != len(VarKeys(constrDict[rel])):
645                #msg += "\nSome (but not all) variables in constraint are not defined:\n\t"
646                #msg += _FormatConstraint(constrDict[rel],fixedList[rel])
647                #msg += '\nNot used: ' + notused + '\n'
648                shortmsg += notused+" not used in constraint "+_FormatConstraint(constrDict[rel],fixedList[rel])
649            elif varied > 0 and varied != len(VarKeys(constrDict[rel])):
650                #msg += "\nNot all variables refined in constraint:\n\t"
651                #msg += _FormatConstraint(constrDict[rel],fixedList[rel])
652                #msg += '\nNot refined: ' + notvaried + '\n'
653                shortmsg += notvaried+" not varied in constraint "+_FormatConstraint(constrDict[rel],fixedList[rel])
654    # if there were errors found, go no farther
655    if shortmsg and SeqHist is not None:
656        if msg:
657            print ' *** ERROR in constraint definitions! ***'
658            print msg
659            raise Exception
660        print '*** Sequential refinement: ignoring constraint definition(s): ***'
661        print shortmsg
662        msg = ''
663    elif shortmsg:
664        msg += shortmsg
665    if msg:
666        print ' *** ERROR in constraint definitions! ***'
667        print msg
668        raise Exception
669               
670    # now process each group and create the relations that are needed to form
671    # a non-singular square matrix
672    # If all are varied and this is a constraint equation, then set VaryFree flag
673    # so that the newly created relationships will be varied
674    for group,varlist in zip(groups,parmlist):
675        if len(varlist) == 1: continue
676        # for constraints, if all included variables are refined,
677        # set the VaryFree flag, and remaining degrees of freedom will be
678        # varied (since consistency was checked, if any one variable is
679        # refined, then assume that all are)
680        varsList = [] # make a list of all the referenced variables as well
681        VaryFree = False
682        for rel in group:
683            varied = 0
684            unused = 0
685            for var in VarKeys(constrDict[rel]):
686                var = translateTable.get(var,var) # replace wildcards
687                if parmDict is not None and var not in parmDict:
688                    unused += 1                   
689                if var not in varsList: varsList.append(var)
690                if var in varyList: varied += 1
691            if fixedList[rel] is not None and varied > 0:
692                VaryFree = True
693        if len(varlist) < len(group): # too many relationships -- no can do
694            msg = 'too many relationships'
695            break
696        # Since we checked before, if any variables are unused, then all must be.
697        # If so, this set of relationships can be ignored
698        if unused:
699            if debug: print('Constraint ignored (all variables undefined)')
700            if debug: print ('    '+_FormatConstraint(constrDict[rel],fixedList[rel]))
701            continue
702        # fill in additional degrees of freedom
703        try:
704            arr = _FillArray(group,constrDict,varlist)
705            _RowEchelon(len(group),arr,varlist)
706            constrArr = _FillArray(group,constrDict,varlist,FillDiagonals=True)
707            GramSchmidtOrtho(constrArr,len(group))
708        except:
709            msg = 'Singular relationships'
710            break
711        mapvar = []
712        group = group[:]
713        # scan through all generated and input relationships, we need to add to the varied list
714        # all the new parameters where VaryFree has been set or where a New Var is varied.
715        #
716        # If a group does not contain any fixed values (constraint equations)
717        # and nothing in the group is varied, drop this group, so that the
718        # dependent parameters can be refined individually.
719        unused = True
720        for i in range(len(varlist)):
721            if len(group) > 0: # get the original equation reference
722                rel = group.pop(0)
723                fixedval = fixedList[rel]
724                varyflag = constrDict[rel].get('_vary',False)
725                varname = constrDict[rel].get('_name','')
726            else: # this relationship has been generated
727                varyflag = False
728                varname = ''
729                fixedval = None
730            if fixedval is None: # this is a new variable, not a constraint
731                if not varname:
732                    varname = paramPrefix + str(consNum) # no assigned name, create one
733                    consNum += 1
734                mapvar.append(varname)
735                # vary the new relationship if it is a degree of freedom in
736                # a set of contraint equations or if a New Var is flagged to be varied.
737                if VaryFree or varyflag: 
738                    unused = False
739                    varyList.append(varname)
740                    # fix (prevent varying) of all the variables inside the constraint group
741                    # (dependent vars)
742                    for var in varsList:
743                        if var in varyList: varyList.remove(var)
744            else:
745                unused = False
746                mapvar.append(fixedval)
747        if unused: # skip over constraints that don't matter (w/o fixed value or any refined variables)
748            if debug: print('Constraint ignored (all variables unrefined)')
749            if debug: print ('   '+_FormatConstraint(constrDict[rel],fixedList[rel]))
750            continue 
751        dependentParmList.append([translateTable.get(var,var) for var in varlist])
752        arrayList.append(constrArr)
753        invarrayList.append(np.linalg.inv(constrArr))
754        indParmList.append(mapvar)
755        symGenList.append(False)
756    if msg:
757        print ' *** ERROR in constraint definitions! ***'
758        print msg
759        print VarRemapShow(varyList)
760        raise Exception
761    # setup dictionary containing the fixed values
762    global fixedDict
763    # key is original ascii string, value is float
764    for fixedval in fixedList:
765        if fixedval:
766            fixedDict[fixedval] = float(fixedval)
767
768    # make list of dependent and independent variables (after dropping unused)
769    global dependentVars
770    global independentVars
771    dependentVars = []
772    independentVars = []
773    for varlist,mapvars in zip(dependentParmList,indParmList):  # process all constraints
774        for mv in mapvars:
775            if mv in dropVarList: continue
776            if mv not in independentVars: independentVars.append(mv)
777        for mv in varlist:
778            if mv in dropVarList: continue
779            if mv not in dependentVars: dependentVars.append(mv)
780    if debug: # on debug, show what is parsed & generated, semi-readable
781        print 50*'-'
782        print VarRemapShow(varyList)
783        print 'Varied: ',varyList
784        print 'Not Varied: ',fixedVarList
785
786def StoreEquivalence(independentVar,dependentList,symGen=True):
787    '''Takes a list of dependent parameter(s) and stores their
788    relationship to a single independent parameter (independentVar)
789
790    :param str independentVar: name of master parameter that will be used to determine the value
791      to set the dependent variables
792
793    :param list dependentList: a list of parameters that will set from
794         independentVar. Each item in the list can be a string with the parameter
795         name or a tuple containing a name and multiplier:
796         ``['::parm1',('::parm2',.5),]``
797
798    '''
799   
800    global dependentParmList,arrayList,invarrayList,indParmList,symGenList
801    mapList = []
802    multlist = []
803    for var in dependentList:
804        if isinstance(var, basestring):
805            mult = 1.0
806        elif len(var) == 2:
807            var,mult = var
808        else:
809            raise Exception("Cannot parse "+repr(var) + " as var or (var,multiplier)")
810        mapList.append(var)
811        multlist.append(tuple((mult,)))
812    # added relationships to stored values
813    arrayList.append(None)
814    invarrayList.append(np.array(multlist))
815    indParmList.append(tuple((independentVar,)))
816    dependentParmList.append(mapList)
817    symGenList.append(symGen)
818    return
819
820def GetDependentVars():
821    '''Return a list of dependent variables: e.g. variables that are
822    constrained in terms of other variables
823
824    :returns: a list of variable names
825
826    '''
827    global dependentVars
828    return dependentVars
829
830def GetIndependentVars():
831    '''Return a list of independent variables: e.g. variables that are
832    created by constraints of other variables
833
834    :returns: a list of variable names
835
836    '''
837    global independentVars
838    return independentVars
839
840def PrintIndependentVars(parmDict,varyList,sigDict,PrintAll=False,pFile=None):
841    '''Print the values and uncertainties on the independent variables'''
842    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict
843    printlist = []
844    mapvars = GetIndependentVars()
845    for i,name in enumerate(mapvars):
846        if name in fixedDict: continue
847        if PrintAll or name in varyList:
848            sig = sigDict.get(name)
849            printlist.append([name,parmDict[name],sig])
850    if len(printlist) == 0: return
851    s1 = ''
852    s2 = ''
853    s3 = ''
854    print >>pFile,130*'-'
855    print >>pFile,"Variables generated by constraints"
856    printlist.append(3*[None])
857    for name,val,esd in printlist:
858        if len(s1) > 120 or name is None:
859            print >>pFile,''
860            print >>pFile,s1
861            print >>pFile,s2
862            print >>pFile,s3
863            s1 = ''
864            if name is None: break
865        if s1 == "":
866            s1 = ' name  :'
867            s2 = ' value :'
868            s3 = ' sig   :'
869        s1 += '%15s' % (name)
870        s2 += '%15.4f' % (val)
871        if esd is None:
872            s3 += '%15s' % ('n/a')
873        else:   
874            s3 += '%15.4f' % (esd)
875
876def ComputeDepESD(covMatrix,varyList,parmDict):
877    '''Compute uncertainties for dependent parameters from independent ones
878    returns a dictionary containing the esd values for dependent parameters
879    '''
880    sigmaDict = {}
881    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
882        #if invmultarr is None: continue # probably not needed
883#        try:
884#            valuelist = [parmDict[var] for var in mapvars]
885#        except KeyError:
886#            continue
887        # get the v-covar matrix for independent parameters
888        vcov = np.zeros((len(mapvars),len(mapvars)))
889        for i1,name1 in enumerate(mapvars):
890            if name1 not in varyList: continue
891            iv1 = varyList.index(name1)
892            for i2,name2 in enumerate(mapvars):
893                if name2 not in varyList: continue
894                iv2 = varyList.index(name2)
895                vcov[i1][i2] = covMatrix[iv1][iv2]
896        # vec is the vector that multiplies each of the independent values
897        for v,vec in zip(varlist,invmultarr):
898            sigmaDict[v] = np.sqrt(np.inner(vec.T,np.inner(vcov,vec)))
899    return sigmaDict
900
901def _FormatConstraint(RelDict,RelVal):
902    '''Formats a Constraint or Function for use in a convenient way'''
903    linelen = 45
904    s = [""]
905    for var,val in RelDict.items():
906        if var.startswith('_'): continue
907        if len(s[-1]) > linelen: s.append(' ')
908        m = val
909        if s[-1] != "" and m >= 0:
910            s[-1] += ' + '
911        elif s[-1] != "":
912            s[-1] += ' - '
913            m = abs(m)
914        s[-1] += '%.3f*%s '%(m,var)
915    if len(s[-1]) > linelen: s.append(' ')
916    if RelVal is None:
917        s[-1] += ' = New variable'
918    else:
919        s[-1] += ' = ' + RelVal
920    s1 = ''
921    for s2 in s:
922        if s1 != '': s1 += '\n\t'
923        s1 += s2
924    return s1
925
926def VarRemapShow(varyList,inputOnly=False):
927    '''List out the saved relationships. This should be done after the constraints have been
928    defined using :func:`StoreEquivalence`, :func:`GroupConstraints` and :func:`GenerateConstraints`.
929
930    :returns: a string containing the details of the contraint relationships
931    '''
932    s = ''
933    if len(fixedVarList) > 0:
934        s += 'Fixed Variables:\n'
935        for v in fixedVarList:
936            s += '    ' + v + '\n'
937    s += 'User-supplied variable mapping relations:\n'
938    symout = ''
939    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict,symGenList
940
941    for varlist,mapvars,multarr,invmultarr,symFlag in zip(
942        dependentParmList,indParmList,arrayList,invarrayList,symGenList):
943        for i,mv in enumerate(mapvars):
944            if multarr is None:
945#                s1 = '  ' + str(mv) + ' is equivalent to parameter(s): '
946                if len(varlist) == 1:
947                    s1 = '   ' + str(mv) + ' is equivalent to '
948                else:
949                    s1 = '   ' + str(mv) + ' is equivalent to parameters: '
950                j = 0
951                for v,m in zip(varlist,invmultarr):
952                    if debug: print 'v,m[0]: ',v,m[0]
953                    if len(s1.split('\n')[-1]) > 75: s1 += '\n        '
954                    if j > 0: s1 += ' & '
955                    j += 1
956                    s1 += str(v)
957                    if m != 1:
958                        s1 += " / " + str(m[0])
959                if symFlag:
960                    symout += s1 + '\n'
961                else:
962                    s += s1 + '\n'
963                continue
964            s += %s = ' % mv
965            j = 0
966            for m,v in zip(multarr[i,:],varlist):
967                if m == 0: continue
968                if j > 0: s += ' + '
969                j += 1
970                s += '(%s * %s)' % (m,v)
971            if mv in varyList: s += ' VARY'
972            s += '\n'
973    if symout:
974        s += 'Symmetry-generated relations:\n' + symout
975    if inputOnly: return s
976    s += 'Inverse variable mapping relations:\n'
977    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
978        for i,mv in enumerate(varlist):
979            s += %s = ' % mv
980            j = 0
981            for m,v in zip(invmultarr[i,:],mapvars):
982                if m == 0: continue
983                if j > 0: s += ' + '
984                j += 1
985                s += '(%s * %s)' % (m,v)
986            s += '\n'
987    return s
988
989def GetSymEquiv():
990    '''Return the automatically generated (equivalence) relationships.
991
992    :returns: a list of strings containing the details of the contraint relationships
993    '''
994    symout = []
995    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict,symGenList
996
997    for varlist,mapvars,multarr,invmultarr,symFlag in zip(
998        dependentParmList,indParmList,arrayList,invarrayList,symGenList):
999        for i,mv in enumerate(mapvars):
1000            if not symFlag: continue
1001            if multarr is None:
1002                s1 = str(mv) + ' = '
1003                j = 0
1004                for v,m in zip(varlist,invmultarr):
1005                    if debug: print 'v,m[0]: ',v,m[0]
1006                    if len(s1.split('\n')[-1]) > 75: s1 += '\n        '
1007                    if j > 0: s1 += ' =  '
1008                    j += 1
1009                    s1 += str(v)
1010                    if m != 1:
1011                        s1 += " / " + str(m[0])
1012                    symout.append(s1)
1013                continue
1014            else:
1015                s = %s = ' % mv
1016                j = 0
1017                for m,v in zip(multarr[i,:],varlist):
1018                    if m == 0: continue
1019                    if j > 0: s += ' + '
1020                    j += 1
1021                    s += '(%s * %s)' % (m,v)
1022                print 'unexpected sym op=',s
1023    return symout
1024
1025def Dict2Deriv(varyList,derivDict,dMdv):
1026    '''Compute derivatives for Independent Parameters from the
1027    derivatives for the original parameters
1028
1029    :param list varyList: a list of parameters names that will be varied
1030
1031    :param dict derivDict: a dict containing derivatives for parameter values keyed by the
1032      parameter names.
1033
1034    :param list dMdv: a Jacobian, as a list of np.array containing derivatives for dependent
1035      parameter computed from derivDict
1036
1037    '''
1038    global dependentParmList,arrayList,invarrayList,indParmList,invarrayList
1039    for varlist,mapvars,multarr,invmultarr in zip(dependentParmList,indParmList,arrayList,invarrayList):
1040        for i,name in enumerate(mapvars):
1041            # grouped variables: need to add in the derv. w/r
1042            # dependent variables to the independent ones
1043            if name not in varyList: continue # skip if independent var not varied
1044            if multarr is None:
1045                for v,m in zip(varlist,invmultarr):
1046                    if debug: print 'start dMdv',dMdv[varyList.index(name)]
1047                    if debug: print 'add derv',v,'/',m[0],'to derv',name,'add=',derivDict[v] / m[0]
1048                    if m == 0: continue
1049                    dMdv[varyList.index(name)] += derivDict[v] / m[0]
1050            else:
1051                for v,m in zip(varlist,multarr[i,:]):
1052                    if debug: print 'start dMdv',dMdv[varyList.index(name)]
1053                    if debug: print 'add derv',v,'*',m,'to derv',name,'add=',m * derivDict[v]
1054                    if m == 0: continue
1055                    dMdv[varyList.index(name)] += m * derivDict[v]
1056
1057def Map2Dict(parmDict,varyList):
1058    '''Create (or update) the Independent Parameters from the original
1059    set of Parameters
1060
1061    Removes dependent variables from the varyList
1062
1063    This should be done once, after the constraints have been
1064    defined using :func:`StoreEquivalence`,
1065    :func:`GroupConstraints` and :func:`GenerateConstraints` and
1066    before any variable refinement is done. This completes the parameter
1067    dictionary by defining independent parameters and it satisfies the
1068    constraint equations in the initial parameters
1069
1070    :param dict parmDict: a dict containing parameter values keyed by the
1071      parameter names.
1072      This will contain updated values for both dependent and independent
1073      parameters after Dict2Map is called. It will also contain some
1074      unexpected entries of every constant value {'0':0.0} & {'1.0':1.0},
1075      which do not cause any problems.
1076
1077    :param list varyList: a list of parameters names that will be varied
1078   
1079
1080    '''
1081    # process the Independent vars: remove dependent ones from varylist
1082    # and then compute values for the Independent ones from their dependents
1083    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict
1084    for varlist,mapvars,multarr in zip(dependentParmList,indParmList,arrayList):
1085        for item in varlist:
1086            try:
1087                varyList.remove(item)
1088            except ValueError:
1089                pass
1090        if multarr is None: continue
1091        valuelist = [parmDict[var] for var in varlist]
1092        parmDict.update(zip(mapvars,
1093                            np.dot(multarr,np.array(valuelist)))
1094                        )
1095    # now remove fixed variables from the varyList
1096    global fixedVarList
1097    for item in fixedVarList:
1098        try:
1099            varyList.remove(item)
1100        except ValueError:
1101            pass
1102    # Set constrained parameters to their fixed values
1103    parmDict.update(fixedDict)
1104
1105def Dict2Map(parmDict,varyList):
1106    '''Applies the constraints defined using :func:`StoreEquivalence`,
1107    :func:`GroupConstraints` and :func:`GenerateConstraints` by changing
1108    values in a dict containing the parameters. This should be
1109    done before the parameters are used for any computations
1110
1111    :param dict parmDict: a dict containing parameter values keyed by the
1112      parameter names.
1113      This will contain updated values for both dependent and independent
1114      parameters after Dict2Map is called. It will also contain some
1115      unexpected entries of every constant value {'0':0.0} & {'1.0':1.0},
1116      which do not cause any problems.
1117
1118    :param list varyList: a list of parameters names that will be varied
1119   
1120    '''
1121    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict
1122    # reset fixed values (should not be needed, but very quick)
1123    # - this seems to update parmDict with {'0':0.0} & {'1.0':1.0} - probably not what was intended
1124    # not needed, but also not a problem - BHT
1125    parmDict.update(fixedDict)
1126    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
1127        #if invmultarr is None: continue
1128        try: 
1129            valuelist = [parmDict[var] for var in mapvars]
1130        except KeyError:
1131            continue
1132        parmDict.update(zip(varlist,np.dot(invmultarr,np.array(valuelist))))
1133
1134#======================================================================
1135# internal routines follow (these routines are unlikely to be called
1136# from outside the module)
1137
1138def GramSchmidtOrtho(a,nkeep=0):
1139    '''Use the Gram-Schmidt process (http://en.wikipedia.org/wiki/Gram-Schmidt) to
1140    find orthonormal unit vectors relative to first row.
1141
1142    If nkeep is non-zero, the first nkeep rows in the array are not changed
1143   
1144    input:
1145       arrayin: a 2-D non-singular square array
1146    returns:
1147       a orthonormal set of unit vectors as a square array
1148    '''
1149    def proj(a,b):
1150        'Projection operator'
1151        return a*(np.dot(a,b)/np.dot(a,a))
1152    for j in range(nkeep,len(a)):
1153        for i in range(j):
1154            a[j] -= proj(a[i],a[j])
1155        if np.allclose(np.linalg.norm(a[j]),0.0):
1156            raise Exception,"Singular input to GramSchmidtOrtho"
1157        a[j] /= np.linalg.norm(a[j])
1158    return a
1159
1160def _FillArray(sel,dict,collist,FillDiagonals=False):
1161    '''Construct a n by n matrix (n = len(collist)
1162    filling in the rows using the relationships defined
1163    in the dictionaries selected by sel
1164
1165    If FillDiagonals is True, diagonal elements in the
1166    array are set to 1.0
1167    '''
1168    n = len(collist)
1169    if FillDiagonals:
1170        arr = np.eye(n)
1171    else:
1172        arr = np.zeros(2*[n,])
1173    # fill the top rows
1174    for i,cnum in enumerate(sel):
1175        for j,var in enumerate(collist):
1176            arr[i,j] = dict[cnum].get(var,0)
1177    return arr
1178
1179def _SwapColumns(i,m,v):
1180    '''Swap columns in matrix m as well as the labels in v
1181    so that element (i,i) is replaced by the first non-zero element in row i after that element
1182
1183    Throws an exception if there are no non-zero elements in that row
1184    '''
1185    for j in range(i+1,len(v)):
1186        if not np.allclose(m[i,j],0):
1187            m[:,(i,j)] = m[:,(j,i)]
1188            v[i],v[j] = v[j],v[i]
1189            return
1190    else:
1191        raise Exception,'Singular input'
1192
1193def _RowEchelon(m,arr,collist):
1194    '''Convert the first m rows in Matrix arr to row-echelon form
1195    exchanging columns in the matrix and collist as needed.
1196
1197    throws an exception if the matrix is singular because
1198    the first m rows are not linearly independent
1199    '''
1200    for i in range(m):
1201        if np.allclose(arr[i,i],0):
1202            _SwapColumns(i,arr,collist)
1203        arr[i,:] /= arr[i,i] # normalize row
1204        # subtract current row from subsequent rows to set values to left of diagonal to 0
1205        for j in range(i+1,m):
1206            arr[j,:] -= arr[i,:] * arr[j,i]
1207
1208if __name__ == "__main__":
1209    parmdict = {}
1210    constrDict = [
1211        {'0:12:Scale': 2.0, '0:11:Scale': 1.0, '0:14:Scale': 4.0, '0:13:Scale': 3.0, '0:0:Scale': 0.5},
1212        {'0:0:eA': 0.0},
1213        {'2::C(10,6,1)': 1.0, '1::C(10,6,1)': 1.0},
1214        {'1::C(10,0,1)': 1.0, '2::C(10,0,1)': 1.0},
1215        {'1::AUiso:0': 1.0, '0::AUiso:0': 1.0},
1216        {'0::A0': 0.0}
1217        ]
1218    fixedList = ['5.0', '0', None, None, '1.0', '0']
1219    StoreEquivalence('2::atomx:3',('2::atomy:3', ('2::atomz:3',2,), ))
1220    #StoreEquivalence('1::atomx:3',('2::atomx:3', ('2::atomz:3',2,), )) # error: dependent & independent vars mixed
1221    #StoreEquivalence('1::atomx:3',('2::atomy:3', ('2::atomz:3',2,), )) # error: dependent vars repeated
1222    #StoreEquivalence('0:1:eA',('0:0:eA',)) # error: equiv & fixed
1223    #StoreEquivalence('0:99:Scale',('0:12:Scale',)) # error: equiv & constrained
1224    #StoreEquivalence('0:12:Scale',('0:99:Scale',)) # error: equiv & constrained
1225    varylist = ['2::atomx:3',
1226                '2::C(10,6,1)', '1::C(10,6,1)',
1227                '2::atomy:3', '2::atomz:3',
1228                '0:12:Scale', '0:11:Scale', '0:14:Scale', '0:13:Scale', '0:0:Scale']
1229#    e,w = CheckConstraints([,
1230#                     [{'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}],
1231#                     ['1.0'])
1232#    if e: print 'error=',e
1233#    if w: print 'error=',w
1234#    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']
1235#    constrDict = [
1236#        {'0::Afrac:4': 24.0, '0::Afrac:1': 16.0, '0::Afrac:3': 24.0, '0::Afrac:2': 16.0},
1237#        {'0::Afrac:1': 1.0, '0::Afrac:2': 1.0},
1238#        {'0::Afrac:4': 1.0, '0::Afrac:3': 1.0}]
1239#    fixedList = ['40.0', '1.0', '1.0']
1240
1241    errmsg, warnmsg = CheckConstraints(varylist,constrDict,fixedList)
1242    if errmsg:
1243        print "*** Error ********************"
1244        print errmsg
1245    if warnmsg:
1246        print "*** Warning ********************"
1247        print warnmsg
1248    if errmsg or warnmsg:
1249        sys.exit()
1250    groups,parmlist = GroupConstraints(constrDict)
1251    GenerateConstraints(groups,parmlist,varylist,constrDict,fixedList)
1252    print VarRemapShow(varylist)
1253    parmdict.update( {
1254        '0:12:Scale': 1.0, '0:11:Scale': 1.0, '0:14:Scale': 1.0, '0:13:Scale': 1.0, '0:0:Scale': 2.0,
1255        '0:0:eA': 0.0,
1256        '2::C(10,6,1)': 0.2, '1::C(10,6,1)': 0.3,
1257        '1::C(10,0,1)': 0.2, '2::C(10,0,1)': 0.3,
1258        '1::AUiso:0': 0.02, '0::AUiso:0': 0.03,
1259        '0::A0': 0.0,
1260        '2::atomx:3':0.23,'2::atomy:3':-.23, '2::atomz:3':-0.11,
1261        })
1262    print 'parmdict start',parmdict
1263    print 'varylist start',varylist
1264    before = parmdict.copy()
1265    Map2Dict(parmdict,varylist)
1266    print 'parmdict before and after Map2Dict'
1267    print '  key / before / after'
1268    for key in sorted(parmdict.keys()):
1269        print '  '+key,'\t',before.get(key),'\t',parmdict[key]
1270    print 'varylist after',varylist
1271    before = parmdict.copy()
1272    Dict2Map(parmdict,varylist)
1273    print 'after Dict2Map'
1274    print '  key / before / after'
1275    for key in sorted(parmdict.keys()):
1276        print '  '+key,'\t',before.get(key),'\t',parmdict[key]
1277#    dMdv = len(varylist)*[0]
1278#    deriv = {}
1279#    for i,v in enumerate(parmdict.keys()): deriv[v]=i
1280#    Dict2Deriv(varylist,deriv,dMdv)
Note: See TracBrowser for help on using the repository browser.