# TODO: revisit SeqRefine and :meth:`GSASIIdataGUI.GSASII.OnSeqRefine` and :func:`GSASIIseqGUI.UpdateSeqResults`

# -*- coding: utf-8 -*-
"""
Classes and routines defined in :mod:`GSASIImapvars` follow. 

Note that parameter names in GSAS-II are strings of form ``<ph#>:<hst#>:<nam>`` or ``<ph#>::<nam>:<at#>``
where ``<ph#>`` is a phase number, ``<hst#>`` is a histogram number and ``<at#>`` is an atom number. 
``<nam>`` is a name that determines the parameter type (see :func:`GSASIIobj.CompileVarDesc`). When 
stored in the data tree, parameters are saved as :class:`GSASIIobj.G2VarObj` objects 
so that they can be resolved if the phase/histogram order changes.
"""
# Note that documentation for GSASIImapvars.py has been moved
# to file docs/source/GSASIImapvars.rst

from __future__ import division, print_function
import copy
import numpy as np
from . import GSASIIpath
from . import GSASIIobj as G2obj 
# data used for constraints; 
debug = False # turns on printing as constraint input is processed

#------------------------------------------------------------------------------------
# Global vars used for storing constraints and equivalences after processing
#   note that new var and constraint equations are stored together in groups,
#   where each constraint group contains those parameters that must be handled 
#   together. Equivalences are also stored in these 

dependentParmList = []
'''a list of lists where each item contains a list of parameters in each constraint group. 
note that parameters listed in dependentParmList should not be refined directly.'''
indParmList = [] # a list of names for the new parameters
'''a list of lists where each item contains a list for each constraint group with 
fixed values for constraint equations and names of generated/New Var parameters.
In the case of equivalences, the name of a single independent parameter is stored.
'''
arrayList = []
'''a list of of relationship matrices that map model parameters in each 
constraint group (in :data:`dependentParmList`) to 
generated (New Var) parameters.
'''
invarrayList = []
'''a list of of inverse-relationship matrices that map constrained values and 
generated (New Var) parameters (in :data:`indParmList`) to model parameters
(in :data:`dependentParmList`). 
'''
symGenList = []
'''A list of flags that if True indicates a constraint was generated by symmetry
'''
holdParmList = []
'''List of parameters that should not be refined ("Hold"s).
Set in :func:`StoreHold`. Initialized in :func:`InitVars`.
'''
holdParmType = {}
'''The reason why a parameter has been marked as "Hold". 
Initialized in :func:`InitVars`; set in :func:`StoreHold`.
'''
constrParms = {'indep-equiv':[], 'indep-constr':[], 'dep-equiv':[], 'dep-constr':[]}
'''A dict with parameters in equivalences, compiled from 
(:data:`dependentParmList`) and (:data:`indParmList`).
Used within :func:`GetIndependentVars` and :func:`GetDependentVars`.
'''
saveVaryList = []
'''A list of the varied parameters that was last supplied when constraints were
processed. This is set in :func:`GenerateConstraints` and updated in 
:func:`Map2Dict`. Used in :func:`VarRemapShow`
'''
droppedSym = []
'''A list of symmetry generated equivalences that have been converted to 
constraint equations in :func:`CheckEquivalences`
'''
#------------------------------------------------------------------------------------
# Global vars set in :func:`GenerateConstraints`. Used for intermediate processing of
# constraints.

constrVarList = []
'List of parameters used in "Constr" and "New Var" constraints'
indepVarList = []
'A list of all independent parameters in equivalences'
depVarList = []
'A list of all dependent parameters in equivalences'

#------------------------------------------------------------------------------------
# global variables used in :func:`getConstrError` to store error and warning information.
# set in CheckEquivalences and in GenerateConstraints
convVarList = []
'parameters in equivalences that were converted to "Const" constraints'
multdepVarList = []
'parameters used as dependents multiple times in equivalences'
unvariedParmsList = []
'parameters used in equivalences that are not varied'
undefinedVars = []
'parameters used in equivalences that are not defined in the parameter dict'
groupErrors = []
'parameters in constraints where parameter grouping and matrix inversion fails'

#------------------------------------------------------------------------------------
paramPrefix = "::constr"
'A prefix for generated parameter names'
consNum = 0
'The number to be assigned to the next constraint to be created'

#------------------------------------------------------------------------------------
class ConstraintException(Exception):
    '''Defines an Exception that is used when an exception is raised processing constraints.
    Raised in :func:`GenerateConstraints` during sequential fits. Possible (but highly unlikely) 
    to be raised in :func:`CheckEquivalences` (called by :func:`GenerateConstraints`) if an 
    infinite loop is detected.
    Also raised in :func:`GramSchmidtOrtho` and :func:`_SwapColumns` but caught 
    within :func:`GenerateConstraints`.
    '''
    pass

def InitVars():
    '''Initializes all constraint information'''
    global dependentParmList,arrayList,invarrayList,indParmList,consNum,symGenList
    dependentParmList = [] # contains a list of parameters in each group
    arrayList = [] # a list of of relationship matrices 
    invarrayList = [] # a list of inverse relationship matrices 
    indParmList = [] # a list of names for the new parameters
    consNum = 0 # number of the next constraint to be created
    symGenList = [] # Flag if constraint is generated by symmetry
    global holdParmList,holdParmType
    holdParmList = []
    holdParmType = {}
    global droppedSym
    droppedSym = []
    global constrVarList,indepVarList,depVarList
    constrVarList = []
    indepVarList = []
    depVarList = []
    global constrParms
    for k in constrParms:
        constrParms[k] = []

def VarKeys(constr):
    """Finds the keys in a constraint that represent parameters
    e.g. eliminates any that start with '_'

    :param dict constr: a single constraint entry of form::

        {'var1': mult1, 'var2': mult2,... '_notVar': val,...}

        (see :func:`GroupConstraints`)
    :returns: a list of keys where any keys beginning with '_' are
      removed.
    """
    return [i for i in constr.keys() if not i.startswith('_')]


def GroupConstraints(constrDict):
    """Divide the constraints into groups that share no parameters.

    :param dict constrDict: a list of dicts defining relationships/constraints

    ::
    
       constrDict = [{<constr1>}, {<constr2>}, ...]

    where {<constr1>} is {'var1': mult1, 'var2': mult2,... }

    :returns: two lists of lists:
    
      * a list of grouped contraints where each constraint grouped containts a list
        of indices for constraint constrDict entries
      * a list containing lists of parameter names contained in each group
      
    """
    assignedlist = [] # relationships that have been used
    groups = [] # contains a list of grouplists
    ParmList = []
    for i,constrI in enumerate(constrDict):
        if i in assignedlist: continue # already in a group, skip
        # starting a new group
        grouplist = [i,]
        assignedlist.append(i)
        groupset = set(VarKeys(constrI))
        changes = True # always loop at least once
        while(changes): # loop until we can't find anything to add to the current group
            changes = False # but don't loop again unless we find something
            for j,constrJ in enumerate(constrDict):
                if j in assignedlist: continue # already in a group, skip
                if len(set(VarKeys(constrJ)) & groupset) > 0: # true if this needs to be added
                    changes = True
                    grouplist.append(j)
                    assignedlist.append(j)
                    groupset = groupset | set(VarKeys(constrJ))
        group = sorted(grouplist)
        varlist = sorted(list(groupset))
        groups.append(group)
        ParmList.append(varlist)
    return groups,ParmList

