2022-06-08 13:51:09 +02:00
|
|
|
# Call with spectra_computation.py <time> <kind>
|
2022-07-04 12:01:37 +02:00
|
|
|
# time = 'ics' for ICs, = 'z=1' for redshift z=1, = 'end' for final results
|
2022-06-08 13:51:09 +02:00
|
|
|
# kind = 'power' for power spectra comparing same resolution, 'cross' for comparing across all resolutions
|
2022-06-08 12:34:59 +02:00
|
|
|
|
|
|
|
import itertools
|
2022-06-07 16:42:58 +02:00
|
|
|
import subprocess
|
2022-12-20 16:08:51 +01:00
|
|
|
from pathlib import Path
|
2022-06-08 12:34:59 +02:00
|
|
|
from sys import argv
|
2022-06-07 16:42:58 +02:00
|
|
|
|
2022-12-20 15:39:26 +01:00
|
|
|
from jobrun.jobrun import jobrun
|
2022-06-07 16:42:58 +02:00
|
|
|
from paths import base_dir, spectra_dir
|
2022-07-21 12:34:57 +02:00
|
|
|
from spectra_plot import waveforms
|
2022-06-07 16:42:58 +02:00
|
|
|
|
2022-12-20 15:39:26 +01:00
|
|
|
vsc = True
|
|
|
|
|
|
|
|
|
|
|
|
def spectra_jobrun(args):
|
|
|
|
if vsc:
|
2022-12-21 15:47:28 +01:00
|
|
|
jobrun(args, time="12:00:00", nodes=1, source=Path("/gpfs/data/fs71636/lwinkler/spack-latest.sh"))
|
2022-12-20 15:39:26 +01:00
|
|
|
else:
|
|
|
|
subprocess.run(args, check=True)
|
|
|
|
|
2022-06-07 16:42:58 +02:00
|
|
|
|
2022-08-10 16:26:30 +02:00
|
|
|
def run_spectra(
|
2022-12-20 15:39:26 +01:00
|
|
|
waveform: str, resolution_1: int, resolution_2: int, Lbox: int, time: str
|
2022-08-10 16:26:30 +02:00
|
|
|
):
|
2022-06-08 13:51:09 +02:00
|
|
|
print("starting")
|
2022-08-10 16:26:30 +02:00
|
|
|
setup_1 = f"{waveform}_{resolution_1}_{Lbox}"
|
|
|
|
setup_2 = f"{waveform}_{resolution_2}_{Lbox}"
|
2022-06-07 16:42:58 +02:00
|
|
|
|
2022-06-08 12:34:59 +02:00
|
|
|
# #For ICs: time == 'ics'
|
2022-08-10 16:26:30 +02:00
|
|
|
if time == "ics":
|
|
|
|
output_file = (
|
2022-12-20 15:39:26 +01:00
|
|
|
base_dir
|
|
|
|
/ f"spectra/{waveform}_{Lbox}/{waveform}_{Lbox}_ics_{resolution_1}_{resolution_2}_cross_spectrum.txt"
|
2022-08-10 16:26:30 +02:00
|
|
|
)
|
2022-12-22 17:22:56 +01:00
|
|
|
output_file.parent.mkdir(exist_ok=True)
|
2022-07-11 12:09:27 +02:00
|
|
|
if output_file.exists():
|
2022-08-10 16:26:30 +02:00
|
|
|
print(f"{output_file} already exists, skipping.")
|
2022-07-11 12:09:27 +02:00
|
|
|
return
|
2022-12-20 15:39:26 +01:00
|
|
|
spectra_jobrun(
|
2022-08-10 16:26:30 +02:00
|
|
|
[
|
|
|
|
str(spectra),
|
|
|
|
"--ngrid",
|
|
|
|
"2048",
|
|
|
|
"--format=4", # This seems to work, but is not as readable
|
|
|
|
"--output",
|
|
|
|
str(
|
|
|
|
base_dir
|
|
|
|
/ f"spectra/{waveform}_{Lbox}/{waveform}_{Lbox}_ics_{resolution_1}_{resolution_2}"
|
|
|
|
),
|
|
|
|
"--input",
|
|
|
|
str(base_dir / f"{setup_1}/ics_{setup_1}.hdf5"),
|
|
|
|
"--input",
|
|
|
|
str(base_dir / f"{setup_2}/ics_{setup_2}.hdf5"),
|
2022-12-20 15:39:26 +01:00
|
|
|
]
|
2022-08-10 16:26:30 +02:00
|
|
|
)
|
|
|
|
|
2022-12-20 15:40:31 +01:00
|
|
|
# #For evaluation of results at redshift z=1: time == 'z=1'
|
2022-08-10 16:26:30 +02:00
|
|
|
elif time == "z=1":
|
|
|
|
output_file = (
|
2022-12-20 15:39:26 +01:00
|
|
|
base_dir
|
|
|
|
/ f"spectra/{waveform}_{Lbox}/{waveform}_{Lbox}_a2_{resolution_1}_{resolution_2}_cross_spectrum.txt"
|
2022-08-10 16:26:30 +02:00
|
|
|
)
|
2022-12-22 17:22:56 +01:00
|
|
|
output_file.parent.mkdir(exist_ok=True)
|
2022-07-11 12:09:27 +02:00
|
|
|
if output_file.exists():
|
2022-08-10 16:26:30 +02:00
|
|
|
print(f"{output_file} already exists, skipping.")
|
2022-07-11 12:09:27 +02:00
|
|
|
return
|
2022-12-20 15:39:26 +01:00
|
|
|
spectra_jobrun(
|
2022-08-10 16:26:30 +02:00
|
|
|
[
|
|
|
|
str(spectra),
|
|
|
|
"--ngrid",
|
2022-12-20 15:40:31 +01:00
|
|
|
"2048",
|
2022-08-10 16:26:30 +02:00
|
|
|
"--format=3",
|
|
|
|
"--output",
|
|
|
|
str(
|
|
|
|
base_dir
|
|
|
|
/ f"spectra/{waveform}_{Lbox}/{waveform}_{Lbox}_a2_{resolution_1}_{resolution_2}"
|
|
|
|
),
|
|
|
|
"--input",
|
|
|
|
str(base_dir / f"{setup_1}/output_0002.hdf5"),
|
|
|
|
"--input",
|
|
|
|
str(base_dir / f"{setup_2}/output_0002.hdf5"),
|
2022-12-20 15:39:26 +01:00
|
|
|
]
|
2022-08-10 16:26:30 +02:00
|
|
|
)
|
2022-06-08 13:51:09 +02:00
|
|
|
|
2022-12-20 15:39:26 +01:00
|
|
|
|
2022-06-08 12:34:59 +02:00
|
|
|
# #For evaluation of final results: time == 'end'
|
2022-08-10 16:26:30 +02:00
|
|
|
elif time == "end":
|
|
|
|
output_file = (
|
2022-12-20 15:39:26 +01:00
|
|
|
base_dir
|
|
|
|
/ f"spectra/{waveform}_{Lbox}/{waveform}_{Lbox}_a4_{resolution_1}_{resolution_2}_cross_spectrum.txt"
|
2022-08-10 16:26:30 +02:00
|
|
|
)
|
2022-12-22 17:22:56 +01:00
|
|
|
output_file.parent.mkdir(exist_ok=True)
|
2022-07-11 12:09:27 +02:00
|
|
|
if output_file.exists():
|
2022-08-10 16:26:30 +02:00
|
|
|
print(f"{output_file} already exists, skipping.")
|
2022-07-11 12:09:27 +02:00
|
|
|
return
|
2022-12-20 15:39:26 +01:00
|
|
|
spectra_jobrun(
|
2022-08-10 16:26:30 +02:00
|
|
|
[
|
|
|
|
str(spectra),
|
|
|
|
"--ngrid",
|
|
|
|
"2048",
|
|
|
|
"--format=3",
|
|
|
|
"--output",
|
|
|
|
str(
|
|
|
|
base_dir
|
|
|
|
/ f"spectra/{waveform}_{Lbox}/{waveform}_{Lbox}_a4_{resolution_1}_{resolution_2}"
|
|
|
|
),
|
|
|
|
"--input",
|
|
|
|
str(base_dir / f"{setup_1}/output_0004.hdf5"),
|
|
|
|
"--input",
|
|
|
|
str(base_dir / f"{setup_2}/output_0004.hdf5"),
|
2022-12-20 15:39:26 +01:00
|
|
|
]
|
2022-08-10 16:26:30 +02:00
|
|
|
)
|
2022-06-08 13:51:09 +02:00
|
|
|
else:
|
|
|
|
raise ValueError(f"invalid time ({time})")
|
|
|
|
|
2022-07-04 12:01:37 +02:00
|
|
|
print("end")
|
2022-06-08 12:34:59 +02:00
|
|
|
|
2022-08-10 16:26:30 +02:00
|
|
|
|
2022-07-21 12:34:57 +02:00
|
|
|
def power_run(resolutions: list, Lbox: int, time: str):
|
2022-06-08 13:51:09 +02:00
|
|
|
args = []
|
2022-06-08 13:29:41 +02:00
|
|
|
for waveform in waveforms:
|
|
|
|
for resolution in resolutions:
|
2022-08-10 16:26:30 +02:00
|
|
|
args.append((waveform, resolution, resolution, Lbox, time))
|
2022-06-08 13:51:09 +02:00
|
|
|
return args
|
|
|
|
|
2022-06-08 12:34:59 +02:00
|
|
|
|
2022-07-21 12:34:57 +02:00
|
|
|
def cross_run(resolutions: list, Lbox: int, time: str):
|
2022-06-08 13:51:09 +02:00
|
|
|
args = []
|
2022-06-08 13:29:41 +02:00
|
|
|
for waveform in waveforms:
|
2022-06-08 13:51:09 +02:00
|
|
|
for res1, res2 in itertools.combinations(resolutions, 2):
|
2022-08-10 16:26:30 +02:00
|
|
|
args.append((waveform, res1, res2, Lbox, time))
|
2022-06-08 13:51:09 +02:00
|
|
|
return args
|
2022-06-08 12:34:59 +02:00
|
|
|
|
|
|
|
|
2022-08-10 16:26:30 +02:00
|
|
|
if __name__ == "__main__":
|
|
|
|
# input("are you sure you want to run this? This might need a large amount of memory")
|
2022-06-08 13:31:13 +02:00
|
|
|
Lbox = 100
|
2022-07-11 12:09:27 +02:00
|
|
|
resolutions = [128, 256, 512, 1024]
|
2022-06-08 12:34:59 +02:00
|
|
|
|
2022-08-10 16:26:30 +02:00
|
|
|
spectra = spectra_dir / "spectra"
|
2022-06-08 13:25:24 +02:00
|
|
|
time = argv[1]
|
2022-06-08 12:34:59 +02:00
|
|
|
|
2022-08-10 16:26:30 +02:00
|
|
|
if argv[2] == "power":
|
2022-07-21 12:34:57 +02:00
|
|
|
args = power_run(resolutions=resolutions, Lbox=Lbox, time=time)
|
2022-06-08 12:34:59 +02:00
|
|
|
|
2022-08-10 16:26:30 +02:00
|
|
|
elif argv[2] == "cross":
|
2022-07-21 12:34:57 +02:00
|
|
|
args = cross_run(resolutions=resolutions, Lbox=Lbox, time=time)
|
2022-06-08 13:51:09 +02:00
|
|
|
else:
|
|
|
|
raise ValueError("missing argv[2] (power|cross)")
|
2022-12-20 16:08:51 +01:00
|
|
|
for a in args:
|
|
|
|
run_spectra(*a)
|
|
|
|
# with Pool(processes=1) as p:
|
|
|
|
# p.starmap(run_spectra, args)
|