maths_fns.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 'iso', 'axial', or 'aniso'. 62 63 diff_params: An array with the diffusion parameters. 64 65 scaling_matrix: A diagonal matrix of scaling factors. 66 67 68 69 Additional layer of equations to simplify the relaxation equations, gradients, and Hessians. 70 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 71 72 The R1 and R2 equations are left alone, while the NOE is calculated from the R1 and 73 sigma_noe values. 74 75 76 The relaxation equations 77 ~~~~~~~~~~~~~~~~~~~~~~~~ 78 79 Data structure: data.ri 80 Dimension: 1D, (relaxation data) 81 Type: Numeric array, Float64 82 Dependencies: data.ri_prime 83 Required by: data.chi2, data.dchi2, data.d2chi2 84 85 86 R1() = R1'() 87 88 89 R2() = R2'() 90 91 gH sigma_noe() 92 NOE() = 1 + -- . ----------- 93 gN R1() 94 95 96 97 The relaxation gradients 98 ~~~~~~~~~~~~~~~~~~~~~~~~ 99 100 Data structure: data.dri 101 Dimension: 2D, (parameters, relaxation data) 102 Type: Numeric array, Float64 103 Dependencies: data.ri_prime, data.dri_prime 104 Required by: data.dchi2, data.d2chi2 105 106 107 dR1() dR1'() 108 ------- = ------- 109 dthetaj dthetaj 110 111 112 dR2() dR2'() 113 ------- = ------- 114 dthetaj dthetaj 115 116 117 dNOE() gH 1 / dsigma_noe() dR1() \ 118 ------- = -- . ------- . | R1() . ------------ - sigma_noe() . ------- | 119 dthetaj gN R1()**2 \ dthetaj dthetaj / 120 121 122 123 The relaxation Hessians 124 ~~~~~~~~~~~~~~~~~~~~~~~ 125 126 Data structure: data.d2ri 127 Dimension: 3D, (parameters, parameters, relaxation data) 128 Type: Numeric array, Float64 129 Dependencies: data.ri_prime, data.dri_prime, data.d2ri_prime 130 Required by: data.d2chi2 131 132 133 d2R1() d2R1'() 134 --------------- = --------------- 135 dthetai.dthetaj dthetai.dthetaj 136 137 138 d2R2() d2R2'() 139 --------------- = --------------- 140 dthetai.dthetaj dthetai.dthetaj 141 142 143 d2NOE() gH 1 / / dR1() dR1() d2R1() \ 144 --------------- = -- . ------- . | sigma_noe() . | 2 . ------- . ------- - R1() . --------------- | 145 dthetai.dthetaj gN R1()**3 \ \ dthetai dthetaj dthetai.dthetaj / 146 147 / dsigma_noe() dR1() dR1() dsigma_noe() d2sigma_noe() \ \ 148 - R1() . | ------------ . ------- + ------- . ------------ - R1() . --------------- | | 149 \ dthetai dthetaj dthetai dthetaj dthetai.dthetaj / / 150 151 152 153 The chi-sqared equation 154 ~~~~~~~~~~~~~~~~~~~~~~~ 155 _n_ 156 \ (Ri - Ri()) ** 2 157 Chi2 = > ---------------- 158 /__ sigma_i ** 2 159 i=1 160 161 where: 162 Ri are the values of the measured relaxation data set. 163 Ri() are the values of the back calculated relaxation data set. 164 sigma_i are the values of the error set. 165 166 167 The chi-sqared gradient 168 ~~~~~~~~~~~~~~~~~~~~~~~ 169 _n_ 170 dChi2 \ / Ri - Ri() dRi() \ 171 ------- = -2 > | ---------- . ------- | 172 dthetaj /__ \ sigma_i**2 dthetaj / 173 i=1 174 175 where: 176 Ri are the values of the measured relaxation data set. 177 Ri() are the values of the back calculated relaxation data set. 178 sigma_i are the values of the error set. 179 180 181 The chi-sqared Hessian 182 ~~~~~~~~~~~~~~~~~~~~~~ 183 _n_ 184 d2chi2 \ 1 / dRi() dRi() d2Ri() \ 185 --------------- = 2 > ---------- | ------- . ------- - (Ri - Ri()) . --------------- | 186 dthetaj.dthetak /__ sigma_i**2 \ dthetaj dthetak dthetaj.dthetak / 187 i=1 188 189 where: 190 Ri are the values of the measured relaxation data set. 191 Ri() are the values of the back calculated relaxation data set. 192 sigma_i are the values of the error set. 193 """ 194 195 # Arguments. 196 self.param_set = param_set 197 self.total_num_params = len(init_params) 198 self.scaling_matrix = scaling_matrix 199 self.num_res = num_res 200 self.params = 1.0 * init_params 201 202 # Data structures for tests set to some random array (in this case all pi). 203 self.func_test = pi * ones(self.total_num_params, Float64) 204 self.grad_test = pi * ones(self.total_num_params, Float64) 205 self.hess_test = pi * ones(self.total_num_params, Float64) 206 207 # Initialise the data class for storing diffusion tensor data. 208 self.diff_data = Data() 209 self.diff_data.type = diff_type 210 self.diff_data.params = diff_params 211 self.init_diff_data(self.diff_data) 212 213 # Total number of ri. 214 self.total_num_ri = 0 215 216 # Set the function for packaging diffusion tensor parameters. 217 if self.diff_data.type == 'iso': 218 self.param_index = 1 219 self.diff_end_index = 1 220 elif self.diff_data.type == 'axial': 221 self.param_index = 4 222 self.diff_end_index = 4 223 elif self.diff_data.type == 'aniso': 224 self.param_index = 6 225 self.diff_end_index = 6 226 if self.param_set != 'all': 227 self.param_index = 0 228 229 # Create the data array used to store data. 230 self.data = [] 231 for i in xrange(self.num_res): 232 # Total number of ri. 233 self.total_num_ri = self.total_num_ri + num_ri[i] 234 235 # Append the class instance Data to the data array. 236 self.data.append(Data()) 237 238 # Number of indecies. 239 self.data[i].num_indecies = self.diff_data.num_indecies 240 241 # Calculate the five frequencies per field strength which cause R1, R2, and NOE relaxation. 242 self.data[i].frq_list = zeros((num_frq[i], 5), Float64) 243 self.data[i].frq_list_ext = zeros((num_frq[i], 5, self.diff_data.num_indecies), Float64) 244 self.data[i].frq_sqrd_list_ext = zeros((num_frq[i], 5, self.diff_data.num_indecies), Float64) 245 for j in xrange(num_frq[i]): 246 frqH = 2.0 * pi * frq[i][j] 247 frqX = frqH / g_ratio 248 self.data[i].frq_list[j, 1] = frqX 249 self.data[i].frq_list[j, 2] = frqH - frqX 250 self.data[i].frq_list[j, 3] = frqH 251 self.data[i].frq_list[j, 4] = frqH + frqX 252 self.data[i].frq_sqrd_list = self.data[i].frq_list ** 2 253 for j in xrange(self.diff_data.num_indecies): 254 self.data[i].frq_list_ext[:, :, j] = self.data[i].frq_list 255 self.data[i].frq_sqrd_list_ext[:, :, j] = self.data[i].frq_sqrd_list 256 257 # Store supplied data in self.data 258 self.data[i].gh = gh 259 self.data[i].gx = gx 260 self.data[i].g_ratio = g_ratio 261 self.data[i].h_bar = h_bar 262 self.data[i].mu0 = mu0 263 self.data[i].equations = equations[i] 264 self.data[i].param_types = param_types[i] 265 self.data[i].relax_data = relax_data[i] 266 self.data[i].errors = errors[i] 267 self.data[i].bond_length = bond_length[i] 268 self.data[i].csa = csa[i] 269 self.data[i].num_ri = num_ri[i] 270 self.data[i].num_frq = num_frq[i] 271 self.data[i].frq = frq[i] 272 self.data[i].remap_table = remap_table[i] 273 self.data[i].noe_r1_table = noe_r1_table[i] 274 self.data[i].ri_labels = ri_labels[i] 275 self.data[i].num_params = num_params[i] 276 self.data[i].xh_unit_vector = vectors[i] 277 278 # Parameter values for minimising soley the diffusion tensor parameters. 279 if self.param_set == 'diff': 280 self.data[i].param_values = param_values[i] 281 282 # Indecies for constructing the global generic model-free gradient and Hessian kite. 283 if i == 0: 284 self.data[i].start_index = self.diff_data.num_params 285 else: 286 self.data[i].start_index = self.data[i-1].end_index 287 self.data[i].end_index = self.data[i].start_index + self.data[i].num_params 288 289 # Total number of parameters. 290 if self.param_set == 'mf' or self.param_set == 'local_tm': 291 self.data[i].total_num_params = self.data[i].num_params 292 elif self.param_set == 'diff': 293 self.data[i].total_num_params = self.diff_data.num_params 294 else: 295 self.data[i].total_num_params = self.data[i].num_params + self.diff_data.num_params 296 297 # Initialise the residue specific data. 298 self.init_res_data(self.data[i], self.diff_data) 299 300 # Setup the residue specific equations. 301 if not self.setup_equations(self.data[i]): 302 raise RelaxError, "The model-free equations could not be setup." 303 304 # Diffusion tensor parameters. 305 if self.param_set == 'local_tm': 306 self.diff_data.params = self.params[0:1] 307 elif self.param_set == 'diff' or self.param_set == 'all': 308 self.diff_data.params = self.params[0:self.diff_end_index] 309 310 # Calculate the correlation times ti. 311 self.diff_data.calc_ti(self.data[i], self.diff_data) 312 313 # ti spectral density components. 314 self.data[i].w_ti_sqrd = self.data[i].frq_sqrd_list_ext * self.data[i].ti ** 2 315 self.data[i].fact_ti = 1.0 / (1.0 + self.data[i].w_ti_sqrd) 316 317 # 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. 318 missing_r1 = 0 319 for j in xrange(self.data[i].num_ri): 320 if self.data[i].ri_labels[j] == 'NOE' and self.data[i].noe_r1_table[j] == None: 321 missing_r1 = 1 322 if missing_r1: 323 self.init_res_r1_data(self.data[i]) 324 325 # Scaling initialisation. 326 if self.scaling_matrix: 327 self.scaling_flag = 1 328 else: 329 self.scaling_flag = 0 330 331 # Initialise the total chi2 value, gradient, and Hessian data structures. 332 self.total_chi2 = 0.0 333 self.total_dchi2 = zeros((self.total_num_params), Float64) 334 self.total_d2chi2 = zeros((self.total_num_params, self.total_num_params), Float64) 335 336 # Initialise the total ri gradient data structure (for Levenberg-Marquardt minimisation). 337 self.total_dri = zeros((self.total_num_params, self.total_num_ri), Float64) 338 339 # Set the functions self.func, self.dfunc, and self.d2func. 340 ########################################################### 341 342 # Functions for minimising model-free parameters for a single residue. 343 if self.param_set == 'mf': 344 self.func = self.func_mf 345 self.dfunc = self.dfunc_mf 346 self.d2func = self.d2func_mf 347 348 # Functions for minimising model-free parameters for a single residue with a local tm. 349 elif self.param_set == 'local_tm': 350 self.func = self.func_local_tm 351 self.dfunc = self.dfunc_local_tm 352 self.d2func = self.d2func_local_tm 353 354 # Functions for minimising diffusion tensor parameters with all model-free parameters fixed. 355 elif self.param_set == 'diff': 356 self.func = self.func_diff 357 self.dfunc = self.dfunc_diff 358 self.d2func = self.d2func_diff 359 360 # Functions for minimising diffusion tensor parameters together with all model-free parameters. 361 elif self.param_set == 'all': 362 self.func = self.func_all 363 self.dfunc = self.dfunc_all 364 self.d2func = self.d2func_all

365 366

367 - def func_mf(self, params):

368 """Function for calculating the chi-squared value. 369 370 Used in the minimisation of model-free parameters for a single residue. 371 """ 372 373 # Store the parameter values in self.func_test for testing. 374 self.func_test = params * 1.0 375 376 # Set self.data[0] to data. 377 data = self.data[0] 378 379 # Scaling. 380 if self.scaling_flag: 381 params = matrixmultiply(params, self.scaling_matrix) 382 383 # Diffusion tensor geometry calculations. 384 if self.diff_data.calc_geom: 385 self.diff_data.calc_geom(data, self.diff_data) 386 387 # Diffusion tensor weight calculations. 388 self.diff_data.calc_ci(data, self.diff_data) 389 390 # Diffusion tensor correlation times. 391 self.diff_data.calc_ti(data, self.diff_data) 392 393 # Calculate the components of the spectral densities. 394 if data.calc_jw_comps: 395 data.calc_jw_comps(data, params) 396 397 # Calculate the spectral density values. 398 data.jw = data.calc_jw(data, params) 399 400 # Calculate the relaxation formula components. 401 data.create_ri_comps(data, params) 402 403 # Calculate the R1, R2, and sigma_noe values. 404 data.ri_prime = data.create_ri_prime(data) 405 406 # Calculate the NOE values. 407 data.ri = data.ri_prime * 1.0 408 for m in xrange(data.num_ri): 409 if data.create_ri[m]: 410 data.create_ri[m](data, m, data.remap_table[m], data.get_r1, params) 411 412 # Calculate the chi-squared value. 413 data.chi2 = chi2(data.relax_data, data.ri, data.errors) 414 415 return data.chi2

416 417

418 - def func_local_tm(self, params):

419 """Function for calculating the chi-squared value. 420 421 Used in the minimisation of model-free parameters for a single residue with a local tm. 422 """ 423 424 # Store the parameter values in self.func_test for testing. 425 self.func_test = params * 1.0 426 427 # Set self.data[0] to data. 428 data = self.data[0] 429 430 # Scaling. 431 if self.scaling_flag: 432 params = matrixmultiply(params, self.scaling_matrix) 433 434 # Diffusion tensor parameters. 435 self.diff_data.params = params[0:1] 436 437 # Diffusion tensor weight calculations. 438 self.diff_data.calc_ci(data, self.diff_data) 439 440 # Diffusion tensor correlation times. 441 self.diff_data.calc_ti(data, self.diff_data) 442 443 # ti spectral density components. 444 data.w_ti_sqrd = data.frq_sqrd_list_ext * data.ti ** 2 445 data.fact_ti = 1.0 / (1.0 + data.w_ti_sqrd) 446 447 # Calculate the components of the spectral densities. 448 if data.calc_jw_comps: 449 data.calc_jw_comps(data, params) 450 451 # Calculate the spectral density values. 452 data.jw = data.calc_jw(data, params) 453 454 # Calculate the relaxation formula components. 455 data.create_ri_comps(data, params) 456 457 # Calculate the R1, R2, and sigma_noe values. 458 data.ri_prime = data.create_ri_prime(data) 459 460 # Calculate the NOE values. 461 data.ri = data.ri_prime * 1.0 462 for m in xrange(data.num_ri): 463 if data.create_ri[m]: 464 data.create_ri[m](data, m, data.remap_table[m], data.get_r1, params) 465 466 # Calculate the chi-squared value. 467 data.chi2 = chi2(data.relax_data, data.ri, data.errors) 468 469 return data.chi2

470 471

472 - def func_diff(self, params):

473 """Function for calculating the chi-squared value. 474 475 Used in the minimisation of diffusion tensor parameters with all model-free parameters 476 fixed. 477 """ 478 479 # Store the parameter values in self.func_test for testing. 480 self.func_test = params * 1.0 481 482 # Scaling. 483 if self.scaling_flag: 484 params = matrixmultiply(params, self.scaling_matrix) 485 486 # Diffusion tensor parameters. 487 self.diff_data.params = params[0:self.diff_end_index] 488 489 # Set the total chi2 to zero. 490 self.total_chi2 = 0.0 491 492 # Loop over the residues. 493 for i in xrange(self.num_res): 494 # Set self.data[i] to data. 495 data = self.data[i] 496 497 # Diffusion tensor geometry calculations. 498 if self.diff_data.calc_geom: 499 self.diff_data.calc_geom(data, self.diff_data) 500 501 # Diffusion tensor weight calculations. 502 self.diff_data.calc_ci(data, self.diff_data) 503 504 # Diffusion tensor correlation times. 505 self.diff_data.calc_ti(data, self.diff_data) 506 507 # ti spectral density components. 508 data.w_ti_sqrd = data.frq_sqrd_list_ext * data.ti ** 2 509 data.fact_ti = 1.0 / (1.0 + data.w_ti_sqrd) 510 511 # Calculate the components of the spectral densities. 512 if data.calc_jw_comps: 513 data.calc_jw_comps(data, data.param_values) 514 515 # Calculate the spectral density values. 516 data.jw = data.calc_jw(data, data.param_values) 517 518 # Calculate the relaxation formula components. 519 data.create_ri_comps(data, data.param_values) 520 521 # Calculate the R1, R2, and sigma_noe values. 522 data.ri_prime = data.create_ri_prime(data) 523 524 # Calculate the NOE values. 525 data.ri = data.ri_prime * 1.0 526 for m in xrange(data.num_ri): 527 if data.create_ri[m]: 528 data.create_ri[m](data, m, data.remap_table[m], data.get_r1, data.param_values) 529 530 # Calculate the chi-squared value. 531 data.chi2 = chi2(data.relax_data, data.ri, data.errors) 532 533 # Add the residue specific chi2 to the total chi2. 534 self.total_chi2 = self.total_chi2 + data.chi2 535 536 return self.total_chi2

537 538

539 - def func_all(self, params):

540 """Function for calculating the chi-squared value. 541 542 Used in the minimisation of diffusion tensor parameters together with all model-free 543 parameters. 544 """ 545 546 # Store the parameter values in self.func_test for testing. 547 self.func_test = params * 1.0 548 549 # Scaling. 550 if self.scaling_flag: 551 params = matrixmultiply(params, self.scaling_matrix) 552 553 # Diffusion tensor parameters. 554 self.diff_data.params = params[0:self.diff_end_index] 555 556 # Set the total chi2 to zero. 557 self.total_chi2 = 0.0 558 559 # Loop over the residues. 560 for i in xrange(self.num_res): 561 # Set self.data[i] to data. 562 data = self.data[i] 563 564 # Diffusion tensor geometry calculations. 565 if self.diff_data.calc_geom: 566 self.diff_data.calc_geom(data, self.diff_data) 567 568 # Diffusion tensor weight calculations. 569 self.diff_data.calc_ci(data, self.diff_data) 570 571 # Diffusion tensor correlation times. 572 self.diff_data.calc_ti(data, self.diff_data) 573 574 # ti spectral density components. 575 data.w_ti_sqrd = data.frq_sqrd_list_ext * data.ti ** 2 576 data.fact_ti = 1.0 / (1.0 + data.w_ti_sqrd) 577 578 # Calculate the components of the spectral densities. 579 if data.calc_jw_comps: 580 data.calc_jw_comps(data, params) 581 582 # Calculate the spectral density values. 583 data.jw = data.calc_jw(data, params) 584 585 # Calculate the relaxation formula components. 586 data.create_ri_comps(data, params) 587 588 # Calculate the R1, R2, and sigma_noe values. 589 data.ri_prime = data.create_ri_prime(data) 590 591 # Calculate the NOE values. 592 data.ri = data.ri_prime * 1.0 593 for m in xrange(data.num_ri): 594 if data.create_ri[m]: 595 data.create_ri[m](data, m, data.remap_table[m], data.get_r1, params) 596 597 # Calculate the chi-squared value. 598 data.chi2 = chi2(data.relax_data, data.ri, data.errors) 599 600 # Add the residue specific chi2 to the total chi2. 601 self.total_chi2 = self.total_chi2 + data.chi2 602 603 return self.total_chi2

604 605

606 - def dfunc_mf(self, params):

607 """Function for calculating the chi-squared gradient. 608 609 Used in the minimisation of model-free parameters for a single residue. 610 """ 611 612 # Test if the function has already been called, otherwise run self.func. 613 if sum(params == self.func_test) != self.total_num_params: 614 self.func(params) 615 616 # Store the parameter values in self.grad_test for testing. 617 self.grad_test = params * 1.0 618 619 # Set self.data[0] to data. 620 data = self.data[0] 621 622 # Scaling. 623 if self.scaling_flag: 624 params = matrixmultiply(params, self.scaling_matrix) 625 626 # Calculate the spectral density gradient components. 627 if data.calc_djw_comps: 628 data.calc_djw_comps(data, params) 629 630 # Loop over the gradient. 631 for j in xrange(data.total_num_params): 632 # Calculate the spectral density gradients. 633 if data.calc_djw[j]: 634 data.djw = data.calc_djw[j](data, params, j) 635 else: 636 data.djw = data.djw * 0.0 637 638 # Calculate the relaxation gradient components. 639 data.create_dri_comps(data, params) 640 641 # Calculate the R1, R2, and sigma_noe gradients. 642 data.dri_prime[j] = data.create_dri_prime[j](data) 643 644 # Loop over the relaxation values and modify the NOE gradients. 645 data.dri[j] = data.dri_prime[j] 646 for m in xrange(data.num_ri): 647 if data.create_dri[m]: 648 data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j) 649 650 # Calculate the chi-squared gradient. 651 data.dchi2[j] = dchi2(data.relax_data, data.ri, data.dri[j], data.errors) 652 653 # Diagonal scaling. 654 if self.scaling_flag: 655 data.dchi2 = matrixmultiply(data.dchi2, self.scaling_matrix) 656 657 # Return a copy of the gradient. 658 return data.dchi2 * 1.0

659 660

661 - def dfunc_local_tm(self, params):

662 """Function for calculating the chi-squared gradient. 663 664 Used in the minimisation of model-free parameters for a single residue with a local tm. 665 """ 666 667 # Test if the function has already been called, otherwise run self.func. 668 if sum(params == self.func_test) != self.total_num_params: 669 self.func(params) 670 671 # Store the parameter values in self.grad_test for testing. 672 self.grad_test = params * 1.0 673 674 # Set self.data[0] to data. 675 data = self.data[0] 676 677 # Scaling. 678 if self.scaling_flag: 679 params = matrixmultiply(params, self.scaling_matrix) 680 681 # Diffusion tensor parameters. 682 self.diff_data.params = params[0:1] 683 684 # Calculate the spectral density gradient components. 685 if data.calc_djw_comps: 686 data.calc_djw_comps(data, params) 687 688 # Diffusion tensor correlation times. 689 self.diff_data.calc_dti(data, self.diff_data) 690 691 # Loop over the gradient. 692 for j in xrange(data.total_num_params): 693 # Calculate the spectral density gradients. 694 if data.calc_djw[j]: 695 data.djw = data.calc_djw[j](data, params, j) 696 else: 697 data.djw = data.djw * 0.0 698 699 # Calculate the relaxation gradient components. 700 data.create_dri_comps(data, params) 701 702 # Calculate the R1, R2, and sigma_noe gradients. 703 data.dri_prime[j] = data.create_dri_prime[j](data) 704 705 # Loop over the relaxation values and modify the NOE gradients. 706 data.dri[j] = data.dri_prime[j] 707 for m in xrange(data.num_ri): 708 if data.create_dri[m]: 709 data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j) 710 711 # Calculate the chi-squared gradient. 712 data.dchi2[j] = dchi2(data.relax_data, data.ri, data.dri[j], data.errors) 713 714 # Diagonal scaling. 715 if self.scaling_flag: 716 data.dchi2 = matrixmultiply(data.dchi2, self.scaling_matrix) 717 718 # Return a copy of the gradient. 719 return data.dchi2 * 1.0

720 721

722 - def dfunc_diff(self, params):

723 """Function for calculating the chi-squared gradient. 724 725 Used in the minimisation of diffusion tensor parameters with all model-free parameters 726 fixed. 727 """ 728 729 # Test if the function has already been called, otherwise run self.func. 730 if sum(params == self.func_test) != self.total_num_params: 731 self.func(params) 732 733 # Store the parameter values in self.grad_test for testing. 734 self.grad_test = params * 1.0 735 736 # Scaling. 737 if self.scaling_flag: 738 params = matrixmultiply(params, self.scaling_matrix) 739 740 # Diffusion tensor parameters. 741 self.diff_data.params = params[0:self.diff_end_index] 742 743 # Set the total chi2 gradient to zero. 744 self.total_dchi2 = self.total_dchi2 * 0.0 745 746 # Loop over the residues. 747 for i in xrange(self.num_res): 748 749 # Set self.data[i] to data. 750 data = self.data[i] 751 752 # Diffusion tensor geometry calculations. 753 if self.diff_data.calc_dgeom: 754 self.diff_data.calc_dgeom(data, self.diff_data) 755 756 # Diffusion tensor weight calculations. 757 if self.diff_data.calc_dci: 758 self.diff_data.calc_dci(data, self.diff_data) 759 760 # Diffusion tensor correlation times. 761 self.diff_data.calc_dti(data, self.diff_data) 762 763 # Calculate the spectral density gradient components. 764 if data.calc_djw_comps: 765 data.calc_djw_comps(data, data.param_values) 766 767 # Loop over the gradient. 768 for j in xrange(data.total_num_params): 769 # Calculate the spectral density gradients. 770 if data.calc_djw[j]: 771 data.djw = data.calc_djw[j](data, data.param_values, j) 772 else: 773 data.djw = data.djw * 0.0 774 775 # Calculate the relaxation gradient components. 776 data.create_dri_comps(data, data.param_values) 777 778 # Calculate the R1, R2, and sigma_noe gradients. 779 data.dri_prime[j] = data.create_dri_prime[j](data) 780 781 # Loop over the relaxation values and modify the NOE gradients. 782 data.dri[j] = data.dri_prime[j] 783 for m in xrange(data.num_ri): 784 if data.create_dri[m]: 785 data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j) 786 787 # Calculate the chi-squared gradient. 788 data.dchi2[j] = dchi2(data.relax_data, data.ri, data.dri[j], data.errors) 789 790 # Index for the construction of the global generic model-free gradient. 791 index = self.diff_data.num_params 792 793 # Diffusion parameter part of the global generic model-free gradient. 794 self.total_dchi2[0:index] = self.total_dchi2[0:index] + data.dchi2[0:index] 795 796 # Diagonal scaling. 797 if self.scaling_flag: 798 self.total_dchi2 = matrixmultiply(self.total_dchi2, self.scaling_matrix) 799 800 # Return a copy of the gradient. 801 return self.total_dchi2 * 1.0

802 803

804 - def dfunc_all(self, params):

805 """Function for calculating the chi-squared gradient. 806 807 Used in the minimisation of diffusion tensor parameters together with all model-free 808 parameters. 809 """ 810 811 # Test if the function has already been called, otherwise run self.func. 812 if sum(params == self.func_test) != self.total_num_params: 813 self.func(params) 814 815 # Store the parameter values in self.grad_test for testing. 816 self.grad_test = params * 1.0 817 818 # Scaling. 819 if self.scaling_flag: 820 params = matrixmultiply(params, self.scaling_matrix) 821 822 # Diffusion tensor parameters. 823 self.diff_data.params = params[0:self.diff_end_index] 824 825 # Set the total chi2 gradient to zero. 826 self.total_dchi2 = self.total_dchi2 * 0.0 827 828 # Loop over the residues. 829 for i in xrange(self.num_res): 830 # Set self.data[i] to data. 831 data = self.data[i] 832 833 # Diffusion tensor geometry calculations. 834 if self.diff_data.calc_dgeom: 835 self.diff_data.calc_dgeom(data, self.diff_data) 836 837 # Diffusion tensor weight calculations. 838 if self.diff_data.calc_dci: 839 self.diff_data.calc_dci(data, self.diff_data) 840 841 # Diffusion tensor correlation times. 842 self.diff_data.calc_dti(data, self.diff_data) 843 844 # Calculate the spectral density gradient components. 845 if data.calc_djw_comps: 846 data.calc_djw_comps(data, params) 847 848 # Loop over the gradient. 849 for j in xrange(data.total_num_params): 850 # Calculate the spectral density gradients. 851 if data.calc_djw[j]: 852 data.djw = data.calc_djw[j](data, params, j) 853 else: 854 data.djw = data.djw * 0.0 855 856 # Calculate the relaxation gradient components. 857 data.create_dri_comps(data, params) 858 859 # Calculate the R1, R2, and sigma_noe gradients. 860 data.dri_prime[j] = data.create_dri_prime[j](data) 861 862 # Loop over the relaxation values and modify the NOE gradients. 863 data.dri[j] = data.dri_prime[j] 864 for m in xrange(data.num_ri): 865 if data.create_dri[m]: 866 data.create_dri[m](data, m, data.remap_table[m], data.get_dr1, params, j) 867 868 # Calculate the chi-squared gradient. 869 data.dchi2[j] = dchi2(data.relax_data, data.ri, data.dri[j], data.errors) 870 871 # Index for the construction of the global generic model-free gradient. 872 index = self.diff_data.num_params 873 874 # Diffusion parameter part of the global generic model-free gradient. 875 self.total_dchi2[0:index] = self.total_dchi2[0:index] + data.dchi2[0:index] 876 877 # Model-free parameter part of the global generic model-free gradient. 878 self.total_dchi2[data.start_index:data.end_index] = self.total_dchi2[data.start_index:data.end_index] + data.dchi2[index:] 879 880 # Diagonal scaling. 881 if self.scaling_flag: 882 self.total_dchi2 = matrixmultiply(self.total_dchi2, self.scaling_matrix) 883 884 # Return a copy of the gradient. 885 return self.total_dchi2 * 1.0

886 887

888 - def d2func_mf(self, params):

889 """Function for calculating the chi-squared Hessian. 890 891 Used in the minimisation of model-free parameters for a single residue. 892 """ 893 894 # Test if the gradient has already been called, otherwise run self.dfunc. 895 if sum(params == self.grad_test) != self.total_num_params: 896 self.dfunc(params) 897 898 # Set self.data[0] to data. 899 data = self.data[0] 900 901 # Scaling. 902 if self.scaling_flag: 903 params = matrixmultiply(params, self.scaling_matrix) 904 905 # Loop over the lower triangle of the Hessian. 906 for j in xrange(data.total_num_params): 907 for k in xrange(j + 1): 908 # Calculate the spectral density Hessians. 909 if data.calc_d2jw[j][k]: 910 data.d2jw = data.calc_d2jw[j][k](data, params, j, k) 911 else: 912 data.d2jw = data.d2jw * 0.0 913 914 # Calculate the relaxation Hessian components. 915 data.create_d2ri_comps(data, params) 916 917 # Calculate the R1, R2, and sigma_noe Hessians. 918 if data.create_d2ri_prime[j][k]: 919 data.d2ri_prime[j, k] = data.create_d2ri_prime[j][k](data) 920 921 # Loop over the relaxation values and modify the NOE Hessians. 922 data.d2ri[j, k] = data.d2ri_prime[j, k] 923 for m in xrange(data.num_ri): 924 if data.create_d2ri[m]: 925 data.create_d2ri[m](data, m, data.remap_table[m], data.get_d2r1, params, j, k) 926 927 # Calculate the chi-squared Hessian. 928 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) 929 930 # Diagonal scaling. 931 if self.scaling_flag: 932 data.d2chi2 = matrixmultiply(self.scaling_matrix, matrixmultiply(data.d2chi2, self.scaling_matrix)) 933 934 # Return a copy of the Hessian. 935 return data.d2chi2 * 1.0

936 937

938 - def d2func_local_tm(self, params):

939 """Function for calculating the chi-squared Hessian. 940 941 Used in the minimisation of model-free parameters for a single residue with a local tm. 942 """ 943 944 # Test if the gradient has already been called, otherwise run self.dfunc. 945 if sum(params == self.grad_test) != self.total_num_params: 946 self.dfunc(params) 947 948 # Set self.data[0] to data. 949 data = self.data[0] 950 951 # Scaling. 952 if self.scaling_flag: 953 params = matrixmultiply(params, self.scaling_matrix) 954 955 # Diffusion tensor parameters. 956 self.diff_data.params = params[0:1] 957 958 # Loop over the lower triangle of the Hessian. 959 for j in xrange(data.total_num_params): 960 for k in xrange(j + 1): 961 # Calculate the spectral density Hessians. 962 if data.calc_d2jw[j][k]: 963 data.d2jw = data.calc_d2jw[j][k](data, params, j, k) 964 else: 965 data.d2jw = data.d2jw * 0.0 966 967 # Calculate the relaxation Hessian components. 968 data.create_d2ri_comps(data, params) 969 970 # Calculate the R1, R2, and sigma_noe Hessians. 971 if data.create_d2ri_prime[j][k]: 972 data.d2ri_prime[j, k] = data.create_d2ri_prime[j][k](data) 973 974 # Loop over the relaxation values and modify the NOE Hessians. 975 data.d2ri[j, k] = data.d2ri_prime[j, k] 976 for m in xrange(data.num_ri): 977 if data.create_d2ri[m]: 978 data.create_d2ri[m](data, m, data.remap_table[m], data.get_d2r1, params, j, k) 979 980 # Calculate the chi-squared Hessian. 981 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) 982 983 # Diagonal scaling. 984 if self.scaling_flag: 985 data.d2chi2 = matrixmultiply(self.scaling_matrix, matrixmultiply(data.d2chi2, self.scaling_matrix)) 986 987 # Return a copy of the Hessian. 988 return data.d2chi2 * 1.0

989 990

991 - def d2func_diff(self, params):

992 """Function for calculating the chi-squared Hessian. 993 994 Used in the minimisation of diffusion tensor parameters with all model-free parameters 995 fixed. 996 """ 997 998 # Test if the gradient has already been called, otherwise run self.dfunc. 999 if sum(params == self.grad_test) != self.total_num_params: 1000 self.dfunc(params) 1001 1002 # Scaling. 1003 if self.scaling_flag: 1004 params = matrixmultiply(params, self.scaling_matrix) 1005 1006 # Diffusion tensor parameters. 1007 self.diff_data.params = params[0:self.diff_end_index] 1008 1009 # Set the total chi2 Hessian to zero. 1010 self.total_d2chi2 = self.total_d2chi2 * 0.0 1011 1012 # Loop over the residues. 1013 for i in xrange(self.num_res): 1014 # Set self.data[i] to data. 1015 data = self.data[i] 1016 1017 # Diffusion tensor geometry calculations. 1018 if self.diff_data.calc_d2geom: 1019 self.diff_data.calc_d2geom(data, self.diff_data) 1020 1021 # Diffusion tensor weight calculations. 1022 if self.diff_data.calc_d2ci: 1023 self.diff_data.calc_d2ci(data, self.diff_data) 1024 1025 # Diffusion tensor correlation times. 1026 if self.diff_data.calc_d2ti: 1027 self.diff_data.calc_d2ti(data, self.diff_data) 1028 1029 # Loop over the lower triangle of the Hessian. 1030 for j in xrange(data.total_num_params): 1031 for k in xrange(j + 1): 1032 # Calculate the spectral density Hessians. 1033 if data.calc_d2jw[j][k]: 1034 data.d2jw = data.calc_d2jw[j][k](data, data.param_values, j, k) 1035 else: 1036 data.d2jw = data.d2jw * 0.0 1037 1038 # Calculate the relaxation Hessian components. 1039 data.create_d2ri_comps(data, data.param_values) 1040 1041 # Calculate the R1, R2, and sigma_noe Hessians. 1042 if data.create_d2ri_prime[j][k]: 1043 data.d2ri_prime[j, k] = data.create_d2ri_prime[j][k](data) 1044 1045 # Loop over the relaxation values and modify the NOE Hessians. 1046 data.d2ri[j, k] = data.d2ri_prime[j, k] 1047 for m in xrange(data.num_ri): 1048 if data.create_d2ri[m]: 1049 data.create_d2ri[m](data, m, data.remap_table[m], data.get_d2r1, params, j, k) 1050 1051 # Calculate the chi-squared Hessian. 1052 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) 1053 1054 # Pure diffusion parameter part of the global generic model-free Hessian. 1055 self.total_d2chi2 = self.total_d2chi2 + data.d2chi2 1056 1057 # Diagonal scaling. 1058 if self.scaling_flag: 1059 self.total_d2chi2 = matrixmultiply(self.scaling_matrix, matrixmultiply(self.total_d2chi2, self.scaling_matrix)) 1060 1061 # Return a copy of the Hessian. 1062 return self.total_d2chi2 * 1.0

1063 1064

1065 - def d2func_all(self, params):

1066 """Function for calculating the chi-squared Hessian. 1067 1068 Used in the minimisation of diffusion tensor parameters together with all model-free 1069 parameters. 1070 """ 1071 1072 # Test if the gradient has already been called, otherwise run self.dfunc. 1073 if sum(params == self.grad_test) != self.total_num_params: 1074 self.dfunc(params) 1075 1076 # Scaling. 1077 if self.scaling_flag: 1078 params = matrixmultiply(params, self.scaling_matrix) 1079 1080 # Diffusion tensor parameters. 1081 self.diff_data.params = params[0:self.diff_end_index] 1082 1083 # Set the total chi2 Hessian to zero. 1084 self.total_d2chi2 = self.total_d2chi2 * 0.0 1085 1086 # Loop over the residues. 1087 for i in xrange(self.num_res): 1088 # Set self.data[i] to data. 1089 data = self.data[i] 1090 1091 # Diffusion tensor geometry calculations. 1092 if self.diff_data.calc_d2geom: 1093 self.diff_data.calc_d2geom(data, self.diff_data) 1094 1095 # Diffusion tensor weight calculations. 1096 if self.diff_data.calc_d2ci: 1097 self.diff_data.calc_d2ci(data, self.diff_data) 1098 1099 # Diffusion tensor correlation times. 1100 if self.diff_data.calc_d2ti: 1101 self.diff_data.calc_d2ti(data, self.diff_data) 1102 1103 # Loop over the lower triangle of the Hessian. 1104 for j in xrange(data.total_num_params): 1105 for k in xrange(j + 1): 1106 # Calculate the spectral density Hessians. 1107 if data.calc_d2jw[j][k]: 1108 data.d2jw = data.calc_d2jw[j][k](data, params, j, k) 1109 else: 1110 data.d2jw = data.d2jw * 0.0 1111 1112 # Calculate the relaxation Hessian components. 1113 data.create_d2ri_comps(data, params) 1114 1115 # Calculate the R1, R2, and sigma_noe Hessians. 1116 if data.create_d2ri_prime[j][k]: 1117 data.d2ri_prime[j, k] = data.create_d2ri_prime[j][k](data) 1118 1119 # Loop over the relaxation values and modify the NOE Hessians. 1120 data.d2ri[j, k] = data.d2ri_prime[j, k] 1121 for m in xrange(data.num_ri): 1122 if data.create_d2ri[m]: 1123 data.create_d2ri[m](data, m, data.remap_table[m], data.get_d2r1, params, j, k) 1124 1125 # Calculate the chi-squared Hessian. 1126 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) 1127 1128 # Index for the construction of the global generic model-free Hessian. 1129 index = self.diff_data.num_params 1130 1131 # Pure diffusion parameter part of the global generic model-free Hessian. 1132 self.total_d2chi2[0:index, 0:index] = self.total_d2chi2[0:index, 0:index] + data.d2chi2[0:index, 0:index] 1133 1134 # Pure model-free parameter part of the global generic model-free Hessian. 1135 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:] 1136 1137 # Off diagonal diffusion and model-free parameter parts of the global generic model-free Hessian. 1138 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:] 1139 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] 1140 1141 # Diagonal scaling. 1142 if self.scaling_flag: 1143 self.total_d2chi2 = matrixmultiply(self.scaling_matrix, matrixmultiply(self.total_d2chi2, self.scaling_matrix)) 1144 1145 # Return a copy of the Hessian. 1146 return self.total_d2chi2 * 1.0

1147 1148

1149 - def calc_ri(self):

1150 """Function for calculating relaxation values.""" 1151 1152 # Function call. 1153 chi2 = self.func_mf(self.params) 1154 1155 # Return the single value. 1156 return self.data[0].ri[0]

1157 1158

1159 - def init_diff_data(self, diff_data):

1160 """Function for the initialisation of diffusion tensor specific data.""" 1161 1162 # Isotropic diffusion. 1163 if diff_data.type == 'iso': 1164 # Number of diffusion parameters. 1165 diff_data.num_params = 1 1166 1167 # Number of indecies in the generic equations. 1168 diff_data.num_indecies = 1 1169 1170 # Geometry function, gradient, and Hessian. 1171 diff_data.calc_geom = None 1172 diff_data.calc_dgeom = None 1173 diff_data.calc_d2geom = None 1174 1175 # Weight function, gradient, and Hessian. 1176 diff_data.calc_ci = calc_iso_ci 1177 diff_data.calc_dci = None 1178 diff_data.calc_d2ci = None 1179 1180 # Global correlation time function, gradient, and Hessian. 1181 diff_data.calc_ti = calc_iso_ti 1182 diff_data.calc_dti = calc_iso_dti 1183 diff_data.calc_d2ti = None 1184 1185 1186 # Axially symmetric diffusion. 1187 elif diff_data.type == 'axial': 1188 # Number of diffusion parameters. 1189 diff_data.num_params = 4 1190 1191 # Number of indecies in the generic equations. 1192 diff_data.num_indecies = 3 1193 1194 # Geometry function, gradient, and Hessian. 1195 diff_data.calc_geom = calc_axial_geom 1196 diff_data.calc_dgeom = calc_axial_dgeom 1197 diff_data.calc_d2geom = calc_axial_d2geom 1198 1199 # Weight function, gradient, and Hessian. 1200 diff_data.calc_ci = calc_axial_ci 1201 diff_data.calc_dci = calc_axial_dci 1202 diff_data.calc_d2ci = calc_axial_d2ci 1203 1204 # Global correlation time function, gradient, and Hessian. 1205 diff_data.calc_ti = calc_axial_ti 1206 diff_data.calc_dti = calc_axial_dti 1207 diff_data.calc_d2ti = calc_axial_d2ti 1208 1209 # Unit vectors. 1210 diff_data.dpar = zeros(3, Float64) 1211 1212 # Unit vector gradients. 1213 diff_data.dpar_dtheta = zeros(3, Float64) 1214 diff_data.dpar_dphi = zeros(3, Float64) 1215 1216 # Unit vector Hessians. 1217 diff_data.dpar_dtheta2 = zeros(3, Float64) 1218 diff_data.dpar_dthetadphi = zeros(3, Float64) 1219 diff_data.dpar_dphi2 = zeros(3, Float64) 1220 1221 # Anisotropic diffusion. 1222 elif diff_data.type == 'aniso': 1223 # Number of diffusion parameters. 1224 diff_data.num_params = 6 1225 1226 # Number of indecies in the generic equations. 1227 diff_data.num_indecies = 5 1228 1229 # Geometry function, gradient, and Hessian. 1230 diff_data.calc_geom = calc_aniso_geom 1231 diff_data.calc_dgeom = calc_aniso_dgeom 1232 diff_data.calc_d2geom = calc_aniso_d2geom 1233 1234 # Weight function, gradient, and Hessian. 1235 diff_data.calc_ci = calc_aniso_ci 1236 diff_data.calc_dci = calc_aniso_dci 1237 diff_data.calc_d2ci = calc_aniso_d2ci 1238 1239 # Global correlation time function, gradient, and Hessian. 1240 diff_data.calc_ti = calc_aniso_ti 1241 diff_data.calc_dti = calc_aniso_dti 1242 diff_data.calc_d2ti = calc_aniso_d2ti 1243 1244 # Unit vectors. 1245 diff_data.dx = zeros(3, Float64) 1246 diff_data.dy = zeros(3, Float64) 1247 diff_data.dz = zeros(3, Float64) 1248 1249 # Unit vector gradients. 1250 diff_data.ddx_dalpha = zeros(3, Float64) 1251 diff_data.ddx_dbeta = zeros(3, Float64) 1252 diff_data.ddx_dgamma = zeros(3, Float64) 1253 1254 diff_data.ddy_dalpha = zeros(3, Float64) 1255 diff_data.ddy_dbeta = zeros(3, Float64) 1256 diff_data.ddy_dgamma = zeros(3, Float64) 1257 1258 diff_data.ddz_dalpha = zeros(3, Float64) 1259 diff_data.ddz_dbeta = zeros(3, Float64) 1260 diff_data.ddz_dgamma = zeros(3, Float64) 1261 1262 # Unit vector Hessians. 1263 diff_data.d2dx_dalpha2 = zeros(3, Float64) 1264 diff_data.d2dx_dalpha_dbeta = zeros(3, Float64) 1265 diff_data.d2dx_dalpha_dgamma = zeros(3, Float64) 1266 diff_data.d2dx_dbeta2 = zeros(3, Float64) 1267 diff_data.d2dx_dbeta_dgamma = zeros(3, Float64) 1268 diff_data.d2dx_dgamma2 = zeros(3, Float64) 1269 1270 diff_data.d2dy_dalpha2 = zeros(3, Float64) 1271 diff_data.d2dy_dalpha_dbeta = zeros(3, Float64) 1272 diff_data.d2dy_dalpha_dgamma = zeros(3, Float64) 1273 diff_data.d2dy_dbeta2 = zeros(3, Float64) 1274 diff_data.d2dy_dbeta_dgamma = zeros(3, Float64) 1275 diff_data.d2dy_dgamma2 = zeros(3, Float64) 1276 1277 diff_data.d2dz_dalpha2 = zeros(3, Float64) 1278 diff_data.d2dz_dalpha_dbeta = zeros(3, Float64) 1279 diff_data.d2dz_dalpha_dgamma = zeros(3, Float64) 1280 diff_data.d2dz_dbeta2 = zeros(3, Float64) 1281 diff_data.d2dz_dbeta_dgamma = zeros(3, Float64) 1282 diff_data.d2dz_dgamma2 = zeros(3, Float64)

1283 1284

1285 - def init_res_data(self, data, diff_data):

1286 """Function for the initialisation of the residue specific data.""" 1287 1288 # Weights and global correlation time values. 1289 data.ci = zeros(diff_data.num_indecies, Float64) 1290 data.ti = zeros(diff_data.num_indecies, Float64) 1291 1292 # Isotropic diffusion. 1293 if self.diff_data.type == 'iso': 1294 # Global correlation time gradient and Hessian. 1295 data.dti = zeros((1, diff_data.num_indecies), Float64) 1296 data.d2ti = zeros((1, 1, diff_data.num_indecies), Float64) 1297 1298 # Axially symmetric diffusion. 1299 elif self.diff_data.type == 'axial': 1300 # Weight gradient and Hessian. 1301 data.dci = zeros((4, diff_data.num_indecies), Float64) 1302 data.d2ci = zeros((4, 4, diff_data.num_indecies), Float64) 1303 1304 # Global correlation time gradient and Hessian. 1305 data.dti = zeros((2, diff_data.num_indecies), Float64) 1306 data.d2ti = zeros((2, 2, diff_data.num_indecies), Float64) 1307 1308 # Dot product. 1309 data.delta = 0 1310 1311 # Dot product gradient. 1312 data.ddelta_dpsi = zeros(2, Float64) 1313 1314 # Dot product Hessian. 1315 data.d2delta_dpsi2 = zeros((2, 2), Float64) 1316 1317 # Anisotropic diffusion. 1318 elif self.diff_data.type == 'aniso': 1319 # Weight gradient and Hessian. 1320 data.dci = zeros((6, diff_data.num_indecies), Float64) 1321 data.d2ci = zeros((6, 6, diff_data.num_indecies), Float64) 1322 1323 # Global correlation time gradient and Hessian. 1324 data.dti = zeros((3, diff_data.num_indecies), Float64) 1325 data.d2ti = zeros((3, 3, diff_data.num_indecies), Float64) 1326 1327 # Dot products. 1328 data.delta_alpha = 0.0 1329 data.delta_beta = 0.0 1330 data.delta_gamma = 0.0 1331 1332 # Dot product gradients. 1333 data.ddelta_alpha_dpsi = zeros(3, Float64) 1334 data.ddelta_beta_dpsi = zeros(3, Float64) 1335 data.ddelta_gamma_dpsi = zeros(3, Float64) 1336 1337 # Dot product Hessians. 1338 data.d2delta_alpha_dpsi2 = zeros((3, 3), Float64) 1339 data.d2delta_beta_dpsi2 = zeros((3, 3), Float64) 1340 data.d2delta_gamma_dpsi2 = zeros((3, 3), Float64) 1341 1342 # Empty spectral density components. 1343 data.w_ti_sqrd = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1344 data.fact_ti = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1345 data.w_te_ti_sqrd = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1346 data.w_tf_ti_sqrd = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1347 data.w_ts_ti_sqrd = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1348 data.inv_te_denom = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1349 data.inv_tf_denom = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1350 data.inv_ts_denom = zeros((data.num_frq, 5, diff_data.num_indecies), Float64) 1351 1352 # Empty spectral density values, gradients, and Hessians. 1353 data.jw = zeros((data.num_frq, 5), Float64) 1354 data.djw = zeros((data.num_frq, 5), Float64) 1355 data.d2jw = zeros((data.num_frq, 5), Float64) 1356 1357 # Calculate the fixed components of the dipolar and CSA constants. 1358 data.csa_const_fixed = zeros(data.num_frq, Float64) 1359 data.dip_const_fixed = None 1360 calc_fixed_csa(data) 1361 calc_fixed_dip(data) 1362 1363 # Dipolar and CSA constants. 1364 data.dip_const_func = 0.0 1365 data.dip_const_grad = 0.0 1366 data.dip_const_hess = 0.0 1367 data.csa_const_func = zeros(data.num_frq, Float64) 1368 data.csa_const_grad = zeros(data.num_frq, Float64) 1369 data.csa_const_hess = zeros(data.num_frq, Float64) 1370 1371 # Components of the transformed relaxation equations. 1372 data.dip_comps_func = zeros(data.num_ri, Float64) 1373 data.csa_comps_func = zeros(data.num_ri, Float64) 1374 data.rex_comps_func = zeros(data.num_ri, Float64) 1375 data.dip_jw_comps_func = zeros(data.num_ri, Float64) 1376 data.csa_jw_comps_func = zeros(data.num_ri, Float64) 1377 1378 # First partial derivative components of the transformed relaxation equations. 1379 data.dip_comps_grad = zeros(data.num_ri, Float64) 1380 data.csa_comps_grad = zeros(data.num_ri, Float64) 1381 data.rex_comps_grad = zeros(data.num_ri, Float64) 1382 data.dip_jw_comps_grad = zeros(data.num_ri, Float64) 1383 data.csa_jw_comps_grad = zeros(data.num_ri, Float64) 1384 1385 # First partial derivative components of the transformed relaxation equations. 1386 data.dip_comps_hess = zeros(data.num_ri, Float64) 1387 data.csa_comps_hess = zeros(data.num_ri, Float64) 1388 data.rex_comps_hess = zeros(data.num_ri, Float64) 1389 data.dip_jw_comps_hess = zeros(data.num_ri, Float64) 1390 data.csa_jw_comps_hess = zeros(data.num_ri, Float64) 1391 1392 # Transformed relaxation values, gradients, and Hessians. 1393 data.ri_prime = zeros((data.num_ri), Float64) 1394 data.dri_prime = zeros((data.total_num_params, data.num_ri), Float64) 1395 data.d2ri_prime = zeros((data.total_num_params, data.total_num_params, data.num_ri), Float64) 1396 1397 # Data structures containing the Ri values. 1398 data.ri = zeros(data.num_ri, Float64) 1399 data.dri = zeros((data.total_num_params, data.num_ri), Float64) 1400 data.d2ri = zeros((data.total_num_params, data.total_num_params, data.num_ri), Float64) 1401 1402 # Data structures containing the R1 values at the position of and corresponding to the NOE. 1403 data.r1 = zeros(data.num_ri, Float64) 1404 data.dr1 = zeros((data.total_num_params, data.num_ri), Float64) 1405 data.d2r1 = zeros((data.total_num_params, data.total_num_params, data.num_ri), Float64) 1406 1407 # Data structures containing the chi-squared values. 1408 data.chi2 = 0.0 1409 data.dchi2 = zeros((data.total_num_params), Float64) 1410 data.d2chi2 = zeros((data.total_num_params, data.total_num_params), Float64)

1411 1412

1413 - def init_res_r1_data(self, data):

1414 """Function for initialisation of the R1 data class. 1415 1416 This data class is only used if an NOE data set is collected but no R1 data set 1417 corresponding to the same frequency exists. 1418 """ 1419 1420 # Initialise an instance of Data. 1421 r1_data = Data() 1422 1423 # Copy a few things from 'data' to 'r1_data' 1424 r1_data.num_frq = data.num_frq 1425 r1_data.dip_const_fixed = data.dip_const_fixed 1426 r1_data.csa_const_fixed = data.csa_const_fixed 1427 1428 # Components of the transformed relaxation equations. 1429 r1_data.dip_comps_func = zeros(data.num_ri, Float64) 1430 r1_data.csa_comps_func = zeros(data.num_ri, Float64) 1431 r1_data.dip_jw_comps_func = zeros(data.num_ri, Float64) 1432 r1_data.csa_jw_comps_func = zeros(data.num_ri, Float64) 1433 1434 # Initialise the first partial derivative components of the transformed relaxation equations. 1435 r1_data.dip_comps_grad = zeros(data.num_ri, Float64) 1436 r1_data.csa_comps_grad = zeros(data.num_ri, Float64) 1437 r1_data.rex_comps_grad = zeros(data.num_ri, Float64) 1438 r1_data.dip_jw_comps_grad = zeros(data.num_ri, Float64) 1439 r1_data.csa_jw_comps_grad = zeros(data.num_ri, Float64) 1440 1441 # Initialise the first partial derivative components of the transformed relaxation equations. 1442 r1_data.dip_comps_hess = zeros(data.num_ri, Float64) 1443 r1_data.csa_comps_hess = zeros(data.num_ri, Float64) 1444 r1_data.rex_comps_hess = zeros(data.num_ri, Float64) 1445 r1_data.dip_jw_comps_hess = zeros(data.num_ri, Float64) 1446 r1_data.csa_jw_comps_hess = zeros(data.num_ri, Float64) 1447 1448 # Initialise the transformed relaxation values, gradients, and Hessians. 1449 r1_data.ri_prime = zeros(data.num_ri, Float64) 1450 r1_data.dri_prime = zeros((data.num_ri, data.total_num_params), Float64) 1451 r1_data.d2ri_prime = zeros((data.num_ri, data.total_num_params, data.total_num_params), Float64) 1452 1453 # 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. 1454 r1_data.create_dri_prime = data.create_dri_prime 1455 r1_data.create_d2ri_prime = data.create_d2ri_prime 1456 1457 # CSA, bond length, and Rex indecies. 1458 r1_data.csa_index = data.csa_index 1459 r1_data.r_index = data.r_index 1460 r1_data.rex_index = data.rex_index 1461 1462 # Place 'r1_data' into 'data'. 1463 data.r1_data = r1_data

1464 1465

1466 - def lm_dri(self):

1467 """Return the function used for Levenberg-Marquardt minimisation.""" 1468 1469 # Create dri. 1470 if self.param_set == 'mf' or self.param_set == 'local_tm': 1471 dri = self.data[0].dri 1472 elif self.param_set == 'diff': 1473 # Set the total dri gradient to zero. 1474 self.total_dri = self.total_dri * 0.0 1475 1476 # Ri indecies. 1477 ri_start_index = 0 1478 ri_end_index = 0 1479 1480 # Loop over the residues. 1481 for i in xrange(self.num_res): 1482 # Set self.data[i] to data. 1483 data = self.data[i] 1484 1485 # Increment Ri end index. 1486 ri_end_index = ri_end_index + data.num_ri 1487 1488 # Diffusion parameter part of the global generic model-free gradient. 1489 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] 1490 1491 # Increment Ri start index. 1492 ri_start_index = ri_start_index + data.num_ri 1493 1494 # dri. 1495 dri = self.total_dri 1496 1497 elif self.param_set == 'all': 1498 # Set the total dri gradient to zero. 1499 self.total_dri = self.total_dri * 0.0 1500 1501 # Ri indecies. 1502 ri_start_index = 0 1503 ri_end_index = 0 1504 1505 # Loop over the residues. 1506 for i in xrange(self.num_res): 1507 # Set self.data[i] to data. 1508 data = self.data[i] 1509 1510 # Increment Ri end index. 1511 ri_end_index = ri_end_index + data.num_ri 1512 1513 # Diffusion parameter part of the global generic model-free gradient. 1514 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] 1515 1516 # Model-free parameter part of the global generic model-free gradient. 1517 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:] 1518 1519 # Increment Ri start index. 1520 ri_start_index = ri_start_index + data.num_ri 1521 1522 # dri. 1523 dri = self.total_dri 1524 1525 # Make the proper Jacobian. 1526 dri = transpose(dri) 1527 1528 # Diagonal scaling. 1529 if self.scaling_flag: 1530 dri = matrixmultiply(dri, self.scaling_matrix) 1531 1532 # Return dri. 1533 return dri

1534 1535

1536 - def setup_equations(self, data):

1537 """Setup all the residue specific equations.""" 1538 1539 # Initialisation. 1540 ################# 1541 1542 # The number of diffusion parameters. 1543 if self.param_set != 'all': 1544 num_diff_params = 0 1545 elif self.diff_data.type == 'iso': 1546 num_diff_params = 1 1547 elif self.diff_data.type == 'axial': 1548 num_diff_params = 4 1549 elif self.diff_data.type == 'aniso': 1550 num_diff_params = 6 1551 1552 # Indecies. 1553 data.tm_index, data.tm_local_index = None, None 1554 data.s2_index, data.s2_local_index = None, None 1555 data.s2f_index, data.s2f_local_index = None, None 1556 data.s2s_index, data.s2s_local_index = None, None 1557 data.te_index, data.te_local_index = None, None 1558 data.tf_index, data.tf_local_index = None, None 1559 data.ts_index, data.ts_local_index = None, None 1560 data.rex_index, data.rex_local_index = None, None 1561 data.r_index, data.r_local_index = None, None 1562 data.csa_index, data.csa_local_index = None, None 1563 1564 # Set up the spectral density functions. 1565 ######################################## 1566 1567 # Create empty spectral density gradient and Hessian function data structures. 1568 data.calc_djw = [] 1569 data.calc_d2jw = [] 1570 for i in xrange(data.total_num_params): 1571 data.calc_djw.append(None) 1572 data.calc_d2jw.append([]) 1573 for j in xrange(data.total_num_params): 1574 data.calc_d2jw[i].append(None) 1575 1576 1577 # The original model-free equations {S2, te, Rex, r, CSA}. 1578 ########################################################## 1579 1580 if data.equations == 'mf_orig': 1581 # Find the indecies of the model-free parameters. 1582 for i in xrange(data.num_params): 1583 if data.param_types[i] == 'S2': 1584 data.s2_local_index = num_diff_params + i 1585 data.s2_index = self.param_index + i 1586 elif data.param_types[i] == 'te': 1587 data.te_local_index = num_diff_params + i 1588 data.te_index = self.param_index + i 1589 elif data.param_types[i] == 'Rex': 1590 data.rex_local_index = num_diff_params + i 1591 data.rex_index = self.param_index + i 1592 elif data.param_types[i] == 'r': 1593 data.r_local_index = num_diff_params + i 1594 data.r_index = self.param_index + i 1595 elif data.param_types[i] == 'CSA': 1596 data.csa_local_index = num_diff_params + i 1597 data.csa_index = self.param_index + i 1598 elif data.param_types[i] == 'tm': 1599 pass 1600 else: 1601 print "Unknown parameter." 1602 return 0 1603 1604 # Increment the parameter index. 1605 self.param_index = self.param_index + data.num_params 1606 1607 # Single residue minimisation with fixed diffusion parameters. 1608 if self.param_set == 'mf': 1609 # No model-free parameters {}. 1610 if data.s2_index == None and data.te_index == None: 1611 # Equation. 1612 data.calc_jw_comps = None 1613 data.calc_jw = calc_jw 1614 1615 # Gradient. 1616 data.calc_djw_comps = None 1617 1618 # Model-free parameters {S2}. 1619 elif data.s2_index != None and data.te_index == None: 1620 # Equation. 1621 data.calc_jw_comps = None 1622 data.calc_jw = calc_S2_jw 1623 1624 # Gradient. 1625 data.calc_djw_comps = None 1626 data.calc_djw[data.s2_local_index] = calc_S2_djw_dS2 1627 1628 # Model-free parameters {S2, te}. 1629 elif data.s2_index != None and data.te_index != None: 1630 # Equation. 1631 data.calc_jw_comps = calc_S2_te_jw_comps 1632 data.calc_jw = calc_S2_te_jw 1633 1634 # Gradient. 1635 data.calc_djw_comps = calc_S2_te_djw_comps 1636 data.calc_djw[data.s2_local_index] = calc_S2_te_djw_dS2 1637 data.calc_djw[data.te_local_index] = calc_S2_te_djw_dte 1638 1639 # Hessian. 1640 data.calc_d2jw[data.s2_local_index][data.te_local_index] = data.calc_d2jw[data.te_local_index][data.s2_local_index] = calc_S2_te_d2jw_dS2dte 1641 data.calc_d2jw[data.te_local_index][data.te_local_index] = calc_S2_te_d2jw_dte2 1642 1643 # Bad parameter combination. 1644 else: 1645 print "Invalid combination of parameters for the extended model-free equation." 1646 return 0 1647 1648 # Minimisation with variable diffusion parameters. 1649 else: 1650 # Diffusion parameters and no model-free parameters {}. 1651 if data.s2_index == None and data.te_index == None: 1652 # Equation. 1653 data.calc_jw_comps = None 1654 data.calc_jw = calc_jw 1655 1656 # Gradient. 1657 data.calc_djw_comps = calc_diff_djw_comps 1658 1659 # Isotropic diffusion. 1660 if self.diff_data.type == 'iso': 1661 # Gradient. 1662 data.calc_djw[0] = calc_diff_djw_dDj 1663 1664 # Hessian. 1665 data.calc_d2jw[0][0] = calc_diff_d2jw_dDjdDk 1666 1667 # Axially symmetric diffusion. 1668 elif self.diff_data.type == 'axial': 1669 # Gradient. 1670 data.calc_djw[0] = data.calc_djw[1] = calc_diff_djw_dDj 1671 data.calc_djw[2] = data.calc_djw[3] = calc_diff_djw_dPsij 1672 1673 # Hessian. 1674 data.calc_d2jw[0][0] = calc_diff_d2jw_dDjdDk 1675 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_d2jw_dDjdDk 1676 data.calc_d2jw[1][1] = calc_diff_d2jw_dDjdDk 1677 1678 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_d2jw_dDjdPsij 1679 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_d2jw_dDjdPsij 1680 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_d2jw_dDjdPsij 1681 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_d2jw_dDjdPsij 1682 1683 data.calc_d2jw[2][2] = calc_diff_d2jw_dPsijdPsik 1684 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_d2jw_dPsijdPsik 1685 data.calc_d2jw[3][3] = calc_diff_d2jw_dPsijdPsik 1686 1687 # Anisotropic diffusion. 1688 elif self.diff_data.type == 'aniso': 1689 # Gradient. 1690 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_aniso_djw_dDj 1691 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_djw_dPsij 1692 1693 # Hessian. 1694 data.calc_d2jw[0][0] = calc_aniso_d2jw_dDjdDk 1695 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_aniso_d2jw_dDjdDk 1696 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_aniso_d2jw_dDjdDk 1697 data.calc_d2jw[1][1] = calc_aniso_d2jw_dDjdDk 1698 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_aniso_d2jw_dDjdDk 1699 data.calc_d2jw[2][2] = calc_aniso_d2jw_dDjdDk 1700 1701 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_aniso_d2jw_dDjdPsij 1702 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_aniso_d2jw_dDjdPsij 1703 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_aniso_d2jw_dDjdPsij 1704 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_aniso_d2jw_dDjdPsij 1705 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_aniso_d2jw_dDjdPsij 1706 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_aniso_d2jw_dDjdPsij 1707 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_aniso_d2jw_dDjdPsij 1708 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_aniso_d2jw_dDjdPsij 1709 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_aniso_d2jw_dDjdPsij 1710 1711 data.calc_d2jw[3][3] = calc_diff_d2jw_dPsijdPsik 1712 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_d2jw_dPsijdPsik 1713 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_d2jw_dPsijdPsik 1714 data.calc_d2jw[4][4] = calc_diff_d2jw_dPsijdPsik 1715 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_d2jw_dPsijdPsik 1716 data.calc_d2jw[5][5] = calc_diff_d2jw_dPsijdPsik 1717 1718 # Diffusion parameters and model-free parameters {S2}. 1719 elif data.s2_index != None and data.te_index == None: 1720 # Equation. 1721 data.calc_jw_comps = None 1722 data.calc_jw = calc_S2_jw 1723 1724 # Gradient. 1725 data.calc_djw_comps = calc_diff_djw_comps 1726 1727 if self.param_set != 'diff': 1728 # Gradient. 1729 data.calc_djw[data.s2_local_index] = calc_S2_djw_dS2 1730 1731 # Isotropic diffusion. 1732 if self.diff_data.type == 'iso': 1733 # Gradient. 1734 data.calc_djw[0] = calc_diff_S2_djw_dDj 1735 1736 # Hessian. 1737 data.calc_d2jw[0][0] = calc_diff_S2_d2jw_dDjdDk 1738 if self.param_set != 'diff': 1739 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_diff_S2_d2jw_dDjdS2 1740 1741 # Axially symmetric diffusion. 1742 elif self.diff_data.type == 'axial': 1743 # Gradient. 1744 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2_djw_dDj 1745 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2_djw_dPsij 1746 1747 # Hessian. 1748 data.calc_d2jw[0][0] = calc_diff_S2_d2jw_dDjdDk 1749 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2_d2jw_dDjdDk 1750 data.calc_d2jw[1][1] = calc_diff_S2_d2jw_dDjdDk 1751 1752 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2_d2jw_dDjdPsij 1753 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2_d2jw_dDjdPsij 1754 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2_d2jw_dDjdPsij 1755 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2_d2jw_dDjdPsij 1756 1757 data.calc_d2jw[2][2] = calc_diff_S2_d2jw_dPsijdPsik 1758 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2_d2jw_dPsijdPsik 1759 data.calc_d2jw[3][3] = calc_diff_S2_d2jw_dPsijdPsik 1760 1761 if self.param_set != 'diff': 1762 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_diff_S2_d2jw_dDjdS2 1763 data.calc_d2jw[1][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][1] = calc_diff_S2_d2jw_dDjdS2 1764 data.calc_d2jw[2][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][2] = calc_diff_S2_d2jw_dPsijdS2 1765 data.calc_d2jw[3][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][3] = calc_diff_S2_d2jw_dPsijdS2 1766 1767 # Anisotropic diffusion. 1768 elif self.diff_data.type == 'aniso': 1769 # Gradient. 1770 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_aniso_S2_djw_dDj 1771 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2_djw_dPsij 1772 1773 # Hessian. 1774 data.calc_d2jw[0][0] = calc_aniso_S2_d2jw_dDjdDk 1775 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_aniso_S2_d2jw_dDjdDk 1776 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_aniso_S2_d2jw_dDjdDk 1777 data.calc_d2jw[1][1] = calc_aniso_S2_d2jw_dDjdDk 1778 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_aniso_S2_d2jw_dDjdDk 1779 data.calc_d2jw[2][2] = calc_aniso_S2_d2jw_dDjdDk 1780 1781 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_aniso_S2_d2jw_dDjdPsij 1782 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_aniso_S2_d2jw_dDjdPsij 1783 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_aniso_S2_d2jw_dDjdPsij 1784 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_aniso_S2_d2jw_dDjdPsij 1785 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_aniso_S2_d2jw_dDjdPsij 1786 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_aniso_S2_d2jw_dDjdPsij 1787 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_aniso_S2_d2jw_dDjdPsij 1788 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_aniso_S2_d2jw_dDjdPsij 1789 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_aniso_S2_d2jw_dDjdPsij 1790 1791 data.calc_d2jw[3][3] = calc_diff_S2_d2jw_dPsijdPsik 1792 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2_d2jw_dPsijdPsik 1793 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2_d2jw_dPsijdPsik 1794 data.calc_d2jw[4][4] = calc_diff_S2_d2jw_dPsijdPsik 1795 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2_d2jw_dPsijdPsik 1796 data.calc_d2jw[5][5] = calc_diff_S2_d2jw_dPsijdPsik 1797 1798 if self.param_set != 'diff': 1799 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_aniso_S2_d2jw_dDjdS2 1800 data.calc_d2jw[1][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][1] = calc_aniso_S2_d2jw_dDjdS2 1801 data.calc_d2jw[2][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][2] = calc_aniso_S2_d2jw_dDjdS2 1802 data.calc_d2jw[3][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][3] = calc_diff_S2_d2jw_dPsijdS2 1803 data.calc_d2jw[4][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][4] = calc_diff_S2_d2jw_dPsijdS2 1804 data.calc_d2jw[5][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][5] = calc_diff_S2_d2jw_dPsijdS2 1805 1806 1807 # Diffusion parameters and model-free parameters {S2, te}. 1808 elif data.s2_index != None and data.te_index != None: 1809 # Equation. 1810 data.calc_jw_comps = calc_S2_te_jw_comps 1811 data.calc_jw = calc_S2_te_jw 1812 1813 # Gradient. 1814 data.calc_djw_comps = calc_diff_S2_te_djw_comps 1815 1816 if self.param_set != 'diff': 1817 # Gradient. 1818 data.calc_djw[data.s2_local_index] = calc_S2_te_djw_dS2 1819 data.calc_djw[data.te_local_index] = calc_S2_te_djw_dte 1820 1821 # Hessian. 1822 data.calc_d2jw[data.s2_local_index][data.te_local_index] = data.calc_d2jw[data.te_local_index][data.s2_local_index] = calc_S2_te_d2jw_dS2dte 1823 data.calc_d2jw[data.te_local_index][data.te_local_index] = calc_S2_te_d2jw_dte2 1824 1825 # Isotropic diffusion. 1826 if self.diff_data.type == 'iso': 1827 # Gradient. 1828 data.calc_djw[0] = calc_diff_S2_te_djw_dDj 1829 1830 # Hessian. 1831 data.calc_d2jw[0][0] = calc_diff_S2_te_d2jw_dDjdDk 1832 if self.param_set != 'diff': 1833 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_diff_S2_te_d2jw_dDjdS2 1834 data.calc_d2jw[0][data.te_local_index] = data.calc_d2jw[data.te_local_index][0] = calc_diff_S2_te_d2jw_dDjdte 1835 1836 # Axially symmetric diffusion. 1837 elif self.diff_data.type == 'axial': 1838 # Gradient. 1839 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2_te_djw_dDj 1840 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2_te_djw_dPsij 1841 1842 # Hessian. 1843 data.calc_d2jw[0][0] = calc_diff_S2_te_d2jw_dDjdDk 1844 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2_te_d2jw_dDjdDk 1845 data.calc_d2jw[1][1] = calc_diff_S2_te_d2jw_dDjdDk 1846 1847 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2_te_d2jw_dDjdPsij 1848 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2_te_d2jw_dDjdPsij 1849 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2_te_d2jw_dDjdPsij 1850 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2_te_d2jw_dDjdPsij 1851 1852 data.calc_d2jw[2][2] = calc_diff_S2_te_d2jw_dPsijdPsik 1853 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2_te_d2jw_dPsijdPsik 1854 data.calc_d2jw[3][3] = calc_diff_S2_te_d2jw_dPsijdPsik 1855 1856 if self.param_set != 'diff': 1857 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_diff_S2_te_d2jw_dDjdS2 1858 data.calc_d2jw[1][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][1] = calc_diff_S2_te_d2jw_dDjdS2 1859 data.calc_d2jw[2][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][2] = calc_diff_S2_te_d2jw_dPsijdS2 1860 data.calc_d2jw[3][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][3] = calc_diff_S2_te_d2jw_dPsijdS2 1861 1862 data.calc_d2jw[0][data.te_local_index] = data.calc_d2jw[data.te_local_index][0] = calc_diff_S2_te_d2jw_dDjdte 1863 data.calc_d2jw[1][data.te_local_index] = data.calc_d2jw[data.te_local_index][1] = calc_diff_S2_te_d2jw_dDjdte 1864 data.calc_d2jw[2][data.te_local_index] = data.calc_d2jw[data.te_local_index][2] = calc_diff_S2_te_d2jw_dPsijdte 1865 data.calc_d2jw[3][data.te_local_index] = data.calc_d2jw[data.te_local_index][3] = calc_diff_S2_te_d2jw_dPsijdte 1866 1867 # Anisotropic diffusion. 1868 elif self.diff_data.type == 'aniso': 1869 # Gradient. 1870 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_aniso_S2_te_djw_dDj 1871 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2_te_djw_dPsij 1872 1873 # Hessian. 1874 data.calc_d2jw[0][0] = calc_aniso_S2_te_d2jw_dDjdDk 1875 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_aniso_S2_te_d2jw_dDjdDk 1876 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_aniso_S2_te_d2jw_dDjdDk 1877 data.calc_d2jw[1][1] = calc_aniso_S2_te_d2jw_dDjdDk 1878 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_aniso_S2_te_d2jw_dDjdDk 1879 data.calc_d2jw[2][2] = calc_aniso_S2_te_d2jw_dDjdDk 1880 1881 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_aniso_S2_te_d2jw_dDjdPsij 1882 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_aniso_S2_te_d2jw_dDjdPsij 1883 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_aniso_S2_te_d2jw_dDjdPsij 1884 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_aniso_S2_te_d2jw_dDjdPsij 1885 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_aniso_S2_te_d2jw_dDjdPsij 1886 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_aniso_S2_te_d2jw_dDjdPsij 1887 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_aniso_S2_te_d2jw_dDjdPsij 1888 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_aniso_S2_te_d2jw_dDjdPsij 1889 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_aniso_S2_te_d2jw_dDjdPsij 1890 1891 data.calc_d2jw[3][3] = calc_diff_S2_te_d2jw_dPsijdPsik 1892 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2_te_d2jw_dPsijdPsik 1893 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2_te_d2jw_dPsijdPsik 1894 data.calc_d2jw[4][4] = calc_diff_S2_te_d2jw_dPsijdPsik 1895 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2_te_d2jw_dPsijdPsik 1896 data.calc_d2jw[5][5] = calc_diff_S2_te_d2jw_dPsijdPsik 1897 1898 if self.param_set != 'diff': 1899 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_aniso_S2_te_d2jw_dDjdS2 1900 data.calc_d2jw[1][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][1] = calc_aniso_S2_te_d2jw_dDjdS2 1901 data.calc_d2jw[2][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][2] = calc_aniso_S2_te_d2jw_dDjdS2 1902 data.calc_d2jw[3][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][3] = calc_diff_S2_te_d2jw_dPsijdS2 1903 data.calc_d2jw[4][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][4] = calc_diff_S2_te_d2jw_dPsijdS2 1904 data.calc_d2jw[5][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][5] = calc_diff_S2_te_d2jw_dPsijdS2 1905 1906 data.calc_d2jw[0][data.te_local_index] = data.calc_d2jw[data.te_local_index][0] = calc_aniso_S2_te_d2jw_dDjdte 1907 data.calc_d2jw[1][data.te_local_index] = data.calc_d2jw[data.te_local_index][1] = calc_aniso_S2_te_d2jw_dDjdte 1908 data.calc_d2jw[2][data.te_local_index] = data.calc_d2jw[data.te_local_index][2] = calc_aniso_S2_te_d2jw_dDjdte 1909 data.calc_d2jw[3][data.te_local_index] = data.calc_d2jw[data.te_local_index][3] = calc_diff_S2_te_d2jw_dPsijdte 1910 data.calc_d2jw[4][data.te_local_index] = data.calc_d2jw[data.te_local_index][4] = calc_diff_S2_te_d2jw_dPsijdte 1911 data.calc_d2jw[5][data.te_local_index] = data.calc_d2jw[data.te_local_index][5] = calc_diff_S2_te_d2jw_dPsijdte 1912 1913 # Bad parameter combination. 1914 else: 1915 print "Invalid combination of parameters for the extended model-free equation." 1916 return 0 1917 1918 1919 1920 # The extended model-free equations {S2f, tf, S2, ts, Rex, r, CSA}. 1921 ################################################################### 1922 1923 elif data.equations == 'mf_ext': 1924 # Find the indecies of the model-free parameters. 1925 for i in xrange(data.num_params): 1926 if data.param_types[i] == 'S2f': 1927 data.s2f_local_index = num_diff_params + i 1928 data.s2f_index = self.param_index + i 1929 elif data.param_types[i] == 'tf': 1930 data.tf_local_index = num_diff_params + i 1931 data.tf_index = self.param_index + i 1932 elif data.param_types[i] == 'S2': 1933 data.s2_local_index = num_diff_params + i 1934 data.s2_index = self.param_index + i 1935 elif data.param_types[i] == 'ts': 1936 data.ts_local_index = num_diff_params + i 1937 data.ts_index = self.param_index + i 1938 elif data.param_types[i] == 'Rex': 1939 data.rex_local_index = num_diff_params + i 1940 data.rex_index = self.param_index + i 1941 elif data.param_types[i] == 'r': 1942 data.r_local_index = num_diff_params + i 1943 data.r_index = self.param_index + i 1944 elif data.param_types[i] == 'CSA': 1945 data.csa_local_index = num_diff_params + i 1946 data.csa_index = self.param_index + i 1947 elif data.param_types[i] == 'tm': 1948 pass 1949 else: 1950 print "Unknown parameter." 1951 return 0 1952 1953 # Increment the parameter index. 1954 self.param_index = self.param_index + data.num_params 1955 1956 # Single residue minimisation with fixed diffusion parameters. 1957 if self.param_set == 'mf': 1958 # Model-free parameters {S2f, S2, ts}. 1959 if data.s2f_index != None and data.tf_index == None and data.s2_index != None and data.ts_index != None: 1960 # Equation. 1961 data.calc_jw_comps = calc_S2f_S2_ts_jw_comps 1962 data.calc_jw = calc_S2f_S2_ts_jw 1963 1964 # Gradient. 1965 data.calc_djw_comps = calc_S2f_S2_ts_djw_comps 1966 data.calc_djw[data.s2f_local_index] = calc_S2f_S2_ts_djw_dS2f 1967 data.calc_djw[data.s2_local_index] = calc_S2f_S2_ts_djw_dS2 1968 data.calc_djw[data.ts_local_index] = calc_S2f_S2_ts_djw_dts 1969 1970 # Hessian. 1971 data.calc_d2jw[data.s2f_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2f_local_index] = calc_S2f_S2_ts_d2jw_dS2fdts 1972 data.calc_d2jw[data.s2_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2_local_index] = calc_S2f_S2_ts_d2jw_dS2dts 1973 data.calc_d2jw[data.ts_local_index][data.ts_local_index] = calc_S2f_S2_ts_d2jw_dts2 1974 1975 # Model-free parameters {S2f, tf, S2, ts}. 1976 elif data.s2f_index != None and data.tf_index != None and data.s2_index != None and data.ts_index != None: 1977 # Equation. 1978 data.calc_jw_comps = calc_S2f_tf_S2_ts_jw_comps 1979 data.calc_jw = calc_S2f_tf_S2_ts_jw 1980 1981 # Gradient. 1982 data.calc_djw_comps = calc_S2f_tf_S2_ts_djw_comps 1983 data.calc_djw[data.s2f_local_index] = calc_S2f_tf_S2_ts_djw_dS2f 1984 data.calc_djw[data.tf_local_index] = calc_S2f_tf_S2_ts_djw_dtf 1985 data.calc_djw[data.s2_local_index] = calc_S2f_S2_ts_djw_dS2 1986 data.calc_djw[data.ts_local_index] = calc_S2f_S2_ts_djw_dts 1987 1988 # Hessian. 1989 data.calc_d2jw[data.s2f_local_index][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][data.s2f_local_index] = calc_S2f_tf_S2_ts_d2jw_dS2fdtf 1990 data.calc_d2jw[data.s2f_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2f_local_index] = calc_S2f_S2_ts_d2jw_dS2fdts 1991 data.calc_d2jw[data.s2_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2_local_index] = calc_S2f_S2_ts_d2jw_dS2dts 1992 data.calc_d2jw[data.tf_local_index][data.tf_local_index] = calc_S2f_tf_S2_ts_d2jw_dtf2 1993 data.calc_d2jw[data.ts_local_index][data.ts_local_index] = calc_S2f_S2_ts_d2jw_dts2 1994 1995 # Bad parameter combination. 1996 else: 1997 print "Invalid combination of parameters for the extended model-free equation." 1998 return 0 1999 2000 # Minimisation with variable diffusion parameters. 2001 else: 2002 # Diffusion parameters and model-free parameters {S2f, S2, ts}. 2003 if data.s2f_index != None and data.tf_index == None and data.s2_index != None and data.ts_index != None: 2004 # Equation. 2005 data.calc_jw_comps = calc_S2f_S2_ts_jw_comps 2006 data.calc_jw = calc_S2f_S2_ts_jw 2007 2008 # Gradient. 2009 data.calc_djw_comps = calc_diff_S2f_S2_ts_djw_comps 2010 2011 if self.param_set != 'diff': 2012 # Gradient. 2013 data.calc_djw[data.s2f_local_index] = calc_S2f_S2_ts_djw_dS2f 2014 data.calc_djw[data.s2_local_index] = calc_S2f_S2_ts_djw_dS2 2015 data.calc_djw[data.ts_local_index] = calc_S2f_S2_ts_djw_dts 2016 2017 # Hessian. 2018 data.calc_d2jw[data.s2f_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2f_local_index] = calc_S2f_S2_ts_d2jw_dS2fdts 2019 data.calc_d2jw[data.s2_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2_local_index] = calc_S2f_S2_ts_d2jw_dS2dts 2020 data.calc_d2jw[data.ts_local_index][data.ts_local_index] = calc_S2f_S2_ts_d2jw_dts2 2021 2022 # Isotropic diffusion. 2023 if self.diff_data.type == 'iso': 2024 # Gradient. 2025 data.calc_djw[0] = calc_diff_S2f_S2_ts_djw_dDj 2026 2027 # Hessian. 2028 data.calc_d2jw[0][0] = calc_diff_S2f_S2_ts_d2jw_dDjdDk 2029 if self.param_set != 'diff': 2030 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2f 2031 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2032 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2033 2034 # Axially symmetric diffusion. 2035 elif self.diff_data.type == 'axial': 2036 # Gradient. 2037 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2f_S2_ts_djw_dDj 2038 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2f_S2_ts_djw_dPsij 2039 2040 # Hessian. 2041 data.calc_d2jw[0][0] = calc_diff_S2f_S2_ts_d2jw_dDjdDk 2042 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2f_S2_ts_d2jw_dDjdDk 2043 data.calc_d2jw[1][1] = calc_diff_S2f_S2_ts_d2jw_dDjdDk 2044 2045 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2f_S2_ts_d2jw_dDjdPsij 2046 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2f_S2_ts_d2jw_dDjdPsij 2047 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2f_S2_ts_d2jw_dDjdPsij 2048 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2f_S2_ts_d2jw_dDjdPsij 2049 2050 data.calc_d2jw[2][2] = calc_diff_S2f_S2_ts_d2jw_dPsijdPsik 2051 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2f_S2_ts_d2jw_dPsijdPsik 2052 data.calc_d2jw[3][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdPsik 2053 2054 if self.param_set != 'diff': 2055 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2f 2056 data.calc_d2jw[1][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][1] = calc_diff_S2f_S2_ts_d2jw_dDjdS2f 2057 data.calc_d2jw[2][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][2] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2f 2058 data.calc_d2jw[3][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2f 2059 2060 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2061 data.calc_d2jw[1][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][1] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2062 data.calc_d2jw[2][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][2] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2063 data.calc_d2jw[3][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2064 2065 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2066 data.calc_d2jw[1][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][1] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2067 data.calc_d2jw[2][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][2] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2068 data.calc_d2jw[3][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2069 2070 # Anisotropic diffusion. 2071 elif self.diff_data.type == 'aniso': 2072 # Gradient. 2073 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_aniso_S2f_S2_ts_djw_dDj 2074 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2f_S2_ts_djw_dPsij 2075 2076 # Hessian. 2077 data.calc_d2jw[0][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdDk 2078 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdDk 2079 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdDk 2080 data.calc_d2jw[1][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdDk 2081 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdDk 2082 data.calc_d2jw[2][2] = calc_aniso_S2f_S2_ts_d2jw_dDjdDk 2083 2084 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdPsij 2085 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdPsij 2086 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdPsij 2087 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdPsij 2088 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdPsij 2089 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdPsij 2090 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_aniso_S2f_S2_ts_d2jw_dDjdPsij 2091 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_aniso_S2f_S2_ts_d2jw_dDjdPsij 2092 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_aniso_S2f_S2_ts_d2jw_dDjdPsij 2093 2094 data.calc_d2jw[3][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdPsik 2095 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdPsik 2096 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdPsik 2097 data.calc_d2jw[4][4] = calc_diff_S2f_S2_ts_d2jw_dPsijdPsik 2098 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2f_S2_ts_d2jw_dPsijdPsik 2099 data.calc_d2jw[5][5] = calc_diff_S2f_S2_ts_d2jw_dPsijdPsik 2100 2101 if self.param_set != 'diff': 2102 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdS2f 2103 data.calc_d2jw[1][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdS2f 2104 data.calc_d2jw[2][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][2] = calc_aniso_S2f_S2_ts_d2jw_dDjdS2f 2105 data.calc_d2jw[3][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2f 2106 data.calc_d2jw[4][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][4] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2f 2107 data.calc_d2jw[5][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][5] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2f 2108 2109 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdS2 2110 data.calc_d2jw[1][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdS2 2111 data.calc_d2jw[2][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][2] = calc_aniso_S2f_S2_ts_d2jw_dDjdS2 2112 data.calc_d2jw[3][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2113 data.calc_d2jw[4][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][4] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2114 data.calc_d2jw[5][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][5] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2115 2116 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdts 2117 data.calc_d2jw[1][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdts 2118 data.calc_d2jw[2][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][2] = calc_aniso_S2f_S2_ts_d2jw_dDjdts 2119 data.calc_d2jw[3][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2120 data.calc_d2jw[4][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][4] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2121 data.calc_d2jw[5][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][5] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2122 2123 2124 # Diffusion parameters and model-free parameters {S2f, tf, S2, ts}. 2125 elif data.s2f_index != None and data.tf_index != None and data.s2_index != None and data.ts_index != None: 2126 # Equation. 2127 data.calc_jw_comps = calc_S2f_tf_S2_ts_jw_comps 2128 data.calc_jw = calc_S2f_tf_S2_ts_jw 2129 2130 # Gradient. 2131 data.calc_djw_comps = calc_diff_S2f_tf_S2_ts_djw_comps 2132 2133 if self.param_set != 'diff': 2134 # Gradient. 2135 data.calc_djw[data.s2f_local_index] = calc_S2f_tf_S2_ts_djw_dS2f 2136 data.calc_djw[data.tf_local_index] = calc_S2f_tf_S2_ts_djw_dtf 2137 data.calc_djw[data.s2_local_index] = calc_S2f_S2_ts_djw_dS2 2138 data.calc_djw[data.ts_local_index] = calc_S2f_S2_ts_djw_dts 2139 2140 # Hessian. 2141 data.calc_d2jw[data.s2f_local_index][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][data.s2f_local_index] = calc_S2f_tf_S2_ts_d2jw_dS2fdtf 2142 data.calc_d2jw[data.s2f_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2f_local_index] = calc_S2f_S2_ts_d2jw_dS2fdts 2143 data.calc_d2jw[data.tf_local_index][data.tf_local_index] = calc_S2f_tf_S2_ts_d2jw_dtf2 2144 data.calc_d2jw[data.s2_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2_local_index] = calc_S2f_S2_ts_d2jw_dS2dts 2145 data.calc_d2jw[data.ts_local_index][data.ts_local_index] = calc_S2f_S2_ts_d2jw_dts2 2146 2147 # Isotropic diffusion. 2148 if self.diff_data.type == 'iso': 2149 # Gradient. 2150 data.calc_djw[0] = calc_diff_S2f_tf_S2_ts_djw_dDj 2151 2152 # Hessian. 2153 data.calc_d2jw[0][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdDk 2154 if self.param_set != 'diff': 2155 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2f 2156 data.calc_d2jw[0][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdtf 2157 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2158 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2159 2160 # Axially symmetric diffusion. 2161 elif self.diff_data.type == 'axial': 2162 # Gradient. 2163 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2f_tf_S2_ts_djw_dDj 2164 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2f_tf_S2_ts_djw_dPsij 2165 2166 # Hessian. 2167 data.calc_d2jw[0][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdDk 2168 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdDk 2169 data.calc_d2jw[1][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdDk 2170 2171 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdPsij 2172 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdPsij 2173 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdPsij 2174 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdPsij 2175 2176 data.calc_d2jw[2][2] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdPsik 2177 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdPsik 2178 data.calc_d2jw[3][3] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdPsik 2179 2180 if self.param_set != 'diff': 2181 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdS2f 2182 data.calc_d2jw[1][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdS2f 2183 data.calc_d2jw[2][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][2] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdS2f 2184 data.calc_d2jw[3][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][3] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdS2f 2185 2186 data.calc_d2jw[0][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][0] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdtf 2187 data.calc_d2jw[1][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][1] = calc_diff_S2f_tf_S2_ts_d2jw_dDjdtf 2188 data.calc_d2jw[2][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][2] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdtf 2189 data.calc_d2jw[3][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][3] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdtf 2190 2191 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2192 data.calc_d2jw[1][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][1] = calc_diff_S2f_S2_ts_d2jw_dDjdS2 2193 data.calc_d2jw[2][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][2] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2194 data.calc_d2jw[3][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2195 2196 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2197 data.calc_d2jw[1][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][1] = calc_diff_S2f_S2_ts_d2jw_dDjdts 2198 data.calc_d2jw[2][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][2] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2199 data.calc_d2jw[3][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2200 2201 # Anisotropic diffusion. 2202 elif self.diff_data.type == 'aniso': 2203 # Gradient. 2204 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_aniso_S2f_tf_S2_ts_djw_dDj 2205 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2f_tf_S2_ts_djw_dPsij 2206 2207 # Hessian. 2208 data.calc_d2jw[0][0] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdDk 2209 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdDk 2210 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdDk 2211 data.calc_d2jw[1][1] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdDk 2212 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdDk 2213 data.calc_d2jw[2][2] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdDk 2214 2215 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdPsij 2216 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdPsij 2217 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdPsij 2218 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdPsij 2219 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdPsij 2220 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdPsij 2221 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdPsij 2222 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdPsij 2223 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdPsij 2224 2225 data.calc_d2jw[3][3] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdPsik 2226 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdPsik 2227 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdPsik 2228 data.calc_d2jw[4][4] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdPsik 2229 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdPsik 2230 data.calc_d2jw[5][5] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdPsik 2231 2232 if self.param_set != 'diff': 2233 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdS2f 2234 data.calc_d2jw[1][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][1] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdS2f 2235 data.calc_d2jw[2][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][2] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdS2f 2236 data.calc_d2jw[3][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][3] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdS2f 2237 data.calc_d2jw[4][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][4] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdS2f 2238 data.calc_d2jw[5][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][5] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdS2f 2239 2240 data.calc_d2jw[0][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][0] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdtf 2241 data.calc_d2jw[1][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][1] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdtf 2242 data.calc_d2jw[2][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][2] = calc_aniso_S2f_tf_S2_ts_d2jw_dDjdtf 2243 data.calc_d2jw[3][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][3] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdtf 2244 data.calc_d2jw[4][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][4] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdtf 2245 data.calc_d2jw[5][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][5] = calc_diff_S2f_tf_S2_ts_d2jw_dPsijdtf 2246 2247 data.calc_d2jw[0][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdS2 2248 data.calc_d2jw[1][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdS2 2249 data.calc_d2jw[2][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][2] = calc_aniso_S2f_S2_ts_d2jw_dDjdS2 2250 data.calc_d2jw[3][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2251 data.calc_d2jw[4][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][4] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2252 data.calc_d2jw[5][data.s2_local_index] = data.calc_d2jw[data.s2_local_index][5] = calc_diff_S2f_S2_ts_d2jw_dPsijdS2 2253 2254 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_aniso_S2f_S2_ts_d2jw_dDjdts 2255 data.calc_d2jw[1][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][1] = calc_aniso_S2f_S2_ts_d2jw_dDjdts 2256 data.calc_d2jw[2][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][2] = calc_aniso_S2f_S2_ts_d2jw_dDjdts 2257 data.calc_d2jw[3][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][3] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2258 data.calc_d2jw[4][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][4] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2259 data.calc_d2jw[5][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][5] = calc_diff_S2f_S2_ts_d2jw_dPsijdts 2260 2261 # Bad parameter combination. 2262 else: 2263 print "Invalid combination of parameters for the extended model-free equation." 2264 return 0 2265 2266 2267 # The extended 2 model-free equations {tm, S2f, tf, S2s, ts, Rex, r, CSA}. 2268 ######################################################################### 2269 2270 elif data.equations == 'mf_ext2': 2271 # Find the indecies of the model-free parameters. 2272 for i in xrange(data.num_params): 2273 if data.param_types[i] == 'S2f': 2274 data.s2f_local_index = num_diff_params + i 2275 data.s2f_index = self.param_index + i 2276 elif data.param_types[i] == 'tf': 2277 data.tf_local_index = num_diff_params + i 2278 data.tf_index = self.param_index + i 2279 elif data.param_types[i] == 'S2s': 2280 data.s2s_local_index = num_diff_params + i 2281 data.s2s_index = self.param_index + i 2282 elif data.param_types[i] == 'ts': 2283 data.ts_local_index = num_diff_params + i 2284 data.ts_index = self.param_index + i 2285 elif data.param_types[i] == 'Rex': 2286 data.rex_local_index = num_diff_params + i 2287 data.rex_index = self.param_index + i 2288 elif data.param_types[i] == 'r': 2289 data.r_local_index = num_diff_params + i 2290 data.r_index = self.param_index + i 2291 elif data.param_types[i] == 'CSA': 2292 data.csa_local_index = num_diff_params + i 2293 data.csa_index = self.param_index + i 2294 elif data.param_types[i] == 'tm': 2295 pass 2296 else: 2297 print "Unknown parameter." 2298 return 0 2299 2300 # Increment the parameter index. 2301 self.param_index = self.param_index + data.num_params 2302 2303 # Single residue minimisation with fixed diffusion parameters. 2304 if self.param_set == 'mf': 2305 # Model-free parameters {S2f, S2s, ts}. 2306 if data.s2f_index != None and data.tf_index == None and data.s2s_index != None and data.ts_index != None: 2307 # Equation. 2308 data.calc_jw_comps = calc_S2f_S2s_ts_jw_comps 2309 data.calc_jw = calc_S2f_S2s_ts_jw 2310 2311 # Gradient. 2312 data.calc_djw_comps = calc_S2f_S2s_ts_djw_comps 2313 data.calc_djw[data.s2f_local_index] = calc_S2f_S2s_ts_djw_dS2f 2314 data.calc_djw[data.s2s_local_index] = calc_S2f_S2s_ts_djw_dS2s 2315 data.calc_djw[data.ts_local_index] = calc_S2f_S2s_ts_djw_dts 2316 2317 # Hessian. 2318 data.calc_d2jw[data.s2f_local_index][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][data.s2f_local_index] = calc_S2f_S2s_ts_d2jw_dS2fdS2s 2319 data.calc_d2jw[data.s2f_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2f_local_index] = calc_S2f_S2s_ts_d2jw_dS2fdts 2320 data.calc_d2jw[data.s2s_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2s_local_index] = calc_S2f_S2s_ts_d2jw_dS2sdts 2321 data.calc_d2jw[data.ts_local_index][data.ts_local_index] = calc_S2f_S2s_ts_d2jw_dts2 2322 2323 # Model-free parameters {S2f, tf, S2s, ts}. 2324 elif data.s2f_index != None and data.tf_index != None and data.s2s_index != None and data.ts_index != None: 2325 # Equation. 2326 data.calc_jw_comps = calc_S2f_tf_S2s_ts_jw_comps 2327 data.calc_jw = calc_S2f_tf_S2s_ts_jw 2328 2329 # Gradient. 2330 data.calc_djw_comps = calc_S2f_tf_S2s_ts_djw_comps 2331 data.calc_djw[data.s2f_local_index] = calc_S2f_tf_S2s_ts_djw_dS2f 2332 data.calc_djw[data.tf_local_index] = calc_S2f_tf_S2s_ts_djw_dtf 2333 data.calc_djw[data.s2s_local_index] = calc_S2f_tf_S2s_ts_djw_dS2s 2334 data.calc_djw[data.ts_local_index] = calc_S2f_tf_S2s_ts_djw_dts 2335 2336 # Hessian. 2337 data.calc_d2jw[data.s2f_local_index][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][data.s2f_local_index] = calc_S2f_S2s_ts_d2jw_dS2fdS2s 2338 data.calc_d2jw[data.s2f_local_index][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][data.s2f_local_index] = calc_S2f_tf_S2s_ts_d2jw_dS2fdtf 2339 data.calc_d2jw[data.s2f_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2f_local_index] = calc_S2f_tf_S2s_ts_d2jw_dS2fdts 2340 data.calc_d2jw[data.s2s_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2s_local_index] = calc_S2f_tf_S2s_ts_d2jw_dS2sdts 2341 data.calc_d2jw[data.tf_local_index][data.tf_local_index] = calc_S2f_tf_S2s_ts_d2jw_dtf2 2342 data.calc_d2jw[data.ts_local_index][data.ts_local_index] = calc_S2f_tf_S2s_ts_d2jw_dts2 2343 2344 # Bad parameter combination. 2345 else: 2346 print "Invalid combination of parameters for the extended model-free equation." 2347 return 0 2348 2349 # Minimisation with variable diffusion parameters. 2350 else: 2351 # Diffusion parameters and model-free parameters {S2f, S2s, ts}. 2352 if data.s2f_index != None and data.tf_index == None and data.s2s_index != None and data.ts_index != None: 2353 # Equation. 2354 data.calc_jw_comps = calc_diff_S2f_S2s_ts_jw_comps 2355 data.calc_jw = calc_S2f_S2s_ts_jw 2356 2357 # Gradient. 2358 data.calc_djw_comps = calc_diff_S2f_S2s_ts_djw_comps 2359 2360 if self.param_set != 'diff': 2361 # Gradient. 2362 data.calc_djw[data.s2f_local_index] = calc_diff_S2f_S2s_ts_djw_dS2f 2363 data.calc_djw[data.s2s_local_index] = calc_diff_S2f_S2s_ts_djw_dS2s 2364 data.calc_djw[data.ts_local_index] = calc_diff_S2f_S2s_ts_djw_dts 2365 2366 # Hessian. 2367 data.calc_d2jw[data.s2f_local_index][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][data.s2f_local_index] = calc_S2f_S2s_ts_d2jw_dS2fdS2s 2368 data.calc_d2jw[data.s2f_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2f_local_index] = calc_diff_S2f_S2s_ts_d2jw_dS2fdts 2369 data.calc_d2jw[data.s2s_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2s_local_index] = calc_diff_S2f_S2s_ts_d2jw_dS2sdts 2370 data.calc_d2jw[data.ts_local_index][data.ts_local_index] = calc_diff_S2f_S2s_ts_d2jw_dts2 2371 2372 # Isotropic diffusion. 2373 if self.diff_data.type == 'iso': 2374 # Gradient. 2375 data.calc_djw[0] = calc_diff_S2f_S2s_ts_djw_dDj 2376 2377 # Hessian. 2378 data.calc_d2jw[0][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdDk 2379 if self.param_set != 'diff': 2380 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2f 2381 data.calc_d2jw[0][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2382 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2383 2384 # Axially symmetric diffusion. 2385 elif self.diff_data.type == 'axial': 2386 # Gradient. 2387 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2f_S2s_ts_djw_dDj 2388 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2f_S2s_ts_djw_dPsij 2389 2390 # Hessian. 2391 data.calc_d2jw[0][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdDk 2392 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdDk 2393 data.calc_d2jw[1][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdDk 2394 2395 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdPsij 2396 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdPsij 2397 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdPsij 2398 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdPsij 2399 2400 data.calc_d2jw[2][2] = calc_diff_S2f_S2s_ts_d2jw_dPsijdPsik 2401 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2f_S2s_ts_d2jw_dPsijdPsik 2402 data.calc_d2jw[3][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdPsik 2403 2404 if self.param_set != 'diff': 2405 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2f 2406 data.calc_d2jw[1][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2f 2407 data.calc_d2jw[2][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][2] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2f 2408 data.calc_d2jw[3][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2f 2409 2410 data.calc_d2jw[0][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2411 data.calc_d2jw[1][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2412 data.calc_d2jw[2][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][2] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2413 data.calc_d2jw[3][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2414 2415 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2416 data.calc_d2jw[1][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2417 data.calc_d2jw[2][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][2] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2418 data.calc_d2jw[3][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2419 2420 # Anisotropic diffusion. 2421 elif self.diff_data.type == 'aniso': 2422 # Gradient. 2423 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_aniso_S2f_S2s_ts_djw_dDj 2424 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2f_S2s_ts_djw_dPsij 2425 2426 # Hessian. 2427 data.calc_d2jw[0][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdDk 2428 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdDk 2429 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdDk 2430 data.calc_d2jw[1][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdDk 2431 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdDk 2432 data.calc_d2jw[2][2] = calc_aniso_S2f_S2s_ts_d2jw_dDjdDk 2433 2434 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdPsij 2435 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdPsij 2436 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdPsij 2437 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdPsij 2438 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdPsij 2439 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdPsij 2440 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_aniso_S2f_S2s_ts_d2jw_dDjdPsij 2441 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_aniso_S2f_S2s_ts_d2jw_dDjdPsij 2442 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_aniso_S2f_S2s_ts_d2jw_dDjdPsij 2443 2444 data.calc_d2jw[3][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdPsik 2445 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdPsik 2446 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdPsik 2447 data.calc_d2jw[4][4] = calc_diff_S2f_S2s_ts_d2jw_dPsijdPsik 2448 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2f_S2s_ts_d2jw_dPsijdPsik 2449 data.calc_d2jw[5][5] = calc_diff_S2f_S2s_ts_d2jw_dPsijdPsik 2450 2451 if self.param_set != 'diff': 2452 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdS2f 2453 data.calc_d2jw[1][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdS2f 2454 data.calc_d2jw[2][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][2] = calc_aniso_S2f_S2s_ts_d2jw_dDjdS2f 2455 data.calc_d2jw[3][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2f 2456 data.calc_d2jw[4][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][4] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2f 2457 data.calc_d2jw[5][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][5] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2f 2458 2459 data.calc_d2jw[0][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdS2s 2460 data.calc_d2jw[1][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdS2s 2461 data.calc_d2jw[2][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][2] = calc_aniso_S2f_S2s_ts_d2jw_dDjdS2s 2462 data.calc_d2jw[3][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2463 data.calc_d2jw[4][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][4] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2464 data.calc_d2jw[5][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][5] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2465 2466 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdts 2467 data.calc_d2jw[1][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdts 2468 data.calc_d2jw[2][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][2] = calc_aniso_S2f_S2s_ts_d2jw_dDjdts 2469 data.calc_d2jw[3][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2470 data.calc_d2jw[4][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][4] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2471 data.calc_d2jw[5][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][5] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2472 2473 # Diffusion parameters and model-free parameters {S2f, tf, S2s, ts}. 2474 elif data.s2f_index != None and data.tf_index != None and data.s2s_index != None and data.ts_index != None: 2475 # Equation. 2476 data.calc_jw_comps = calc_diff_S2f_tf_S2s_ts_jw_comps 2477 data.calc_jw = calc_S2f_tf_S2s_ts_jw 2478 2479 # Gradient. 2480 data.calc_djw_comps = calc_diff_S2f_tf_S2s_ts_djw_comps 2481 2482 if self.param_set != 'diff': 2483 # Gradient. 2484 data.calc_djw[data.s2f_local_index] = calc_diff_S2f_tf_S2s_ts_djw_dS2f 2485 data.calc_djw[data.tf_local_index] = calc_diff_S2f_tf_S2s_ts_djw_dtf 2486 data.calc_djw[data.s2s_local_index] = calc_diff_S2f_tf_S2s_ts_djw_dS2s 2487 data.calc_djw[data.ts_local_index] = calc_diff_S2f_tf_S2s_ts_djw_dts 2488 2489 # Hessian. 2490 data.calc_d2jw[data.s2f_local_index][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][data.s2f_local_index] = calc_S2f_S2s_ts_d2jw_dS2fdS2s 2491 data.calc_d2jw[data.s2f_local_index][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][data.s2f_local_index] = calc_diff_S2f_tf_S2s_ts_d2jw_dS2fdtf 2492 data.calc_d2jw[data.s2f_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2f_local_index] = calc_diff_S2f_tf_S2s_ts_d2jw_dS2fdts 2493 data.calc_d2jw[data.tf_local_index][data.tf_local_index] = calc_diff_S2f_tf_S2s_ts_d2jw_dtf2 2494 data.calc_d2jw[data.s2s_local_index][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][data.s2s_local_index] = calc_diff_S2f_tf_S2s_ts_d2jw_dS2sdts 2495 data.calc_d2jw[data.ts_local_index][data.ts_local_index] = calc_diff_S2f_tf_S2s_ts_d2jw_dts2 2496 2497 # Isotropic diffusion. 2498 if self.diff_data.type == 'iso': 2499 # Gradient. 2500 data.calc_djw[0] = calc_diff_S2f_tf_S2s_ts_djw_dDj 2501 2502 # Hessian. 2503 data.calc_d2jw[0][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdDk 2504 if self.param_set != 'diff': 2505 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdS2f 2506 data.calc_d2jw[0][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdtf 2507 data.calc_d2jw[0][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdS2s 2508 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdts 2509 2510 # Axially symmetric diffusion. 2511 elif self.diff_data.type == 'axial': 2512 # Gradient. 2513 data.calc_djw[0] = data.calc_djw[1] = calc_diff_S2f_tf_S2s_ts_djw_dDj 2514 data.calc_djw[2] = data.calc_djw[3] = calc_diff_S2f_tf_S2s_ts_djw_dPsij 2515 2516 # Hessian. 2517 data.calc_d2jw[0][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdDk 2518 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdDk 2519 data.calc_d2jw[1][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdDk 2520 2521 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdPsij 2522 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdPsij 2523 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdPsij 2524 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdPsij 2525 2526 data.calc_d2jw[2][2] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdPsik 2527 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdPsik 2528 data.calc_d2jw[3][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdPsik 2529 2530 if self.param_set != 'diff': 2531 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdS2f 2532 data.calc_d2jw[1][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdS2f 2533 data.calc_d2jw[2][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][2] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdS2f 2534 data.calc_d2jw[3][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdS2f 2535 2536 data.calc_d2jw[0][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][0] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdtf 2537 data.calc_d2jw[1][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][1] = calc_diff_S2f_tf_S2s_ts_d2jw_dDjdtf 2538 data.calc_d2jw[2][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][2] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdtf 2539 data.calc_d2jw[3][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdtf 2540 2541 data.calc_d2jw[0][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2542 data.calc_d2jw[1][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdS2s 2543 data.calc_d2jw[2][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][2] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2544 data.calc_d2jw[3][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2545 2546 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2547 data.calc_d2jw[1][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][1] = calc_diff_S2f_S2s_ts_d2jw_dDjdts 2548 data.calc_d2jw[2][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][2] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2549 data.calc_d2jw[3][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2550 2551 # Anisotropic diffusion. 2552 elif self.diff_data.type == 'aniso': 2553 # Gradient. 2554 data.calc_djw[0] = data.calc_djw[1] = data.calc_djw[2] = calc_aniso_S2f_tf_S2s_ts_djw_dDj 2555 data.calc_djw[3] = data.calc_djw[4] = data.calc_djw[5] = calc_diff_S2f_tf_S2s_ts_djw_dPsij 2556 2557 # Hessian. 2558 data.calc_d2jw[0][0] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdDk 2559 data.calc_d2jw[0][1] = data.calc_d2jw[1][0] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdDk 2560 data.calc_d2jw[0][2] = data.calc_d2jw[2][0] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdDk 2561 data.calc_d2jw[1][1] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdDk 2562 data.calc_d2jw[1][2] = data.calc_d2jw[2][1] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdDk 2563 data.calc_d2jw[2][2] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdDk 2564 2565 data.calc_d2jw[0][3] = data.calc_d2jw[3][0] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdPsij 2566 data.calc_d2jw[0][4] = data.calc_d2jw[4][0] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdPsij 2567 data.calc_d2jw[0][5] = data.calc_d2jw[5][0] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdPsij 2568 data.calc_d2jw[1][3] = data.calc_d2jw[3][1] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdPsij 2569 data.calc_d2jw[1][4] = data.calc_d2jw[4][1] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdPsij 2570 data.calc_d2jw[1][5] = data.calc_d2jw[5][1] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdPsij 2571 data.calc_d2jw[2][3] = data.calc_d2jw[3][2] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdPsij 2572 data.calc_d2jw[2][4] = data.calc_d2jw[4][2] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdPsij 2573 data.calc_d2jw[2][5] = data.calc_d2jw[5][2] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdPsij 2574 2575 data.calc_d2jw[3][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdPsik 2576 data.calc_d2jw[3][4] = data.calc_d2jw[4][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdPsik 2577 data.calc_d2jw[3][5] = data.calc_d2jw[5][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdPsik 2578 data.calc_d2jw[4][4] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdPsik 2579 data.calc_d2jw[4][5] = data.calc_d2jw[5][4] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdPsik 2580 data.calc_d2jw[5][5] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdPsik 2581 2582 if self.param_set != 'diff': 2583 data.calc_d2jw[0][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][0] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdS2f 2584 data.calc_d2jw[1][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][1] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdS2f 2585 data.calc_d2jw[2][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][2] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdS2f 2586 data.calc_d2jw[3][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdS2f 2587 data.calc_d2jw[4][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][4] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdS2f 2588 data.calc_d2jw[5][data.s2f_local_index] = data.calc_d2jw[data.s2f_local_index][5] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdS2f 2589 2590 data.calc_d2jw[0][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][0] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdtf 2591 data.calc_d2jw[1][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][1] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdtf 2592 data.calc_d2jw[2][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][2] = calc_aniso_S2f_tf_S2s_ts_d2jw_dDjdtf 2593 data.calc_d2jw[3][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][3] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdtf 2594 data.calc_d2jw[4][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][4] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdtf 2595 data.calc_d2jw[5][data.tf_local_index] = data.calc_d2jw[data.tf_local_index][5] = calc_diff_S2f_tf_S2s_ts_d2jw_dPsijdtf 2596 2597 data.calc_d2jw[0][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdS2s 2598 data.calc_d2jw[1][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdS2s 2599 data.calc_d2jw[2][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][2] = calc_aniso_S2f_S2s_ts_d2jw_dDjdS2s 2600 data.calc_d2jw[3][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2601 data.calc_d2jw[4][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][4] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2602 data.calc_d2jw[5][data.s2s_local_index] = data.calc_d2jw[data.s2s_local_index][5] = calc_diff_S2f_S2s_ts_d2jw_dPsijdS2s 2603 2604 data.calc_d2jw[0][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][0] = calc_aniso_S2f_S2s_ts_d2jw_dDjdts 2605 data.calc_d2jw[1][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][1] = calc_aniso_S2f_S2s_ts_d2jw_dDjdts 2606 data.calc_d2jw[2][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][2] = calc_aniso_S2f_S2s_ts_d2jw_dDjdts 2607 data.calc_d2jw[3][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][3] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2608 data.calc_d2jw[4][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][4] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2609 data.calc_d2jw[5][data.ts_local_index] = data.calc_d2jw[data.ts_local_index][5] = calc_diff_S2f_S2s_ts_d2jw_dPsijdts 2610 2611 # Bad parameter combination. 2612 else: 2613 print "Invalid combination of parameters for the extended model-free equation." 2614 return 0 2615 2616 # Unknown model-free equation. 2617 else: 2618 print "Unknown model-free equation." 2619 return 0 2620 2621 2622 # Initialise function data structures. 2623 ###################################### 2624 2625 # Relaxation equation components. 2626 data.create_dip_func, data.create_dip_grad, data.create_dip_hess = [], [], [] 2627 data.create_csa_func, data.create_csa_grad, data.create_csa_hess = [], [], [] 2628 data.create_rex_func, data.create_rex_grad = [], [] 2629 2630 data.create_dip_jw_func, data.create_dip_jw_grad, data.create_dip_jw_hess = [], [], [] 2631 data.create_csa_jw_func, data.create_csa_jw_grad, data.create_csa_jw_hess = [], [], [] 2632 2633 # Ri' 2634 data.create_ri_prime = None 2635 data.create_dri_prime = [] 2636 data.create_d2ri_prime = [] 2637 2638 # Ri 2639 data.create_ri, data.create_dri, data.create_d2ri = [], [], [] 2640 data.get_r1, data.get_dr1, data.get_d2r1 = [], [], [] 2641 2642 # Fill the structures with None. 2643 for i in xrange(data.num_ri): 2644 data.create_dip_func.append(None) 2645 data.create_dip_grad.append(None) 2646 data.create_dip_hess.append(None) 2647 data.create_csa_func.append(None) 2648 data.create_csa_grad.append(None) 2649 data.create_csa_hess.append(None) 2650 data.create_rex_func.append(None) 2651 data.create_rex_grad.append(None) 2652 data.create_dip_jw_func.append(None) 2653 data.create_dip_jw_grad.append(None) 2654 data.create_dip_jw_hess.append(None) 2655 data.create_csa_jw_func.append(None) 2656 data.create_csa_jw_grad.append(None) 2657 data.create_csa_jw_hess.append(None) 2658 data.create_ri.append(None) 2659 data.create_dri.append(None) 2660 data.create_d2ri.append(None) 2661 data.get_r1.append(None) 2662 data.get_dr1.append(None) 2663 data.get_d2r1.append(None) 2664 2665 2666 # Select the functions for the calculation of ri_prime, dri_prime, and d2ri_prime components. 2667 ############################################################################################# 2668 2669 for i in xrange(data.num_ri): 2670 # The R1 equations. 2671 if data.ri_labels[i] == 'R1': 2672 data.create_csa_func[i] = comp_r1_csa_const 2673 data.create_csa_grad[i] = comp_r1_csa_const 2674 data.create_csa_hess[i] = comp_r1_csa_const 2675 data.create_dip_jw_func[i] = comp_r1_dip_jw 2676 data.create_dip_jw_grad[i] = comp_r1_dip_jw 2677 data.create_dip_jw_hess[i] = comp_r1_dip_jw 2678 data.create_csa_jw_func[i] = comp_r1_csa_jw 2679 data.create_csa_jw_grad[i] = comp_r1_csa_jw 2680 data.create_csa_jw_hess[i] = comp_r1_csa_jw 2681 2682 # The R2 equations. 2683 elif data.ri_labels[i] == 'R2': 2684 data.create_dip_func[i] = comp_r2_dip_const 2685 data.create_dip_grad[i] = comp_r2_dip_const 2686 data.create_dip_hess[i] = comp_r2_dip_const 2687 data.create_csa_func[i] = comp_r2_csa_const 2688 data.create_csa_grad[i] = comp_r2_csa_const 2689 data.create_csa_hess[i] = comp_r2_csa_const 2690 data.create_rex_func[i] = comp_rex_const_func 2691 data.create_rex_grad[i] = comp_rex_const_grad 2692 data.create_dip_jw_func[i] = comp_r2_dip_jw 2693 data.create_dip_jw_grad[i] = comp_r2_dip_jw 2694 data.create_dip_jw_hess[i] = comp_r2_dip_jw 2695 data.create_csa_jw_func[i] = comp_r2_csa_jw 2696 data.create_csa_jw_grad[i] = comp_r2_csa_jw 2697 data.create_csa_jw_hess[i] = comp_r2_csa_jw 2698 2699 # The NOE equations. 2700 elif data.ri_labels[i] == 'NOE': 2701 data.create_dip_jw_func[i] = comp_sigma_noe_dip_jw 2702 data.create_dip_jw_grad[i] = comp_sigma_noe_dip_jw 2703 data.create_dip_jw_hess[i] = comp_sigma_noe_dip_jw 2704 data.create_ri[i] = calc_noe 2705 data.create_dri[i] = calc_dnoe 2706 data.create_d2ri[i] = calc_d2noe 2707 if data.noe_r1_table[i] == None: 2708 data.get_r1[i] = calc_r1 2709 data.get_dr1[i] = calc_dr1 2710 data.get_d2r1[i] = calc_d2r1 2711 else: 2712 data.get_r1[i] = extract_r1 2713 data.get_dr1[i] = extract_dr1 2714 data.get_d2r1[i] = extract_d2r1 2715 2716 2717 # Select the functions for the calculation of ri_prime, dri_prime, and d2ri_prime. 2718 ################################################################################## 2719 2720 # ri_prime. 2721 if data.rex_index == None: 2722 data.create_ri_prime = func_ri_prime 2723 else: 2724 data.create_ri_prime = func_ri_prime_rex 2725 2726 # dri_prime and d2ri_prime. 2727 for i in xrange(data.total_num_params): 2728 # Diffusion tensor parameters are the only parameters. 2729 if self.param_set == 'diff': 2730 # Gradient. 2731 data.create_dri_prime.append(func_dri_djw_prime) 2732 2733 # Hessian. 2734 data.create_d2ri_prime.append([]) 2735 for j in xrange(data.total_num_params): 2736 data.create_d2ri_prime[i].append(func_d2ri_djwidjwj_prime) 2737 2738 # Skip to the next parameter index. 2739 continue 2740 2741 # Residue specific parameter index. 2742 index = i - num_diff_params 2743 if index < 0: 2744 index = None 2745 2746 # Rex. 2747 if index != None and data.param_types[index] == 'Rex': 2748 # Gradient. 2749 data.create_dri_prime.append(func_dri_drex_prime) 2750 2751 # Hessian. 2752 data.create_d2ri_prime.append([]) 2753 for j in xrange(data.total_num_params): 2754 # Residue specific parameter index. 2755 index2 = j - num_diff_params 2756 if index2 < 0: 2757 index2 = None 2758 2759 # Rex. 2760 if index2 != None and data.param_types[index2] == 'Rex': 2761 data.create_d2ri_prime[i].append(None) 2762 2763 # Bond length. 2764 elif index2 != None and data.param_types[index2] == 'r': 2765 data.create_d2ri_prime[i].append(None) 2766 2767 # CSA. 2768 elif index2 != None and data.param_types[index2] == 'CSA': 2769 data.create_d2ri_prime[i].append(None) 2770 2771 # Any other parameter. 2772 else: 2773 data.create_d2ri_prime[i].append(None) 2774 2775 # Bond length. 2776 elif index != None and data.param_types[index] == 'r': 2777 # Gradient. 2778 data.create_dri_prime.append(func_dri_dr_prime) 2779 2780 # Hessian. 2781 data.create_d2ri_prime.append([]) 2782 for j in xrange(data.total_num_params): 2783 # Residue specific parameter index. 2784 index2 = j - num_diff_params 2785 if index2 < 0: 2786 index2 = None 2787 2788 # Rex. 2789 if index2 != None and data.param_types[index2] == 'Rex': 2790 data.create_d2ri_prime[i].append(None) 2791 2792 # Bond length. 2793 elif index2 != None and data.param_types[index2] == 'r': 2794 data.create_d2ri_prime[i].append(func_d2ri_dr2_prime) 2795 2796 # CSA. 2797 elif index2 != None and data.param_types[index2] == 'CSA': 2798 data.create_d2ri_prime[i].append(None) 2799 2800 # Any other parameter. 2801 else: 2802 data.create_d2ri_prime[i].append(func_d2ri_drdjw_prime) 2803 2804 # CSA. 2805 elif index != None and data.param_types[index] == 'CSA': 2806 # Gradient. 2807 data.create_dri_prime.append(func_dri_dcsa_prime) 2808 2809 # Hessian. 2810 data.create_d2ri_prime.append([]) 2811 for j in xrange(data.total_num_params): 2812 # Residue specific parameter index. 2813 index2 = j - num_diff_params 2814 if index2 < 0: 2815 index2 = None 2816 2817 # Rex. 2818 if index2 != None and data.param_types[index2] == 'Rex': 2819 data.create_d2ri_prime[i].append(None) 2820 2821 # Bond length. 2822 elif index2 != None and data.param_types[index2] == 'r': 2823 data.create_d2ri_prime[i].append(None) 2824 2825 # CSA. 2826 elif index2 != None and data.param_types[index2] == 'CSA': 2827 data.create_d2ri_prime[i].append(func_d2ri_dcsa2_prime) 2828 2829 # Any other parameter. 2830 else: 2831 data.create_d2ri_prime[i].append(func_d2ri_dcsadjw_prime) 2832 2833 # Any other parameter. 2834 else: 2835 # Gradient. 2836 data.create_dri_prime.append(func_dri_djw_prime) 2837 2838 # Hessian. 2839 data.create_d2ri_prime.append([]) 2840 for j in xrange(data.total_num_params): 2841 # Residue specific parameter index. 2842 index2 = j - num_diff_params 2843 if index2 < 0: 2844 index2 = None 2845 2846 # Rex. 2847 if index2 != None and data.param_types[index2] == 'Rex': 2848 data.create_d2ri_prime[i].append(None) 2849 2850 # Bond length. 2851 elif index2 != None and data.param_types[index2] == 'r': 2852 data.create_d2ri_prime[i].append(func_d2ri_djwdr_prime) 2853 2854 # CSA. 2855 elif index2 != None and data.param_types[index2] == 'CSA': 2856 data.create_d2ri_prime[i].append(func_d2ri_djwdcsa_prime) 2857 2858 # Any other parameter. 2859 else: 2860 data.create_d2ri_prime[i].append(func_d2ri_djwidjwj_prime) 2861 2862 2863 # Both the bond length and CSA are fixed {}. 2864 ############################################ 2865 2866 if data.r_index == None and data.csa_index == None: 2867 # The main ri component functions 2868 if data.rex_index == None: 2869 data.create_ri_comps = ri_comps 2870 data.create_dri_comps = dri_comps 2871 data.create_d2ri_comps = d2ri_comps 2872 else: 2873 data.create_ri_comps = ri_comps_rex 2874 data.create_dri_comps = dri_comps_rex 2875 data.create_d2ri_comps = d2ri_comps 2876 2877 # Calculate the dipolar and CSA constant components. 2878 comp_dip_const_func(data, data.bond_length) 2879 comp_csa_const_func(data, data.csa) 2880 for i in xrange(data.num_ri): 2881 data.dip_comps_func[i] = data.dip_const_func 2882 if data.create_dip_func[i]: 2883 data.dip_comps_func[i] = data.create_dip_func[i](data.dip_const_func) 2884 if data.create_csa_func[i]: 2885 data.csa_comps_func[i] = data.create_csa_func[i](data.csa_const_func[data.remap_table[i]]) 2886 2887 2888 # The bond length is a parameter {r}. 2889 ##################################### 2890 2891 elif data.r_index != None and data.csa_index == None: 2892 # The main ri component functions 2893 if data.rex_index == None: 2894 data.create_ri_comps = ri_comps_r 2895 data.create_dri_comps = dri_comps_r 2896 data.create_d2ri_comps = d2ri_comps_r 2897 else: 2898 data.create_ri_comps = ri_comps_r_rex 2899 data.create_dri_comps = dri_comps_r_rex 2900 data.create_d2ri_comps = d2ri_comps_r 2901 2902 # Calculate the CSA constant. 2903 comp_csa_const_func(data, data.csa) 2904 for i in xrange(data.num_ri): 2905 if data.create_csa_func[i]: 2906 data.csa_comps_func[i] = data.create_csa_func[i](data.csa_const_func[data.remap_table[i]]) 2907 2908 2909 # The CSA is a parameter {CSA}. 2910 ############################### 2911 2912 elif data.r_index == None and data.csa_index != None: 2913 # The main ri component functions 2914 if data.rex_index == None: 2915 data.create_ri_comps = ri_comps_csa 2916 data.create_dri_comps = dri_comps_csa 2917 data.create_d2ri_comps = d2ri_comps_csa 2918 else: 2919 data.create_ri_comps = ri_comps_csa_rex 2920 data.create_dri_comps = dri_comps_csa_rex 2921 data.create_d2ri_comps = d2ri_comps_csa 2922 2923 # Calculate the dipolar constant. 2924 comp_dip_const_func(data, data.bond_length) 2925 for i in xrange(data.num_ri): 2926 data.dip_comps_func[i] = data.dip_const_func 2927 if data.create_dip_func[i]: 2928 data.dip_comps_func[i] = data.create_dip_func[i](data.dip_const_func) 2929 2930 2931 # Both the bond length and CSA are parameters {r, CSA}. 2932 ####################################################### 2933 2934 elif data.r_index != None and data.csa_index != None: 2935 # The main ri component functions 2936 if data.rex_index == None: 2937 data.create_ri_comps = ri_comps_r_csa 2938 data.create_dri_comps = dri_comps_r_csa 2939 data.create_d2ri_comps = d2ri_comps_r_csa 2940 else: 2941 data.create_ri_comps = ri_comps_r_csa_rex 2942 data.create_dri_comps = dri_comps_r_csa_rex 2943 data.create_d2ri_comps = d2ri_comps_r_csa 2944 2945 2946 # Invalid combination of parameters. 2947 #################################### 2948 2949 else: 2950 print "Invalid combination of parameters for the model-free equations." 2951 return 0 2952 2953 return 1

Source Code for Module maths_fns.mf