def GenerateConstraints(varyList,constrDict,fixedList,parmDict=None,
                        seqHistNum=None,raiseException=False):
    '''Takes a list of relationship entries that have been stored by 
    :func:`ProcessConstraints` into lists ``constrDict`` and ``fixedList``

    This routine then calls :func:`CheckEquivalences` for internal 
    consistency. This includes converting equivalenced variables into 
    constraints when a variable is used in both. 

    Once checked, parameters are grouped so that any parameters that are used in 
    more than one constraint are grouped together. This allows checking for incompatible
    logic (for example, when four constraints are specified for three variables). 

    If parmDict is not None, the parameter groups are checked for constraints where 
    some parameters are varied, but not others. If so, the value for that unvaried 
    parameter is subtracted from the constant in the constraint. 

    Once all checks are complete, the constraints are then
    converted to the form used to apply them, saving them as global variables within 
    this module. 

    :param list varyList: a list of parameters names (strings of form
      ``<ph>:<hst>:<nam>``) that will be varied. Note that this is changed 
      here unless set to None. None is used to indicate that all constraints
      should be generated.
    
    :param dict constrDict: a list of dicts defining relationships/constraints
      (as described in :func:`GroupConstraints`)

    :param list fixedList: a list of values specifying a fixed value for each
      dict in constrDict. Values are either strings that can be converted to
      floats, float values or None if the constraint defines a new parameter.
      
    :param dict parmDict: a dict containing all parameters defined in current
      refinement.

    :param int seqHistNum: the hId number of the current histogram in a sequential
      fit. None (default) otherwise. 
      
    :param bool raiseException: When True, generation of an error causes
       an exception to be raised (used in sequential fits)

    :returns: errmsg,warning,groups,parmlist

      **errmsg**
        Is an error message or empty if no errors were found
      **warning**
        Is a warning message about constraints that have been ignored or changed
      **groups**
        Lists parameter groups
      **parmlist**
        Lists parameters in each parameter groups
    '''
    warninfo = {'msg':'', 'shown':{}}
    def warn(msg,cdict=None,val=None,prefix=None):
        if prefix is None: prefix = '\nProblem with '
        if cdict is not None and cdict != warninfo['shown']:
            warninfo['shown'] = cdict
            if warninfo['msg']: warninfo['msg'] += '\n'
            if '_vary' in cdict:
                warninfo['msg'] += prefix + 'new var expression: ' + _FormatConstraint(cdict,cdict.get('_name','New Var'))
            else:
                warninfo['msg'] += prefix + 'constraint equation: ' + _FormatConstraint(cdict,val)
        if warninfo['msg']: warninfo['msg'] += '\n'
        warninfo['msg'] += '  ' + msg
    
    global dependentParmList,arrayList,invarrayList,indParmList,consNum
    # lists of parameters used for error reporting
    global undefinedVars # parameters that are used in equivalences but are not defined
    undefinedVars = []
    global groupErrors # parameters in constraints that cause grouping errors
    groupErrors = []
    global saveVaryList
    saveVaryList = copy.copy(varyList)
    
    errmsg = ''  # save error messages here. If non-blank, constraints cannot be used.
    warning = '' # save informational text messages here.

    # find parameters used in constraint equations & new var assignments (all are dependent)
    global constrVarList 
    constrVarList = []
    for cnum,(cdict,fixVal) in enumerate(zip(constrDict,fixedList)):
        constrVarList += [i for i in cdict if i not in constrVarList and not i.startswith('_')]

    # Process the equivalences; If there are conflicting parameters, move them into constraints
    warning = CheckEquivalences(constrDict,varyList,fixedList,parmDict,seqHistNum=seqHistNum)
    

    # look through "Constr" and "New Var" constraints looking for zero multipliers and
    # Hold, Unvaried & Undefined parameters
    skipList = []
    invalidParms = []
    for cnum,(cdict,fixVal) in enumerate(zip(constrDict,fixedList)):
        #constrVarList += [i for i in cdict if i not in constrVarList and not i.startswith('_')]
        valid = 0          # count of good parameters
        # error reporting
        zeroList = []      # parameters with zero multipliers
        holdList = []      # parameters with "Hold"'s
        noVaryList = []    # parameters not varied
        noWildcardList = [] # wildcard parameters in non-sequential fit
        notDefList = []    # parameters not defined
        # processing to be done
        problem = False    # constraint must be dropped 
        dropList = []      # parameters to remove
        for var in VarKeys(cdict):   # assemble warning info
            if cdict[var] == 0:    # zero multiplier
                if var not in zeroList: zeroList.append(var)
                if var not in dropList: dropList.append(var)
            elif var in holdParmList:   # hold invalid in New Var, drop from constraint eqn
                holdList.append(var)
                if fixVal is None:   # hold in a newvar is not allowed
                    problem = True
                else:
                    if var not in dropList: dropList.append(var)
            elif ':*:' in var :  # wildcard still present should be treated as undefined
                if var not in undefinedVars: undefinedVars.append(var)
                noWildcardList.append(var)
                problem = True
            elif parmDict is not None and var not in parmDict: # not defined, constraint will not be used
                if var not in undefinedVars: undefinedVars.append(var)
                notDefList.append(var)
                if seqHistNum is None:
                    if ':dAx:' in var or ':dAy:' in var or ':dAz:' in var: # coordinates from undefined atoms 
                        if fixVal is None:
                            problem = True  # invalid in New Var
                        else:
                            if var not in dropList: dropList.append(var) # ignore in constraint eqn
                    else:
                        problem = True
            elif varyList is None:
                valid += 1                
            elif var not in varyList and fixVal is not None:  # unvaried, constraint eq. only
                if var not in unvariedParmsList: unvariedParmsList.append(var)
                noVaryList.append(var)
                dropList.append(var)
            else:
                valid += 1
        if seqHistNum is not None and len(notDefList) > 0 and valid == 0:
            # for sequential ref can quietly ignore constraints with all undefined vars
            notDefList = []
        
        if noVaryList and cdict.get('_vary',True): # true for constr eq & varied New Var
            msg = "parameter(s) not varied: "
            for i,v in enumerate(noVaryList):
                if i != 0: msg += ', '
                msg += v
            warn(msg,cdict,fixVal,prefix='\nUnused ')
        for l,m in ((zeroList,"have zero multipliers"), # show warning
                      (holdList,'set as "Hold"'),
                      #(noVaryList,"not varied"),
                      (noWildcardList,"wildcard in non-sequential fit"),
                      (notDefList,"not defined")):
            if l and cdict.get('_vary',True): # true for constr eq & varied New Var
                msg = "parameter(s) " + m + ': '
                for i,v in enumerate(l):
                    if i != 0: msg += ', '
                    msg += v
                warn(msg,cdict,fixVal)
        if valid == 0: # no valid entries
            if seqHistNum is None:
                warn('Ignoring constraint, contains no usable parameters',cdict,prefix='\nUnused ')
            skipList.append(cnum)
        elif problem: # mix of valid & refined and undefined items, cannot use this
            if cdict.get('_vary',True): # true for constr eq & varied New Var
                warn('New Var constraint will be ignored',cdict,prefix='\nUnused ')
                skipList.append(cnum)
            invalidParms += VarKeys(cdict)
        elif len(dropList) > 0: # mix of valid and problematic items, drop problem vars, but keep rest
            if GSASIIpath.GetConfigValue('debug'): 
                msg = ''
                for v in dropList:
                    if msg: msg += ' ,'
                    msg += v
                warn('removing: '+msg,cdict)
            value = fixedList[cnum]
            for var in dropList:   # do cleanup
                # NB expressions in constraint multipliers have already been evaluated
                if ':dAx:' in var or ':dAy:' in var or ':dAz:' in var: 
                    pass # treat delta coords as 0; no change in total is needed
                elif cdict[var] != 0: 
                    value = float(value) - cdict[var]*parmDict[var]
                del cdict[var]
            if float(value) != float(fixedList[cnum]): fixedList[cnum] = float(np.round(value,12))
            if GSASIIpath.GetConfigValue('debug'):
                warn('revised as: '+_FormatConstraint(constrDict[cnum],fixedList[cnum]))
    for i in list(range(len(constrDict)-1,-1,-1)): # remove the dropped constraints
        if i in skipList:
            del constrDict[i]
            del fixedList[i]
            
    for i in invalidParms: StoreHold(i,"Used in invalid constraint")
    if warning: warning += '\n'
    warning += warninfo['msg']
            
    groups,parmlist = GroupConstraints(constrDict)

    # now process each group and create the relations that are needed to form
    # a non-singular square matrix
    # Now check that all parameters are varied (probably do not need to do this
    # any more). For constraint equations, if all are varied, set VaryFree to True
    # and all newly created relationships will be varied. For NewVar constraints,
    # vary if the vary flag was set. 
    for group,depPrmList in zip(groups,parmlist):
        if len(depPrmList) < len(group): # too many relationships -- no can do
            if errmsg: errmsg += '\n'
            errmsg += "Over-constrained input. "
            errmsg += "There are more constraints (" + str(len(group))
            errmsg += ") than parameters (" + str(len(depPrmList)) + ")\nin these constraints:"
            for rel in group:
                errmsg += '\n\t'+ _FormatConstraint(constrDict[rel],fixedList[rel])
            groupErrors += depPrmList
            continue # go on to next group

        try:
            constrArr = _FillArray(group,constrDict,depPrmList)
        except Exception as err:
            if errmsg: errmsg += '\n'
            if 'Initial' in str(err): 
                errmsg += "\nSingular input. "
                errmsg += "There are internal inconsistencies in these constraints:"
            else:
                errmsg += "\nError expanding matrix with these constraints:"
                for rel in group:
                    errmsg += '\n\t' + _FormatConstraint(constrDict[rel],fixedList[rel])
            groupErrors += depPrmList
            continue

        try:
            GramSchmidtOrtho(constrArr,len(group))
        except:
            if errmsg: errmsg += '\n'
            errmsg += "\nUnexpected singularity with constraints group (in Gram-Schmidt)"
            for rel in group:
                errmsg += '\n\t' + _FormatConstraint(constrDict[rel],fixedList[rel])
            groupErrors += depPrmList
            continue
        
        try: 
            invConstrArr = np.linalg.inv(constrArr)
        except:
            if errmsg: errmsg += '\n'
            errmsg += "\nSingular input. "
            errmsg += "The following constraints are not "
            errmsg += "linearly independent\nor do not "
            errmsg += "allow for generation of a non-singular set.\n"
            errmsg += 'This is unexpected. Please report this (toby@anl.gov)'
            for rel in group:
                errmsg += '\n\t' + _FormatConstraint(constrDict[rel],fixedList[rel])
            groupErrors += depPrmList
            continue

        # scan through current group looking for new var assignments
        hasNewVar = False
        for rel in group:
            if fixedList[rel] is None: 
                hasNewVar = True # there a New Var relationship in this group
                break
        else: # only constraint equations, check for unvaried parameters in one
            unvaried = False
            # this should not happen as they should have been removed
            for var in depPrmList:
                if varyList is not None and var not in varyList: 
                    unvaried = True
                    break
            if unvaried: # something is not varied: skip group & remove all parameters from varyList
                for var in depPrmList: StoreHold(var,'Unexpected: mixed use')
                if GSASIIpath.GetConfigValue('debug'): 
                    print('Unexpected: Constraint group ignored (some parameters unvaried)')
                    for rel in group:
                        print ('   '+_FormatConstraint(constrDict[rel],fixedList[rel]))                    
                continue
        maplist = []   # value or param name mapped to each row in expression matrix            
        if not hasNewVar: # constraint equations; all remaining entries varied, vary generated 
            for i in range(len(depPrmList)):
                if i >= len(group): # tag generated degrees of freedom with names and vary them
                    varname = paramPrefix + str(consNum) # assign a unique name
                    consNum += 1
                    maplist.append(varname)
                    if varyList is not None: varyList.append(varname)
                else:
                    rel = group[i]
                    maplist.append(fixedList[rel])
        else:   # ------------------------- groups with new var assignments, vary only NV's w/flags set
            for i,rel in enumerate(group):
                if fixedList[rel] is None:
                    varname = constrDict[rel].get('_name','::?????')
                    maplist.append(varname)
                    if  constrDict[rel].get('_vary',False):
                        if varyList is not None: varyList.append(varname)
                else:
                   maplist.append(fixedList[rel])   # constraint equation
            for i in range(len(depPrmList)):
                if i >= len(group): # name generated degrees of freedom
                    varname = paramPrefix + str(consNum) # assign a unique name
                    consNum += 1
                    maplist.append(varname)
            for var in depPrmList: StoreHold(var,'New Var use')
        # keep this group
        dependentParmList.append(depPrmList)
        arrayList.append(constrArr)
        invarrayList.append(invConstrArr)
        indParmList.append(maplist)
        symGenList.append(False)
            
    if errmsg and raiseException:
        print (' *** ERROR in constraint definitions! ***')
        print (errmsg)
        if warning:
            print (' also note warnings in constraint processing:')
            print (warning)
        raise ConstraintException
    elif errmsg:
        return errmsg,warning,None,None

    # Make list of dependent and independent variables for all constraints
    constrParms['dep-equiv'] = []
    constrParms['dep-constr'] = []
    constrParms['indep-equiv'] = []
    constrParms['indep-constr'] = []
    for varlist,mapvars,multarr in zip(dependentParmList,indParmList,arrayList):  # process all constraints
        for mv in mapvars:
            if type(mv) is float or type(mv) is int: continue
            if multarr is None and mv not in constrParms['indep-equiv']: 
                constrParms['indep-equiv'].append(mv)
            elif mv not in constrParms['indep-constr']: 
                constrParms['indep-constr'].append(mv)
        for mv in varlist:
            if multarr is None and mv not in constrParms['dep-equiv']: 
                constrParms['dep-equiv'].append(mv)
            elif mv not in constrParms['dep-constr']:
                constrParms['dep-constr'].append(mv)
            StoreHold(mv,'dependent param')
    saveVaryList = copy.copy(varyList)  # save varyList so it can be used within module
    if varyList is None: saveVaryList = []

    # if equivMoved:
    #     print(60*'=')
    #     print('Constraints were reclassified to avoid conflicts, as below:')
    #     print(mvMsg)
    #     print('New constraints are:')
    #     print (VarRemapShow(varyList,True))
    #     print(60*'=')
    return errmsg,warning,groups,parmlist # saved for sequential fits
    
