Package maths_fns :: Module mf
[hide private]
[frames] | no frames]

Source Code for Module maths_fns.mf

   1  ############################################################################### 
   2  #                                                                             # 
   3  # Copyright (C) 2003-2005 Edward d'Auvergne                                   # 
   4  #                                                                             # 
   5  # This file is part of the program relax.                                     # 
   6  #                                                                             # 
   7  # relax is free software; you can redistribute it and/or modify               # 
   8  # it under the terms of the GNU General Public License as published by        # 
   9  # the Free Software Foundation; either version 2 of the License, or           # 
  10  # (at your option) any later version.                                         # 
  11  #                                                                             # 
  12  # relax is distributed in the hope that it will be useful,                    # 
  13  # but WITHOUT ANY WARRANTY; without even the implied warranty of              # 
  14  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               # 
  15  # GNU General Public License for more details.                                # 
  16  #                                                                             # 
  17  # You should have received a copy of the GNU General Public License           # 
  18  # along with relax; if not, write to the Free Software                        # 
  19  # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA   # 
  20  #                                                                             # 
  21  ############################################################################### 
  22   
  23   
  24  from Numeric import Float64, matrixmultiply, ones, sum, transpose, zeros 
  25  from math import pi 
  26  import sys 
  27   
  28  from direction_cosine import * 
  29  from weights import * 
  30  from correlation_time import * 
  31  from jw_mf_comps import * 
  32  from jw_mf import * 
  33  from ri_comps import * 
  34  from ri_prime import * 
  35  from ri import * 
  36  from chi2 import * 
  37   
  38   
