Python function cis.read_data() equivalent in Julia

Specific Domains Data
, ,
1

I want to write Python function cis.read_data() equivalent in Julia.

def read_data(data_list, read_function):
    """
    Wrapper for calling an HDF reading function for each dataset, and then concatenating the result.
    :param list data_list: A list of data objects to read
    :param callable or str read_function: A function for reading the data, or 'SD' or 'VD' for default reading routines.
    :return: A single numpy array of concatenated data values.
    """
    if callable(read_function):
        out = utils.concatenate([read_function(i) for i in data_list])
    elif read_function == 'VD':
        out = utils.concatenate([hdf_vd.get_data(i) for i in data_list])
    elif read_function == 'SD':
        out = utils.concatenate([hdf_sd.get_data(i) for i in data_list])
    else:
        raise ValueError("Invalid read-function: {}, please supply a callable read "
                         "function, 'VD' or 'SD' only".format(read_function))
    return out

I have to read a part of data from a hdf5 file and i will use a config file for selecting that part of data.

image
image840Ă—494 43.2 KB

and the config file of Athena++ is :point_down:

sane00.athinput 
<job>
problem_id = sane00

<time>
cfl_number = 0.25
nlim       = -1
tlim       = 4.0E4

<output1>
file_type = hdf5
variable  = prim
id        = prim
dt        = 10.0
xdmf      = false

<output2>
file_type = hdf5
variable  = uov
id        = user
dt        = 100.0
xdmf      = false

<output3>
file_type = rst
dt        = 1.E4

<mesh>
nx1    = 72
x1min  = 1.6
x1max  = 1200.0
x1rat  = 1.0963050292915717
ix1_bc = user
ox1_bc = user

nx2    = 32
x2min  = 0.0
x2max  = 3.141592653589793
ix2_bc = polar
ox2_bc = polar

nx3    = 16
x3min  = 0.0
x3max  = 6.283185307179586
ix3_bc = periodic
ox3_bc = periodic

refinement = static

<meshblock>
nx1 = 18
nx2 = 4
nx3 = 16

<refinement1>
level = 1
x1min = 1.6
x1max = 200.0
x2min = 0.41
x2max = 2.731592653589793
x3min = 0.0
x3max = 6.283185307179586

<refinement2>
level = 2
x1min = 1.6
x1max = 30.0
x2min = 0.61
x2max = 2.531592653589793
x3min = 0.0
x3max = 6.283185307179586

<refinement3>
level = 1
x1min = 8.5
x1max = 1200.0
x2min = 0.0
x2max = 0.380
x3min = 0.0
x3max = 6.283185307179586

<refinement4>
level = 1
x1min = 8.5
x1max = 1200.0
x2min = 2.761592653589793
x2max = 3.141592653589793
x3min = 0.0
x3max = 6.283185307179586

<refinement5>
level = 2
x1min = 44.4
x1max = 1200.0
x2min = 0.0
x2max = 0.18
x3min = 0.0
x3max = 6.283185307179586

<refinement6>
level = 2
x1min = 44.4
x1max = 1200.0
x2min = 2.961592653589793
x2max = 3.141592653589793
x3min = 0.0
x3max = 6.283185307179586

<coord>
m = 1.0
a = 0.0

<hydro>
gamma     = 1.3333333333333333
dfloor    = 1.0e-6
pfloor    = 1.0e-8
rho_min   = 1.0e-2
rho_pow   = -1.5
pgas_min  = 1.0e-2
pgas_pow  = -2.5
sigma_max = 100.0
beta_min  = 1.0e-3
gamma_max = 50.0

<problem>
k_adi   = 1.0
r_edge  = 40.5
r_peak  = 80.0
l       = 0.0
rho_max = 1.0

field_config         = loops
potential_cutoff     = 0.2
potential_start      = 25.0
potential_end        = 550.0
potential_wavelength = 3.75
beta_min             = 0.05833

potential_r_pow      = 0.0
potential_rho_pow    = 1.0 

sample_n_r     = 288
sample_n_theta = 128
sample_r_rat   = 1.0232525877105085
sample_cutoff  = 0.2