def CheckEquivalences(constrDict,varyList,fixedList,parmDict=None,seqHistNum=None):
    '''Process equivalence constraints, looking for conflicts such as 
    where a parameter is used in both an equivalence and a constraint expression
    or where chaining is done (A->B and B->C). 

    Removes equivalences or parameters from equivalences or converts equivalences to
    constraints as described for :ref:`Equivalence Checking and Reorganization <CheckEquivalences>`.

    :param dict constrDict: a list of dicts defining relationships/constraints
    :param list varyList: list of varied parameters (defined during refinements only)
    :param list fixedList: a list of values specifying a fixed value for each
       dict in constrDict. Values are either strings that can be converted to
       floats or ``None`` if the constraint defines a new parameter rather
       than a constant.
    :param dict parmDict: a dict containing defined parameters and their values. Used to find 
       equivalences where a parameter is has been removed from a refinement. 
    :param int seqHistNum: the hId number of the current histogram in a sequential
      fit. None (default) otherwise. 

    :returns: warning messages about changes that need to be made to equivalences 
    '''
    
    warninfo = {'msg':'', 'shown':-1}
    def warnEqv(msg,cnum=None):
        if cnum is not None and cnum != warninfo['shown']:
            warninfo['shown'] = cnum
            if warninfo['msg']: warninfo['msg'] += '\n'
            warninfo['msg'] += '\nProblem with equivalence: ' + _showEquiv(
                dependentParmList[cnum],indParmList[cnum],invarrayList[cnum])
        if warninfo['msg']: warninfo['msg'] += '\n'
        warninfo['msg'] += '  ' + msg

    global depVarList # parameters used in equivalences as dependent parameters
    global indepVarList # parameters used in equivalences as independent parameters
    global constrVarList  # parameters used in other constraints
             
    # lists of parameters used for error reporting
    global undefinedVars # parameters that are used in equivalences but are not defined
    global convVarList # parameters in equivalences that will be converted to constraints
    convVarList = [] # parameters in equivalences to be made into constraints
    global multdepVarList
    multdepVarList = [] # list of dependent parameters used in more than one equivalence

    # local vars
    dropVarList = []  # parameters that can be removed from equivalences 
    removeList = []  # equivalences that are not needed
    convertList = [] # equivalences that should be converted to "Const" constraints
    
    # tabulate parameters used in equivalences by type
    depVarList = []  # list of all dependent parameters in equivalences
    indepVarList = []  # list of all independent parameters in equivalences
    for cnum,(varlist,mapvars,multarr,invmultarr) in enumerate(zip(
            dependentParmList,indParmList,arrayList,invarrayList)):
        #if multarr is not None: continue # equivalence
        indepVarList += [mv for mv in mapvars if mv not in indepVarList]
        depVarList += [v for v in varlist if v not in depVarList]

    # process equivalences: make a list of dependent and independent vars
    #    and check for repeated uses (repetition of a parameter as an
    #    independent var is OK)
    # look for parameters in equivalences that are used more than once as dependent parameters
    seenOnce = []
    for cnum,(varlist,multarr) in enumerate(zip(dependentParmList,arrayList)):
        if multarr is not None: continue # equivalences only
        for v in varlist:
            if v not in seenOnce:
                seenOnce.append(v)
            elif v not in multdepVarList:
                multdepVarList.append(v)                

    # scan through equivalences looking for other "dual uses". Stop when no new ones are found
    changed = True
    count = 0
    while changed:
        changed = False
        count += 1
        if count > 1000:
            raise ConstraintException("Too many loops in CheckEquivalences")
        
        # look for repeated dependent vars
        convVarList = [] # parameters in equivalences to be made into constraints
        for cnum,(varlist,mapvars,multarr,invmultarr) in enumerate(zip(
            dependentParmList,indParmList,arrayList,invarrayList)):
            if multarr is not None: continue # equivalences only
            if cnum in convertList:
                convVarList += [v for v in mapvars+varlist if v not in convVarList and type(v) is not float]
                continue
            
            # identify equivalences that need to be converted to constraints.
            #  Where parameters:
            #     are used both in equivalences & constraints,
            #     are used as dependent multiple times or
            #     where are used as both dependent and independent (chained)
            msg = False
            for v in mapvars:
                if v in constrVarList+convVarList:
                    changed = True
                    msg = True
                    warnEqv("Independent parameter "+str(v)+' used in constraint',cnum)
                    if cnum not in convertList: convertList.append(cnum)
            for v in varlist:
                if v in multdepVarList:
                    changed = True
                    msg = True
                    warnEqv("Dependent parameter "+str(v)+' repeated',cnum)
                    if cnum not in convertList: convertList.append(cnum)
                elif v in indepVarList:
                    changed = True
                    msg = True
                    warnEqv("Dependent parameter "+str(v)+' used elsewhere as independent',cnum)
                    if cnum not in convertList: convertList.append(cnum)
                elif v in constrVarList+convVarList:
                    changed = True
                    msg = True
                    warnEqv("Dependent parameter "+str(v)+' used in constraint',cnum)
                    if cnum not in convertList: convertList.append(cnum)
            if msg:
                warnEqv('Converting to "Constr"',cnum)
            if warninfo['msg'] and seqHistNum is not None:
                # sequential fit -- print the recasts but don't stop with a warning
                print(warninfo['msg'])
                warninfo = {'msg':'', 'shown':-1}
                
    global unvariedParmsList
    unvariedParmsList = []  # parameters in equivalences that are not varied
    # scan equivalences: look for holds
    for cnum,(varlist,mapvars,multarr,invmultarr) in enumerate(zip(
        dependentParmList,indParmList,arrayList,invarrayList)):
        if multarr is not None: continue # not an equivalence
        if cnum in convertList: continue

        # look for holds
        gotHold = False
        holdList = []
        for v in varlist+mapvars:
            if v in holdParmList:
                gotHold = True
            elif type(v) is not float:
                holdList.append(v)
        if gotHold:
            if holdList:
                msg = '  Some parameters set as "Hold"; setting remainder as "Hold": '
                for i,var in enumerate(holdList):
                    if i != 0: msg += ", "
                    msg += var
                    StoreHold(var,'Equiv fixed')
            else:
                msg = '  All parameters set as "Hold" '
            msg += "\n   Will ignore equivalence"
            warnEqv(msg,cnum)
            removeList.append(cnum)
            continue
        
        # look for unvaried parameters
        gotVary = False
        gotNotVary = False
        holdList = []
        for v in varlist+mapvars:
            if varyList is None or v in varyList:
                gotVary = True
                holdList.append(v)
            elif type(v) is not float:
                gotNotVary = True
                if v not in unvariedParmsList: unvariedParmsList.append(v)
        if gotNotVary:  # at least some unvaried parameters
            if gotVary:  # mix of varied and unvaried parameters
                msg = '\nSome parameters not varied; setting remainder as "Hold": '
                for i,var in enumerate(holdList):
                    if i != 0: msg += ", "
                    msg += var
                    StoreHold(var,'Equiv fixed')
            elif seqHistNum is not None: # don't need to warn for sequential fit
                removeList.append(cnum)
                continue
            else:
                msg = 'No parameters varied '
            msg += "\n  Will ignore equivalence"
            warnEqv(msg,cnum)
            removeList.append(cnum)
            continue
            
        # look for undefined or zero multipliers
        holdList = []
        drop = 0
        for v,m in zip(varlist,invmultarr):    
            if parmDict is not None and v not in parmDict:
                if v not in undefinedVars: undefinedVars.append(v)
                if v not in dropVarList: dropVarList.append(v)
                drop += 1
            elif m == 0:
                warnEqv("Parameter "+str(v)+" has a zero multiplier, dropping",cnum)
                if v not in dropVarList: dropVarList.append(v)
                drop += 1
            else:
                holdList.append(v)
        if drop == len(varlist):
            warnEqv("No dependent parameters defined, ignoring equivalence",cnum)
            removeList.append(cnum)
            continue
        for mv in mapvars:
            if type(mv) is float: continue
            if parmDict is not None and mv not in parmDict:
                # independent parameter is undefined, but some dependent parameters are defined
                # hold them
                if mv not in undefinedVars: undefinedVars.append(mv)
                msg = "Parameter(s) "+str(mv)
                for v in varlist:
                    if v in dropVarList:
                        msg += ', ' + v
                msg += "not defined in this refinement\n"
                msg = "Setting holds for: "
                for i,var in enumerate(holdList):
                    if i != 0: msg += ", "
                    msg += var
                warnEqv(msg,cnum)
                drop += 1
        if drop: # independent var and at least one dependent variable is defined
            msg = "Dropping undefined parameter(s) "
            i = 0
            for v in varlist:
                if v in dropVarList:
                    if i != 0: msg += ', '
                    i += 1
                    msg += v
            warnEqv(msg,cnum)
            msg = "Some parameters not defined. Setting holds for: "
            for i,var in enumerate(holdList):
                if i != 0: msg += ", "
                msg += var
                StoreHold(var,'Equiv fixed')
            warnEqv(msg,cnum)
    
    # Convert equivalences where noted
    for cnum,varlist in enumerate(dependentParmList):
        if cnum not in convertList: continue
        indvar = indParmList[cnum][0]
        # msg = '\nChanging equivalence:\n    ' + _showEquiv(
        #     dependentParmList[cnum],indParmList[cnum],invarrayList[cnum])
        for dep,mult in zip(dependentParmList[cnum],invarrayList[cnum]):
            constrDict += [{indvar:-1.,dep:1./mult[0]}]
            fixedList += ['0.0']
        #msg += '\n  to constraint(s):'
        #msg += '\n    ' + _FormatConstraint(constrDict[-1],fixedList[-1])
        removeList.append(cnum)
    # Drop equivalences where noted
    global droppedSym
    droppedSym = []
    if removeList:
        for i in set(removeList):
            if symGenList[i]: 
                droppedSym.append((dependentParmList[i],indParmList[i],arrayList[i],invarrayList[i]))
        for i in sorted(set(removeList),reverse=True):
            del dependentParmList[i],indParmList[i],arrayList[i],invarrayList[i],symGenList[i]
    # Drop variables from remaining equivalences
    for cnum,varlist in enumerate(dependentParmList):
        for j,v in enumerate(varlist):
            drop = []
            if v in dropVarList:
                drop.append(j)
        if drop:
            for j in sorted(drop,reverse=True):
                del indParmList[cnum][j]
    return warninfo['msg']

