source: trunk/GSASIImapvars.py @ 4379

Last change on this file since 4379 was 4379, checked in by toby, 19 months ago

fix isodistort constraint processing

  • Property svn:eol-style set to native
  • Property svn:keywords set to Date Author Revision URL Id
File size: 78.8 KB
Line 
1# -*- coding: utf-8 -*-
2########### SVN repository information ###################
3# $Date: 2020-03-23 01:53:24 +0000 (Mon, 23 Mar 2020) $
4# $Author: toby $
5# $Revision: 4379 $
6# $URL: trunk/GSASIImapvars.py $
7# $Id: GSASIImapvars.py 4379 2020-03-23 01:53:24Z toby $
8########### SVN repository information ###################
9"""
10*GSASIImapvars: Parameter constraints*
11======================================
12
13Module to implements algebraic contraints, parameter redefinition
14and parameter simplification contraints.
15
16Types of constraints
17--------------------
18
19There are four ways to specify constraints, as listed below.
20Note that parameters are initially stored in the
21main section of the GSAS-II data tree under heading ``Constraints``.
22This dict has four keys, 'Hist', 'HAP', 'Global', and 'Phase',
23each containing a list of constraints. An additional set of constraints
24are generated for each phase based on symmetry considerations by calling
25:func:`GSASIIstrIO.GetPhaseData`.
26
27Note that in the constraints, as stored in the GSAS-II data tree, parameters
28are stored as :class:`GSASIIobj.G2VarObj` objects, as these objects allow for
29changes in numbering of phases, histograms and atoms.
30When they are interpreted (in :func:`GSASIIstrIO.ProcessConstraints`),
31references
32to numbered objects are resolved using the appropriate random ids and the
33parameter object is converted to a string of form ``ph:hst:VARNAM:at``.
34
35Alternate parameters (New Var)
36^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
37
38Parameter redefinition ("New Var" constraints)
39is done by creating an expression that relates several
40parameters:
41
42::
43
44   Mx1 * Px + My1 * Py +...
45   Mx2 * Px + Mz2 * Pz + ...
46
47where Pj is a GSAS-II parameter name and Mjk is a constant (float) multiplier.
48Alternately, multipliers Mjk can contain a formula (str) that will be evaluated prior
49to the start of the refinement. In a formula, GSAS-II parameters will be replaced by the
50value of the parameter before the formula is evaluated, so ``'np.cos(0::Ax:2)'`` is a valid
51multiplier. At present, only phase (atom/cell) parameters are available for use in
52a formula, but this can be expanded if needed.
53
54This type of constraint describes an alternate
55degree of freedom where parameter Px and Py, etc. are varied to keep
56their ratio
57fixed according the expression. A new variable parameter is assigned to each degree of
58freedom when refined. An example where this can be valuable is when
59two parameters, P1 and P2, have similar values and are highly correlated. It is often better to create a new variable, Ps = P1 + P2, and refine Ps.
60In the later stages of refinement, a second
61variable, Pd = P1 - P2, can be defined and it can be seen if refining Pd is
62supported by the data. Another use will be to define parameters that
63express "irrep modes" in terms of the fundamental structural parameters.
64
65These "New Var" constraints are stored as described for type "f" in the
66:ref:`constraint definitions table <Constraint_definitions_table>`.
67
68Constrained parameters (Const)
69^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
70
71A constraint on a set of variables can be supplied in the form of a
72linear algebraic equation: ::
73
74  Nx1 * Px + Ny1 * Py +... = C1
75
76where Cn is a constant (float), where Pj is a GSAS-II parameter name,
77and where Njk is a constant multiplier (float) or a formula (str) that will be evaluated prior
78to the start of the refinement. In a formula, GSAS-II parameters will be replaced by the
79value of the parameter before the formula is evaluated, so ``'np.cos(0::Ax:2)'`` is a valid
80multiplier. At present, only phase (atom/cell) parameters are available for use in
81a formula, but this can be expanded if needed.
82
83These equations set an interdependence between parameters.
84Common uses of parameter constraints are to set rules that decrease the number of parameters,
85such as restricting the sum of fractional occupancies for atoms that share
86a site to sum to unity, thus reducing the effective number of variables by one.
87Likewise, the Uiso value for a H atom "riding" on a C, N or O atom
88can be related by a fixed offset (the so called B+1 "rule").
89
90A "Const" constraint is stored as described for type "c" in the
91:ref:`constraint definitions table <Constraint_definitions_table>`.
92
93Equivalenced parameters (Equiv)
94^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
95
96A simplifed way to set up a constraint equation is to define an equivalence,
97which can be of form: ::
98
99  C1 * P1 = C2 * Py
100
101or::
102
103  C1 * P1 = C2 * P2 = C3 * P3 = ...
104
105where Cn is a constant (float), where Pj is a GSAS-II parameter name. This
106means that parameters Py (or P2 and P3) are determined from (or "slaved" to)
107parameter P1. Alternately, equivalences can be created with :func:`StoreEquivalence`
108where the multipliers can be a formula (str) that will be evaluated prior to the start of
109the refinement. In a formula, GSAS-II parameters will be replaced by the value of the
110parameter before the formula is evaluate, so ``'np.cos(0::Ax:2)'`` is a valid multiplier.
111At present, only phase (atom/cell) parameters are available for use in
112a formula, but this can be expanded if needed.
113Note that
114the latter constraint expression is conceptually identical to
115defining constraint equations. In practice, however,
116equivalenced parameters are processed in a
117different and more direct manner than constraint equations. The previous
118set of equalities could also be written in this way as a set of constraint
119equations: ::
120
121  C1 * P1 - C2 * P2 = 0
122  C1 * P1 - C3 * P3 = 0
123  ...
124
125
126The first parameter (P1 above)
127is considered the independent variable
128and the remaining parameters are dependent variables. The dependent variables
129are set from the independent variable.
130An example of how this may be used woul be if, for example,
131a material has a number of O atoms, all in fairly similar bonding environments
132and the diffraction data are sparse, one my reduce the complexity of the model
133by defining Uiso for the first O atoms to be identical to the remaining atoms.
134The results of this refinement will be simpler to understand than if a set of
135constraint equations is used because the refined parameter will be the independent variable, which will be as ph::Uiso:n, corresponding to the first O atom.
136
137A parameter can be used in multiple
138equivalences as independent variable,
139but if parameter is used as both a dependent and independent variable or
140a parameter is used in equivalences and in "New Var" or "Const" constraints,
141this create conflicts that cannot be resolved within the equivalences implementation
142but can be handled as constraint equations.
143The equivalences that violate this are discovered in :func:`CheckEquivalences`
144and then :func:`MoveConfEquiv` is used to change these equivalences to "Const" equations.
145
146Equivalenced parameters ("EQUIV" constraints), when defined by users,
147are stored as described for type "e" in the
148:ref:`constraint definitions table <Constraint_definitions_table>`.
149Other equvalences are generated by symmetry prior
150to display or refinement in :func:`GSASIIstrIO.GetPhaseData`.
151These are not stored.
152
153Fixed parameters (Hold)
154^^^^^^^^^^^^^^^^^^^^^^^^
155
156When parameters are refined where a single refinement flag determines that several variables
157are refined at the same time (examples are: cell parameters, atom positions, anisotropic
158displacement parameters, magnetic moments,...) it can be useful to specify that a
159specific parameter should not be varied. These will most commonly be generated due to symmetry,
160but under specific conditions, there may be other good reasons to constrain a parameter.
161
162A "Hold" constraint is stored as described for type "h" in the
163:ref:`constraint definitions table <Constraint_definitions_table>`.
164
165Constraint Processing
166---------------------
167
168When constraints will be used or edited, they are processed using a series of
169calls:
170
171* First all of the stored constraints are appended into a single list. They are
172  initially stored in separate lists only to improve their creation and display
173  in the GUI.
174
175* Then :func:`InitVars` is used to initialize the global variables in
176  this module (:mod:`GSASIImapvars`).
177
178* Then :func:`GSASIIstrIO.ProcessConstraints` is used to initially process the
179  constraints, as described below.
180
181* Symmetry-generated equivalences are then created in
182  :func:`GSASIIstrIO.GetPhaseData`, which also calls
183  :func:`GSASIIstrIO.cellVary` and for Pawley refinements
184  :func:`GSASIIstrIO.GetPawleyConstr`. These are entered directly into this
185  module's globals using :func:`StoreEquivalence`.
186* Constraints/equivalences are then checked for possible conflicts with
187  :func:`CheckConstraints`, this requires grouping the constraints,
188  as described below.
189* In refinements, :func:`GenerateConstraints` is then called to
190  create the constraints that will be used, see below for
191* For debugging constraints, :func:`VarRemapShow` can be called after
192  :func:`GenerateConstraints` to display the generated constraints.
193
194Constraint Reorganization (:func:`~GSASIIstrIO.ProcessConstraints`)
195^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
196
197:func:`GSASIIstrIO.ProcessConstraints` is used to initially process the
198constraints. This does these things:
199
2001. The "Hold", "Const" and "New Var" expressions are split between two paired lists,
201   :data:`constDictList` and :data:`fixedList` which are set:
202 
203   * For "Hold" entries a dict with a single entry is placed in constDictList where
204     the key is the parameter name (associated value is 0.0) and fixedList gets a
205     value of 0.0.
206   * For "Const" entries, a dict with multiple entries is placed in constDictList where
207     the key is the parameter name and the value is the multiplier for the parameter,
208     while fixedList gets a string value corresponding to the constant value for
209     the expression.
210   * For "New Var" entries, a dict with multiple entries is placed in constDictList
211     where the key is the parameter name and the value is the multiplier
212     for the parameter; an additional key "_vary" is given the value of True or False
213     depending on the refinement flag setting. The corresponding entry in
214     fixedList is None
215
216   The output from this will have this form where the first entry is a "Const", the
217   second is a "New Var" and the final is a "Hold".
218
219  .. code-block:: python
220
221    constrDict = [
222         {'0:12:Scale': 2.0, '0:14:Scale': 4.0, '0:13:Scale': 3.0, '0:0:Scale': 0.5},
223         {'2::C(10,6,1)': 1.0, '1::C(10,6,1)': 1.0, '_vary':True},
224         {'0::A0': 0.0}]
225    fixedList = ['5.0', None, '0']
226
227
2282. Equivalences are stored using :func:`StoreEquivalence`
229   into this module's globals (:data:`~GSASIImapvars.arrayList`,
230   :data:`~GSASIImapvars.invarrayList`,    :data:`~GSASIImapvars.indParmList`,
231   :data:`~GSASIImapvars.dependentParmList` and  :data:`~GSASIImapvars.symGenList`)
232
233
234Parameter Grouping (:func:`GenerateConstraints`)
235^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
236
237Functions :func:`CheckConstraints` and :func:`GenerateConstraints` are used to
238process the parameter equivalences and constraint lists created in
239:func:`~GSASIIstrIO.ProcessConstraints`. The former is used to generate
240error messages and the latter to generate the internal information used to
241apply the constraints.
242
243Initially, in both a list of parameters that are fixed and those used in
244constraint relations are tabulated in :func:`CheckEquivalences`. The equivalence
245relations are the scanned for the following potential problems:
246
2471. a parameter is used as a dependent variable in more than one
248   equivalence relation
2492. a parameter is fixed and used in in an equivalence relation either
250   as a dependent or independent variable
2513. a parameter is used as a dependent variable in one equivalence relation and
252   as a independent variable in another
2534. a parameter is used in in an equivalence relation (either
254   as a dependent or independent variable) and is used in a constraint
255   expression
2565. a parameter is not defined in a particular refinement, but is used in an
257   equivalence relation
2586. a parameter uses a wildcard for the histogram number (sequential refinements)
259
260Cases 1 & 2 above cannot be corrected, and result in errors. Cases 3 & 4 are
261potentially corrected with :func:`MoveConfEquiv`, as described below.
262Case 5 causes the equivalence to
263be dropped. Case 6 causes the current histogram number to be substituted for
264the wildcard.
265
266For cases 3 & 4, :func:`MoveConfEquiv` is used to change these equivalences
267into "Const" equations. This can potentially mean that additional
268equivalences will be problematic, so if there are changes made by
269:func:`MoveConfEquiv`, :func:`CheckEquivalences` is repeated. If any problem
270cases are noted, the refinement cannot be performed.
271
272Constraint expressions ("Const" and "New Var") are sorted into
273groups so that each group contains the minimum number of entries that
274ensures each parameter is referenced in only one group in
275:func:`GroupConstraints`. This is done by scanning the
276list of dicts in :data:`constDictList` one by one and making a list
277of parameters used in that constraint expression. Any expression that contains
278a parameter in is in that list is added to the current group and those
279parameters are added to this list of parameters. The list of ungrouped
280expressions is then scanned again until no more expressions are added to the
281current group. This process is repeated until every expression has been
282placed in a group. Function :func:`GroupConstraints` returns two lists of lists.
283The first has, for each group, a list of the indices in :data:`constDictList`
284that comprise the group (there can be only one). The second list contains,
285for each group, the unique parameter names in that group.
286
287Each constraint group is then processed. First, wildcard parameters are
288renamed (in a sequential refinement). Any fixed parameters that are used
289in constraints are noted as errors. The number of refined parameters and
290the number of parameters that are not defined in the current refinement are
291also noted. It is fine if all parameters in a group are not defined or all are
292not varied, but if some are defined and others not or some are varied and
293others not, this constitutes an error.
294
295The contents of each group is then examined. Groups with a single
296parameter (holds) are ignored. Then for each group, the number
297of parameters in the group (Np) and the number of expressions in the
298group (Nc) are counted and for each expression. If Nc > Np, then the constraint
299is overdetermined, which also constitutes an error.
300
301The parameter multipliers for each expression are then assembled:
302
303::
304
305   M1a * P1 + M2a * P2 +... Mka * Pk
306   M1b * P1 + M2b * P2 +... Mkb * Pk
307   ...
308   M1j * P1 + M2j * P2 +... Mkj * Pk
309
310
311From this it becomes possible to create a Nc x Np matrix, which is
312called the constraint matrix:
313
314 .. math::
315
316   \\left( \\begin{matrix}
317   M_{1a}  & M_{2a} &... & M_{ka} \\\\
318   M_{1b}  & M_{2b} &... & M_{kb} \\\\
319   ...  \\\\
320   M_{1j}  & M_{2j}  &... & M_{kj}
321   \\end{matrix}\\right)
322
323When Nc<Np, then additional rows need to be added to the matrix and to
324the vector that contains the value for each row (:data:`fixedList`) where
325values are ``None`` for New Vars and a constant for fixed values.
326This then can describe a system of Np simultaneous equations:
327
328 .. math::
329
330   \\left( \\begin{matrix}
331   M_{1a}  & M_{2a} &... & M_{ka} \\\\
332   M_{1b}  & M_{2b} &... & M_{kb} \\\\
333   ...  \\\\
334   M_{1j}  & M_{2j}  &... & M_{kj}
335   \\end{matrix}\\right)
336   \\left( \\begin{matrix}
337   P_{1} \\\\
338   P_{2} \\\\
339   ...  \\\\
340   P_{k}
341   \\end{matrix}\\right)
342   =
343   \\left( \\begin{matrix}
344   C_{1} & C_{2} &  ... & C_{k}
345   \\end{matrix}\\right)
346
347This is done by creating a square matrix from the group using
348:func:`_FillArray` with parameter ``FillDiagonals=False`` (the default). Any
349unspecified rows are left as all zero. The first Nc rows in the
350array are then coverted to row-echelon form using :func:`_RowEchelon`. This
351will create an Exception if any two rows are linearly dependent (which means
352that no matter what values are used for the remaining rows, that the matrix
353will be singular). :func:`_FillArray` is then called with parameter
354``FillDiagonals=True``, which again creates a square matrix but where
355unspecified rows are zero except for the diagonal elements. The 
356`Gram-Schmidt process <http://en.wikipedia.org/wiki/Gram-Schmidt>`_,
357implemented  in :func:`GramSchmidtOrtho`, is used to find orthonormal unit
358vectors for the remaining Np-Nc rows of the matrix. This will fail with
359a ``ConstraintException`` if this is not possible (singular matrix) or
360the result is placed in :data:`constrArr` as a numpy array.
361
362Rows in the matrix corresponding to "New Var" constraints and those that
363were generated by the Gram-Schmidt process are provided with parameter names
364(this can be specified if a "New Var" entry by using a ``"_name"`` element
365in the constraint dict, but at present this is not implemented.) Names are
366generated using :data:`paramPrefix` which is set to ``"::constr"``, plus a
367number to make the new parameter name unique. Global dict :data:`genVarLookup`
368provides a lookup table, where the names of the parameters related to this new
369parameter can be looked up easily.
370
371Finally the parameters used as input to the constraint are placed in
372this module's globals
373:data:`~GSASIImapvars.dependentParmList` and the constraint matrix is
374placed in into  :data:`~GSASIImapvars.arrayList`. This can be used to compute
375the initial values for "New Var" parameters. The inverse of the
376constraint matrix is placed in :data:`~GSASIImapvars.invarrayList` and a list
377of the "New Var" parameters and a list of the fixed values (as str's)
378is placed in :data:`~GSASIImapvars.indParmList`. A lookup table for
379fixed values as floats is placed in :data:`~GSASIImapvars.fixedDict`.
380Finally the appropriate entry in :data:`~GSASIImapvars.symGenList` is set to
381False to indicate that this is not a symmetry generated constraint.
382
383
384*Externally-Accessible Routines*
385---------------------------------
386
387To define a set of constrained and unconstrained relations, one
388defines a list of dictionary defining constraint parameters and their
389values, a list of fixed values for each constraint and a list of
390parameters to be varied. In addition, one uses
391:func:`StoreEquivalence` to define parameters that are equivalent.
392Use :func:`EvaluateMultipliers` to convert formula-based constraint/equivalence
393multipliers to numbers and then
394use :func:`CheckConstraints` to check that the input is
395internally consistent and finally :func:`GroupConstraints` and
396:func:`GenerateConstraints` to generate the internally used
397tables. Routine :func:`Map2Dict` is used to initialize the parameter
398dictionary and routine :func:`Dict2Map`, :func:`Dict2Deriv`, and
399:func:`ComputeDepESD` are used to apply constraints. Routine
400:func:`VarRemapShow` is used to print out the constraint information,
401as stored by :func:`GenerateConstraints`. Further information on each routine
402is below:
403
404:func:`InitVars`
405  This is optionally used to clear out all defined previously defined constraint information
406 
407:func:`StoreEquivalence`
408  To implement parameter redefinition, one calls StoreEquivalence. This should be called for every set of
409  equivalence relationships. There is no harm in using StoreEquivalence with the same independent variable:
410
411  .. code-block:: python
412
413       StoreEquivalence('x',('y',))
414       StoreEquivalence('x',('z',))
415
416  or equivalently
417
418  .. code-block:: python
419
420       StoreEquivalence('x',('y','z'))
421
422  The latter will run more efficiently. Note that mixing independent
423  and dependent variables would require assignments, such as
424
425  .. code-block:: python
426
427        StoreEquivalence('x',('y',))
428        StoreEquivalence('y',('z',))
429
430  would require that equivalences be applied in a particular order and
431  thus is implemented as a constraint equation rather than an equivalence.
432       
433  Use StoreEquivalence before calling GenerateConstraints or CheckConstraints
434
435:func:`CheckConstraints`
436   check that input in internally consistent
437
438:func:`GenerateConstraints`
439   generate the internally used tables from constraints and equivalences
440
441:func:`EvaluateMultipliers`
442   Convert any string-specified (formula-based) multipliers to numbers. Call this before
443   using :func:`CheckConstraints` or :func:`GenerateConstraints`.
444   At present, the code may pass only the dict for phase (atom/cell) parameters,
445   but this could be expanded if needed.
446
447:func:`Map2Dict`
448   To determine values for the parameters created in this module, one
449   calls Map2Dict. This will not apply contraints.
450
451:func:`Dict2Map`
452   To take values from the new independent parameters and constraints,
453   one calls Dict2Map and set the parameter array, thus appling contraints.
454
455:func:`Dict2Deriv`
456   Use Dict2Deriv to determine derivatives on independent parameters
457   from those on dependent ones.
458
459:func:`ComputeDepESD`     
460   Use ComputeDepESD to compute uncertainties on dependent variables.
461
462:func:`VarRemapShow`
463   To show a summary of the parameter remapping, one calls VarRemapShow.
464
465*Global Variables*
466------------------
467
468dependentParmList:
469   a list containing group of lists of
470   parameters used in the group. Note that parameters listed in
471   dependentParmList should not be refined as they will not affect
472   the model
473
474indParmList:
475     a list containing groups of Independent parameters defined in
476     each group. This contains both parameters used in parameter
477     redefinitions as well as names of generated new parameters.
478
479arrayList:
480     a list containing group of relationship matrices to relate
481     parameters in dependentParmList to those in indParmList. Unlikely
482     to be used externally.
483
484invarrayList:
485     a list containing group of relationship matrices to relate
486     parameters in indParmList to those in dependentParmList. Unlikely
487     to be used externally.
488
489fixedVarList:
490     a list of parameters that have been 'fixed'
491     by defining them as equal to a constant (::var: = 0). Note that
492     the constant value is ignored at present. These parameters are
493     later removed from varyList which prevents them from being refined.
494     Unlikely to be used externally.
495
496fixedDict:
497     a dictionary containing the fixed values corresponding
498     to parameter equations.  The dict key is an ascii string, but the
499     dict value is a float.  Unlikely to be used externally.
500
501symGenList:
502     a list of boolean values that will be True to indicate that a constraint
503     (only equivalences) is generated by symmetry (or Pawley overlap)
504     
505problemVars:
506     a list containing parameters that show up in constraints producing errors
507
508
509
510*Routines/variables*
511---------------------
512
513Note that parameter names in GSAS-II are strings of form ``<ph#>:<hst#>:<nam>`` or ``<ph#>::<nam>:<at#>``.
514"""
515
516from __future__ import division, print_function
517import numpy as np
518import sys
519import GSASIIpath
520GSASIIpath.SetVersionNumber("$Revision: 4379 $")
521# data used for constraints;
522debug = False # turns on printing as constraint input is processed
523
524# note that constraints are stored listed by contraint groups,
525# where each constraint
526# group contains those parameters that must be handled together
527dependentParmList = []
528'''a list of lists where each item contains a list of parameters in each constraint group.
529note that parameters listed in dependentParmList should not be refined directly.'''
530indParmList = [] # a list of names for the new parameters
531'''a list of lists where each item contains a list for each constraint group with
532fixed values for constraint equations and names of generated (New Var) parameters.
533'''
534arrayList = []
535'''a list of of relationship matrices that map model parameters in each
536constraint group (in :data:`dependentParmList`) to
537generated (New Var) parameters.
538'''
539invarrayList = []
540'''a list of of inverse-relationship matrices that map constrained values and
541generated (New Var) parameters (in :data:`indParmList`) to model parameters
542(in :data:`dependentParmList`).
543'''
544fixedDict = {}
545'''A dict lookup-table containing the fixed values corresponding
546to defined parameter equations. Note the key is the original ascii string
547and the value in the dict is a float.
548'''
549fixedVarList = []
550'''List of parameters that should not be refined.
551'''
552symGenList = []
553'''A list of flags that if True indicates a constraint was generated by symmetry
554'''
555problemVars = []
556'''a list of parameters causing errors
557'''
558dependentVars = []
559'A list of dependent variables, taken from (:data:`dependentParmList`).'
560independentVars = []
561'A list of dependent variables, taken from (:data:`indParmList`).'
562genVarLookup = {}
563'provides a list of parameters that are related to each generated parameter'
564paramPrefix = "::constr"
565'A prefix for generated parameter names'
566consNum = 0
567'The number to be assigned to the next constraint to be created'
568
569class ConstraintException(Exception):
570    '''Defines an Exception that is used when an exception is raised processing constraints
571    '''
572    pass
573
574def InitVars():
575    '''Initializes all constraint information'''
576    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict,consNum,symGenList
577    dependentParmList = [] # contains a list of parameters in each group
578    arrayList = [] # a list of of relationship matrices
579    invarrayList = [] # a list of inverse relationship matrices
580    indParmList = [] # a list of names for the new parameters
581    fixedDict = {} # a dictionary containing the fixed values corresponding to defined parameter equations
582    symGenList = [] # Flag if constraint is generated by symmetry
583    consNum = 0 # number of the next constraint to be created
584    global genVarLookup
585    genVarLookup = {}
586
587def VarKeys(constr):
588    """Finds the keys in a constraint that represent parameters
589    e.g. eliminates any that start with '_'
590
591    :param dict constr: a single constraint entry of form::
592
593        {'var1': mult1, 'var2': mult2,... '_notVar': val,...}
594
595        (see :func:`GroupConstraints`)
596    :returns: a list of keys where any keys beginning with '_' are
597      removed.
598    """
599    return [i for i in constr.keys() if not i.startswith('_')]
600
601
602def GroupConstraints(constrDict):
603    """divide the constraints into groups that share no parameters.
604
605    :param dict constrDict: a list of dicts defining relationships/constraints
606
607    ::
608   
609       constrDict = [{<constr1>}, {<constr2>}, ...]
610
611    where {<constr1>} is {'var1': mult1, 'var2': mult2,... }
612
613    :returns: two lists of lists:
614   
615      * a list of grouped contraints where each constraint grouped containts a list
616        of indices for constraint constrDict entries
617      * a list containing lists of parameter names contained in each group
618     
619      """
620    assignedlist = [] # relationships that have been used
621    groups = [] # contains a list of grouplists
622    ParmList = []
623    for i,consi in enumerate(constrDict):
624        if i in assignedlist: continue # already in a group, skip
625        # starting a new group
626        grouplist = [i,]
627        assignedlist.append(i)
628        groupset = set(VarKeys(consi))
629        changes = True # always loop at least once
630        while(changes): # loop until we can't find anything to add to the current group
631            changes = False # but don't loop again unless we find something
632            for j,consj in enumerate(constrDict):
633                if j in assignedlist: continue # already in a group, skip
634                if len(set(VarKeys(consj)) & groupset) > 0: # true if this needs to be added
635                    changes = True
636                    grouplist.append(j)
637                    assignedlist.append(j)
638                    groupset = groupset | set(VarKeys(consj))
639        group = sorted(grouplist)
640        varlist = sorted(list(groupset))
641        groups.append(group)
642        ParmList.append(varlist)
643    return groups,ParmList
644
645def CheckConstraints(varyList,constrDict,fixedList):
646    '''Takes a list of relationship entries comprising a group of
647    constraints and checks for inconsistencies such as conflicts in
648    parameter/variable definitions and or inconsistently varied parameters.
649
650    :param list varyList: a list of parameters names that will be varied
651
652    :param dict constrDict: a list of dicts defining relationships/constraints
653      (as created in :func:`GSASIIstrIO.ProcessConstraints` and
654      documented in :func:`GroupConstraints`)
655
656    :param list fixedList: a list of values specifying a fixed value for each
657      dict in constrDict. Values are either strings that can be converted to
658      floats or ``None`` if the constraint defines a new parameter rather
659      than a constant.
660
661    :returns: two strings:
662
663      * the first lists conflicts internal to the specified constraints
664      * the second lists conflicts where the varyList specifies some
665        parameters in a constraint, but not all
666       
667      If there are no errors, both strings will be empty
668    '''
669    import re
670    global dependentParmList,arrayList,invarrayList,indParmList,consNum
671    global problemVars
672    # Process the equivalences
673    #    If there are conflicting parameters, move them into constraints. This
674    #    may create new conflicts, requiring movement of other parameters, so
675    #    loop until there are no more changes to make.
676    parmsChanged = True
677    while parmsChanged:
678        parmsChanged = 0
679        errmsg,warnmsg,fixVlist,dropVarList,translateTable = CheckEquivalences(
680            constrDict,varyList)
681        #print('debug: using MoveConfEquiv to address =',errmsg)
682        if problemVars: parmsChanged = MoveConfEquiv(constrDict,fixedList)
683#    GSASIIpath.IPyBreak()
684
685    groups,parmlist = GroupConstraints(constrDict)
686    # scan through parameters in each relationship. Are all varied? If only some are
687    # varied, create a warning message.
688    for group,varlist in zip(groups,parmlist):
689        if len(varlist) == 1:   # process fixed (held) variables
690            var = varlist[0]
691            if var not in fixedVarList:
692                fixedVarList.append(var)
693            continue
694        for rel in group:
695            varied = 0
696            notvaried = ''
697            for var in constrDict[rel]:
698                if var.startswith('_'): continue
699                if not re.match('[0-9]*:[0-9\*]*:',var):
700                    warnmsg += "\nParameter "+str(var)+" does not begin with a ':'"
701                if var in varyList:
702                    varied += 1
703                else:
704                    if notvaried: notvaried += ', '
705                    notvaried += var
706                if var in fixVlist:
707                    errmsg += '\nParameter '+var+" is Fixed and used in a constraint:\n\t"
708                    if var not in problemVars: problemVars.append(var)
709                    errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
710            if varied > 0 and varied != len(VarKeys(constrDict[rel])):
711                warnmsg += "\nNot all parameters refined in constraint:\n\t"
712                warnmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
713                warnmsg += '\nNot refined: ' + notvaried + '\n'
714    if errmsg or warnmsg:
715        return errmsg,warnmsg
716
717    # now look for process each group and create the relations that are needed to form
718    # non-singular square matrix
719    for group,varlist in zip(groups,parmlist):
720        if len(varlist) == 1: continue # a constraint group with a single parameter can be ignored
721        if len(varlist) < len(group): # too many relationships -- no can do
722            errmsg += "\nOver-constrained input. "
723            errmsg += "There are more constraints " + str(len(group))
724            errmsg += "\n\tthan parameters " + str(len(varlist)) + "\n"
725            for rel in group:
726                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
727                errmsg += "\n"
728                continue
729        try:
730            multarr = _FillArray(group,constrDict,varlist)
731            _RowEchelon(len(group),multarr,varlist)
732        except:
733            errmsg += "\nSingular input. "
734            errmsg += "There are internal inconsistencies in these constraints\n"
735            for rel in group:
736                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
737                errmsg += "\n"
738            continue
739        try:
740            multarr = _FillArray(group,constrDict,varlist,FillDiagonals=True)
741            GramSchmidtOrtho(multarr,len(group))
742        except:
743            errmsg += "\nUnexpected singularity with constraints (in Gram-Schmidt)\n"
744            for rel in group:
745                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
746                errmsg += "\n"
747            continue
748        mapvar = []
749        group = group[:]
750        # scan through all generated and input parameters
751        # Check again for inconsistent parameter use
752        # for new parameters -- where varied and unvaried parameters get grouped
753        # together. I don't think this can happen when not flagged before, but
754        # it does not hurt to check again.
755        for i in range(len(varlist)):
756            varied = 0
757            notvaried = ''
758            if len(group) > 0:
759                rel = group.pop(0)
760                fixedval = fixedList[rel]
761                for var in VarKeys(constrDict[rel]):
762                    if var in varyList:
763                        varied += 1
764                    else:
765                        if notvaried: notvaried += ', '
766                        notvaried += var
767            else:
768                fixedval = None
769            if fixedval is None:
770                varname = paramPrefix + str(consNum) # assign a name to a parameter
771                mapvar.append(varname)
772                consNum += 1
773            else:
774                mapvar.append(fixedval)
775            if varied > 0 and notvaried != '':
776                warnmsg += "\nNot all parameters refined in generated constraint"
777                warnmsg += '\nPlease report this unexpected error\n'
778                for rel in group:
779                    warnmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
780                    warnmsg += "\n"
781                warnmsg += '\n\tNot refined: ' + notvaried + '\n'
782        try:
783            np.linalg.inv(multarr)           
784        except:
785            errmsg += "\nSingular input. "
786            errmsg += "The following constraints are not "
787            errmsg += "linearly independent\n\tor do not "
788            errmsg += "allow for generation of a non-singular set\n"
789            errmsg += 'Please report this unexpected error\n'
790            for rel in group:
791                errmsg += _FormatConstraint(constrDict[rel],fixedList[rel])
792                errmsg += "\n"
793    _setVarLists([])
794    return errmsg,warnmsg
795
796def GenerateConstraints(varyList,constrDict,fixedList,parmDict=None,SeqHist=None):
797    '''Takes a list of relationship entries comprising a group of
798    constraints and builds the relationship lists and their inverse
799    and stores them in global parameters Also checks for internal
800    conflicts or inconsistencies in parameter/variable definitions.
801
802    :param list varyList: a list of parameters names (strings of form
803      ``<ph>:<hst>:<nam>``) that will be varied. Note that this is changed here.
804   
805    :param dict constrDict: a list of dicts defining relationships/constraints
806      (as defined in :func:`GroupConstraints`)
807
808    :param list fixedList: a list of values specifying a fixed value for each
809      dict in constrDict. Values are either strings that can be converted to
810      floats, float values or None if the constraint defines a new parameter.
811     
812    :param dict parmDict: a dict containing all parameters defined in current
813      refinement.
814
815    :param int SeqHist: number of current histogram, when used in a sequential
816      refinement. None (default) otherwise. Wildcard parameter names are
817      set to the current histogram, when found if not None.
818    '''
819    global dependentParmList,arrayList,invarrayList,indParmList,consNum
820    global genVarLookup
821    msg = ''
822    shortmsg = ''
823    changed = False
824   
825    # Process the equivalences
826    #    If there are conflicting parameters, move them into constraints. This
827    #    may create new conflicts, requiring movement of other parameters, so
828    #    loop until there are no more changes to make.
829    parmsChanged = True
830    while parmsChanged:
831        parmsChanged = 0
832        errmsg,warnmsg,fixVlist,dropVarList,translateTable = CheckEquivalences(
833            constrDict,varyList,parmDict,SeqHist)
834        if problemVars:
835            parmsChanged = MoveConfEquiv(constrDict,fixedList)
836            changed = True
837    if errmsg:
838        msg = errmsg
839    if warnmsg:
840        if msg: msg += '\n'
841        msg += warnmsg
842
843    # scan through parameters in each relationship. Are all varied? If only some are
844    # varied, create an error message.
845    groups,parmlist = GroupConstraints(constrDict)
846    for group,varlist in zip(groups,parmlist):
847        if len(varlist) == 1:   # process fixed (held) variables
848            var = varlist[0]
849            if var not in fixedVarList:
850                fixedVarList.append(var)
851            continue
852        for rel in group:
853            varied = 0
854            notvaried = ''
855            unused = 0
856            notused = ''
857            for var in constrDict[rel]:
858                if var.startswith('_'): continue
859                if var.split(':')[1] == '*' and SeqHist is not None:
860                    # convert wildcard var to reference current histogram; save translation in table
861                    sv = var.split(':')
862                    sv[1] = str(SeqHist)
863                    translateTable[var] = ':'.join(sv)
864                    var = translateTable[var]
865                if parmDict is not None and var not in parmDict:
866                    unused += 1
867                    if notvaried: notused += ', '
868                    notused += var
869                if var in varyList:
870                    varied += 1
871                else:
872                    if notvaried: notvaried += ', '
873                    notvaried += var
874                if var in fixedVarList:
875                    msg += '\nError: parameter '+var+" is Fixed and used in a constraint:\n\t"
876                    msg += _FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
877            #if unused > 0:# and unused != len(VarKeys(constrDict[rel])):
878            if unused > 0 and unused != len(VarKeys(constrDict[rel])):
879                #msg += "\nSome (but not all) parameters in constraint are not defined:\n\t"
880                #msg += _FormatConstraint(constrDict[rel],fixedList[rel])
881                #msg += '\nNot used: ' + notused + '\n'
882                shortmsg += notused+" not used in constraint "+_FormatConstraint(constrDict[rel],fixedList[rel])
883            elif varied > 0 and varied != len(VarKeys(constrDict[rel])):
884                #msg += "\nNot all parameters refined in constraint:\n\t"
885                #msg += _FormatConstraint(constrDict[rel],fixedList[rel])
886                #msg += '\nNot refined: ' + notvaried + '\n'
887                shortmsg += notvaried+" not varied in constraint "+_FormatConstraint(constrDict[rel],fixedList[rel])
888    # if there were errors found, go no farther
889    if shortmsg and SeqHist is not None:
890        if msg:
891            print (' *** ERROR in constraint definitions! ***')
892            print (msg)
893            raise ConstraintException
894        print ('*** Sequential refinement: ignoring constraint definition(s): ***')
895        print (shortmsg)
896        msg = ''
897    elif shortmsg:
898        msg += shortmsg
899    if msg:
900        print (' *** ERROR in constraint definitions! ***')
901        print (msg)
902        raise ConstraintException
903               
904    # now process each group and create the relations that are needed to form
905    # a non-singular square matrix
906    # If all are varied and this is a constraint equation, then set VaryFree flag
907    # so that the newly created relationships will be varied
908    for group,varlist in zip(groups,parmlist):
909        if len(varlist) == 1: continue
910        # for constraints, if all included parameters are refined,
911        # set the VaryFree flag, and remaining degrees of freedom will be
912        # varied (since consistency was checked, if any one parameter is
913        # refined, then assume that all are)
914        varsList = [] # make a list of all the referenced parameters as well
915        VaryFree = False
916        for rel in group:
917            varied = 0
918            unused = 0
919            for var in VarKeys(constrDict[rel]):
920                var = translateTable.get(var,var) # replace wildcards
921                if parmDict is not None and var not in parmDict:
922                    unused += 1                   
923                if var not in varsList: varsList.append(var)
924                if var in varyList: varied += 1
925            if fixedList[rel] is not None and varied > 0:
926                VaryFree = True
927        if len(varlist) < len(group): # too many relationships -- no can do
928            msg = 'too many relationships'
929            break
930        # Since we checked before, if any parameters are unused, then all must be.
931        # If so, this set of relationships can be ignored
932        if unused:
933            if debug: print('Constraint ignored (all parameters undefined)')
934            if debug: print ('    '+_FormatConstraint(constrDict[rel],fixedList[rel]))
935            continue
936        # fill in additional degrees of freedom
937        try:
938            arr = _FillArray(group,constrDict,varlist)
939            _RowEchelon(len(group),arr,varlist)
940            constrArr = _FillArray(group,constrDict,varlist,FillDiagonals=True)
941            GramSchmidtOrtho(constrArr,len(group))
942        except:
943            msg = 'Singular relationships'
944            break
945        mapvar = []
946        group = group[:]
947        # scan through all generated and input relationships, we need to add to the varied list
948        # all the new parameters where VaryFree has been set or where a New Var is varied.
949        #
950        # If a group does not contain any fixed values (constraint equations)
951        # and nothing in the group is varied, drop this group, so that the
952        # dependent parameters can be refined individually.
953        unused = True
954        for i in range(len(varlist)):
955            if len(group) > 0: # get the original equation reference
956                rel = group.pop(0)
957                fixedval = fixedList[rel]
958                varyflag = constrDict[rel].get('_vary',False)
959                varname = constrDict[rel].get('_name','')
960            else: # this relationship has been generated
961                varyflag = False
962                varname = ''
963                fixedval = None
964            if fixedval is None: # this is a new parameter, not a constraint
965                if not varname:
966                    varname = paramPrefix + str(consNum) # no assigned name, create one
967                    consNum += 1
968                mapvar.append(varname)
969                genVarLookup[varname] = varlist # save list of parameters related to this new var
970                # vary the new relationship if it is a degree of freedom in
971                # a set of contraint equations or if a New Var is flagged to be varied.
972                if VaryFree or varyflag: 
973                    unused = False
974                    varyList.append(varname)
975                    # fix (prevent varying) of all the parameters inside the constraint group
976                    # (dependent vars)
977                    for var in varsList:
978                        if var in varyList: varyList.remove(var)
979            else:
980                unused = False
981                mapvar.append(fixedval)
982        if unused: # skip over constraints that don't matter (w/o fixed value or any refined parameters)
983            if debug: print('Constraint ignored (all parameters unrefined)')
984            if debug: print ('   '+_FormatConstraint(constrDict[rel],fixedList[rel]))
985            continue 
986        dependentParmList.append([translateTable.get(var,var) for var in varlist])
987        arrayList.append(constrArr)
988        invarrayList.append(np.linalg.inv(constrArr))
989        indParmList.append(mapvar)
990        symGenList.append(False)
991    if msg:
992        print (' *** ERROR in constraint definitions! ***')
993        print (msg)
994        print (VarRemapShow(varyList))
995        raise ConstraintException
996    # setup dictionary containing the fixed values
997    global fixedDict
998    # key is original ascii string, value is float
999    for fixedval in fixedList:
1000        if fixedval:
1001            fixedDict[fixedval] = float(fixedval)
1002    _setVarLists(dropVarList)
1003    if changed:
1004        print(60*'=')
1005        print('Constraints were reclassified to avoid conflicts, as below:')
1006        print (VarRemapShow(varyList,True))
1007        print(60*'=')
1008    return groups,parmlist # saved for sequential fits
1009   
1010def _setVarLists(dropVarList):
1011    '''Make list of dependent and independent variables (after dropping unused vars in dropVarList)
1012    '''
1013    global dependentParmList,indParmList
1014    global dependentVars
1015    global independentVars
1016    dependentVars = []
1017    independentVars = []
1018    for varlist,mapvars in zip(dependentParmList,indParmList):  # process all constraints
1019        for mv in mapvars:
1020            if mv in dropVarList: continue
1021            if mv not in independentVars: independentVars.append(mv)
1022        for mv in varlist:
1023            if mv in dropVarList: continue
1024            if mv not in dependentVars: dependentVars.append(mv)
1025    if debug: # on debug, show what is parsed & generated, semi-readable
1026        print (50*'-')
1027        #print (VarRemapShow(varyList))
1028        #print ('Varied: ',varyList)
1029        print ('Not Varied: ',fixedVarList)
1030
1031# def CheckEquivalences(constrDict,varyList):
1032#     global dependentParmList,arrayList,invarrayList,indParmList,consNum
1033#     global problemVars
1034#     warnmsg = ''
1035#     errmsg = ''
1036#     problemVars = []
1037#     # process fixed variables (holds)
1038#     fixVlist = [] # list of Fixed vars
1039#     constrVars = [] # list of vars in constraint expressions
1040#     for cdict in constrDict:
1041#         # N.B. No "_" names in holds
1042#         if len(cdict) == 1:
1043#             fixVlist.append(list(cdict.keys())[0])
1044#         else:
1045#             constrVars += cdict.keys() # this will include _vary (not a problem)
1046#     # process equivalences: make a list of dependent and independent vars
1047#     #    and check for repeated uses (repetition of a parameter as an
1048#     #    independent var is OK)
1049#     indepVarList = []
1050#     depVarList = []
1051#     multdepVarList = []
1052#     for varlist,mapvars,multarr,invmultarr in zip(
1053#         dependentParmList,indParmList,arrayList,invarrayList):
1054#         if multarr is None: # an equivalence
1055#             zeromult = False
1056#             for mv in mapvars:
1057#                 varied = 0
1058#                 notvaried = ''
1059#                 if mv in varyList:
1060#                     varied += 1
1061#                 else:
1062#                     if notvaried: notvaried += ', '
1063#                     notvaried += mv
1064#                 if mv not in indepVarList: indepVarList.append(mv)
1065#                 for v,m in zip(varlist,invmultarr):
1066#                     if v in indepVarList:
1067#                         errmsg += '\nVariable '+v+' is used to set values in a constraint before its value is set in another constraint\n'
1068#                         if v not in problemVars: problemVars.append(v)
1069#                     if m == 0: zeromult = True
1070#                     if v in varyList:
1071#                         varied += 1
1072#                     else:
1073#                         if notvaried: notvaried += ', '
1074#                         notvaried += v
1075#                     if v in depVarList:
1076#                         multdepVarList.append(v)
1077#                     else:
1078#                         depVarList.append(v)
1079#             if varied > 0 and varied != len(varlist)+1:
1080#                 warnmsg += "\nNot all variables refined in equivalence:\n\t"
1081#                 s = ""
1082#                 for v in varlist:
1083#                     if s != "": s+= " & "
1084#                     s += str(v)           
1085#                 warnmsg += str(mv) + " => " + s
1086#                 warnmsg += '\nNot refined: ' + notvaried + '\n'
1087#             if zeromult:
1088#                 errmsg += "\nZero multiplier is invalid in equivalence:\n\t"
1089#                 s = ""
1090#                 for v in varlist:
1091#                     if s != "": s+= " & "
1092#                     s += str(v)           
1093#                 errmsg += str(mv) + " => " + s + '\n'
1094#     # check for errors:
1095#     if len(multdepVarList) > 0:
1096#         errmsg += "\nThe following parameters(s) are used in conflicting Equivalence relations as dependent variables:\n"
1097#         s = ''
1098#         for var in sorted(set(multdepVarList)):
1099#             if v not in problemVars: problemVars.append(v)
1100#             if s != "": s+= ", "
1101#             s += str(var)           
1102#         errmsg += '\t'+ s + '\n'
1103#     equivVarList = list(set(indepVarList).union(set(depVarList)))
1104#     if debug: print ('equivVarList',equivVarList)
1105#     # check for parameters that are both fixed and in an equivalence (not likely)
1106#     inboth = set(fixVlist).intersection(set(equivVarList))
1107#     if len(inboth) > 0:
1108#         errmsg += "\nThe following parameter(s) are used in both Equivalence and Fixed constraints:\n"
1109#         s = ''
1110#         for var in sorted(inboth):
1111#             if var not in problemVars: problemVars.append(var)
1112#             if s != "": s+= ", "
1113#             s += str(var)
1114#         errmsg += '\t'+ s + '\n'
1115#     # check for parameters that in both an equivalence and a constraint expression
1116#     inboth = set(constrVars).intersection(set(equivVarList))
1117#     if len(inboth) > 0:
1118#         errmsg += "\nThe following parameter(s) are used in both Equivalence and Equiv or new var constraints:\n"
1119#         s = ''
1120#         for var in sorted(inboth):
1121#             if var not in problemVars: problemVars.append(var)
1122#             if s != "": s+= ", "
1123#             s += str(var)
1124#         errmsg += '\t'+ s + '\n'
1125#     return errmsg,warnmsg,fixVlist
1126
1127def CheckEquivalences(constrDict,varyList,parmDict=None,SeqHist=None):
1128    '''Process equivalence constraints, looking for conflicts such as
1129    where a parameter is used in both an equivalence and a constraint expression
1130    or where chaining is done (A->B and B->C).
1131    When called during refinements, parmDict is defined, and for sequential refinement
1132    SeqHist ia also defined.
1133
1134      * parmDict is used to remove equivalences where a parameter is not present
1135        in a refinement
1136      * SeqHist is used to rename wild-card parameter names in sequential
1137        refinements to use the current histogram.
1138    '''
1139    global dependentParmList,arrayList,invarrayList,indParmList,consNum
1140    global problemVars
1141    warnmsg = ''
1142    errmsg = ''
1143    problemVars = []
1144    # process fixed parameters (holds)
1145    fixVlist = [] # list of Fixed vars
1146    constrVars = [] # list of vars in constraint expressions
1147    for cdict in constrDict:
1148        # N.B. No "_" names in holds
1149        if len(cdict) == 1:
1150            fixVlist.append(list(cdict.keys())[0])
1151        else:
1152            constrVars += cdict.keys() # this will include _vary (not a problem)
1153    # process equivalences: make a list of dependent and independent vars
1154    #    and check for repeated uses (repetition of a parameter as an
1155    #    independent var is OK)
1156    indepVarList = []
1157    depVarList = []
1158    multdepVarList = []
1159    dropVarList = []
1160    translateTable = {} # lookup table for wildcard referenced parameters
1161    for varlist,mapvars,multarr,invmultarr in zip(
1162        dependentParmList,indParmList,arrayList,invarrayList):
1163        if multarr is None: # an equivalence
1164            zeromult = False
1165            for i,mv in enumerate(mapvars):
1166                if mv.split(':')[1] == '*' and SeqHist is not None:
1167                    # convert wildcard var to reference current histogram; save translation in table
1168                    sv = mv.split(':')
1169                    sv[1] = str(SeqHist)
1170                    mv = translateTable[mv] = ':'.join(sv)
1171                    mapvars[i] = mv
1172                varied = 0
1173                notvaried = ''
1174                if mv in varyList:
1175                    varied += 1
1176                else:
1177                    if notvaried: notvaried += ', '
1178                    notvaried += mv
1179                if parmDict is not None and mv not in parmDict:
1180                    print ("Dropping equivalence for parameter "+str(mv)+". Not defined in this refinement")
1181                    if mv not in dropVarList: dropVarList.append(mv)
1182                if mv not in indepVarList: indepVarList.append(mv)
1183            for i,(v,m) in enumerate(zip(varlist,invmultarr)):
1184                if v.split(':')[1] == '*' and SeqHist is not None:
1185                    # convert wildcard var to reference current histogram; save translation in table
1186                    sv = v.split(':')
1187                    sv[1] = str(SeqHist)
1188                    varlist[i] = v = translateTable[v] = ':'.join(sv)
1189                if parmDict is not None and v not in parmDict:
1190                    print ("Dropping equivalence for dep. variable "+str(v)+". Not defined in this refinement")
1191                    if v not in dropVarList: dropVarList.append(v)
1192                    continue
1193                if m == 0: zeromult = True
1194                if v in varyList:
1195                    varied += 1
1196                else:
1197                    if notvaried: notvaried += ', '
1198                    notvaried += v
1199                if v in indepVarList:
1200                    errmsg += '\nParameter '+v+' is used to set values in a constraint before its value is set in another constraint\n'
1201                    if v not in problemVars: problemVars.append(v)
1202                if v in depVarList:
1203                    multdepVarList.append(v)
1204                else:
1205                    depVarList.append(v)
1206            if varied > 0 and varied != len(varlist)+1:
1207                warnmsg += "\nNot all parameters refined in equivalence:\n\t"
1208                s = ""
1209                for v in varlist:
1210                    if s != "": s+= " & "
1211                    s += str(v)           
1212                warnmsg += str(mv) + " => " + s
1213                warnmsg += '\nNot refined: ' + notvaried + '\n'
1214            if zeromult:
1215                errmsg += "\nZero multiplier is invalid in equivalence:\n\t"
1216                s = ""
1217                for v in varlist:
1218                    if s != "": s+= " & "
1219                    s += str(v)           
1220                errmsg += str(mv) + " => " + s + '\n'
1221    # check for errors:
1222    if len(multdepVarList) > 0:
1223        errmsg += "\nThe following parameters(s) are used in conflicting Equivalence relations as dependent variables:\n"
1224        s = ''
1225        for var in sorted(set(multdepVarList)):
1226            if v not in problemVars: problemVars.append(v)
1227            if s != "": s+= ", "
1228            s += str(var)           
1229        errmsg += '\t'+ s + '\n'
1230    equivVarList = list(set(indepVarList).union(set(depVarList)))
1231    if debug: print ('equivVarList',equivVarList)
1232    # check for parameters that are both fixed and in an equivalence (not likely)
1233    inboth = set(fixVlist).intersection(set(equivVarList))
1234    if len(inboth) > 0:
1235        errmsg += "\nThe following parameter(s) are used in both Equivalence and Fixed constraints:\n"
1236        s = ''
1237        for var in sorted(inboth):
1238            if var not in problemVars: problemVars.append(var)
1239            if s != "": s+= ", "
1240            s += str(var)
1241        errmsg += '\t'+ s + '\n'
1242    # check for parameters that in both an equivalence and a constraint expression
1243    inboth = set(constrVars).intersection(set(equivVarList))
1244    if len(inboth) > 0:
1245        errmsg += "\nThe following parameter(s) are used in both Equivalence and Equiv or new var constraints:\n"
1246        s = ''
1247        for var in sorted(inboth):
1248            if var not in problemVars: problemVars.append(var)
1249            if s != "": s+= ", "
1250            s += str(var)
1251        errmsg += '\t'+ s + '\n'
1252    return errmsg,warnmsg,fixVlist,dropVarList,translateTable
1253
1254def MoveConfEquiv(constrDict,fixedList):
1255    '''Address conflicts in Equivalence constraints by creating an constraint equation
1256    that has the same action as the equivalence and removing the Equivalence
1257    '''
1258    global dependentParmList,arrayList,invarrayList,indParmList,consNum
1259    global problemVars
1260    parmsChanged = 0
1261    for i,(varlist,mapvars) in enumerate(zip(dependentParmList,indParmList)):
1262        conf = False
1263        for mv in mapvars:
1264            if mv in problemVars:
1265                conf = True
1266                break
1267        for v in varlist:
1268            if v in problemVars:
1269                conf = True
1270                break
1271        if conf:
1272            parmsChanged += 1
1273            indvar = indParmList[i][0]
1274            for dep,mult in zip(dependentParmList[i],invarrayList[i]):
1275                #print('debug replacing equiv with constraint equation 0=',{indvar:-1.,dep:mult[0]})
1276                constrDict += [{indvar:-1.,dep:mult[0]}]
1277                fixedList += ['0.0']
1278            dependentParmList[i] = None
1279    if parmsChanged:
1280        for i in range(len(dependentParmList)-1,-1,-1):
1281            if dependentParmList[i] is None:
1282                del dependentParmList[i],indParmList[i],arrayList[i],invarrayList[i]
1283    return parmsChanged
1284
1285def StoreEquivalence(independentVar,dependentList,symGen=True):
1286    '''Takes a list of dependent parameter(s) and stores their
1287    relationship to a single independent parameter (independentVar).
1288
1289    Called with user-supplied constraints by :func:`GSASIIstrIO.ProcessConstraints,
1290    with Pawley constraints from :func:`GSASIIstrIO.GetPawleyConstr`,
1291    with Unit Cell constraints from :func:`GSASIIstrIO.cellVary`
1292    with symmetry-generated atom constraints from :func:`GSASIIstrIO.GetPhaseData`
1293
1294    :param str independentVar: name of master parameter that will be used to determine the value
1295      to set the dependent variables
1296
1297    :param list dependentList: a list of parameters that will set from
1298         independentVar. Each item in the list can be a string with the parameter
1299         name or a tuple containing a name and multiplier:
1300         ``['::parm1',('::parm2',.5),]``
1301
1302    '''
1303   
1304    global dependentParmList,arrayList,invarrayList,indParmList,symGenList
1305    mapList = []
1306    multlist = []
1307    allfloat = True
1308    for var in dependentList:
1309        if isinstance(var, str):
1310            mult = 1.0
1311        elif len(var) == 2:
1312            var,mult = var
1313        else:
1314            raise Exception("Cannot parse "+repr(var) + " as var or (var,multiplier)")
1315        mapList.append(var)
1316        try:           
1317            multlist.append(tuple((float(mult),)))
1318        except:
1319            allfloat = False
1320            multlist.append(tuple((mult,)))
1321    # added relationships to stored values
1322    arrayList.append(None)
1323    if allfloat:
1324        invarrayList.append(np.array(multlist))
1325    else:
1326        invarrayList.append(multlist)
1327    indParmList.append(list((independentVar,)))
1328    dependentParmList.append(mapList)
1329    symGenList.append(symGen)
1330    return
1331
1332def EvaluateMultipliers(constList,*dicts):
1333    '''Convert multipliers for constraints and equivalences that are specified
1334    as strings into values. The strings can specify values in the parameter dicts as
1335    well as normal Python functions, such as "2*np.cos(0::Ax:2/2.)"
1336   
1337    :param list constList: a list of dicts containing constraint expressions
1338    :param \*dicts: one or more dicts containing GSAS-II parameters and their values
1339       can be specified
1340    :returns: an empty string if there were no errors, or an error message listing
1341       the strings that could not be converted.
1342    '''
1343    def SubfromParmDict(s,prmDict):
1344        for key in prmDict:
1345            if key in s:
1346                s = s.replace(key,str(prmDict[key]))
1347        return eval(s)
1348    prmDict = {}
1349    for d in dicts: prmDict.update(d) # combine all passed parameter dicts
1350    problemList = ""
1351    # loop through multipliers in contraint expressions
1352    for const in constList:
1353        for key in const:
1354            if key.startswith('_'): continue
1355            try: # is this already a float, etc?
1356                1+const[key]
1357                continue
1358            except:
1359                pass
1360            try:
1361                newval = SubfromParmDict(const[key][:],prmDict)
1362                if GSASIIpath.GetConfigValue('debug'):
1363                    print('Changing ',const[key],'to',newval)
1364                const[key] = newval               
1365            except:
1366                if problemList: problemList += ", "
1367                problemList += const[key]
1368    # loop through multipliers in equivalences
1369    global arrayList,invarrayList
1370    for i,(a,valList) in enumerate(zip(arrayList,invarrayList)):
1371        if a is not None: continue # ignore if not equiv
1372        try:
1373            valList.shape
1374            continue # ignore if already a numpy array
1375        except:
1376            pass
1377        repList = []
1378        for v in valList:
1379            try:
1380                1+v[0]
1381                repList.append(tuple((v[0],)))
1382                continue
1383            except:
1384                pass
1385            try:
1386                newval = SubfromParmDict(v[0][:],prmDict)
1387                if GSASIIpath.GetConfigValue('debug'):
1388                    print('Changing ',v[0],'to',newval)
1389                repList.append(tuple((newval,)))
1390            except:
1391                if problemList: problemList += ", "
1392                problemList += v[0]
1393                repList.append(tuple(('error',)))
1394        invarrayList[i] = np.array(repList)
1395    return problemList
1396
1397def GetDependentVars():
1398    '''Return a list of dependent variables: e.g. parameters that are
1399    constrained in terms of other parameters
1400
1401    :returns: a list of parameter names
1402
1403    '''
1404    global dependentVars
1405    return dependentVars
1406
1407def GetIndependentVars():
1408    '''Return a list of independent variables: e.g. parameters that are
1409    slaved to other parameters by constraints
1410
1411    :returns: a list of parameter names
1412
1413    '''
1414    global independentVars
1415    return independentVars
1416
1417def PrintIndependentVars(parmDict,varyList,sigDict,PrintAll=False,pFile=None):
1418    '''Print the values and uncertainties on the independent parameters'''
1419    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict
1420    printlist = []
1421    mapvars = GetIndependentVars()
1422    for i,name in enumerate(mapvars):
1423        if name in fixedDict: continue
1424        if PrintAll or name in varyList:
1425            sig = sigDict.get(name)
1426            printlist.append([name,parmDict[name],sig])
1427    if len(printlist) == 0: return
1428    s1 = ''
1429    s2 = ''
1430    s3 = ''
1431    pFile.write(130*'-'+'\n')
1432    pFile.write("Parameters generated by constraints\n")
1433    printlist.append(3*[None])
1434    for name,val,esd in printlist:
1435        if len(s1) > 120 or name is None:
1436            pFile.write(''+'\n')
1437            pFile.write(s1+'\n')
1438            pFile.write(s2+'\n')
1439            pFile.write(s3+'\n')
1440            s1 = ''
1441            if name is None: break
1442        if s1 == "":
1443            s1 = ' name  :'
1444            s2 = ' value :'
1445            s3 = ' sig   :'
1446        s1 += '%15s' % (name)
1447        s2 += '%15.4f' % (val)
1448        if esd is None:
1449            s3 += '%15s' % ('n/a')
1450        else:   
1451            s3 += '%15.4f' % (esd)
1452
1453def ComputeDepESD(covMatrix,varyList,parmDict):
1454    '''Compute uncertainties for dependent parameters from independent ones
1455    returns a dictionary containing the esd values for dependent parameters
1456    '''
1457    sigmaDict = {}
1458    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
1459        #if invmultarr is None: continue # probably not needed
1460#        try:
1461#            valuelist = [parmDict[var] for var in mapvars]
1462#        except KeyError:
1463#            continue
1464        # get the v-covar matrix for independent parameters
1465        vcov = np.zeros((len(mapvars),len(mapvars)))
1466        for i1,name1 in enumerate(mapvars):
1467            if name1 not in varyList: continue
1468            iv1 = varyList.index(name1)
1469            for i2,name2 in enumerate(mapvars):
1470                if name2 not in varyList: continue
1471                iv2 = varyList.index(name2)
1472                vcov[i1][i2] = covMatrix[iv1][iv2]
1473        # vec is the vector that multiplies each of the independent values
1474        for v,vec in zip(varlist,invmultarr):
1475            sigmaDict[v] = np.sqrt(np.inner(vec.T,np.inner(vcov,vec)))
1476    return sigmaDict
1477
1478def _FormatConstraint(RelDict,RelVal):
1479    '''Formats a Constraint or Function for use in a convenient way'''
1480    linelen = 45
1481    s = [""]
1482    for var,val in RelDict.items():
1483        if var.startswith('_'): continue
1484        if len(s[-1]) > linelen: s.append(' ')
1485        m = val
1486        if s[-1] != "" and m >= 0:
1487            s[-1] += ' + '
1488        elif s[-1] != "":
1489            s[-1] += ' - '
1490            m = abs(m)
1491        s[-1] += '%.3f*%s '%(m,var)
1492    if len(s[-1]) > linelen: s.append(' ')
1493    if RelVal is None:
1494        s[-1] += ' = New variable'
1495    else:
1496        s[-1] += ' = ' + RelVal
1497    s1 = ''
1498    for s2 in s:
1499        if s1 != '': s1 += '\n\t'
1500        s1 += s2
1501    return s1
1502
1503def VarRemapShow(varyList,inputOnly=False):
1504    '''List out the saved relationships. This should be done after the constraints have been
1505    defined using :func:`StoreEquivalence`, :func:`GroupConstraints` and :func:`GenerateConstraints`.
1506
1507    :returns: a string containing the details of the contraint relationships
1508    '''
1509    s = ''
1510    if len(fixedVarList) > 0:
1511        s += 'Fixed Parameters:\n'
1512        for v in fixedVarList:
1513            s += '    ' + v + '\n'
1514    if not inputOnly:
1515        s += 'User-supplied parameter mapping relations:\n'
1516    symout = ''
1517    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict,symGenList
1518
1519    for varlist,mapvars,multarr,invmultarr,symFlag in zip(
1520        dependentParmList,indParmList,arrayList,invarrayList,symGenList):
1521        for i,mv in enumerate(mapvars):
1522            if multarr is None:
1523#                s1 = '  ' + str(mv) + ' is equivalent to parameter(s): '
1524                if len(varlist) == 1:
1525                    s1 = '   ' + str(mv) + ' is equivalent to '
1526                else:
1527                    s1 = '   ' + str(mv) + ' is equivalent to parameters: '
1528                j = 0
1529                for v,m in zip(varlist,invmultarr):
1530                    if debug: print ('v,m[0]: ',v,m[0])
1531                    if len(s1.split('\n')[-1]) > 75: s1 += '\n        '
1532                    if j > 0: s1 += ' & '
1533                    j += 1
1534                    s1 += str(v)
1535                    if m != 1:
1536                        s1 += " / " + str(m[0])
1537                if symFlag:
1538                    symout += s1 + '\n'
1539                else:
1540                    s += s1 + '\n'
1541                continue
1542            s += %s = ' % mv
1543            j = 0
1544            for m,v in zip(multarr[i,:],varlist):
1545                if m == 0: continue
1546                if j > 0: s += ' + '
1547                j += 1
1548                s += '(%s * %s)' % (m,v)
1549            if mv in varyList: s += ' VARY'
1550            s += '\n'
1551    if symout:
1552        s += 'Symmetry-generated relations:\n' + symout
1553    if inputOnly: return s
1554    s += 'Inverse parameter mapping relations:\n'
1555    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
1556        for i,mv in enumerate(varlist):
1557            s += %s = ' % mv
1558            j = 0
1559            for m,v in zip(invmultarr[i,:],mapvars):
1560                if m == 0: continue
1561                if j > 0: s += ' + '
1562                j += 1
1563                s += '(%s * %s)' % (m,v)
1564            s += '\n'
1565    return s
1566
1567def GetSymEquiv():
1568    '''Return the automatically generated (equivalence) relationships.
1569
1570    :returns: a list of strings containing the details of the contraint relationships
1571    '''
1572    symout = []
1573    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict,symGenList
1574
1575    for varlist,mapvars,multarr,invmultarr,symFlag in zip(
1576        dependentParmList,indParmList,arrayList,invarrayList,symGenList):
1577        for i,mv in enumerate(mapvars):
1578            if not symFlag: continue
1579            if multarr is None:
1580                #s1 = str(mv) + ' = '
1581                s1 = ''
1582                s2 = ' = ' + str(mv)
1583                j = 0
1584                for v,m in zip(varlist,invmultarr):
1585                    if debug: print ('v,m[0]: ',v,m[0])
1586                    if len(s1.split('\n')[-1]) > 75: s1 += '\n        '
1587                    if j > 0: s1 += ' =  '
1588                    j += 1
1589                    s1 += str(v)
1590                    if m != 1:
1591                        s1 += " / " + str(m[0])
1592                symout.append(s1+s2)
1593                continue
1594            else:
1595                s = %s = ' % mv
1596                j = 0
1597                for m,v in zip(multarr[i,:],varlist):
1598                    if m == 0: continue
1599                    if j > 0: s += ' + '
1600                    j += 1
1601                    s += '(%s * %s)' % (m,v)
1602                print ('unexpected sym op='+s)
1603    return symout
1604
1605def Dict2Deriv(varyList,derivDict,dMdv):
1606    '''Compute derivatives for Independent Parameters from the
1607    derivatives for the original parameters
1608
1609    :param list varyList: a list of parameters names that will be varied
1610
1611    :param dict derivDict: a dict containing derivatives for parameter values keyed by the
1612      parameter names.
1613
1614    :param list dMdv: a Jacobian, as a list of np.array containing derivatives for dependent
1615      parameter computed from derivDict
1616
1617    '''
1618    global dependentParmList,arrayList,invarrayList,indParmList,invarrayList
1619    for varlist,mapvars,multarr,invmultarr in zip(dependentParmList,indParmList,arrayList,invarrayList):
1620        for i,name in enumerate(mapvars):
1621            # grouped parameters: need to add in the derv. w/r
1622            # dependent variables to the independent ones
1623            if name not in varyList: continue # skip if independent var not varied
1624            if multarr is None:
1625                for v,m in zip(varlist,invmultarr):
1626                    if debug: print ('start dMdv',dMdv[varyList.index(name)])
1627                    if debug: print ('add derv',v,'/',m[0],'to derv',name,'add=',derivDict[v] / m[0])
1628                    if m == 0: continue
1629                    dMdv[varyList.index(name)] += derivDict[v] / m[0]
1630            else:
1631                for v,m in zip(varlist,multarr[i,:]):
1632                    if debug: print ('start dMdv',dMdv[varyList.index(name)])
1633                    if debug: print ('add derv',v,'*',m,'to derv',name,'add=',m * derivDict[v])
1634                    if m == 0: continue
1635                    dMdv[varyList.index(name)] += m * derivDict[v]
1636
1637def Map2Dict(parmDict,varyList):
1638    '''Create (or update) the Independent Parameters from the original
1639    set of Parameters
1640
1641    Removes dependent variables from the varyList
1642
1643    This should be done once, after the constraints have been
1644    defined using :func:`StoreEquivalence`,
1645    :func:`GroupConstraints` and :func:`GenerateConstraints` and
1646    before any parameter refinement is done. This completes the parameter
1647    dictionary by defining independent parameters and it satisfies the
1648    constraint equations in the initial parameters
1649
1650    :param dict parmDict: a dict containing parameter values keyed by the
1651      parameter names.
1652      This will contain updated values for both dependent and independent
1653      parameters after Dict2Map is called. It will also contain some
1654      unexpected entries of every constant value {'0':0.0} & {'1.0':1.0},
1655      which do not cause any problems.
1656
1657    :param list varyList: a list of parameters names that will be varied
1658   
1659
1660    '''
1661    # process the Independent vars: remove dependent ones from varylist
1662    # and then compute values for the Independent ones from their dependents
1663    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict
1664    for varlist,mapvars,multarr in zip(dependentParmList,indParmList,arrayList):
1665        for item in varlist:
1666            try:
1667                varyList.remove(item)
1668            except ValueError:
1669                pass
1670        if multarr is None: continue
1671        valuelist = [parmDict[var] for var in varlist]
1672        parmDict.update(zip(mapvars,
1673                            np.dot(multarr,np.array(valuelist)))
1674                        )
1675    # now remove fixed parameters from the varyList
1676    global fixedVarList
1677    for item in fixedVarList:
1678        try:
1679            varyList.remove(item)
1680        except ValueError:
1681            pass
1682    # Set constrained parameters to their fixed values
1683    parmDict.update(fixedDict)
1684
1685def Dict2Map(parmDict,varyList):
1686    '''Applies the constraints defined using :func:`StoreEquivalence`,
1687    :func:`GroupConstraints` and :func:`GenerateConstraints` by changing
1688    values in a dict containing the parameters. This should be
1689    done before the parameters are used for any computations
1690
1691    :param dict parmDict: a dict containing parameter values keyed by the
1692      parameter names.
1693      This will contain updated values for both dependent and independent
1694      parameters after Dict2Map is called. It will also contain some
1695      unexpected entries of every constant value {'0':0.0} & {'1.0':1.0},
1696      which do not cause any problems.
1697
1698    :param list varyList: a list of parameters names that will be varied
1699   
1700    '''
1701    global dependentParmList,arrayList,invarrayList,indParmList,fixedDict
1702    # reset fixed values (should not be needed, but very quick)
1703    # - this seems to update parmDict with {'0':0.0} & {'1.0':1.0} - probably not what was intended
1704    # not needed, but also not a problem - BHT
1705    parmDict.update(fixedDict)
1706    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
1707        #if invmultarr is None: continue
1708        try: 
1709            valuelist = [parmDict[var] for var in mapvars]
1710        except KeyError:
1711            continue
1712        parmDict.update(zip(varlist,np.dot(invmultarr,np.array(valuelist))))
1713
1714#======================================================================
1715# internal routines follow (these routines are unlikely to be called
1716# from outside the module)
1717
1718def GramSchmidtOrtho(a,nkeep=0):
1719    '''Use the Gram-Schmidt process (http://en.wikipedia.org/wiki/Gram-Schmidt) to
1720    find orthonormal unit vectors relative to first row.
1721
1722    If nkeep is non-zero, the first nkeep rows in the array are not changed
1723   
1724    input:
1725       arrayin: a 2-D non-singular square array
1726    returns:
1727       a orthonormal set of unit vectors as a square array
1728    '''
1729    def proj(a,b):
1730        'Projection operator'
1731        return a*(np.dot(a,b)/np.dot(a,a))
1732    for j in range(nkeep,len(a)):
1733        for i in range(j):
1734            a[j] -= proj(a[i],a[j])
1735        if np.allclose(np.linalg.norm(a[j]),0.0):
1736            raise ConstraintException("Singular input to GramSchmidtOrtho")
1737        a[j] /= np.linalg.norm(a[j])
1738    return a
1739
1740def _FillArray(sel,dict,collist,FillDiagonals=False):
1741    '''Construct a n by n matrix (n = len(collist)
1742    filling in the rows using the relationships defined
1743    in the dictionaries selected by sel
1744
1745    If FillDiagonals is True, diagonal elements in the
1746    array are set to 1.0
1747    '''
1748    n = len(collist)
1749    if FillDiagonals:
1750        arr = np.eye(n)
1751    else:
1752        arr = np.zeros(2*[n,])
1753    # fill the top rows
1754    for i,cnum in enumerate(sel):
1755        for j,var in enumerate(collist):
1756            arr[i,j] = dict[cnum].get(var,0)
1757    return arr
1758
1759def _SwapColumns(i,m,v):
1760    '''Swap columns in matrix m as well as the labels in v
1761    so that element (i,i) is replaced by the first non-zero element in row i after that element
1762
1763    Throws an exception if there are no non-zero elements in that row
1764    '''
1765    for j in range(i+1,len(v)):
1766        if not np.allclose(m[i,j],0):
1767            m[:,(i,j)] = m[:,(j,i)]
1768            v[i],v[j] = v[j],v[i]
1769            return
1770    else:
1771        raise ConstraintException('Singular input')
1772
1773def _RowEchelon(m,arr,collist):
1774    '''Convert the first m rows in Matrix arr to row-echelon form
1775    exchanging columns in the matrix and collist as needed.
1776
1777    throws an exception if the matrix is singular because
1778    the first m rows are not linearly independent
1779    '''
1780    for i in range(m):
1781        if np.allclose(arr[i,i],0):
1782            _SwapColumns(i,arr,collist)
1783        arr[i,:] /= arr[i,i] # normalize row
1784        # subtract current row from subsequent rows to set values to left of diagonal to 0
1785        for j in range(i+1,m):
1786            arr[j,:] -= arr[i,:] * arr[j,i]
1787
1788if __name__ == "__main__":
1789    parmdict = {}
1790    constrDict = [
1791        {'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},
1792        {'0:0:eA': 0.0},
1793        {'2::C(10,6,1)': 1.0, '1::C(10,6,1)': 1.0},
1794        {'1::C(10,0,1)': 1.0, '2::C(10,0,1)': 1.0},
1795        {'1::AUiso:0': 1.0, '0::AUiso:0': 1.0},
1796        {'0::A0': 0.0}
1797        ]
1798    fixedList = ['5.0', '0', None, None, '1.0', '0']
1799    StoreEquivalence('2::atomx:3',('2::atomy:3', ('2::atomz:3',2,), ))
1800    #StoreEquivalence('1::atomx:3',('2::atomx:3', ('2::atomz:3',2,), )) # error: dependent & independent vars mixed
1801    #StoreEquivalence('1::atomx:3',('2::atomy:3', ('2::atomz:3',2,), )) # error: dependent vars repeated
1802    #StoreEquivalence('0:1:eA',('0:0:eA',)) # error: equiv & fixed
1803    #StoreEquivalence('0:99:Scale',('0:12:Scale',)) # error: equiv & constrained
1804    #StoreEquivalence('0:12:Scale',('0:99:Scale',)) # error: equiv & constrained
1805    varylist = ['2::atomx:3',
1806                '2::C(10,6,1)', '1::C(10,6,1)',
1807                '2::atomy:3', '2::atomz:3',
1808                '0:12:Scale', '0:11:Scale', '0:14:Scale', '0:13:Scale', '0:0:Scale']
1809#    e,w = CheckConstraints([,
1810#                     [{'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}],
1811#                     ['1.0'])
1812#    if e: print 'error=',e
1813#    if w: print 'error=',w
1814#    varyList = ['0::A0', '0::AUiso:0', '0::Afrac:1', '0::Afrac:2', '0::Afrac:3', '0::Afrac:4',
1815#       '0::dAx:5', '0::dAy:5', '0::dAz:5', '0::AUiso:5', ':0:Back;0', ':0:Back;1', ':0:Back;2', ':0:Back;3',
1816#       ':0:Back;4', ':0:Back;5', ':0:Back;6', ':0:Back;7', ':0:Back;8', ':0:Back;9', ':0:Back;10', ':0:Back;11'
1817#       :0:U', ':0:V', ':0:W', ':0:X', ':0:Y', ':0:Scale', ':0:DisplaceX', ':0:DisplaceY']
1818#    constrDict = [
1819#        {'0::Afrac:4': 24.0, '0::Afrac:1': 16.0, '0::Afrac:3': 24.0, '0::Afrac:2': 16.0},
1820#        {'0::Afrac:1': 1.0, '0::Afrac:2': 1.0},
1821#        {'0::Afrac:4': 1.0, '0::Afrac:3': 1.0}]
1822#    fixedList = ['40.0', '1.0', '1.0']
1823
1824    errmsg, warnmsg = CheckConstraints(varylist,constrDict,fixedList)
1825    if errmsg:
1826        print ("*** Error ********************")
1827        print (errmsg)
1828    if warnmsg:
1829        print ("*** Warning ********************")
1830        print (warnmsg)
1831    if errmsg or warnmsg:
1832        sys.exit()
1833    groups,parmlist = GroupConstraints(constrDict)
1834    GenerateConstraints(groups,parmlist,varylist,constrDict,fixedList)
1835    print (VarRemapShow(varylist))
1836    parmdict.update( {
1837        '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,
1838        '0:0:eA': 0.0,
1839        '2::C(10,6,1)': 0.2, '1::C(10,6,1)': 0.3,
1840        '1::C(10,0,1)': 0.2, '2::C(10,0,1)': 0.3,
1841        '1::AUiso:0': 0.02, '0::AUiso:0': 0.03,
1842        '0::A0': 0.0,
1843        '2::atomx:3':0.23,'2::atomy:3':-.23, '2::atomz:3':-0.11,
1844        })
1845    print ('parmdict start',parmdict)
1846    print ('varylist start',varylist)
1847    before = parmdict.copy()
1848    Map2Dict(parmdict,varylist)
1849    print ('parmdict before and after Map2Dict')
1850    print ('  key / before / after')
1851    for key in sorted(list(parmdict.keys())):
1852        print ('  '+key,'\t',before.get(key),'\t',parmdict[key])
1853    print ('varylist after',varylist)
1854    before = parmdict.copy()
1855    Dict2Map(parmdict,varylist)
1856    print ('after Dict2Map')
1857    print ('  key / before / after')
1858    for key in sorted(list(parmdict.keys())):
1859        print ('  '+key,'\t',before.get(key),'\t',parmdict[key])
1860#    dMdv = len(varylist)*[0]
1861#    deriv = {}
1862#    for i,v in enumerate(parmdict.keys()): deriv[v]=i
1863#    Dict2Deriv(varylist,deriv,dMdv)
Note: See TracBrowser for help on using the repository browser.