guilhem_lavaux/cosmotool
0
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Working ICs with ramses

Browse Source
This commit is contained in:
Guilhem Lavaux 2014-06-03 18:50:04 +02:00
parent 8f582707da
commit a82688b0b1
3 changed files with 41 additions and 19 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
python_sample/icgen/borgicgen.py
View File

@ -71,7 +71,7 @@ def run_generation(input_borg, a_borg, a_ic, cosmo, supersample=1, do_lpt2=True)
# Compute LPT scaling coefficient
D1 = cgrowth.D(a_ic)
D1_0 = D1/cgrowth.D(a_borg)
velmul = cgrowth.compute_velmul(a_ic)*D1_0
velmul = cgrowth.compute_velmul(a_ic)
D2 = -3./7 * D1_0**2
@ -87,6 +87,8 @@ def run_generation(input_borg, a_borg, a_ic, cosmo, supersample=1, do_lpt2=True)
# Generate vel
vel.append((psi*velmul).astype(np.float32))
print("velmul=%lg" % (cosmo['h']*velmul))
density = np.fft.irfftn(density_hat*D1_0)*(N/L)**3
return posx,vel,density,N,L,a_ic
@ -97,7 +99,7 @@ def write_icfiles(*generated_ic, **cosmo):
ct.simpleWriteGadget("borg.gad", posx, velocities=vel, boxsize=L, Hubble=cosmo['h'], Omega_M=cosmo['omega_M_0'], time=a_ic)
for i,c in enumerate(['x','y','z']):
ct.writeGrafic("ic_velc%s" % c, vel[i].reshape((N,N,N)), L, a_ic, **cosmo)
ct.writeGrafic("ic_velc%s" % c, vel[i].reshape((N,N,N)).transpose(), L, a_ic, **cosmo)
ct.writeGrafic("ic_deltab", density, L, a_ic, **cosmo)
ct.writeWhitePhase("white.dat", density)

20
python_sample/icgen/gen_ic_from_borg.py
View File

@ -8,23 +8,9 @@ cosmo['omega_k_0'] = 0
cosmo['omega_B_0']=0.049
cosmo['SIGMA8']=0.8344
TestCase=True
zstart=10
zstart=200
astart=1/(1.+zstart)
if TestCase:
pos,_,density,N,L,_ = bic.run_generation("initial_condition_borg.dat", 0.001, astart, cosmo, supersample=1, do_lpt2=True)
dcic = ct.cicParticles(pos, L, N)
dcic /= np.average(np.average(np.average(dcic, axis=0), axis=0), axis=0)
dcic -= 1
dcic_hat = np.fft.rfftn(dcic)*(L/N)**3
Pcic, bcic = bic.bin_power(np.abs(dcic_hat)**2/L**3, L, bins=50)
borg_evolved = ct.read_borg_vol("final_density_1380.dat")
snapid=1
if __name__=="__main__":
if not TestCase:
bic.write_icfiles(*bic.run_generation("initial_condition_borg.dat", 0.001, astart, cosmo, do_lpt2=True), **cosmo)
bic.write_icfiles(*bic.run_generation("initial_condition_borg.dat", 0.001, astart, cosmo, do_lpt2=False), **cosmo)

34
python_sample/icgen/test_ic_from_borg.py Normal file
View File

@ -0,0 +1,34 @@
import numpy as np
import cosmotool as ct
import borgicgen as bic
cosmo={'omega_M_0':0.3175, 'h':0.6711}
cosmo['omega_lambda_0']=1-cosmo['omega_M_0']
cosmo['omega_k_0'] = 0
cosmo['omega_B_0']=0.049
cosmo['SIGMA8']=0.8344
snap_id=11
s = ct.loadRamsesAll("/nethome/lavaux/remote2/borgsim/", snap_id, doublePrecision=True)
astart=s.getTime()
pos,_,density,N,L,_ = bic.run_generation("initial_condition_borg.dat", 0.001, astart, cosmo, supersample=1, do_lpt2=True)
dcic = ct.cicParticles(pos, L, N)
dcic /= np.average(np.average(np.average(dcic, axis=0), axis=0), axis=0)
dcic -= 1
dsim = ct.cicParticles(s.getPositions(), L, N)
dsim /= np.average(np.average(np.average(dsim, axis=0), axis=0), axis=0)
dsim -= 1
dcic_hat = np.fft.rfftn(dcic)*(L/N)**3
dsim_hat = np.fft.rfftn(dsim)*(L/N)**3
Pcic, bcic = bic.bin_power(np.abs(dcic_hat)**2/L**3, L, bins=50)
Psim, bsim = bic.bin_power(np.abs(dsim_hat)**2/L**3, L, bins=50)
borg_evolved = ct.read_borg_vol("final_density_1380.dat")