def ProcessConstraints(constList,seqmode='use-all',seqhst=None):
    """Interpret the constraints in the constList input into a dictionary, etc.
    All :class:`GSASIIobj.G2VarObj` objects are mapped to the appropriate
    phase/hist/atoms based on the object internals (random Ids). If this can't be
    done (if a phase has been deleted, etc.), the variable is ignored.
    If the constraint cannot be used due to too many dropped variables,
    it is counted as ignored. In the case of sequential refinements, 
    the current histogram number is substituted for a histogram number of "*".

    NB: this processing does not include symmetry imposed constraints
    
    :param list constList: a list of lists where each item in the outer list
      specifies a constraint of some form, as described in the :mod:`GSASIIobj`
      :ref:`Constraint definitions <Constraint_definitions_table>`.
    :param str seqmode: one of 'use-all', 'wildcards-only' or 'auto-wildcard'. 
       When seqmode=='wildcards-only' then any constraint with a numerical 
       histogram number is skipped. With seqmode=='auto-wildcard',
       any non-null constraint number is set to the selected histogram.
    :param int seqhst: number for current histogram (used for 
      'wildcards-only' or 'auto-wildcard' only). Should be None for 
      non-sequential fits.

    :returns:  a tuple of (constrDict,fixedList,ignored) where:
      
      * constrDict (list of dicts) contains the constraint relationships
      * fixedList (list) contains the fixed values for each type
        of constraint.
      * ignored (int) counts the number of invalid constraint items
        (should always be zero!)
    """
    constrDict = []
    fixedList = []
    ignored = 0
    namedVarList = []
    for constr in constList:
        terms = copy.deepcopy(constr[:-3]) # don't change the tree contents
        # deal with wildcards in sequential fits
        if seqmode == 'wildcards-only' and seqhst is not None:
            skip = False
            for term in terms:
                if term[1].histogram == '*':
                    term[1] = term[1].varname(seqhst)
                elif term[1].histogram:
                    skip = True
            if skip: continue
        elif seqmode == 'auto-wildcard' and seqhst is not None:
            for term in terms:
                term[1] = term[1].varname(seqhst)
        elif seqhst is not None:
            for term in terms:
                if term[1].histogram == '*':
                    term[1] = term[1].varname(seqhst)
                # else:
                #     term[1] = term[1].varname()   # does this change anything???
        # separate processing by constraint type
        if constr[-1] == 'h':
            # process a hold
            var = str(terms[0][1])
            if '?' not in var:
                StoreHold(var,'User supplied')
            else:
                ignored += 1
        elif constr[-1] == 'f':
            # process a new variable
            fixedList.append(None)
            D = {}
            varyFlag = constr[-2]
            varname = constr[-3]
            for term in terms:
                var = str(term[1])
                if '?' not in var:
                    D[var] = term[0]
            # add extra dict terms for input variable name and vary flag
            if varname is None: # no assigned name, create one
                global consNum
                varname = str(consNum) 
                consNum += 1
            else:
                varname = str(varname) # in case this is a G2VarObj
            if '::' in varname:
                D['_name'] = varname.replace('::','::nv-')
            else:
                D['_name'] = '::nv-' + varname
            D['_name'] = G2obj.MakeUniqueLabel(D['_name'],namedVarList)
            D['_vary'] = varyFlag == True # force to bool
            constrDict.append(D)
        elif constr[-1] == 'c': 
            # process a constraint equation
            D = {}
            for term in terms:
                var = str(term[1])
                if '?' not in var:
                    D[var] = term[0]
            if len(D) >= 1:
                fixedList.append(float(constr[-3]))
                constrDict.append(D)
            else:
                ignored += 1
        elif constr[-1] == 'e':
            # process an equivalence
            firstmult = None
            eqlist = []
            for term in terms:
                if term[0] == 0: term[0] = 1.0
                var = str(term[1])
                if '?' in var: continue
                if firstmult is None:
                    firstmult = term[0]
                    firstvar = var
                else:
                    eqlist.append([var,firstmult/term[0]])
            if len(eqlist) > 0:
                StoreEquivalence(firstvar,eqlist,False)
            else:
                ignored += 1
        else:
            ignored += 1
    return constrDict,fixedList,ignored

def StoreHold(var,holdType=None):
    '''Takes a variable name and prepares it to be removed from the 
    refined variables.

    Called with user-supplied constraints by :func:`ProcessConstraints`.
    At present symGen is not used, but could be set up to track Holds generated
    by symmetry.
    '''
    global holdParmList,holdParmType
    if var not in holdParmList:
        holdParmList.append(var)
        if holdType: holdParmType[var] = holdType
    
def StoreEquivalence(independentVar,dependentList,symGen=True):
    '''Takes a list of dependent parameter(s) and stores their
    relationship to a single independent parameter (independentVar).

    Called with user-supplied constraints by :func:`ProcessConstraints`,
    with Pawley constraints from :func:`GSASIIstrIO.GetPawleyConstr`, 
    with Unit Cell constraints from :func:`GSASIIstrIO.cellVary`
    with symmetry-generated atom constraints from :func:`GSASIIstrIO.GetPhaseData`

      There is no harm in using StoreEquivalence with the same independent variable::

       StoreEquivalence('x',('y',))
       StoreEquivalence('x',('z',))

      but the same outcome can be obtained with a single call::

       StoreEquivalence('x',('y','z'))

      The latter will run more efficiently. 


      Note that mixing independent and dependent variables, such as:: 

        StoreEquivalence('x',('y',))
        StoreEquivalence('y',('z',))

      is a poor choice. The module will attempt to fix this by transforming the equivalence to a 
      "Const" constraint.

    :param str independentVar: name of master parameter that will be used to determine the value
      to set the dependent variables

    :param list dependentList: a list of parameters that will set from
         independentVar. Each item in the list can be a string with the parameter
         name or a tuple containing a name and multiplier:
         ``['::parm1',('::parm2',.5),]``

    '''
    
    global dependentParmList,arrayList,invarrayList,indParmList,symGenList
    mapList = []
    multlist = []
    allfloat = True
    for var in dependentList:
        if isinstance(var, str):
            mult = 1.0
        elif len(var) == 2:
            var,mult = var
        else:
            raise Exception("Cannot parse "+repr(var) + " as var or (var,multiplier)")
        mapList.append(var)
        try:            
            multlist.append(tuple((float(mult),)))
        except:
            allfloat = False
            multlist.append(tuple((mult,)))
    # added relationships to stored values
    arrayList.append(None)
    if allfloat:
        invarrayList.append(np.array(multlist))
    else:
        invarrayList.append(multlist)
    indParmList.append(list((independentVar,)))
    dependentParmList.append(mapList)
    symGenList.append(symGen)
    return

def SubfromParmDict(s,prmDict):
    '''Process a string as a multiplier and convert it to a float value. This
    is done by subsituting any GSAS-II parameter names that appear in the 
    string that have associated values in the parameter dict with the value 
    for that parameter. 

    :param str s: a string to be converted to a value
    :param dict prmDict: a dictionary with keys as GSAS-II parameter names
      and values the corresponding parameter value.
    :returns: the evaluated expression as a float.
    '''
    # TODO: perhaps SubfromParmDict should be called to convert the
    # fixed-val in constraint equations from strings to values.
    for key in prmDict:
        if key in s:
            s = s.replace(key,str(prmDict[key]))
    return eval(s)

def EvaluateMultipliers(constList,*dicts):
    '''Convert multipliers for constraints and equivalences that are specified
    as strings into values. The strings can specify values in the parameter dicts as 
    well as normal Python functions, such as "2*np.cos(0::Ax:2/2.)"
    
    :param list constList: a list of dicts containing constraint expressions
    :param \\*dicts: one or more dicts containing GSAS-II parameters and their values 
       can be specified
    :returns: an empty string if there were no errors, or an error message listing
       the strings that could not be converted.
    '''
    prmDict = {}
    for d in dicts: prmDict.update(d) # combine all passed parameter dicts
    problemList = ""
    # loop through multipliers in contraint expressions
    for const in constList:
        for key in const:
            if key.startswith('_'): continue
            try: # is this already a float, etc? 
                1+const[key]
                continue
            except:
                pass
            try:
                newval = SubfromParmDict(const[key][:],prmDict)
                if GSASIIpath.GetConfigValue('debug'):
                    print('Changing ',const[key],'to',newval)
                const[key] = newval                
            except:
                if problemList: problemList += ", "
                problemList += const[key]
    # loop through multipliers in equivalences
    global arrayList,invarrayList
    for i,(a,valList) in enumerate(zip(arrayList,invarrayList)):
        if a is not None: continue # ignore if not equiv
        try:
            valList.shape
            continue # ignore if already a numpy array
        except:
            pass
        repList = []
        for v in valList:
            try:
                1+v[0]
                repList.append(tuple((v[0],)))
                continue
            except:
                pass
            try:
                newval = SubfromParmDict(v[0][:],prmDict)
                if GSASIIpath.GetConfigValue('debug'):
                    print('Changing ',v[0],'to',newval)
                repList.append(tuple((newval,)))
            except:
                if problemList: problemList += ", "
                problemList += v[0]
                repList.append(tuple(('error',)))
        invarrayList[i] = np.array(repList)
    return problemList

