Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

add Gadget2Particle #10

Open
wants to merge 20 commits into
base: master
Choose a base branch
from
Open

Conversation

elehcim
Copy link
Contributor

@elehcim elehcim commented Nov 25, 2021

  • Many routines translated from pynbody to read gadget2 files (format 1 and 2).
  • Drastically reduced code repetition.

TODO:

  • A config file needed for collecting gadget user defined fields and relative units.
  • StructArray data type (e.g. Float32 of Float64) can't be instantiated at runtime. I think it should be possible, but still I need to find a way.
  • Rebase with latest modifications upstream
  • Implement write operations with new block based mechanics
  • Implement Mass block header dependant read/write operation

I pushed because I saw you're very nicely going forward. I did this work in the past two weeks but newly arrived commitments hindered me to finishing this up. By having it here, you can have a look and probably improve or pick up something for the codebase. I'll try to keep this going in any case.

The main architecture change is the introduction of a Gadget2Block type which detects and reads 1D and 3D blocks and knows about the datatype (Float32, Float64, Int) of the data block itself. This uses some routines from pynbody and reduces a lot of code repetition.
Then a newly created Gadget2Particle struct is added containing fields commonly used in Gadget2 files. This type is then used to load a StructArray and fill it up with a common file-reading machinery on blocks.

- Many routines translated from pynbody to read gadget2 files (format 1 and 2).
- Drastically reduced code repetition.

TODO:
- A config file needed for collect gadget user defined fields and relative units.
- StructArray data type (e.g. Float32 of Float64) can't be instantiated at runtime. I think it should be possible, but still I need to find a way.
- Rebase with latest modifications upstream
@islent islent self-assigned this Nov 26, 2021
@islent
Copy link
Member

islent commented Nov 26, 2021

Thanks! I will review the code next week (soon after I have updated documentations and release packages of JuliaAstroSim). Before that, AstroIO will release a minor patch to solve some compat problems.

StructArray data type (e.g. Float32 of Float64) can't be instantiated at runtime. I think it should be possible, but still I need to find a way.

My clumsy solution is to construct a small StructArray and then empty it:

AstroIO.jl/src/Gadget.jl

Lines 285 to 289 in 901d337

data = StructArray(Star(units))
empty!(data)
for k in 1:6
append!(data, StructArray(Star(units, collection = GadgetTypes[k]) for i = 1:header.npart[k]))
end

gadget user defined fields and relative units.

I decide to force user to determine what units they are using while reading snapshots, for example,

header, data = read_gadget2("snapshot.gadget2", uAstro)

where argument uAstro is required thereafter. And the units in snapshot file can be change by

header, data = read_gadget2("snapshot.gadget2", uAstro, uGadget2) #! There will be unit conversions

write_gadget2("output.Gadget2", header, data, uGadget2)

where uGadget2 is the default units for I/O.

More elegant way, not having to append!
(which was causing a 'type does not have a definite number of fields' error).
Also probably more efficient because memory is contiguous from the beginning.
@elehcim
Copy link
Contributor Author

elehcim commented Nov 29, 2021

My clumsy solution is to construct a small StructArray and then empty it:

I think I found a better solution by initializing a full-length StructArray and then write the proper Collection value in the correct positions. See last commit.
This way there's no need to append different sub-StructArray.

The append! was also creating a 'type does not have a definite number of fields' error with the newly defined Gadget2Particle with explicit parametric types.

Now I'm updating the code with the possibility to use other units like you introduced recently.

Last step is to make the choice of fields more flexible via a configuration file and probably some metaprogramming in order to create on-demand additional fields in the StructArray.
For example if the user needs to read a file with a custom field, the user should be able to enter the block name, the mapped name and its fileunits, e.g. POS, :Position, "kpc".
Accordingly, the now hardcoded fields and units should be moved in a configuration file in the future.
For reference, the idea of a configuration file comes from pynbody with the .pynbodyrc conf file.

@islent
Copy link
Member

islent commented Nov 29, 2021

initializing a full-length StructArray and then write the proper Collection value in the correct positions

Wow! Great idea!

As for the conf file, ConfParser.jl should be helpful. I had thought about it once before, however it's hard to make the whole interface user-friendly. pynbody might be a good reference.

My idea is that, users can choose to

  1. write a conf file (AstroIO will construct a handler (struct, Dict or Array) that handles all information. POS, :Position, "kpc" are packed in a struct)
  2. construct the handler manually

@elehcim
Copy link
Contributor Author

elehcim commented Nov 29, 2021

Ok, thanks. I'll have a look at ConfParser.jl

@islent
Copy link
Member

islent commented Nov 29, 2021

Are you pushing to https://github.com/elehcim/AstroIO.jl? Where can I contribute some help?

@elehcim
Copy link
Contributor Author

elehcim commented Nov 29, 2021

Now I have a problem regarding unit conversion which seems to change the datatype.

For example:

source_units = getuVel(uGadget2);
target_units = getuVel(uAstro);
v = PVector(3.f0, 3.f0, 4.f0, source_units);
println(typeof(uconvert(target_units, v.x)))
     Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(Gyr^-1, kpc), 𝐋 𝐓^-1, nothing}}

but I'd like it to stay a Float32.

I don't know now the solution. Do you think a type stable uconvert for PVector would be a solution? Any idea?

I can open an issue in PhysicalParticles

@elehcim
Copy link
Contributor Author

elehcim commented Nov 29, 2021

Are you pushing to https://github.com/elehcim/AstroIO.jl? Where can I contribute some help?

Yes, I'm pushing there.
You're very welcome to contrbute, don't know exactly what is the best way :D. This PR allows maintainers to edit.

