1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24 """The model-free analysis optimisation functions."""
25
26
27 from minfx.generic import generic_minimise
28 from minfx.grid import grid, grid_point_array
29 from numpy import array, dot, float64
30 import sys
31
32
33 import lib.arg_check
34 from lib.errors import RelaxError, RelaxInfError, RelaxMultiVectorError, RelaxNaNError
35 from lib.float import isNaN, isInf
36 from lib.periodic_table import periodic_table
37 from lib.text.sectioning import subsection
38 from multi import Memo, Result_command, Slave_command
39 from pipe_control import pipes
40 from pipe_control.interatomic import return_interatom_list
41 from pipe_control.mol_res_spin import return_spin, return_spin_from_index, spin_loop
42 from specific_analyses.model_free.parameters import assemble_param_vector, disassemble_param_vector
43 from target_functions.mf import Mf
44
45
46 -def disassemble_result(param_vector=None, func=None, iter=None, fc=None, gc=None, hc=None, warning=None, spin=None, sim_index=None, model_type=None, scaling_matrix=None):
47 """Disassemble the optimisation results.
48
49 @keyword param_vector: The model-free parameter vector.
50 @type param_vector: numpy array
51 @keyword func: The optimised chi-squared value.
52 @type func: float
53 @keyword iter: The number of optimisation steps required to find the minimum.
54 @type iter: int
55 @keyword fc: The function count.
56 @type fc: int
57 @keyword gc: The gradient count.
58 @type gc: int
59 @keyword hc: The Hessian count.
60 @type hc: int
61 @keyword warning: Any optimisation warnings.
62 @type warning: str or None
63 @keyword spin: The spin container.
64 @type spin: SpinContainer instance or None
65 @keyword sim_index: The Monte Carlo simulation index.
66 @type sim_index: int or None
67 @keyword model_type: The model-free model type, one of 'mf', 'local_tm', 'diff', or
68 'all'.
69 @type model_type: str
70 @keyword scaling_matrix: The diagonal, square scaling matrix.
71 @type scaling_matrix: numpy diagonal matrix
72 """
73
74
75 if param_vector == None:
76 return
77
78
79 cdp = pipes.get_pipe()
80
81
82 if isInf(func):
83 raise RelaxInfError('chi-squared')
84
85
86 if isNaN(func):
87 raise RelaxNaNError('chi-squared')
88
89
90 if scaling_matrix != None:
91 param_vector = dot(scaling_matrix, param_vector)
92
93
94 if sim_index == None:
95
96 chi2 = None
97 if (model_type == 'mf' or model_type == 'local_tm') and hasattr(cdp, 'chi2'):
98 chi2 = spin.chi2
99 if (model_type == 'diff' or model_type == 'all') and hasattr(cdp, 'chi2'):
100 chi2 = cdp.chi2
101
102
103 spin_text = ''
104 if spin != None and hasattr(spin, '_spin_ids') and len(spin._spin_ids):
105 spin_text = " for the spin '%s'" % spin._spin_ids[0]
106
107
108 if chi2 != None and func >= chi2:
109 print("Discarding the optimisation results%s, the optimised chi-squared value is higher than the current value (%s >= %s)." % (spin_text, func, chi2))
110
111
112 return
113
114
115 else:
116 print("Storing the optimisation results%s, the optimised chi-squared value is lower than the current value (%s < %s)." % (spin_text, func, chi2))
117
118
119 disassemble_param_vector(model_type, param_vector=param_vector, spin=spin, sim_index=sim_index)
120
121
122 if sim_index != None:
123
124 if model_type == 'mf' or model_type == 'local_tm':
125
126
127 spin.chi2_sim[sim_index] = func
128
129
130 spin.iter_sim[sim_index] = iter
131
132
133 spin.f_count_sim[sim_index] = fc
134
135
136 spin.g_count_sim[sim_index] = gc
137
138
139 spin.h_count_sim[sim_index] = hc
140
141
142 spin.warning_sim[sim_index] = warning
143
144
145 elif model_type == 'diff' or model_type == 'all':
146
147 cdp.chi2_sim[sim_index] = func
148
149
150 cdp.iter_sim[sim_index] = iter
151
152
153 cdp.f_count_sim[sim_index] = fc
154
155
156 cdp.g_count_sim[sim_index] = gc
157
158
159 cdp.h_count_sim[sim_index] = hc
160
161
162 cdp.warning_sim[sim_index] = warning
163
164
165 else:
166
167 if model_type == 'mf' or model_type == 'local_tm':
168
169 spin.chi2 = func
170
171
172 spin.iter = iter
173
174
175 spin.f_count = fc
176
177
178 spin.g_count = gc
179
180
181 spin.h_count = hc
182
183
184 spin.warning = warning
185
186
187 elif model_type == 'diff' or model_type == 'all':
188
189 cdp.chi2 = func
190
191
192 cdp.iter = iter
193
194
195 cdp.f_count = fc
196
197
198 cdp.g_count = gc
199
200
201 cdp.h_count = hc
202
203
204 cdp.warning = warning
205
206
207 -def minimise_data_setup(data_store, min_algor, num_data_sets, min_options, spin=None, sim_index=None):
208 """Set up all the data required for minimisation.
209
210 @param data_store: A data storage container.
211 @type data_store: class instance
212 @param min_algor: The minimisation algorithm to use.
213 @type min_algor: str
214 @param num_data_sets: The number of data sets.
215 @type num_data_sets: int
216 @param min_options: The minimisation options array.
217 @type min_options: list
218 @keyword spin: The spin data container.
219 @type spin: SpinContainer instance
220 @keyword sim_index: The optional MC simulation index.
221 @type sim_index: int
222 @return: An insane tuple. The full tuple is (ri_data, ri_data_err, equations, param_types, param_values, r, csa, num_frq, frq, num_ri, remap_table, noe_r1_table, ri_types, num_params, xh_unit_vectors, diff_type, diff_params)
223 @rtype: tuple
224 """
225
226
227 data_store.ri_data = []
228 data_store.ri_data_err = []
229 data_store.equations = []
230 data_store.param_types = []
231 data_store.param_values = None
232 data_store.r = []
233 data_store.csa = []
234 data_store.num_frq = []
235 data_store.frq = []
236 data_store.num_ri = []
237 data_store.remap_table = []
238 data_store.noe_r1_table = []
239 data_store.ri_types = []
240 data_store.gx = []
241 data_store.gh = []
242 data_store.num_params = []
243 data_store.xh_unit_vectors = []
244 if data_store.model_type == 'local_tm':
245 data_store.mf_params = []
246 elif data_store.model_type == 'diff':
247 data_store.param_values = []
248
249
250 if min_algor == 'back_calc':
251
252 data_store.ri_data = [0.0]
253 data_store.ri_data_err = [0.000001]
254 data_store.equations = [spin.equation]
255 data_store.param_types = [spin.params]
256 data_store.csa = [spin.csa]
257 data_store.num_frq = [1]
258 data_store.frq = [[min_options[3]]]
259 data_store.num_ri = [1]
260 data_store.remap_table = [[0]]
261 data_store.noe_r1_table = [[None]]
262 data_store.ri_types = [[min_options[2]]]
263 data_store.gx = [periodic_table.gyromagnetic_ratio(spin.isotope)]
264
265
266 interatoms = return_interatom_list(data_store.spin_id)
267 for i in range(len(interatoms)):
268
269 if not interatoms[i].dipole_pair:
270 continue
271
272
273 if data_store.spin_id != interatoms[i].spin_id1:
274 spin_id2 = interatoms[i].spin_id1
275 else:
276 spin_id2 = interatoms[i].spin_id2
277 spin2 = return_spin(spin_id2)
278
279
280 data_store.r = [interatoms[i].r]
281 data_store.gh = [periodic_table.gyromagnetic_ratio(spin2.isotope)]
282 if data_store.model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere':
283 data_store.xh_unit_vectors = [interatoms[i].vector]
284 else:
285 data_store.xh_unit_vectors = [None]
286
287
288 data_store.num_params = [len(spin.params)]
289
290
291 for j in range(num_data_sets):
292
293 if data_store.model_type == 'diff' or data_store.model_type == 'all':
294 spin_index = j
295 spin, data_store.spin_id = return_spin_from_index(global_index=spin_index, return_spin_id=True)
296
297
298 if not spin.select:
299 continue
300
301
302 if not hasattr(spin, 'ri_data') or not hasattr(spin, 'ri_data_err'):
303 continue
304
305
306 for ri_id in cdp.ri_ids:
307
308 if not ri_id in spin.ri_data_err:
309 continue
310
311
312 err = spin.ri_data_err[ri_id]
313
314
315 if err != None and err == 0.0:
316 raise RelaxError("Zero error for spin '%s' for the relaxation data ID '%s', minimisation not possible." % (data_store.spin_id, ri_id))
317 elif err != None and err < 0.0:
318 raise RelaxError("Negative error of %s for spin '%s' for the relaxation data ID '%s', minimisation not possible." % (err, data_store.spin_id, ri_id))
319
320
321 data = relax_data_opt_structs(spin, sim_index=sim_index)
322
323
324 data_store.ri_data.append(data[0])
325 data_store.ri_data_err.append(data[1])
326 data_store.num_frq.append(data[2])
327 data_store.num_ri.append(data[3])
328 data_store.ri_types.append(data[4])
329 data_store.frq.append(data[5])
330 data_store.remap_table.append(data[6])
331 data_store.noe_r1_table.append(data[7])
332 if sim_index == None or data_store.model_type == 'diff':
333 data_store.csa.append(spin.csa)
334 else:
335 data_store.csa.append(spin.csa_sim[sim_index])
336
337
338 data_store.equations.append(spin.equation)
339 data_store.param_types.append(spin.params)
340 data_store.gx.append(periodic_table.gyromagnetic_ratio(spin.isotope))
341
342
343 interatoms = return_interatom_list(data_store.spin_id)
344 for i in range(len(interatoms)):
345
346 if not interatoms[i].dipole_pair:
347 continue
348
349
350 if data_store.spin_id != interatoms[i].spin_id1:
351 spin_id2 = interatoms[i].spin_id1
352 else:
353 spin_id2 = interatoms[i].spin_id2
354 spin2 = return_spin(spin_id2)
355
356
357 data_store.gh.append(periodic_table.gyromagnetic_ratio(spin2.isotope))
358 if sim_index == None or data_store.model_type == 'diff' or not hasattr(interatoms[i], 'r_sim'):
359 data_store.r.append(interatoms[i].r)
360 else:
361 data_store.r.append(interatoms[i].r_sim[sim_index])
362
363
364 if data_store.model_type != 'local_tm' and cdp.diff_tensor.type != 'sphere':
365
366 if lib.arg_check.is_num_list(interatoms[i].vector[0], raise_error=False):
367 raise RelaxMultiVectorError(data_store.spin_id)
368
369
370 data_store.xh_unit_vectors.append(interatoms[i].vector)
371
372
373 else:
374 data_store.xh_unit_vectors.append(None)
375
376
377 break
378
379
380 if data_store.model_type == 'local_tm':
381 pass
382
383
384 data_store.num_params.append(len(spin.params))
385
386
387 if data_store.model_type == 'diff':
388 data_store.param_values.append(assemble_param_vector(model_type='mf'))
389
390
391 for k in range(len(data_store.ri_data)):
392 data_store.ri_data[k] = array(data_store.ri_data[k], float64)
393 data_store.ri_data_err[k] = array(data_store.ri_data_err[k], float64)
394
395
396 if data_store.model_type == 'local_tm':
397 data_store.diff_type = 'sphere'
398 else:
399 data_store.diff_type = cdp.diff_tensor.type
400
401
402 data_store.diff_params = None
403 if data_store.model_type == 'mf':
404
405 if data_store.diff_type == 'sphere':
406 data_store.diff_params = [cdp.diff_tensor.tm]
407
408
409 elif data_store.diff_type == 'spheroid':
410 data_store.diff_params = [cdp.diff_tensor.tm, cdp.diff_tensor.Da, cdp.diff_tensor.theta, cdp.diff_tensor.phi]
411
412
413 elif data_store.diff_type == 'ellipsoid':
414 data_store.diff_params = [cdp.diff_tensor.tm, cdp.diff_tensor.Da, cdp.diff_tensor.Dr, cdp.diff_tensor.alpha, cdp.diff_tensor.beta, cdp.diff_tensor.gamma]
415 elif min_algor == 'back_calc' and data_store.model_type == 'local_tm':
416
417 data_store.diff_params = [spin.local_tm]
418
419
421 """Package the relaxation data into the data structures used for optimisation.
422
423 @param spin: The spin container to extract the data from.
424 @type spin: SpinContainer instance
425 @keyword sim_index: The optional MC simulation index.
426 @type sim_index: int
427 @return: The structures ri_data, ri_data_err, num_frq, num_ri, ri_ids, frq, remap_table, noe_r1_table.
428 @rtype: tuple
429 """
430
431
432 ri_data = []
433 ri_data_err = []
434 ri_labels = []
435 frq = []
436 remap_table = []
437 noe_r1_table = []
438
439
440 for ri_id in cdp.ri_ids:
441
442 if ri_id not in spin.ri_data:
443 continue
444
445
446 if sim_index == None:
447 data = spin.ri_data[ri_id]
448 else:
449 data = spin.ri_data_sim[ri_id][sim_index]
450
451
452 err = spin.ri_data_err[ri_id]
453
454
455 if data == None or err == None:
456 continue
457
458
459 ri_data.append(data)
460 ri_data_err.append(err)
461
462
463 ri_labels.append(cdp.ri_type[ri_id])
464
465
466 if cdp.spectrometer_frq[ri_id] not in frq:
467 frq.append(cdp.spectrometer_frq[ri_id])
468
469
470 remap_table.append(frq.index(cdp.spectrometer_frq[ri_id]))
471
472
473 noe_r1_table.append(None)
474
475
476 num_ri = len(ri_data)
477
478
479 for i in range(num_ri):
480
481 if cdp.ri_type[cdp.ri_ids[i]] == 'NOE':
482 for j in range(num_ri):
483 if cdp.ri_type[cdp.ri_ids[j]] == 'R1' and cdp.spectrometer_frq[cdp.ri_ids[i]] == cdp.spectrometer_frq[cdp.ri_ids[j]]:
484 noe_r1_table[i] = j
485
486
487 return ri_data, ri_data_err, len(frq), num_ri, ri_labels, frq, remap_table, noe_r1_table
488
489
490
492 """The model-free memo class.
493
494 Not quite a momento so a memo.
495 """
496
497 - def __init__(self, model_free=None, model_type=None, spin=None, sim_index=None, scaling_matrix=None):
498 """Initialise the model-free memo class.
499
500 This memo stores the model-free class instance so that the disassemble_result() method can be called to store the optimisation results. The other args are those required by this method but not generated through optimisation.
501
502 @keyword model_free: The model-free class instance.
503 @type model_free: specific_analyses.model_free.Model_free instance
504 @keyword spin: The spin data container. If this argument is supplied, then the spin_id argument will be ignored.
505 @type spin: SpinContainer instance
506 @keyword sim_index: The optional MC simulation index.
507 @type sim_index: int
508 @keyword scaling_matrix: The diagonal, square scaling matrix.
509 @type scaling_matrix: numpy diagonal matrix
510 """
511
512
513 super(MF_memo, self).__init__()
514
515
516 self.model_free = model_free
517 self.model_type = model_type
518 self.spin = spin
519 self.sim_index = sim_index
520 self.scaling_matrix = scaling_matrix
521
522
523
525 """Command class for standard model-free minimisation."""
526
532
533
535 """Model-free optimisation.
536
537 @return: The optimisation results consisting of the parameter vector, function value, iteration count, function count, gradient count, Hessian count, and warnings.
538 @rtype: tuple of numpy array, float, int, int, int, int, str
539 """
540
541
542 results = generic_minimise(func=self.mf.func, dfunc=self.mf.dfunc, d2func=self.mf.d2func, args=(), x0=self.opt_params.param_vector, min_algor=self.opt_params.min_algor, min_options=self.opt_params.min_options, func_tol=self.opt_params.func_tol, grad_tol=self.opt_params.grad_tol, maxiter=self.opt_params.max_iterations, A=self.opt_params.A, b=self.opt_params.b, full_output=True, print_flag=self.opt_params.verbosity)
543
544
545 return results
546
547
548 - def run(self, processor, completed):
549 """Setup and perform the model-free optimisation."""
550
551
552 self.mf = Mf(init_params=self.opt_params.param_vector, model_type=self.data.model_type, diff_type=self.data.diff_type, diff_params=self.data.diff_params, scaling_matrix=self.data.scaling_matrix, num_spins=self.data.num_spins, equations=self.data.equations, param_types=self.data.param_types, param_values=self.data.param_values, relax_data=self.data.ri_data, errors=self.data.ri_data_err, bond_length=self.data.r, csa=self.data.csa, num_frq=self.data.num_frq, frq=self.data.frq, num_ri=self.data.num_ri, remap_table=self.data.remap_table, noe_r1_table=self.data.noe_r1_table, ri_labels=self.data.ri_types, gx=self.data.gx, gh=self.data.gh, h_bar=self.data.h_bar, mu0=self.data.mu0, num_params=self.data.num_params, vectors=self.data.xh_unit_vectors)
553
554
555 if self.opt_params.verbosity >= 1 and (self.data.model_type == 'mf' or self.data.model_type == 'local_tm'):
556 subsection(file=sys.stdout, text="Optimisation: Spin '%s'" % self.data.spin_id, prespace=2, postspace=0)
557
558
559 results = self.optimise()
560
561
562 param_vector, func, iter, fc, gc, hc, warning = results
563
564 processor.return_object(MF_result_command(processor, self.memo_id, param_vector, func, iter, fc, gc, hc, warning, completed=False))
565
566
568 """Store all the data required for model-free optimisation.
569
570 @param data: The data used to initialise the model-free target function class.
571 @type data: class instance
572 @param opt_params: The parameters and data required for optimisation using minfx.
573 @type opt_params: class instance
574 """
575
576
577 self.data = data
578 self.opt_params = opt_params
579
580
581
583 """Command class for the model-free grid search."""
584
590
591
593 """Model-free grid search.
594
595 @return: The optimisation results consisting of the parameter vector, function value, iteration count, function count, gradient count, Hessian count, and warnings.
596 @rtype: tuple of numpy array, float, int, int, int, int, str
597 """
598
599
600 if not hasattr(self.opt_params, 'subdivision'):
601 results = grid(func=self.mf.func, args=(), num_incs=self.opt_params.inc, lower=self.opt_params.lower, upper=self.opt_params.upper, A=self.opt_params.A, b=self.opt_params.b, verbosity=self.opt_params.verbosity)
602
603
604 else:
605 results = grid_point_array(func=self.mf.func, args=(), points=self.opt_params.subdivision, verbosity=self.opt_params.verbosity)
606
607
608 param_vector, func, iter, warning = results
609 fc = iter
610 gc = 0.0
611 hc = 0.0
612
613
614 return param_vector, func, iter, fc, gc, hc, warning
615
616
617
619 """Class for processing the model-free results."""
620
621 - def __init__(self, processor, memo_id, param_vector, func, iter, fc, gc, hc, warning, completed):
636
637
638 - def run(self, processor, memo):
639 """Disassemble the model-free optimisation results.
640
641 @param processor: Unused!
642 @type processor: None
643 @param memo: The model-free memo.
644 @type memo: memo
645 """
646
647
648 disassemble_result(param_vector=self.param_vector, func=self.func, iter=self.iter, fc=self.fc, gc=self.gc, hc=self.hc, warning=self.warning, spin=memo.spin, sim_index=memo.sim_index, model_type=memo.model_type, scaling_matrix=memo.scaling_matrix)
649