def GetDependentVars(opt=None):
    '''Return a list of dependent variables: e.g. parameters that are
    constrained in terms of other parameters

    :param str opt: type of dependent variables. 
       'equiv': from equivalences,
       'constr': from constraints 
       None (default): all
    :returns: a list of parameter names
    '''
    if opt == 'equiv':
        return constrParms['dep-equiv']
    elif opt == 'constr':
        return constrParms['dep-constr']
    else:
        return constrParms['dep-equiv']+constrParms['dep-constr']

def GetIndependentVars():
    '''Return a list of independent variables: e.g. parameters that are
    slaved to other parameters by constraints

    :returns: a list of parameter names

    '''
    return constrParms['indep-equiv']+constrParms['indep-constr']

def getConstrError(constrLst,seqmode,seqhst):
    '''This is used to display error messages for constraints and 
    equivalence relations

    :parm list constrLst: a single constraint or equivalence as saved 
      in the data tree
      (see :ref:`constraint definitions <Constraint_definitions_table>`).
    :param str seqmode: one of 'use-all', 'wildcards-only' or 'auto-wildcard'
    :param int seqhst: number for current histogram (used for 
      'wildcards-only' or 'auto-wildcard' only). Should be None for 
      non-sequential fits.

    :returns: error, msg where error (bool) is True if the 
      constraint/equivalence creates an error, msg (str) can be a warning 
      or an error
    '''
    msg = ''
    note = ''
    terms = copy.deepcopy(constrLst[:-3])
    if seqmode == 'wildcards-only' and seqhst is not None:
        if constrLst[-1] == 'e':
            msg = 'equivalence'
        else:
            msg = 'constraint'
        for term in terms:
            if term[1].histogram == '*':
                term[1] = term[1].varname(seqhst)
            elif term[1].histogram:
                return False,"Ignoring non-wildcard "+msg, "Ignore"
    elif seqmode == 'auto-wildcard' and seqhst is not None:
        for term in terms:
            term[1] = term[1].varname(seqhst)
    else:
        for term in terms:
            if term[1].histogram == '*':
                if seqhst is None:
                    msg = "Parameter "+str(terms[0][1])+" contains a wildcard, which are used only sequential refinements. Constraint ignored."
                    return False,msg, "Ignored"
                else:
                    term[1] = term[1].varname(seqhst)
            else:
                term[1] = term[1].varname()
    if constrLst[-1] == 'e':
        # conflicting uses
        if terms[0][1] in constrVarList+convVarList:
            msg = "Parameter "+str(terms[0][1])+" used in constraint. To be recast as constraint."
            return False,msg, "Recast as constraint."
        varList = []
        for m,v in terms:
            if v in constrVarList+convVarList:
                varList.append(str(v))
        if varList:
            msg = "Parameter(s) used in constraint: "
            for i,v in enumerate(varList):
                if i != 0: msg += ', '
                msg += v
            return False,msg,"Recast as constraint."
        varList = []
        for m,v in terms[1:]:
            if v in indepVarList:
                varList.append(str(v))
        if varList:
            msg = "Parameter(s) used elsewhere as independent: "
            for i,v in enumerate(varList):
                if i != 0: msg += ', '
                msg += v
            return False,msg,"Recast as constraint."
        varList = []
        for m,v in terms[1:]:
            if v in multdepVarList:
                varList.append(str(v))
        if varList:
            msg += "Parameter(s) repeated as dependent: "
            for i,v in enumerate(varList):
                if i != 0: msg += ', '
                msg += v
            return False,msg,"Recast as constraint."

        # zero multiplier
        varList = []
        valid = 0
        for m,v in terms:
            if m == 0:
                varList.append(str(v))
            else:
                valid += 1
        if varList and valid > 1:
            msg += "Parameter(s) with zero multipliers: "
            for i,v in enumerate(varList):
                if i != 0: msg += ', '
                msg += v
            msg += " will be ignored"
        elif varList:
            msg += "Parameter(s) with zero multipliers:"
            for i,v in enumerate(varList):
                if i != 0: msg += ', '
                msg += v
            return False,msg,"Ignored: not varied"
        
        # hold parameters
        s = ''
        for m,v in terms:
            if v in holdParmList and (
                    "User supplied" in holdParmType.get(v,'') or
                    "symmetry" in holdParmType.get(v,'') or
                    "rigid body" in holdParmType.get(v,'')):
                if s: s += ', '
                s += str(v)
        if s:
            if msg: msg += '; '
            msg += "Parameter(s) set as Hold: "+s
            msg += "\nThis equivalence will be ignored; all parameters are held."
            return False,msg,'Has holds: Not varied.'
                
        # unrefined parameters
        gotVary = False
        gotNotVary = False
        s = ''
        for m,v in terms:
            if v in unvariedParmsList:
                gotNotVary = True
                if s: s += ', '
                s += str(v)
            else:
                gotVary = True
        if gotNotVary and gotVary:  # mix of varied and unvaried parameters
            if msg: msg += '. '
            msg += 'Unvaried parameter(s): '+s+"; remainder set Hold. All parameters fixed."
            return False,msg, 'Ignored: not all varied'
        elif gotNotVary:
            if msg: msg += '. '
            msg += 'All parameters not varied. Equivalence Ignored.'
            return False,msg, 'Ignored: not varied'
            
        # undefined parameters
        undef = 0
        s = ''
        for m,v in terms[1:]:
            if v in undefinedVars:
                undef += 1
                if s: s += ', '
                s += str(v)
        if undef == len(terms[1:]):
            msg += 'Ignored: None of the dependent parameters are defined'
        elif terms[0][1] in undefinedVars:
            if s:
                s = terms[0][1] + ', ' + s
            else:
                s = terms[0][1]
            msg += 'Undefined parameter(s): '+s+'. Remainder will be fixed'
        elif undef:
            msg += 'Undefined parameter(s): '+s+' will be dropped'
    elif constrLst[-1] == 'h':
        v = terms[0][1]
        if v in undefinedVars: return False,"Parameter is undefined","Ignored"
        if v in unvariedParmsList: return False,"Parameter is not refined","Ignored"
    else:
        # check for post-grouping errors in new var & constr. eq. groups
        for m,v in terms:
            if v in groupErrors:
                return True,'Constraint singularity: see error listing','Singular'
        zeroList = []
        toBeUsed = []
        undef = []
        unvar = []
        hold = []
        for m,v in terms: # check for zero multiplier, undefined, unvaried or hold
            if constrLst[-1] == 'f':    # new var expressions; 
                if m == 0:
                    zeroList.append(str(v))
                elif v in undefinedVars:
                    undef.append(str(v))
                elif (v in holdParmList and constrLst[-2] and
                        "User supplied" in holdParmType.get(v,'') or
                        "symmetry" in holdParmType.get(v,'') or
                        "rigid body" in holdParmType.get(v,'')):                        
                    hold.append(str(v))
            else:                        # constraint equation
                if m == 0:
                    zeroList.append(str(v))
                elif v in undefinedVars:
                    undef.append(str(v))
                elif v in unvariedParmsList:
                    unvar.append(str(v))
                elif v in holdParmList and holdParmType.get(v,'') != 'dependent param':
                    hold.append(str(v))
                else:
                    toBeUsed.append(str(v))

        s = ''
        for v in zeroList:
            if s: s += ', '
            s += str(v)
        if s:
            if msg: msg += '; '
            msg += "Parameter(s) with zero multipliers: "+s
        s = ''
        for v in undef:
            if s: s += ', '
            s += str(v)
        if s:
            if msg: msg += '; '
            msg += "Undefined parameter(s): "+s
        s = ''
        for v in unvar:
            if s: s += ', '
            s += str(v)
        if s:
            if msg: msg += '; '
            msg += "Unrefined parameter(s) will be dropped: "+s
        s = ''
        for v in hold:
            if s: s += ', '
            s += str(v)
        if s:
            if msg: msg += '; '
            msg += '"Hold" parameter(s): '+s
        if hold and constrLst[-1] == 'f':
            if msg: msg += '; '
            msg += "\nNew var with holds cannot be processed; constraint ignored."
            note = 'Ignored'
        elif undef and toBeUsed:
            s = ''
            for v in toBeUsed:
                if s: s += ', '
                s += str(v)
            if msg: msg += '; '
            msg += "Adding Holds on "+s+"; Constraint Ignored."
            note = 'Ignored'
        elif undef or (len(toBeUsed) == 0 and (zeroList or unvar or hold)):
            if msg: msg += '; '
            msg += "No parameters remain. Constraint Ignored."
            note = 'Ignored'
        elif len(toBeUsed) == 1 and (zeroList or unvar or hold):
            if msg: msg += '; '
            msg += "One parameter is retained; converted to fixed value."
            note = 'Converted'
        elif zeroList or unvar or hold:
            if msg: msg += ': '
            msg += '\nConstraint adjusted for fixed values and retained.'
            note = 'Parameter(s) removed'
    return False,msg,note
        
