diff --git a/sbmy_control/low_level.py b/sbmy_control/low_level.py index 26054bf..b5382c0 100644 --- a/sbmy_control/low_level.py +++ b/sbmy_control/low_level.py @@ -64,10 +64,14 @@ def stdout_redirector(stream): try: # Create a temporary file and redirect stdout to it tfile = tempfile.TemporaryFile(mode="w+b") + _redirect_stdout(tfile.fileno()) + # Yield to caller, then redirect stdout back to the saved fd yield + _redirect_stdout(saved_stdout_fd) + # Copy contents of temporary file to the given stream tfile.flush() tfile.seek(0, io.SEEK_SET) @@ -215,9 +219,9 @@ def get_progress_from_logfile(filename): with open(filename, "r") as f: lines = f.readlines() for line in lines: - if " operation " in line: + if "Step " in line: try: - splitted_line = line.split(" operation ")[1] + splitted_line = line.split("Step ")[1] splitted_line = splitted_line.split(":")[0] current_operation = int(splitted_line.split("/")[0]) total_operations = int(splitted_line.split("/")[1]) @@ -237,7 +241,7 @@ def progress_bar_from_logfile(filename:str, desc:str="", verbose:int=1, **kwargs from tqdm import tqdm from time import sleep k=0 - limit=600 + limit=6000 update_interval=0.2 wait_until_file_exists(filename, verbose=verbose, limit=limit) current_operation, total_operations = get_progress_from_logfile(filename) diff --git a/sbmy_control/main.py b/sbmy_control/main.py index 8bb7ff0..0dd9afe 100644 --- a/sbmy_control/main.py +++ b/sbmy_control/main.py @@ -199,25 +199,26 @@ def main(parsed_args): def check_consistency(card_dict, mode): - ## Check consistency of EvolutionMode and ModulePMCOLA - if mode == "PM" or mode == "allPM": - if card_dict["EvolutionMode"] != 1: - raise ValueError(f"EvolutionMode is not 1: EvolutionMode={card_dict["EvolutionMode"]}") - if card_dict["ModulePMCOLA"] != 1: - raise ValueError(f"ModulePMCOLA is not 1: ModulePMCOLA={card_dict["ModulePMCOLA"]}") - elif mode == "tCOLA" or mode == "alltCOLA": - if card_dict["EvolutionMode"] != 2: - raise ValueError(f"EvolutionMode is not 2: EvolutionMode={card_dict["EvolutionMode"]}") - if card_dict["ModulePMCOLA"] != 1: - raise ValueError(f"ModulePMCOLA is not 1: ModulePMCOLA={card_dict["ModulePMCOLA"]}") - elif mode == "LPT": - if card_dict["ModulePMCOLA"] !=0: - raise ValueError(f"ModulePMCOLA is not 0: ModulePMCOLA={card_dict["ModulePMCOLA"]}") - elif mode == "pre_sCOLA" or mode == "post_sCOLA" or mode == "sCOLA" or mode=="allsCOLA": - if card_dict["EvolutionMode"] != 3: - raise ValueError(f"EvolutionMode is not 3: EvolutionMode={card_dict["EvolutionMode"]}") - if card_dict["ModulePMCOLA"] != 1: - raise ValueError(f"ModulePMCOLA is not 1: ModulePMCOLA={card_dict["ModulePMCOLA"]}") + # ## Check consistency of EvolutionMode and ModulePMCOLA + # if mode == "PM" or mode == "allPM": + # if card_dict["EvolutionMode"] != 1: + # raise ValueError(f"EvolutionMode is not 1: EvolutionMode={card_dict["EvolutionMode"]}") + # if card_dict["ModulePMCOLA"] != 1: + # raise ValueError(f"ModulePMCOLA is not 1: ModulePMCOLA={card_dict["ModulePMCOLA"]}") + # elif mode == "tCOLA" or mode == "alltCOLA": + # if card_dict["EvolutionMode"] != 2: + # raise ValueError(f"EvolutionMode is not 2: EvolutionMode={card_dict["EvolutionMode"]}") + # if card_dict["ModulePMCOLA"] != 1: + # raise ValueError(f"ModulePMCOLA is not 1: ModulePMCOLA={card_dict["ModulePMCOLA"]}") + # elif mode == "LPT": + # if card_dict["ModulePMCOLA"] !=0: + # raise ValueError(f"ModulePMCOLA is not 0: ModulePMCOLA={card_dict["ModulePMCOLA"]}") + # elif mode == "pre_sCOLA" or mode == "post_sCOLA" or mode == "sCOLA" or mode=="allsCOLA": + # if card_dict["EvolutionMode"] != 3: + # raise ValueError(f"EvolutionMode is not 3: EvolutionMode={card_dict["EvolutionMode"]}") + # if card_dict["ModulePMCOLA"] != 1: + # raise ValueError(f"ModulePMCOLA is not 1: ModulePMCOLA={card_dict["ModulePMCOLA"]}") + pass def console_main(): diff --git a/sbmy_control/parameters_card.py b/sbmy_control/parameters_card.py index bbdc1a6..cfeda44 100644 --- a/sbmy_control/parameters_card.py +++ b/sbmy_control/parameters_card.py @@ -47,6 +47,15 @@ def register_arguments_card(parser:ArgumentParser): parser.add_argument("--OutputFCsDensity", type=str, default=None, help="Output FCs density file.") parser.add_argument("--OutputFCsSnapshot", type=str, default=None, help="Output FCs snapshot file.") parser.add_argument("--OutputRngStateLPT", type=str, default=None, help="Output RNG state file.") + ## Tests with phiBCs and density + parser.add_argument("--WriteGravPot", type=bool, default=False, help="Write gravitational potential.") + parser.add_argument("--OutputGravitationalPotentialBase", type=str, default=None, help="Output gravitational potential base.") + parser.add_argument("--MeshGravPot", type=int, default=None, help="Mesh for gravitational potential.") + parser.add_argument("--WriteDensity", type=bool, default=False, help="Write density.") + parser.add_argument("--OutputDensityBase", type=str, default=None, help="Output density base.") + parser.add_argument("--MeshDensity", type=int, default=None, help="Mesh for density.") + parser.add_argument("--LoadPhiBCs", type=bool, default=False, help="Load phiBCs.") + parser.add_argument("--InputPhiBCsBase", type=str, default=None, help="Input phiBCs file base.") parser.add_argument("--WriteReferenceFrame", type=bool, default=False, help="Write reference frame (COCA).") parser.add_argument("--OutputMomentaBase", type=str, default=None, help="Output momenta base (COCA).") parser.add_argument("--ReadReferenceFrame", type=bool, default=False, help="Read reference frame (COCA).") @@ -126,6 +135,16 @@ def parse_arguments_card(parsed_args): OutputFCsDensity=parsed_args.OutputFCsDensity, OutputFCsSnapshot=parsed_args.OutputFCsSnapshot, OutputRngStateLPT=parsed_args.OutputRngStateLPT, + ## Tests with phiBCs and density + WriteGravPot=parsed_args.WriteGravPot, + OutputGravitationalPotentialBase=parsed_args.OutputGravitationalPotentialBase, + MeshGravPot=parsed_args.MeshGravPot, + WriteDensity=parsed_args.WriteDensity, + OutputDensityBase=parsed_args.OutputDensityBase, + MeshDensity=parsed_args.MeshDensity, + LoadPhiBCs=parsed_args.LoadPhiBCs, + InputPhiBCsBase=parsed_args.InputPhiBCsBase, + ## Cosmological parameters h=cosmo_dict["h"], Omega_m=cosmo_dict["Omega_m"], Omega_b=cosmo_dict["Omega_b"], @@ -211,6 +230,17 @@ def parse_arguments_card(parsed_args): card_dict["OutputFCsSnapshot"] = main_dict["resultdir"]+"final_particles_"+main_dict["simname"]+".gadget3" if card_dict["OutputRngStateLPT"] is None: card_dict["OutputRngStateLPT"] = main_dict["workdir"]+"rng_state.h5" + ## Tests with phiBCs and density + if card_dict["OutputGravitationalPotentialBase"] is None: + card_dict["OutputGravitationalPotentialBase"] = main_dict["workdir"]+"gravpot_"+main_dict["simname"] + if card_dict["MeshGravPot"] is None: + card_dict["MeshGravPot"] = card_dict["N_PM_mesh"] + if card_dict["OutputDensityBase"] is None: + card_dict["OutputDensityBase"] = main_dict["workdir"]+"density_"+main_dict["simname"] + if card_dict["MeshDensity"] is None: + card_dict["MeshDensity"] = card_dict["N_PM_mesh"] + if card_dict["InputPhiBCsBase"] is None: + card_dict["InputPhiBCsBase"] = main_dict["workdir"]+"gravpot_tCOLA" if card_dict["OutputMomentaBase"] is None: card_dict["OutputMomentaBase"] = main_dict["workdir"]+"momenta_"+main_dict["simname"]+"_" if card_dict["InputMomentaBase"] is None: @@ -333,6 +363,15 @@ def create_parameter_card_dict( OutputFCsDensity:str = 'fcs_density.h5', OutputFCsSnapshot:str = 'fcs_particles.gadget3', + ## Tests with phiBCs and density + WriteGravPot:bool = True, + OutputGravitationalPotentialBase:str = 'gravitational_potential.h5', + MeshGravPot:int = 128, + WriteDensity:bool = False, + OutputDensityBase:str = 'density.h5', + MeshDensity:int = 128, + LoadPhiBCs:bool = False, + InputPhiBCsBase:str = 'gravitational_potential.h5', ## COCA parameters WriteReferenceFrame:bool = False, OutputMomentaBase:str = 'momenta_', @@ -403,6 +442,46 @@ def create_parameter_card_dict( OutputLPTPotential1=OutputLPTPotential1, OutputLPTPotential2=OutputLPTPotential2, OutputTilesBase=OutputTilesBase, + + # Tests with phiBCs and density + WriteGravPot=int(WriteGravPot), + OutputGravitationalPotentialBase=OutputGravitationalPotentialBase, + MeshGravPot=MeshGravPot, + WriteDensity=int(WriteDensity), + OutputDensityBase=OutputDensityBase, + MeshDensity=MeshDensity, + LoadPhiBCs=int(LoadPhiBCs), + InputPhiBCsBase=InputPhiBCsBase, + + + # P3M default parameters + n_Tiles = 60, + # Number of tiles in each dimension for the P3M code + nPairsForceDiagnostic =0, + # Number of particle-hole pairs to be used for the force diagnostic + nBinsForceDiagnostic =0, + maxTrialsForceDiagnostic =0, + OutputForceDiagnostic ='', + PrintOutputTimestepsLog = 0, + OutputTimestepsLog = "", + fac_dyn_custom = 0.03, + da_max_early_custom = 0.0013, + fac_H_custom = 0.05, + fac_bend = 0.3125, + sub_bend1 = 0.012, + sub_bend2 = 0.014, + fac_p3m_fit = 0.15, + da_max_late_custom = 0.02, + # Polynomial coefficients for the P3M time step limiter. Only used if use_p3m_fit=1 + use_p3m_fit = 1, + p3m_fit_coeff0 = -3.19544528, + p3m_fit_coeff1 = -0.01948223, + p3m_fit_coeff2 = 0.5180241, + p3m_fit_coeff3 = -1.08372045, + p3m_fit_coeff4 = -0.71808255, + p3m_fit_coeff5 = -0.14370312, + RunForceDiagnostic = 0, + h=h, Omega_m=Omega_m, Omega_b=Omega_b,