Hi,
1) Move matrix_exponential_rankN, to lib.dispersion.matrix_exponential. OK. Done.
This is good, as the function is dispersion specific and is not general.
1a) Split into matrix_exponential_rank6, and matrix_exponential_rank7. Add unit tests. OK. Done.
This should be better for speed.
2) In matrices, use the float64 type in the array() function, or drop the '.0' from all numbers. ? I could skip the ".0". The "multiply.outer" with float 64, takes care of the conversion. No need to do float 64 as type.
This is for the fixed matrices. So for example, you could change: m_r1rho_prime = array([ [-1.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, -1.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, -1.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, -1.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]]) to: m_r1rho_prime = array([ [-1, 0, 0, 0, 0, 0], [ 0, -1, 0, 0, 0, 0], [ 0, 0, 0, 0, 0, 0], [ 0, 0, 0, -1, 0, 0], [ 0, 0, 0, 0, -1, 0], [ 0, 0, 0, 0, 0, 0]], float64) This also introduces double spacing, but that is ok when it is for better formatting and hence easier reading. It is much easier to see in the second version which numbers in the matrix are non-zero.
Where should the matrices "live" ?
The important thing is to fully understand the scope and namespace concepts for different objects in Python (for example see https://docs.python.org/2/tutorial/classes.html). The above mentioned m_r1rho_prime never changes - this is a constant. Therefore it only needs to ever be created once. This structure is currently defined in the rr1rho_3d_2site_rankN() function of the lib/dispersion/ns_matrices.py module. You have two choices for the namespace here. You can either use the function namespace or the module namespace. So you should know the differences. The function namespace comes into existence when the function is called, and then it is destroyed or garbage collected once you exit the function. It is temporary. The module namespace is very different. It comes into existence when the module is imported for the first time, and it is only destroyed once the program exists. Therefore if you define m_r1rho_prime in the module namespace rather than in the current function namespace, it will only be created once when relax starts up and is only destroyed when relax exits. This is better than creating m_r1rho_prime every time rr1rho_3d_2site_rankN() is called, once per target function call. This is exactly what you are using for the 'eta_scale' variable in the lib.dispersion.cr72 module.
3) Move rcpmg_star_rankN() from /lib/dispersion/ns_matrices.py into /lib/dispersion/ns_cpmg_2site_3d_matrices.py ? Several possibilities: a) Stay as it is. Good to keep all matrices in one place. b) Move all matrices into each dispersion file. Also very good. This would eliminate ns_matrices.py file. c) Move to separate files. I do not prefer this.
With the change in point 2), shifting of m_r1rho_prime, m_wA, and all m_* structures out of the function namespace and into the module namespace, it will be clearer to avoid a). Option b) is also good, as it is cleaner than having multiple files. Note that if modules get too big in relax, what I normally do is to convert them into a new package. So for example the directory lib/dispersion/ns_cpmg_2site_3d is created and all the functions and variables are shifted into new modules such as lib.dispersion.ns_cpmg_2site_3d.r2eff, lib.dispersion.ns_cpmg_2site_3d.matrices, etc. Just something to keep in mind, I'm not suggesting that solution for now.
4) Shift all the m_* star matrices out. They are "not changing" matrices. Will do.
This is the same as point 2).
5) Shift "dA=omega - offset", and "theta=atan2(spin_lock_fields_i, dA)" out of loops. They dont change: Yep. That is true.
I hope the numpy.atan2 ufunc handles this. If so, this will speed up the model.
6) Unit tests of ns matrices are the same for low and high dimensionality. Yep, this is good. First determine where to store low and high dimensional matrices. Where is the ns matrices located?
I don't understand the question? You can initialise the numpy arrays to pass into the functions in the unit test method, then directly call the functions of lib.dispersion.ns_matrices, and finally compare the results.
7) The complex64 in ns_mmq_2site_mq. Change to complex128, and change systemtests.
What was wrong with the system tests? As these are slow, the system tests are deliberately of low quality and hence change a lot between systems and if the functions are slightly changed. Is the difference significant? Using numpy.complex64 would be better for speed however, as float operations are always much faster than double operations in BLAS/LAPACK. That's why they are defined in different functions in BLAS/LAPACK.
8) Profiling scripts for the 3-site dispersion model. Same as the 2-site models and where states B and C are identical Coming up.
Cheers!
9) Handling numpy v 1.6. Relevant for all systemtests of NS models. How? Both system and unit tests. Under numpy 1.6, no NS models is possible !
If you look at the system tests, you will see in the __init__() method that I have blacklisted a number of tests if the C modules are not compiled (the tests requiring curve-fitting). The same can be done for numpy using the version number: if float(numpy.version.version[:3]) < 1.6: blacklist_the_test Just copy the logic of the C module blacklisted tests in __init__(), and append the tests requiring numpy.eigsum to status.skipped_tests.
10) Move M1, M2, M1* and M2* out of for loops. Coming up.
This may be complicated, but it should be doable and will be great for speed. Cheers, Edward