def ComputeDepESD(covMatrix,varyList,noSym=False):
    '''Compute uncertainties for dependent parameters from independent ones
    returns a dictionary containing the esd values for dependent parameters
    
    :param np.array covMatrix: the full covariance matrix
    :param list varyList: the names of the variables matching the columns 
      and rows in covMatrix
    :param bool noSym: When True symmetry generated parameters are 
      not included. Do this so that redundant s.u.'s eare not shown.
      When False (default) s.u. values for all dependent 
      parameters are placed in the returned dict. 
    '''
    sigmaDict = {}
    for varlist,mapvars,multarr,invmultarr,symgen in zip(
            dependentParmList,indParmList,arrayList,invarrayList,symGenList):
        if symgen and noSym: continue # skip symmetry generted 
        varied = 0
        # get the v-covar matrix for independent parameters 
        vcov = np.zeros((len(mapvars),len(mapvars)))
        for i1,name1 in enumerate(mapvars):
            if name1 not in varyList: continue
            varied += 1
            iv1 = varyList.index(name1)
            for i2,name2 in enumerate(mapvars):
                if name2 not in varyList: continue
                iv2 = varyList.index(name2)
                vcov[i1][i2] = covMatrix[iv1][iv2]
        # vec is the vector that multiplies each of the independent values
        for i,(v,vec) in enumerate(zip(varlist,invmultarr)):
            #if i == varied: break # this limits the number of generated params
            # to match the number varied. Not sure why I did this. 
            sigmaDict[v] = np.sqrt(np.inner(vec.T,np.inner(vcov,vec)))
    return sigmaDict

def _FormatConstraint(RelDict,RelVal):
    '''Formats a Constraint or Function for use in a convenient way'''
    linelen = 65
    s = [""]
    for var,val in RelDict.items():
        if var.startswith('_'): continue
        if len(s[-1]) > linelen: s.append(' ')
        m = val
        if s[-1] != "" and m >= 0:
            s[-1] += ' + '
        elif s[-1] != "":
            s[-1] += ' - '
            m = abs(m)
        if m == 1:
            s[-1] += '%s '%var
        else:
            s[-1] += '%.3f*%s '%(m,var)
    if len(s[-1]) > linelen: s.append(' ')
    if RelVal is None:
        s[-1] += ' = New variable'
    else:
        s[-1] += ' = ' + str(RelVal)
    s1 = ''
    for s2 in s:
        if s1 != '': s1 += '\n  '
        s1 += s2
    return s1

def _showEquiv(varlist,mapvars,invmultarr,longmsg=False):
    '''Format an equivalence relationship, note that varlist, mapvars, invmultarr 
    are elements of dependentParmList, indParmList, invarrayList
    '''
    for i,mv in enumerate(mapvars):
        s1 = str(mv)
        if not longmsg:
            s1 += ' ==> '
        elif len(varlist) == 1:
            s1 += ' is equivalent to '
        else:
            s1 += ' is equivalent to parameters: '
        j = 0
        for v,m in zip(varlist,invmultarr):
            if debug: print ('v,m[0]: ',v,m[0])
            if len(s1.split('\n')[-1]) > 60: s1 += '\n        '
            if j > 0: s1 += ' & '
            j += 1
            s1 += str(v)
            if m != 1:
                s1 += " / " + str(m[0])
    return s1

def VarRemapSumm():
    '''Summarize the constraints in a single line

    :returns: a string summarizing the contraint relationships. May be blank
    '''
    varyList = saveVaryList
    s = ''
    depVars = constrParms['dep-equiv']+constrParms['dep-constr']
        
    userOut = 0
    symOut = 0
    consOut = 0
    varOut = 0
    freeOut = 0
    variedOut = 0
    global dependentParmList,arrayList,invarrayList,indParmList,symGenList
    for varlist,mapvars,multarr,invmultarr,symFlag in zip(
        dependentParmList,indParmList,arrayList,invarrayList,symGenList):
        for i,mv in enumerate(mapvars):
            if multarr is None:
                if symFlag:
                    symOut += 1
                else:
                    userOut += 1
                continue            
            if mv in varyList: 
                variedOut += 1
            if type(mv) is float:
                consOut += 1
            elif '::nv-' in mv:
                varOut += 1
            else:
                freeOut += 1
    numHold = len([i for i in holdParmList if i not in depVars])
    if userOut:
        if s: s += '; '
        s += f'Equivalences: {userOut}'
    if numHold > 0:
        if s: s += '; '
        s += f'Hold param: {numHold}'
#    if symOut:  # don't need to log this
#        if s: s += '; '
#        s += f'Sym gen param: {symOut}'
    if consOut + freeOut:
        if s: s += '; '
        s += f'Constr Eqns: {consOut + freeOut}'
    if freeOut:
        s += f' ({freeOut} generated)'
    if varOut:
        if s: s += '; '
        s += f'New Var Eqns: {varOut}'
    if variedOut:
        s += f' ({freeOut} varied)'
    if len(depVars) > 0:
        if s: s += '. '
        s += f'Depend. param: {len(depVars)}'
    return s
    
def VarRemapShow(varyList=None,inputOnly=False,linelen=60):
    '''List out the saved relationships. This should be done after the constraints have been
    defined using :func:`StoreEquivalence`, :func:`GroupConstraints` and :func:`GenerateConstraints`.

    :returns: a string containing the details of the contraint relationships
    '''
    if varyList is None:
        varyList = saveVaryList
    s = ''
    depVars = constrParms['dep-equiv']+constrParms['dep-constr']
    if len(depVars) > 0:
        s += '\nDependent parameters from constraints & equivalences:\n'
        for v in sorted(set(depVars)):
            s += '    ' + v + '\n'
    first = True
    for v in sorted(holdParmList):
        if v in depVars: continue
        if v in holdParmType: v += '\t'+holdParmType[v]
        if first:
            s += '\nParameters set as Hold:\n'
        first = False
        s += '    ' + v + '\n'
        
    userOut = ''
    symOut = ''
    consOut = ''
    varOut = ''
    freeOut = ''
    global dependentParmList,arrayList,invarrayList,indParmList,symGenList

    for varlist,mapvars,multarr,invmultarr,symFlag in zip(
        dependentParmList,indParmList,arrayList,invarrayList,symGenList):
        for i,mv in enumerate(mapvars):
            if multarr is None:
#                s1 = '  ' + str(mv) + ' is equivalent to parameter(s): '
                if len(varlist) == 1:
                    s1 = '   ' + str(mv) + ' is equivalent to '
                else:
                    s1 = '   ' + str(mv) + ' is equivalent to parameters: '
                j = 0
                for v,m in zip(varlist,invmultarr):
                    if debug: print ('v,m[0]: ',v,m[0])
                    if j > 0: s1 += ' & '
                    j += 1
                    s1 += str(v)
                    if m != 1:
                        s1 += " / " + '{:.4f}'.format(m[0])
                    #if len(s1.split('\n')[-1]) > 70: 
                    #    s1 = ' \n          &'.join(s1.rsplit('&',1))
                if symFlag:
                    symOut += s1 + '\n'
                else:
                    userOut += s1 + '\n'
                continue
            ln = ''
            # if mv in varyList: 
            #     lineOut = '  (V)* {} = '.format(mv)
            # else:
            #     lineOut = '  {} = '.format(mv)
            j = 0 
            lineOut = '  {} = '.format(mv)
            for (m,v) in zip(multarr[i,:],varlist):
                if m == 0: continue
                if m < 0:
                    lineOut += ' - '
                    m *= -1
                elif j != 0:
                    lineOut += ' + '
                j += 1
                if len(lineOut) > linelen:
                    ln += lineOut
                    lineOut = '\n  '
                if m == 1:
                    lineOut += '{}'.format(v)
                else:
                    lineOut += '({:.4g} * {})'.format(m,v)
            if mv in varyList: 
                lineOut += '\t *VARIED*'
            ln += lineOut
            
            if type(mv) is float:
                consOut += ln + '\n'
            elif '::nv-' in mv:
                varOut += ln + '\n'
            else:
                freeOut += ln + '\n'
    if userOut:
        s += '\nUser-supplied equivalences:\n' + userOut
    if droppedSym:
        for varlist,mapvars,multarr,invmultarr in droppedSym:
            for i,mv in enumerate(mapvars):
                if len(varlist) == 1:
                    s1 = '   ' + str(mv) + ' is equivalent to '
                else:
                    s1 = '   ' + str(mv) + ' is equivalent to parameters: '
                j = 0
                for v,m in zip(varlist,invmultarr):
                    if debug: print ('v,m[0]: ',v,m[0])
                    if j > 0: s1 += ' & '
                    j += 1
                    s1 += str(v)
                    if m != 1:
                        s1 += " / " + '{:.4f}'.format(m[0])
                symOut += s1 + '\t *recast as equation*\n'        
    if symOut:
        s += '\nSymmetry-generated equivalences:\n' + symOut
    if consOut:
        s += '\nConstraint Equations:\n' + consOut
    if varOut:
        s += '\nNew Variable assignments:\n' + varOut
    if freeOut:
        s += '\nGenerated constraint equations:\n' + freeOut
    if not (userOut or consOut or varOut or symOut):
        return s + '\nNo constraints or equivalences in use'
    elif inputOnly:
        return s
        
    s += '\nInverse parameter mapping relations:\n'
    lineDict = {} # store so we can sort them
    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
        varied = False
        for i,mv in enumerate(varlist):
            constrVal = 0    # offset to constraint from constant terms
            lineOut = '  {} = '.format(mv)
            j = 0 # number of terms shown
            for m,v in zip(invmultarr[i,:],mapvars):
                if np.isclose(m,0): continue
                try:
                    if np.isclose(float(v),0): continue
                except:
                    pass
                if v in varyList: varied = True
                try:
                    constrVal += m*v  # v is a constant; add to offset
                    continue
                except:
                    pass
                if m < 0:
                    lineOut += ' - '
                    m *= -1
                elif j != 0:
                    lineOut += ' + '
                j += 1
                if len(lineOut) > linelen:
                    s += lineOut
                    lineOut = '\n  '
                if m == 1:
                    lineOut += '{}'.format(v)
                else:
                    lineOut += '({:.4g} * {})'.format(m,v)
            if constrVal < 0:
                lineOut += ' - {:.4g}'.format(-constrVal)
            elif constrVal > 0:
                lineOut += ' + {:.4g}'.format(constrVal)
            elif j == 0:
                lineOut += '0'  # no terms, no constants: var fixed at zero
            if varied: lineOut += '\t *VARIED*'
            lineDict[mv] = lineOut
    for key in sorted(lineDict):
        s += lineDict[key] + '\n'
    return s