39 -class Mf:
40 - def __init__(self, init_params=None, param_set=None, diff_type=None, diff_params=None, scaling_matrix=None, num_res=None, equations=None, param_types=None, param_values=None, relax_data=None, errors=None, bond_length=None, csa=None, num_frq=0, frq=None, num_ri=None, remap_table=None, noe_r1_table=None, ri_labels=None, gx=0, gh=0, g_ratio=0, h_bar=0, mu0=0, num_params=None, vectors=None):
41 """The model-free minimisation class. 42 43 This class should be initialised before every calculation. 44 45 46 Arguments 47 ~~~~~~~~~ 48 49 equation: The model-free equation string which should be either 'mf_orig' or 'mf_ext'. 50 51 param_types: An array of the parameter types used in minimisation. 52 53 relax_data: An array containing the experimental relaxation values. 54 55 errors: An array containing the experimental errors. 56 57 bond_length: The fixed bond length in meters. 58 59 csa: The fixed CSA value. 60 61 diff_type: The diffusion tensor string which should be either 'sphere', 'spheroid', or 62 'ellipsoid'. 63 64 diff_params: An array with the diffusion parameters. 65 66 scaling_matrix: A diagonal matrix of scaling factors. 67 68 69 70 Additional layer of equations to simplify the relaxation equations, gradients, and Hessians. 71 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 72 73 The R1 and R2 equations are left alone, while the NOE is calculated from the R1 and 74 sigma_noe values. 75 76 77 The relaxation equations 78 ~~~~~~~~~~~~~~~~~~~~~~~~ 79 80 Data structure: data.ri 81 Dimension: 1D, (relaxation data) 82 Type: Numeric array, Float64 83 Dependencies: data.ri_prime 84 Required by: data.chi2, data.dchi2, data.d2chi2 85 86 87 R1() = R1'() 88 89 90 R2() = R2'() 91 92 gH sigma_noe() 93 NOE() = 1 + -- . ----------- 94 gN R1() 95 96 97 98 The relaxation gradients 99 ~~~~~~~~~~~~~~~~~~~~~~~~ 100 101 Data structure: data.dri 102 Dimension: 2D, (parameters, relaxation data) 103 Type: Numeric array, Float64 104 Dependencies: data.ri_prime, data.dri_prime 105 Required by: data.dchi2, data.d2chi2 106 107 108 dR1() dR1'() 109 ------- = ------- 110 dthetaj dthetaj 111 112 113 dR2() dR2'() 114 ------- = ------- 115 dthetaj dthetaj 116 117 118 dNOE() gH 1 / dsigma_noe() dR1() \ 119 ------- = -- . ------- . | R1() . ------------ - sigma_noe() . ------- | 120 dthetaj gN R1()**2 \ dthetaj dthetaj / 121 122 123 124 The relaxation Hessians 125 ~~~~~~~~~~~~~~~~~~~~~~~ 126 127 Data structure: data.d2ri 128 Dimension: 3D, (parameters, parameters, relaxation data) 129 Type: Numeric array, Float64 130 Dependencies: data.ri_prime, data.dri_prime, data.d2ri_prime 131 Required by: data.d2chi2 132 133 134 d2R1() d2R1'() 135 --------------- = --------------- 136 dthetai.dthetaj dthetai.dthetaj 137 138 139 d2R2() d2R2'() 140 --------------- = --------------- 141 dthetai.dthetaj dthetai.dthetaj 142 143 144 d2NOE() gH 1 / / dR1() dR1() d2R1() \ 145 --------------- = -- . ------- . | sigma_noe() . | 2 . ------- . ------- - R1() . --------------- | 146 dthetai.dthetaj gN R1()**3 \ \ dthetai dthetaj dthetai.dthetaj / 147 148 / dsigma_noe() dR1() dR1() dsigma_noe() d2sigma_noe() \ \ 149 - R1() . | ------------ . ------- + ------- . ------------ - R1() . --------------- | | 150 \ dthetai dthetaj dthetai dthetaj dthetai.dthetaj / / 151 152 153 154 The chi-sqared equation 155 ~~~~~~~~~~~~~~~~~~~~~~~ 156 _n_ 157 \ (Ri - Ri()) ** 2 158 Chi2 = > ---------------- 159 /__ sigma_i ** 2 160 i=1 161 162 where: 163 Ri are the values of the measured relaxation data set. 164 Ri() are the values of the back calculated relaxation data set. 165 sigma_i are the values of the error set. 166 167 168 The chi-sqared gradient 169 ~~~~~~~~~~~~~~~~~~~~~~~ 170 _n_ 171 dChi2 \ / Ri - Ri() dRi() \ 172 ------- = -2 > | ---------- . ------- | 173 dthetaj /__ \ sigma_i**2 dthetaj / 174 i=1 175 176 where: 177 Ri are the values of the measured relaxation data set. 178 Ri() are the values of the back calculated relaxation data set. 179 sigma_i are the values of the error set. 180 181 182 The chi-sqared Hessian 183 ~~~~~~~~~~~~~~~~~~~~~~ 184 _n_ 185 d2chi2 \ 1 / dRi() dRi() d2Ri() \ 186 --------------- = 2 > ---------- | ------- . ------- - (Ri - Ri()) . --------------- | 187 dthetaj.dthetak /__ sigma_i**2 \ dthetaj dthetak dthetaj.dthetak / 188 i=1 189 190 where: 191 Ri are the values of the measured relaxation data set. 192 Ri() are the values of the back calculated relaxation data set. 193 sigma_i are the values of the error set. 194 """ 195 196 # Arguments. 197 self.param_set = param_set 198 self.total_num_params = len(init_params) 199 self.scaling_matrix = scaling_matrix 200 self.num_res = num_res 201 self.params = 1.0 * init_params 202 203 # Data structures for tests set to some random array (in this case all pi). 204 self.func_test = pi * ones(self.total_num_params, Float64) 205 self.grad_test = pi * ones(self.total_num_params, Float64) 206 self.hess_test = pi * ones(self.total_num_params, Float64) 207 208 # Initialise the data class for storing diffusion tensor data. 209 self.diff_data = Data() 210 self.diff_data.type = diff_type 211 self.diff_data.params = diff_params 212 self.init_diff_data(self.diff_data) 213 214 # Total number of ri. 215 self.total_num_ri = 0 216 217 # Set the function for packaging diffusion tensor parameters. 218 if self.diff_data.type == 'sphere': 219 self.param_index = 1 220 self.diff_end_index = 1 221 elif self.diff_data.type == 'spheroid': 222 self.param_index = 4 223 self.diff_end_index = 4 224 elif self.diff_data.type == 'ellipsoid': 225 self.param_index = 6 226 self.diff_end_index = 6 227 if self.param_set != 'all': 228 self.param_index = 0 229 230 # Create the data array used to store data. 231 self.data = [] 232 for i in xrange(self.num_res): 233 # Total number of ri. 234 self.total_num_ri = self.total_num_ri + num_ri[i] 235 236 # Append the class instance Data to the data array. 237 self.data.append(Data()) 238 239 # Number of indecies. 240 self.data[i].num_indecies = self.diff_data.num_indecies 241 242 # Calculate the five frequencies per field strength which cause R1, R2, and NOE relaxation. 243 self.data[i].frq_list = zeros((num_frq[i], 5), Float64) 244 self.data[i].frq_list_ext = zeros((num_frq[i], 5, self.diff_data.num_indecies), Float64) 245 self.data[i].frq_sqrd_list_ext = zeros((num_frq[i], 5, self.diff_data.num_indecies), Float64) 246 for j in xrange(num_frq[i]): 247 frqH = 2.0 * pi * frq[i][j] 248 frqX = frqH / g_ratio 249 self.data[i].frq_list[j, 1] = frqX 250 self.data[i].frq_list[j, 2] = frqH - frqX 251 self.data[i].frq_list[j, 3] = frqH 252 self.data[i].frq_list[j, 4] = frqH + frqX 253 self.data[i].frq_sqrd_list = self.data[i].frq_list ** 2 254 for j in xrange(self.diff_data.num_indecies): 255 self.data[i].frq_list_ext[:, :, j] = self.data[i].frq_list 256 self.data[i].frq_sqrd_list_ext[:, :, j] = self.data[i].frq_sqrd_list 257 258 # Store supplied data in self.data 259 self.data[i].gh = gh 260 self.data[i].gx = gx 261 self.data[i].g_ratio = g_ratio 262 self.data[i].h_bar = h_bar 263 self.data[i].mu0 = mu0 264 self.data[i].equations = equations[i] 265 self.data[i].param_types = param_types[i] 266 self.data[i].relax_data = relax_data[i] 267 self.data[i].errors = errors[i] 268 self.data[i].bond_length = bond_length[i] 269 self.data[i].csa = csa[i] 270 self.data[i].num_ri = num_ri[i] 271 self.data[i].num_frq = num_frq[i] 272 self.data[i].frq = frq[i] 273 self.data[i].remap_table = remap_table[i] 274 self.data[i].noe_r1_table = noe_r1_table[i] 275 self.data[i].ri_labels = ri_labels[i] 276 self.data[i].num_params = num_params[i] 277 self.data[i].xh_unit_vector = vectors[i] 278 279 # Parameter values for minimising soley the diffusion tensor parameters. 280 if self.param_set == 'diff': 281 self.data[i].param_values = param_values[i] 282 283 # Indecies for constructing the global generic model-free gradient and Hessian kite. 284 if i == 0: 285 self.data[i].start_index = self.diff_data.num_params 286 else: 287 self.data[i].start_index = self.data[i-1].end_index 288 self.data[i].end_index = self.data[i].start_index + self.data[i].num_params 289 290 # Total number of parameters. 291 if self.param_set == 'mf' or self.param_set == 'local_tm': 292 self.data[i].total_num_params = self.data[i].num_params 293 elif self.param_set == 'diff': 294 self.data[i].total_num_params = self.diff_data.num_params 295 else: 296 self.data[i].total_num_params = self.data[i].num_params + self.diff_data.num_params 297 298 # Initialise the residue specific data. 299 self.init_res_data(self.data[i], self.diff_data) 300 301 # Setup the residue specific equations. 302 if not self.setup_equations(self.data[i]): 303 raise RelaxError, "The model-free equations could not be setup." 304 305 # Diffusion tensor parameters. 306 if self.param_set == 'local_tm': 307 self.diff_data.params = self.params[0:1] 308 elif self.param_set == 'diff' or self.param_set == 'all': 309 self.diff_data.params = self.params[0:self.diff_end_index] 310 311 # Calculate the correlation times ti. 312 self.diff_data.calc_ti(self.data[i], self.diff_data) 313 314 # ti spectral density components. 315 self.data[i].w_ti_sqrd = self.data[i].frq_sqrd_list_ext * self.data[i].ti ** 2 316 self.data[i].fact_ti = 1.0 / (1.0 + self.data[i].w_ti_sqrd) 317 318 # Initialise the R1 data class. This is used only if an NOE data set is collected but the R1 data of the same frequency has not. 319 missing_r1 = 0 320 for j in xrange(self.data[i].num_ri): 321 if self.data[i].ri_labels[j] == 'NOE' and self.data[i].noe_r1_table[j] == None: 322 missing_r1 = 1 323 if missing_r1: 324 self.init_res_r1_data(self.data[i]) 325 326 # Scaling initialisation. 327 if self.scaling_matrix: 328 self.scaling_flag = 1 329 else: 330 self.scaling_flag = 0 331 332 # Initialise the total chi2 value, gradient, and Hessian data structures. 333 self.total_chi2 = 0.0 334 self.total_dchi2 = zeros((self.total_num_params), Float64) 335 self.total_d2chi2 = zeros((self.total_num_params, self.total_num_params), Float64) 336 337 # Initialise the total ri gradient data structure (for Levenberg-Marquardt minimisation). 338 self.total_dri = zeros((self.total_num_params, self.total_num_ri), Float64) 339 340 # Set the functions self.func, self.dfunc, and self.d2func. 341 ########################################################### 342 343 # Functions for minimising model-free parameters for a single residue. 344 if self.param_set == 'mf': 345 self.func = self.func_mf 346 self.dfunc = self.dfunc_mf 347 self.d2func = self.d2func_mf 348 349 # Functions for minimising model-free parameters for a single residue with a local tm. 350 elif self.param_set == 'local_tm': 351 self.func = self.func_local_tm 352 self.dfunc = self.dfunc_local_tm 353 self.d2func = self.d2func_local_tm 354 355 # Functions for minimising diffusion tensor parameters with all model-free parameters fixed. 356 elif self.param_set == 'diff': 357 self.func = self.func_diff 358 self.dfunc = self.dfunc_diff 359 self.d2func = self.d2func_diff 360 361 # Functions for minimising diffusion tensor parameters together with all model-free parameters. 362 elif self.param_set == 'all': 363 self.func = self.func_all 364 self.dfunc = self.dfunc_all 365 self.d2func = self.d2func_all
366 367
368 - def func_mf(self, params):
369 """Function for calculating the chi-squared value. 370 371 Used in the minimisation of model-free parameters for a single residue. 372 """ 373 374 # Store the parameter values in self.func_test for testing. 375 self.func_test = params * 1.0 376 377 # Set self.data[0] to data. 378 data = self.data[0] 379 380 # Scaling. 381 if self.scaling_flag: 382 params = matrixmultiply(params, self.scaling_matrix) 383 384 # Direction cosine calculations. 385 if self.diff_data.calc_di: 386 self.diff_data.calc_di(data, self.diff_data) 387 388 # Diffusion tensor weight calculations. 389 self.diff_data.calc_ci(data, self.diff_data) 390 391 # Diffusion tensor correlation times. 392 self.diff_data.calc_ti(data, self.diff_data) 393 394 # Calculate the components of the spectral densities. 395 if data.calc_jw_comps: 396 data.calc_jw_comps(data, params) 397 398 # Calculate the spectral density values. 399 data.jw = data.calc_jw(data, params) 400 401 # Calculate the relaxation formula components. 402 data.create_ri_comps(data, params) 403 404 # Calculate the R1, R2, and sigma_noe values. 405 data.ri_prime = data.create_ri_prime(data) 406 407 # Calculate the NOE values. 408 data.ri = data.ri_prime * 1.0 409 for m in xrange(data.num_ri): 410 if data.create_ri[m]: 411 data.create_ri[m](data, m, data.remap_table[m], data.get_r1, params) 412 413 # Calculate the chi-squared value. 414 data.chi2 = chi2(data.relax_data, data.ri, data.errors) 415 416 return data.chi2
417 418
419 - def func_local_tm(self, params):
420 """Function for calculating the chi-squared value. 421 422 Used in the minimisation of model-free parameters for a single residue with a local tm. 423 """ 424 425 # Store the parameter values in self.func_test for testing. 426 self.func_test = params * 1.0 427 428 # Set self.data[0] to data. 429 data = self.data[0] 430 431 # Scaling. 432 if self.scaling_flag: 433 params = matrixmultiply(params, self.scaling_matrix) 434 435 # Diffusion tensor parameters. 436 self.diff_data.params = params[0:1] 437 438 # Diffusion tensor weight calculations. 439 self.diff_data.calc_ci(data, self.diff_data) 440 441 # Diffusion tensor correlation times. 442 self.diff_data.calc_ti(data, self.diff_data) 443 444 # ti spectral density components. 445 data.w_ti_sqrd = data.frq_sqrd_list_ext * data.ti ** 2 446 data.fact_ti = 1.0 / (1.0 + data.w_ti_sqrd) 447 448 # Calculate the components of the spectral densities. 449 if data.calc_jw_comps: 450 data.calc_jw_comps(data, params) 451 452 # Calculate the spectral density values. 453 data.jw = data.calc_jw(data, params) 454 455 # Calculate the relaxation formula components. 456 data.create_ri_comps(data, params) 457 458 # Calculate the R1, R2, and sigma_noe values. 459 data.ri_prime = data.create_ri_prime(data) 460 461 # Calculate the NOE values. 462 data.ri = data.ri_prime * 1.0 463 for m in xrange(data.num_ri): 464 if data.create_ri[m]: 465 data.create_ri[m](data, m, data.remap_table[m], data.get_r1, params) 466 467 # Calculate the chi-squared value. 468 data.chi2 = chi2(data.relax_data, data.ri, data.errors) 469 470 return data.chi2
471 472
473 - def func_diff(self, params):
474 """Function for calculating the chi-squared value. 475 476 Used in the minimisation of diffusion tensor parameters with all model-free parameters 477 fixed. 478 """ 479 480 # Store the parameter values in self.func_test for testing. 481 self.func_test = params * 1.0 482 483 # Scaling. 484 if self.scaling_flag: 485 params = matrixmultiply(params, self.scaling_matrix) 486 487 # Diffusion tensor parameters. 488 self.diff_data.params = params[0:self.diff_end_index] 489 490 # Set the total chi2 to zero. 491 self.total_chi2 = 0.0 492 493 # Loop over the residues. 494 for i in xrange(self.num_res): 495 # Set self.data[i] to data. 496 data = self.data[i] 497 498 # Direction cosine calculations. 499 if self.diff_data.calc_di: 500 self.diff_data.calc_di(data, self.diff_data) 501 502 # Diffusion tensor weight calculations. 503 self.diff_data.calc_ci(data, self.diff_data) 504 505 # Diffusion tensor correlation times. 506 self.diff_data.calc_ti(data, self.diff_data) 507 508 # ti spectral density components. 509 data.w_ti_sqrd = data.frq_sqrd_list_ext * data.ti ** 2 510 data.fact_ti = 1.0 / (1.0 + data.w_ti_sqrd) 511 512 # Calculate the components of the spectral densities. 513 if data.calc_jw_comps: 514 data.calc_jw_comps(data, data.param_values) 515 516 # Calculate the spectral density values. 517 data.jw = data.calc_jw(data, data.param_values) 518 519 # Calculate the relaxation formula components. 520 data.create_ri_comps(data, data.param_values) 521 522 # Calculate the R1, R2, and sigma_noe values. 523 data.ri_prime = data.create_ri_prime(data) 524 525 # Calculate the NOE values. 526 data.ri = data.ri_prime * 1.0 527 for m in xrange(data.num_ri): 528 if data.create_ri[m]: 529 data.create_ri[m](data, m, data.remap_table[m], data.get_r1, data.param_values) 530 531 # Calculate the chi-squared value. 532 data.chi2 = chi2(data.relax_data, data.ri, data.errors) 533 534 # Add the residue specific chi2 to the total chi2. 535 self.total_chi2 = self.total_chi2 + data.chi2 536 537 return self.total_chi2
538 539
540 - def func_all(self, params):
541 """Function for calculating the chi-squared value. 542 543 Used in the minimisation of diffusion tensor parameters together with all model-free 544 parameters. 545 """ 546 547 # Store the parameter values in self.func_test for testing. 548 self.func_test = params * 1.0 549 550 # Scaling. 551 if self.scaling_flag: 552 params = matrixmultiply(params, self.scaling_matrix) 553 554 # Diffusion tensor parameters. 555 self.diff_data.params = params[0:self.diff_end_index] 556 557 # Set the total chi2 to zero. 558 self.total_chi2 = 0.0 559 560 # Loop over the residues. 561 for i in xrange(self.num_res): 562 # Set self.data[i] to data. 563 data = self.data[i] 564 565 # Direction cosine calculations. 566 if self.diff_data.calc_di: 567 self.diff_data.calc_di(data, self.diff_data) 568 569 # Diffusion tensor weight calculations. 570 self.diff_data.calc_ci(data, self.diff_data) 571 572 # Diffusion tensor correlation times. 573 self.diff_data.calc_ti(data, self.diff_data) 574 575 # ti spectral density components. 576 data.w_ti_sqrd = data.frq_sqrd_list_ext * data.ti ** 2 577 data.fact_ti = 1.0 / (1.0 + data.w_ti_sqrd) 578 579 # Calculate the components of the spectral densities. 580 if data.calc_jw_comps: 581 data.calc_jw_comps(data, params) 582 583 # Calculate the spectral density values. 584 data.jw = data.calc_jw(data, params) 585 586 # Calculate the relaxation formula components. 587 data.create_ri_comps(data, params) 588 589 # Calculate the R1, R2, and sigma_noe values. 590 data.ri_prime = data.create_ri_prime(data) 591 592 # Calculate the NOE values. 593 data.ri = data.ri_prime * 1.0 594 for m in xrange(data.num_ri): 595 if data.create_ri[m]: 596 data.create_ri[m](data, m, data.remap_table[m], data.get_r1, params) 597 598 # Calculate the chi-squared value. 599 data.chi2 = chi2(data.relax_data, data.ri, data.errors) 600 601 # Add the residue specific chi2 to the total chi2. 602 self.total_chi2 = self.total_chi2 + data.chi2 603 604 return self.total_chi2
605 606
607 - def dfunc_mf(self, params):
608 """Function for calculating the chi-squared gradient. 609 610 Used in the minimisation of model-free parameters for a single residue. 611 """ 612 613 # Test if the function has already been called, otherwise run self.func. 614 if sum(params == self.func_test) != self.total_num_params: 615 self.func(params) 616 617 # Store the parameter values in self.grad_test for testing. 618 self.grad_test = params * 1.0 619 620 # Set self.data[0] to data. 621 data = self.data[0] 622 623 # Scaling. 624 if self.scaling_flag: 625 params = matrixmultiply(params, self.scaling_matrix) 626 627 # Calculate the spectral density gradient components. 628 if data.calc_djw_comps: 629 data.calc_djw_comps(data, params) 630 631 # Loop over the gradient. 632 for j in xrange(data.total_num_params): 633 # Calculate the spectral density gradients. 634 if data.calc_djw[j]: 635 data.djw = data.calc_djw[j](data, params, j) 636 else: 637 data.djw = data.djw * 0.0 638 639 # Calculate the relaxation gradient components. 640 data.create_dri_comps(data, params) 641 642 # Calculate the R1, R2, and sigma_noe gradients. 643 data.dri_prime[j] = data.create_dri_prime[j](data) 644 645 # Loop over the relaxation values and modify the NOE gradients. 646 data.dri[j] = data.dri_prime[j] 647 for m in xrange(data.num_ri): 648 if data.create_dri[m]: 649 data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j) 650 651 # Calculate the chi-squared gradient. 652 data.dchi2[j] = dchi2(data.relax_data, data.ri, data.dri[j], data.errors) 653 654 # Diagonal scaling. 655 if self.scaling_flag: 656 data.dchi2 = matrixmultiply(data.dchi2, self.scaling_matrix) 657 658 # Return a copy of the gradient. 659 return data.dchi2 * 1.0
660 661
662 - def dfunc_local_tm(self, params):
663 """Function for calculating the chi-squared gradient. 664 665 Used in the minimisation of model-free parameters for a single residue with a local tm. 666 """ 667 668 # Test if the function has already been called, otherwise run self.func. 669 if sum(params == self.func_test) != self.total_num_params: 670 self.func(params) 671 672 # Store the parameter values in self.grad_test for testing. 673 self.grad_test = params * 1.0 674 675 # Set self.data[0] to data. 676 data = self.data[0] 677 678 # Scaling. 679 if self.scaling_flag: 680 params = matrixmultiply(params, self.scaling_matrix) 681 682 # Diffusion tensor parameters. 683 self.diff_data.params = params[0:1] 684 685 # Calculate the spectral density gradient components. 686 if data.calc_djw_comps: 687 data.calc_djw_comps(data, params) 688 689 # Diffusion tensor correlation times. 690 self.diff_data.calc_dti(data, self.diff_data) 691 692 # Loop over the gradient. 693 for j in xrange(data.total_num_params): 694 # Calculate the spectral density gradients. 695 if data.calc_djw[j]: 696 data.djw = data.calc_djw[j](data, params, j) 697 else: 698 data.djw = data.djw * 0.0 699 700 # Calculate the relaxation gradient components. 701 data.create_dri_comps(data, params) 702 703 # Calculate the R1, R2, and sigma_noe gradients. 704 data.dri_prime[j] = data.create_dri_prime[j](data) 705 706 # Loop over the relaxation values and modify the NOE gradients. 707 data.dri[j] = data.dri_prime[j] 708 for m in xrange(data.num_ri): 709 if data.create_dri[m]: 710 data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j) 711 712 # Calculate the chi-squared gradient. 713 data.dchi2[j] = dchi2(data.relax_data, data.ri, data.dri[j], data.errors) 714 715 # Diagonal scaling. 716 if self.scaling_flag: 717 data.dchi2 = matrixmultiply(data.dchi2, self.scaling_matrix) 718 719 # Return a copy of the gradient. 720 return data.dchi2 * 1.0
721 722
723 - def dfunc_diff(self, params):
724 """Function for calculating the chi-squared gradient. 725 726 Used in the minimisation of diffusion tensor parameters with all model-free parameters 727 fixed. 728 """ 729 730 # Test if the function has already been called, otherwise run self.func. 731 if sum(params == self.func_test) != self.total_num_params: 732 self.func(params) 733 734 # Store the parameter values in self.grad_test for testing. 735 self.grad_test = params * 1.0 736 737 # Scaling. 738 if self.scaling_flag: 739 params = matrixmultiply(params, self.scaling_matrix) 740 741 # Diffusion tensor parameters. 742 self.diff_data.params = params[0:self.diff_end_index] 743 744 # Set the total chi2 gradient to zero. 745 self.total_dchi2 = self.total_dchi2 * 0.0 746 747 # Loop over the residues. 748 for i in xrange(self.num_res): 749 750 # Set self.data[i] to data. 751 data = self.data[i] 752 753 # Direction cosine calculations. 754 if self.diff_data.calc_ddi: 755 self.diff_data.calc_ddi(data, self.diff_data) 756 757 # Diffusion tensor weight calculations. 758 if self.diff_data.calc_dci: 759 self.diff_data.calc_dci(data, self.diff_data) 760 761 # Diffusion tensor correlation times. 762 self.diff_data.calc_dti(data, self.diff_data) 763 764 # Calculate the spectral density gradient components. 765 if data.calc_djw_comps: 766 data.calc_djw_comps(data, data.param_values) 767 768 # Loop over the gradient. 769 for j in xrange(data.total_num_params): 770 # Calculate the spectral density gradients. 771 if data.calc_djw[j]: 772 data.djw = data.calc_djw[j](data, data.param_values, j) 773 else: 774 data.djw = data.djw * 0.0 775 776 # Calculate the relaxation gradient components. 777 data.create_dri_comps(data, data.param_values) 778 779 # Calculate the R1, R2, and sigma_noe gradients. 780 data.dri_prime[j] = data.create_dri_prime[j](data) 781 782 # Loop over the relaxation values and modify the NOE gradients. 783 data.dri[j] = data.dri_prime[j] 784 for m in xrange(data.num_ri): 785 if data.create_dri[m]: 786 data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j) 787 788 # Calculate the chi-squared gradient. 789 data.dchi2[j] = dchi2(data.relax_data, data.ri, data.dri[j], data.errors) 790 791 # Index for the construction of the global generic model-free gradient. 792 index = self.diff_data.num_params 793 794 # Diffusion parameter part of the global generic model-free gradient. 795 self.total_dchi2[0:index] = self.total_dchi2[0:index] + data.dchi2[0:index] 796 797 # Diagonal scaling. 798 if self.scaling_flag: 799 self.total_dchi2 = matrixmultiply(self.total_dchi2, self.scaling_matrix) 800 801 # Return a copy of the gradient. 802 return self.total_dchi2 * 1.0
803 804
805 - def dfunc_all(self, params):
806 """Function for calculating the chi-squared gradient. 807 808 Used in the minimisation of diffusion tensor parameters together with all model-free 809 parameters. 810 """ 811 812 # Test if the function has already been called, otherwise run self.func. 813 if sum(params == self.func_test) != self.total_num_params: 814 self.func(params) 815 816 # Store the parameter values in self.grad_test for testing. 817 self.grad_test = params * 1.0 818 819 # Scaling. 820 if self.scaling_flag: 821 params = matrixmultiply(params, self.scaling_matrix) 822 823 # Diffusion tensor parameters. 824 self.diff_data.params = params[0:self.diff_end_index] 825 826 # Set the total chi2 gradient to zero. 827 self.total_dchi2 = self.total_dchi2 * 0.0 828 829 # Loop over the residues. 830 for i in xrange(self.num_res): 831 # Set self.data[i] to data. 832 data = self.data[i] 833 834 # Direction cosine calculations. 835 if self.diff_data.calc_ddi: 836 self.diff_data.calc_ddi(data, self.diff_data) 837 838 # Diffusion tensor weight calculations. 839 if self.diff_data.calc_dci: 840 self.diff_data.calc_dci(data, self.diff_data) 841 842 # Diffusion tensor correlation times. 843 self.diff_data.calc_dti(data, self.diff_data) 844 845 # Calculate the spectral density gradient components. 846 if data.calc_djw_comps: 847 data.calc_djw_comps(data, params) 848 849 # Loop over the gradient. 850 for j in xrange(data.total_num_params): 851 # Calculate the spectral density gradients. 852 if data.calc_djw[j]: 853 data.djw = data.calc_djw[j](data, params, j) 854 else: 855 data.djw = data.djw * 0.0 856 857 # Calculate the relaxation gradient components. 858 data.create_dri_comps(data, params) 859 860 # Calculate the R1, R2, and sigma_noe gradients. 861 data.dri_prime[j] = data.create_dri_prime[j](data) 862 863 # Loop over the relaxation values and modify the NOE gradients. 864 data.dri[j] = data.dri_prime[j] 865 for m in xrange(data.num_ri): 866 if data.create_dri[m]: 867 data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j) 868 869 # Calculate the chi-squared gradient. 870 data.dchi2[j] = dchi2(data.relax_data, data.ri, data.dri[j], data.errors) 871 872 # Index for the construction of the global generic model-free gradient. 873 index = self.diff_data.num_params 874 875 # Diffusion parameter part of the global generic model-free gradient. 876 self.total_dchi2[0:index] = self.total_dchi2[0:index] + data.dchi2[0:index] 877 878 # Model-free parameter part of the global generic model-free gradient. 879 self.total_dchi2[data.start_index:data.end_index] = self.total_dchi2[data.start_index:data.end_index] + data.dchi2[index:] 880 881 # Diagonal scaling. 882 if self.scaling_flag: 883 self.total_dchi2 = matrixmultiply(self.total_dchi2, self.scaling_matrix) 884 885 # Return a copy of the gradient. 886 return self.total_dchi2 * 1.0
887 888
889 - def d2func_mf(self, params):
890 """Function for calculating the chi-squared Hessian. 891 892 Used in the minimisation of model-free parameters for a single residue. 893 """ 894 895 # Test if the gradient has already been called, otherwise run self.dfunc. 896 if sum(params == self.grad_test) != self.total_num_params: 897 self.dfunc(params) 898 899 # Set self.data[0] to data. 900 data = self.data[0] 901 902 # Scaling. 903 if self.scaling_flag: 904 params = matrixmultiply(params, self.scaling_matrix) 905 906 # Loop over the lower triangle of the Hessian. 907 for j in xrange(data.total_num_params): 908 for k in xrange(j + 1): 909 # Calculate the spectral density Hessians. 910 if data.calc_d2jw[j][k]: 911 data.d2jw = data.calc_d2jw[j][k](data, params, j, k) 912 else: 913 data.d2jw = data.d2jw * 0.0 914 915 # Calculate the relaxation Hessian components. 916 data.create_d2ri_comps(data, params) 917 918 # Calculate the R1, R2, and sigma_noe Hessians. 919 if data.create_d2ri_prime[j][k]: 920 data.d2ri_prime[j, k] = data.create_d2ri_prime[j][k](data) 921 922 # Loop over the relaxation values and modify the NOE Hessians. 923 data.d2ri[j, k] = data.d2ri_prime[j, k] 924 for m in xrange(data.num_ri): 925 if data.create_d2ri[m]: 926 data.create_d2ri[m](data, m, data.remap_table[m], data.get_d2r1, params, j, k) 927 928 # Calculate the chi-squared Hessian. 929 data.d2chi2[j, k] = data.d2chi2[k, j] = d2chi2(data.relax_data, data.ri, data.dri[j], data.dri[k], data.d2ri[j, k], data.errors) 930 931 # Diagonal scaling. 932 if self.scaling_flag: 933 data.d2chi2 = matrixmultiply(self.scaling_matrix, matrixmultiply(data.d2chi2, self.scaling_matrix)) 934 935 # Return a copy of the Hessian. 936 return data.d2chi2 * 1.0
937 938
939 - def d2func_local_tm(self, params):
940 """Function for calculating the chi-squared Hessian. 941 942 Used in the minimisation of model-free parameters for a single residue with a local tm. 943 """ 944 945 # Test if the gradient has already been called, otherwise run self.dfunc. 946 if sum(params == self.grad_test) != self.total_num_params: 947 self.dfunc(params) 948 949 # Set self.data[0] to data. 950 data = self.data[0] 951 952 # Scaling. 953 if self.scaling_flag: 954 params = matrixmultiply(params, self.scaling_matrix) 955 956 # Diffusion tensor parameters. 957 self.diff_data.params = params[0:1] 958 959 # Loop over the lower triangle of the Hessian. 960 for j in xrange(data.total_num_params): 961 for k in xrange(j + 1): 962 # Calculate the spectral density Hessians. 963 if data.calc_d2jw[j][k]: 964 data.d2jw = data.calc_d2jw[j][k](data, params, j, k) 965 else: 966 data.d2jw = data.d2jw * 0.0 967 968 # Calculate the relaxation Hessian components. 969 data.create_d2ri_comps(data, params) 970 971 # Calculate the R1, R2, and sigma_noe Hessians. 972 if data.create_d2ri_prime[j][k]: 973 data.d2ri_prime[j, k] = data.create_d2ri_prime[j][k](data) 974 975 # Loop over the relaxation values and modify the NOE Hessians. 976 data.d2ri[j, k] = data.d2ri_prime[j, k] 977 for m in xrange(data.num_ri): 978 if data.create_d2ri[m]: 979 data.create_d2ri[m](data, m, data.remap_table[m], data.get_d2r1, params, j, k) 980 981 # Calculate the chi-squared Hessian. 982 data.d2chi2[j, k] = data.d2chi2[k, j] = d2chi2(data.relax_data, data.ri, data.dri[j], data.dri[k], data.d2ri[j, k], data.errors) 983 984 # Diagonal scaling. 985 if self.scaling_flag: 986 data.d2chi2 = matrixmultiply(self.scaling_matrix, matrixmultiply(data.d2chi2, self.scaling_matrix)) 987 988 # Return a copy of the Hessian. 989 return data.d2chi2 * 1.0
990 991
992 - def d2func_diff(self, params):
993 """Function for calculating the chi-squared Hessian. 994 995 Used in the minimisation of diffusion tensor parameters with all model-free parameters 996 fixed. 997 """ 998 999 # Test if the gradient has already been called, otherwise run self.dfunc. 1000 if sum(params == self.grad_test) != self.total_num_params: 1001 self.dfunc(params) 1002 1003 # Scaling. 1004 if self.scaling_flag: 1005 params = matrixmultiply(params, self.scaling_matrix) 1006 1007 # Diffusion tensor parameters. 1008 self.diff_data.params = params[0:self.diff_end_index] 1009 1010 # Set the total chi2 Hessian to zero. 1011 self.total_d2chi2 = self.total_d2chi2 * 0.0 1012 1013 # Loop over the residues. 1014 for i in xrange(self.num_res): 1015 # Set self.data[i] to data. 1016 data = self.data[i] 1017 1018 # Direction cosine calculations. 1019 if self.diff_data.calc_d2di: 1020 self.diff_data.calc_d2di(data, self.diff_data) 1021 1022 # Diffusion tensor weight calculations. 1023 if self.diff_data.calc_d2ci: 1024 self.diff_data.calc_d2ci(data, self.diff_data) 1025 1026 # Diffusion tensor correlation times. 1027 if self.diff_data.calc_d2ti: 1028 self.diff_data.calc_d2ti(data, self.diff_data) 1029 1030 # Loop over the lower triangle of the Hessian. 1031 for j in xrange(data.total_num_params): 1032 for k in xrange(j + 1): 1033 # Calculate the spectral density Hessians. 1034 if data.calc_d2jw[j][k]: 1035 data.d2jw = data.calc_d2jw[j][k](data, data.param_values, j, k) 1036 else: 1037 data.d2jw = data.d2jw * 0.0 1038 1039 # Calculate the relaxation Hessian components. 1040 data.create_d2ri_comps(data, data.param_values) 1041 1042 # Calculate the R1, R2, and sigma_noe Hessians. 1043 if data.create_d2ri_prime[j][k]: 1044 data.d2ri_prime[j, k] = data.create_d2ri_prime[j][k](data) 1045 1046 # Loop over the relaxation values and modify the NOE Hessians. 1047 data.d2ri[j, k] = data.d2ri_prime[j, k] 1048 for m in xrange(data.num_ri): 1049 if data.create_d2ri[m]: 1050 data.create_d2ri[m](data, m, data.remap_table[m], data.get_d2r1, params, j, k) 1051 1052 # Calculate the chi-squared Hessian. 1053 data.d2chi2[j, k] = data.d2chi2[k, j] = d2chi2(data.relax_data, data.ri, data.dri[j], data.dri[k], data.d2ri[j, k], data.errors) 1054 1055 # Pure diffusion parameter part of the global generic model-free Hessian. 1056 self.total_d2chi2 = self.total_d2chi2 + data.d2chi2 1057 1058 # Diagonal scaling. 1059 if self.scaling_flag: 1060 self.total_d2chi2 = matrixmultiply(self.scaling_matrix, matrixmultiply(self.total_d2chi2, self.scaling_matrix)) 1061 1062 # Return a copy of the Hessian. 1063 return self.total_d2chi2 * 1.0
1064 1065
1066 - def d2func_all(self, params):
1067 """Function for calculating the chi-squared Hessian. 1068 1069 Used in the minimisation of diffusion tensor parameters together with all model-free 1070 parameters. 1071 """ 1072 1073 # Test if the gradient has already been called, otherwise run self.dfunc. 1074 if sum(params == self.grad_test) != self.total_num_params: 1075 self.dfunc(params) 1076 1077 # Scaling. 1078 if self.scaling_flag: 1079 params = matrixmultiply(params, self.scaling_matrix) 1080 1081 # Diffusion tensor parameters. 1082 self.diff_data.params = params[0:self.diff_end_index] 1083 1084 # Set the total chi2 Hessian to zero. 1085 self.total_d2chi2 = self.total_d2chi2 * 0.0 1086 1087 # Loop over the residues. 1088 for i in xrange(self.num_res): 1089 # Set self.data[i] to data. 1090 data = self.data[i] 1091 1092 # Direction cosine calculations. 1093 if self.diff_data.calc_d2di: 1094 self.diff_data.calc_d2di(data, self.diff_data) 1095 1096 # Diffusion tensor weight calculations. 1097 if self.diff_data.calc_d2ci: 1098 self.diff_data.calc_d2ci(data, self.diff_data) 1099 1100 # Diffusion tensor correlation times. 1101 if self.diff_data.calc_d2ti: 1102 self.diff_data.calc_d2ti(data, self.diff_data) 1103 1104 # Loop over the lower triangle of the Hessian. 1105 for j in xrange(data.total_num_params): 1106 for k in xrange(j + 1): 1107 # Calculate the spectral density Hessians. 1108 if data.calc_d2jw[j][k]: 1109 data.d2jw = data.calc_d2jw[j][k](data, params, j, k) 1110 else: 1111 data.d2jw = data.d2jw * 0.0 1112 1113 # Calculate the relaxation Hessian components. 1114 data.create_d2ri_comps(data, params) 1115 1116 # Calculate the R1, R2, and sigma_noe Hessians. 1117 if data.create_d2ri_prime[j][k]: 1118 data.d2ri_prime[j, k] = data.create_d2ri_prime[j][k](data) 1119 1120 # Loop over the relaxation values and modify the NOE Hessians. 1121 data.d2ri[j, k] = data.d2ri_prime[j, k] 1122 for m in xrange(data.num_ri): 1123 if data.create_d2ri[m]: 1124 data.create_d2ri[m](data, m, data.remap_table[m], data.get_d2r1, params, j, k) 1125 1126 # Calculate the chi-squared Hessian. 1127 data.d2chi2[j, k] = data.d2chi2[k, j] = d2chi2(data.relax_data, data.ri, data.dri[j], data.dri[k], data.d2ri[j, k], data.errors) 1128 1129 # Index for the construction of the global generic model-free Hessian. 1130 index = self.diff_data.num_params 1131 1132 # Pure diffusion parameter part of the global generic model-free Hessian. 1133 self.total_d2chi2[0:index, 0:index] = self.total_d2chi2[0:index, 0:index] + data.d2chi2[0:index, 0:index] 1134 1135 # Pure model-free parameter part of the global generic model-free Hessian. 1136 self.total_d2chi2[data.start_index:data.end_index, data.start_index:data.end_index] = self.total_d2chi2[data.start_index:data.end_index, data.start_index:data.end_index] + data.d2chi2[index:, index:] 1137 1138 # Off diagonal diffusion and model-free parameter parts of the global generic model-free Hessian. 1139 self.total_d2chi2[0:index, data.start_index:data.end_index] = self.total_d2chi2[0:index, data.start_index:data.end_index] + data.d2chi2[0:index, index:] 1140 self.total_d2chi2[data.start_index:data.end_index, 0:index] = self.total_d2chi2[data.start_index:data.end_index, 0:index] + data.d2chi2[index:, 0:index] 1141 1142 # Diagonal scaling. 1143 if self.scaling_flag: 1144 self.total_d2chi2 = matrixmultiply(self.scaling_matrix, matrixmultiply(self.total_d2chi2, self.scaling_matrix)) 1145 1146 # Return a copy of the Hessian. 1147 return self.total_d2chi2 * 1.0
1148 1149
1150 - def calc_ri(self):
1151 """Function for calculating relaxation values.""" 1152 1153 # Function call. 1154 chi2 = self.func_mf(self.params) 1155 1156 # Return the single value. 1157 return self.data[0].ri[0]
1158 1159
1160 - def init_diff_data(self, diff_data):
1161 """Function for the initialisation of diffusion tensor specific data.""" 1162 1163 # Diffusion as an sphere. 1164 if diff_data.type == 'sphere': 1165 # Number of diffusion parameters. 1166 diff_data.num_params = 1 1167 1168 # Number of indecies in the generic equations. 1169 diff_data.num_indecies = 1 1170 1171 # Direction cosine function, gradient, and Hessian. 1172 diff_data.calc_di = None 1173 diff_data.calc_ddi = None 1174 diff_data.calc_d2di = None 1175 1176 # Weight function, gradient, and Hessian. 1177 diff_data.calc_ci = calc_sphere_ci 1178 diff_data.calc_dci = None 1179 diff_data.calc_d2ci = None 1180 1181 # Global correlation time function, gradient, and Hessian. 1182 diff_data.calc_ti = calc_sphere_ti 1183 diff_data.calc_dti = calc_sphere_dti 1184 diff_data.calc_d2ti = None 1185 1186 1187 # Diffusion as an spheroid. 1188 elif diff_data.type == 'spheroid': 1189 # Number of diffusion parameters. 1190 diff_data.num_params = 4 1191 1192 # Number of indecies in the generic equations. 1193 diff_data.num_indecies = 3 1194 1195 # Direction cosine function, gradient, and Hessian. 1196 diff_data.calc_di = calc_spheroid_di 1197 diff_data.calc_ddi = calc_spheroid_ddi 1198 diff_data.calc_d2di = calc_spheroid_d2di 1199 1200 # Weight function, gradient, and Hessian. 1201 diff_data.calc_ci = calc_spheroid_ci 1202 diff_data.calc_dci = calc_spheroid_dci 1203 diff_data.calc_d2ci = calc_spheroid_d2ci 1204 1205 # Global correlation time function, gradient, and Hessian. 1206 diff_data.calc_ti = calc_spheroid_ti 1207 diff_data.calc_dti = calc_spheroid_dti 1208 diff_data.calc_d2ti = calc_spheroid_d2ti 1209 1210 # Unit vectors. 1211 diff_data.dpar = zeros(3, Float64) 1212 1213 # Unit vector gradients. 1214 diff_data.dpar_dtheta = zeros(3, Float64) 1215 diff_data.dpar_dphi = zeros(3, Float64) 1216 1217 # Unit vector Hessians. 1218 diff_data.dpar_dtheta2 = zeros(3, Float64) 1219 diff_data.dpar_dthetadphi = zeros(3, Float64) 1220 diff_data.dpar_dphi2 = zeros(3, Float64) 1221 1222 # Diffusion as an ellipsoid. 1223 elif diff_data.type == 'ellipsoid': 1224 # Number of diffusion parameters. 1225 diff_data.num_params = 6 1226 1227 # Number of indecies in the generic equations. 1228 diff_data.num_indecies = 5 1229 1230 # Direction cosine function, gradient, and Hessian. 1231 diff_data.calc_di = calc_ellipsoid_di 1232 diff_data.calc_ddi = calc_ellipsoid_ddi 1233 diff_data.calc_d2di = calc_ellipsoid_d2di 1234 1235 # Weight function, gradient, and Hessian. 1236 diff_data.calc_ci = calc_ellipsoid_ci 1237 diff_data.calc_dci = calc_ellipsoid_dci 1238 diff_data.calc_d2ci = calc_ellipsoid_d2ci 1239 1240 # Global correlation time function, gradient, and Hessian. 1241 diff_data.calc_ti = calc_ellipsoid_ti 1242 diff_data.calc_dti = calc_ellipsoid_dti 1243 diff_data.calc_d2ti = calc_ellipsoid_d2ti 1244 1245 # Unit vectors. 1246 diff_data.dx = zeros(3, Float64) 1247 diff_data.dy = zeros(3, Float64) 1248 diff_data.dz = zeros(3, Float64) 1249 1250 # Unit vector gradients. 1251 diff_data.ddx_dalpha = zeros(3, Float64) 1252 diff_data.ddx_dbeta = zeros(3, Float64) 1253 diff_data.ddx_dgamma = zeros(3, Float64) 1254 1255 diff_data.ddy_dalpha = zeros(3, Float64) 1256 diff_data.ddy_dbeta = zeros(3, Float64) 1257 diff_data.ddy_dgamma = zeros(3, Float64) 1258 1259 diff_data.ddz_dalpha = zeros(3, Float64) 1260 diff_data.ddz_dbeta = zeros(3, Float64) 1261 diff_data.ddz_dgamma = zeros(3, Float64) 1262 1263 # Unit vector Hessians. 1264 diff_data.d2dx_dalpha2 = zeros(3, Float64) 1265 diff_data.d2dx_dalpha_dbeta = zeros(3, Float64) 1266 diff_data.d2dx_dalpha_dgamma = zeros(3, Float64) 1267 diff_data.d2dx_dbeta2 = zeros(3, Float64) 1268 diff_data.d2dx_dbeta_dgamma = zeros(3, Float64) 1269 diff_data.d2dx_dgamma2 = zeros(3, Float64) 1270 1271 diff_data.d2dy_dalpha2 = zeros(3, Float64) 1272 diff_data.d2dy_dalpha_dbeta = zeros(3, Float64) 1273 diff_data.d2dy_dalpha_dgamma = zeros(3, Float64) 1274 diff_data.d2dy_dbeta2 = zeros(3, Float64) 1275 diff_data.d2dy_dbeta_dgamma = zeros(3, Float64) 1276 diff_data.d2dy_dgamma2 = zeros(3, Float64) 1277 1278 diff_data.d2dz_dalpha2 = zeros(3, Float64) 1279 diff_data.d2dz_dalpha_dbeta = zeros(3, Float64) 1280 diff_data.d2dz_dalpha_dgamma = zeros(3, Float64) 1281 diff_data.d2dz_dbeta2 = zeros(3, Float64) 1282 diff_data.d2dz_dbeta_dgamma = zeros(3, Float64) 1283 diff_data.d2dz_dgamma2 = zeros(3, Float64)
1284 1285
1286 - def init_res_data(self, data, diff_data):
1287 """Function for the initialisation of the residue specific data.""" 1288 1289 # Correlation times. 1290 data.ci = zeros(diff_data.num_indecies, Float64) 1291 data.ci_comps = zeros(diff_data.num_indecies, Float64) 1292 1293 # Weights. 1294 data.ti = zeros(diff_data.num_indecies, Float64) 1295 data.tau_comps = zeros(diff_data.num_indecies, Float64) 1296 data.tau_comps_sqrd = zeros(diff_data.num_indecies, Float64) 1297 data.tau_comps_cubed = zeros(diff_data.num_indecies, Float64) 1298 data.tau_scale = zeros(diff_data.num_indecies, Float64) 1299 1300 # Diffusion as a sphere. 1301 if self.diff_data.type == 'sphere': 1302 # Global correlation time gradient and Hessian. 1303 data.dti = zeros((1, diff_data.num_indecies), Float64) 1304 data.d2ti = zeros((1, 1, diff_data.num_indecies), Float64) 1305 1306 # Diffusion as a spheroid. 1307 elif self.diff_data.type == 'spheroid': 1308 # Weight gradient and Hessian. 1309 data.dci = zeros((4, diff_data.num_indecies), Float64) 1310 data.d2ci = zeros((4, 4, diff_data.num_indecies), Float64) 1311 1312 # Global correlation time gradient and Hessian. 1313 data.dti = zeros((2, diff_data.num_indecies), Float64) 1314 data.d2ti = zeros((2, 2, diff_data.num_indecies), Float64) 1315 1316 # Dot product. 1317 data.dz = 0 1318 1319 # Dot product gradient. 1320 data.ddz_dO = zeros(2, Float64) 1321 1322 # Dot product Hessian. 1323 data.d2dz_dO2 = zeros((2, 2), Float64) 1324 1325 # Diffusion as an ellipsoid. 1326 elif self.diff_data.type == 'ellipsoid': 1327 # Weight gradient and Hessian. 1328 data.dci = zeros((6, diff_data.num_indecies), Float64) 1329 data.d2ci = zeros((6, 6, diff_data.num_indecies), Float64) 1330 1331 # Global correlation time gradient and Hessian. 1332 data.dti = zeros((3, diff_data.num_indecies), Float64) 1333 data.d2ti = zeros((3, 3, diff_data.num_indecies), Float64) 1334 1335 # Dot products. 1336 data.dx = 0.0 1337 data.dy = 0.0 1338 data.dz = 0.0 1339 1340 # Dot product gradients. 1341 data.ddx_dO = zeros(3, Float64) 1342 data.ddy_dO = zeros(3, Float64) 1343 data.ddz_dO = zeros(3, Float64) 1344 1345 # Dot product Hessians. 1346 data.d2dx_dO2 = zeros((3, 3), Float64) 1347 data.d2dy_dO2 = zeros((3, 3), Float64) 1348 data.d2dz_dO2 = zeros((3, 3), Float64) 1349 1350 # Empty spectral density components. 1351 data.w_ti_sqrd = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1352 data.fact_ti = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1353 data.w_te_ti_sqrd = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1354 data.w_tf_ti_sqrd = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1355 data.w_ts_ti_sqrd = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1356 data.inv_te_denom = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1357 data.inv_tf_denom = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1358 data.inv_ts_denom = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1359 1360 # Empty spectral density values, gradients, and Hessians. 1361 data.jw = zeros((data.num_frq, 5), Float64) 1362 data.djw = zeros((data.num_frq, 5), Float64) 1363 data.d2jw = zeros((data.num_frq, 5), Float64) 1364 1365 # Calculate the fixed components of the dipolar and CSA constants. 1366 data.csa_const_fixed = zeros(data.num_frq, Float64) 1367 data.dip_const_fixed = None 1368 calc_fixed_csa(data) 1369 calc_fixed_dip(data) 1370 1371 # Dipolar and CSA constants. 1372 data.dip_const_func = 0.0 1373 data.dip_const_grad = 0.0 1374 data.dip_const_hess = 0.0 1375 data.csa_const_func = zeros(data.num_frq, Float64) 1376 data.csa_const_grad = zeros(data.num_frq, Float64) 1377 data.csa_const_hess = zeros(data.num_frq, Float64) 1378 1379 # Components of the transformed relaxation equations. 1380 data.dip_comps_func = zeros(data.num_ri, Float64) 1381 data.csa_comps_func = zeros(data.num_ri, Float64) 1382 data.rex_comps_func = zeros(data.num_ri, Float64) 1383 data.dip_jw_comps_func = zeros(data.num_ri, Float64) 1384 data.csa_jw_comps_func = zeros(data.num_ri, Float64) 1385 1386 # First partial derivative components of the transformed relaxation equations. 1387 data.dip_comps_grad = zeros(data.num_ri, Float64) 1388 data.csa_comps_grad = zeros(data.num_ri, Float64) 1389 data.rex_comps_grad = zeros(data.num_ri, Float64) 1390 data.dip_jw_comps_grad = zeros(data.num_ri, Float64) 1391 data.csa_jw_comps_grad = zeros(data.num_ri, Float64) 1392 1393 # First partial derivative components of the transformed relaxation equations. 1394 data.dip_comps_hess = zeros(data.num_ri, Float64) 1395 data.csa_comps_hess = zeros(data.num_ri, Float64) 1396 data.rex_comps_hess = zeros(data.num_ri, Float64) 1397 data.dip_jw_comps_hess = zeros(data.num_ri, Float64) 1398 data.csa_jw_comps_hess = zeros(data.num_ri, Float64) 1399 1400 # Transformed relaxation values, gradients, and Hessians. 1401 data.ri_prime = zeros((data.num_ri), Float64) 1402 data.dri_prime = zeros((data.total_num_params, data.num_ri), Float64) 1403 data.d2ri_prime = zeros((data.total_num_params, data.total_num_params, data.num_ri), Float64) 1404 1405 # Data structures containing the Ri values. 1406 data.ri = zeros(data.num_ri, Float64) 1407 data.dri = zeros((data.total_num_params, data.num_ri), Float64) 1408 data.d2ri = zeros((data.total_num_params, data.total_num_params, data.num_ri), Float64) 1409 1410 # Data structures containing the R1 values at the position of and corresponding to the NOE. 1411 data.r1 = zeros(data.num_ri, Float64) 1412 data.dr1 = zeros((data.total_num_params, data.num_ri), Float64) 1413 data.d2r1 = zeros((data.total_num_params, data.total_num_params, data.num_ri), Float64) 1414 1415 # Data structures containing the chi-squared values. 1416 data.chi2 = 0.0 1417 data.dchi2 = zeros((data.total_num_params), Float64) 1418 data.d2chi2 = zeros((data.total_num_params, data.total_num_params), Float64)
1419 1420
1421 - def init_res_r1_data(self, data):
1422 """Function for initialisation of the R1 data class. 1423 1424 This data class is only used if an NOE data set is collected but no R1 data set 1425 corresponding to the same frequency exists. 1426 """ 1427 1428 # Initialise an instance of Data. 1429 r1_data = Data() 1430 1431 # Copy a few things from 'data' to 'r1_data' 1432 r1_data.num_frq = data.num_frq 1433 r1_data.dip_const_fixed = data.dip_const_fixed 1434 r1_data.csa_const_fixed = data.csa_const_fixed 1435 1436 # Components of the transformed relaxation equations. 1437 r1_data.dip_comps_func = zeros(data.num_ri, Float64) 1438 r1_data.csa_comps_func = zeros(data.num_ri, Float64) 1439 r1_data.dip_jw_comps_func = zeros(data.num_ri, Float64) 1440 r1_data.csa_jw_comps_func = zeros(data.num_ri, Float64) 1441 1442 # Initialise the first partial derivative components of the transformed relaxation equations. 1443 r1_data.dip_comps_grad = zeros(data.num_ri, Float64) 1444 r1_data.csa_comps_grad = zeros(data.num_ri, Float64) 1445 r1_data.rex_comps_grad = zeros(data.num_ri, Float64) 1446 r1_data.dip_jw_comps_grad = zeros(data.num_ri, Float64) 1447 r1_data.csa_jw_comps_grad = zeros(data.num_ri, Float64) 1448 1449 # Initialise the first partial derivative components of the transformed relaxation equations. 1450 r1_data.dip_comps_hess = zeros(data.num_ri, Float64) 1451 r1_data.csa_comps_hess = zeros(data.num_ri, Float64) 1452 r1_data.rex_comps_hess = zeros(data.num_ri, Float64) 1453 r1_data.dip_jw_comps_hess = zeros(data.num_ri, Float64) 1454 r1_data.csa_jw_comps_hess = zeros(data.num_ri, Float64) 1455 1456 # Initialise the transformed relaxation values, gradients, and Hessians. 1457 r1_data.ri_prime = zeros(data.num_ri, Float64) 1458 r1_data.dri_prime = zeros((data.num_ri, data.total_num_params), Float64) 1459 r1_data.d2ri_prime = zeros((data.num_ri, data.total_num_params, data.total_num_params), Float64) 1460 1461 # Place a few function arrays in the data class for the calculation of the R1 value when an NOE data set exists but the R1 set does not. 1462 r1_data.create_dri_prime = data.create_dri_prime 1463 r1_data.create_d2ri_prime = data.create_d2ri_prime 1464 1465 # CSA, bond length, and Rex indecies. 1466 r1_data.csa_i = data.csa_i 1467 r1_data.r_i = data.r_i 1468 r1_data.rex_i = data.rex_i 1469 1470 # Place 'r1_data' into 'data'. 1471 data.r1_data = r1_data
1472 1473
1474 - def lm_dri(self):
1475 """Return the function used for Levenberg-Marquardt minimisation.""" 1476 1477 # Create dri. 1478 if self.param_set == 'mf' or self.param_set == 'local_tm': 1479 dri = self.data[0].dri 1480 elif self.param_set == 'diff': 1481 # Set the total dri gradient to zero. 1482 self.total_dri = self.total_dri * 0.0 1483 1484 # Ri indecies. 1485 ri_start_index = 0 1486 ri_end_index = 0 1487 1488 # Loop over the residues. 1489 for i in xrange(self.num_res): 1490 # Set self.data[i] to data. 1491 data = self.data[i] 1492 1493 # Increment Ri end index. 1494 ri_end_index = ri_end_index + data.num_ri 1495 1496 # Diffusion parameter part of the global generic model-free gradient. 1497 self.total_dri[0:self.diff_data.num_params, ri_start_index:ri_end_index] = self.total_dri[0:self.diff_data.num_params, ri_start_index:ri_end_index] + data.dri[0:self.diff_data.num_params] 1498 1499 # Increment Ri start index. 1500 ri_start_index = ri_start_index + data.num_ri 1501 1502 # dri. 1503 dri = self.total_dri 1504 1505 elif self.param_set == 'all': 1506 # Set the total dri gradient to zero. 1507 self.total_dri = self.total_dri * 0.0 1508 1509 # Ri indecies. 1510 ri_start_index = 0 1511 ri_end_index = 0 1512 1513 # Loop over the residues. 1514 for i in xrange(self.num_res): 1515 # Set self.data[i] to data. 1516 data = self.data[i] 1517 1518 # Increment Ri end index. 1519 ri_end_index = ri_end_index + data.num_ri 1520 1521 # Diffusion parameter part of the global generic model-free gradient. 1522 self.total_dri[0:self.diff_data.num_params, ri_start_index:ri_end_index] = self.total_dri[0:self.diff_data.num_params, ri_start_index:ri_end_index] + data.dri[0:self.diff_data.num_params] 1523 1524 # Model-free parameter part of the global generic model-free gradient. 1525 self.total_dri[data.start_index:data.end_index, ri_start_index:ri_end_index] = self.total_dri[data.start_index:data.end_index, ri_start_index:ri_end_index] + data.dri[self.diff_data.num_params:] 1526 1527 # Increment Ri start index. 1528 ri_start_index = ri_start_index + data.num_ri 1529 1530 # dri. 1531 dri = self.total_dri 1532 1533 # Make the proper Jacobian. 1534 dri = transpose(dri) 1535 1536 # Diagonal scaling. 1537 if self.scaling_flag: 1538 dri = matrixmultiply(dri, self.scaling_matrix) 1539 1540 # Return dri. 1541 return dri
1542 1543
1544 - def setup_equations(self, data):
1545 """Setup all the residue specific equations.""" 1546 1547 # Initialisation. 1548 ################# 1549 1550 # The number of diffusion parameters. 1551 if self.param_set != 'all': 1552 num_diff_params = 0 1553 elif self.diff_data.type == 'sphere': 1554 num_diff_params = 1 1555 elif self.diff_data.type == 'spheroid': 1556 num_diff_params = 4 1557 elif self.diff_data.type == 'ellipsoid': 1558 num_diff_params = 6 1559 1560 # Indecies. 1561 data.tm_i, data.tm_li = None, None 1562 data.s2_i, data.s2_li = None, None 1563 data.s2f_i, data.s2f_li = None, None 1564 data.s2s_i, data.s2s_li = None, None 1565 data.te_i, data.te_li = None, None 1566 data.tf_i, data.tf_li = None, None 1567 data.ts_i, data.ts_li = None, None 1568 data.rex_i, data.rex_li = None, None 1569 data.r_i, data.r_li = None, None 1570 data.csa_i, data.csa_li = None, None 1571 1572 # Set up the spectral density functions. 1573 ######################################## 1574 1575 # Create empty spectral density gradient and Hessian function data structures. 1576 data.calc_djw = [] 1577 data.calc_d2jw = [] 1578 for i in xrange(data.total_num_params): 1579 data.calc_djw.append(None) 1580 data.calc_d2jw.append([]) 1581 for j in xrange(data.total_num_params): 1582 data.calc_d2jw[i].append(None) 1583 1584 1585 # The original model-free equations {S2, te, Rex, r, CSA}. 1586 ########################################################## 1587 1588 if data.equations == 'mf_orig': 1589 # Find the indecies of the model-free parameters. 1590 for i in xrange(data.num_params): 1591 if data.param_types[i] == 'S2': 1592 data.s2_li = num_diff_params + i 1593 data.s2_i = self.param_index + i 1594 elif data.param_types[i] == 'te': 1595 data.te_li = num_diff_params + i 1596 data.te_i = self.param_index + i 1597 elif data.param_types[i] == 'Rex': 1598 data.rex_li = num_diff_params + i 1599 data.rex_i = self.param_index + i 1600 elif data.param_types[i] == 'r': 1601 data.r_li = num_diff_params + i 1602 data.r_i = self.param_index + i 1603 elif data.param_types[i] == 'CSA': 1604 data.csa_li = num_diff_params + i 1605 data.csa_i = self.param_index + i 1606 elif data.param_types[i] == 'tm': 1607 pass 1608 else: 1609 print "Unknown parameter." 1610 return 0 1611 1612 # Increment the parameter index. 1613 self.param_index = self.param_index + data.num_params 1614 1615 # Single residue minimisation with fixed diffusion parameters. 1616 if self.param_set == 'mf': 1617 # No model-free parameters {}. 1618 if data.s2_i == None and data.te_i == None: 1619 # Equation. 1620 data.calc_jw_comps = None 1621 data.calc_jw = calc_jw 1622 1623 # Gradient. 1624 data.calc_djw_comps = None 1625 1626 # Model-free parameters {S2}. 1627 elif data.s2_i != None and data.te_i == None: 1628 # Equation. 1629 data.calc_jw_comps = None 1630 data.calc_jw = calc_S2_jw 1631 1632 # Gradient. 1633 data.calc_djw_comps = None 1634 data.calc_djw[data.s2_li] = calc_S2_djw_dS2 1635 1636 # Model-free parameters {S2, te}. 1637 elif data.s2_i != None and data.te_i != None: 1638 # Equation. 1639 data.calc_jw_comps = calc_S2_te_jw_comps 1640 data.calc_jw = calc_S2_te_jw 1641 1642 # Gradient. 1643 data.calc_djw_comps = calc_S2_te_djw_comps 1644 data.calc_djw[data.s2_li] = calc_S2_te_djw_dS2 1645 data.calc_djw[data.te_li] = calc_S2_te_djw_dte 1646 1647 # Hessian. 1648 data.calc_d2jw[data.s2_li][data.te_li] = data.calc_d2jw[data.te_li][data.s2_li] = calc_S2_te_d2jw_dS2dte 1649 data.calc_d2jw[data.te_li][data.te_li] = calc_S2_te_d2jw_dte2 1650 1651 # Bad parameter combination. 1652 else: 1653 print "Invalid combination of parameters for the extended model-free equation." 1654 return 0 1655 1656 # Minimisation with variable diffusion parameters. 1657 else: 1658 # Diffusion parameters and no model-free parameters {}. 1659 if data.s2_i == None and data.te_i == None: 1660 # Equation. 1661 data.calc_jw_comps = None 1662 data.calc_jw = calc_jw 1663 1664 # Gradient. 1665 data.calc_djw_comps = calc_diff_djw_comps 1666 1667 # Diffusion as a sphere. 1668 if self.diff_data.type == 'sphere': 1669 # Gradient. 1670 data.calc_djw[0] = calc_diff_djw_dDj 1671 1672 # Hessian. 1673 data.calc_d2jw[0][0] = calc_diff_d2jw_dDjdDk 1674 1675 # Diffusion as a spheroid. 1676 elif self.diff_data.type == 'spheroid': 1677 # Gradient. 1678 data.calc_djw[0] = data.calc_djw[1] = calc_diff_djw_dDj 1679 data.calc_djw[2] = data.calc_djw[3] = calc_diff_djw_dOj 1680 1681 # Hessian. 1682 data.calc_d2jw[0][0] = calc_diff_d2jw_dDjdDk 1683 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_d2jw_dDjdDk 1684 data.calc_d2jw[1][1] = calc_diff_d2jw_dDjdDk 1685 1686 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_d2jw_dDjdOj 1687 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_d2jw_dDjdOj 1688 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_d2jw_dDjdOj 1689 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_d2jw_dDjdOj 1690 1691 data.calc_d2jw[2][2] = calc_diff_d2jw_dOjdOk 1692 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_d2jw_dOjdOk 1693 data.calc_d2jw[3][3] = calc_diff_d2jw_dOjdOk 1694 1695 # Diffusion as an ellipsoid. 1696 elif self.diff_data.type == 'ellipsoid': 1697 # Gradient. 1698 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_ellipsoid_djw_dDj 1699 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_djw_dOj 1700 1701 # Hessian. 1702 data.calc_d2jw[0][0] = calc_ellipsoid_d2jw_dDjdDk 1703 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_ellipsoid_d2jw_dDjdDk 1704 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_ellipsoid_d2jw_dDjdDk 1705 data.calc_d2jw[1][1] = calc_ellipsoid_d2jw_dDjdDk 1706 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_ellipsoid_d2jw_dDjdDk 1707 data.calc_d2jw[2][2] = calc_ellipsoid_d2jw_dDjdDk 1708 1709 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_ellipsoid_d2jw_dDjdOj 1710 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_ellipsoid_d2jw_dDjdOj 1711 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_ellipsoid_d2jw_dDjdOj 1712 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_ellipsoid_d2jw_dDjdOj 1713 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_ellipsoid_d2jw_dDjdOj 1714 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_ellipsoid_d2jw_dDjdOj 1715 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_ellipsoid_d2jw_dDjdOj 1716 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_ellipsoid_d2jw_dDjdOj 1717 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_ellipsoid_d2jw_dDjdOj 1718 1719 data.calc_d2jw[3][3] = calc_diff_d2jw_dOjdOk 1720 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_d2jw_dOjdOk 1721 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_d2jw_dOjdOk 1722 data.calc_d2jw[4][4] = calc_diff_d2jw_dOjdOk 1723 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_d2jw_dOjdOk 1724 data.calc_d2jw[5][5] = calc_diff_d2jw_dOjdOk 1725 1726 # Diffusion parameters and model-free parameters {S2}. 1727 elif data.s2_i != None and data.te_i == None: 1728 # Equation. 1729 data.calc_jw_comps = None 1730 data.calc_jw = calc_S2_jw 1731 1732 # Gradient. 1733 data.calc_djw_comps = calc_diff_djw_comps 1734 1735 if self.param_set != 'diff': 1736 # Gradient. 1737 data.calc_djw[data.s2_li] = calc_S2_djw_dS2 1738 1739 # Diffusion as a sphere. 1740 if self.diff_data.type == 'sphere': 1741 # Gradient. 1742 data.calc_djw[0] = calc_diff_S2_djw_dDj 1743 1744 # Hessian. 1745 data.calc_d2jw[0][0] = calc_diff_S2_d2jw_dDjdDk 1746 if self.param_set != 'diff': 1747 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_diff_S2_d2jw_dDjdS2 1748 1749 # Diffusion as a spheroid. 1750 elif self.diff_data.type == 'spheroid': 1751 # Gradient. 1752 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2_djw_dDj 1753 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2_djw_dOj 1754 1755 # Hessian. 1756 data.calc_d2jw[0][0] = calc_diff_S2_d2jw_dDjdDk 1757 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2_d2jw_dDjdDk 1758 data.calc_d2jw[1][1] = calc_diff_S2_d2jw_dDjdDk 1759 1760 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2_d2jw_dDjdOj 1761 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2_d2jw_dDjdOj 1762 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2_d2jw_dDjdOj 1763 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2_d2jw_dDjdOj 1764 1765 data.calc_d2jw[2][2] = calc_diff_S2_d2jw_dOjdOk 1766 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2_d2jw_dOjdOk 1767 data.calc_d2jw[3][3] = calc_diff_S2_d2jw_dOjdOk 1768 1769 if self.param_set != 'diff': 1770 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_diff_S2_d2jw_dDjdS2 1771 data.calc_d2jw[1][data.s2_li] = data.calc_d2jw[data.s2_li][1] = calc_diff_S2_d2jw_dDjdS2 1772 data.calc_d2jw[2][data.s2_li] = data.calc_d2jw[data.s2_li][2] = calc_diff_S2_d2jw_dOjdS2 1773 data.calc_d2jw[3][data.s2_li] = data.calc_d2jw[data.s2_li][3] = calc_diff_S2_d2jw_dOjdS2 1774 1775 # Diffusion as an ellipsoid. 1776 elif self.diff_data.type == 'ellipsoid': 1777 # Gradient. 1778 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_ellipsoid_S2_djw_dDj 1779 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2_djw_dOj 1780 1781 # Hessian. 1782 data.calc_d2jw[0][0] = calc_ellipsoid_S2_d2jw_dDjdDk 1783 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_ellipsoid_S2_d2jw_dDjdDk 1784 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_ellipsoid_S2_d2jw_dDjdDk 1785 data.calc_d2jw[1][1] = calc_ellipsoid_S2_d2jw_dDjdDk 1786 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_ellipsoid_S2_d2jw_dDjdDk 1787 data.calc_d2jw[2][2] = calc_ellipsoid_S2_d2jw_dDjdDk 1788 1789 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_ellipsoid_S2_d2jw_dDjdOj 1790 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_ellipsoid_S2_d2jw_dDjdOj 1791 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_ellipsoid_S2_d2jw_dDjdOj 1792 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_ellipsoid_S2_d2jw_dDjdOj 1793 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_ellipsoid_S2_d2jw_dDjdOj 1794 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_ellipsoid_S2_d2jw_dDjdOj 1795 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_ellipsoid_S2_d2jw_dDjdOj 1796 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_ellipsoid_S2_d2jw_dDjdOj 1797 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_ellipsoid_S2_d2jw_dDjdOj 1798 1799 data.calc_d2jw[3][3] = calc_diff_S2_d2jw_dOjdOk 1800 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2_d2jw_dOjdOk 1801 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2_d2jw_dOjdOk 1802 data.calc_d2jw[4][4] = calc_diff_S2_d2jw_dOjdOk 1803 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2_d2jw_dOjdOk 1804 data.calc_d2jw[5][5] = calc_diff_S2_d2jw_dOjdOk 1805 1806 if self.param_set != 'diff': 1807 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_ellipsoid_S2_d2jw_dDjdS2 1808 data.calc_d2jw[1][data.s2_li] = data.calc_d2jw[data.s2_li][1] = calc_ellipsoid_S2_d2jw_dDjdS2 1809 data.calc_d2jw[2][data.s2_li] = data.calc_d2jw[data.s2_li][2] = calc_ellipsoid_S2_d2jw_dDjdS2 1810 data.calc_d2jw[3][data.s2_li] = data.calc_d2jw[data.s2_li][3] = calc_diff_S2_d2jw_dOjdS2 1811 data.calc_d2jw[4][data.s2_li] = data.calc_d2jw[data.s2_li][4] = calc_diff_S2_d2jw_dOjdS2 1812 data.calc_d2jw[5][data.s2_li] = data.calc_d2jw[data.s2_li][5] = calc_diff_S2_d2jw_dOjdS2 1813 1814 1815 # Diffusion parameters and model-free parameters {S2, te}. 1816 elif data.s2_i != None and data.te_i != None: 1817 # Equation. 1818 data.calc_jw_comps = calc_S2_te_jw_comps 1819 data.calc_jw = calc_S2_te_jw 1820 1821 # Gradient. 1822 data.calc_djw_comps = calc_diff_S2_te_djw_comps 1823 1824 if self.param_set != 'diff': 1825 # Gradient. 1826 data.calc_djw[data.s2_li] = calc_S2_te_djw_dS2 1827 data.calc_djw[data.te_li] = calc_S2_te_djw_dte 1828 1829 # Hessian. 1830 data.calc_d2jw[data.s2_li][data.te_li] = data.calc_d2jw[data.te_li][data.s2_li] = calc_S2_te_d2jw_dS2dte 1831 data.calc_d2jw[data.te_li][data.te_li] = calc_S2_te_d2jw_dte2 1832 1833 # Diffusion as a sphere. 1834 if self.diff_data.type == 'sphere': 1835 # Gradient. 1836 data.calc_djw[0] = calc_diff_S2_te_djw_dDj 1837 1838 # Hessian. 1839 data.calc_d2jw[0][0] = calc_diff_S2_te_d2jw_dDjdDk 1840 if self.param_set != 'diff': 1841 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_diff_S2_te_d2jw_dDjdS2 1842 data.calc_d2jw[0][data.te_li] = data.calc_d2jw[data.te_li][0] = calc_diff_S2_te_d2jw_dDjdte 1843 1844 # Diffusion as a spheroid. 1845 elif self.diff_data.type == 'spheroid': 1846 # Gradient. 1847 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2_te_djw_dDj 1848 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2_te_djw_dOj 1849 1850 # Hessian. 1851 data.calc_d2jw[0][0] = calc_diff_S2_te_d2jw_dDjdDk 1852 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2_te_d2jw_dDjdDk 1853 data.calc_d2jw[1][1] = calc_diff_S2_te_d2jw_dDjdDk 1854 1855 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2_te_d2jw_dDjdOj 1856 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2_te_d2jw_dDjdOj 1857 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2_te_d2jw_dDjdOj 1858 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2_te_d2jw_dDjdOj 1859 1860 data.calc_d2jw[2][2] = calc_diff_S2_te_d2jw_dOjdOk 1861 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2_te_d2jw_dOjdOk 1862 data.calc_d2jw[3][3] = calc_diff_S2_te_d2jw_dOjdOk 1863 1864 if self.param_set != 'diff': 1865 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_diff_S2_te_d2jw_dDjdS2 1866 data.calc_d2jw[1][data.s2_li] = data.calc_d2jw[data.s2_li][1] = calc_diff_S2_te_d2jw_dDjdS2 1867 data.calc_d2jw[2][data.s2_li] = data.calc_d2jw[data.s2_li][2] = calc_diff_S2_te_d2jw_dOjdS2 1868 data.calc_d2jw[3][data.s2_li] = data.calc_d2jw[data.s2_li][3] = calc_diff_S2_te_d2jw_dOjdS2 1869 1870 data.calc_d2jw[0][data.te_li] = data.calc_d2jw[data.te_li][0] = calc_diff_S2_te_d2jw_dDjdte 1871 data.calc_d2jw[1][data.te_li] = data.calc_d2jw[data.te_li][1] = calc_diff_S2_te_d2jw_dDjdte 1872 data.calc_d2jw[2][data.te_li] = data.calc_d2jw[data.te_li][2] = calc_diff_S2_te_d2jw_dOjdte 1873 data.calc_d2jw[3][data.te_li] = data.calc_d2jw[data.te_li][3] = calc_diff_S2_te_d2jw_dOjdte 1874 1875 # Diffusion as an ellipsoid. 1876 elif self.diff_data.type == 'ellipsoid': 1877 # Gradient. 1878 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_ellipsoid_S2_te_djw_dDj 1879 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2_te_djw_dOj 1880 1881 # Hessian. 1882 data.calc_d2jw[0][0] = calc_ellipsoid_S2_te_d2jw_dDjdDk 1883 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_ellipsoid_S2_te_d2jw_dDjdDk 1884 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_ellipsoid_S2_te_d2jw_dDjdDk 1885 data.calc_d2jw[1][1] = calc_ellipsoid_S2_te_d2jw_dDjdDk 1886 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_ellipsoid_S2_te_d2jw_dDjdDk 1887 data.calc_d2jw[2][2] = calc_ellipsoid_S2_te_d2jw_dDjdDk 1888 1889 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_ellipsoid_S2_te_d2jw_dDjdOj 1890 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_ellipsoid_S2_te_d2jw_dDjdOj 1891 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_ellipsoid_S2_te_d2jw_dDjdOj 1892 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_ellipsoid_S2_te_d2jw_dDjdOj 1893 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_ellipsoid_S2_te_d2jw_dDjdOj 1894 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_ellipsoid_S2_te_d2jw_dDjdOj 1895 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_ellipsoid_S2_te_d2jw_dDjdOj 1896 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_ellipsoid_S2_te_d2jw_dDjdOj 1897 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_ellipsoid_S2_te_d2jw_dDjdOj 1898 1899 data.calc_d2jw[3][3] = calc_diff_S2_te_d2jw_dOjdOk 1900 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2_te_d2jw_dOjdOk 1901 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2_te_d2jw_dOjdOk 1902 data.calc_d2jw[4][4] = calc_diff_S2_te_d2jw_dOjdOk 1903 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2_te_d2jw_dOjdOk 1904 data.calc_d2jw[5][5] = calc_diff_S2_te_d2jw_dOjdOk 1905 1906 if self.param_set != 'diff': 1907 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_ellipsoid_S2_te_d2jw_dDjdS2 1908 data.calc_d2jw[1][data.s2_li] = data.calc_d2jw[data.s2_li][1] = calc_ellipsoid_S2_te_d2jw_dDjdS2 1909 data.calc_d2jw[2][data.s2_li] = data.calc_d2jw[data.s2_li][2] = calc_ellipsoid_S2_te_d2jw_dDjdS2 1910 data.calc_d2jw[3][data.s2_li] = data.calc_d2jw[data.s2_li][3] = calc_diff_S2_te_d2jw_dOjdS2 1911 data.calc_d2jw[4][data.s2_li] = data.calc_d2jw[data.s2_li][4] = calc_diff_S2_te_d2jw_dOjdS2 1912 data.calc_d2jw[5][data.s2_li] = data.calc_d2jw[data.s2_li][5] = calc_diff_S2_te_d2jw_dOjdS2 1913 1914 data.calc_d2jw[0][data.te_li] = data.calc_d2jw[data.te_li][0] = calc_ellipsoid_S2_te_d2jw_dDjdte 1915 data.calc_d2jw[1][data.te_li] = data.calc_d2jw[data.te_li][1] = calc_ellipsoid_S2_te_d2jw_dDjdte 1916 data.calc_d2jw[2][data.te_li] = data.calc_d2jw[data.te_li][2] = calc_ellipsoid_S2_te_d2jw_dDjdte 1917 data.calc_d2jw[3][data.te_li] = data.calc_d2jw[data.te_li][3] = calc_diff_S2_te_d2jw_dOjdte 1918 data.calc_d2jw[4][data.te_li] = data.calc_d2jw[data.te_li][4] = calc_diff_S2_te_d2jw_dOjdte 1919 data.calc_d2jw[5][data.te_li] = data.calc_d2jw[data.te_li][5] = calc_diff_S2_te_d2jw_dOjdte 1920 1921 # Bad parameter combination. 1922 else: 1923 print "Invalid combination of parameters for the extended model-free equation." 1924 return 0 1925 1926 1927 1928 # The extended model-free equations {S2f, tf, S2, ts, Rex, r, CSA}. 1929 ################################################################### 1930 1931 elif data.equations == 'mf_ext': 1932 # Find the indecies of the model-free parameters. 1933 for i in xrange(data.num_params): 1934 if data.param_types[i] == 'S2f': 1935 data.s2f_li = num_diff_params + i 1936 data.s2f_i = self.param_index + i 1937 elif data.param_types[i] == 'tf': 1938 data.tf_li = num_diff_params + i 1939 data.tf_i = self.param_index + i 1940 elif data.param_types[i] == 'S2': 1941 data.s2_li = num_diff_params + i 1942 data.s2_i = self.param_index + i 1943 elif data.param_types[i] == 'ts': 1944 data.ts_li = num_diff_params + i 1945 data.ts_i = self.param_index + i 1946 elif data.param_types[i] == 'Rex': 1947 data.rex_li = num_diff_params + i 1948 data.rex_i = self.param_index + i 1949 elif data.param_types[i] == 'r': 1950 data.r_li = num_diff_params + i 1951 data.r_i = self.param_index + i 1952 elif data.param_types[i] == 'CSA': 1953 data.csa_li = num_diff_params + i 1954 data.csa_i = self.param_index + i 1955 elif data.param_types[i] == 'tm': 1956 pass 1957 else: 1958 print "Unknown parameter." 1959 return 0 1960 1961 # Increment the parameter index. 1962 self.param_index = self.param_index + data.num_params 1963 1964 # Single residue minimisation with fixed diffusion parameters. 1965 if self.param_set == 'mf': 1966 # Model-free parameters {S2f, S2, ts}. 1967 if data.s2f_i != None and data.tf_i == None and data.s2_i != None and data.ts_i != None: 1968 # Equation. 1969 data.calc_jw_comps = calc_S2f_S2_ts_jw_comps 1970 data.calc_jw = calc_S2f_S2_ts_jw 1971 1972 # Gradient. 1973 data.calc_djw_comps = calc_S2f_S2_ts_djw_comps 1974 data.calc_djw[data.s2f_li] = calc_S2f_S2_ts_djw_dS2f 1975 data.calc_djw[data.s2_li] = calc_S2f_S2_ts_djw_dS2 1976 data.calc_djw[data.ts_li] = calc_S2f_S2_ts_djw_dts 1977 1978 # Hessian. 1979 data.calc_d2jw[data.s2f_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2f_li] = calc_S2f_S2_ts_d2jw_dS2fdts 1980 data.calc_d2jw[data.s2_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2_li] = calc_S2f_S2_ts_d2jw_dS2dts 1981 data.calc_d2jw[data.ts_li][data.ts_li] = calc_S2f_S2_ts_d2jw_dts2 1982 1983 # Model-free parameters {S2f, tf, S2, ts}. 1984 elif data.s2f_i != None and data.tf_i != None and data.s2_i != None and data.ts_i != None: 1985 # Equation. 1986 data.calc_jw_comps = calc_S2f_tf_S2_ts_jw_comps 1987 data.calc_jw = calc_S2f_tf_S2_ts_jw 1988 1989 # Gradient. 1990 data.calc_djw_comps = calc_S2f_tf_S2_ts_djw_comps 1991 data.calc_djw[data.s2f_li] = calc_S2f_tf_S2_ts_djw_dS2f 1992 data.calc_djw[data.tf_li] = calc_S2f_tf_S2_ts_djw_dtf 1993 data.calc_djw[data.s2_li] = calc_S2f_S2_ts_djw_dS2 1994 data.calc_djw[data.ts_li] = calc_S2f_S2_ts_djw_dts 1995 1996 # Hessian. 1997 data.calc_d2jw[data.s2f_li][data.tf_li] = data.calc_d2jw[data.tf_li][data.s2f_li] = calc_S2f_tf_S2_ts_d2jw_dS2fdtf 1998 data.calc_d2jw[data.s2f_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2f_li] = calc_S2f_S2_ts_d2jw_dS2fdts 1999 data.calc_d2jw[data.s2_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2_li] = calc_S2f_S2_ts_d2jw_dS2dts 2000 data.calc_d2jw[data.tf_li][data.tf_li] = calc_S2f_tf_S2_ts_d2jw_dtf2 2001 data.calc_d2jw[data.ts_li][data.ts_li] = calc_S2f_S2_ts_d2jw_dts2 2002 2003 # Bad parameter combination. 2004 else: 2005 print "Invalid combination of parameters for the extended model-free equation." 2006 return 0 2007 2008 # Minimisation with variable diffusion parameters. 2009 else: 2010 # Diffusion parameters and model-free parameters {S2f, S2, ts}. 2011 if data.s2f_i != None and data.tf_i == None and data.s2_i != None and data.ts_i != None: 2012 # Equation. 2013 data.calc_jw_comps = calc_S2f_S2_ts_jw_comps 2014 data.calc_jw = calc_S2f_S2_ts_jw 2015 2016 # Gradient. 2017 data.calc_djw_comps = calc_diff_S2f_S2_ts_djw_comps 2018 2019 if self.param_set != 'diff': 2020 # Gradient. 2021 data.calc_djw[data.s2f_li] = calc_S2f_S2_ts_djw_dS2f 2022 data.calc_djw[data.s2_li] = calc_S2f_S2_ts_djw_dS2 2023 data.calc_djw[data.ts_li] = calc_S2f_S2_ts_djw_dts 2024 2025 # Hessian. 2026 data.calc_d2jw[data.s2f_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2f_li] = calc_S2f_S2_ts_d2jw_dS2fdts 2027 data.calc_d2jw[data.s2_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2_li] = calc_S2f_S2_ts_d2jw_dS2dts 2028 data.calc_d2jw[data.ts_li][data.ts_li] = calc_S2f_S2_ts_d2jw_dts2 2029 2030 # Diffusion as a sphere. 2031 if self.diff_data.type == 'sphere': 2032 # Gradient. 2033 data.calc_djw[0] = calc_diff_S2f_S2_ts_djw_dDj 2034 2035 # Hessian. 2036 data.calc_d2jw[0][0] = calc_diff_S2f_S2_ts_d2jw_dDjdDk 2037 if self.param_set != 'diff': 2038 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2f 2039 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2040 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2041 2042 # Diffusion as a spheroid. 2043 elif self.diff_data.type == 'spheroid': 2044 # Gradient. 2045 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2f_S2_ts_djw_dDj 2046 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2f_S2_ts_djw_dOj 2047 2048 # Hessian. 2049 data.calc_d2jw[0][0] = calc_diff_S2f_S2_ts_d2jw_dDjdDk 2050 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2f_S2_ts_d2jw_dDjdDk 2051 data.calc_d2jw[1][1] = calc_diff_S2f_S2_ts_d2jw_dDjdDk 2052 2053 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2f_S2_ts_d2jw_dDjdOj 2054 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2f_S2_ts_d2jw_dDjdOj 2055 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2f_S2_ts_d2jw_dDjdOj 2056 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2f_S2_ts_d2jw_dDjdOj 2057 2058 data.calc_d2jw[2][2] = calc_diff_S2f_S2_ts_d2jw_dOjdOk 2059 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2f_S2_ts_d2jw_dOjdOk 2060 data.calc_d2jw[3][3] = calc_diff_S2f_S2_ts_d2jw_dOjdOk 2061 2062 if self.param_set != 'diff': 2063 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2f 2064 data.calc_d2jw[1][data.s2f_li] = data.calc_d2jw[data.s2f_li][1] = calc_diff_S2f_S2_ts_d2jw_dDjdS2f 2065 data.calc_d2jw[2][data.s2f_li] = data.calc_d2jw[data.s2f_li][2] = calc_diff_S2f_S2_ts_d2jw_dOjdS2f 2066 data.calc_d2jw[3][data.s2f_li] = data.calc_d2jw[data.s2f_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdS2f 2067 2068 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2069 data.calc_d2jw[1][data.s2_li] = data.calc_d2jw[data.s2_li][1] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2070 data.calc_d2jw[2][data.s2_li] = data.calc_d2jw[data.s2_li][2] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2071 data.calc_d2jw[3][data.s2_li] = data.calc_d2jw[data.s2_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2072 2073 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2074 data.calc_d2jw[1][data.ts_li] = data.calc_d2jw[data.ts_li][1] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2075 data.calc_d2jw[2][data.ts_li] = data.calc_d2jw[data.ts_li][2] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2076 data.calc_d2jw[3][data.ts_li] = data.calc_d2jw[data.ts_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2077 2078 # Diffusion as an ellipsoid. 2079 elif self.diff_data.type == 'ellipsoid': 2080 # Gradient. 2081 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_ellipsoid_S2f_S2_ts_djw_dDj 2082 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2f_S2_ts_djw_dOj 2083 2084 # Hessian. 2085 data.calc_d2jw[0][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdDk 2086 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdDk 2087 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdDk 2088 data.calc_d2jw[1][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdDk 2089 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdDk 2090 data.calc_d2jw[2][2] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdDk 2091 2092 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdOj 2093 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdOj 2094 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdOj 2095 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdOj 2096 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdOj 2097 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdOj 2098 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdOj 2099 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdOj 2100 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdOj 2101 2102 data.calc_d2jw[3][3] = calc_diff_S2f_S2_ts_d2jw_dOjdOk 2103 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2f_S2_ts_d2jw_dOjdOk 2104 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2f_S2_ts_d2jw_dOjdOk 2105 data.calc_d2jw[4][4] = calc_diff_S2f_S2_ts_d2jw_dOjdOk 2106 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2f_S2_ts_d2jw_dOjdOk 2107 data.calc_d2jw[5][5] = calc_diff_S2f_S2_ts_d2jw_dOjdOk 2108 2109 if self.param_set != 'diff': 2110 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdS2f 2111 data.calc_d2jw[1][data.s2f_li] = data.calc_d2jw[data.s2f_li][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdS2f 2112 data.calc_d2jw[2][data.s2f_li] = data.calc_d2jw[data.s2f_li][2] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdS2f 2113 data.calc_d2jw[3][data.s2f_li] = data.calc_d2jw[data.s2f_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdS2f 2114 data.calc_d2jw[4][data.s2f_li] = data.calc_d2jw[data.s2f_li][4] = calc_diff_S2f_S2_ts_d2jw_dOjdS2f 2115 data.calc_d2jw[5][data.s2f_li] = data.calc_d2jw[data.s2f_li][5] = calc_diff_S2f_S2_ts_d2jw_dOjdS2f 2116 2117 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdS2 2118 data.calc_d2jw[1][data.s2_li] = data.calc_d2jw[data.s2_li][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdS2 2119 data.calc_d2jw[2][data.s2_li] = data.calc_d2jw[data.s2_li][2] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdS2 2120 data.calc_d2jw[3][data.s2_li] = data.calc_d2jw[data.s2_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2121 data.calc_d2jw[4][data.s2_li] = data.calc_d2jw[data.s2_li][4] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2122 data.calc_d2jw[5][data.s2_li] = data.calc_d2jw[data.s2_li][5] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2123 2124 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdts 2125 data.calc_d2jw[1][data.ts_li] = data.calc_d2jw[data.ts_li][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdts 2126 data.calc_d2jw[2][data.ts_li] = data.calc_d2jw[data.ts_li][2] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdts 2127 data.calc_d2jw[3][data.ts_li] = data.calc_d2jw[data.ts_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2128 data.calc_d2jw[4][data.ts_li] = data.calc_d2jw[data.ts_li][4] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2129 data.calc_d2jw[5][data.ts_li] = data.calc_d2jw[data.ts_li][5] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2130 2131 2132 # Diffusion parameters and model-free parameters {S2f, tf, S2, ts}. 2133 elif data.s2f_i != None and data.tf_i != None and data.s2_i != None and data.ts_i != None: 2134 # Equation. 2135 data.calc_jw_comps = calc_S2f_tf_S2_ts_jw_comps 2136 data.calc_jw = calc_S2f_tf_S2_ts_jw 2137 2138 # Gradient. 2139 data.calc_djw_comps = calc_diff_S2f_tf_S2_ts_djw_comps 2140 2141 if self.param_set != 'diff': 2142 # Gradient. 2143 data.calc_djw[data.s2f_li] = calc_S2f_tf_S2_ts_djw_dS2f 2144 data.calc_djw[data.tf_li] = calc_S2f_tf_S2_ts_djw_dtf 2145 data.calc_djw[data.s2_li] = calc_S2f_S2_ts_djw_dS2 2146 data.calc_djw[data.ts_li] = calc_S2f_S2_ts_djw_dts 2147 2148 # Hessian. 2149 data.calc_d2jw[data.s2f_li][data.tf_li] = data.calc_d2jw[data.tf_li][data.s2f_li] = calc_S2f_tf_S2_ts_d2jw_dS2fdtf 2150 data.calc_d2jw[data.s2f_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2f_li] = calc_S2f_S2_ts_d2jw_dS2fdts 2151 data.calc_d2jw[data.s2_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2_li] = calc_S2f_S2_ts_d2jw_dS2dts 2152 data.calc_d2jw[data.tf_li][data.tf_li] = calc_S2f_tf_S2_ts_d2jw_dtf2 2153 data.calc_d2jw[data.ts_li][data.ts_li] = calc_S2f_S2_ts_d2jw_dts2 2154 2155 # Diffusion as a sphere. 2156 if self.diff_data.type == 'sphere': 2157 # Gradient. 2158 data.calc_djw[0] = calc_diff_S2f_tf_S2_ts_djw_dDj 2159 2160 # Hessian. 2161 data.calc_d2jw[0][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdDk 2162 if self.param_set != 'diff': 2163 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdS2f 2164 data.calc_d2jw[0][data.tf_li] = data.calc_d2jw[data.tf_li][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdtf 2165 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2166 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2167 2168 # Diffusion as a spheroid. 2169 elif self.diff_data.type == 'spheroid': 2170 # Gradient. 2171 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2f_tf_S2_ts_djw_dDj 2172 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2f_tf_S2_ts_djw_dOj 2173 2174 # Hessian. 2175 data.calc_d2jw[0][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdDk 2176 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdDk 2177 data.calc_d2jw[1][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdDk 2178 2179 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdOj 2180 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdOj 2181 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdOj 2182 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdOj 2183 2184 data.calc_d2jw[2][2] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdOk 2185 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdOk 2186 data.calc_d2jw[3][3] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdOk 2187 2188 if self.param_set != 'diff': 2189 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdS2f 2190 data.calc_d2jw[1][data.s2f_li] = data.calc_d2jw[data.s2f_li][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdS2f 2191 data.calc_d2jw[2][data.s2f_li] = data.calc_d2jw[data.s2f_li][2] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdS2f 2192 data.calc_d2jw[3][data.s2f_li] = data.calc_d2jw[data.s2f_li][3] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdS2f 2193 2194 data.calc_d2jw[0][data.tf_li] = data.calc_d2jw[data.tf_li][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdtf 2195 data.calc_d2jw[1][data.tf_li] = data.calc_d2jw[data.tf_li][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdtf 2196 data.calc_d2jw[2][data.tf_li] = data.calc_d2jw[data.tf_li][2] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdtf 2197 data.calc_d2jw[3][data.tf_li] = data.calc_d2jw[data.tf_li][3] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdtf 2198 2199 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2200 data.calc_d2jw[1][data.s2_li] = data.calc_d2jw[data.s2_li][1] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2201 data.calc_d2jw[2][data.s2_li] = data.calc_d2jw[data.s2_li][2] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2202 data.calc_d2jw[3][data.s2_li] = data.calc_d2jw[data.s2_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2203 2204 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2205 data.calc_d2jw[1][data.ts_li] = data.calc_d2jw[data.ts_li][1] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2206 data.calc_d2jw[2][data.ts_li] = data.calc_d2jw[data.ts_li][2] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2207 data.calc_d2jw[3][data.ts_li] = data.calc_d2jw[data.ts_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2208 2209 # Diffusion as an ellipsoid. 2210 elif self.diff_data.type == 'ellipsoid': 2211 # Gradient. 2212 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_ellipsoid_S2f_tf_S2_ts_djw_dDj 2213 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2f_tf_S2_ts_djw_dOj 2214 2215 # Hessian. 2216 data.calc_d2jw[0][0] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdDk 2217 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdDk 2218 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdDk 2219 data.calc_d2jw[1][1] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdDk 2220 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdDk 2221 data.calc_d2jw[2][2] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdDk 2222 2223 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdOj 2224 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdOj 2225 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdOj 2226 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdOj 2227 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdOj 2228 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdOj 2229 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdOj 2230 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdOj 2231 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdOj 2232 2233 data.calc_d2jw[3][3] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdOk 2234 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdOk 2235 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdOk 2236 data.calc_d2jw[4][4] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdOk 2237 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdOk 2238 data.calc_d2jw[5][5] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdOk 2239 2240 if self.param_set != 'diff': 2241 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdS2f 2242 data.calc_d2jw[1][data.s2f_li] = data.calc_d2jw[data.s2f_li][1] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdS2f 2243 data.calc_d2jw[2][data.s2f_li] = data.calc_d2jw[data.s2f_li][2] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdS2f 2244 data.calc_d2jw[3][data.s2f_li] = data.calc_d2jw[data.s2f_li][3] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdS2f 2245 data.calc_d2jw[4][data.s2f_li] = data.calc_d2jw[data.s2f_li][4] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdS2f 2246 data.calc_d2jw[5][data.s2f_li] = data.calc_d2jw[data.s2f_li][5] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdS2f 2247 2248 data.calc_d2jw[0][data.tf_li] = data.calc_d2jw[data.tf_li][0] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdtf 2249 data.calc_d2jw[1][data.tf_li] = data.calc_d2jw[data.tf_li][1] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdtf 2250 data.calc_d2jw[2][data.tf_li] = data.calc_d2jw[data.tf_li][2] = calc_ellipsoid_S2f_tf_S2_ts_d2jw_dDjdtf 2251 data.calc_d2jw[3][data.tf_li] = data.calc_d2jw[data.tf_li][3] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdtf 2252 data.calc_d2jw[4][data.tf_li] = data.calc_d2jw[data.tf_li][4] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdtf 2253 data.calc_d2jw[5][data.tf_li] = data.calc_d2jw[data.tf_li][5] = calc_diff_S2f_tf_S2_ts_d2jw_dOjdtf 2254 2255 data.calc_d2jw[0][data.s2_li] = data.calc_d2jw[data.s2_li][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdS2 2256 data.calc_d2jw[1][data.s2_li] = data.calc_d2jw[data.s2_li][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdS2 2257 data.calc_d2jw[2][data.s2_li] = data.calc_d2jw[data.s2_li][2] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdS2 2258 data.calc_d2jw[3][data.s2_li] = data.calc_d2jw[data.s2_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2259 data.calc_d2jw[4][data.s2_li] = data.calc_d2jw[data.s2_li][4] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2260 data.calc_d2jw[5][data.s2_li] = data.calc_d2jw[data.s2_li][5] = calc_diff_S2f_S2_ts_d2jw_dOjdS2 2261 2262 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdts 2263 data.calc_d2jw[1][data.ts_li] = data.calc_d2jw[data.ts_li][1] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdts 2264 data.calc_d2jw[2][data.ts_li] = data.calc_d2jw[data.ts_li][2] = calc_ellipsoid_S2f_S2_ts_d2jw_dDjdts 2265 data.calc_d2jw[3][data.ts_li] = data.calc_d2jw[data.ts_li][3] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2266 data.calc_d2jw[4][data.ts_li] = data.calc_d2jw[data.ts_li][4] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2267 data.calc_d2jw[5][data.ts_li] = data.calc_d2jw[data.ts_li][5] = calc_diff_S2f_S2_ts_d2jw_dOjdts 2268 2269 # Bad parameter combination. 2270 else: 2271 print "Invalid combination of parameters for the extended model-free equation." 2272 return 0 2273 2274 2275 # The extended 2 model-free equations {tm, S2f, tf, S2s, ts, Rex, r, CSA}. 2276 ######################################################################### 2277 2278 elif data.equations == 'mf_ext2': 2279 # Find the indecies of the model-free parameters. 2280 for i in xrange(data.num_params): 2281 if data.param_types[i] == 'S2f': 2282 data.s2f_li = num_diff_params + i 2283 data.s2f_i = self.param_index + i 2284 elif data.param_types[i] == 'tf': 2285 data.tf_li = num_diff_params + i 2286 data.tf_i = self.param_index + i 2287 elif data.param_types[i] == 'S2s': 2288 data.s2s_li = num_diff_params + i 2289 data.s2s_i = self.param_index + i 2290 elif data.param_types[i] == 'ts': 2291 data.ts_li = num_diff_params + i 2292 data.ts_i = self.param_index + i 2293 elif data.param_types[i] == 'Rex': 2294 data.rex_li = num_diff_params + i 2295 data.rex_i = self.param_index + i 2296 elif data.param_types[i] == 'r': 2297 data.r_li = num_diff_params + i 2298 data.r_i = self.param_index + i 2299 elif data.param_types[i] == 'CSA': 2300 data.csa_li = num_diff_params + i 2301 data.csa_i = self.param_index + i 2302 elif data.param_types[i] == 'tm': 2303 pass 2304 else: 2305 print "Unknown parameter." 2306 return 0 2307 2308 # Increment the parameter index. 2309 self.param_index = self.param_index + data.num_params 2310 2311 # Single residue minimisation with fixed diffusion parameters. 2312 if self.param_set == 'mf': 2313 # Model-free parameters {S2f, S2s, ts}. 2314 if data.s2f_i != None and data.tf_i == None and data.s2s_i != None and data.ts_i != None: 2315 # Equation. 2316 data.calc_jw_comps = calc_S2f_S2s_ts_jw_comps 2317 data.calc_jw = calc_S2f_S2s_ts_jw 2318 2319 # Gradient. 2320 data.calc_djw_comps = calc_S2f_S2s_ts_djw_comps 2321 data.calc_djw[data.s2f_li] = calc_S2f_S2s_ts_djw_dS2f 2322 data.calc_djw[data.s2s_li] = calc_S2f_S2s_ts_djw_dS2s 2323 data.calc_djw[data.ts_li] = calc_S2f_S2s_ts_djw_dts 2324 2325 # Hessian. 2326 data.calc_d2jw[data.s2f_li][data.s2s_li] = data.calc_d2jw[data.s2s_li][data.s2f_li] = calc_S2f_S2s_ts_d2jw_dS2fdS2s 2327 data.calc_d2jw[data.s2f_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2f_li] = calc_S2f_S2s_ts_d2jw_dS2fdts 2328 data.calc_d2jw[data.s2s_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2s_li] = calc_S2f_S2s_ts_d2jw_dS2sdts 2329 data.calc_d2jw[data.ts_li][data.ts_li] = calc_S2f_S2s_ts_d2jw_dts2 2330 2331 # Model-free parameters {S2f, tf, S2s, ts}. 2332 elif data.s2f_i != None and data.tf_i != None and data.s2s_i != None and data.ts_i != None: 2333 # Equation. 2334 data.calc_jw_comps = calc_S2f_tf_S2s_ts_jw_comps 2335 data.calc_jw = calc_S2f_tf_S2s_ts_jw 2336 2337 # Gradient. 2338 data.calc_djw_comps = calc_S2f_tf_S2s_ts_djw_comps 2339 data.calc_djw[data.s2f_li] = calc_S2f_tf_S2s_ts_djw_dS2f 2340 data.calc_djw[data.tf_li] = calc_S2f_tf_S2s_ts_djw_dtf 2341 data.calc_djw[data.s2s_li] = calc_S2f_tf_S2s_ts_djw_dS2s 2342 data.calc_djw[data.ts_li] = calc_S2f_tf_S2s_ts_djw_dts 2343 2344 # Hessian. 2345 data.calc_d2jw[data.s2f_li][data.s2s_li] = data.calc_d2jw[data.s2s_li][data.s2f_li] = calc_S2f_S2s_ts_d2jw_dS2fdS2s 2346 data.calc_d2jw[data.s2f_li][data.tf_li] = data.calc_d2jw[data.tf_li][data.s2f_li] = calc_S2f_tf_S2s_ts_d2jw_dS2fdtf 2347 data.calc_d2jw[data.s2f_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2f_li] = calc_S2f_tf_S2s_ts_d2jw_dS2fdts 2348 data.calc_d2jw[data.s2s_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2s_li] = calc_S2f_tf_S2s_ts_d2jw_dS2sdts 2349 data.calc_d2jw[data.tf_li][data.tf_li] = calc_S2f_tf_S2s_ts_d2jw_dtf2 2350 data.calc_d2jw[data.ts_li][data.ts_li] = calc_S2f_tf_S2s_ts_d2jw_dts2 2351 2352 # Bad parameter combination. 2353 else: 2354 print "Invalid combination of parameters for the extended model-free equation." 2355 return 0 2356 2357 # Minimisation with variable diffusion parameters. 2358 else: 2359 # Diffusion parameters and model-free parameters {S2f, S2s, ts}. 2360 if data.s2f_i != None and data.tf_i == None and data.s2s_i != None and data.ts_i != None: 2361 # Equation. 2362 data.calc_jw_comps = calc_diff_S2f_S2s_ts_jw_comps 2363 data.calc_jw = calc_S2f_S2s_ts_jw 2364 2365 # Gradient. 2366 data.calc_djw_comps = calc_diff_S2f_S2s_ts_djw_comps 2367 2368 if self.param_set != 'diff': 2369 # Gradient. 2370 data.calc_djw[data.s2f_li] = calc_diff_S2f_S2s_ts_djw_dS2f 2371 data.calc_djw[data.s2s_li] = calc_diff_S2f_S2s_ts_djw_dS2s 2372 data.calc_djw[data.ts_li] = calc_diff_S2f_S2s_ts_djw_dts 2373 2374 # Hessian. 2375 data.calc_d2jw[data.s2f_li][data.s2s_li] = data.calc_d2jw[data.s2s_li][data.s2f_li] = calc_S2f_S2s_ts_d2jw_dS2fdS2s 2376 data.calc_d2jw[data.s2f_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2f_li] = calc_diff_S2f_S2s_ts_d2jw_dS2fdts 2377 data.calc_d2jw[data.s2s_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2s_li] = calc_diff_S2f_S2s_ts_d2jw_dS2sdts 2378 data.calc_d2jw[data.ts_li][data.ts_li] = calc_diff_S2f_S2s_ts_d2jw_dts2 2379 2380 # Diffusion as a sphere. 2381 if self.diff_data.type == 'sphere': 2382 # Gradient. 2383 data.calc_djw[0] = calc_diff_S2f_S2s_ts_djw_dDj 2384 2385 # Hessian. 2386 data.calc_d2jw[0][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdDk 2387 if self.param_set != 'diff': 2388 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2f 2389 data.calc_d2jw[0][data.s2s_li] = data.calc_d2jw[data.s2s_li][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2390 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2391 2392 # Diffusion as a spheroid. 2393 elif self.diff_data.type == 'spheroid': 2394 # Gradient. 2395 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2f_S2s_ts_djw_dDj 2396 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2f_S2s_ts_djw_dOj 2397 2398 # Hessian. 2399 data.calc_d2jw[0][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdDk 2400 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdDk 2401 data.calc_d2jw[1][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdDk 2402 2403 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdOj 2404 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdOj 2405 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdOj 2406 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdOj 2407 2408 data.calc_d2jw[2][2] = calc_diff_S2f_S2s_ts_d2jw_dOjdOk 2409 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2f_S2s_ts_d2jw_dOjdOk 2410 data.calc_d2jw[3][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdOk 2411 2412 if self.param_set != 'diff': 2413 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2f 2414 data.calc_d2jw[1][data.s2f_li] = data.calc_d2jw[data.s2f_li][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2f 2415 data.calc_d2jw[2][data.s2f_li] = data.calc_d2jw[data.s2f_li][2] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2f 2416 data.calc_d2jw[3][data.s2f_li] = data.calc_d2jw[data.s2f_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2f 2417 2418 data.calc_d2jw[0][data.s2s_li] = data.calc_d2jw[data.s2s_li][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2419 data.calc_d2jw[1][data.s2s_li] = data.calc_d2jw[data.s2s_li][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2420 data.calc_d2jw[2][data.s2s_li] = data.calc_d2jw[data.s2s_li][2] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2421 data.calc_d2jw[3][data.s2s_li] = data.calc_d2jw[data.s2s_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2422 2423 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2424 data.calc_d2jw[1][data.ts_li] = data.calc_d2jw[data.ts_li][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2425 data.calc_d2jw[2][data.ts_li] = data.calc_d2jw[data.ts_li][2] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2426 data.calc_d2jw[3][data.ts_li] = data.calc_d2jw[data.ts_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2427 2428 # Diffusion as an ellipsoid. 2429 elif self.diff_data.type == 'ellipsoid': 2430 # Gradient. 2431 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_ellipsoid_S2f_S2s_ts_djw_dDj 2432 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2f_S2s_ts_djw_dOj 2433 2434 # Hessian. 2435 data.calc_d2jw[0][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdDk 2436 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdDk 2437 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdDk 2438 data.calc_d2jw[1][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdDk 2439 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdDk 2440 data.calc_d2jw[2][2] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdDk 2441 2442 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdOj 2443 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdOj 2444 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdOj 2445 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdOj 2446 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdOj 2447 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdOj 2448 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdOj 2449 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdOj 2450 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdOj 2451 2452 data.calc_d2jw[3][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdOk 2453 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdOk 2454 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdOk 2455 data.calc_d2jw[4][4] = calc_diff_S2f_S2s_ts_d2jw_dOjdOk 2456 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2f_S2s_ts_d2jw_dOjdOk 2457 data.calc_d2jw[5][5] = calc_diff_S2f_S2s_ts_d2jw_dOjdOk 2458 2459 if self.param_set != 'diff': 2460 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdS2f 2461 data.calc_d2jw[1][data.s2f_li] = data.calc_d2jw[data.s2f_li][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdS2f 2462 data.calc_d2jw[2][data.s2f_li] = data.calc_d2jw[data.s2f_li][2] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdS2f 2463 data.calc_d2jw[3][data.s2f_li] = data.calc_d2jw[data.s2f_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2f 2464 data.calc_d2jw[4][data.s2f_li] = data.calc_d2jw[data.s2f_li][4] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2f 2465 data.calc_d2jw[5][data.s2f_li] = data.calc_d2jw[data.s2f_li][5] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2f 2466 2467 data.calc_d2jw[0][data.s2s_li] = data.calc_d2jw[data.s2s_li][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdS2s 2468 data.calc_d2jw[1][data.s2s_li] = data.calc_d2jw[data.s2s_li][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdS2s 2469 data.calc_d2jw[2][data.s2s_li] = data.calc_d2jw[data.s2s_li][2] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdS2s 2470 data.calc_d2jw[3][data.s2s_li] = data.calc_d2jw[data.s2s_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2471 data.calc_d2jw[4][data.s2s_li] = data.calc_d2jw[data.s2s_li][4] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2472 data.calc_d2jw[5][data.s2s_li] = data.calc_d2jw[data.s2s_li][5] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2473 2474 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdts 2475 data.calc_d2jw[1][data.ts_li] = data.calc_d2jw[data.ts_li][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdts 2476 data.calc_d2jw[2][data.ts_li] = data.calc_d2jw[data.ts_li][2] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdts 2477 data.calc_d2jw[3][data.ts_li] = data.calc_d2jw[data.ts_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2478 data.calc_d2jw[4][data.ts_li] = data.calc_d2jw[data.ts_li][4] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2479 data.calc_d2jw[5][data.ts_li] = data.calc_d2jw[data.ts_li][5] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2480 2481 # Diffusion parameters and model-free parameters {S2f, tf, S2s, ts}. 2482 elif data.s2f_i != None and data.tf_i != None and data.s2s_i != None and data.ts_i != None: 2483 # Equation. 2484 data.calc_jw_comps = calc_diff_S2f_tf_S2s_ts_jw_comps 2485 data.calc_jw = calc_S2f_tf_S2s_ts_jw 2486 2487 # Gradient. 2488 data.calc_djw_comps = calc_diff_S2f_tf_S2s_ts_djw_comps 2489 2490 if self.param_set != 'diff': 2491 # Gradient. 2492 data.calc_djw[data.s2f_li] = calc_diff_S2f_tf_S2s_ts_djw_dS2f 2493 data.calc_djw[data.tf_li] = calc_diff_S2f_tf_S2s_ts_djw_dtf 2494 data.calc_djw[data.s2s_li] = calc_diff_S2f_tf_S2s_ts_djw_dS2s 2495 data.calc_djw[data.ts_li] = calc_diff_S2f_tf_S2s_ts_djw_dts 2496 2497 # Hessian. 2498 data.calc_d2jw[data.s2f_li][data.s2s_li] = data.calc_d2jw[data.s2s_li][data.s2f_li] = calc_S2f_S2s_ts_d2jw_dS2fdS2s 2499 data.calc_d2jw[data.s2f_li][data.tf_li] = data.calc_d2jw[data.tf_li][data.s2f_li] = calc_diff_S2f_tf_S2s_ts_d2jw_dS2fdtf 2500 data.calc_d2jw[data.s2f_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2f_li] = calc_diff_S2f_tf_S2s_ts_d2jw_dS2fdts 2501 data.calc_d2jw[data.tf_li][data.tf_li] = calc_diff_S2f_tf_S2s_ts_d2jw_dtf2 2502 data.calc_d2jw[data.s2s_li][data.ts_li] = data.calc_d2jw[data.ts_li][data.s2s_li] = calc_diff_S2f_tf_S2s_ts_d2jw_dS2sdts 2503 data.calc_d2jw[data.ts_li][data.ts_li] = calc_diff_S2f_tf_S2s_ts_d2jw_dts2 2504 2505 # Diffusion as a sphere. 2506 if self.diff_data.type == 'sphere': 2507 # Gradient. 2508 data.calc_djw[0] = calc_diff_S2f_tf_S2s_ts_djw_dDj 2509 2510 # Hessian. 2511 data.calc_d2jw[0][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdDk 2512 if self.param_set != 'diff': 2513 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdS2f 2514 data.calc_d2jw[0][data.tf_li] = data.calc_d2jw[data.tf_li][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdtf 2515 data.calc_d2jw[0][data.s2s_li] = data.calc_d2jw[data.s2s_li][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdS2s 2516 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdts 2517 2518 # Diffusion as a spheroid. 2519 elif self.diff_data.type == 'spheroid': 2520 # Gradient. 2521 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2f_tf_S2s_ts_djw_dDj 2522 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2f_tf_S2s_ts_djw_dOj 2523 2524 # Hessian. 2525 data.calc_d2jw[0][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdDk 2526 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdDk 2527 data.calc_d2jw[1][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdDk 2528 2529 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdOj 2530 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdOj 2531 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdOj 2532 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdOj 2533 2534 data.calc_d2jw[2][2] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdOk 2535 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdOk 2536 data.calc_d2jw[3][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdOk 2537 2538 if self.param_set != 'diff': 2539 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdS2f 2540 data.calc_d2jw[1][data.s2f_li] = data.calc_d2jw[data.s2f_li][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdS2f 2541 data.calc_d2jw[2][data.s2f_li] = data.calc_d2jw[data.s2f_li][2] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdS2f 2542 data.calc_d2jw[3][data.s2f_li] = data.calc_d2jw[data.s2f_li][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdS2f 2543 2544 data.calc_d2jw[0][data.tf_li] = data.calc_d2jw[data.tf_li][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdtf 2545 data.calc_d2jw[1][data.tf_li] = data.calc_d2jw[data.tf_li][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdtf 2546 data.calc_d2jw[2][data.tf_li] = data.calc_d2jw[data.tf_li][2] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdtf 2547 data.calc_d2jw[3][data.tf_li] = data.calc_d2jw[data.tf_li][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdtf 2548 2549 data.calc_d2jw[0][data.s2s_li] = data.calc_d2jw[data.s2s_li][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2550 data.calc_d2jw[1][data.s2s_li] = data.calc_d2jw[data.s2s_li][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2551 data.calc_d2jw[2][data.s2s_li] = data.calc_d2jw[data.s2s_li][2] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2552 data.calc_d2jw[3][data.s2s_li] = data.calc_d2jw[data.s2s_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2553 2554 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2555 data.calc_d2jw[1][data.ts_li] = data.calc_d2jw[data.ts_li][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2556 data.calc_d2jw[2][data.ts_li] = data.calc_d2jw[data.ts_li][2] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2557 data.calc_d2jw[3][data.ts_li] = data.calc_d2jw[data.ts_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2558 2559 # Diffusion as an ellipsoid. 2560 elif self.diff_data.type == 'ellipsoid': 2561 # Gradient. 2562 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_ellipsoid_S2f_tf_S2s_ts_djw_dDj 2563 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2f_tf_S2s_ts_djw_dOj 2564 2565 # Hessian. 2566 data.calc_d2jw[0][0] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdDk 2567 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdDk 2568 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdDk 2569 data.calc_d2jw[1][1] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdDk 2570 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdDk 2571 data.calc_d2jw[2][2] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdDk 2572 2573 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdOj 2574 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdOj 2575 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdOj 2576 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdOj 2577 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdOj 2578 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdOj 2579 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdOj 2580 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdOj 2581 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdOj 2582 2583 data.calc_d2jw[3][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdOk 2584 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdOk 2585 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdOk 2586 data.calc_d2jw[4][4] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdOk 2587 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdOk 2588 data.calc_d2jw[5][5] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdOk 2589 2590 if self.param_set != 'diff': 2591 data.calc_d2jw[0][data.s2f_li] = data.calc_d2jw[data.s2f_li][0] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdS2f 2592 data.calc_d2jw[1][data.s2f_li] = data.calc_d2jw[data.s2f_li][1] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdS2f 2593 data.calc_d2jw[2][data.s2f_li] = data.calc_d2jw[data.s2f_li][2] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdS2f 2594 data.calc_d2jw[3][data.s2f_li] = data.calc_d2jw[data.s2f_li][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdS2f 2595 data.calc_d2jw[4][data.s2f_li] = data.calc_d2jw[data.s2f_li][4] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdS2f 2596 data.calc_d2jw[5][data.s2f_li] = data.calc_d2jw[data.s2f_li][5] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdS2f 2597 2598 data.calc_d2jw[0][data.tf_li] = data.calc_d2jw[data.tf_li][0] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdtf 2599 data.calc_d2jw[1][data.tf_li] = data.calc_d2jw[data.tf_li][1] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdtf 2600 data.calc_d2jw[2][data.tf_li] = data.calc_d2jw[data.tf_li][2] = calc_ellipsoid_S2f_tf_S2s_ts_d2jw_dDjdtf 2601 data.calc_d2jw[3][data.tf_li] = data.calc_d2jw[data.tf_li][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdtf 2602 data.calc_d2jw[4][data.tf_li] = data.calc_d2jw[data.tf_li][4] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdtf 2603 data.calc_d2jw[5][data.tf_li] = data.calc_d2jw[data.tf_li][5] = calc_diff_S2f_tf_S2s_ts_d2jw_dOjdtf 2604 2605 data.calc_d2jw[0][data.s2s_li] = data.calc_d2jw[data.s2s_li][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdS2s 2606 data.calc_d2jw[1][data.s2s_li] = data.calc_d2jw[data.s2s_li][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdS2s 2607 data.calc_d2jw[2][data.s2s_li] = data.calc_d2jw[data.s2s_li][2] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdS2s 2608 data.calc_d2jw[3][data.s2s_li] = data.calc_d2jw[data.s2s_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2609 data.calc_d2jw[4][data.s2s_li] = data.calc_d2jw[data.s2s_li][4] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2610 data.calc_d2jw[5][data.s2s_li] = data.calc_d2jw[data.s2s_li][5] = calc_diff_S2f_S2s_ts_d2jw_dOjdS2s 2611 2612 data.calc_d2jw[0][data.ts_li] = data.calc_d2jw[data.ts_li][0] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdts 2613 data.calc_d2jw[1][data.ts_li] = data.calc_d2jw[data.ts_li][1] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdts 2614 data.calc_d2jw[2][data.ts_li] = data.calc_d2jw[data.ts_li][2] = calc_ellipsoid_S2f_S2s_ts_d2jw_dDjdts 2615 data.calc_d2jw[3][data.ts_li] = data.calc_d2jw[data.ts_li][3] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2616 data.calc_d2jw[4][data.ts_li] = data.calc_d2jw[data.ts_li][4] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2617 data.calc_d2jw[5][data.ts_li] = data.calc_d2jw[data.ts_li][5] = calc_diff_S2f_S2s_ts_d2jw_dOjdts 2618 2619 # Bad parameter combination. 2620 else: 2621 print "Invalid combination of parameters for the extended model-free equation." 2622 return 0 2623 2624 # Unknown model-free equation. 2625 else: 2626 print "Unknown model-free equation." 2627 return 0 2628 2629 2630 # Initialise function data structures. 2631 ###################################### 2632 2633 # Relaxation equation components. 2634 data.create_dip_func, data.create_dip_grad, data.create_dip_hess = [], [], [] 2635 data.create_csa_func, data.create_csa_grad, data.create_csa_hess = [], [], [] 2636 data.create_rex_func, data.create_rex_grad = [], [] 2637 2638 data.create_dip_jw_func, data.create_dip_jw_grad, data.create_dip_jw_hess = [], [], [] 2639 data.create_csa_jw_func, data.create_csa_jw_grad, data.create_csa_jw_hess = [], [], [] 2640 2641 # Ri' 2642 data.create_ri_prime = None 2643 data.create_dri_prime = [] 2644 data.create_d2ri_prime = [] 2645 2646 # Ri 2647 data.create_ri, data.create_dri, data.create_d2ri = [], [], [] 2648 data.get_r1, data.get_dr1, data.get_d2r1 = [], [], [] 2649 2650 # Fill the structures with None. 2651 for i in xrange(data.num_ri): 2652 data.create_dip_func.append(None) 2653 data.create_dip_grad.append(None) 2654 data.create_dip_hess.append(None) 2655 data.create_csa_func.append(None) 2656 data.create_csa_grad.append(None) 2657 data.create_csa_hess.append(None) 2658 data.create_rex_func.append(None) 2659 data.create_rex_grad.append(None) 2660 data.create_dip_jw_func.append(None) 2661 data.create_dip_jw_grad.append(None) 2662 data.create_dip_jw_hess.append(None) 2663 data.create_csa_jw_func.append(None) 2664 data.create_csa_jw_grad.append(None) 2665 data.create_csa_jw_hess.append(None) 2666 data.create_ri.append(None) 2667 data.create_dri.append(None) 2668 data.create_d2ri.append(None) 2669 data.get_r1.append(None) 2670 data.get_dr1.append(None) 2671 data.get_d2r1.append(None) 2672 2673 2674 # Select the functions for the calculation of ri_prime, dri_prime, and d2ri_prime components. 2675 ############################################################################################# 2676 2677 for i in xrange(data.num_ri): 2678 # The R1 equations. 2679 if data.ri_labels[i] == 'R1': 2680 data.create_csa_func[i] = comp_r1_csa_const 2681 data.create_csa_grad[i] = comp_r1_csa_const 2682 data.create_csa_hess[i] = comp_r1_csa_const 2683 data.create_dip_jw_func[i] = comp_r1_dip_jw 2684 data.create_dip_jw_grad[i] = comp_r1_dip_jw 2685 data.create_dip_jw_hess[i] = comp_r1_dip_jw 2686 data.create_csa_jw_func[i] = comp_r1_csa_jw 2687 data.create_csa_jw_grad[i] = comp_r1_csa_jw 2688 data.create_csa_jw_hess[i] = comp_r1_csa_jw 2689 2690 # The R2 equations. 2691 elif data.ri_labels[i] == 'R2': 2692 data.create_dip_func[i] = comp_r2_dip_const 2693 data.create_dip_grad[i] = comp_r2_dip_const 2694 data.create_dip_hess[i] = comp_r2_dip_const 2695 data.create_csa_func[i] = comp_r2_csa_const 2696 data.create_csa_grad[i] = comp_r2_csa_const 2697 data.create_csa_hess[i] = comp_r2_csa_const 2698 data.create_rex_func[i] = comp_rex_const_func 2699 data.create_rex_grad[i] = comp_rex_const_grad 2700 data.create_dip_jw_func[i] = comp_r2_dip_jw 2701 data.create_dip_jw_grad[i] = comp_r2_dip_jw 2702 data.create_dip_jw_hess[i] = comp_r2_dip_jw 2703 data.create_csa_jw_func[i] = comp_r2_csa_jw 2704 data.create_csa_jw_grad[i] = comp_r2_csa_jw 2705 data.create_csa_jw_hess[i] = comp_r2_csa_jw 2706 2707 # The NOE equations. 2708 elif data.ri_labels[i] == 'NOE': 2709 data.create_dip_jw_func[i] = comp_sigma_noe_dip_jw 2710 data.create_dip_jw_grad[i] = comp_sigma_noe_dip_jw 2711 data.create_dip_jw_hess[i] = comp_sigma_noe_dip_jw 2712 data.create_ri[i] = calc_noe 2713 data.create_dri[i] = calc_dnoe 2714 data.create_d2ri[i] = calc_d2noe 2715 if data.noe_r1_table[i] == None: 2716 data.get_r1[i] = calc_r1 2717 data.get_dr1[i] = calc_dr1 2718 data.get_d2r1[i] = calc_d2r1 2719 else: 2720 data.get_r1[i] = extract_r1 2721 data.get_dr1[i] = extract_dr1 2722 data.get_d2r1[i] = extract_d2r1 2723 2724 2725 # Select the functions for the calculation of ri_prime, dri_prime, and d2ri_prime. 2726 ################################################################################## 2727 2728 # ri_prime. 2729 if data.rex_i == None: 2730 data.create_ri_prime = func_ri_prime 2731 else: 2732 data.create_ri_prime = func_ri_prime_rex 2733 2734 # dri_prime and d2ri_prime. 2735 for i in xrange(data.total_num_params): 2736 # Diffusion tensor parameters are the only parameters. 2737 if self.param_set == 'diff': 2738 # Gradient. 2739 data.create_dri_prime.append(func_dri_djw_prime) 2740 2741 # Hessian. 2742 data.create_d2ri_prime.append([]) 2743 for j in xrange(data.total_num_params): 2744 data.create_d2ri_prime[i].append(func_d2ri_djwidjwj_prime) 2745 2746 # Skip to the next parameter index. 2747 continue 2748 2749 # Residue specific parameter index. 2750 index = i - num_diff_params 2751 if index < 0: 2752 index = None 2753 2754 # Rex. 2755 if index != None and data.param_types[index] == 'Rex': 2756 # Gradient. 2757 data.create_dri_prime.append(func_dri_drex_prime) 2758 2759 # Hessian. 2760 data.create_d2ri_prime.append([]) 2761 for j in xrange(data.total_num_params): 2762 # Residue specific parameter index. 2763 index2 = j - num_diff_params 2764 if index2 < 0: 2765 index2 = None 2766 2767 # Rex. 2768 if index2 != None and data.param_types[index2] == 'Rex': 2769 data.create_d2ri_prime[i].append(None) 2770 2771 # Bond length. 2772 elif index2 != None and data.param_types[index2] == 'r': 2773 data.create_d2ri_prime[i].append(None) 2774 2775 # CSA. 2776 elif index2 != None and data.param_types[index2] == 'CSA': 2777 data.create_d2ri_prime[i].append(None) 2778 2779 # Any other parameter. 2780 else: 2781 data.create_d2ri_prime[i].append(None) 2782 2783 # Bond length. 2784 elif index != None and data.param_types[index] == 'r': 2785 # Gradient. 2786 data.create_dri_prime.append(func_dri_dr_prime) 2787 2788 # Hessian. 2789 data.create_d2ri_prime.append([]) 2790 for j in xrange(data.total_num_params): 2791 # Residue specific parameter index. 2792 index2 = j - num_diff_params 2793 if index2 < 0: 2794 index2 = None 2795 2796 # Rex. 2797 if index2 != None and data.param_types[index2] == 'Rex': 2798 data.create_d2ri_prime[i].append(None) 2799 2800 # Bond length. 2801 elif index2 != None and data.param_types[index2] == 'r': 2802 data.create_d2ri_prime[i].append(func_d2ri_dr2_prime) 2803 2804 # CSA. 2805 elif index2 != None and data.param_types[index2] == 'CSA': 2806 data.create_d2ri_prime[i].append(None) 2807 2808 # Any other parameter. 2809 else: 2810 data.create_d2ri_prime[i].append(func_d2ri_djwdr_prime) 2811 2812 # CSA. 2813 elif index != None and data.param_types[index] == 'CSA': 2814 # Gradient. 2815 data.create_dri_prime.append(func_dri_dcsa_prime) 2816 2817 # Hessian. 2818 data.create_d2ri_prime.append([]) 2819 for j in xrange(data.total_num_params): 2820 # Residue specific parameter index. 2821 index2 = j - num_diff_params 2822 if index2 < 0: 2823 index2 = None 2824 2825 # Rex. 2826 if index2 != None and data.param_types[index2] == 'Rex': 2827 data.create_d2ri_prime[i].append(None) 2828 2829 # Bond length. 2830 elif index2 != None and data.param_types[index2] == 'r': 2831 data.create_d2ri_prime[i].append(None) 2832 2833 # CSA. 2834 elif index2 != None and data.param_types[index2] == 'CSA': 2835 data.create_d2ri_prime[i].append(func_d2ri_dcsa2_prime) 2836 2837 # Any other parameter. 2838 else: 2839 data.create_d2ri_prime[i].append(func_d2ri_djwdcsa_prime) 2840 2841 # Any other parameter. 2842 else: 2843 # Gradient. 2844 data.create_dri_prime.append(func_dri_djw_prime) 2845 2846 # Hessian. 2847 data.create_d2ri_prime.append([]) 2848 for j in xrange(data.total_num_params): 2849 # Residue specific parameter index. 2850 index2 = j - num_diff_params 2851 if index2 < 0: 2852 index2 = None 2853 2854 # Rex. 2855 if index2 != None and data.param_types[index2] == 'Rex': 2856 data.create_d2ri_prime[i].append(None) 2857 2858 # Bond length. 2859 elif index2 != None and data.param_types[index2] == 'r': 2860 data.create_d2ri_prime[i].append(func_d2ri_djwdr_prime) 2861 2862 # CSA. 2863 elif index2 != None and data.param_types[index2] == 'CSA': 2864 data.create_d2ri_prime[i].append(func_d2ri_djwdcsa_prime) 2865 2866 # Any other parameter. 2867 else: 2868 data.create_d2ri_prime[i].append(func_d2ri_djwidjwj_prime) 2869 2870 2871 # Both the bond length and CSA are fixed {}. 2872 ############################################ 2873 2874 if data.r_i == None and data.csa_i == None: 2875 # The main ri component functions 2876 if data.rex_i == None: 2877 data.create_ri_comps = ri_comps 2878 data.create_dri_comps = dri_comps 2879 data.create_d2ri_comps = d2ri_comps 2880 else: 2881 data.create_ri_comps = ri_comps_rex 2882 data.create_dri_comps = dri_comps_rex 2883 data.create_d2ri_comps = d2ri_comps 2884 2885 # Calculate the dipolar and CSA constant components. 2886 comp_dip_const_func(data, data.bond_length) 2887 comp_csa_const_func(data, data.csa) 2888 for i in xrange(data.num_ri): 2889 data.dip_comps_func[i] = data.dip_const_func 2890 if data.create_dip_func[i]: 2891 data.dip_comps_func[i] = data.create_dip_func[i](data.dip_const_func) 2892 if data.create_csa_func[i]: 2893 data.csa_comps_func[i] = data.create_csa_func[i](data.csa_const_func[data.remap_table[i]]) 2894 2895 2896 # The bond length is a parameter {r}. 2897 ##################################### 2898 2899 elif data.r_i != None and data.csa_i == None: 2900 # The main ri component functions 2901 if data.rex_i == None: 2902 data.create_ri_comps = ri_comps_r 2903 data.create_dri_comps = dri_comps_r 2904 data.create_d2ri_comps = d2ri_comps_r 2905 else: 2906 data.create_ri_comps = ri_comps_r_rex 2907 data.create_dri_comps = dri_comps_r_rex 2908 data.create_d2ri_comps = d2ri_comps_r 2909 2910 # Calculate the CSA constant. 2911 comp_csa_const_func(data, data.csa) 2912 for i in xrange(data.num_ri): 2913 if data.create_csa_func[i]: 2914 data.csa_comps_func[i] = data.create_csa_func[i](data.csa_const_func[data.remap_table[i]]) 2915 2916 2917 # The CSA is a parameter {CSA}. 2918 ############################### 2919 2920 elif data.r_i == None and data.csa_i != None: 2921 # The main ri component functions 2922 if data.rex_i == None: 2923 data.create_ri_comps = ri_comps_csa 2924 data.create_dri_comps = dri_comps_csa 2925 data.create_d2ri_comps = d2ri_comps_csa 2926 else: 2927 data.create_ri_comps = ri_comps_csa_rex 2928 data.create_dri_comps = dri_comps_csa_rex 2929 data.create_d2ri_comps = d2ri_comps_csa 2930 2931 # Calculate the dipolar constant. 2932 comp_dip_const_func(data, data.bond_length) 2933 for i in xrange(data.num_ri): 2934 data.dip_comps_func[i] = data.dip_const_func 2935 if data.create_dip_func[i]: 2936 data.dip_comps_func[i] = data.create_dip_func[i](data.dip_const_func) 2937 2938 2939 # Both the bond length and CSA are parameters {r, CSA}. 2940 ####################################################### 2941 2942 elif data.r_i != None and data.csa_i != None: 2943 # The main ri component functions 2944 if data.rex_i == None: 2945 data.create_ri_comps = ri_comps_r_csa 2946 data.create_dri_comps = dri_comps_r_csa 2947 data.create_d2ri_comps = d2ri_comps_r_csa 2948 else: 2949 data.create_ri_comps = ri_comps_r_csa_rex 2950 data.create_dri_comps = dri_comps_r_csa_rex 2951 data.create_d2ri_comps = d2ri_comps_r_csa 2952 2953 2954 # Invalid combination of parameters. 2955 #################################### 2956 2957 else: 2958 print "Invalid combination of parameters for the model-free equations." 2959 return 0 2960 2961 return 1
2962 2963
2964 -class Data:
2965 - def __init__(self):
2966 """Empty container for storing data."""
2967