Changeset 570 for trunk/GSASIImapvars.py


Ignore:
Timestamp:
Apr 24, 2012 12:19:19 PM (9 years ago)
Author:
toby
Message:

Add checks on constraints: if problems display error window when displayed, added/edited & on Refine; handle GPX w/o doPawley flag

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/GSASIImapvars.py

    r567 r570  
    5858parameters, there must exist N-C new parameters, where C is the number
    5959of contraint equations in the group. Routine GenerateConstraints takes
    60 the output from InputParse and GroupConstraints generates the
     60the output from GroupConstraints and generates the
    6161"missing" relationships and saves that information in the module's
    6262global variables. Each generated parameter is named sequentially using paramPrefix.
     
    7676External Routines:
    7777   To define a set of constrained and unconstrained relations, one
    78      calls InputParse, GroupConstraints and GenerateConstraints. It
     78     calls InitVars, GroupConstraints and GenerateConstraints. It
    7979     is best to supply all relations/equations in one call to these
    8080     routines so that grouping is done correctly.
     
    8888     (though StoreEquivalence('x',('y','z')) will run more
    8989     efficiently) but mixing independent and dependent variables is
    90      problematic. This may not work properly:
     90     problematic. This is not allowed:
    9191        StoreEquivalence('x',('y',))
    9292        StoreEquivalence('y',('z',))
     93   Use StoreEquivalence before calling GenerateConstraints or
     94      CheckConstraints
     95
     96   To check that input in internally consistent, use CheckConstraints
    9397
    9498   To show a summary of the parameter remapping, one calls
     
    100104   To take values from the new independent parameters and constraints,
    101105      one calls Dict2Map. This will apply contraints.
     106
     107   Use Dict2Deriv to determine derivatives on independent parameters
     108      from those on dependent ones
     109     
     110   Use ComputeDepESD to compute uncertainties on dependent variables
    102111
    103112Global Variables:
     
    159168    Input
    160169       constrDict: a list of dicts defining relationships/constraints
    161        (see InputParse)
     170       constrDict = [{<constr1>}, {<constr2>}, ...]
     171       where {<constr1>} is {'param1': mult1, 'param2': mult2,...}
    162172    Returns two lists of lists:
    163       a list of relationship list indices for each group and
     173      a list of relationship list indices for each group pointing to
     174        elements in constrDict and
    164175      a list containing the parameters used in each group
    165176      """
     
    189200    return groups,ParmList
    190201
     202def CheckConstraints(varyList,constrDict,fixedList):
     203    '''Takes a list of relationship entries comprising a group of
     204    constraints and checks for inconsistencies such as conflicts in
     205    parameter/variable definitions and or inconsistently varied parameters.
     206    Input: see GenerateConstraints
     207    Output: returns two strings
     208      the first lists conflicts internal to the specified constraints
     209      the second lists conflicts where the varyList specifies some
     210        parameters in a constraint, but not all
     211      If there are no errors, both strings will be empty
     212    '''
     213    global dependentParmList,arrayList,invarrayList,indParmList,consNum
     214    errmsg = ''
     215    warnmsg = ''
     216    fixVlist = []
     217    # process fixed (held) variables
     218    for cdict in constrDict:
     219        if len(cdict) == 1:
     220            fixVlist.append(cdict.keys()[0])
     221   
     222    # process equivalences: make a list of dependent and independent vars
     223    #    and check for repeated uses (repetition of a parameter as an
     224    #    independent var is OK)
     225    indepVarList = []
     226    depVarList = []
     227    multdepVarList = []
     228    for varlist,mapvars,multarr,invmultarr in zip(
     229        dependentParmList,indParmList,arrayList,invarrayList):
     230        if multarr is None: # an equivalence
     231            zeromult = False
     232            for mv in mapvars:
     233                varied = 0
     234                notvaried = ''
     235                if mv in varyList:
     236                    varied += 1
     237                else:
     238                    if notvaried: notvaried += ', '
     239                    notvaried += mv
     240                if mv not in indepVarList: indepVarList.append(mv)
     241                for v,m in zip(varlist,invmultarr):
     242                    if m == 0: zeromult = True
     243                    if v in varyList:
     244                        varied += 1
     245                    else:
     246                        if notvaried: notvaried += ', '
     247                        notvaried += v
     248                    if v in depVarList:
     249                        multdepVarList.append(v)
     250                    else:
     251                        depVarList.append(v)
     252            if varied > 0 and varied != len(varlist)+1:
     253                warnmsg += "\nNot all variables refined in equivalence:\n\t"
     254                s = ""
     255                for v in varlist:
     256                    if s != "": s+= " & "
     257                    s += str(v)           
     258                warnmsg += str(mv) + " => " + s
     259                warnmsg += '\nNot refined: ' + notvaried + '\n'
     260            if zeromult:
     261                errmsg += "\nZero multiplier is invalid in equivalence:\n\t"
     262                s = ""
     263                for v in varlist:
     264                    if s != "": s+= " & "
     265                    s += str(v)           
     266                errmsg += str(mv) + " => " + s + '\n'
     267
     268    #print 'indepVarList',indepVarList
     269    #print 'depVarList',depVarList
     270    # check for errors:
     271    inboth = set(indepVarList).intersection(set(depVarList))
     272    if len(inboth) > 0:
     273        errmsg += "\nThe following parameters(s) are used as both dependent and independent variables in Equivalence relations:\n"
     274        s = ''
     275        for var in inboth:
     276            if s != "": s+= ", "
     277            s += str(v)
     278        errmsg += '\t'+ s + '\n'
     279    if len(multdepVarList) > 0:
     280        errmsg += "\nThe following parameters(s) are used in multiple Equivalence relations as dependent variables:\n"
     281        s = ''
     282        for var in multdepVarList:
     283            if s != "": s+= ", "
     284            s += str(v)           
     285        errmsg += '\t'+ s + '\n'
     286    equivVarList = list(set(indepVarList).union(set(depVarList)))
     287    #print 'equivVarList',equivVarList
     288    inboth = set(fixVlist).intersection(set(equivVarList))
     289    if len(inboth) > 0:
     290        errmsg += "\nThe following parameter(s) are used in both Equivalence and Fixed constraints:\n"
     291        s = ''
     292        for var in inboth:
     293            if s != "": s+= ", "
     294            s += str(v)
     295        errmsg += '\t'+ s + '\n'
     296
     297    groups,parmlist = GroupConstraints(constrDict)
     298    # scan through parameters in each relationship. Are all varied? If only some are
     299    # varied, create a warning message.
     300    for group,varlist in zip(groups,parmlist):
     301        if len(varlist) == 1: continue
     302        VaryFree = False
     303        for rel in group:
     304            varied = 0
     305            notvaried = ''
     306            for var in constrDict[rel]:
     307                if var in varyList:
     308                    varied += 1
     309                else:
     310                    if notvaried: notvaried += ', '
     311                    notvaried += var
     312                if var in fixVlist:
     313                    errmsg += '\nParameter '+var+" is Fixed and used in a constraint:\n\t"
     314                    errmsg += FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
     315                if var in equivVarList:
     316                    errmsg += '\nParameter '+var+" is Equivalenced and used in a constraint:\n\t"
     317                    errmsg += FormatConstraint(constrDict[rel],fixedList[rel])+"\n"
     318            if varied > 0 and varied != len(constrDict[rel]):
     319                warnmsg += "\nNot all variables refined in constraint:\n\t"
     320                warnmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     321                warnmsg += '\nNot refined: ' + notvaried + '\n'
     322    if errmsg or warnmsg: return errmsg,warnmsg
     323
     324    # now look for process each group and create the relations that are needed to form
     325    # non-singular square matrix
     326    for group,varlist in zip(groups,parmlist):
     327        if len(varlist) == 1: continue
     328        try:
     329            arr = MakeArray(constrDict, group, varlist)
     330        except:
     331            errmsg += "\nOver-constrained input. "
     332            errmsg += "There are more constraints" + str(len(group))
     333            errmsg += "\n\tthan variables" + str(len(varlist)) + "\n"
     334            for rel in group:
     335                errmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     336                errmsg += "\n"
     337        try:
     338            constrArr = FillInMissingRelations(arr,len(group))
     339        except Exception,msg:
     340            if msg.message.startswith('VectorProduct'):
     341                errmsg += "\nSingular input. "
     342                errmsg += "This set of constraints is not linearly independent:\n\n"
     343            else:
     344                errmsg += "\nInconsistent input. "
     345                errmsg += "Cound not generate a set of non-singular constraints"
     346                errmsg += "\n\tfor this set of input:\n\n"
     347            for rel in group:
     348                errmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     349                errmsg += "\n"
     350            #import traceback
     351            #print traceback.format_exc()
     352            continue
     353       
     354        mapvar = []
     355        group = group[:]
     356        # scan through all generated and input variables
     357        # Check again for inconsistent variable use
     358        # for new variables -- where varied and unvaried parameters get grouped
     359        # together. I don't think this can happen when not flagged before, but
     360        # it does not hurt to check again.
     361        for i in range(len(varlist)):
     362            varied = 0
     363            notvaried = ''
     364            if len(group) > 0:
     365                rel = group.pop(0)
     366                fixedval = fixedList[rel]
     367                for var in constrDict[rel]:
     368                    if var in varyList:
     369                        varied += 1
     370                    else:
     371                        if notvaried: notvaried += ', '
     372                        notvaried += var
     373            else:
     374                fixedval = None
     375            if fixedval is None:
     376                varname = paramPrefix + str(consNum)
     377                mapvar.append(varname)
     378                consNum += 1
     379                if VaryFree or varied > 0:
     380                    varyList.append(varname)
     381            else:
     382                mapvar.append(fixedval)
     383            if varied > 0 and notvaried != '':
     384                warnmsg += "\nNot all variables refined in generated constraint"
     385                warnmsg += '\nPlease report this unexpected error\n'
     386                for rel in group:
     387                    warnmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     388                    warnmsg += "\n"
     389                warnmsg += '\n\tNot refined: ' + notvaried + '\n'
     390        try:
     391            np.linalg.inv(constrArr)
     392        except:
     393            errmsg += "\nSingular input. "
     394            errmsg += "The following constraints are not "
     395            errmsg += "linearly independent\n\tor do not "
     396            errmsg += "allow for generation of a non-singular set\n"
     397            errmsg += 'Please report this unexpected error\n'
     398            for rel in group:
     399                errmsg += FormatConstraint(constrDict[rel],fixedList[rel])
     400                errmsg += "\n"
     401    return errmsg,warnmsg
     402
    191403def GenerateConstraints(groups,parmlist,varyList,constrDict,fixedList):
    192404    '''Takes a list of relationship entries comprising a group of
     
    194406    and stores them in global variables Also checks for internal
    195407    conflicts or inconsistencies in parameter/variable definitions.
     408    Input:
     409       groups,parmlist: see GroupConstraints
     410       varyList: a list of parameters that will be varied
     411       constrDict: a list of dicts defining relationships/constraints
     412         constrDict = [{<constr1>}, {<constr2>}, ...]
     413         where {<constr1>} is {'param1': mult1, 'param2': mult2,...}
     414       fixedList: a list of values for each constraint/variable in
     415          constrDict, value is either a float (contraint) or None (for
     416          a new variable).
    196417    '''
    197418    global dependentParmList,arrayList,invarrayList,indParmList,consNum
     
    413634    return
    414635
    415 # def SetVaryFlags(varyList,fixedList):
    416 #     '''Adds independent variables to the varyList provided that all
    417 #     dependent variables are being varied.
    418 #     Ignores independent variables where no dependent variables are
    419 #     being varied.
    420 #     Returns a non-empty text message where some but not all dependent
    421 #     variables are being varied.
    422 #     '''
    423 #     global dependentParmList,arrayList,invarrayList,indParmList,fixedDict
    424 #     msg = ""
    425 #     for varlist,mapvars,multarr in zip(dependentParmList,indParmList,arrayList):
    426 #         for mv in mapvars:
    427 #             varied = []
    428 #             notvaried = []
    429 #             if mv in fixedDict: continue
    430 #             if multarr is None:
    431 #                 if mv in varyList:
    432 #                     varied.append(mv)
    433 #                 else:
    434 #                     notvaried.append(mv)
    435 #             for v in varlist:
    436 #                 if v in varyList:
    437 #                     varied.append(v)
    438 #                 else:
    439 #                     notvaried.append(v)
    440 #             if len(varied) > 0 and len(notvaried) > 0:
    441 #                 if msg != "": msg += "\n"
    442 #                 msg += "Error: inconsistent use of parameter " + mv
    443 #                 msg += "\n  varied: "
    444 #                 for v in varied: msg += v
    445 #                 msg += "\n  not varied: "
    446 #                 for v in notvaried: msg += v   
    447 #             #elif len(varied) > 0 and multarr is not None:
    448 #             #    varyList.append(mv)
    449 #     return msg
    450 
    451636def GetDependentVars():
    452637    '''Return a list of dependent variables: e.g. variables that are
     
    460645def GetIndependentVars():
    461646    '''Return a list of independent variables: e.g. variables that are
    462     created by constrains of other variables'''
     647    created by constraints of other variables'''
    463648    independentVars = []
    464649    global indParmList,fixedDict
     
    507692    returns a dictionary containing the esd values for dependent parameters
    508693    '''
    509     # I think this fails for equivalencies: the ESD for the master (independent variable)
    510     # needs to be adjusted.
    511694    sigmaDict = {}
    512695    for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList):
     
    529712def FormatConstraint(RelDict,RelVal):
    530713    '''Formats a Constraint or Function for use in a convenient way'''
     714    linelen = 45
    531715    s = [""]
    532716    for var,val in RelDict.items():
    533         if len(s[-1]) > 60: s.append(' ')
     717        if len(s[-1]) > linelen: s.append(' ')
    534718        m = val
    535719        if s[-1] != "" and m >= 0:
     
    539723            m = abs(m)
    540724        s[-1] += '%.3f*%s '%(m,var)
    541     if len(s[-1]) > 60: s.append(' ')
     725    if len(s[-1]) > linelen: s.append(' ')
    542726    if RelVal is None:
    543727        s[-1] += ' = New variable'
     
    8561040                '2::atomy:3', '2::atomz:3',
    8571041                '0:12:Scale', '0:11:Scale', '0:14:Scale', '0:13:Scale', '0:0:Scale']
     1042    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']
     1043    constrDict = [
     1044        {'0::Afrac:4': 24.0, '0::Afrac:1': 16.0, '0::Afrac:3': 24.0, '0::Afrac:2': 16.0},
     1045        {'0::Afrac:1': 1.0, '0::Afrac:2': 1.0},
     1046        {'0::Afrac:4': 1.0, '0::Afrac:3': 1.0}]
     1047    fixedList = ['40.0', '1.0', '1.0']
     1048
     1049    errmsg, warnmsg = CheckConstraints(varylist,constrDict,fixedList)
     1050    if errmsg:
     1051        print "*** Error ********************"
     1052        print errmsg
     1053    if warnmsg:
     1054        print "*** Warning ********************"
     1055        print warnmsg
     1056    if errmsg or warnmsg:
     1057        sys.exit()
    8581058    groups,parmlist = GroupConstraints(constrDict)
    8591059    GenerateConstraints(groups,parmlist,varylist,constrDict,fixedList)
Note: See TracChangeset for help on using the changeset viewer.