replaced jacobi() with an open-source version #360
Merged
Conversation
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
…:Diagonalize() from "jacobi_pd.h"
When I submitted this pull request earlier tonight, it seems that I was making edits to the LAMMPS version of this file which is older (instead of starting from the Colvars/colvars repository version). Fortunately, thew two versions of the file are almost the same, and it was easy to make an update.
When I submitted this pull request earlier tonight, it seems that I was making edits to the LAMMPS version of this file which is older (instead of starting from the Colvars/colvars repository version). Fortunately, thew two versions of the file are almost the same, and it was easy to make an update.
a typo (due to trying to edit the file directly on the github web page). sorry about that.
Minor issue: #ifdef should be consistent with the name of the file ("jacobi_pd.h")
11 tasks
|
Many thanks @jewettaij ! Will try to test this asap. |
|
@jewettaij Thanks for setting this up. Indeed the combination of licensing and portability issues of most libraries is very annoying. The general principle here is that compiling a code like LAMMPS with Colvars should not be any more complicated than it already is without it. Basically, we felt that minimizing effort by sysadmins was more important than resolving the NR issue. At least, in the way that we include the Jacobi routine here, it is not immediately usable by others without modifying the Colvars code, at which point people would rather get the function from somewhere else. Looking closer at this code, there are two issues:
Build LAMMPS and then run the following in
|
|
(Closed accidentally by child interference.) |
Thanks. That was what I needed. I think I fixed this, but please take a look at the change and verify that I didn't introduce a memory leak or something else bad. (I allocated some arrays, but did not bother to explicitly deallocate them. Hopefully the destructors for cvm::vector1d and cvm::matrix2d will take care of that.)
The good news is that I don't think I'm using any fancy C 11 specific features that aren't available in C 03. (Other than using "nullptr" instead of "NULL".) I take a quick stab at converting this code to C 03 and if I can get that working within the hour, I'll make another post. |
|
I made a couple changes to "jacobi_pd.h" which hopefully remove the dependence on C 11. I tested this code on my computer and it compiles and runs. However my compiler is C 11, so it's possible there could be some C 11-specific stuff hiding in jacobi_pd.h somewhere. So please let me know if you run into problems. (The LAMMPS developers finally made the switch to C 11) Also, (as I mentioned earlier), please look at the previous commit I made (e77a259) and check that I didn't introduce a memory leak into colvartypes.cpp. Cheers. (When you eventually move to C 11, I suggest uncommenting the move constructor. It is currently commented out with this pattern: //C 11//.) |
|
@jewettaij Thanks. I converted those commented-out bits to The code does build with C 03, at least just the library files. However, I am getting a ton of warnings with the Sun/Oracle compiler about hidden members (see below). These are probably harmless, but they will definitely pop up again when we submit this code to John Stone (lead VMD dev). I have not yet run numerical tests, but will do next. It would be really good if you could integrate these local changes in your original repository, and potentially break up your original As far as LAMMPS is concerned, I'll submit a PR to update Colvars using the code that is currently on |
Hopefully I fixed that in my most recent commit.
I admit that the unit test code at https://github.com/jewettaij/jacobi_pd is pretty messy. (It is full of #ifdefs.) So if it helps, I posted some cleaner testing code along with some results comparing the new code with the old LAMMPS code. (I compared against the "jacobi()" function that was in one of the LAMMPS source files, "math_extra.h". That function was derived from the Numerical Recipes code and it works only with 3x3 matrices.)
Do you want me to split up math_eigen.h or jacobi_pd.h? The "math_eigen.h" file contains two different eigensolvers: Jacobi and a more advanced and much more general solver ("LambdaLanczos") which works on large, sparse matrices. (The "jacobi_pd.h" file only contains Jacobi.) Do you want me to split "math_eigen.h" into two different files (containing Jacobi and LambdaLanczos)? (LambdaLanczos is in "math_eigen.h" for historical reasons. Jacob Gissinger and I originally used the LambdaLanczos code when we were working together on "fix_bond_react.cpp". Later I created "jacobi_pd" and we switched to that. If I draw attention to the fact that we are no longer using the LambdaLanczos code in LAMMPS, I worry Axel will ask me to remove it from "math_eigen.h". That would make me sad. LambdaLanczos is amazing. It's only a couple hundred lines, and it does so much. Somebody will find it very useful eventually.) The current math_eigen.h file is self-contained.
By "original", do you mean the "math_eigen.h" file in the LAMMPS github repository (containing Jacobi LambdaLanczos), or the jacobi_pd repository? I was actually thinking of creating a new repository on github containing the "math_eigen.h" file which includes both eigensolvers (Jacobi and LambdaLanczos). Together these two algorithms cover a wide range of uses. LAMMPS and Colvars could link to this code.
I am hoping there won't be many changes to these files in the future. Perhaps that's naive. But I agree that it would be nice to minimize redundant (in "colvars/src/jacobi_pd.h" and "lammps/src/math_eigen.h"). To make that possible, we would have to persuade Axel to change the LAMMPS makefiles similar to the way your Makefile includes Lepton. Being lazy and having two copies of this code is the first step. I'm not certain whether Axel will accept this new code, but if he does, we can bring this up eventually.
It sounds good to me. |
|
I meant splitting Jacobi from LambdaLanczos (and from any other classes that may be used independently) in Several warnings about hidden members remain. I reduced the number of their repeated occurrences by moving the inclusion of your header to one file only, but that doesn't remove them. One thing I noticed while trying to silence the remaining warnings, is the use of STL math functions, which were problematic in some older versions of Visual Studio (another VMD issue). Asking you to fix these to appease a compiler that may not even be needed any more is definitely not reasonable. See this comment and the parent issue for more details about our portability troubles: I'll keep following the LAMMPS PR to see what is going on and we can update this one accordingly. |
|
Regarding LambdaLanczos, I don't know much about it, but generically Lanczos seems more efficient for larger matrices. The only diagonalization here is 4x4, so in this case any version that is optimized at compile time would actually be even better. |
Sounds good to me.
Agreed. |
|
In case the LAMMPS developers are reading this conversation... Although it doesn't effect compilation speed, I am happy to split "math_eigen.h" file and remove LambdaLanczos as Giacomo suggested. In that case, I'd prefer to keep that code as a separate file if that's okay. I do think it is likely that a future user will eventually want to implement a fix that requires diagonalizing larger matrices, and that code will be really helpful. It's only a couple hundred lines of extra code. As Giacomo mentioned, the licenses used by popular eigensolver libraries are a headache. |
|
I am a little curious about the license issue. The Eigen C library uses Mozilla Public License 2.0 and it should be fine with LGPL, which is used by Colvars. Is including their header files an issue for Colvars? |
Great question. Eigen is almost 9 mb in size. According to the FAQ, I think you are right about the license. Here's an excerpt from something Axel said recently |
@giacomofiorin just discussed what we want to do in LAMMPS with the LOTL and he was happy with the approach I proposed. |
|
Great, thank you for the posting an update here! The most reasonable thing to do is probably adapting the code in this PR so that it conditionally compiles when LAMMPS is being built, and falling back to the NR version of We could probably use your advice on how to best detect through the preprocessor that LAMMPS is being built. Can we e.g. do this? |
nope. |
need to think about this one for a bit. Is this for compiling the colvars library or the USER-COLVARS code? |
|
Right! I saw the snippet in Unfortunately the Jacobi routine is needed in the library code, so yeah, that's not trivial at all without including LAMMPS headers there. We can definitely make the check at runtime without losing performance, but that most likely means carrying the NR version into LAMMPS without a conditional compilation. The simplest option would be to copy @jewettaij's code also into NAMD and VMD, but even if he was okay with that, it is unlikely to make it into VMD without significant changes. |
|
One simpler option may be defining a LAMMPS-specific processor macro just for the Colvars library, by editing the CMake recipe and its Makefile from the LAMMPS distribution. VMD and NAMD would be unaffected. |
Yes. adding something like |
It's okay with me. Feel free to make any changes you like. |
|
Just implemented this in a small-fixes-collection branch that will soon become a pull request. akohlmey/lammps@b965121 |
|
@akohlmey Thanks for your Makefile tweaks: I pushed one further commit to your LAMMPS branch to allow Colvars to include @jewettaij's header directly from the LAMMPS source directory without need to create a copy. If this is not okay please amend it in your branch. There is one instantiation of the Jacobi solver (for 4x4 dimensionality), but it may be replaced with linking a I build and ran the regtest for LAMMPS, NAMD and VMD and they all work (using |
giacomofiorin
added a commit
to giacomofiorin/lammps
that referenced
this pull request
Sep 18, 2020
Small changes Implement global map of components (@HanatoK) Colvars/colvars#363 Format code examples with colored background (@giacomofiorin) Colvars/colvars#361 replaced jacobi() with an open-source version (@jewettaij) Colvars/colvars#360
7 tasks
giacomofiorin
added a commit
to giacomofiorin/lammps
that referenced
this pull request
Sep 18, 2020
Small changes Implement global map of components (@HanatoK) Colvars/colvars#363 Format code examples with colored background (@giacomofiorin) Colvars/colvars#361 replaced jacobi() with an open-source version (@jewettaij) Colvars/colvars#360
shelllbw
pushed a commit
to nufeb/NUFEB-2
that referenced
this pull request
Aug 27, 2021
Small changes Implement global map of components (@HanatoK) Colvars/colvars#363 Format code examples with colored background (@giacomofiorin) Colvars/colvars#361 replaced jacobi() with an open-source version (@jewettaij) Colvars/colvars#360
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
I replaced the matrix diagonalization code with a new file ("jacobi_pd.h").
Why?
This solves a legal problem with the current version of "colvartypes.cpp". Currently that file contains code from "Numerical Recipes in C", which is proprietary and is incompatible with LGPL or GPL.
The new matrix diagonalization code has been extensively tested and benchmarked. (For 4x4 matrices, the new code is actually slightly faster than the current Numerical Recipes code.) The new uses a the CC0 public domain license which is compatible with, LGPL and MIT.
This new code nearly identical to the code from the "math_eigen.h" file distributed with LAMMPS.
(Incidentally, LAMMPS also now includes superpose3d_h, which implements the RMSD-superposition problem. But it looks like you have implemented this already.)
Status: compiled but untested
While the "jacobi_pd" code has been carefully tested, it is possible (if unlikely) that I made a mistake in the way I combined it with your code. I verified that LAMMPS is able to compile using this new code. But I have not tested the code to see whether it is behaving correctly. (I have not invested the time to learn how to use fix_colvars in LAMMPS yet.)
Feel free to ignore this pull request. Either way, I am on a quest to clean out the Numerical Recipes code from LAMMPS.
Cheers
-Andrew