def CountUserConstraints():
    '''Count the number of user-supplied constraints in use for the 
    current refinement. Symmetry generated constraints are not counted
    as well as constraints that involve variables that are not varied.

    This is used for CIF reporting.

    :returns: an integer count
    '''
    varyList = saveVaryList
    global dependentParmList,arrayList,invarrayList,indParmList,symGenList

    count = 0
    for varlist,mapvars,multarr,invmultarr,symFlag in zip(
        dependentParmList,indParmList,arrayList,invarrayList,symGenList):
        if symFlag: continue
        if multarr is None:
            for i,mv in enumerate(mapvars):
                if mv not in varyList: continue
                count += len(varlist)
        else:
            # constraints = original # vars - number still generated
            # = number of fixed items in constraint
            count += sum([1 for mv in mapvars if mv not in varyList])
    return count

def getInvConstraintEq(var,varyList):
    '''For a dependent variable, find the constraint that 
    defines the dependent variable in terms of varied independent variables.
    This works for constraint equations (via new var or generated parameters)
    or equivalences. For equivalences the result will lists of length 1
    
    :param str var: named of refined variable (e.g. 0:0:Scale)
    :param list varyList: list of refined variables
    :returns: vList,mList where vList is a list of variables and 
      mList is a list of multipliers for that variable (floats)
    '''
    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
        if var not in varlist: continue
        i = varlist.index(var)
        vList = []
        mList = []
        for m,v in zip(invmultarr[i,:],mapvars):
            if v not in varyList: continue
            if m == 0: continue
            if v == 0: continue
            vList.append(v)
            mList.append(m)
        return vList,mList
    if GSASIIpath.GetConfigValue('debug'): print('getInvConstraintEq: not found: ',var)
    return [],[]   # unexpected -- not an independent parameter
    
def GetSymEquiv(seqmode,seqhistnum):
    '''Return the automatically generated (equivalence) relationships.

    :returns: a list of strings containing the details of the contraint relationships
    '''
    symout = []
    symerr = []
    symhelp = []
    global dependentParmList,arrayList,invarrayList,indParmList,symGenList

    for varlist,mapvars,multarr,invmultarr,symFlag in zip(
        dependentParmList,indParmList,arrayList,invarrayList,symGenList):
        if not symFlag: continue
        for i,mv in enumerate(mapvars):
            cnstr = [[1,G2obj.G2VarObj(mv)]]
            if multarr is None:
                s1 = ''
                s2 = ' = ' + str(mv)
                j = 0
                helptext = 'Variable {:} '.format(mv) + " ("+ G2obj.fmtVarDescr(mv) + ")"
                if len(varlist) == 1:
                    cnstr.append([invmultarr[0][0],G2obj.G2VarObj(varlist[0])])
                    # format the way Bob prefers
                    if invmultarr[0][0] == 1: 
                        s1 = str(varlist[0]) + ' = ' + str(mv)
                    else:
                        s1 = str(varlist[0]) + ' = ' + str(
                            invmultarr[0][0]) + ' * '+ str(mv)
                    s2 = ''
                    
                    m = 1./invmultarr[0][0]
                    var1 = str(varlist[0])
                    helptext += "\n\nis equivalent to "
                    if m == 1:
                        helptext += '\n  {:} '.format(var1) + " ("+ G2obj.fmtVarDescr(var1) + ")"
                    else:
                        helptext += '\n  {:3g} * {:} '.format(m,var1) + " ("+ G2obj.fmtVarDescr(var1) + ")"
                else:
                    helptext += "\n\nis equivalent to the following:"
                    for v,m in zip(varlist,invmultarr):
                        cnstr.append([m,G2obj.G2VarObj(v)])
                        #if debug: print ('v,m[0]: ',v,m[0])
                        if len(s1.split('\n')[-1]) > 75: s1 += '\n        '
                        if j > 0: s1 += ' =  '
                        j += 1
                        s1 += str(v)
                        if m != 1:
                            s1 += " / " + str(m[0])
                            helptext += '\n  {:3g} * {:} '.format(m,v) + " ("+ G2obj.fmtVarDescr(v) + ")"
                        else:
                            helptext += '\n  {:} '.format(v) + " ("+ G2obj.fmtVarDescr(v) + ")"
                err,msg,note = getConstrError(cnstr+[None,None,'e'],seqmode,seqhistnum)
                symerr.append([msg,note])
                symout.append(s1+s2)
                symhelp.append(helptext)
            else:
                s = '  %s = ' % mv
                j = 0
                for m,v in zip(multarr[i,:],varlist):
                    if m == 0: continue
                    if j > 0: s += ' + '
                    j += 1
                    s += '(%s * %s)' % (m,v)
                print ('unexpected sym op='+s)
    return symout,symerr,symhelp

def GetDroppedSym(seqmode,seqhistnum):
    '''Return automatically generated (equivalence) relationships that were
    converted to constraint equations

    :returns: a list of strings containing the details of the equivalences
    '''
    symout = []
    symerr = []
    symhelp = []
    global droppedSym
    for varlist,mapvars,multarr,invmultarr in droppedSym:
        for i,mv in enumerate(mapvars):
            cnstr = [[1,G2obj.G2VarObj(mv)]]
            if multarr is None:
                s1 = ''
                s2 = ' = ' + str(mv)
                j = 0
                helptext = 'Variable {:} '.format(mv) + " ("+ G2obj.fmtVarDescr(mv) + ")"
                if len(varlist) == 1:
                    cnstr.append([invmultarr[0][0],G2obj.G2VarObj(varlist[0])])
                    # format the way Bob prefers
                    if invmultarr[0][0] == 1: 
                        s1 = str(varlist[0]) + ' = ' + str(mv)
                    else:
                        s1 = str(varlist[0]) + ' = ' + str(
                            invmultarr[0][0]) + ' * '+ str(mv)
                    s2 = ''
                    
                    m = 1./invmultarr[0][0]
                    var1 = str(varlist[0])
                    helptext += "\n\nis equivalent to "
                    if m == 1:
                        helptext += '\n  {:} '.format(var1) + " ("+ G2obj.fmtVarDescr(var1) + ")"
                    else:
                        helptext += '\n  {:3g} * {:} '.format(m,var1) + " ("+ G2obj.fmtVarDescr(var1) + ")"
                else:
                    helptext += "\n\nis equivalent to the following:"
                    for v,m in zip(varlist,invmultarr):
                        cnstr.append([m,G2obj.G2VarObj(v)])
                        #if debug: print ('v,m[0]: ',v,m[0])
                        if len(s1.split('\n')[-1]) > 75: s1 += '\n        '
                        if j > 0: s1 += ' =  '
                        j += 1
                        s1 += str(v)
                        if m != 1:
                            s1 += " / " + str(m[0])
                            helptext += '\n  {:3g} * {:} '.format(m,v) + " ("+ G2obj.fmtVarDescr(v) + ")"
                        else:
                            helptext += '\n  {:} '.format(v) + " ("+ G2obj.fmtVarDescr(v) + ")"
                err,msg,note = getConstrError(cnstr+[None,None,'e'],seqmode,seqhistnum)
                symerr.append([msg,note])
                symout.append(s1+s2)
                symhelp.append(helptext)
            else:
                s = '  %s = ' % mv
                j = 0
                for m,v in zip(multarr[i,:],varlist):
                    if m == 0: continue
                    if j > 0: s += ' + '
                    j += 1
                    s += '(%s * %s)' % (m,v)
                print ('unexpected sym op='+s)
    return symout,symerr,symhelp

def Dict2Deriv(varyList,derivDict,dMdv):
    '''Compute derivatives for Independent Parameters from the
    derivatives for the original parameters

    :param list varyList: a list of parameters names that will be varied

    :param dict derivDict: a dict containing derivatives for parameter values keyed by the
      parameter names.

    :param list dMdv: a Jacobian, as a list of np.array containing derivatives for dependent
      parameter computed from derivDict

    '''
    global dependentParmList,arrayList,invarrayList,indParmList,invarrayList
    for varlist,mapvars,multarr,invmultarr in zip(dependentParmList,indParmList,arrayList,invarrayList):
        for i,name in enumerate(mapvars):
            # grouped parameters: need to add in the derv. w/r
            # dependent variables to the independent ones
            if name not in varyList: continue # skip if independent var not varied
            if multarr is None:
                if debug: print ('start dMdv for',name,dMdv[varyList.index(name)])
                for v,m in zip(varlist,invmultarr):
                    if m[0] == 0: continue
                    dMdv[varyList.index(name)] += derivDict[v]/ m[0] 
            else:
                for v,m in zip(varlist,invmultarr[:,i]):
                    if m == 0: continue
                    dMdv[varyList.index(name)] += m * derivDict[v]

def Map2Dict(parmDict,varyList):
    '''Updates the parameter dictionary and the varyList using the 
    equivalence and constraint input. This should be called at least once, after 
    the constraints have been defined using :func:`StoreEquivalence`, 
    :func:`GroupConstraints` and :func:`GenerateConstraints` and before any 
    parameter refinement is done. 

    This completes the parameter dictionary by defining values for parameters 
    created by constraints based on the constraints that define them 
    using the values for the current parameters. It also removes all dependent 
    variables from the varyList

    :param dict parmDict: a dict containing parameter values keyed by the
      parameter names. For new variables created by constraints, entries
      will be added to the dictionary, if not alreay present, or the 
      values will be recomputed.

    :param list varyList: a list of parameters names. Will be modified.
    '''
    # remove fixed parameters from the varyList
    for item in holdParmList:
        if varyList is not None and item in varyList: varyList.remove(item)
            
    # process the independent parameters:
    # * remove dependent ones from varylist
    # * for equivalences apply the independent parameters onto dependent variables
    global dependentParmList,arrayList,invarrayList,indParmList
    for varlist,mapvars,multarr,invmultarr in zip(dependentParmList,indParmList,arrayList,invarrayList):
        for item in varlist: # TODO: is this still needed?
            if varyList is not None and item in varyList: varyList.remove(item)
        if multarr is None:
            #for v,val in zip(  # shows values to be set
            #    varlist,
            #    np.dot(invmultarr,np.array([parmDict[var] for var in mapvars]))
            #    ): print('parmDict set',v,':',val)
            parmDict.update(zip(varlist,np.dot(invmultarr,np.array([parmDict[var] for var in mapvars]))))

    # * for the created parameters, compute them from their dependents
    for varlist,mapvars,multarr in zip(dependentParmList,indParmList,arrayList):
        if multarr is None: continue
        # evaluate constraints in the forward direction
        A = np.array([parmDict[var] for var in varlist])
        z = zip(mapvars,np.dot(multarr,A))
        # add/replace in parameter dict
        parmDict.update([i for i in z if type(i[0]) is not float and ':' in i[0]])
