include RNSID in ET

Issue #2077 closed
Roland Haas created an issue

Pull request is here: https://bitbucket.org/einsteintoolkit/einsteininitialdata/pull-requests/3/new-initial-data-thorn-for-grhydro/diff

It generate ID for rotating star either with uniform rotational profile or describe by differentialy rotating j-law.

    testsuite par file provided and tested
    example of usage parfile in par subdirectory as perl script that generate par file with resolution dx=0.75 suited for workstation tests.
    Hydro_RNSID/par/speed.rpar generate par file usefull to peform speed test
    two stand alone utility are generate at compile time to [RNS] generate initial models (to be imported in Cactus) amd to read [RNS_read] the generated initial model

Keyword: Hydro_RNSID

Comments (10)

  1. Jonah Miller
    • removed comment

    Hi all, Roland Haas asked me to take a look this thorn. I've played with the examples and the documentation and poked my head in the code. It's very nice to have this tool packaged with the Toolkit. Thank you for that! I've included most of my comments inline, but I'll summarize here:

    As far as I can tell, the code is clean and readable and looks like it should be maintainable. I have a few minor complaints---such as missing include guards in header files, but mostly this is nice and clean. My major complaint here is that it would be nice if there was better encapsulation of the original RNSID code. RNS.c and rnsid.c are basically identical except that rnsid.c contains the appropriate interpolation to a 3D grid at the very end. I'd rather there be a backend file which does all the numerics, then separate files for the CLI and ET frontends which simply read in parameter commands and output a solution or interpolate to the grid respectively.

    The tests and examples are really complete and nice. I also really like that there is a single rpar file to generate a wide variety of initial data. That said, a little explanation at the top of each rpar file explaining how to use it and which scenarios it covers might be nice.

    Other than some really minor things (like the inclusion of a nonexistent header file) the Cactus files (like schedule.ccl) all look good.

    I'd really like to see more documentation. The current documentation is essentially a list of required literature. As with all scientific code, a complete description of the method, which already exists in the literature, would be very excessive. That said, I found it frustrating to dig through said literature to look up simple things---things like the ansatz for the differential rotation, the exact form of the ansatz metric, the compactification scheme, and the precise meaning of named variables such as s, mu, and axes_ratio. The documentation should provide a reference document for these details so that the educated yet forgetful reader doesn't need to dig through the literature to find them. At the very least the documentation should give equation numbers. Algorithmic details, such as the type of iteratiion method used (Newton-Raphson?) and the type of derivatives used (a particularly clever application of finite differences!) would also be nice.

    There's a lot of commented-out code. This doesn't hurt maintainability, but it does hurt readability. Since the code is under version control, is there a reason this code can't be removed?

    That about covers it. Overall, great to see this!

  2. Jonah Miller
    • removed comment

    A few more comments.

    I would really like to see a reader tool for data generated with RNSID, just to show the user how it works (in case they want to output for some reason.) It doesn't have to be fancy. For example, something like this (with maybe a bit more boiler plate) would be cool:

    import  numpy as  np
    import matplotlib as mpl
    from matplotlib import pyplot as plt
    import h5py 
    import sys
    
    filename = sys.argv[1]
    
    with h5py.File(filename,'r') as f:
        K=f.attrs['poly K'][0]
        Gamma = f.attrs['poly Gamma'][0]
        axis_ratio = f.attrs['axis_ratio'][0]
        re = f.attrs['re'][0]
        rhoc = f.attrs['rhoc'][0]
        A_diff = f.attrs['A diff'][0]
        Omega = f.attrs['Omega'][0]
        Omega_e = f.attrs['Omega_e'][0]
        SDIV=f.attrs['SDIV'][0]
        MDIV=f.attrs['MDIV'][0]
        rtype = f.attrs['rotation_type'][0]
        eos_type=f.attrs['eos_type'][0]
        eos_file = f.attrs['eos_file'][0]
        Omega_diff = f['Omega_diff'][:]
        alpha = f['alpha'][:]
        energy = f['energy'][:]
        enthalpy = f['enthalpy'][:]
        gama = f['gama'][:]
        mu = f['mu'][:]
        omega = f['omega'][:]
        pressure = f['pressure'][:]
        rho_potential = f['rho_potential'][:]
        s = f['s_qp'][:]
        v2 = f['velocity_sq'][:]
    
    S,MU = np.meshgrid(s,mu,indexing='ij')
    NU = np.sqrt(1 - MU**2)
    XC = S*MU
    ZC = S*NU
    
    plt.pcolor(XC,ZC,np.log(pressure))
    plt.savefig('star_compactified.png',bbox_inches='tight')
    plt.clf()
    
    R = -re*S/(S-1)
    rmax = 2*re
    rmask = R <= rmax
    X = R*MU
    Y = R*NU
    
    plt.pcolor(X[rmask],Y[rmask],np.log(pressure[rmask]))
    plt.savefig('star_spherical_polar.png',bbox_inches='tight')
    
    
    The tabulated EOS that RNSID expects is different from the tabulated EOS that GRHydro expects. The former is a text file, the latter an hdf5 file. But it's not clear how to correctly extract the former from the latter. Is there a utility to do this? If so, can it be packaged? And if not can it be provided?
    
    In a similar vein---the format for the tabulated EOS should be documented in the thorndoc. What does RNSID need/expect? In what units does it expect the table? Etc.
    
    Compilation generates lots of warnings about unused and uninstantiated variables. These warnings should be fixed if possible.
    
  3. Roberto De Pietri
    • removed comment

    Updated the documentation a answered most of the comments. A read of the 2d generated file is produced as cactus utility as well as the stand alone version to pre-generate initial data before submitting a job to HPC machine. A 2d binary file reader that read generate an hdf5 file of interpolated data on a 3d cartesian grid is also provided.

  4. Log in to comment