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