#        parmDict.update([i for i in zip(mapvars,np.dot(multarr,A)) if ':' in i[0]])
    global saveVaryList
    if varyList is not None:
        saveVaryList = copy.copy(varyList)
    else:
        saveVaryList = []
    
def normParms(parmDict):
    '''Attempt to put parameters into the right ballpark by scaling to 
    enforce constraint equations
    '''
    for varlist,mapvars,multarr,invmultarr in zip(
            dependentParmList,indParmList,arrayList,invarrayList):
        if multarr is None or invmultarr is None: continue # unexpected
        for i,s in enumerate(mapvars):
            try:
                s = float(s)
            except: 
                continue
            sumcons = 0.
            for var,m in zip(varlist,multarr[i]):
                if var not in parmDict: 
                    print('normParms error: Parameter',var,'not in parmDict')
                    break
                sumcons += parmDict[var] * m
            if sumcons != s and sumcons != 0 and s != 0:
                for var in varlist:
                    parmDict[var] *= s/sumcons  
            
def Dict2Map(parmDict):
    '''Applies the constraints defined using :func:`StoreEquivalence`,
    :func:`GroupConstraints` and :func:`GenerateConstraints` by changing
    values in a dict containing the parameters. This should be
    done after refinement and before the parameters are used for 
    any computations

    :param dict parmDict: a dict containing parameter values keyed by the
      parameter names. After this is called, all the dependent variables
      will be updated based on constraints and equivalences.
    '''
    global dependentParmList,arrayList,invarrayList,indParmList
    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
        if invmultarr is None:  # is this needed?
            if GSASIIpath.GetConfigValue('debug'): 
                print('Why does this constraint have None for invmultarr?',varlist,mapvars)
            continue
        valslist = np.array([float(parmDict.get(var,var)) for var in mapvars])
        #for v,val in zip(varlist,np.dot(invmultarr,np.array(valslist))): print(v,val) # shows what is being set
        parmDict.update(zip(varlist,np.dot(invmultarr,valslist)))
        
#======================================================================
# internal routines follow (these routines are unlikely to be called
# from outside the module)

def GramSchmidtOrtho(a,nkeep=0):
    '''Use the Gram-Schmidt process (http://en.wikipedia.org/wiki/Gram-Schmidt) to
    find orthonormal unit vectors relative to first row.

    If nkeep is non-zero, the first nkeep rows in the array are not changed
    
    input: 
       arrayin: a 2-D non-singular square array
    returns:
       a orthonormal set of unit vectors as a square array
    '''
    def proj(a,b):
        'Projection operator'
        return a*(np.dot(a,b)/np.dot(a,a))
    for j in range(nkeep,len(a)):
        for i in range(j):
            a[j] -= proj(a[i],a[j])
        if np.allclose(np.linalg.norm(a[j]),0.0):
            raise ConstraintException("Singular input to GramSchmidtOrtho")
        a[j] /= np.linalg.norm(a[j])
    return a

def _FillArray(sel,d,collist):
    '''Construct a n by n matrix [n = len(collist)]
    with the initial m rows [m = len(sel)] using the 
    relationships defined in the expressions dict, d.
    Since m may be smaller than n, the remaining rows
    are filled with rows that are tested to not create
    a singular matrix. 

    :param list sel:     a list of indices in dict d
    :param list d:       a list of dict's where each dict describes an 
      expression from a constraint equation or a new var
    :param list collist: a list parameter names.
    :returns: an n by n numpy.array matrix
    '''
    n = len(collist)
    m = len(sel)
    arr = np.zeros(2*[n,])
    # fill the top rows
    for i,cnum in enumerate(sel):
        for j,var in enumerate(collist):
            arr[i,j] = d[cnum].get(var,0)
    try:
        _RowEchelon(m,copy.copy(arr),collist[:])
    except:
        raise Exception('Initial constraints singular')
    for i in range(m,n):
        arr[i][i] = 1   # add a diagonal element
        try:
            _RowEchelon(i+1,copy.copy(arr),collist[:])
            continue
        except:
            pass
        for j in range(n):
            if j == i: continue
            arr[i][j] = 1   # add another element
            try:
                _RowEchelon(i+1,copy.copy(arr),collist[:])
                break
            except:
                pass
            arr[i][j] = -1   # try a different valuefor this element
            try:
                _RowEchelon(i+1,copy.copy(arr),collist[:])
                break
            except:
                arr[i][j] = 0   # reset to add another element
        else:
            raise Exception('Unable to create non-singular matrix')
    return arr

def _SwapColumns(i,m,v):
    '''Swap columns in matrix m as well as the labels in v 
    so that element (i,i) is replaced by the first non-zero element in row i after that element

    Throws an exception if there are no non-zero elements in that row
    '''
    for j in range(i+1,len(v)):
        if not np.allclose(m[i,j],0):
            m[:,(i,j)] = m[:,(j,i)]
            v[i],v[j] = v[j],v[i]
            return
    else:
        raise ConstraintException('Singular input')

def _RowEchelon(m,arr,collist):
    '''Convert the first m rows in Matrix arr to row-echelon form
    exchanging columns in the matrix and collist as needed.

    throws an exception if the matrix is singular because
    the first m rows are not linearly independent
    '''
    for i in range(m):
        if np.allclose(arr[i,i],0):
            _SwapColumns(i,arr,collist)
        arr[i,:] /= arr[i,i] # normalize row
        # subtract current row from subsequent rows to set values to left of diagonal to 0
        for j in range(i+1,m):
            arr[j,:] -= arr[i,:] * arr[j,i]

if __name__ == "__main__":
    d = [{'a': 0, 'b': 1.5,'d':0}, {'d': -1}]
    lbls = ['a','c','b','d']
    sel = (0,1)
    try:
        arr2 = _FillArray(sel,d,lbls)
    except Exception as err:
        if 'Initial' in str(err): 
            print('initial error')
        else:
            print('unable to extend matrix error')
        import sys; sys.exit()
    print(arr2)

    d = [{'a': 1, 'b': 1,}]
    lbls = ['a','b']
    sel = (0,)
    arr1 = _FillArray(sel,d,lbls)
    print(arr1)
    print(GramSchmidtOrtho(arr1,1))
    
    import sys; sys.exit()
    
    parmdict = {}
    constrDict = [
        {'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},
        {'0:0:eA': 0.0},
        {'2::C(10,6,1)': 1.0, '1::C(10,6,1)': 1.0},
        {'1::C(10,0,1)': 1.0, '2::C(10,0,1)': 1.0},
        {'1::AUiso:0': 1.0, '0::AUiso:0': 1.0},
        {'0::A0': 0.0}
        ]
    fixedList = ['5.0', '0', None, None, '1.0', '0']
    StoreEquivalence('2::atomx:3',('2::atomy:3', ('2::atomz:3',2,), ))
    #StoreEquivalence('1::atomx:3',('2::atomx:3', ('2::atomz:3',2,), )) # error: dependent & independent vars mixed
    #StoreEquivalence('1::atomx:3',('2::atomy:3', ('2::atomz:3',2,), )) # error: dependent vars repeated
    #StoreEquivalence('0:1:eA',('0:0:eA',)) # error: equiv & fixed
    #StoreEquivalence('0:99:Scale',('0:12:Scale',)) # error: equiv & constrained
    #StoreEquivalence('0:12:Scale',('0:99:Scale',)) # error: equiv & constrained
    varyList = ['2::atomx:3',
                '2::C(10,6,1)', '1::C(10,6,1)',
                '2::atomy:3', '2::atomz:3',
                '0:12:Scale', '0:11:Scale', '0:14:Scale', '0:13:Scale', '0:0:Scale']
#    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']
#    constrDict = [
#        {'0::Afrac:4': 24.0, '0::Afrac:1': 16.0, '0::Afrac:3': 24.0, '0::Afrac:2': 16.0},
#        {'0::Afrac:1': 1.0, '0::Afrac:2': 1.0},
#        {'0::Afrac:4': 1.0, '0::Afrac:3': 1.0}]
#    fixedList = ['40.0', '1.0', '1.0']

    msg,warning,groups,parmlist = GenerateConstraints(varyList,constrDict,fixedList,parmdict)
    print (VarRemapShow(varyList))
    parmdict.update( {
        '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,
        '0:0:eA': 0.0,
        '2::C(10,6,1)': 0.2, '1::C(10,6,1)': 0.3,
        '1::C(10,0,1)': 0.2, '2::C(10,0,1)': 0.3,
        '1::AUiso:0': 0.02, '0::AUiso:0': 0.03,
        '0::A0': 0.0,
        '2::atomx:3':0.23,'2::atomy:3':-.23, '2::atomz:3':-0.11,
        })
    print ('parmdict start',parmdict)
    print ('varylist start',varyList)
    before = parmdict.copy()
    Map2Dict(parmdict,varyList)
    print ('parmdict before and after Map2Dict')
    print ('  key / before / after')
    for key in sorted(list(parmdict.keys())):
        print ('  '+key,'\t',before.get(key),'\t',parmdict[key])
    print ('varylist after',varyList)
    before = parmdict.copy()
    Dict2Map(parmdict)
    print ('after Dict2Map')
    print ('  key / before / after')
    for key in sorted(list(parmdict.keys())):
        print ('  '+key,'\t',before.get(key),'\t',parmdict[key])
#    dMdv = len(varylist)*[0]
#    deriv = {}
#    for i,v in enumerate(parmdict.keys()): deriv[v]=i
#    Dict2Deriv(varylist,deriv,dMdv)
