source: trunk/GSASIImapvars.py @ 3508

Last change on this file since 3508 was 3371, checked in by toby, 7 years ago

constraint fixes: implement wild-card equivalences, fix use of all for phase selection in atom vars; prevent use of Edit Constr menu items for sym-generated equivalences

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