PrincetonUniversity / SPEC

The Stepped-Pressure Equilibrium Code, an advanced MRxMHD equilibrium solver.
https://princetonuniversity.github.io/SPEC/
GNU General Public License v3.0
25 stars 6 forks source link

Problem with utilities with new HDF5 output #87

Closed jloizu closed 4 years ago

jloizu commented 4 years ago

@jonathanschilling I realized that the simplest reading routine read_hdf5 does not work when reading output from results generated using newest master branch. I get this error

Undefined function or variable 'filename'.
Error in read_hdf5>getGroup (line 58)
data_info = h5info(filename,root);

when calling read_hdf5(filename) with filename being the hdf5 output file name. Any clue?

jonathanschilling commented 4 years ago

@jloizu This is expected, since read_hdf5.m was the old reading routine for the HDF5 part of the many output files. The corresponding new routine is read_spec.m, which replaces read_hdf5.m, read_spec_field.m, read_spec_iota.m, read_spec_poincare.m and read_spec_grid.m. The output format of the read_spec.m is basically the unification of the outputs of the separate routines. Therefore, we have to adjust all other matlab routines which directly use the output format, since the new format has e.g. the outputs in data.output.Rbc instead of data.Rbc and so on for a data struct generated by calling data=read_spec('G3V02L1Fi.001.sp.h5');. As a starting point, I adjusted plot_spec_poincare.m accordingly. Sorry for the inconvenience.

Edit: I pushed some small hotfix changes to master which handle the case when no Poincare data is present in the output file.

jloizu commented 4 years ago

@jonathanschilling Thanks, this makes total sense. I get it now :) and I will try to adapt the few matlab routines that need it. Shall we remove then the routines read_hdf5.m, read_spec_field.m, read_spec_iota.m, read_spec_poincare.m and read_spec_grid.m?

jonathanschilling commented 4 years ago

@jloizu Cool that it works for you now :-) I left the other reading routines still there to have them at hand when converting the other matlab routines, but I agree that they shall be deleted after the conversion has taken place.

jloizu commented 4 years ago

@jonathanschilling How about the Hessian matrix, is it written somewhere? The info used to be stored by SPEC in the HDF5 file and I would read it with read_spec_hessian.m

jloizu commented 4 years ago

@jonathanschilling by the way, Antoine and I are creating now a branch from master called matlab_update in which we will work on all modifications required to bring these tools to a working level with the newest version of SPEC.

jonathanschilling commented 4 years ago

@jloizu The hessian is still written in the old output format, so one has to use read_spec_hessian.m for the moment. I plan to include it in the "final" HDF5 layout, which I am currently working on.

@jloizu @abaillod Thanks for taking care of the matlab routines! Please let me know if you need any help regarding e.g. which quantities from the old output files ended up where in the new output file.

jloizu commented 4 years ago

The hessian is still written in the old output format, so one has to use read_spec_hessian.m for the moment. I plan to include it in the "final" HDF5 layout, which I am currently working on.

Well there are really "two hessians": the original one (which is an approximate hessian useful for the newton iterations) written by Stuart in binary. And another one (which is more exact and useful for stability calculations) written by me in the HDF5 file when the input flag LHmatrix=T. This last one I cannot find anymore.

jonathanschilling commented 4 years ago

Hmm... I was indeed referring to Stuart's binary force-gradient matrix in my last comment. When I branched off for the unified HDF5 file, I included the original HDF5 writing routine hdfint.f90 as subroutine hdfint in the file sphdf5.f90. Right now, I cannot find the hessian part in the old hdfint.h. Have you added this feature after I started working on the issue68 branch?

As far as I can see, the flag LHmatrix is only used (at the time of branching off) in hesian.h, line 398, which should be the binary force-gradient matrix. This one is written in .<ext>.GF.ma and is therefore hidden on Unix systems.

jloizu commented 4 years ago

OK I think that regarding the LHmatrix flag (which I created some time ago, even though I cannot go enough back in time in the commit history, it's like we lost memory of what SPEC use to be, say, at initial import, and thereafter...any clue?) it is fine, right now the Hessian is still written in a binary file .GF.ma (not HDF5 as I had said, that was never the case). So we need to work on this, namely, to have this written in the HDF5 as well, or to keep it as binary (if potentially too large) and adapt the routine in matlab that can read it.

On the other hand, I think that because of the time you branched off, the variables rpoland rtorthat I added are still here but not written in the HDF5 output anymore. We should add that in sphdf5.

jloizu commented 4 years ago

On the other hand, I think that because of the time you branched off, the variables rpol and rtor that I added are still here but not written in the HDF5 output anymore. We should add that in sphdf5.

I just did that on the matlab_update branch, so it will come to master as we merge.

jloizu commented 4 years ago

So we need to work on this, namely, to have this written in the HDF5 as well, or to keep it as binary (if potentially too large) and adapt the routine in matlab that can read it.

I adapted the routine read_spec_hessian.m to read the binary containing the hessian matrix. We can work with this for now. Maybe it is the best case anyway since the hessian can become a pretty large matrix.