Changeset 570 for trunk/GSASIImapvars.py
 Timestamp:
 Apr 24, 2012 12:19:19 PM (9 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

trunk/GSASIImapvars.py
r567 r570 58 58 parameters, there must exist NC new parameters, where C is the number 59 59 of contraint equations in the group. Routine GenerateConstraints takes 60 the output from InputParse and GroupConstraintsgenerates the60 the output from GroupConstraints and generates the 61 61 "missing" relationships and saves that information in the module's 62 62 global variables. Each generated parameter is named sequentially using paramPrefix. … … 76 76 External Routines: 77 77 To define a set of constrained and unconstrained relations, one 78 calls In putParse, GroupConstraints and GenerateConstraints. It78 calls InitVars, GroupConstraints and GenerateConstraints. It 79 79 is best to supply all relations/equations in one call to these 80 80 routines so that grouping is done correctly. … … 88 88 (though StoreEquivalence('x',('y','z')) will run more 89 89 efficiently) but mixing independent and dependent variables is 90 problematic. This may not work properly:90 problematic. This is not allowed: 91 91 StoreEquivalence('x',('y',)) 92 92 StoreEquivalence('y',('z',)) 93 Use StoreEquivalence before calling GenerateConstraints or 94 CheckConstraints 95 96 To check that input in internally consistent, use CheckConstraints 93 97 94 98 To show a summary of the parameter remapping, one calls … … 100 104 To take values from the new independent parameters and constraints, 101 105 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 102 111 103 112 Global Variables: … … 159 168 Input 160 169 constrDict: a list of dicts defining relationships/constraints 161 (see InputParse) 170 constrDict = [{<constr1>}, {<constr2>}, ...] 171 where {<constr1>} is {'param1': mult1, 'param2': mult2,...} 162 172 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 164 175 a list containing the parameters used in each group 165 176 """ … … 189 200 return groups,ParmList 190 201 202 def 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 # nonsingular 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 += "\nOverconstrained 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 nonsingular 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 nonsingular 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 191 403 def GenerateConstraints(groups,parmlist,varyList,constrDict,fixedList): 192 404 '''Takes a list of relationship entries comprising a group of … … 194 406 and stores them in global variables Also checks for internal 195 407 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). 196 417 ''' 197 418 global dependentParmList,arrayList,invarrayList,indParmList,consNum … … 413 634 return 414 635 415 # def SetVaryFlags(varyList,fixedList):416 # '''Adds independent variables to the varyList provided that all417 # dependent variables are being varied.418 # Ignores independent variables where no dependent variables are419 # being varied.420 # Returns a nonempty text message where some but not all dependent421 # variables are being varied.422 # '''423 # global dependentParmList,arrayList,invarrayList,indParmList,fixedDict424 # 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: continue430 # 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 " + mv443 # msg += "\n varied: "444 # for v in varied: msg += v445 # msg += "\n not varied: "446 # for v in notvaried: msg += v447 # #elif len(varied) > 0 and multarr is not None:448 # # varyList.append(mv)449 # return msg450 451 636 def GetDependentVars(): 452 637 '''Return a list of dependent variables: e.g. variables that are … … 460 645 def GetIndependentVars(): 461 646 '''Return a list of independent variables: e.g. variables that are 462 created by constrain s of other variables'''647 created by constraints of other variables''' 463 648 independentVars = [] 464 649 global indParmList,fixedDict … … 507 692 returns a dictionary containing the esd values for dependent parameters 508 693 ''' 509 # I think this fails for equivalencies: the ESD for the master (independent variable)510 # needs to be adjusted.511 694 sigmaDict = {} 512 695 for varlist,mapvars,invmultarr in zip(dependentParmList,indParmList,invarrayList): … … 529 712 def FormatConstraint(RelDict,RelVal): 530 713 '''Formats a Constraint or Function for use in a convenient way''' 714 linelen = 45 531 715 s = [""] 532 716 for var,val in RelDict.items(): 533 if len(s[1]) > 60: s.append(' ')717 if len(s[1]) > linelen: s.append(' ') 534 718 m = val 535 719 if s[1] != "" and m >= 0: … … 539 723 m = abs(m) 540 724 s[1] += '%.3f*%s '%(m,var) 541 if len(s[1]) > 60: s.append(' ')725 if len(s[1]) > linelen: s.append(' ') 542 726 if RelVal is None: 543 727 s[1] += ' = New variable' … … 856 1040 '2::atomy:3', '2::atomz:3', 857 1041 '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() 858 1058 groups,parmlist = GroupConstraints(constrDict) 859 1059 GenerateConstraints(groups,parmlist,varylist,constrDict,fixedList)
Note: See TracChangeset
for help on using the changeset viewer.