@codecov
Copy link

codecov bot commented Dec 5, 2021

Welcome to Codecov 🎉

Once you merge this PR into your default branch, you're all set! Codecov will compare coverage reports and display results in all future pull requests.

Thanks for integrating Codecov - We've got you covered ☂️

@elehcim elehcim marked this pull request as ready for review December 5, 2021 21:23
@elehcim elehcim changed the title add Gadget2Particle as it is. WIP - request for comments add Gadget2Particle Dec 5, 2021
@elehcim
Copy link
Contributor Author

elehcim commented Dec 6, 2021

How can I infer the numeric type from a PVector?
Example:

PVector{Quantity{Float32, 𝐋, Unitful.FreeUnits{(kpc,), 𝐋, nothing}}}  --> Float32

to be used like:

julia> infer_numeric_type(eltype(data.Pos))
Float32

@elehcim
Copy link
Contributor Author

elehcim commented Dec 6, 2021

How can I infer the numeric type from a PVector? Example:

PVector{Quantity{Float32, 𝐋, Unitful.FreeUnits{(kpc,), 𝐋, nothing}}}  --> Float32

to be used like:

julia> infer_numeric_type(eltype(data.Pos))
Float32

For now I'll do:

eltype(d.Pos).parameters[1].parameters[1]

@islent
Copy link
Member

islent commented Dec 7, 2021

eltype(d.Pos).parameters[1].parameters[1]

I get a better solution (in the latest master branch of PhysicalParticles):
https://github.com/JuliaAstroSim/PhysicalParticles.jl/blob/4114ddc0e4ac959fa8c6cfe9c7275f5662ed6983/src/PVector.jl#L198-L205

Have fun with it :D

julia> using PhysicalParticles

julia> a = PVector(1.0f0, 2.0f0, 3.0f0)
PVector{Float32}(1.0f0, 2.0f0, 3.0f0)

julia> numeric_type(a)
Float32

julia> numeric_type([a])
Float32

@islent
Copy link
Member

islent commented Dec 16, 2021

Merged through #12. (Sorry that I was not familiar with editing PRs before merging. I should have merged this PR first and then commit to master branch)

A config file needed for collecting gadget user defined fields and relative units.

I have implemented a config file loading procedure in ConfParser.jl:

  1. load a config file with field [io]
  2. define IOConfig
  3. load snapshot using IOConfig

There is an example of loading JLD2 file:

AstroIO.jl/test/runtests.jl

Lines 142 to 145 in 1bb9949

@testset "ConfParser" begin
ioconfig = loadconfig("conf-jld2.ini")
data = loadfromconfig(ioconfig)
@test length(data) == 10

How to define blocks from iofields (which is a Vector of 4-Char strings, for example, ["POS ", "VEL "...])?

struct IOConfigGadget2 <: IOConfig
filename::String
format2::Bool
iofields::Vector{String}
end

@islent
Copy link
Member

islent commented Dec 16, 2021

I'm also planing to remove SPH fields from PhysicalParticles.Star and migrate Gadget2Particle to PhysicalParticles. Will these cause any error?

@elehcim
Copy link
Contributor Author

elehcim commented Dec 20, 2021

Hi!

Merged through #12. (Sorry that I was not familiar with editing PRs before merging. I should have merged this PR first and then commit to master branch)

No problem, thanks for merging. I think I'll close this one and open a new one dedicated to the write operations. And also for pushing forward the config file.
I'm thinking to use the StructTypes.jl model.

A config file needed for collecting gadget user defined fields and relative units.

I have implemented a config file loading procedure in ConfParser.jl:

  1. load a config file with field [io]
  2. define IOConfig
  3. load snapshot using IOConfig

Thanks, I'll look into it.

There is an example of loading JLD2 file:

AstroIO.jl/test/runtests.jl

Lines 142 to 145 in 1bb9949

@testset "ConfParser" begin
ioconfig = loadconfig("conf-jld2.ini")
data = loadfromconfig(ioconfig)
@test length(data) == 10

How to define blocks from iofields (which is a Vector of 4-Char strings, for example, ["POS ", "VEL "...])?

struct IOConfigGadget2 <: IOConfig
filename::String
format2::Bool
iofields::Vector{String}
end

@elehcim
Copy link
Contributor Author

elehcim commented Dec 20, 2021

I'm also planing to remove SPH fields from PhysicalParticles.Star and migrate Gadget2Particle to PhysicalParticles. Will these cause any error?

I don't think so. But I'm thinking that in PhysicalParticles.jl there should be a struct which is your interface to your simulation code.
I see Gadget2Particle as a struct for post processing simulation data.

It may contain what is there in the gadget2 file (guided by the config file).

If we convene to use StructTypes.jl, it actually may contain many custom fields generated on the fly (if I understand well the contructfrom method).
Also StructTypes.jl has a nice architecture to select serializable fields.

Let me know what do you think.

@elehcim
Copy link
Contributor Author

elehcim commented Dec 20, 2021

Ah, after commit b655ee3, I think I can continue here.

@islent
Copy link
Member

islent commented Dec 21, 2021

But I'm thinking that in PhysicalParticles.jl there should be a struct which is your interface to your simulation code.

Yep. I'm using Star for gravity-only simulations. So SPH info is redundant for this.

There are only two restrictions to user-defined particle types for a N-body force solver:

  1. Be a subtype of AbstractParticle
  2. Contain at least the same data fields in Star

So It's ok to leave Gadget2Particle to AstroIO. And the situation is same with StructTypes.jl (just to emphasize it, use immutable struct for higher performance)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants