Refactor GRHD project

Major rework of GRHD project to include new broadcasting interface. The goal is to get a TOV evolution with the AV method.

PR contains various refactors:

• initialdata/: We promote this to a separate env in order to add CairoMakie for plotting.
• refactored and improve initialdata/SphericalAccretion.jl solver so that it agrees with various references: solver from bugner's thesis (extract from bamps), plots from two papers from Papadoulos.
• refactored initialdata/TOV.jl to not depend on dg1d.jl anymore, instead use OrindaryDiffEq and CairoMakie for plots
• copy of cons2prim from SRHD. these two versions should be merged eventually, but will do in a separate PR
• implement spherical1d formulation: this is a rewrite of the evolution GRHD equations to work with non-zero radial shift vector. this is needed to evolve the spherical accretion test in cowling approximation where there should be no dynamics.
• implement rescaled_spherical1d formulation: this is a rewrite of the previous formulation based on the generalized Valencia formulation (cf. Montero, Baumgarte paper from 2017); the auxiliary metric is chosen such that coordinate singularities at r=0 in the source terms are cancelled; doesn't work atm
• allow the spherical1d to include r=0 in the domain by using symmetry conditions to evaluate source terms at r=0 (workaround for rescaled_spherical1d not working right now)
• new ref tests: bondi accretion and TOV star
• implement a differentiate! method on Mesh in order to compute derivatives by hand, if needed

examples/grhd_bondi/fv.toml

+[EquationOfState]
+eos = "polytrope"
+polytrope_k = 0.02143872868864732
+polytrope_gamma = "$(4/3)"
+
+[GRHD]
+bc = "from_id"
+id = "bondi_accretion"
+id_filename = """$(joinpath(ROOTDIR,"initialdata","bondi_accretion.h5"))"""
+
+[Mesh]
+range = [ 1.8, 10.0 ]
+k = 1024
+scheme = "FV"
+periodic = false
+
+[Output]
+every_iteration = 1
+variables       = [ "D", "Sr", "τ", "p", "ρ", "ϵ", "vr", "rhs_D", "rhs_Sr", "rhs_τ", "src_D", "src_Sr", "src_τ" ]
+# aligned_ts      = "$(collect(range(0.01,0.6,step=0.01)))"
+enable1d        = true
+
+[Evolution]
+cfl = 1.0
+tend = 20.0