source: trunk/GSASIImapvars.py @ 3057

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

more work on constraints: prevent use of perviously constrained variables in new equivalences; include sym. generated constraints in checks; do not allow error-generating constraints to be added

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