pert_amp = 1.00E-1
pert_kr  = 6.28E-2
pert_kz  = 6.28E-2
2

You could do:

"""
    read_data(read_function, data_list)
Wrapper for calling an HDF reading function for each dataset, and then concatenating the result.
 - data_list: A list of data objects to read
 - read_function: A function for reading the data, or 'SD' or 'VD' for default reading routines.
Return single array of concatenated data values.
"""
function read_data(read_function::Function, data_list)
    # note I change the order of args because it's more Julian this way
    return stack(read_function, data_list)
end

# if you really want that string argument as well you could define something like
function read_data(functionname::String, data_list)
    # now determine the function to call and just call 
    return read_data(myfun, data_list)
end
1 Like
3

You wrote incomplete code. Is read_function defined in Julia? Please tell in detail.

4

read_function in my code is as well defined as it is in your code - it is simply the function’s argument. I don’t know what function your are referring to and just rewrote the function that you gave in your post as a Julia function. The connect to HDF5 below is unclear.

5

See cis.read_data()
I have to convert following python code to julia.

from pyfiles import * 
runname='./' #define start dir and file, uncomment below
sane = 1#0 for MAD98, 1 for SANE00
if sane == 1:
    inputfile=runname+'sane00.athinput'
datafile='sane00.prim.01800.athdf'
data=read_data(datafile,athinput=inputfile)

so i think read_function in my case is config file named
sane00.athinput

6

You link the very same code you posted here that I translated to Julia for you.

The script you posted contains the call read_data(datafile,athinput=inputfile) where both arguments are strings which is not compatible with the function definition you gave me.

If you want to replicate the functionality of that Python project, then the much easier route is to just use PythonCall.jl and call it directly.

7

I have to write efficient code as it will be used for plotting astrophysical numerical simulations outputs of Athena++ so it would be better if i use julia instead of using PythonCall.jl.

8

I have to agree with @abraemer here – from the snippet you posted, his code is the translation of that.

Also this thread sounds quite like “Please program this for me” – I think you will not really find someone doing that for you for free.

PythonCall.jl would also be a good idea to do this step by step. Inside the functions provided above as a start you could first call the corresponding python functions and then turn that into Julia code step by step. That way you could always check whether it still works.

It also looks like the python package you link has a certain maturity already, so if you “just want to reprogram it in Julia” – I doubt you will hit their efficiency “on first try”. Then PythonCall might still be a better start there.

Finally, since you said you want to write it:

  • What have you tried?
  • Where are you stuck?
  • What is your current code?
1 Like
9

In addition to @kellertuer’s questions I have and additional one:
Are you sure this effort is even necessary/worth it? Plotting is usually not the limiting factor in large-scale numerical simulations.

10

I tried it like :point_down: but it didn’t work.

function read_data(data_list::Vector{Any}, read_function::Union{Function, String})
    if typeof(read_function) == Function
        out = vcat([read_function(i) for i in data_list]...)
    elseif read_function == "VD"
        out = vcat([hdf_vd_get_data(i) for i in data_list]...)
    elseif read_function == "SD"
        out = vcat([hdf_sd_get_data(i) for i in data_list]...)
    else
        throw(ArgumentError("Invalid read-function: $read_function, please supply a callable read function, 'VD' or 'SD' only"))
    end
    return out
end

I am unable to find alternative of “VD” and “SD” in Julia.

Honestly i don’t know how much speed will change by doing it but i think it would be better if i keep everything in Julia. Don’t Julia have some alternative of that cis.read_data() function ?

11

To me that seems require a lot of effort for uncertain gains and it is not even clear whether these gains translate to something noticeable in the real world. I’d recommend you focus your effort elsewhere for now and only if it becomes clear that this part is a bottleneck, then you come back and try to improve it.

Not even Python has this function. It belongs to CIS which is a command-line tool that can be used as a Python library. It is a multi-year open source project. You won’t be able to replicate all of this just on your own in any reasonable amount of time. You really should consider just using this tool as is and not try to replicate its functionality.

2 Likes