diff --git a/.git-blame-ignore-revs b/.git-blame-ignore-revs index 107cb9ef..0c984a95 100644 --- a/.git-blame-ignore-revs +++ b/.git-blame-ignore-revs @@ -9,3 +9,6 @@ # apply inifix-format to all .ini files ac260f4af48221229aef4fc1f93ffbf9df450739 1ca7676fee7a578bd64e64ec9b6bdfbceb921775 + +# apply ruff format to all .py files +8d0437bf23d7d8007355fa35b9f094287f7de987 diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 9df5af9d..ddd79ae6 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -30,8 +30,10 @@ repos: - repo: https://github.com/astral-sh/ruff-pre-commit rev: v0.15.8 hooks: + - id: ruff-format + args: [--target-version=py310] - id: ruff-check - args: [--quiet, --fix, --show-fixes] + args: [--target-version=py310, --quiet, --fix, --show-fixes] - repo: https://github.com/la-niche/inifix-pre-commit rev: v1.0.0 diff --git a/doc/source/conf.py b/doc/source/conf.py index f47f6659..1c772f49 100644 --- a/doc/source/conf.py +++ b/doc/source/conf.py @@ -18,13 +18,12 @@ # -- Project information ----------------------------------------------------- -project = 'Idefix' -copyright = 'Geoffroy Lesur et al.' -author = 'Geoffroy Lesur' +project = "Idefix" +copyright = "Geoffroy Lesur et al." +author = "Geoffroy Lesur" # The full version, including alpha/beta/rc tags -release = '2.3.0' - +release = "2.3.0" # -- General configuration --------------------------------------------------- @@ -34,20 +33,24 @@ # ones. extensions = [ "sphinx_rtd_theme", - 'sphinx_git', - 'matplotlib.sphinxext.plot_directive', + "sphinx_git", + "matplotlib.sphinxext.plot_directive", "breathe", "exhale", "m2r2", - "sphinx_copybutton" - ] + "sphinx_copybutton", +] source_suffix = [".rst", ".md"] -cpp_id_attributes=["KOKKOS_FORCEINLINE_FUNCTION","KOKKOS_INLINE_FUNCTION","KOKKOS_RESTRICT"] +cpp_id_attributes = [ + "KOKKOS_FORCEINLINE_FUNCTION", + "KOKKOS_INLINE_FUNCTION", + "KOKKOS_RESTRICT", +] # Add any paths that contain templates here, relative to this directory. -templates_path = ['_templates'] +templates_path = ["_templates"] # List of patterns, relative to source directory, that match files and # directories to ignore when looking for source files. @@ -60,45 +63,47 @@ # The theme to use for HTML and HTML Help pages. See the documentation for # a list of builtin themes. # -html_theme = 'sphinx_rtd_theme' +html_theme = "sphinx_rtd_theme" # Add any paths that contain custom static files (such as style sheets) here, # relative to this directory. They are copied after the builtin static files, # so a file named "default.css" will overwrite the builtin "default.css". -html_static_path = ['_static'] +html_static_path = ["_static"] + def setup(app): - app.add_css_file('my_theme.css') + app.add_css_file("my_theme.css") + ############## # Breathe/Exhale options # Setup the breathe extension breathe_projects = { - "Idefix": "./xml" + "Idefix": "./xml", } breathe_default_project = "Idefix" # Setup the exhale extension exhale_args = { # These arguments are required - "containmentFolder": "./api", - "rootFileName": "library_root.rst", - "rootFileTitle": "Idefix API", + "containmentFolder": "./api", + "rootFileName": "library_root.rst", + "rootFileTitle": "Idefix API", # Suggested optional arguments - "createTreeView": False, + "createTreeView": False, # TIP: if using the sphinx-bootstrap-theme, you need # "treeViewIsBootstrap": True, "doxygenStripFromPath": "../../src", "exhaleExecutesDoxygen": True, - "exhaleDoxygenStdin": textwrap.dedent(''' + "exhaleDoxygenStdin": textwrap.dedent(""" INPUT = ../../src EXCLUDE = ../../src/kokkos - ''') + """), } # Tell sphinx what the primary language being documented is. -primary_domain = 'cpp' +primary_domain = "cpp" # Tell sphinx what the pygments highlight language should be. -highlight_language = 'cpp' +highlight_language = "cpp" diff --git a/doc/source/plot_idefix_bench.py b/doc/source/plot_idefix_bench.py index baa74ac8..7810751e 100644 --- a/doc/source/plot_idefix_bench.py +++ b/doc/source/plot_idefix_bench.py @@ -4,24 +4,24 @@ def do_plot(title, bench_file, gpumodels): - with open(bench_file, 'r') as f: + with open(bench_file, "r") as f: benches = json.load(f) plt.figure() - xmax=0 - ymax=0 + xmax = 0 + ymax = 0 for gpumodel in gpumodels: - select = [bench for bench in benches if bench['gpumodel'] == gpumodel][-1] + select = [bench for bench in benches if bench["gpumodel"] == gpumodel][-1] - xs = [r['nbgpu'] for r in select['results']] - ys = [r['cell_updates'] for r in select['results']] - plt.plot(xs, ys,'o-',label=gpumodel) - xmax=max(xmax,max(xs)) - ymax=max(ymax,max(ys)) + xs = [r["nbgpu"] for r in select["results"]] + ys = [r["cell_updates"] for r in select["results"]] + plt.plot(xs, ys, "o-", label=gpumodel) + xmax = max(xmax, max(xs)) + ymax = max(ymax, max(ys)) plt.xscale("log", base=2) - plt.ylim(0,ymax*1.1) - plt.xlim(1,xmax*1.1) + plt.ylim(0, ymax * 1.1) + plt.xlim(1, xmax * 1.1) plt.legend() plt.xlabel("Number of GPUs/GCDs") plt.ylabel("Performance (cells / second / GPU)") diff --git a/make_tarballs.py b/make_tarballs.py index dfc1b55f..df57b94b 100644 --- a/make_tarballs.py +++ b/make_tarballs.py @@ -89,7 +89,7 @@ def _make_self_contained_tarball(output_dir, *, suffix="", exclude_list=None): os.replace(os.path.join(work_dir, tarball), final_tarfile) filesize = os.path.getsize(final_tarfile) - print("Done ! Final file size is %.1f MB" % (filesize / 1024 ** 2)) + print("Done ! Final file size is %.1f MB" % (filesize / 1024**2)) def main(argv=None): diff --git a/pytools/dump_io.py b/pytools/dump_io.py index 7d0ec580..6c833842 100644 --- a/pytools/dump_io.py +++ b/pytools/dump_io.py @@ -4,6 +4,7 @@ @author: lesurg """ + import os import re import struct @@ -119,6 +120,7 @@ def _read_fields(self, fh): def __repr__(self): return "DumpDataset('%s')" % self.filename + # public API def readDump(filename): return DumpDataset(filename) diff --git a/pytools/idfx_io.py b/pytools/idfx_io.py index 0e53a39a..e530117a 100644 --- a/pytools/idfx_io.py +++ b/pytools/idfx_io.py @@ -13,9 +13,9 @@ HEADER_SIZE = 128 + # There is one .idfx file per processor class IdfxFileField(object): - def __init__(self, fh, byteorder="little"): # read entry name q = fh.read(NAME_SIZE) @@ -32,6 +32,7 @@ def __init__(self, fh, byteorder="little"): raw = struct.unpack(str(ntot) + "d", fh.read(DOUBLE_SIZE * ntot)) self.array = np.asarray(raw).reshape(dims[::-1]) + class IdfxFileDataset(object): def __init__(self, filename): # identical to DumpDataset @@ -45,12 +46,11 @@ def _read_header(self, fh): # could easily be identical to DumpDataset headerSize = 128 q = fh.read(headerSize) - n = q.index(b'\x00') - self.header = q[:n].decode('utf-8') + n = q.index(b"\x00") + self.header = q[:n].decode("utf-8") self.metadata["byteorder"] = "little" - def _read_field(self, fh): # identical to DumpDataset if self.metadata["byteorder"] is None: diff --git a/pytools/idfx_test.py b/pytools/idfx_test.py index ee3fbee3..ac72d6dd 100644 --- a/pytools/idfx_test.py +++ b/pytools/idfx_test.py @@ -13,648 +13,744 @@ class bcolors: - HEADER = '\033[95m' - OKBLUE = '\033[94m' - OKCYAN = '\033[96m' - OKGREEN = '\033[92m' - WARNING = '\033[93m' - FAIL = '\033[91m' - ENDC = '\033[0m' - BOLD = '\033[1m' - UNDERLINE = '\033[4m' + HEADER = "\033[95m" + OKBLUE = "\033[94m" + OKCYAN = "\033[96m" + OKGREEN = "\033[92m" + WARNING = "\033[93m" + FAIL = "\033[91m" + ENDC = "\033[0m" + BOLD = "\033[1m" + UNDERLINE = "\033[4m" -class idfxTest: - def __init__ (self, current_test_file, name=""): - parser = argparse.ArgumentParser() - - idefix_dir_env = os.getenv("IDEFIX_DIR") - parser.add_argument("-noplot", - dest="noplot", - help="disable plotting in standard tests", - action="store_false") - - parser.add_argument("-ploterr", - help="Enable plotting on error in regression tests", - action="store_true") - - parser.add_argument("-cmake", - default=[], - help="CMake options", - nargs='+') - - parser.add_argument("-definitions", - default="", - help="definitions.hpp file") - - parser.add_argument("-dec", - default="", - help="MPI domain decomposition", - nargs='+') - - parser.add_argument("-check", - help="Only perform regression tests without compilation", - action="store_true") - - parser.add_argument("-cuda", - help="Test on Nvidia GPU using CUDA", - action="store_true") - - parser.add_argument("-intel", - help="Test compiling with Intel OneAPI", - action="store_true") - - parser.add_argument("-hip", - help="Test on AMD GPU using HIP", - action="store_true") - - parser.add_argument("-single", - help="Enable single precision", - action="store_true") - - parser.add_argument("-vectPot", - help="Enable vector potential formulation", - action="store_true") - - parser.add_argument("-reconstruction", - type=int, - default=2, - help="set reconstruction scheme (2=PLM, 3=LimO3, 4=PPM)") - - parser.add_argument("-idefixDir", - default=idefix_dir_env, - help="Set directory for idefix source files (default $IDEFIX_DIR)") - - parser.add_argument("-mpi", - help="Enable MPI", - action="store_true") - - parser.add_argument("-all", - help="Do all test suite (otherwise, just do the test with the current configuration)", - action="store_true") - - parser.add_argument("-init", - help="Reinit reference files for non-regression tests (dangerous!)", - action="store_true") - - parser.add_argument("-Werror", - help="Consider warnings as errors", - action="store_true") - - parser.add_argument("-ccache", - help="Use ccache to reduce the build time over multiple run of the test suite.", - action="store_true") - - parser.add_argument("-restart", - help="Enable creating a restart from a checkpoint.", - action='store_true') - - parser.add_argument("-v", "--verbose", - help="Enable verbose mode, by not capturing the output.", - action="store_true") - - parser.add_argument("--help-pytest", - help="Display the options you can transmit directly to pytest in addition to the specific to idefix tests.", - action="store_true") - - parser.add_argument("-fake", - help="Make a fake run by just logging the actions to validate that we generate same command over refactoring.", - action="store_true") - - # this option is not used directly by direct users of idfxTest but by idx_test_gen - parser.add_argument("-subdir", - default="./test", - help="Select the test in the given subdir not to run all.", - type=str) - - args, unknown=parser.parse_known_args() - - # transform all arguments from args into attributes of this instance - self.__dict__.update(vars(args)) - # store the full path of problem directory - self.currentTestFile = current_test_file - self.currentTestName = name - self.problemDir=os.path.dirname(current_test_file) - self.referenceDirectory = os.path.join(idefix_dir_env,"reference") - # current directory relative to $IDEFIX_DIR/test (used to retrieve the path ot reference files) - self.testDir=os.path.relpath(self.problemDir,os.path.join(idefix_dir_env,"test")) - # build directory, currently inside the test named build-test - self.buildDir=os.path.join(self.problemDir,"build-test") - # remind what build we dit last - self.lastCmakeCmd="" - # subdir - self.filterSubdir=args.subdir - # save - self.cmdArgs = vars(args) - self.cmdArgs.update({ - "restart_no_overwrite": [], - }) - self.log=[] - # when making a restart we should not overrite those files (will be checked) - self.restart_no_overwrite=[] - - # forward args for pytest - if args.verbose: - unknown.append("--capture=no") - if args.help_pytest: - unknown.append("--help") - # remaining args - self.remainingArgs=unknown - - def addLog(self, entry): - if self.fake: - self.log.append(entry) - with open(os.path.join(self.problemDir,"testsuite.log.json"), "w+") as fp: - json.dump(self.log, fp, indent='\t') - - def applyConfig(self, config: dict | None = None): - # check args - if config is None: - config = {} - for key, value in config.items(): - if key not in ['ini', 'testfile', 'testname', 'dumpname', 'check_file_produced']: - assert key in self.cmdArgs, f"The given configuration overriding try to set an invalid paramater : {key}={value}" - - # override options - newArgs = {} - newArgs.update(self.cmdArgs) - newArgs.update(config) - - # replace in dict - self.__dict__.update(newArgs) - - def _genCmakeCommand(self,definitionFile=""): - comm=["cmake"] - # add source directory - comm.append(self.idefixDir) - - # we will build in ./build-test so problem dir is parent - comm.append("-DIdefix_PROBLEM_DIR="+self.problemDir) - - # add specific options - for opt in self.cmake: - comm.append("-D"+opt) - - if self.cuda: - comm.append("-DKokkos_ENABLE_CUDA=ON") - # disable fmad operations on Cuda to make it compatible with CPU arithmetics - comm.append("-DIdefix_CXX_FLAGS=--fmad=false") - # disable Async cuda malloc for tests performed on old UCX implementations - comm.append("-DKokkos_ENABLE_IMPL_CUDA_MALLOC_ASYNC=OFF") - - if self.intel: - # disable fmad operations on Cuda to make it compatible with CPU arithmetics - comm.append("-DIdefix_CXX_FLAGS=-fp-model=strict") - comm.append("-DCMAKE_CXX_COMPILER=icpx") - comm.append("-DCMAKE_C_COMPILER=icx") - - if self.hip: - comm.append("-DKokkos_ENABLE_HIP=ON") - # disable fmad operations on HIP to make it compatible with CPU arithmetics - comm.append("-DIdefix_CXX_FLAGS=-ffp-contract=off") - - #if we use single precision - if(self.single): - comm.append("-DIdefix_PRECISION=Single") - else: - comm.append("-DIdefix_PRECISION=Double") - - - if(self.vectPot): - comm.append("-DIdefix_EVOLVE_VECTOR_POTENTIAL=ON") - else: - comm.append("-DIdefix_EVOLVE_VECTOR_POTENTIAL=OFF") - - if(self.Werror): - comm.append("-DIdefix_WERROR=ON") - - # add a definition file if provided - if(definitionFile): - self.definitions=definitionFile - else: - self.definitions="definitions.hpp" - - comm.append("-DIdefix_DEFS="+os.path.join(self.problemDir, self.definitions)) - - if(self.mpi): - comm.append("-DIdefix_MPI=ON") - else: - comm.append("-DIdefix_MPI=OFF") - - if(self.reconstruction == 2): - comm.append("-DIdefix_RECONSTRUCTION=Linear") - elif(self.reconstruction==3): - comm.append("-DIdefix_RECONSTRUCTION=LimO3") - elif(self.reconstruction==4): - comm.append("-DIdefix_RECONSTRUCTION=Parabolic") - - # export ccache env - if self.ccache: - comm.append("-DCMAKE_CXX_COMPILER_LAUNCHER=ccache") - - # ok - return comm - - - def configure(self,definitionFile="", reuse_last_same_build=True, override: dict | None = None): - if override is None: - override = {} - # log - self.addLog({"call": "configure", "args":{ - 'definitionFile': definitionFile, - 'reuse_last_same_build': reuse_last_same_build, - 'override': override - }}) - - # gen command - comm = self._genCmakeCommand(definitionFile) - - # log - print("***************** CALLING CMAKE *******************") - print(f"mkdir -p {self.buildDir}") - print(f"cd {self.buildDir}") - print(' '.join(comm)) - print("***************************************************") - - # not action needed if same command or force no rebuild - if reuse_last_same_build == True and self.lastCmakeCmd == ' '.join(comm): - print("SKIP ALREADY DONE") - return - - # run - try: - # do cleanup - self.clean() +class idfxTest: + def __init__(self, current_test_file, name=""): + parser = argparse.ArgumentParser() + + idefix_dir_env = os.getenv("IDEFIX_DIR") + + parser.add_argument( + "-noplot", + dest="noplot", + help="disable plotting in standard tests", + action="store_false", + ) + + parser.add_argument( + "-ploterr", + help="Enable plotting on error in regression tests", + action="store_true", + ) + + parser.add_argument("-cmake", default=[], help="CMake options", nargs="+") + + parser.add_argument("-definitions", default="", help="definitions.hpp file") + + parser.add_argument( + "-dec", default="", help="MPI domain decomposition", nargs="+" + ) + + parser.add_argument( + "-check", + help="Only perform regression tests without compilation", + action="store_true", + ) + + parser.add_argument( + "-cuda", help="Test on Nvidia GPU using CUDA", action="store_true" + ) + + parser.add_argument( + "-intel", help="Test compiling with Intel OneAPI", action="store_true" + ) + + parser.add_argument( + "-hip", help="Test on AMD GPU using HIP", action="store_true" + ) + + parser.add_argument( + "-single", help="Enable single precision", action="store_true" + ) + + parser.add_argument( + "-vectPot", help="Enable vector potential formulation", action="store_true" + ) + + parser.add_argument( + "-reconstruction", + type=int, + default=2, + help="set reconstruction scheme (2=PLM, 3=LimO3, 4=PPM)", + ) + + parser.add_argument( + "-idefixDir", + default=idefix_dir_env, + help="Set directory for idefix source files (default $IDEFIX_DIR)", + ) + + parser.add_argument("-mpi", help="Enable MPI", action="store_true") + + parser.add_argument( + "-all", + help="Do all test suite (otherwise, just do the test with the current configuration)", + action="store_true", + ) + + parser.add_argument( + "-init", + help="Reinit reference files for non-regression tests (dangerous!)", + action="store_true", + ) + + parser.add_argument( + "-Werror", help="Consider warnings as errors", action="store_true" + ) + + parser.add_argument( + "-ccache", + help="Use ccache to reduce the build time over multiple run of the test suite.", + action="store_true", + ) + + parser.add_argument( + "-restart", + help="Enable creating a restart from a checkpoint.", + action="store_true", + ) + + parser.add_argument( + "-v", + "--verbose", + help="Enable verbose mode, by not capturing the output.", + action="store_true", + ) + + parser.add_argument( + "--help-pytest", + help="Display the options you can transmit directly to pytest in addition to the specific to idefix tests.", + action="store_true", + ) + + parser.add_argument( + "-fake", + help="Make a fake run by just logging the actions to validate that we generate same command over refactoring.", + action="store_true", + ) + + # this option is not used directly by direct users of idfxTest but by idx_test_gen + parser.add_argument( + "-subdir", + default="./test", + help="Select the test in the given subdir not to run all.", + type=str, + ) + + args, unknown = parser.parse_known_args() + + # transform all arguments from args into attributes of this instance + self.__dict__.update(vars(args)) + # store the full path of problem directory + self.currentTestFile = current_test_file + self.currentTestName = name + self.problemDir = os.path.dirname(current_test_file) + self.referenceDirectory = os.path.join(idefix_dir_env, "reference") + # current directory relative to $IDEFIX_DIR/test (used to retrieve the path ot reference files) + self.testDir = os.path.relpath( + self.problemDir, os.path.join(idefix_dir_env, "test") + ) + # build directory, currently inside the test named build-test + self.buildDir = os.path.join(self.problemDir, "build-test") + # remind what build we dit last + self.lastCmakeCmd = "" + # subdir + self.filterSubdir = args.subdir + # save + self.cmdArgs = vars(args) + self.cmdArgs.update({"restart_no_overwrite": []}) + self.log = [] + # when making a restart we should not overrite those files (will be checked) + self.restart_no_overwrite = [] + + # forward args for pytest + if args.verbose: + unknown.append("--capture=no") + if args.help_pytest: + unknown.append("--help") + # remaining args + self.remainingArgs = unknown + + def addLog(self, entry): + if self.fake: + self.log.append(entry) + with open(os.path.join(self.problemDir, "testsuite.log.json"), "w+") as fp: + json.dump(self.log, fp, indent="\t") + + def applyConfig(self, config: dict | None = None): + # check args + if config is None: + config = {} + for key, value in config.items(): + if key not in [ + "ini", + "testfile", + "testname", + "dumpname", + "check_file_produced", + ]: + assert key in self.cmdArgs, ( + f"The given configuration overriding try to set an invalid paramater : {key}={value}" + ) + + # override options + newArgs = {} + newArgs.update(self.cmdArgs) + newArgs.update(config) + + # replace in dict + self.__dict__.update(newArgs) + + def _genCmakeCommand(self, definitionFile=""): + comm = ["cmake"] + # add source directory + comm.append(self.idefixDir) + + # we will build in ./build-test so problem dir is parent + comm.append("-DIdefix_PROBLEM_DIR=" + self.problemDir) + + # add specific options + for opt in self.cmake: + comm.append("-D" + opt) + + if self.cuda: + comm.append("-DKokkos_ENABLE_CUDA=ON") + # disable fmad operations on Cuda to make it compatible with CPU arithmetics + comm.append("-DIdefix_CXX_FLAGS=--fmad=false") + # disable Async cuda malloc for tests performed on old UCX implementations + comm.append("-DKokkos_ENABLE_IMPL_CUDA_MALLOC_ASYNC=OFF") + + if self.intel: + # disable fmad operations on Cuda to make it compatible with CPU arithmetics + comm.append("-DIdefix_CXX_FLAGS=-fp-model=strict") + comm.append("-DCMAKE_CXX_COMPILER=icpx") + comm.append("-DCMAKE_C_COMPILER=icx") + + if self.hip: + comm.append("-DKokkos_ENABLE_HIP=ON") + # disable fmad operations on HIP to make it compatible with CPU arithmetics + comm.append("-DIdefix_CXX_FLAGS=-ffp-contract=off") + + # if we use single precision + if self.single: + comm.append("-DIdefix_PRECISION=Single") + else: + comm.append("-DIdefix_PRECISION=Double") + + if self.vectPot: + comm.append("-DIdefix_EVOLVE_VECTOR_POTENTIAL=ON") + else: + comm.append("-DIdefix_EVOLVE_VECTOR_POTENTIAL=OFF") + + if self.Werror: + comm.append("-DIdefix_WERROR=ON") + + # add a definition file if provided + if definitionFile: + self.definitions = definitionFile + else: + self.definitions = "definitions.hpp" + + comm.append("-DIdefix_DEFS=" + os.path.join(self.problemDir, self.definitions)) + + if self.mpi: + comm.append("-DIdefix_MPI=ON") + else: + comm.append("-DIdefix_MPI=OFF") + + if self.reconstruction == 2: + comm.append("-DIdefix_RECONSTRUCTION=Linear") + elif self.reconstruction == 3: + comm.append("-DIdefix_RECONSTRUCTION=LimO3") + elif self.reconstruction == 4: + comm.append("-DIdefix_RECONSTRUCTION=Parabolic") + + # export ccache env + if self.ccache: + comm.append("-DCMAKE_CXX_COMPILER_LAUNCHER=ccache") + + # ok + return comm + + def configure( + self, + definitionFile="", + reuse_last_same_build=True, + override: dict | None = None, + ): + if override is None: + override = {} + # log + self.addLog( + { + "call": "configure", + "args": { + "definitionFile": definitionFile, + "reuse_last_same_build": reuse_last_same_build, + "override": override, + }, + } + ) + + # gen command + comm = self._genCmakeCommand(definitionFile) + + # log + print("***************** CALLING CMAKE *******************") + print(f"mkdir -p {self.buildDir}") + print(f"cd {self.buildDir}") + print(" ".join(comm)) + print("***************************************************") + # not action needed if same command or force no rebuild + if reuse_last_same_build == True and self.lastCmakeCmd == " ".join(comm): + print("SKIP ALREADY DONE") + return + + # run + try: + # do cleanup + self.clean() + + # log and in fake mode we do not execute + self.addLog({"command": comm}) + + # call cmake + if not self.fake: + cmake = subprocess.run(comm, cwd=os.path.abspath(self.buildDir)) + cmake.check_returncode() + except subprocess.CalledProcessError as e: + print(bcolors.FAIL + "***************************************************") + print("Cmake failed") + print("***************************************************" + bcolors.ENDC) + raise e + finally: + # remind for next time + self.lastCmakeCmd = " ".join(comm) + + def clean(self): # log and in fake mode we do not execute - self.addLog({"command": comm}) - - # call cmake - if not self.fake: - cmake=subprocess.run(comm, cwd=os.path.abspath(self.buildDir)) - cmake.check_returncode() - except subprocess.CalledProcessError as e: - print(bcolors.FAIL+"***************************************************") - print("Cmake failed") - print("***************************************************"+bcolors.ENDC) - raise e - finally: - # remind for next time - self.lastCmakeCmd = ' '.join(comm) - - def clean(self): - # log and in fake mode we do not execute - self.addLog({"call": "clean", "args":{}}) - if self.fake: - return - - # remove the build directory before re-creating it - if os.path.exists(self.buildDir): - shutil.rmtree(self.buildDir) - - # recreate - os.makedirs(self.buildDir, exist_ok=False) - - def compile(self,jobs=8): - self.addLog({"call": "compile", "args":{ - 'jobs': jobs, - }}) - - try: - comm = ["make","-j"+str(jobs)] - self.addLog({"command": comm}) + self.addLog({"call": "clean", "args": {}}) + if self.fake: + return + + # remove the build directory before re-creating it + if os.path.exists(self.buildDir): + shutil.rmtree(self.buildDir) + + # recreate + os.makedirs(self.buildDir, exist_ok=False) + + def compile(self, jobs=8): + self.addLog({"call": "compile", "args": {"jobs": jobs}}) + + try: + comm = ["make", "-j" + str(jobs)] + self.addLog({"command": comm}) + + print("***************************************************") + print(f"cd {os.getcwd()}") + print(" ".join(comm)) + print("***************************************************") + + if not self.fake: + make = subprocess.run(comm, cwd=os.path.abspath(self.buildDir)) + make.check_returncode() + except subprocess.CalledProcessError as e: + print(bcolors.FAIL + "***************************************************") + print("Compilation failed") + print("***************************************************" + bcolors.ENDC) + raise e + + def run(self, inputFile="", np=2, nowrite=False, restart=-1): + # log + self.addLog( + { + "call": "run", + "args": { + "np": np, + "nowrite": nowrite, + "restart": restart, + }, + } + ) + + comm = [os.path.join(self.buildDir, "idefix")] + if inputFile: + comm.append("-i") + comm.append(inputFile) + if self.mpi: + if self.dec: + np = 1 + for n in range(len(self.dec)): + np = np * int(self.dec[n]) + + comm.insert(0, "mpirun") + comm.insert(1, "-np") + comm.insert(2, str(np)) + if self.dec: + comm.append("-dec") + for n in range(len(self.dec)): + comm.append(str(self.dec[n])) + + if nowrite: + comm.append("-nowrite") + + if self.Werror: + comm.append("-Werror") + + if restart >= 0: + comm.append("-restart") + comm.append(str(restart)) print("***************************************************") print(f"cd {os.getcwd()}") - print(' '.join(comm)) + print(" ".join(comm)) print("***************************************************") - if not self.fake: - make=subprocess.run(comm, cwd=os.path.abspath(self.buildDir)) - make.check_returncode() - except subprocess.CalledProcessError as e: - print(bcolors.FAIL+"***************************************************") - print("Compilation failed") - print("***************************************************"+bcolors.ENDC) - raise e - - def run(self, inputFile="", np=2, nowrite=False, restart=-1): - # log - self.addLog({"call": "run", "args":{ - 'np': np, - 'nowrite': nowrite, - 'restart': restart, - }}) - - comm=[os.path.join(self.buildDir,"idefix")] - if inputFile: - comm.append("-i") - comm.append(inputFile) - if self.mpi: - if self.dec: - np=1 - for n in range(len(self.dec)): - np=np*int(self.dec[n]) - - comm.insert(0,"mpirun") - comm.insert(1,"-np") - comm.insert(2,str(np)) - if self.dec: - comm.append("-dec") - for n in range(len(self.dec)): - comm.append(str(self.dec[n])) - - if nowrite: - comm.append("-nowrite") - - if self.Werror: - comm.append("-Werror") - - if restart>=0: - comm.append("-restart") - comm.append(str(restart)) - - print("***************************************************") - print(f"cd {os.getcwd()}") - print(' '.join(comm)) - print("***************************************************") - - try: - self.addLog({"command": comm}) - if not self.fake: - make=subprocess.run(comm, cwd=self.problemDir) - make.check_returncode() - except subprocess.CalledProcessError as e: - print(bcolors.FAIL+"***************************************************") - print("Execution failed") - print("***************************************************"+bcolors.ENDC) - raise e - - self._readLog() - - def _readLog(self): - if self.fake: - return - - if not os.path.exists('./idefix.0.log'): - # When no idefix file is produced, we leave - return - - with open('./idefix.0.log','r') as file: - log = file.read() - - if "SINGLE PRECISION" in log: - self.single = True - else: - self.single = False - - if "Kokkos CUDA target ENABLED" in log: - self.cuda = True - else: - self.cuda = False - - if "Kokkos HIP target ENABLED" in log: - self.hip = True - else: - self.hip = False - - - self.reconstruction = 2 - if "3rd order (LimO3)" in log: - self.reconstruction = 3 - - if "4th order (PPM)" in log: - self.reconstruction = 4 - - self.mpi=False - if "MPI ENABLED" in log: - self.mpi=True - - - # Get input file from log - line = re.search('(?<=Input Parameters using input file )(.*)', log) - self.inifile=line.group(0)[:-1] - - # Get performances from log - line = re.search('Main: Perfs are (.*) cell', log) - self.perf=float(line.group(1)) - - def checkOnly(self, filename, tolerance=0): - # log - self.addLog({"call": "checkOnly", "args":{ - 'filename': filename, - 'tolerance': tolerance, - }}) - - # Assumes the code has been run manually using some configuration, so we simply - # do the test suite witout configure/compile/run - self._readLog() - self._showConfig() - if self.cuda or self.hip: - print(bcolors.WARNING+"***********************************************") - print("WARNING: Idefix guarantees floating point arithmetic accuracy") - print("ONLY when fmad instruction are explicitely disabled at compilation.") - print("Otheriwse, this test will likely fail.") - print("***********************************************"+bcolors.ENDC) - self.standardTest() - self.nonRegressionTest(filename, tolerance) - - def standardTest(self): - # log and in fake mode do not execute. - self.addLog({"call": "standardTest", "args":{}}) - if self.fake: - return - - if os.path.exists(os.path.join(self.problemDir, 'python', 'testidefix.py')): - oldPwd = os.getcwd() - os.chdir(os.path.join(self.problemDir, "python")) - comm = [sys.executable, "testidefix.py"] - if self.noplot: - comm.append("-noplot") - - print(bcolors.OKCYAN+"Running standard test...") - try: - make=subprocess.run(comm) - make.check_returncode() - except subprocess.CalledProcessError as e: - print(bcolors.FAIL+"***************************************************") - print("Standard test execution failed") - print("***************************************************"+bcolors.ENDC) - raise e - print(bcolors.OKCYAN+"Standard test succeeded"+bcolors.ENDC) - os.chdir(oldPwd) - else: - print(bcolors.WARNING+"No standard testidefix.py for this test"+bcolors.ENDC) - sys.stdout.flush() - - def nonRegressionTest(self, filename,tolerance=0): - # log and in fake mode do not execute. - self.addLog({"call": "nonRegressionTest", "args":{ - "filename": filename, - "tolerance": tolerance - }}) - if self.fake: - return - - fileref=os.path.join(self.referenceDirectory, self.testDir, self._getReferenceFilename()) - if not(os.path.exists(fileref)): - raise Exception("Reference file "+fileref+ " doesn't exist") - - filetest=filename - if not(os.path.exists(filetest)): - raise Exception("Test file "+filetest+ " doesn't exist") - - Vref=readDump(fileref) - Vtest=readDump(filetest) - error=self._computeError(Vref,Vtest) - if error > tolerance: - print(bcolors.FAIL+"Non-Regression test failed!") - self._showConfig() - print(bcolors.ENDC) - if self.ploterr: - self._plotDiff(Vref,Vtest) - assert error <= tolerance, bcolors.FAIL+"Error (%e) above tolerance (%e)"%(error,tolerance)+bcolors.ENDC - print(bcolors.OKGREEN+"Non-regression test succeeded with error=%e"%error+bcolors.ENDC) - sys.stdout.flush() - - def compareDump(self, file1, file2,tolerance=0): - self.addLog({"call": "compareDump", "args":{ - "file1": file1, - "file2": file2, - "tolerance": tolerance, - }}) - if self.fake: - return - - Vref=readDump(file1) - Vtest=readDump(file2) - error=self._computeError(Vref,Vtest) - if error > tolerance: - print(bcolors.FAIL+"Files are different !") - print(bcolors.ENDC) - - self._plotDiff(Vref,Vtest) - assert error <= tolerance, bcolors.FAIL+"Error (%e) above tolerance (%e)"%(error,tolerance)+bcolors.ENDC - print(bcolors.OKGREEN+"Files are identical up to error=%e"%error+bcolors.ENDC) - sys.stdout.flush() - - - def makeReference(self,filename): - # log and in fake mode do not execute. - self.addLog({"call": "compareDump", "args":{ - "filename": filename, - }}) - if self.fake: - return - - self._readLog() - targetDir = os.path.join(self.referenceDirectory,self.testDir) - if not os.path.exists(targetDir): - print("Creating reference directory") - os.makedirs(targetDir, exist_ok=True) - fileout = os.path.join(targetDir, self._getReferenceFilename()) - if(os.path.exists(fileout)): - ans=input(bcolors.WARNING+"This will overwrite already existing reference file:\n"+fileout+"\nDo you confirm? (type yes to continue): "+bcolors.ENDC) - if(ans != "yes"): - print(bcolors.WARNING+"Reference creation aborpted"+bcolors.ENDC) - return - - shutil.copy(filename,fileout) - print(bcolors.OKGREEN+"Reference file "+fileout+" created"+bcolors.ENDC) - sys.stdout.flush() - - def _showConfig(self): - print("**************************************************************") - if self.cuda: - print("Nvidia Cuda enabled.") - if self.hip: - print("AMD HIP enabled.") - print("CMake Opts: " +" ".join(self.cmake)) - print("Definitions file:"+self.definitions) - print("Input File: "+self.inifile) - if(self.single): - print("Precision: Single") - else: - print("Precision: Double") - if(self.reconstruction==2): - print("Reconstruction: PLM") - elif(self.reconstruction==3): - print("Reconstruction: LimO3") - elif(self.reconstruction==4): - print("Reconstruction: PPM") - if(self.vectPot): - print("Vector Potential: ON") - else: - print("Vector Potential: OFF") - if self.mpi: - print("MPI: ON") - else: - print("MPI: OFF") - - print("**************************************************************") - - def _getReferenceFilename(self): - strReconstruction="plm" - if self.reconstruction == 3: - strReconstruction = "limo3" - if self.reconstruction == 4: - strReconstruction= "ppm" - - strPrecision="double" - if self.single: - strPrecision="single" - - fileref='dump.ref.'+strPrecision+"."+strReconstruction+"."+self.inifile - if self.vectPot: - fileref=fileref+".vectPot" - - fileref=fileref+'.dmp' - return(fileref) - - def _computeError(self,Vref,Vtest): - ntested=0 - error=0 - for fld in Vtest.data.keys(): - if(Vtest.data[fld].ndim==3): - if fld in Vref.data.keys(): - #print("error in "+fld+" = "+str(np.sqrt(np.mean((Vref.data[fld]-Vtest.data[fld])**2)))) - error = error+np.sqrt(np.mean((Vref.data[fld]-Vtest.data[fld])**2)) - ntested=ntested+1 - - if ntested==0: - raise Exception(bcolors.FAIL+"There is no common field between the reference and current file"+bcolors.ENDC) - - error=error/ntested - return(error) - - def _plotDiff(self,Vref,Vtest): - - for fld in Vtest.data.keys(): - if(Vtest.data[fld].ndim==3): - if fld in Vref.data.keys(): - plt.figure() - plt.title(fld) - x1=Vref.x1 - if Vref.data[fld].shape[0] == Vref.x1.size+1: - x1=np.zeros(Vref.data[fld].shape[0]) - x1[:-1]=Vref.x1l - x1[-1]=Vref.x1r[-1] - x2=Vref.x2 - if Vref.data[fld].shape[1] == Vref.x2.size+1: - x2=np.zeros(Vref.data[fld].shape[1]) - x2[:-1]=Vref.x2l - x2[-1]=Vref.x2r[-1] - x3=Vref.x3 - if Vref.data[fld].shape[2] == Vref.x3.size+1: - x3=np.zeros(Vref.data[fld].shape[2]) - x3[:-1]=Vref.x3l - x3[-1]=Vref.x3r[-1] - if Vref.data[fld].shape[1]>1: - plt.pcolor(x1, x2, Vref.data[fld][:,:,0].T-Vtest.data[fld][:,:,0].T,cmap='seismic') - plt.xlabel("x1") - plt.ylabel("x2") - plt.colorbar() - else: - plt.plot(x1, Vref.data[fld][:,0,0]-Vtest.data[fld][:,0,0]) - plt.xlabel("x1") - plt.show() + try: + self.addLog({"command": comm}) + if not self.fake: + make = subprocess.run(comm, cwd=self.problemDir) + make.check_returncode() + except subprocess.CalledProcessError as e: + print(bcolors.FAIL + "***************************************************") + print("Execution failed") + print("***************************************************" + bcolors.ENDC) + raise e + + self._readLog() + + def _readLog(self): + if self.fake: + return + + if not os.path.exists("./idefix.0.log"): + # When no idefix file is produced, we leave + return + + with open("./idefix.0.log", "r") as file: + log = file.read() + + if "SINGLE PRECISION" in log: + self.single = True + else: + self.single = False + + if "Kokkos CUDA target ENABLED" in log: + self.cuda = True + else: + self.cuda = False + + if "Kokkos HIP target ENABLED" in log: + self.hip = True + else: + self.hip = False + + self.reconstruction = 2 + if "3rd order (LimO3)" in log: + self.reconstruction = 3 + + if "4th order (PPM)" in log: + self.reconstruction = 4 + + self.mpi = False + if "MPI ENABLED" in log: + self.mpi = True + + # Get input file from log + line = re.search("(?<=Input Parameters using input file )(.*)", log) + self.inifile = line.group(0)[:-1] + + # Get performances from log + line = re.search("Main: Perfs are (.*) cell", log) + self.perf = float(line.group(1)) + + def checkOnly(self, filename, tolerance=0): + # log + self.addLog( + { + "call": "checkOnly", + "args": { + "filename": filename, + "tolerance": tolerance, + }, + } + ) + + # Assumes the code has been run manually using some configuration, so we simply + # do the test suite witout configure/compile/run + self._readLog() + self._showConfig() + if self.cuda or self.hip: + print(bcolors.WARNING + "***********************************************") + print("WARNING: Idefix guarantees floating point arithmetic accuracy") + print("ONLY when fmad instruction are explicitely disabled at compilation.") + print("Otheriwse, this test will likely fail.") + print("***********************************************" + bcolors.ENDC) + self.standardTest() + self.nonRegressionTest(filename, tolerance) + + def standardTest(self): + # log and in fake mode do not execute. + self.addLog({"call": "standardTest", "args": {}}) + if self.fake: + return + + if os.path.exists(os.path.join(self.problemDir, "python", "testidefix.py")): + oldPwd = os.getcwd() + os.chdir(os.path.join(self.problemDir, "python")) + comm = [sys.executable, "testidefix.py"] + if self.noplot: + comm.append("-noplot") + + print(bcolors.OKCYAN + "Running standard test...") + try: + make = subprocess.run(comm) + make.check_returncode() + except subprocess.CalledProcessError as e: + print( + bcolors.FAIL + "***************************************************" + ) + print("Standard test execution failed") + print( + "***************************************************" + bcolors.ENDC + ) + raise e + print(bcolors.OKCYAN + "Standard test succeeded" + bcolors.ENDC) + os.chdir(oldPwd) + else: + print( + bcolors.WARNING + + "No standard testidefix.py for this test" + + bcolors.ENDC + ) + sys.stdout.flush() + + def nonRegressionTest(self, filename, tolerance=0): + # log and in fake mode do not execute. + self.addLog( + { + "call": "nonRegressionTest", + "args": {"filename": filename, "tolerance": tolerance}, + } + ) + if self.fake: + return + + fileref = os.path.join( + self.referenceDirectory, self.testDir, self._getReferenceFilename() + ) + if not (os.path.exists(fileref)): + raise Exception("Reference file " + fileref + " doesn't exist") + + filetest = filename + if not (os.path.exists(filetest)): + raise Exception("Test file " + filetest + " doesn't exist") + + Vref = readDump(fileref) + Vtest = readDump(filetest) + error = self._computeError(Vref, Vtest) + if error > tolerance: + print(bcolors.FAIL + "Non-Regression test failed!") + self._showConfig() + print(bcolors.ENDC) + if self.ploterr: + self._plotDiff(Vref, Vtest) + assert error <= tolerance, ( + bcolors.FAIL + + "Error (%e) above tolerance (%e)" % (error, tolerance) + + bcolors.ENDC + ) + print( + bcolors.OKGREEN + + "Non-regression test succeeded with error=%e" % error + + bcolors.ENDC + ) + sys.stdout.flush() + + def compareDump(self, file1, file2, tolerance=0): + self.addLog( + { + "call": "compareDump", + "args": { + "file1": file1, + "file2": file2, + "tolerance": tolerance, + }, + } + ) + if self.fake: + return + + Vref = readDump(file1) + Vtest = readDump(file2) + error = self._computeError(Vref, Vtest) + if error > tolerance: + print(bcolors.FAIL + "Files are different !") + print(bcolors.ENDC) + + self._plotDiff(Vref, Vtest) + assert error <= tolerance, ( + bcolors.FAIL + + "Error (%e) above tolerance (%e)" % (error, tolerance) + + bcolors.ENDC + ) + print( + bcolors.OKGREEN + + "Files are identical up to error=%e" % error + + bcolors.ENDC + ) + sys.stdout.flush() + + def makeReference(self, filename): + # log and in fake mode do not execute. + self.addLog( + { + "call": "compareDump", + "args": {"filename": filename}, + } + ) + if self.fake: + return + + self._readLog() + targetDir = os.path.join(self.referenceDirectory, self.testDir) + if not os.path.exists(targetDir): + print("Creating reference directory") + os.makedirs(targetDir, exist_ok=True) + fileout = os.path.join(targetDir, self._getReferenceFilename()) + if os.path.exists(fileout): + ans = input( + bcolors.WARNING + + "This will overwrite already existing reference file:\n" + + fileout + + "\nDo you confirm? (type yes to continue): " + + bcolors.ENDC + ) + if ans != "yes": + print(bcolors.WARNING + "Reference creation aborpted" + bcolors.ENDC) + return + + shutil.copy(filename, fileout) + print(bcolors.OKGREEN + "Reference file " + fileout + " created" + bcolors.ENDC) + sys.stdout.flush() + + def _showConfig(self): + print("**************************************************************") + if self.cuda: + print("Nvidia Cuda enabled.") + if self.hip: + print("AMD HIP enabled.") + print("CMake Opts: " + " ".join(self.cmake)) + print("Definitions file:" + self.definitions) + print("Input File: " + self.inifile) + if self.single: + print("Precision: Single") + else: + print("Precision: Double") + if self.reconstruction == 2: + print("Reconstruction: PLM") + elif self.reconstruction == 3: + print("Reconstruction: LimO3") + elif self.reconstruction == 4: + print("Reconstruction: PPM") + if self.vectPot: + print("Vector Potential: ON") + else: + print("Vector Potential: OFF") + if self.mpi: + print("MPI: ON") + else: + print("MPI: OFF") + + print("**************************************************************") + + def _getReferenceFilename(self): + strReconstruction = "plm" + if self.reconstruction == 3: + strReconstruction = "limo3" + if self.reconstruction == 4: + strReconstruction = "ppm" + + strPrecision = "double" + if self.single: + strPrecision = "single" + + fileref = ( + "dump.ref." + strPrecision + "." + strReconstruction + "." + self.inifile + ) + if self.vectPot: + fileref = fileref + ".vectPot" + + fileref = fileref + ".dmp" + return fileref + + def _computeError(self, Vref, Vtest): + ntested = 0 + error = 0 + for fld in Vtest.data.keys(): + if Vtest.data[fld].ndim == 3: + if fld in Vref.data.keys(): + # print("error in "+fld+" = "+str(np.sqrt(np.mean((Vref.data[fld]-Vtest.data[fld])**2)))) + error = error + np.sqrt( + np.mean((Vref.data[fld] - Vtest.data[fld]) ** 2) + ) + ntested = ntested + 1 + + if ntested == 0: + raise Exception( + bcolors.FAIL + + "There is no common field between the reference and current file" + + bcolors.ENDC + ) + + error = error / ntested + return error + + def _plotDiff(self, Vref, Vtest): + + for fld in Vtest.data.keys(): + if Vtest.data[fld].ndim == 3: + if fld in Vref.data.keys(): + plt.figure() + plt.title(fld) + x1 = Vref.x1 + if Vref.data[fld].shape[0] == Vref.x1.size + 1: + x1 = np.zeros(Vref.data[fld].shape[0]) + x1[:-1] = Vref.x1l + x1[-1] = Vref.x1r[-1] + x2 = Vref.x2 + if Vref.data[fld].shape[1] == Vref.x2.size + 1: + x2 = np.zeros(Vref.data[fld].shape[1]) + x2[:-1] = Vref.x2l + x2[-1] = Vref.x2r[-1] + x3 = Vref.x3 + if Vref.data[fld].shape[2] == Vref.x3.size + 1: + x3 = np.zeros(Vref.data[fld].shape[2]) + x3[:-1] = Vref.x3l + x3[-1] = Vref.x3r[-1] + if Vref.data[fld].shape[1] > 1: + plt.pcolor( + x1, + x2, + Vref.data[fld][:, :, 0].T - Vtest.data[fld][:, :, 0].T, + cmap="seismic", + ) + plt.xlabel("x1") + plt.ylabel("x2") + plt.colorbar() + else: + plt.plot(x1, Vref.data[fld][:, 0, 0] - Vtest.data[fld][:, 0, 0]) + plt.xlabel("x1") + plt.show() diff --git a/pytools/idfx_test_gen.py b/pytools/idfx_test_gen.py index ada8b962..eaee8291 100644 --- a/pytools/idfx_test_gen.py +++ b/pytools/idfx_test_gen.py @@ -9,256 +9,271 @@ import pytest -DO_NOT_LOOP_ON = ['restart_no_overwrite', "dec", "multirun", "check_file_produced"] +DO_NOT_LOOP_ON = ["restart_no_overwrite", "dec", "multirun", "check_file_produced"] + class IdefixDirTestGenerator: - ''' - Class used to generate the various configuration to run by parsing the - files `testme.py` found in the hierarchy of the /test directory of Idefix. - ''' - - def __init__(self, currentTestFile: str, name: str = ""): - ''' - Constructor or the class. - - Args: - currentTestFile (str): - Path the the current python file. Normally you - simply pass __file__ to this parameter. - name (str): - Define a name for the test, can be empty. - ''' - - self.currentTestFile = currentTestFile - self.currentTestName = name - - # generate the list of configs to run - def genTestConfigs(self, names:str, params, whenClauses = None, defaultConfig: dict | None = None) -> list: - ''' - Generate the the list of configurations as pytest parameters. - It will unpack the configuration set by looping on all combinations defined - by the given sets. - - Args: - names (str): - Comma separated list of variables do consider to build the - name of the file. - params (dict|list): - A configuration set as a dictionnary or a list of - configuration set. - whenCaluses (dict|list): - Provide a set of clauses to apply after unpacking - the configuration so we can patch some values depending on some others. - defaultConfig (dict): - The default configuration on top of which to apply the variants. - - Returns: - A list of pytest.param() ready to be fiven to parametrized pytest functions. - ''' - if defaultConfig is None: - defaultConfig = {} - if whenClauses is None: - whenClauses = {} - # get name ordering list - nameList = names.split(',') - if '' in nameList: - nameList.remove('') - - # gen list of complete configs - all_configs = [] - if isinstance(params, dict): - all_configs += self._genOneConfigSeries(names, params, defaultConfig=defaultConfig) - elif isinstance(params, list): - for p in params: - all_configs += self._genOneConfigSeries(names, p, defaultConfig=defaultConfig) - else: - raise Exception("Should never be called !") - - # convert as parametrize with nice name - result = [] - for config in all_configs: - # append the file - config['testfile'] = self.currentTestFile - config['testname'] = self.currentTestName - # gen name - nameParts = [self.currentTestName] - for name in nameList: - if isinstance(config[name], bool): - if config[name]: - nameParts.append(name) - elif isinstance(config[name], str): - nameParts.append(str(config[name])) + """ + Class used to generate the various configuration to run by parsing the + files `testme.py` found in the hierarchy of the /test directory of Idefix. + """ + + def __init__(self, currentTestFile: str, name: str = ""): + """ + Constructor or the class. + + Args: + currentTestFile (str): + Path the the current python file. Normally you + simply pass __file__ to this parameter. + name (str): + Define a name for the test, can be empty. + """ + + self.currentTestFile = currentTestFile + self.currentTestName = name + + # generate the list of configs to run + def genTestConfigs( + self, names: str, params, whenClauses=None, defaultConfig: dict | None = None + ) -> list: + """ + Generate the the list of configurations as pytest parameters. + It will unpack the configuration set by looping on all combinations defined + by the given sets. + + Args: + names (str): + Comma separated list of variables do consider to build the + name of the file. + params (dict|list): + A configuration set as a dictionnary or a list of + configuration set. + whenCaluses (dict|list): + Provide a set of clauses to apply after unpacking + the configuration so we can patch some values depending on some others. + defaultConfig (dict): + The default configuration on top of which to apply the variants. + + Returns: + A list of pytest.param() ready to be fiven to parametrized pytest functions. + """ + if defaultConfig is None: + defaultConfig = {} + if whenClauses is None: + whenClauses = {} + # get name ordering list + nameList = names.split(",") + if "" in nameList: + nameList.remove("") + + # gen list of complete configs + all_configs = [] + if isinstance(params, dict): + all_configs += self._genOneConfigSeries( + names, params, defaultConfig=defaultConfig + ) + elif isinstance(params, list): + for p in params: + all_configs += self._genOneConfigSeries( + names, p, defaultConfig=defaultConfig + ) + else: + raise Exception("Should never be called !") + + # convert as parametrize with nice name + result = [] + for config in all_configs: + # append the file + config["testfile"] = self.currentTestFile + config["testname"] = self.currentTestName + # gen name + nameParts = [self.currentTestName] + for name in nameList: + if isinstance(config[name], bool): + if config[name]: + nameParts.append(name) + elif isinstance(config[name], str): + nameParts.append(str(config[name])) + else: + nameParts.append(f"{name}-{config[name]}") + confName = "-".join(nameParts) + + # apply when clause + config = self._applyWhen(config, whenClauses) + + result.append(pytest.param(config, id=confName)) + + # ok + return result + + def extractNamingParameters(self, params) -> list: + """ + Loop on the parameters and check automatically what are the list of variable parameters. + + Args: + params (list|dict): + The list of configuration sets to scan or a single set. + + Returns: + The list of names of the variable parameters. + """ + + # if not a list make a list + if not isinstance(params, list): + params = [params] + + # see params + seen = {} + variables = [] + + # loop + for param_set in params: + for key, value in param_set.items(): + if key in variables: + pass + elif key in DO_NOT_LOOP_ON: + pass + elif isinstance(value, list): + variables.append(key) + elif key in seen and seen[key] != value: + variables.append(key) + elif key not in seen: + seen[key] = value + + # by default sort by alphabetic order about var names + # TODO make something better by assigning a priority to vars + variables.sort() + + # ok + return ",".join(variables) + + def _genNextLevelCombinations( + self, input: list, paramName: str, paramValues: list + ) -> list: + """ + Take an input case list and unpack the given parameter to build the new combinations. + + Args: + input (list): + The incoming list of combinations already unpacked before this call. + paramName (str): + Name of the parameter to unpack. + paramValues (list): + The list of values to unpack and to build combinations for. + + Returns: + The updated list of run sets. + """ + result = [] + for entry in input: + for value in paramValues: + v = copy.deepcopy(entry) + v[paramName] = value + result.append(v) + return result + + def _genOneConfigSeries( + self, names: str, config: dict, defaultConfig: dict | None = None + ) -> list: + """ + Generate the the list of configurations as pytest parameters. + It will unpack the configuration set by looping on all combinations defined + by the given sets. + + Args: + names (str): + Comma separated list of variables do consider to build the + name of the file. + config (dict): + A configuration set as a dictionnary. + defaultConfig (dict): + The default configuration to overload with the variant part. + + Returns: + A list of pytest.param() ready to be fiven to parametrized pytest functions. + """ + if defaultConfig is None: + defaultConfig = {} + # get name ordering list + nameList = names.split(",") + + # if there is ini in the list we put it at the end + loopOrder = nameList.copy() + if "ini" in loopOrder: + loopOrder.remove("ini") + loopOrder.append("ini") + if "" in loopOrder: + loopOrder.remove("") + + # init core with everything not a list + core = copy.deepcopy(defaultConfig) + for key, value in config.items(): + if isinstance(value, list) and key not in DO_NOT_LOOP_ON: + assert key in nameList, ( + f"All variable parameteres should be ordered in the names list, '{key}' is not." + ) + else: + core[key] = copy.deepcopy(value) + + # at start we have only default core + result = [core] + + # loop + for key in loopOrder: + value = config[key] + assert isinstance(value, list), ( + f"This parameter is marked as a list but is not a list : {key}={value} !" + ) + result = self._genNextLevelCombinations(result, key, value) + + # ok + return result + + def _matchWhenClause(self, config: dict, when_clause: dict) -> bool: + """ + Check the matching of a given when clause on the given configuration. + + Args: + config (dict): The configuration. + when_clause (dict): The when clause to check. + + Returns: + True if the clause, match, False otherwise. + """ + for key, value in when_clause.items(): + if config[key] != value: + return False + return True + + def _applyWhen(self, config: dict, when: dict) -> dict: + """ + Check if the when clause applies and apply if if any. + + Args: + config (dict): The configuration. + when_clause (dict): The when clause to check and apply. + + Returns: + The fixed config. + """ + # nothing to do + if when == {}: + return config + + # clone + result = copy.deepcopy(config) + + # loop on when + if isinstance(when, list): + for clause in when: + if self._matchWhenClause(config, clause["conditions"]): + result.update(clause["apply"]) + elif isinstance(when, dict): + if self._matchWhenClause(config, when["conditions"]): + result.update(when["apply"]) else: - nameParts.append(f"{name}-{config[name]}") - confName = "-".join(nameParts) - - # apply when clause - config = self._applyWhen(config, whenClauses) - - result.append(pytest.param(config, id=confName)) - - # ok - return result - - def extractNamingParameters(self, params) -> list: - ''' - Loop on the parameters and check automatically what are the list of variable parameters. - - Args: - params (list|dict): - The list of configuration sets to scan or a single set. - - Returns: - The list of names of the variable parameters. - ''' - - # if not a list make a list - if not isinstance(params, list): - params = [params] - - # see params - seen = {} - variables = [] - - # loop - for param_set in params: - for key, value in param_set.items(): - if key in variables: - pass - elif key in DO_NOT_LOOP_ON: - pass - elif isinstance(value, list): - variables.append(key) - elif key in seen and seen[key] != value: - variables.append(key) - elif key not in seen: - seen[key] = value - - # by default sort by alphabetic order about var names - # TODO make something better by assigning a priority to vars - variables.sort() - - # ok - return ','.join(variables) - - def _genNextLevelCombinations(self, input: list, paramName: str, paramValues: list) -> list: - ''' - Take an input case list and unpack the given parameter to build the new combinations. - - Args: - input (list): - The incoming list of combinations already unpacked before this call. - paramName (str): - Name of the parameter to unpack. - paramValues (list): - The list of values to unpack and to build combinations for. - - Returns: - The updated list of run sets. - ''' - result = [] - for entry in input: - for value in paramValues: - v = copy.deepcopy(entry) - v[paramName] = value - result.append(v) - return result - - def _genOneConfigSeries(self, names: str, config: dict, defaultConfig: dict | None = None) -> list: - ''' - Generate the the list of configurations as pytest parameters. - It will unpack the configuration set by looping on all combinations defined - by the given sets. - - Args: - names (str): - Comma separated list of variables do consider to build the - name of the file. - config (dict): - A configuration set as a dictionnary. - defaultConfig (dict): - The default configuration to overload with the variant part. - - Returns: - A list of pytest.param() ready to be fiven to parametrized pytest functions. - ''' - if defaultConfig is None: - defaultConfig = {} - # get name ordering list - nameList = names.split(',') - - # if there is ini in the list we put it at the end - loopOrder = nameList.copy() - if 'ini' in loopOrder: - loopOrder.remove('ini') - loopOrder.append('ini') - if '' in loopOrder: - loopOrder.remove('') - - # init core with everything not a list - core = copy.deepcopy(defaultConfig) - for key, value in config.items(): - if isinstance(value, list) and key not in DO_NOT_LOOP_ON: - assert key in nameList, f"All variable parameteres should be ordered in the names list, '{key}' is not." - else: - core[key] = copy.deepcopy(value) - - # at start we have only default core - result = [core] - - # loop - for key in loopOrder: - value = config[key] - assert isinstance(value, list), f"This parameter is marked as a list but is not a list : {key}={value} !" - result = self._genNextLevelCombinations(result, key, value) - - # ok - return result - - def _matchWhenClause(self, config: dict, when_clause: dict) -> bool: - ''' - Check the matching of a given when clause on the given configuration. - - Args: - config (dict): The configuration. - when_clause (dict): The when clause to check. - - Returns: - True if the clause, match, False otherwise. - ''' - for key, value in when_clause.items(): - if config[key] != value: - return False - return True - - def _applyWhen(self, config: dict, when: dict) -> dict: - ''' - Check if the when clause applies and apply if if any. - - Args: - config (dict): The configuration. - when_clause (dict): The when clause to check and apply. - - Returns: - The fixed config. - ''' - # nothing to do - if when == {}: - return config - - # clone - result = copy.deepcopy(config) - - # loop on when - if isinstance(when, list): - for clause in when: - if self._matchWhenClause(config, clause['conditions']): - result.update(clause['apply']) - elif isinstance(when, dict): - if self._matchWhenClause(config, when['conditions']): - result.update(when['apply']) - else: - raise Exception(f"Invalid type for 'when' : {when}") - - # ok - return result + raise Exception(f"Invalid type for 'when' : {when}") + + # ok + return result diff --git a/pytools/idfx_test_run.py b/pytools/idfx_test_run.py index b96ce840..25071612 100644 --- a/pytools/idfx_test_run.py +++ b/pytools/idfx_test_run.py @@ -21,10 +21,10 @@ @contextmanager def moveInDir(path): - ''' + """ Change current working dir for the given directoy for what is inside the with statement. - ''' + """ oldpwd = os.getcwd() os.chdir(path) try: @@ -32,236 +32,353 @@ def moveInDir(path): finally: os.chdir(oldpwd) + class IdexPytestRunner: - ''' - Implement the Idefix pytest runner to scan and run all the tests described by the files - testme.json into the /test directory of Idefix sources. - ''' - - def __init__(self, parentScriptFile: str): - self.currentTestFile="" - self.currentTestRunner: tst.idfxTest=None - self.parentScritFile=parentScriptFile - self.filterSubdir = os.environ.get("IDEFIX_TEST_FILTER_SUBDIR", "./test/") - - def _validateNaming(self, namings: str, autoExtracted: str, file: str): - names = namings.split(',') - for n in autoExtracted.split(','): - if not n in names: - raise Exception(f"Naming parameter list not match the auto-detectection, some are missinge. You gave : '{names}', detected '{autoExtracted}' into {file}") - - def _makeVariableArgAsList(self, namings: str, variants) -> list: - # make as a list - if isinstance(variants, list) == False: - variants = [variants] - - # split namings - namings = namings.split(",") - - # loop on each - for variant in variants: - for name in namings: - if name in variant and isinstance(variant[name], list) == False: - variant[name] = [variant[name]] - - def genTests(self) -> list: - sourceDir = os.path.dirname(self.parentScritFile) - - # loop over all tests - result = [] - with moveInDir(sourceDir): - # if missing / - if self.filterSubdir != "" and self.filterSubdir[-1] != "/": - self.filterSubdir += "/" - - # walk in test dir to find the tests & sort by name - testfiles = glob.glob(self.filterSubdir + "**/testme.json", recursive=True) - testfiles.sort() - - # loop over each - for testfile in testfiles: - try: - # calc some paths - testfileRelPath = os.path.relpath(testfile, os.path.join(sourceDir, 'test')) - testfilePath = os.path.abspath(os.path.join('test',testfileRelPath)) - testfileDir = os.path.dirname(testfileRelPath) - - # load json & build the inner test combinations - with open(testfilePath, 'r') as fp: - # load - test = json.load(fp) - idefixTestGenerator=IdefixDirTestGenerator(testfilePath, testfileDir) - - # extract default config & when & variants - defaultConfig = test.get('default', {}) - whenClauses = test.get('when', {}) - variants = test.get('variants', {}) - - if 'namings' in test: - namings = test['namings'] - autoExtracted = idefixTestGenerator.extractNamingParameters(variants) - self._validateNaming(namings, autoExtracted, testfilePath) - else: - namings = idefixTestGenerator.extractNamingParameters(variants) - - # required to simplify the algos later, if var is listed as variable, we need to loop over it. - self._makeVariableArgAsList(namings, variants) - - # gen - result += idefixTestGenerator.genTestConfigs(namings, variants, whenClauses=whenClauses, defaultConfig=defaultConfig) - except Exception as e: - raise Exception(f"Fail to generate tests from {testfileRelPath}") from e - - # ok - return result - - def run(self, config: dict) -> None: - # clone before modify to not modity for caller - config = copy.deepcopy(config) - - # print config - print("***************************************************") - print(json.dumps(config, indent='\t')) - print("***************************************************") - - # extract some infos for local usage - testfile = config["testfile"] - dumpname = config['dumpname'] - testname = config["testname"] - tolerance = config.get("tolerance", 0) - definitionFile = config.get("definitionFile", "") - standardTest = config.get("standardTest", True) - nonRegressionTest = config.get("nonRegressionTest", True) - nonRegressionTestIni = config.get("nonRegressionTestIni", None) - check_file_produced = config.get("check_file_produced", []) - problemDir = os.path.dirname(testfile) - - # cleanup some keyword not handled at the - # level of idx_test so we don't perturbate it - del config['dumpname'] - if 'definitionFile' in config: - del config['definitionFile'] - if 'tolerance' in config: - del config['tolerance'] - if 'standardTest' in config: - del config['standardTest'] - if 'nonRegressionTest' in config: - del config['nonRegressionTest'] - if 'nonRegressionTestIni' in config: - del config['nonRegressionTestIni'] - - # if switch from test, rebuild the runner (a runner make for one dir) - if self.currentTestFile != testfile: - self.currentTestRunner = tst.idfxTest(testfile, name=testname) - self.currentTestFile = testfile - - # run - with moveInDir(problemDir): - self._runNonRegression(dumpname, config['ini'], config, tolerance=tolerance, definitionFile=definitionFile, standardTest=standardTest, nonReg=nonRegressionTest, nonRegIni=nonRegressionTestIni) - - # check produced - for file in check_file_produced: - if not os.path.exists(file) and not self.currentTestRunner.fake: - raise Exception(f"Don't find expected file to be produced by the run : {file} !") - - def _runNonRegression(self, dumpname, ini, config_override, tolerance=0, definitionFile="", nonReg=True, nonRegIni=None, standardTest=True, first_run_ini=None,first_run_dumpname=None,configure_and_compile=True): - if 'multirun' in config_override: - self._runNonRegMultirun(dumpname, ini, config_override, tolerance=tolerance, nonReg=nonReg, nonRegIni=nonRegIni, standardTest=standardTest, configure_and_compile=configure_and_compile, definitionFile=definitionFile, first_run_ini=first_run_ini, first_run_dumpname=first_run_dumpname) - else: - # single basic run - self._runNonRegSingleRun(dumpname, ini, config_override, tolerance=tolerance, nonReg=nonReg, standardTest=standardTest, configure_and_compile=configure_and_compile, definitionFile=definitionFile, first_run_ini=first_run_ini, first_run_dumpname=first_run_dumpname) - - def _runNonRegMultirun(self, dumpname, ini, config_override, tolerance=0, definitionFile="", nonReg=True, nonRegIni=None, standardTest=True, first_run_ini=None,first_run_dumpname=None,configure_and_compile=True): - # check - assert 'multirun' in config_override - - # loop over runs - for run in config_override['multirun']: - # copy config - run_config = copy.deepcopy(config_override) - - # patch a bit - del run_config['multirun'] - run_config.update(run) - nonReg=run_config.get('nonRegressionTest', nonReg) - dumpname=run_config.get('dumpname', dumpname) - if 'nonRegressionTest' in run_config: - del run_config['nonRegressionTest'] - standardTest=run_config.get('standardTest', standardTest) - if 'standardTest' in run_config: - del run_config['standardTest'] - - # make single run - self._runNonRegSingleRun(dumpname, run_config['ini'], run_config, definitionFile=definitionFile, tolerance=tolerance, nonReg=nonReg, nonRegIni=nonRegIni, standardTest=standardTest) - - def _runNonRegSingleRun(self, dumpname, ini, config_override, tolerance=0, definitionFile="", nonReg=True, nonRegIni=None, standardTest=True, first_run_ini=None,first_run_dumpname=None,configure_and_compile=True): - # build the runner - idefixTest = self.currentTestRunner - - # handle special override which should not go into - # idefixTest because it is not supported by the idfx_test layer. - config_override = copy.deepcopy(config_override) - nonRegIni = config_override.get("nonRegressionTestIni", nonRegIni) - if 'nonRegressionTestIni' in config_override: - del config_override['nonRegressionTestIni'] - - # apply config - idefixTest.applyConfig(config_override) - - # recompile if needed - if configure_and_compile: - idefixTest.configure(override=config_override, definitionFile=definitionFile) - idefixTest.compile() - - if first_run_ini: - idefixTest.run(inputFile=first_run_ini) - if nonReg: - idefixTest.nonRegressionTest(filename=first_run_dumpname, tolerance=tolerance) - - # Test the restart option - file_mtime={} - if not idefixTest.fake: - for file in idefixTest.restart_no_overwrite: - file_mtime[file] = os.path.getmtime(file) - - # restart - if idefixTest.restart: - restart = 1 - else: - restart = -1 - - # run - idefixTest.run(inputFile=ini, restart=restart) - - # regen ref if needed - if idefixTest.init: - idefixTest.makeReference(filename=dumpname) - - # check outputs - if standardTest: - idefixTest.standardTest() - if nonReg: - if nonRegIni: - idefixTest.inifile = nonRegIni - idefixTest.nonRegressionTest(filename=dumpname, tolerance=tolerance) - - # check that we didn't overrite the file during the restart - if not idefixTest.fake: - for file in idefixTest.restart_no_overwrite: - assert file_mtime[file] == os.path.getmtime(file), f"Dump file {file} was overwritten on restart" - - def main(self, all: bool = False): - if all: - sys.argv.append("-all") - idefixTest = tst.idfxTest(self.parentScritFile, name="main") - os.environ["IDEFIX_TEST_FILTER_SUBDIR"] = idefixTest.filterSubdir - - if idefixTest.all: - pytest.main(['-v', '--no-header', '--junit-xml=idefix-tests.junit.xml', '--tb=short'] + idefixTest.remainingArgs + [self.parentScritFile]) - else: - raise NotImplementedError("Not yet supported !") - #elif self.check: - # idefixTest.checkOnly(filename=dumpname, tolerance=tolerance) - #else: - # for ini in ini_list: - # self.runNonRegression(dumpname, ini, {}, tolerance=tolerance) + """ + Implement the Idefix pytest runner to scan and run all the tests described by the files + testme.json into the /test directory of Idefix sources. + """ + + def __init__(self, parentScriptFile: str): + self.currentTestFile = "" + self.currentTestRunner: tst.idfxTest = None + self.parentScritFile = parentScriptFile + self.filterSubdir = os.environ.get("IDEFIX_TEST_FILTER_SUBDIR", "./test/") + + def _validateNaming(self, namings: str, autoExtracted: str, file: str): + names = namings.split(",") + for n in autoExtracted.split(","): + if not n in names: + raise Exception( + f"Naming parameter list not match the auto-detectection, some are missinge. You gave : '{names}', detected '{autoExtracted}' into {file}" + ) + + def _makeVariableArgAsList(self, namings: str, variants) -> list: + # make as a list + if isinstance(variants, list) == False: + variants = [variants] + + # split namings + namings = namings.split(",") + + # loop on each + for variant in variants: + for name in namings: + if name in variant and isinstance(variant[name], list) == False: + variant[name] = [variant[name]] + + def genTests(self) -> list: + sourceDir = os.path.dirname(self.parentScritFile) + + # loop over all tests + result = [] + with moveInDir(sourceDir): + # if missing / + if self.filterSubdir != "" and self.filterSubdir[-1] != "/": + self.filterSubdir += "/" + + # walk in test dir to find the tests & sort by name + testfiles = glob.glob(self.filterSubdir + "**/testme.json", recursive=True) + testfiles.sort() + + # loop over each + for testfile in testfiles: + try: + # calc some paths + testfileRelPath = os.path.relpath( + testfile, os.path.join(sourceDir, "test") + ) + testfilePath = os.path.abspath( + os.path.join("test", testfileRelPath) + ) + testfileDir = os.path.dirname(testfileRelPath) + + # load json & build the inner test combinations + with open(testfilePath, "r") as fp: + # load + test = json.load(fp) + idefixTestGenerator = IdefixDirTestGenerator( + testfilePath, testfileDir + ) + + # extract default config & when & variants + defaultConfig = test.get("default", {}) + whenClauses = test.get("when", {}) + variants = test.get("variants", {}) + + if "namings" in test: + namings = test["namings"] + autoExtracted = idefixTestGenerator.extractNamingParameters( + variants + ) + self._validateNaming(namings, autoExtracted, testfilePath) + else: + namings = idefixTestGenerator.extractNamingParameters( + variants + ) + + # required to simplify the algos later, if var is listed as variable, we need to loop over it. + self._makeVariableArgAsList(namings, variants) + + # gen + result += idefixTestGenerator.genTestConfigs( + namings, + variants, + whenClauses=whenClauses, + defaultConfig=defaultConfig, + ) + except Exception as e: + raise Exception( + f"Fail to generate tests from {testfileRelPath}" + ) from e + + # ok + return result + + def run(self, config: dict) -> None: + # clone before modify to not modity for caller + config = copy.deepcopy(config) + + # print config + print("***************************************************") + print(json.dumps(config, indent="\t")) + print("***************************************************") + + # extract some infos for local usage + testfile = config["testfile"] + dumpname = config["dumpname"] + testname = config["testname"] + tolerance = config.get("tolerance", 0) + definitionFile = config.get("definitionFile", "") + standardTest = config.get("standardTest", True) + nonRegressionTest = config.get("nonRegressionTest", True) + nonRegressionTestIni = config.get("nonRegressionTestIni", None) + check_file_produced = config.get("check_file_produced", []) + problemDir = os.path.dirname(testfile) + + # cleanup some keyword not handled at the + # level of idx_test so we don't perturbate it + del config["dumpname"] + if "definitionFile" in config: + del config["definitionFile"] + if "tolerance" in config: + del config["tolerance"] + if "standardTest" in config: + del config["standardTest"] + if "nonRegressionTest" in config: + del config["nonRegressionTest"] + if "nonRegressionTestIni" in config: + del config["nonRegressionTestIni"] + + # if switch from test, rebuild the runner (a runner make for one dir) + if self.currentTestFile != testfile: + self.currentTestRunner = tst.idfxTest(testfile, name=testname) + self.currentTestFile = testfile + + # run + with moveInDir(problemDir): + self._runNonRegression( + dumpname, + config["ini"], + config, + tolerance=tolerance, + definitionFile=definitionFile, + standardTest=standardTest, + nonReg=nonRegressionTest, + nonRegIni=nonRegressionTestIni, + ) + + # check produced + for file in check_file_produced: + if not os.path.exists(file) and not self.currentTestRunner.fake: + raise Exception( + f"Don't find expected file to be produced by the run : {file} !" + ) + + def _runNonRegression( + self, + dumpname, + ini, + config_override, + tolerance=0, + definitionFile="", + nonReg=True, + nonRegIni=None, + standardTest=True, + first_run_ini=None, + first_run_dumpname=None, + configure_and_compile=True, + ): + if "multirun" in config_override: + self._runNonRegMultirun( + dumpname, + ini, + config_override, + tolerance=tolerance, + nonReg=nonReg, + nonRegIni=nonRegIni, + standardTest=standardTest, + configure_and_compile=configure_and_compile, + definitionFile=definitionFile, + first_run_ini=first_run_ini, + first_run_dumpname=first_run_dumpname, + ) + else: + # single basic run + self._runNonRegSingleRun( + dumpname, + ini, + config_override, + tolerance=tolerance, + nonReg=nonReg, + standardTest=standardTest, + configure_and_compile=configure_and_compile, + definitionFile=definitionFile, + first_run_ini=first_run_ini, + first_run_dumpname=first_run_dumpname, + ) + + def _runNonRegMultirun( + self, + dumpname, + ini, + config_override, + tolerance=0, + definitionFile="", + nonReg=True, + nonRegIni=None, + standardTest=True, + first_run_ini=None, + first_run_dumpname=None, + configure_and_compile=True, + ): + # check + assert "multirun" in config_override + + # loop over runs + for run in config_override["multirun"]: + # copy config + run_config = copy.deepcopy(config_override) + + # patch a bit + del run_config["multirun"] + run_config.update(run) + nonReg = run_config.get("nonRegressionTest", nonReg) + dumpname = run_config.get("dumpname", dumpname) + if "nonRegressionTest" in run_config: + del run_config["nonRegressionTest"] + standardTest = run_config.get("standardTest", standardTest) + if "standardTest" in run_config: + del run_config["standardTest"] + + # make single run + self._runNonRegSingleRun( + dumpname, + run_config["ini"], + run_config, + definitionFile=definitionFile, + tolerance=tolerance, + nonReg=nonReg, + nonRegIni=nonRegIni, + standardTest=standardTest, + ) + + def _runNonRegSingleRun( + self, + dumpname, + ini, + config_override, + tolerance=0, + definitionFile="", + nonReg=True, + nonRegIni=None, + standardTest=True, + first_run_ini=None, + first_run_dumpname=None, + configure_and_compile=True, + ): + # build the runner + idefixTest = self.currentTestRunner + + # handle special override which should not go into + # idefixTest because it is not supported by the idfx_test layer. + config_override = copy.deepcopy(config_override) + nonRegIni = config_override.get("nonRegressionTestIni", nonRegIni) + if "nonRegressionTestIni" in config_override: + del config_override["nonRegressionTestIni"] + + # apply config + idefixTest.applyConfig(config_override) + + # recompile if needed + if configure_and_compile: + idefixTest.configure( + override=config_override, definitionFile=definitionFile + ) + idefixTest.compile() + + if first_run_ini: + idefixTest.run(inputFile=first_run_ini) + if nonReg: + idefixTest.nonRegressionTest( + filename=first_run_dumpname, tolerance=tolerance + ) + + # Test the restart option + file_mtime = {} + if not idefixTest.fake: + for file in idefixTest.restart_no_overwrite: + file_mtime[file] = os.path.getmtime(file) + + # restart + if idefixTest.restart: + restart = 1 + else: + restart = -1 + + # run + idefixTest.run(inputFile=ini, restart=restart) + + # regen ref if needed + if idefixTest.init: + idefixTest.makeReference(filename=dumpname) + + # check outputs + if standardTest: + idefixTest.standardTest() + if nonReg: + if nonRegIni: + idefixTest.inifile = nonRegIni + idefixTest.nonRegressionTest(filename=dumpname, tolerance=tolerance) + + # check that we didn't overrite the file during the restart + if not idefixTest.fake: + for file in idefixTest.restart_no_overwrite: + assert file_mtime[file] == os.path.getmtime(file), ( + f"Dump file {file} was overwritten on restart" + ) + + def main(self, all: bool = False): + if all: + sys.argv.append("-all") + idefixTest = tst.idfxTest(self.parentScritFile, name="main") + os.environ["IDEFIX_TEST_FILTER_SUBDIR"] = idefixTest.filterSubdir + + if idefixTest.all: + pytest.main( + [ + "-v", + "--no-header", + "--junit-xml=idefix-tests.junit.xml", + "--tb=short", + ] + + idefixTest.remainingArgs + + [self.parentScritFile] + ) + else: + raise NotImplementedError("Not yet supported !") + # elif self.check: + # idefixTest.checkOnly(filename=dumpname, tolerance=tolerance) + # else: + # for ini in ini_list: + # self.runNonRegression(dumpname, ini, {}, tolerance=tolerance) diff --git a/pytools/sod.py b/pytools/sod.py index ba0e1631..3639825f 100644 --- a/pytools/sod.py +++ b/pytools/sod.py @@ -11,7 +11,7 @@ import scipy.optimize -def sound_speed(gamma, pressure, density, dustFrac=0.): +def sound_speed(gamma, pressure, density, dustFrac=0.0): """ Calculate sound speed, scaled by the dust fraction according to: @@ -21,27 +21,28 @@ def sound_speed(gamma, pressure, density, dustFrac=0.): Where :math:`\epsilon` is the dustFrac """ scale = np.sqrt(1 - dustFrac) - return np.sqrt(gamma * pressure/ density) * scale + return np.sqrt(gamma * pressure / density) * scale -def shock_tube_function(p4, p1, p5, rho1, rho5, gamma, dustFrac=0.): + +def shock_tube_function(p4, p1, p5, rho1, rho5, gamma, dustFrac=0.0): """ Shock tube equation """ - z = (p4 / p5 - 1.) + z = p4 / p5 - 1.0 c1 = sound_speed(gamma, p1, rho1, dustFrac) c5 = sound_speed(gamma, p5, rho5, dustFrac) - gm1 = gamma - 1. - gp1 = gamma + 1. - g2 = 2. * gamma + gm1 = gamma - 1.0 + gp1 = gamma + 1.0 + g2 = 2.0 * gamma - fact = gm1 / g2 * (c5 / c1) * z / np.sqrt(1. + gp1 / g2 * z) - fact = (1. - fact) ** (g2 / gm1) + fact = gm1 / g2 * (c5 / c1) * z / np.sqrt(1.0 + gp1 / g2 * z) + fact = (1.0 - fact) ** (g2 / gm1) return p1 * fact - p4 -def calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma=1.4, dustFrac=0.): +def calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma=1.4, dustFrac=0.0): """ Compute regions :rtype : tuple @@ -68,18 +69,18 @@ def calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma=1.4, dustFrac=0.): p4 = scipy.optimize.fsolve(shock_tube_function, p1, (p1, p5, rho1, rho5, gamma))[0] # compute post-shock density and velocity - z = (p4 / p5 - 1.) + z = p4 / p5 - 1.0 c5 = sound_speed(gamma, p5, rho5, dustFrac) - gm1 = gamma - 1. - gp1 = gamma + 1. + gm1 = gamma - 1.0 + gp1 = gamma + 1.0 gmfac1 = 0.5 * gm1 / gamma gmfac2 = 0.5 * gp1 / gamma - fact = np.sqrt(1. + gmfac2 * z) + fact = np.sqrt(1.0 + gmfac2 * z) u4 = c5 * z / (gamma * fact) - rho4 = rho5 * (1. + gmfac2 * z) / (1. + gmfac1 * z) + rho4 = rho5 * (1.0 + gmfac2 * z) / (1.0 + gmfac1 * z) # shock speed w = c5 * fact @@ -87,11 +88,11 @@ def calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma=1.4, dustFrac=0.): # compute values at foot of rarefaction p3 = p4 u3 = u4 - rho3 = rho1 * (p3 / p1)**(1. / gamma) + rho3 = rho1 * (p3 / p1) ** (1.0 / gamma) return (p1, rho1, u1), (p3, rho3, u3), (p4, rho4, u4), (p5, rho5, u5), w -def calc_positions(pl, pr, region1, region3, w, xi, t, gamma, dustFrac=0.): +def calc_positions(pl, pr, region1, region3, w, xi, t, gamma, dustFrac=0.0): """ :return: tuple of positions in the following order -> Head of Rarefaction: xhd, Foot of Rarefaction: xft, @@ -123,21 +124,39 @@ def region_states(pl, pr, region1, region3, region4, region5): where the value is a string, obviously """ if pl > pr: - return {'Region 1': region1, - 'Region 2': 'RAREFACTION', - 'Region 3': region3, - 'Region 4': region4, - 'Region 5': region5} + return { + "Region 1": region1, + "Region 2": "RAREFACTION", + "Region 3": region3, + "Region 4": region4, + "Region 5": region5, + } else: - return {'Region 1': region5, - 'Region 2': region4, - 'Region 3': region3, - 'Region 4': 'RAREFACTION', - 'Region 5': region1} - - -def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, - npts, gamma, t, xi, dustFrac=0.): + return { + "Region 1": region5, + "Region 2": region4, + "Region 3": region3, + "Region 4": "RAREFACTION", + "Region 5": region1, + } + + +def create_arrays( + pl, + pr, + xl, + xr, + positions, + state1, + state3, + state4, + state5, + npts, + gamma, + t, + xi, + dustFrac=0.0, +): """ :return: tuple of x, p, rho and u values across the domain of interest """ @@ -146,8 +165,8 @@ def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, p3, rho3, u3 = state3 p4, rho4, u4 = state4 p5, rho5, u5 = state5 - gm1 = gamma - 1. - gp1 = gamma + 1. + gm1 = gamma - 1.0 + gp1 = gamma + 1.0 x_arr = np.linspace(xl, xr, npts) rho = np.zeros(npts, dtype=float) @@ -161,10 +180,10 @@ def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, p[i] = p1 u[i] = u1 elif x < xft: - u[i] = 2. / gp1 * (c1 + (x - xi) / t) - fact = 1. - 0.5 * gm1 * u[i] / c1 - rho[i] = rho1 * fact ** (2. / gm1) - p[i] = p1 * fact ** (2. * gamma / gm1) + u[i] = 2.0 / gp1 * (c1 + (x - xi) / t) + fact = 1.0 - 0.5 * gm1 * u[i] / c1 + rho[i] = rho1 * fact ** (2.0 / gm1) + p[i] = p1 * fact ** (2.0 * gamma / gm1) elif x < xcd: rho[i] = rho3 p[i] = p3 @@ -192,10 +211,10 @@ def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, p[i] = p3 u[i] = -u3 elif x < xhd: - u[i] = -2. / gp1 * (c1 + (xi - x) / t) - fact = 1. + 0.5 * gm1 * u[i] / c1 - rho[i] = rho1 * fact ** (2. / gm1) - p[i] = p1 * fact ** (2. * gamma / gm1) + u[i] = -2.0 / gp1 * (c1 + (xi - x) / t) + fact = 1.0 + 0.5 * gm1 * u[i] / c1 + rho[i] = rho1 * fact ** (2.0 / gm1) + p[i] = p1 * fact ** (2.0 * gamma / gm1) else: rho[i] = rho1 p[i] = p1 @@ -204,8 +223,7 @@ def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, return x_arr, p, rho, u -def solve(left_state, right_state, geometry, t, gamma=1.4, npts=500, - dustFrac=0.): +def solve(left_state, right_state, geometry, t, gamma=1.4, npts=500, dustFrac=0.0): """ Solves the Sod shock tube problem (i.e. riemann problem) of discontinuity across an interface. @@ -246,34 +264,57 @@ def solve(left_state, right_state, geometry, t, gamma=1.4, npts=500, # basic checking if xl >= xr: - print('xl has to be less than xr!') + print("xl has to be less than xr!") exit() if xi >= xr or xi <= xl: - print('xi has in between xl and xr!') + print("xi has in between xl and xr!") exit() # calculate regions - region1, region3, region4, region5, w = \ - calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma, dustFrac) + region1, region3, region4, region5, w = calculate_regions( + pl, ul, rhol, pr, ur, rhor, gamma, dustFrac + ) regions = region_states(pl, pr, region1, region3, region4, region5) # calculate positions - x_positions = calc_positions(pl, pr, region1, region3, w, xi, t, gamma, - dustFrac) - - pos_description = ('Head of Rarefaction', 'Foot of Rarefaction', - 'Contact Discontinuity', 'Shock') + x_positions = calc_positions(pl, pr, region1, region3, w, xi, t, gamma, dustFrac) + + pos_description = ( + "Head of Rarefaction", + "Foot of Rarefaction", + "Contact Discontinuity", + "Shock", + ) positions = dict(zip(pos_description, x_positions, strict=True)) # create arrays - x, p, rho, u = create_arrays(pl, pr, xl, xr, x_positions, - region1, region3, region4, region5, - npts, gamma, t, xi, dustFrac) - - energy = p/(rho * (gamma - 1.0)) - rho_total = rho/(1.0 - dustFrac) - val_dict = {'x':x, 'p':p, 'rho':rho, 'u':u, 'energy':energy, - 'rho_total':rho_total} + x, p, rho, u = create_arrays( + pl, + pr, + xl, + xr, + x_positions, + region1, + region3, + region4, + region5, + npts, + gamma, + t, + xi, + dustFrac, + ) + + energy = p / (rho * (gamma - 1.0)) + rho_total = rho / (1.0 - dustFrac) + val_dict = { + "x": x, + "p": p, + "rho": rho, + "u": u, + "energy": energy, + "rho_total": rho_total, + } return positions, regions, val_dict diff --git a/pytools/tests/test_idfx_test_gen.py b/pytools/tests/test_idfx_test_gen.py index d92a7d48..11a9603a 100644 --- a/pytools/tests/test_idfx_test_gen.py +++ b/pytools/tests/test_idfx_test_gen.py @@ -11,200 +11,200 @@ def test_extractNamingParameters(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") - - lst = gen.extractNamingParameters([ - { - "dec": [1,2,3], - "a": 10, - "b": 11, - "c": [True, False] - }, - { - "dec": [1,2,4], - "a": 10, - "b": 12, - "d": [True, False] - }, - ]) - - # valid - assert lst == 'b,c,d' + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + lst = gen.extractNamingParameters( + [ + { + "dec": [1, 2, 3], + "a": 10, + "b": 11, + "c": [True, False], + }, + { + "dec": [1, 2, 4], + "a": 10, + "b": 12, + "d": [True, False], + }, + ] + ) + + # valid + assert lst == "b,c,d" + def test_genNextLevelCombinations(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") - - # init with single element - core = {} - init = [core] - - # gen first set of combinations - lst1 = gen._genNextLevelCombinations(init, "mpi", [True, False]) - assert lst1 == [ - {"mpi": True}, - {"mpi": False}, - ] - - # gen first second set of combinations - lst3 = gen._genNextLevelCombinations(lst1, "name", ["a", "b"]) - assert lst3 == [ - {"mpi": True, "name": "a"}, - {"mpi": True, "name": "b"}, - {"mpi": False, "name": "a"}, - {"mpi": False, "name": "b"}, - ] + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + # init with single element + core = {} + init = [core] + + # gen first set of combinations + lst1 = gen._genNextLevelCombinations(init, "mpi", [True, False]) + assert lst1 == [{"mpi": True}, {"mpi": False}] + + # gen first second set of combinations + lst3 = gen._genNextLevelCombinations(lst1, "name", ["a", "b"]) + assert lst3 == [ + {"mpi": True, "name": "a"}, + {"mpi": True, "name": "b"}, + {"mpi": False, "name": "a"}, + {"mpi": False, "name": "b"}, + ] + def test_genOneConfigSeries(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") - - # gen - result = gen._genOneConfigSeries("mpi,name", { - "mpi": [True, False], - "name": ["a", "b"] - }) - - # check - assert result == [ - { - 'mpi': True, - 'name': 'a', - },{ - 'mpi': True, - 'name': 'b', - },{ - 'mpi': False, - 'name': 'a', - },{ - 'mpi': False, - 'name': 'b', - }, - ] + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + # gen + result = gen._genOneConfigSeries( + "mpi,name", {"mpi": [True, False], "name": ["a", "b"]} + ) + + # check + assert result == [ + {"mpi": True, "name": "a"}, + {"mpi": True, "name": "b"}, + {"mpi": False, "name": "a"}, + {"mpi": False, "name": "b"}, + ] + def test_matchWhenClause_single(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + # calls + assert gen._matchWhenClause({"a": 10, "b": 11}, {"a": 10}) == True + assert gen._matchWhenClause({"a": 10, "b": 11}, {"a": 11}) == False + assert gen._matchWhenClause({"a": 10, "b": 11}, {"b": 11}) == True + assert gen._matchWhenClause({"a": 10, "b": 11}, {"b": 10}) == False - # calls - assert gen._matchWhenClause({"a":10, "b": 11}, {"a":10}) == True - assert gen._matchWhenClause({"a":10, "b": 11}, {"a":11}) == False - assert gen._matchWhenClause({"a":10, "b": 11}, {"b":11}) == True - assert gen._matchWhenClause({"a":10, "b": 11}, {"b":10}) == False def test_matchWhenClause_and(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + # calls + assert gen._matchWhenClause({"a": 10, "b": 11}, {"a": 10, "b": 11}) == True + assert gen._matchWhenClause({"a": 10, "b": 11}, {"a": 11, "b": 11}) == False - # calls - assert gen._matchWhenClause({"a":10, "b": 11}, {"a":10, "b":11}) == True - assert gen._matchWhenClause({"a":10, "b": 11}, {"a":11, "b":11}) == False def test_applyWhen_none(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + # call when clause + res = gen._applyWhen({"a": 10, "b": 11}, when={}) - # call when clause - res = gen._applyWhen({ - "a": 10, - "b": 11 - }, when={}) + # check result + assert res == {"a": 10, "b": 11} - # check result - assert res == {"a": 10, "b": 11} def test_applyWhen_single_apply_yes(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") - - # call when clause - res = gen._applyWhen({ - "a": 10, - "b": 11 - }, when={ - "conditions": { - "a": 10, - }, - "apply": { - "b": 12 - } - }) - - # check result - assert res == {"a": 10, "b": 12} + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + # call when clause + res = gen._applyWhen( + {"a": 10, "b": 11}, + when={ + "conditions": {"a": 10}, + "apply": {"b": 12}, + }, + ) + + # check result + assert res == {"a": 10, "b": 12} + def test_applyWhen_single_apply_no(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") - - # call when clause - res = gen._applyWhen({ - "a": 10, - "b": 11 - }, when={ - "conditions": { - "a": 9, - }, - "apply": { - "b": 12 - } - }) - - # check result - assert res == {"a": 10, "b": 11} + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + # call when clause + res = gen._applyWhen( + {"a": 10, "b": 11}, + when={ + "conditions": {"a": 9}, + "apply": {"b": 12}, + }, + ) + + # check result + assert res == {"a": 10, "b": 11} + def test_applyWhen_list_apply_yes(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") - - # call when clause - res = gen._applyWhen({ - "a": 10, - "b": 11 - }, when=[ - { - "conditions": { - "a": 10, - }, - "apply": { - "b": 12 - } - }, - { - "conditions": { - "a": 13, - }, - "apply": { - "b": 13 - } - } - ]) - - # check result - assert res == {"a": 10, "b": 12} + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + # call when clause + res = gen._applyWhen( + {"a": 10, "b": 11}, + when=[ + { + "conditions": {"a": 10}, + "apply": {"b": 12}, + }, + { + "conditions": {"a": 13}, + "apply": {"b": 13}, + }, + ], + ) + + # check result + assert res == {"a": 10, "b": 12} + def test_gen_full(): - # build generator - gen = IdefixDirTestGenerator(__file__, "unit-test") - - # gen - result = gen.genTestConfigs("mpi,name", { - "mpi": [True, False], - "name": ["a", "b"] - }, { - "conditions": { - "mpi": True - }, - "apply": { - "extra": 10, - } - }) - - # check - assert result == [ - pytest.param({'mpi': True, 'name': 'a', 'testfile': __file__, 'testname': 'unit-test', 'extra': 10}, id='unit-test-mpi-a'), - pytest.param({'mpi': True, 'name': 'b', 'testfile': __file__, 'testname': 'unit-test', 'extra': 10}, id='unit-test-mpi-b'), - pytest.param({'mpi': False, 'name': 'a', 'testfile': __file__, 'testname': 'unit-test'}, id='unit-test-a'), - pytest.param({'mpi': False, 'name': 'b', 'testfile': __file__, 'testname': 'unit-test'}, id='unit-test-b'), - ] + # build generator + gen = IdefixDirTestGenerator(__file__, "unit-test") + + # gen + result = gen.genTestConfigs( + "mpi,name", + {"mpi": [True, False], "name": ["a", "b"]}, + { + "conditions": {"mpi": True}, + "apply": {"extra": 10}, + }, + ) + + # check + assert result == [ + pytest.param( + { + "mpi": True, + "name": "a", + "testfile": __file__, + "testname": "unit-test", + "extra": 10, + }, + id="unit-test-mpi-a", + ), + pytest.param( + { + "mpi": True, + "name": "b", + "testfile": __file__, + "testname": "unit-test", + "extra": 10, + }, + id="unit-test-mpi-b", + ), + pytest.param( + {"mpi": False, "name": "a", "testfile": __file__, "testname": "unit-test"}, + id="unit-test-a", + ), + pytest.param( + {"mpi": False, "name": "b", "testfile": __file__, "testname": "unit-test"}, + id="unit-test-b", + ), + ] diff --git a/pytools/tests/test_idfx_test_run.py b/pytools/tests/test_idfx_test_run.py index 58dc31eb..dfe32ac5 100644 --- a/pytools/tests/test_idfx_test_run.py +++ b/pytools/tests/test_idfx_test_run.py @@ -13,79 +13,91 @@ def test_genTests(): - # build runner - runner = IdexPytestRunner(__file__) + # build runner + runner = IdexPytestRunner(__file__) - # dir - dir = os.path.dirname(__file__) + # dir + dir = os.path.dirname(__file__) - # generate - result = runner.genTests() - assert result == [ - pytest.param( - { - 'dumpname': 'dump.0001.dmp', - 'noplot': True, - 'reconstruction': 2, - 'tolerance': 1e-14, - 'ini': 'idefix.ini', - 'testfile': dir + '/test/pb1/testme.json', - 'testname': 'pb1' - }, marks=(), id='pb1-idefix.ini' - ), - pytest.param( - { - 'dumpname': 'dump.0001.dmp', - 'noplot': True, - 'reconstruction': 2, - 'tolerance': 1e-14, - 'ini': 'idefix-implicit.ini', - 'testfile': dir + '/test/pb1/testme.json', - 'testname': 'pb1' - }, marks=(), id='pb1-idefix-implicit.ini' - ), - pytest.param( - { - 'dumpname': 'dump.0001.dmp', - 'noplot': True, - 'reconstruction': 2, - 'tolerance': 1e-14, - 'ini': 'idefix.ini', - 'testfile': dir + '/test/pb2/testme.json', - 'testname': 'pb2' - }, marks=(), id='pb2-idefix.ini-noplot' - ), - pytest.param( - { - 'dumpname': 'dump.0001.dmp', - 'noplot': True, - 'reconstruction': 2, - 'tolerance': 1e-14, - 'ini': 'idefix-implicit.ini', - 'testfile': dir + '/test/pb2/testme.json', - 'testname': 'pb2' - }, marks=(), id='pb2-idefix-implicit.ini-noplot' - ), - pytest.param( - { - 'dumpname': 'dump.0001.dmp', - 'noplot': False, - 'reconstruction': 2, - 'tolerance': 1e-14, - 'ini': 'idefix.ini', - 'testfile': dir + '/test/pb2/testme.json', - 'testname': 'pb2' - }, marks=(), id='pb2-idefix.ini' - ), - pytest.param( - { - 'dumpname': 'dump.0001.dmp', - 'noplot': False, - 'reconstruction': 2, - 'tolerance': 1e-14, - 'ini': 'idefix-implicit.ini', - 'testfile': dir + '/test/pb2/testme.json', - 'testname': 'pb2' - }, marks=(), id='pb2-idefix-implicit.ini' - ), - ] + # generate + result = runner.genTests() + assert result == [ + pytest.param( + { + "dumpname": "dump.0001.dmp", + "noplot": True, + "reconstruction": 2, + "tolerance": 1e-14, + "ini": "idefix.ini", + "testfile": dir + "/test/pb1/testme.json", + "testname": "pb1", + }, + marks=(), + id="pb1-idefix.ini", + ), + pytest.param( + { + "dumpname": "dump.0001.dmp", + "noplot": True, + "reconstruction": 2, + "tolerance": 1e-14, + "ini": "idefix-implicit.ini", + "testfile": dir + "/test/pb1/testme.json", + "testname": "pb1", + }, + marks=(), + id="pb1-idefix-implicit.ini", + ), + pytest.param( + { + "dumpname": "dump.0001.dmp", + "noplot": True, + "reconstruction": 2, + "tolerance": 1e-14, + "ini": "idefix.ini", + "testfile": dir + "/test/pb2/testme.json", + "testname": "pb2", + }, + marks=(), + id="pb2-idefix.ini-noplot", + ), + pytest.param( + { + "dumpname": "dump.0001.dmp", + "noplot": True, + "reconstruction": 2, + "tolerance": 1e-14, + "ini": "idefix-implicit.ini", + "testfile": dir + "/test/pb2/testme.json", + "testname": "pb2", + }, + marks=(), + id="pb2-idefix-implicit.ini-noplot", + ), + pytest.param( + { + "dumpname": "dump.0001.dmp", + "noplot": False, + "reconstruction": 2, + "tolerance": 1e-14, + "ini": "idefix.ini", + "testfile": dir + "/test/pb2/testme.json", + "testname": "pb2", + }, + marks=(), + id="pb2-idefix.ini", + ), + pytest.param( + { + "dumpname": "dump.0001.dmp", + "noplot": False, + "reconstruction": 2, + "tolerance": 1e-14, + "ini": "idefix-implicit.ini", + "testfile": dir + "/test/pb2/testme.json", + "testname": "pb2", + }, + marks=(), + id="pb2-idefix-implicit.ini", + ), + ] diff --git a/pytools/vtk_io.py b/pytools/vtk_io.py index ac4551b2..6d804109 100644 --- a/pytools/vtk_io.py +++ b/pytools/vtk_io.py @@ -97,7 +97,9 @@ def _load_header(self, fh, geometry=None): self.periodicity = np.fromfile(fh, dtype=dint, count=3).astype(bool) elif NATIVE_COORDINATE_REGEXP.match(d): native_name, _ncomp, native_dim, _dtype = d.split() - self.native_coordinates[native_name] = np.fromfile(fh, dtype=dt, count=int(native_dim)) + self.native_coordinates[native_name] = np.fromfile( + fh, dtype=dt, count=int(native_dim) + ) else: warnings.warn("Found unknown field %s" % d, stacklevel=3) fh.readline() # skip extra linefeed (empty line) @@ -211,7 +213,6 @@ def _load_hydro(self, fh): # Reconstruct the polar coordinate system if self.geometry == "polar": - r = np.sqrt(xcart[:, 0, 0] ** 2 + ycart[:, 0, 0] ** 2) theta = np.unwrap(np.arctan2(ycart[0, :, 0], xcart[0, :, 0])) z = zcart[0, 0, :] @@ -377,6 +378,7 @@ def _setup_coordinates_from_native(self): def __repr__(self): return "VTKDataset('%s')" % self.filename + # ////// public API ////// def readVTK(filename, geometry=None): r"""Read a VTK file for any geometry. diff --git a/test.py b/test.py index 081d3166..1bb6768f 100755 --- a/test.py +++ b/test.py @@ -17,11 +17,13 @@ # to build for each run if we just changed the ini file and run options. gblIdefixPytestRunner = IdexPytestRunner(__file__) + # define the pytest test @pytest.mark.parametrize("config", gblIdefixPytestRunner.genTests()) def test_idefix_build_run_check(config): - gblIdefixPytestRunner.run(config) + gblIdefixPytestRunner.run(config) + # if called directly as a script if __name__ == "__main__": - gblIdefixPytestRunner.main(all=True) + gblIdefixPytestRunner.main(all=True) diff --git a/test/Dust/DustEnergy/python/testidefix.py b/test/Dust/DustEnergy/python/testidefix.py index 14b64e00..24112123 100755 --- a/test/Dust/DustEnergy/python/testidefix.py +++ b/test/Dust/DustEnergy/python/testidefix.py @@ -5,6 +5,7 @@ @author: lesurg """ + import argparse import sys @@ -12,47 +13,43 @@ import numpy as np parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() # solve eq. A6 in Riols & Lesur (2018) - - - # load the dat file produced by the setup -raw=np.loadtxt('../timevol.dat',skiprows=1) -t=raw[:,0] -Ekg=raw[:,1] -Ekd=raw[:,2] -Eth=raw[:,3] - -etot=Ekg+Ekd+Eth - -if not(args.noplot): - plt.figure() - plt.plot(t,Ekg,label="Ek_g") - plt.plot(t,Ekd,label="Ek_d") - plt.plot(t,Eth,label="Eth") - plt.plot(t,etot,'--',label="Etot") - plt.legend() - plt.xlabel("t") - plt.ylabel("Flow energy") - plt.show() +raw = np.loadtxt("../timevol.dat", skiprows=1) +t = raw[:, 0] +Ekg = raw[:, 1] +Ekd = raw[:, 2] +Eth = raw[:, 3] + +etot = Ekg + Ekd + Eth + +if not (args.noplot): + plt.figure() + plt.plot(t, Ekg, label="Ek_g") + plt.plot(t, Ekd, label="Ek_d") + plt.plot(t, Eth, label="Eth") + plt.plot(t, etot, "--", label="Etot") + plt.legend() + plt.xlabel("t") + plt.ylabel("Flow energy") + plt.show() # Compute relative evolution of total energy -error=abs((etot[-1]-etot[0])/etot[0]) +error = abs((etot[-1] - etot[0]) / etot[0]) -print("error=%f"%error) -if(error<1e-4): - print("Success!") +print("error=%f" % error) +if error < 1e-4: + print("Success!") else: - print("Failure!") - sys.exit(1) + print("Failure!") + sys.exit(1) diff --git a/test/Dust/DustEnergy/testme.py b/test/Dust/DustEnergy/testme.py index 1d997a4b..473f39db 100755 --- a/test/Dust/DustEnergy/testme.py +++ b/test/Dust/DustEnergy/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,30 +12,31 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-implicit.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-implicit.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init: - test.makeReference(filename=name) - test.standardTest() - test.nonRegressionTest(filename=name) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init: + test.makeReference(filename=name) + test.standardTest() + test.nonRegressionTest(filename=name) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename=name) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name) + else: + testMe(test) else: - test.noplot = True - test.reconstruction=2 - testMe(test) + test.noplot = True + test.reconstruction = 2 + testMe(test) diff --git a/test/Dust/DustyShock/python/testidefix.py b/test/Dust/DustyShock/python/testidefix.py index df575bfe..b87708be 100755 --- a/test/Dust/DustyShock/python/testidefix.py +++ b/test/Dust/DustyShock/python/testidefix.py @@ -19,80 +19,95 @@ from pytools.vtk_io import readVTK parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() -V=readVTK('../data.0001.vtk', geometry='cartesian') +V = readVTK("../data.0001.vtk", geometry="cartesian") # Find where the shock starts -istart = np.argwhere(V.data['RHO'][:,0,0]>2.0)[0][0]-1 +istart = np.argwhere(V.data["RHO"][:, 0, 0] > 2.0)[0][0] - 1 # Compute analytical solution -M=2 # Mach number -K=np.asarray([1,3,5])/M # Drag coefficients +M = 2 # Mach number +K = np.asarray([1, 3, 5]) / M # Drag coefficients + # (56) from Benitez-Llambay 2018 def get_wg(Mach, wi): - coeff=[1,np.sum(wi-1)-1/Mach**2-1,1/Mach**2] - solution=np.roots(coeff) - return(np.min(solution)) + coeff = [1, np.sum(wi - 1) - 1 / Mach**2 - 1, 1 / Mach**2] + solution = np.roots(coeff) + return np.min(solution) + # (55) from Benitez-Llambay 2018 -def rhs(t,wi,Mach,K): - wg = get_wg(Mach,wi) - return K*(wg-wi) +def rhs(t, wi, Mach, K): + wg = get_wg(Mach, wi) + return K * (wg - wi) + -wi0=np.asarray([1.0,1.0,1.0]) -x=V.x[istart:] -arguments=M,K +wi0 = np.asarray([1.0, 1.0, 1.0]) +x = V.x[istart:] +arguments = M, K -sol = solve_ivp(rhs, t_span=[x[0],x[-1]], y0=wi0, t_eval=x,args=arguments) +sol = solve_ivp(rhs, t_span=[x[0], x[-1]], y0=wi0, t_eval=x, args=arguments) # compute gas velocity -wg=np.zeros(len(x)) -for i in range(0,len(wg)): - wg[i]=get_wg(M,sol.y[:,i]) +wg = np.zeros(len(x)) +for i in range(0, len(wg)): + wg[i] = get_wg(M, sol.y[:, i]) -if(not args.noplot): +if not args.noplot: plt.figure(1) - plt.plot(V.x,V.data['RHO'][:,0,0],'o',markersize=4,color='tab:purple') - plt.plot(x,1/wg,color='tab:purple') - plt.plot(V.x,V.data['Dust0_RHO'][:,0,0],'o',markersize=4,color='tab:orange') - plt.plot(x,1/sol.y[0,:],color='tab:orange') - plt.plot(V.x,V.data['Dust1_RHO'][:,0,0],'o',markersize=4,color='tab:blue') - plt.plot(x,1/sol.y[1,:],color='tab:blue') - plt.plot(V.x,V.data['Dust2_RHO'][:,0,0],'o',markersize=4,color='tab:green') - plt.plot(x,1/sol.y[2,:],color='tab:green') - plt.xlim([0,10]) - plt.title('Density') + plt.plot(V.x, V.data["RHO"][:, 0, 0], "o", markersize=4, color="tab:purple") + plt.plot(x, 1 / wg, color="tab:purple") + plt.plot(V.x, V.data["Dust0_RHO"][:, 0, 0], "o", markersize=4, color="tab:orange") + plt.plot(x, 1 / sol.y[0, :], color="tab:orange") + plt.plot(V.x, V.data["Dust1_RHO"][:, 0, 0], "o", markersize=4, color="tab:blue") + plt.plot(x, 1 / sol.y[1, :], color="tab:blue") + plt.plot(V.x, V.data["Dust2_RHO"][:, 0, 0], "o", markersize=4, color="tab:green") + plt.plot(x, 1 / sol.y[2, :], color="tab:green") + plt.xlim([0, 10]) + plt.title("Density") plt.figure(2) - plt.plot(V.x,V.data['VX1'][:,0,0],'o',markersize=4,color='tab:purple') - plt.plot(x,M*wg,color='tab:purple') - plt.plot(V.x[::4],V.data['Dust0_VX1'][::4,0,0],'o',markersize=4,color='tab:orange') - plt.plot(x,M*sol.y[0,:],color='tab:orange') - plt.plot(V.x[::4],V.data['Dust1_VX1'][::4,0,0],'o',markersize=4,color='tab:blue') - plt.plot(x,M*sol.y[1,:],color='tab:blue') - plt.plot(V.x[::4],V.data['Dust2_VX1'][::4,0,0],'o',markersize=4,color='tab:green') - plt.plot(x,M*sol.y[2,:],color='tab:green') - - plt.xlim([0,10]) - plt.title('Velocity') + plt.plot(V.x, V.data["VX1"][:, 0, 0], "o", markersize=4, color="tab:purple") + plt.plot(x, M * wg, color="tab:purple") + plt.plot( + V.x[::4], V.data["Dust0_VX1"][::4, 0, 0], "o", markersize=4, color="tab:orange" + ) + plt.plot(x, M * sol.y[0, :], color="tab:orange") + plt.plot( + V.x[::4], V.data["Dust1_VX1"][::4, 0, 0], "o", markersize=4, color="tab:blue" + ) + plt.plot(x, M * sol.y[1, :], color="tab:blue") + plt.plot( + V.x[::4], V.data["Dust2_VX1"][::4, 0, 0], "o", markersize=4, color="tab:green" + ) + plt.plot(x, M * sol.y[2, :], color="tab:green") + + plt.xlim([0, 10]) + plt.title("Velocity") plt.ioff() plt.show() # Compute the error on the dust densities -error=np.mean(np.fabs(V.data['Dust0_RHO'][istart:,0,0]-1/sol.y[0,:]) + np.fabs(V.data['Dust1_RHO'][istart:,0,0]-1/sol.y[1,:]) + np.fabs(V.data['Dust2_RHO'][istart:,0,0]-1/sol.y[2,:])) / 3 -print("Error=%e"%error) -if error<3.6e-2: +error = ( + np.mean( + np.fabs(V.data["Dust0_RHO"][istart:, 0, 0] - 1 / sol.y[0, :]) + + np.fabs(V.data["Dust1_RHO"][istart:, 0, 0] - 1 / sol.y[1, :]) + + np.fabs(V.data["Dust2_RHO"][istart:, 0, 0] - 1 / sol.y[2, :]) + ) + / 3 +) +print("Error=%e" % error) +if error < 3.6e-2: print("SUCCESS!") sys.exit(0) else: diff --git a/test/Dust/DustyShock/testme.py b/test/Dust/DustyShock/testme.py index b3cfab71..2aecf025 100755 --- a/test/Dust/DustyShock/testme.py +++ b/test/Dust/DustyShock/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,30 +12,31 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-implicit.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-implicit.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init: - test.makeReference(filename=name) - test.standardTest() - test.nonRegressionTest(filename=name,tolerance=1e-14) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init: + test.makeReference(filename=name) + test.standardTest() + test.nonRegressionTest(filename=name, tolerance=1e-14) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename=name) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name) + else: + testMe(test) else: - test.noplot = True - test.reconstruction=2 - testMe(test) + test.noplot = True + test.reconstruction = 2 + testMe(test) diff --git a/test/Dust/DustyWave/python/testidefix.py b/test/Dust/DustyWave/python/testidefix.py index 1c386e6b..9e05128b 100755 --- a/test/Dust/DustyWave/python/testidefix.py +++ b/test/Dust/DustyWave/python/testidefix.py @@ -5,6 +5,7 @@ @author: lesurg """ + import argparse import sys @@ -12,55 +13,56 @@ import numpy as np parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() # solve eq. A6 in Riols & Lesur (2018) -cs=1.0 -chi=1.0 -kx=2.0*np.pi -taus=1.0 -poly=np.zeros(4,dtype=complex) -poly[0]=1 -poly[1]=1j/taus*(1+chi) -poly[2]=-kx**2*cs**2 -poly[3]=-1j*kx**2*cs**2/taus -sol=np.roots(poly) +cs = 1.0 +chi = 1.0 +kx = 2.0 * np.pi +taus = 1.0 +poly = np.zeros(4, dtype=complex) +poly[0] = 1 +poly[1] = 1j / taus * (1 + chi) +poly[2] = -(kx**2) * cs**2 +poly[3] = -1j * kx**2 * cs**2 / taus +sol = np.roots(poly) # Get the minimal decay rate (this should be the one that pops up) -tau=np.amax(np.imag(sol)) +tau = np.amax(np.imag(sol)) # load the dat file produced by the setup -raw=np.loadtxt('../timevol.dat',skiprows=1) -t=raw[:,0] -vx2=raw[:,1] -rho2=raw[:,2] +raw = np.loadtxt("../timevol.dat", skiprows=1) +t = raw[:, 0] +vx2 = raw[:, 1] +rho2 = raw[:, 2] -etot=vx2+rho2 +etot = vx2 + rho2 -if not(args.noplot): - plt.figure() - plt.semilogy(t,etot,label="idefix") - plt.semilogy(t,np.exp(2*tau*t)*etot[10],'--',label="theoretical decay rate") - plt.legend() - plt.xlabel("t") - plt.ylabel("Wave energy") - plt.show() +if not (args.noplot): + plt.figure() + plt.semilogy(t, etot, label="idefix") + plt.semilogy( + t, np.exp(2 * tau * t) * etot[10], "--", label="theoretical decay rate" + ) + plt.legend() + plt.xlabel("t") + plt.ylabel("Wave energy") + plt.show() # Compute decay rate: -tau_measured=t[-1]/(2*np.log(etot[-1]/etot[0])) +tau_measured = t[-1] / (2 * np.log(etot[-1] / etot[0])) # error on the decay rate: -error=np.abs((tau_measured-tau)/tau) +error = np.abs((tau_measured - tau) / tau) -print("error=%f"%error) -if(error<0.07): - print("Success!") +print("error=%f" % error) +if error < 0.07: + print("Success!") else: - print("Failure!") - sys.exit(1) + print("Failure!") + sys.exit(1) diff --git a/test/Dust/DustyWave/testme.py b/test/Dust/DustyWave/testme.py index b3cfab71..2aecf025 100755 --- a/test/Dust/DustyWave/testme.py +++ b/test/Dust/DustyWave/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,30 +12,31 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-implicit.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-implicit.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init: - test.makeReference(filename=name) - test.standardTest() - test.nonRegressionTest(filename=name,tolerance=1e-14) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init: + test.makeReference(filename=name) + test.standardTest() + test.nonRegressionTest(filename=name, tolerance=1e-14) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename=name) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name) + else: + testMe(test) else: - test.noplot = True - test.reconstruction=2 - testMe(test) + test.noplot = True + test.reconstruction = 2 + testMe(test) diff --git a/test/HD/FargoPlanet/testme.py b/test/HD/FargoPlanet/testme.py index 5a036444..7be82a4b 100755 --- a/test/HD/FargoPlanet/testme.py +++ b/test/HD/FargoPlanet/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,39 +12,41 @@ import pytools.idfx_test as tst -tolerance=1e-13 -def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-rkl.ini"] - mytol=tolerance - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename="dump.0001.dmp") - if ini == "idefix-rkl.ini": - mytol = 1e-12 - else: - mytol = tolerance - test.standardTest() - test.nonRegressionTest(filename="dump.0001.dmp",tolerance=mytol) +tolerance = 1e-13 -test=tst.idfxTest(__file__) +def testMe(test): + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-rkl.ini"] + mytol = tolerance + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename="dump.0001.dmp") + if ini == "idefix-rkl.ini": + mytol = 1e-12 + else: + mytol = tolerance + test.standardTest() + test.nonRegressionTest(filename="dump.0001.dmp", tolerance=mytol) + + +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','2'] + test.dec = ["2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp",tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp", tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - test.mpi=True - testMe(test) + test.noplot = True + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + test.mpi = True + testMe(test) diff --git a/test/HD/MachReflection/testme.py b/test/HD/MachReflection/testme.py index cc8cf72e..7007c6e3 100755 --- a/test/HD/MachReflection/testme.py +++ b/test/HD/MachReflection/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -13,33 +14,33 @@ def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-hll.ini","idefix-hllc.ini","idefix-tvdlf.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-hll.ini", "idefix-hllc.ini", "idefix-tvdlf.ini"] - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename="dump.0001.dmp") - test.standardTest() - test.nonRegressionTest(filename="dump.0001.dmp") + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename="dump.0001.dmp") + test.standardTest() + test.nonRegressionTest(filename="dump.0001.dmp") -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','2'] + test.dec = ["2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp") - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp") + else: + testMe(test) else: - test.noplot = True - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - test.mpi=True - testMe(test) + test.noplot = True + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + test.mpi = True + testMe(test) diff --git a/test/HD/SedovBlastWave/python/testidefix.py b/test/HD/SedovBlastWave/python/testidefix.py index 40d39a79..e926a19d 100755 --- a/test/HD/SedovBlastWave/python/testidefix.py +++ b/test/HD/SedovBlastWave/python/testidefix.py @@ -19,63 +19,61 @@ from pytools.vtk_io import readVTK parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() - -V=readVTK('../data.0001.vtk') -R=readVTK('1Dsolution/1Dsolution.vtk') +args, unknown = parser.parse_known_args() +V = readVTK("../data.0001.vtk") +R = readVTK("1Dsolution/1Dsolution.vtk") # Finally, let's plot solutions -p = R.data['PRS'][:,0,0] -x= R.r +p = R.data["PRS"][:, 0, 0] +x = R.r -if V.geometry=="cartesian": - [x3D, y3D, z3D]= np.meshgrid(V.x,V.y,V.z,indexing='ij') +if V.geometry == "cartesian": + [x3D, y3D, z3D] = np.meshgrid(V.x, V.y, V.z, indexing="ij") -elif V.geometry=="spherical": - [r, theta, phi] = np.meshgrid(V.r,V.theta,V.phi,indexing='ij') - x0=0.1 - y0=0 - z0=1.0 +elif V.geometry == "spherical": + [r, theta, phi] = np.meshgrid(V.r, V.theta, V.phi, indexing="ij") + x0 = 0.1 + y0 = 0 + z0 = 1.0 - x3D=r*np.sin(theta)*np.cos(phi)-x0 - y3D=r*np.sin(theta)*np.sin(phi)-y0 - z3D=r*np.cos(theta)-z0 + x3D = r * np.sin(theta) * np.cos(phi) - x0 + y3D = r * np.sin(theta) * np.sin(phi) - y0 + z3D = r * np.cos(theta) - z0 else: - print("Unknown geometry "+V.geometry) - sys.exit(1) + print("Unknown geometry " + V.geometry) + sys.exit(1) -r3D=np.sqrt(x3D**2+y3D**2+z3D**2).flatten() -p3D=V.data['PRS'].flatten() +r3D = np.sqrt(x3D**2 + y3D**2 + z3D**2).flatten() +p3D = V.data["PRS"].flatten() -index=np.argwhere((r3D>0.4) & (r3D<0.458)) +index = np.argwhere((r3D > 0.4) & (r3D < 0.458)) -solinterp=interp1d(x,p) +solinterp = interp1d(x, p) -if(not args.noplot): +if not args.noplot: plt.figure(1) - plt.plot(r3D.flatten(),p3D.flatten(),'+',markersize=2) - plt.plot(x,p) - plt.title('Pressure') + plt.plot(r3D.flatten(), p3D.flatten(), "+", markersize=2) + plt.plot(x, p) + plt.title("Pressure") plt.ioff() plt.show() -error=np.mean(np.fabs(p3D[index]-solinterp(r3D[index]))) -print("Error=%e"%error) +error = np.mean(np.fabs(p3D[index] - solinterp(r3D[index]))) +print("Error=%e" % error) erroref = 7.7e-2 if V.geometry == "spherical": erroref = 1.9e-1 -if errorV.x.size): - x=V.y - bxSim=V.data['BX2'].flatten() - bySim=V.data['BX3'].flatten() -elif(V.z.size>V.x.size): - x=V.z - bxSim=V.data['BX3'].flatten() - bySim=V.data['BX1'].flatten() +V = readVTK("../data.0001.vtk", geometry="cartesian") + +DSim = V.data["RHO"].flatten() +if V.y.size > V.x.size: + x = V.y + bxSim = V.data["BX2"].flatten() + bySim = V.data["BX3"].flatten() +elif V.z.size > V.x.size: + x = V.z + bxSim = V.data["BX3"].flatten() + bySim = V.data["BX1"].flatten() else: - x=V.x - bxSim=V.data['BX1'].flatten() - bySim=V.data['BX2'].flatten() + x = V.x + bxSim = V.data["BX1"].flatten() + bySim = V.data["BX2"].flatten() -bSim=bySim/np.sqrt(bySim[0]**2+bxSim[0]**2) +bSim = bySim / np.sqrt(bySim[0] ** 2 + bxSim[0] ** 2) # Index where the two solutions are assumed to match -iref=np.argwhere(DSim>5.0)[0][0] +iref = np.argwhere(DSim > 5.0)[0][0] -#spatial index where we cut the solution -iend=np.argwhere(x>xcut)[0][0] +# spatial index where we cut the solution +iend = np.argwhere(x > xcut)[0][0] # compute analytical solution -r=ode(f).set_integrator('vode',rtol=1e-6) -r.set_initial_value(DSim[iref],x[iref]) +r = ode(f).set_integrator("vode", rtol=1e-6) +r.set_initial_value(DSim[iref], x[iref]) -Dth=np.zeros(x.shape) +Dth = np.zeros(x.shape) for i in range(x.size): - r.set_initial_value(DSim[iref],x[iref]) - Dth[i]=r.integrate(x[i]).item() - #print("Dth[%d]=%g"%(i,Dth[i])) + r.set_initial_value(DSim[iref], x[iref]) + Dth[i] = r.integrate(x[i]).item() + # print("Dth[%d]=%g"%(i,Dth[i])) -bTh=np.sqrt(b0**2+2*A**2*(Dth-1)*(1/Dth-1/M**2)) +bTh = np.sqrt(b0**2 + 2 * A**2 * (Dth - 1) * (1 / Dth - 1 / M**2)) -errb=(bSim[:iend]-bTh[:iend])/np.amax(bTh[:iend]) -errD=(DSim[:iend]-Dth[:iend])/np.amax(Dth[:iend]) +errb = (bSim[:iend] - bTh[:iend]) / np.amax(bTh[:iend]) +errD = (DSim[:iend] - Dth[:iend]) / np.amax(Dth[:iend]) -if(not args.noplot): - plt.close('all') +if not args.noplot: + plt.close("all") plt.figure() plt.subplot(211) - plt.plot(x[:iend]/L,bSim[:iend],'--',label='Sim') - plt.plot(x[:iend]/L,bTh[:iend],label='Theoretical') - plt.xlabel('x/L') - plt.ylabel('b') + plt.plot(x[:iend] / L, bSim[:iend], "--", label="Sim") + plt.plot(x[:iend] / L, bTh[:iend], label="Theoretical") + plt.xlabel("x/L") + plt.ylabel("b") plt.subplot(212) - plt.plot(x[:iend]/L,errb,'--',label='Error') - plt.xlabel('x/L') - plt.ylabel('b') + plt.plot(x[:iend] / L, errb, "--", label="Error") + plt.xlabel("x/L") + plt.ylabel("b") plt.legend() plt.figure() plt.subplot(211) - plt.plot(x[:iend]/L,DSim[:iend],label='Sim') - plt.plot(x[:iend]/L,Dth[:iend],'--',label='Theoretical') - plt.xlabel('x/L') - plt.ylabel('D') + plt.plot(x[:iend] / L, DSim[:iend], label="Sim") + plt.plot(x[:iend] / L, Dth[:iend], "--", label="Theoretical") + plt.xlabel("x/L") + plt.ylabel("D") plt.subplot(212) - plt.plot(x[:iend]/L,errD,'--',label='Error') - plt.xlabel('x/L') - plt.ylabel('D') + plt.plot(x[:iend] / L, errD, "--", label="Error") + plt.xlabel("x/L") + plt.ylabel("D") plt.legend() plt.ioff() plt.show() -err=np.mean(0.5*np.sqrt(errb**2+errD**2)) -print("Error total=%e"%err) -if(err<5e-3): +err = np.mean(0.5 * np.sqrt(errb**2 + errD**2)) +print("Error total=%e" % err) +if err < 5e-3: print("Success!") sys.exit(0) else: diff --git a/test/MHD/AmbipolarCshock/testme.py b/test/MHD/AmbipolarCshock/testme.py index b8bd6945..1b3f843f 100755 --- a/test/MHD/AmbipolarCshock/testme.py +++ b/test/MHD/AmbipolarCshock/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,30 +12,31 @@ import pytools.idfx_test as tst -tolerance=0 +tolerance = 0 + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-rkl.ini"] - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename="dump.0001.dmp") - test.standardTest() - test.nonRegressionTest(filename="dump.0001.dmp",tolerance=tolerance) + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-rkl.ini"] + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename="dump.0001.dmp") + test.standardTest() + test.nonRegressionTest(filename="dump.0001.dmp", tolerance=tolerance) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp",tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp", tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) + test.noplot = True + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) diff --git a/test/MHD/AmbipolarCshock3D/python/testidefix.py b/test/MHD/AmbipolarCshock3D/python/testidefix.py index c151dd78..2e62c09d 100755 --- a/test/MHD/AmbipolarCshock3D/python/testidefix.py +++ b/test/MHD/AmbipolarCshock3D/python/testidefix.py @@ -5,6 +5,7 @@ @author: lesurg """ + import os import sys @@ -18,106 +19,107 @@ from pytools.vtk_io import readVTK parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() # Parameters from Maclow+ 1995 -theta=np.pi/4 -M=50 -b0=np.sin(theta) -A=10 -L=1 +theta = np.pi / 4 +M = 50 +b0 = np.sin(theta) +A = 10 +L = 1 + +xcut = 40.0 + -xcut=40.0 +def f(x, D): + b = np.sqrt(b0**2 + 2 * A**2 * (D - 1) * (1 / D - 1 / M**2)) + diff = (b - b0) / A**2 * np.cos(theta) ** 2 + np.sin(theta) + diff = (b - D * diff) / (b**2 + np.cos(theta) ** 2) + diff = (b / A) * diff / (L * (1 / D**2 - 1 / M**2)) -def f(x,D): - b=np.sqrt(b0**2+2*A**2*(D-1)*(1/D-1/M**2)) - diff=(b-b0)/A**2*np.cos(theta)**2+np.sin(theta) - diff=(b-D*diff)/(b**2+np.cos(theta)**2) - diff=(b/A)*diff/(L*(1/D**2-1/M**2)) + return diff - return(diff) # load solution -V=readVTK('../data.0001.vtk', geometry='cartesian') - - -if(V.y.size>V.x.size): - x=V.y - DSim=V.data['RHO'][0,:,0] - bxSim=V.data['BX2'][0,:,0] - bySim=V.data['BX3'][0,:,0] -elif(V.z.size>V.x.size): - x=V.z - DSim=V.data['RHO'][0,0,:] - bxSim=V.data['BX3'][0,0,:] - bySim=V.data['BX1'][0,0,:] +V = readVTK("../data.0001.vtk", geometry="cartesian") + + +if V.y.size > V.x.size: + x = V.y + DSim = V.data["RHO"][0, :, 0] + bxSim = V.data["BX2"][0, :, 0] + bySim = V.data["BX3"][0, :, 0] +elif V.z.size > V.x.size: + x = V.z + DSim = V.data["RHO"][0, 0, :] + bxSim = V.data["BX3"][0, 0, :] + bySim = V.data["BX1"][0, 0, :] else: - x=V.x - DSim=V.data['RHO'][:,0,0] - bxSim=V.data['BX1'][:,0,0] - bySim=V.data['BX2'][:,0,0] + x = V.x + DSim = V.data["RHO"][:, 0, 0] + bxSim = V.data["BX1"][:, 0, 0] + bySim = V.data["BX2"][:, 0, 0] -bSim=bySim/np.sqrt(bySim[0]**2+bxSim[0]**2) +bSim = bySim / np.sqrt(bySim[0] ** 2 + bxSim[0] ** 2) # Index where the two solutions are assumed to match -iref=np.argwhere(DSim>5.0)[0][0] +iref = np.argwhere(DSim > 5.0)[0][0] -#spatial index where we cut the solution -iend=np.argwhere(x>xcut)[0][0] +# spatial index where we cut the solution +iend = np.argwhere(x > xcut)[0][0] # compute analytical solution -r=ode(f).set_integrator('vode',rtol=1e-6) -r.set_initial_value(DSim[iref],x[iref]) +r = ode(f).set_integrator("vode", rtol=1e-6) +r.set_initial_value(DSim[iref], x[iref]) -Dth=np.zeros(x.shape) +Dth = np.zeros(x.shape) for i in range(x.size): - r.set_initial_value(DSim[iref],x[iref]) - Dth[i]=r.integrate(x[i]).item() - #print("Dth[%d]=%g"%(i,Dth[i])) + r.set_initial_value(DSim[iref], x[iref]) + Dth[i] = r.integrate(x[i]).item() + # print("Dth[%d]=%g"%(i,Dth[i])) -bTh=np.sqrt(b0**2+2*A**2*(Dth-1)*(1/Dth-1/M**2)) +bTh = np.sqrt(b0**2 + 2 * A**2 * (Dth - 1) * (1 / Dth - 1 / M**2)) -errb=(bSim[:iend]-bTh[:iend])/np.amax(bTh[:iend]) -errD=(DSim[:iend]-Dth[:iend])/np.amax(Dth[:iend]) +errb = (bSim[:iend] - bTh[:iend]) / np.amax(bTh[:iend]) +errD = (DSim[:iend] - Dth[:iend]) / np.amax(Dth[:iend]) -if(not args.noplot): - plt.close('all') +if not args.noplot: + plt.close("all") plt.figure() plt.subplot(211) - plt.plot(x[:iend]/L,bSim[:iend],'--',label='Sim') - plt.plot(x[:iend]/L,bTh[:iend],label='Theoretical') - plt.xlabel('x/L') - plt.ylabel('b') + plt.plot(x[:iend] / L, bSim[:iend], "--", label="Sim") + plt.plot(x[:iend] / L, bTh[:iend], label="Theoretical") + plt.xlabel("x/L") + plt.ylabel("b") plt.subplot(212) - plt.plot(x[:iend]/L,errb,'--',label='Error') - plt.xlabel('x/L') - plt.ylabel('b') + plt.plot(x[:iend] / L, errb, "--", label="Error") + plt.xlabel("x/L") + plt.ylabel("b") plt.legend() plt.figure() plt.subplot(211) - plt.plot(x[:iend]/L,DSim[:iend],label='Sim') - plt.plot(x[:iend]/L,Dth[:iend],'--',label='Theoretical') - plt.xlabel('x/L') - plt.ylabel('D') + plt.plot(x[:iend] / L, DSim[:iend], label="Sim") + plt.plot(x[:iend] / L, Dth[:iend], "--", label="Theoretical") + plt.xlabel("x/L") + plt.ylabel("D") plt.subplot(212) - plt.plot(x[:iend]/L,errD,'--',label='Error') - plt.xlabel('x/L') - plt.ylabel('D') + plt.plot(x[:iend] / L, errD, "--", label="Error") + plt.xlabel("x/L") + plt.ylabel("D") plt.legend() plt.ioff() plt.show() -err=np.mean(0.5*np.sqrt(errb**2+errD**2)) -print("Error total=%e"%err) -if(err<5.9e-3): +err = np.mean(0.5 * np.sqrt(errb**2 + errD**2)) +print("Error total=%e" % err) +if err < 5.9e-3: print("Success!") sys.exit(0) else: diff --git a/test/MHD/AmbipolarCshock3D/testme.py b/test/MHD/AmbipolarCshock3D/testme.py index 981fbff3..37363d5b 100755 --- a/test/MHD/AmbipolarCshock3D/testme.py +++ b/test/MHD/AmbipolarCshock3D/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,43 +12,45 @@ import pytools.idfx_test as tst -tolerance=3e-14 +tolerance = 3e-14 + + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-rkl.ini"] - for ini in inifiles: - mytol=tolerance - - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename="dump.0001.dmp") - test.standardTest() - # When using RKL, except larger error due to B reconstruction and RKL # of substeps - if test.mpi and ini=="idefix-rkl.ini": - mytol=2e-10 - test.nonRegressionTest(filename="dump.0001.dmp",tolerance=mytol) - - -test=tst.idfxTest(__file__) + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-rkl.ini"] + for ini in inifiles: + mytol = tolerance + + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename="dump.0001.dmp") + test.standardTest() + # When using RKL, except larger error due to B reconstruction and RKL # of substeps + if test.mpi and ini == "idefix-rkl.ini": + mytol = 2e-10 + test.nonRegressionTest(filename="dump.0001.dmp", tolerance=mytol) + + +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','1','1'] + test.dec = ["2", "1", "1"] if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp",tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp", tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - test.vectPot=True - testMe(test) - - test.vectPot=False - test.mpi=True - testMe(test) + test.noplot = True + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + test.vectPot = True + testMe(test) + + test.vectPot = False + test.mpi = True + testMe(test) diff --git a/test/MHD/AmbipolarShearingBox/python/testidefix.py b/test/MHD/AmbipolarShearingBox/python/testidefix.py index 6aec86a6..779ecfae 100755 --- a/test/MHD/AmbipolarShearingBox/python/testidefix.py +++ b/test/MHD/AmbipolarShearingBox/python/testidefix.py @@ -5,6 +5,7 @@ @author: lesurg """ + import os import sys @@ -14,27 +15,27 @@ sys.path.append(os.getenv("IDEFIX_DIR")) from pytools.vtk_io import readVTK -rep='../' -nend=1000 +rep = "../" +nend = 1000 -n=0 -dt=0.1 +n = 0 +dt = 0.1 -Bx=np.zeros(nend) +Bx = np.zeros(nend) -t=dt*np.arange(0,nend) +t = dt * np.arange(0, nend) for n in range(nend): - V=readVTK(rep+'/data.'+'%0*d'%(4,n)+'.vtk', geometry='cartesian') - Bx[n]=np.sqrt(np.mean(np.mean(V.data['BX1']**2,axis=2),axis=0)) + V = readVTK(rep + "/data." + "%0*d" % (4, n) + ".vtk", geometry="cartesian") + Bx[n] = np.sqrt(np.mean(np.mean(V.data["BX1"] ** 2, axis=2), axis=0)) plt.figure() -plt.semilogy(t,Bx) -plt.semilogy(t,Bx[-1]*np.exp(0.171*(t-t[-1]))) +plt.semilogy(t, Bx) +plt.semilogy(t, Bx[-1] * np.exp(0.171 * (t - t[-1]))) -gr=np.log(Bx[-1]/Bx[-700])/(t[-1]-t[-700]) -print('growth rate=%g'%gr) -if np.abs(gr-0.171)/0.171<0.05: +gr = np.log(Bx[-1] / Bx[-700]) / (t[-1] - t[-700]) +print("growth rate=%g" % gr) +if np.abs(gr - 0.171) / 0.171 < 0.05: print("SUCCESS!") else: print("FAILED!") diff --git a/test/MHD/AxisFluxTube/python/checkAxisBounds.py b/test/MHD/AxisFluxTube/python/checkAxisBounds.py index 77406988..cb313a44 100755 --- a/test/MHD/AxisFluxTube/python/checkAxisBounds.py +++ b/test/MHD/AxisFluxTube/python/checkAxisBounds.py @@ -5,6 +5,7 @@ @author: lesurg """ + import os import sys @@ -13,59 +14,59 @@ from pytools.idfx_io import readIdfxFile -rep="../" +rep = "../" -idfx=readIdfxFile(rep+"analysis.99.0.idfx") -d=idfx.data -Bx1s=d['Vs'][0,:-1,:-1,:] -Bx2s=d['Vs'][1,:-1,:,:-1] -Bx3s=d['Vs'][2,:,:-1,:-1] +idfx = readIdfxFile(rep + "analysis.99.0.idfx") +d = idfx.data +Bx1s = d["Vs"][0, :-1, :-1, :] +Bx2s = d["Vs"][1, :-1, :, :-1] +Bx3s = d["Vs"][2, :, :-1, :-1] -nxhalf=10 -Bx1slice=Bx1s[:,:,nxhalf] -Bx2slice=Bx2s[:,:,nxhalf] -Bx3slice=Bx3s[:,:,nxhalf] +nxhalf = 10 +Bx1slice = Bx1s[:, :, nxhalf] +Bx2slice = Bx2s[:, :, nxhalf] +Bx3slice = Bx3s[:, :, nxhalf] -plt.matshow(Bx1slice,cmap='RdBu_r') +plt.matshow(Bx1slice, cmap="RdBu_r") plt.colorbar() plt.title("Bx1") -plt.matshow(Bx2slice,cmap='RdBu_r') +plt.matshow(Bx2slice, cmap="RdBu_r") plt.colorbar() plt.title("Bx2") -plt.matshow(Bx3slice,cmap='RdBu_r') +plt.matshow(Bx3slice, cmap="RdBu_r") plt.colorbar() plt.title("Bx3") # inner boundary -Bx1slice=Bx1s[:,:5,nxhalf] -Bx2slice=Bx2s[:,:5,nxhalf] -Bx3slice=Bx3s[:,:5,nxhalf] +Bx1slice = Bx1s[:, :5, nxhalf] +Bx2slice = Bx2s[:, :5, nxhalf] +Bx3slice = Bx3s[:, :5, nxhalf] -Vx1slice=d['Vc'][1,:,:5,nxhalf] -Vx2slice=d['Vc'][2,:,:5,nxhalf] -Vx3slice=d['Vc'][3,:,:5,nxhalf] +Vx1slice = d["Vc"][1, :, :5, nxhalf] +Vx2slice = d["Vc"][2, :, :5, nxhalf] +Vx3slice = d["Vc"][3, :, :5, nxhalf] -plt.matshow(Bx1slice,cmap='RdBu_r') +plt.matshow(Bx1slice, cmap="RdBu_r") plt.colorbar() plt.title("Bx1") -plt.matshow(Bx2slice,cmap='RdBu_r') +plt.matshow(Bx2slice, cmap="RdBu_r") plt.colorbar() plt.title("Bx2") -plt.matshow(Bx3slice,cmap='RdBu_r') +plt.matshow(Bx3slice, cmap="RdBu_r") plt.colorbar() plt.title("Bx3") -plt.matshow(Vx1slice,cmap='RdBu_r') +plt.matshow(Vx1slice, cmap="RdBu_r") plt.colorbar() plt.title("Vx1") -plt.matshow(Vx2slice,cmap='RdBu_r') +plt.matshow(Vx2slice, cmap="RdBu_r") plt.colorbar() plt.title("Vx2") -plt.matshow(Vx3slice,cmap='RdBu_r') +plt.matshow(Vx3slice, cmap="RdBu_r") plt.colorbar() plt.title("Vx3") diff --git a/test/MHD/AxisFluxTube/testme.py b/test/MHD/AxisFluxTube/testme.py index c5db6cfa..f027b24b 100755 --- a/test/MHD/AxisFluxTube/testme.py +++ b/test/MHD/AxisFluxTube/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,39 +12,40 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + +tolerance = 1e-14 -tolerance=1e-14 def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-coarsening.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-coarsening.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename=name) - test.nonRegressionTest(filename=name,tolerance=tolerance) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename=name) + test.nonRegressionTest(filename=name, tolerance=tolerance) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename=name,tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name, tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True + test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) - # test with MPI - test.mpi=True - testMe(test) + # test with MPI + test.mpi = True + testMe(test) diff --git a/test/MHD/Coarsening/python/testidefix.py b/test/MHD/Coarsening/python/testidefix.py index 206c0160..d901d805 100755 --- a/test/MHD/Coarsening/python/testidefix.py +++ b/test/MHD/Coarsening/python/testidefix.py @@ -19,31 +19,30 @@ from pytools.vtk_io import readVTK parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() # find which ini file was used -with open('../idefix.0.log','r') as fp: - lines = fp.readlines() - for line in lines: - if line.find(".ini") != -1: - s=line.split() - inputFile=s[5][:-1] - -input = inifix.load('../'+inputFile) +with open("../idefix.0.log", "r") as fp: + lines = fp.readlines() + for line in lines: + if line.find(".ini") != -1: + s = line.split() + inputFile = s[5][:-1] + +input = inifix.load("../" + inputFile) componentDir = input["Setup"]["direction"][0] spatialDir = input["Setup"]["direction"][1] # Convert componentDir -componentName=["BX1","BX2","BX3"] +componentName = ["BX1", "BX2", "BX3"] -V=readVTK('../data.0001.vtk', geometry='cartesian') +V = readVTK("../data.0001.vtk", geometry="cartesian") # left_state and right_state set p, rho and u # geometry sets left boundary on 0., right boundary on 1 and initial @@ -52,48 +51,46 @@ # gamma denotes specific heat # note that gamma and npts are default parameters (1.4 and 500) in solve function -x=V.x -if spatialDir==0: - qside=V.data[componentName[componentDir]][:,4,0] - qcenter=V.data[componentName[componentDir]][:,32,0] - x=V.x -if spatialDir==1: - qside=V.data[componentName[componentDir]][4,:,0] - qcenter=V.data[componentName[componentDir]][32,:,0] - x=V.y -if spatialDir==2: - qside=V.data[componentName[componentDir]][4,0,:] - qcenter=V.data[componentName[componentDir]][32,0,:] - x=V.z - -#initial condition in the code -q0=0.*x+0.1*(np.abs(x)<0.1) +x = V.x +if spatialDir == 0: + qside = V.data[componentName[componentDir]][:, 4, 0] + qcenter = V.data[componentName[componentDir]][:, 32, 0] + x = V.x +if spatialDir == 1: + qside = V.data[componentName[componentDir]][4, :, 0] + qcenter = V.data[componentName[componentDir]][32, :, 0] + x = V.y +if spatialDir == 2: + qside = V.data[componentName[componentDir]][4, 0, :] + qcenter = V.data[componentName[componentDir]][32, 0, :] + x = V.z + +# initial condition in the code +q0 = 0.0 * x + 0.1 * (np.abs(x) < 0.1) # compute solution to diffusion equation using Fourier transforms -qf0=np.fft.fft(q0) -k2=(2.0*np.pi*q0.size*np.fft.fftfreq(q0.size))**2 -qf=qf0*np.exp(-k2*V.t*0.05) -qtheory=np.real(np.fft.ifft(qf)) +qf0 = np.fft.fft(q0) +k2 = (2.0 * np.pi * q0.size * np.fft.fftfreq(q0.size)) ** 2 +qf = qf0 * np.exp(-k2 * V.t * 0.05) +qtheory = np.real(np.fft.ifft(qf)) -if(not args.noplot): +if not args.noplot: plt.figure(1) - plt.plot(x,qside,'+',markersize=6,label='side') - plt.plot(x,qcenter,'+',markersize=6,label='center') + plt.plot(x, qside, "+", markersize=6, label="side") + plt.plot(x, qcenter, "+", markersize=6, label="center") - #plt.plot(x,q0) - plt.plot(x,qtheory,label='theory') + # plt.plot(x,q0) + plt.plot(x, qtheory, label="theory") plt.legend() plt.title(componentName[componentDir]) - - plt.ioff() plt.show() -errorCenter=np.mean((qcenter-qtheory)**2) -errorSide=np.mean((qside-qtheory)**2) -print("Error Side=%e, Error Center=%e"%(errorSide,errorCenter)) -assert(errorCenter<2e-5) -assert(errorSide<2e-8) +errorCenter = np.mean((qcenter - qtheory) ** 2) +errorSide = np.mean((qside - qtheory) ** 2) +print("Error Side=%e, Error Center=%e" % (errorSide, errorCenter)) +assert errorCenter < 2e-5 +assert errorSide < 2e-8 print("SUCCESS!") diff --git a/test/MHD/Coarsening/testme.py b/test/MHD/Coarsening/testme.py index e82277ee..c1e63ac8 100755 --- a/test/MHD/Coarsening/testme.py +++ b/test/MHD/Coarsening/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,33 +12,34 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-rkl.ini","idefix-x2.ini","idefix-x3.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-rkl.ini", "idefix-x2.ini", "idefix-x3.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename=name) - test.standardTest() - test.nonRegressionTest(filename=name) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename=name) + test.standardTest() + test.nonRegressionTest(filename=name) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename=name) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name) + else: + testMe(test) else: - test.noplot = True - test.mpi=False - testMe(test) + test.noplot = True + test.mpi = False + testMe(test) - test.mpi=True - testMe(test) + test.mpi = True + testMe(test) diff --git a/test/MHD/FargoMHDSpherical/testme.py b/test/MHD/FargoMHDSpherical/testme.py index ee74c24d..5f003ed5 100755 --- a/test/MHD/FargoMHDSpherical/testme.py +++ b/test/MHD/FargoMHDSpherical/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,50 +12,51 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + +tolerance = 1e-14 -tolerance=1e-14 def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename=name) - test.nonRegressionTest(filename=name,tolerance=tolerance) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename=name) + test.nonRegressionTest(filename=name, tolerance=tolerance) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','2','2'] + test.dec = ["2", "2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename=name,tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name, tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - #test vector potential formulation - test.vectPot=True - testMe(test) - - # test with MPI - test.vectPot=False - test.mpi=True - testMe(test) - - test.vectPot=True - test.mpi=True - testMe(test) + test.noplot = True + + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + # test vector potential formulation + test.vectPot = True + testMe(test) + + # test with MPI + test.vectPot = False + test.mpi = True + testMe(test) + + test.vectPot = True + test.mpi = True + testMe(test) diff --git a/test/MHD/HallWhistler/python/testidefix.py b/test/MHD/HallWhistler/python/testidefix.py index 2ca5d6bf..2e08f387 100755 --- a/test/MHD/HallWhistler/python/testidefix.py +++ b/test/MHD/HallWhistler/python/testidefix.py @@ -5,6 +5,7 @@ @author: lesurg """ + #!/usr/bin/env python3 # -*- coding: utf-8 -*- """ @@ -20,69 +21,67 @@ from scipy.signal import argrelextrema parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() # values for the field strength to compute theoretical eigenfrequency values -lH=1.0 -B=1.0 -Va=1.0 -k=2.0*np.pi +lH = 1.0 +B = 1.0 +Va = 1.0 +k = 2.0 * np.pi # load the dat file produced by the setup -raw=np.loadtxt('../timevol.dat',skiprows=1) -t=raw[:,0] -by=raw[:,1] +raw = np.loadtxt("../timevol.dat", skiprows=1) +t = raw[:, 0] +by = raw[:, 1] # Theoretical speedes -f_w=Va/(2*np.pi)*k*(np.sqrt(1+(k*lH/2)**2)+k*lH/2) -f_i=Va/(2*np.pi)*k*(np.sqrt(1+(k*lH/2)**2)-k*lH/2) -f_A=k*Va/(2*np.pi) - +f_w = Va / (2 * np.pi) * k * (np.sqrt(1 + (k * lH / 2) ** 2) + k * lH / 2) +f_i = Va / (2 * np.pi) * k * (np.sqrt(1 + (k * lH / 2) ** 2) - k * lH / 2) +f_A = k * Va / (2 * np.pi) # Compute temporal spectrum and get frequency -#by=By[:,0] -sp=np.abs(np.fft.rfft(by)) -f=np.arange(sp.size)/(t[-1]-t[0]) - -r=10**np.arange(np.floor(np.log10(np.amin(sp))),np.ceil(np.log10(np.amax(sp))),0.1) - - - +# by=By[:,0] +sp = np.abs(np.fft.rfft(by)) +f = np.arange(sp.size) / (t[-1] - t[0]) -if(not args.noplot): +r = 10 ** np.arange( + np.floor(np.log10(np.amin(sp))), np.ceil(np.log10(np.amax(sp))), 0.1 +) - plt.close('all') +if not args.noplot: + plt.close("all") plt.figure() - plt.loglog(f,sp,label='signal') - plt.loglog(f_w+0*r,r,'--',label='whistler') - plt.loglog(f_i+0*r,r,'-.',label='ion cyclotron') - plt.loglog(f_A+0*r,r,':',label='Alfven') - plt.xlabel('frequency') - plt.ylabel('By amplitude') + plt.loglog(f, sp, label="signal") + plt.loglog(f_w + 0 * r, r, "--", label="whistler") + plt.loglog(f_i + 0 * r, r, "-.", label="ion cyclotron") + plt.loglog(f_A + 0 * r, r, ":", label="Alfven") + plt.xlabel("frequency") + plt.ylabel("By amplitude") plt.legend() plt.ioff() plt.show() -imax=argrelextrema(sp,np.greater,order=2)[0][0] -f_wnum=f[imax] -error=np.fabs(f_w-f_wnum)/f_w -print("Theoretical whistler frequency=%g, numerical=%g, error=%g"%(f_w,f_wnum,error)) +imax = argrelextrema(sp, np.greater, order=2)[0][0] +f_wnum = f[imax] +error = np.fabs(f_w - f_wnum) / f_w +print( + "Theoretical whistler frequency=%g, numerical=%g, error=%g" % (f_w, f_wnum, error) +) -if(error<0.06): +if error < 0.06: print("SUCCESS") sys.exit(0) else: diff --git a/test/MHD/HallWhistler/testme.py b/test/MHD/HallWhistler/testme.py index ce50a7ae..2e0b0abc 100755 --- a/test/MHD/HallWhistler/testme.py +++ b/test/MHD/HallWhistler/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,34 +12,36 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" -tolerance=1e-15 +name = "dump.0001.dmp" +tolerance = 1e-15 + + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - test.standardTest(); - if test.init: - test.makeReference(filename=name) - test.nonRegressionTest(filename=name,tolerance=tolerance) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + test.standardTest() + if test.init: + test.makeReference(filename=name) + test.nonRegressionTest(filename=name, tolerance=tolerance) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename=name,tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name, tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True + test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) diff --git a/test/MHD/LinearWaveTest/python/testidefix.py b/test/MHD/LinearWaveTest/python/testidefix.py index 8603ee8e..f1555991 100755 --- a/test/MHD/LinearWaveTest/python/testidefix.py +++ b/test/MHD/LinearWaveTest/python/testidefix.py @@ -15,54 +15,62 @@ from pytools.dump_io import readDump parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() def getError(rep): try: - V=readDump(rep+'/dump.0000.dmp') - U=readDump(rep+'/dump.0001.dmp') + V = readDump(rep + "/dump.0000.dmp") + U = readDump(rep + "/dump.0001.dmp") except: return math.nan - keylist=[] - keylist=['Vc-RHO','Vc-VX1','Vc-VX2','Vc-VX3','Vs-BX1s','Vs-BX2s','Vs-BX3s','Vc-PRS'] + keylist = [] + keylist = [ + "Vc-RHO", + "Vc-VX1", + "Vc-VX2", + "Vc-VX3", + "Vs-BX1s", + "Vs-BX2s", + "Vs-BX3s", + "Vc-PRS", + ] err = 0.0 for key in keylist: - Q1=V.data[key] - Q2=U.data[key] - - Q1=Q1-np.mean(Q1) - Q2=Q2-np.mean(Q2) - print(key+" error=%e"%((np.mean(np.abs(Q1-Q2)))*1e6)) - err=err + (np.mean(np.abs(Q1-Q2)))**2 - if(not args.noplot): - plt.figure() - plt.contourf(V.data[key][:,:,16]-U.data[key][:,:,16]) - plt.title(key) - - - err = err/len(keylist) + Q1 = V.data[key] + Q2 = U.data[key] + + Q1 = Q1 - np.mean(Q1) + Q2 = Q2 - np.mean(Q2) + print(key + " error=%e" % ((np.mean(np.abs(Q1 - Q2))) * 1e6)) + err = err + (np.mean(np.abs(Q1 - Q2))) ** 2 + if not args.noplot: + plt.figure() + plt.contourf(V.data[key][:, :, 16] - U.data[key][:, :, 16]) + plt.title(key) + + err = err / len(keylist) err = np.sqrt(err) return err + error = getError("..") -print("Error=%e"%(error*1e6)) +print("Error=%e" % (error * 1e6)) -if(not args.noplot): - plt.show() +if not args.noplot: + plt.show() -if error*1e6<3e-2: - print("SUCCESS") - sys.exit(0) +if error * 1e6 < 3e-2: + print("SUCCESS") + sys.exit(0) else: - print("Failed") - sys.exit(1) + print("Failed") + sys.exit(1) diff --git a/test/MHD/LinearWaveTest/testme.py b/test/MHD/LinearWaveTest/testme.py index 5b4aed87..4e6ae749 100755 --- a/test/MHD/LinearWaveTest/testme.py +++ b/test/MHD/LinearWaveTest/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,43 +12,50 @@ import pytools.idfx_test as tst -tolerance=2e-13 -def testMe(test): - test.configure() - test.compile() - inifiles=["idefix-fast.ini","idefix-slow.ini","idefix-alfven.ini","idefix-entropy.ini"] - for ini in inifiles: - mytol=tolerance - - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename="dump.0001.dmp") - test.standardTest() - test.nonRegressionTest(filename="dump.0001.dmp",tolerance=mytol) +tolerance = 2e-13 -test=tst.idfxTest(__file__) +def testMe(test): + test.configure() + test.compile() + inifiles = [ + "idefix-fast.ini", + "idefix-slow.ini", + "idefix-alfven.ini", + "idefix-entropy.ini", + ] + for ini in inifiles: + mytol = tolerance + + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename="dump.0001.dmp") + test.standardTest() + test.nonRegressionTest(filename="dump.0001.dmp", tolerance=mytol) + + +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','2','2'] + test.dec = ["2", "2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp",tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp", tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) + test.noplot = True + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) - test.mpi=True - testMe(test) + test.mpi = True + testMe(test) - test.mpi=False - test.reconstruction=3 - testMe(test) + test.mpi = False + test.reconstruction = 3 + testMe(test) - test.reconstruction=4 - testMe(test) + test.reconstruction = 4 + testMe(test) diff --git a/test/MHD/MTI/python/testidefix.py b/test/MHD/MTI/python/testidefix.py index 18481692..28c72161 100644 --- a/test/MHD/MTI/python/testidefix.py +++ b/test/MHD/MTI/python/testidefix.py @@ -14,60 +14,71 @@ import inifix import numpy as np -n = 1 # fondamental mode +n = 1 # fondamental mode conf = inifix.load("../idefix.ini") -Pr=conf["Setup"]["pr"] -ksi=conf["Setup"]["ksi"] +Pr = conf["Setup"]["pr"] +ksi = conf["Setup"]["ksi"] L = 0.1 -H = 3. -wbuoy = 1/np.sqrt(3) -gamma = 5/3 +H = 3.0 +wbuoy = 1 / np.sqrt(3) +gamma = 5 / 3 -kn = 2*np.pi*n*np.sqrt(2)/L +kn = 2 * np.pi * n * np.sqrt(2) / L -fid=open('../average.dat',"r") +fid = open("../average.dat", "r") # read the first line to get data names -varnames=fid.readline().split() +varnames = fid.readline().split() fid.close() # load the bulk of the file -data=np.loadtxt('../average.dat',skiprows=1) +data = np.loadtxt("../average.dat", skiprows=1) # store this in our data structure -V={} -i=0 +V = {} +i = 0 for name in varnames: - V[name]=data[:,i] - i=i+1 + V[name] = data[:, i] + i = i + 1 + +sigma_vx = ( + 0.5 * (np.log(V["kinx"][20]) - np.log(V["kinx"][15])) / (V["t"][20] - V["t"][15]) +) +sigma_vy = ( + 0.5 * (np.log(V["kiny"][20]) - np.log(V["kiny"][15])) / (V["t"][20] - V["t"][15]) +) +sigma_by = ( + 0.5 * (np.log(V["by_2"][20]) - np.log(V["by_2"][15])) / (V["t"][20] - V["t"][15]) +) +sigma_num = (sigma_vx + sigma_vy + sigma_by) / 3.0 -sigma_vx = 0.5*(np.log(V['kinx'][20]) - np.log(V['kinx'][15]))/(V['t'][20] - V['t'][15]) -sigma_vy = 0.5*(np.log(V['kiny'][20]) - np.log(V['kiny'][15]))/(V['t'][20] - V['t'][15]) -sigma_by = 0.5*(np.log(V['by_2'][20]) - np.log(V['by_2'][15]))/(V['t'][20] - V['t'][15]) -sigma_num = (sigma_vx + sigma_vy + sigma_by)/3. def dispMTI(wcond, k, Pr): - cosk_2 = 0.5 # here cosk_2 = sqrt(bhat) dot khat i.e. xhat dot khat + cosk_2 = 0.5 # here cosk_2 = sqrt(bhat) dot khat i.e. xhat dot khat k_2 = k**2 - ksi = 2.5/k_2/cosk_2*wcond - kx_2 = k_2*cosk_2 - K = kx_2/k_2 + ksi = 2.5 / k_2 / cosk_2 * wcond + kx_2 = k_2 * cosk_2 + K = kx_2 / k_2 - bk_2 = k_2*cosk_2 + bk_2 = k_2 * cosk_2 bkhat_2 = cosk_2 - nu = Pr*ksi - wvisc = 3*nu*bk_2 - v = wvisc*(1 - bkhat_2) - N_2 = wbuoy**2/gamma*((H-1)*gamma-H) - return np.array([1, wcond+v, wcond*v+N_2*kx_2/k_2, K*wcond*-wbuoy**2]) + nu = Pr * ksi + wvisc = 3 * nu * bk_2 + v = wvisc * (1 - bkhat_2) + N_2 = wbuoy**2 / gamma * ((H - 1) * gamma - H) + return np.array( + [1, wcond + v, wcond * v + N_2 * kx_2 / k_2, K * wcond * -(wbuoy**2)] + ) + def sigmaMTI(wcond, k, Pr): - roots = np.roots(dispMTI(wcond,k,Pr)) + roots = np.roots(dispMTI(wcond, k, Pr)) return np.real(roots[np.real(roots) > 0])[0] -wcd = 0.4*ksi*kn**2/2 + +wcd = 0.4 * ksi * kn**2 / 2 sigma_analytic = sigmaMTI(wcd, kn, Pr) -if np.abs(np.log(sigma_analytic/sigma_num)/np.log(sigma_num)) >= 0.1: +if np.abs(np.log(sigma_analytic / sigma_num) / np.log(sigma_num)) >= 0.1: print("Failed") sys.exit(1) else: diff --git a/test/MHD/MTI/testme.py b/test/MHD/MTI/testme.py index baaed630..cb145f3a 100755 --- a/test/MHD/MTI/testme.py +++ b/test/MHD/MTI/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -13,32 +14,32 @@ def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-rkl.ini","idefix-sl.ini"] -# inifiles=["idefix-rkl.ini","idefix-sl.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-rkl.ini", "idefix-sl.ini"] + # inifiles=["idefix-rkl.ini","idefix-sl.ini"] - # loop on all the ini files for this test - name="dump.0001.dmp" - for ini in inifiles: - test.run(inputFile=ini) - test.standardTest() - if test.init: - test.makeReference(filename=name) - test.nonRegressionTest(filename=name) + # loop on all the ini files for this test + name = "dump.0001.dmp" + for ini in inifiles: + test.run(inputFile=ini) + test.standardTest() + if test.init: + test.makeReference(filename=name) + test.nonRegressionTest(filename=name) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp") - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp") + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) diff --git a/test/MHD/OrszagTang/testme.py b/test/MHD/OrszagTang/testme.py index 7f9688eb..3c0e56bc 100755 --- a/test/MHD/OrszagTang/testme.py +++ b/test/MHD/OrszagTang/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,58 +12,64 @@ import pytools.idfx_test as tst -tolerance=1e-12 +tolerance = 1e-12 + + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini", - "idefix-hll.ini","idefix-hlld-arithmetic.ini", - "idefix-hlld-hll.ini", - "idefix-hlld-hlld.ini", - "idefix-hlld-uct0.ini", - "idefix-hlld.ini","idefix-tvdlf.ini"] + test.configure() + test.compile() + inifiles = [ + "idefix.ini", + "idefix-hll.ini", + "idefix-hlld-arithmetic.ini", + "idefix-hlld-hll.ini", + "idefix-hlld-hlld.ini", + "idefix-hlld-uct0.ini", + "idefix-hlld.ini", + "idefix-tvdlf.ini", + ] - for ini in inifiles: - mytol=tolerance + for ini in inifiles: + mytol = tolerance - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename="dump.0001.dmp") + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename="dump.0001.dmp") - if(test.single): - mytol=1e-5 + if test.single: + mytol = 1e-5 - test.nonRegressionTest(filename="dump.0001.dmp",tolerance=mytol) + test.nonRegressionTest(filename="dump.0001.dmp", tolerance=mytol) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','2'] + test.dec = ["2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp",tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp", tolerance=tolerance) + else: + testMe(test) else: - for rec in [2,3,4]: - test.noplot = True - test.vectPot = False - test.single=False - test.reconstruction=rec - test.mpi=False - testMe(test) - test.mpi=True + for rec in [2, 3, 4]: + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = rec + test.mpi = False + testMe(test) + test.mpi = True + testMe(test) + + # single precision validation + test.reconstruction = 2 + test.mpi = False + test.single = True testMe(test) - # single precision validation - test.reconstruction=2 - test.mpi=False - test.single=True - testMe(test) - - # Vector potential validation - test.single=False - test.mpi=False - test.vectPot=True - testMe(test) + # Vector potential validation + test.single = False + test.mpi = False + test.vectPot = True + testMe(test) diff --git a/test/MHD/OrszagTang3D/testme.py b/test/MHD/OrszagTang3D/testme.py index 662ba406..885da1de 100755 --- a/test/MHD/OrszagTang3D/testme.py +++ b/test/MHD/OrszagTang3D/testme.py @@ -3,6 +3,7 @@ @author: glesur """ + import os import sys @@ -12,59 +13,59 @@ # Whether we should reset our reference run (only do that on purpose!) -tolerance=1e-13 +tolerance = 1e-13 + def testMe(test): - test.configure() - test.compile() - tol=tolerance - if test.single: - tol=1e-6 + test.configure() + test.compile() + tol = tolerance + if test.single: + tol = 1e-6 - # default with idefix.ini - test.run() - if test.init: - if not test.mpi: - test.makeReference(filename="dump.0001.dmp") - test.nonRegressionTest(filename="dump.0001.dmp",tolerance=tol) + # default with idefix.ini + test.run() + if test.init: + if not test.mpi: + test.makeReference(filename="dump.0001.dmp") + test.nonRegressionTest(filename="dump.0001.dmp", tolerance=tol) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) # if no decomposition is specified, use that one if not test.dec: - test.dec=["2","2","2"] + test.dec = ["2", "2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp",tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp", tolerance=tolerance) + else: + testMe(test) else: - test.vectPot=False - test.reconstruction=2 - test.mpi=False - testMe(test) - # test in MPI mode - test.mpi=True - testMe(test) - - - # test with vector potential - test.mpi=False - test.vectPot=True - test.reconstruction=2 - testMe(test) - - test.mpi=True - testMe(test) - - # test with other precision - test.single=True - test.vectPot=False - test.reconstruction=2 - test.mpi=False - testMe(test) - # test in MPI mode - test.mpi=True - testMe(test) + test.vectPot = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + # test in MPI mode + test.mpi = True + testMe(test) + + # test with vector potential + test.mpi = False + test.vectPot = True + test.reconstruction = 2 + testMe(test) + + test.mpi = True + testMe(test) + + # test with other precision + test.single = True + test.vectPot = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + # test in MPI mode + test.mpi = True + testMe(test) diff --git a/test/MHD/ResistiveAlfvenWave/python/testidefix.py b/test/MHD/ResistiveAlfvenWave/python/testidefix.py index 52b11c8c..c221544d 100644 --- a/test/MHD/ResistiveAlfvenWave/python/testidefix.py +++ b/test/MHD/ResistiveAlfvenWave/python/testidefix.py @@ -5,41 +5,39 @@ import numpy as np parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() # values for the field strength to compute theoretical decay rate -eta=0.05 -k=2.0*np.pi -B0=1.0 +eta = 0.05 +k = 2.0 * np.pi +B0 = 1.0 # decay rate of the energy, assuming k^2eta^2 < 2B0^2 -gamma=-k**2*eta/2 +gamma = -(k**2) * eta / 2 -raw=np.loadtxt('../timevol.dat',skiprows=1) -t=raw[:,0] -e=raw[:,1] +raw = np.loadtxt("../timevol.dat", skiprows=1) +t = raw[:, 0] +e = raw[:, 1] -eth=e[0]*np.exp(2*gamma*t) -error=np.mean(e/eth-1.0) +eth = e[0] * np.exp(2 * gamma * t) +error = np.mean(e / eth - 1.0) -if(not args.noplot): - - plt.close('all') +if not args.noplot: + plt.close("all") plt.figure() - plt.semilogy(t,e) - plt.semilogy(t,eth) + plt.semilogy(t, e) + plt.semilogy(t, eth) plt.ioff() plt.show() -if(error<0.03): +if error < 0.03: print("SUCCESS") sys.exit(0) else: diff --git a/test/MHD/ResistiveAlfvenWave/testme.py b/test/MHD/ResistiveAlfvenWave/testme.py index 72014c7c..b7ce4a73 100755 --- a/test/MHD/ResistiveAlfvenWave/testme.py +++ b/test/MHD/ResistiveAlfvenWave/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,35 +12,36 @@ import pytools.idfx_test as tst -tolerance=1e-14 +tolerance = 1e-14 + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-rkl.ini"] - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename="dump.0001.dmp") - test.standardTest() - mytol=tolerance - if ini=="idefix-rkl.ini": - mytol=1e-10 - test.nonRegressionTest(filename="dump.0001.dmp",tolerance=mytol) - - -test=tst.idfxTest(__file__) + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-rkl.ini"] + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename="dump.0001.dmp") + test.standardTest() + mytol = tolerance + if ini == "idefix-rkl.ini": + mytol = 1e-10 + test.nonRegressionTest(filename="dump.0001.dmp", tolerance=mytol) + + +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp",tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp", tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - test.mpi=True - testMe(test) + test.noplot = True + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + test.mpi = True + testMe(test) diff --git a/test/MHD/ShearingBox/python/testidefix.py b/test/MHD/ShearingBox/python/testidefix.py index 1ab5b293..6f75c6af 100755 --- a/test/MHD/ShearingBox/python/testidefix.py +++ b/test/MHD/ShearingBox/python/testidefix.py @@ -5,6 +5,7 @@ @author: lesurg """ + import argparse import sys @@ -13,13 +14,13 @@ from scipy.integrate import solve_ivp -#compute theoretical solution (from Balbus & Hawley 2006, using notations from Lesur 2021) +# compute theoretical solution (from Balbus & Hawley 2006, using notations from Lesur 2021) def rhs(t, y, Omega, q, B0y, B0z, k0x, k0y, k0z): - kx = k0x + q*Omega * t * k0y + kx = k0x + q * Omega * t * k0y ky = k0y kz = k0z - k2 = kx*kx+ky*ky+kz*kz + k2 = kx * kx + ky * ky + kz * kz vx = y[0] vy = y[1] vz = y[2] @@ -27,100 +28,114 @@ def rhs(t, y, Omega, q, B0y, B0z, k0x, k0y, k0z): by = y[4] bz = y[5] - gxx = kx*kx/k2 - gxy = kx*ky/k2 - gyy = ky*ky/k2 - gyz = ky*kz/k2 - gxz = kx*kz/k2 - #gzz = kz*kz/k2 + gxx = kx * kx / k2 + gxy = kx * ky / k2 + gyy = ky * ky / k2 + gyz = ky * kz / k2 + gxz = kx * kz / k2 + # gzz = kz*kz/k2 + + kdotB = ky * B0y + kz * B0z - kdotB = ky*B0y+kz*B0z + # we assume b=1j*db hence db=-1j*b - #we assume b=1j*db hence db=-1j*b + dvx = kdotB * bx + 2 * Omega * vy * (1 - gxx) + 2 * (1 - q) * Omega * gxy * vx + dvy = ( + kdotB * by + - q * Omega * vx * gyy + - (2 - q) * Omega * (1 - gyy) * vx + - 2 * Omega * vy * gxy + ) + dvz = kdotB * bz + 2 * (1 - q) * Omega * vx * gyz - 2 * Omega * vy * gxz - dvx = kdotB*bx + 2*Omega*vy*(1-gxx) + 2*(1-q)*Omega*gxy*vx - dvy = kdotB*by - q*Omega*vx*gyy - (2-q)*Omega*(1-gyy)*vx - 2*Omega*vy*gxy - dvz= kdotB*bz +2*(1-q)*Omega*vx*gyz-2*Omega*vy*gxz + dbx = -kdotB * vx + dby = -kdotB * vy - q * Omega * bx + dbz = -kdotB * vz - dbx = - kdotB*vx - dby = - kdotB*vy - q*Omega*bx - dbz = - kdotB*vz + return np.asarray([dvx, dvy, dvz, dbx, dby, dbz]) - return np.asarray([dvx,dvy,dvz,dbx,dby,dbz]) parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() -#initial condition: vr=1, rest is 0, mode initial is nx=0, ny=1, nz=5) -y=solve_ivp(rhs, [0,30], [1, 0, 0, 0, 0, 0], args=(1, 1.5, 0.02, 0.05, 0, 2.0*np.pi, 4*np.pi), dense_output=True) +# initial condition: vr=1, rest is 0, mode initial is nx=0, ny=1, nz=5) +y = solve_ivp( + rhs, + [0, 30], + [1, 0, 0, 0, 0, 0], + args=(1, 1.5, 0.02, 0.05, 0, 2.0 * np.pi, 4 * np.pi), + dense_output=True, +) # read timevol file rep = "../" filename = "timevol.dat" -fid=open(rep+filename,"r") +fid = open(rep + filename, "r") -gamma=5/3 +gamma = 5 / 3 # read the first line to get data names -varnames=fid.readline().split() +varnames = fid.readline().split() fid.close() # load the bulk of the file -data=np.loadtxt(rep+filename,skiprows=1) +data = np.loadtxt(rep + filename, skiprows=1) # store this in our data structure -V={} +V = {} -i=0 +i = 0 for name in varnames: - V[name]=data[:,i] - i=i+1 + V[name] = data[:, i] + i = i + 1 # velocity normalisation -v0=V['vx'][0] +v0 = V["vx"][0] # compute L2 error norm -error=np.sqrt((V['vx']/v0-y.sol(V["t"])[0,:])**2 + (V['vy']/v0-y.sol(V["t"])[1,:])**2+ + (V['vz']/v0-y.sol(V["t"])[2,:])**2) - - - -if(not args.noplot): - # Comparison with exact solution - plt.rc('text', usetex=True) - plt.rc('font', family='serif') - plt.rc('font', size=16) - plt.close('all') - plt.figure(1) - plt.plot(V["t"],V['vx']/v0,'r-',label=r'$u_{R}$') - plt.plot(V["t"],y.sol(V["t"])[0,:],'r--') - plt.plot(V["t"],V['vy']/v0,'b-',label=r'$u_{\varphi}$') - plt.plot(V["t"],y.sol(V["t"])[1,:],'b--') - plt.plot(V["t"],V['vz']/v0,'g-',label=r'$u_{z}$') - plt.plot(V["t"],y.sol(V["t"])[2,:],'g--') - plt.legend() - plt.xlabel('t') - - # plot error - plt.figure() - plt.semilogy(V["t"],error) - plt.xlim([0,10]) - plt.ylim([1e-3,1]) - - plt.ioff() - plt.show() - -err=np.mean(error) -print("Error=",err) - -if(err<0.03): - print("SUCCESS") - sys.exit(0) +error = np.sqrt( + (V["vx"] / v0 - y.sol(V["t"])[0, :]) ** 2 + + (V["vy"] / v0 - y.sol(V["t"])[1, :]) ** 2 + + (V["vz"] / v0 - y.sol(V["t"])[2, :]) ** 2 +) + + +if not args.noplot: + # Comparison with exact solution + plt.rc("text", usetex=True) + plt.rc("font", family="serif") + plt.rc("font", size=16) + plt.close("all") + plt.figure(1) + plt.plot(V["t"], V["vx"] / v0, "r-", label=r"$u_{R}$") + plt.plot(V["t"], y.sol(V["t"])[0, :], "r--") + plt.plot(V["t"], V["vy"] / v0, "b-", label=r"$u_{\varphi}$") + plt.plot(V["t"], y.sol(V["t"])[1, :], "b--") + plt.plot(V["t"], V["vz"] / v0, "g-", label=r"$u_{z}$") + plt.plot(V["t"], y.sol(V["t"])[2, :], "g--") + plt.legend() + plt.xlabel("t") + + # plot error + plt.figure() + plt.semilogy(V["t"], error) + plt.xlim([0, 10]) + plt.ylim([1e-3, 1]) + + plt.ioff() + plt.show() + +err = np.mean(error) +print("Error=", err) + +if err < 0.03: + print("SUCCESS") + sys.exit(0) else: - print("FAILED") + print("FAILED") diff --git a/test/MHD/ShearingBox/testme.py b/test/MHD/ShearingBox/testme.py index 1e30e174..0d782c1c 100755 --- a/test/MHD/ShearingBox/testme.py +++ b/test/MHD/ShearingBox/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,38 +12,39 @@ import pytools.idfx_test as tst -tolerance=1e-14 +tolerance = 1e-14 + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-fargo.ini"] - for ini in inifiles: - mytol=tolerance + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-fargo.ini"] + for ini in inifiles: + mytol = tolerance - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename="dump.0001.dmp") - test.standardTest() - # When using RKL, except larger error due to B reconstruction and RKL # of substeps - test.nonRegressionTest(filename="dump.0001.dmp",tolerance=mytol) + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename="dump.0001.dmp") + test.standardTest() + # When using RKL, except larger error due to B reconstruction and RKL # of substeps + test.nonRegressionTest(filename="dump.0001.dmp", tolerance=mytol) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','1','2'] + test.dec = ["2", "1", "2"] if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp",tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp", tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - test.mpi=True - testMe(test) + test.noplot = True + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + test.mpi = True + testMe(test) diff --git a/test/MHD/clessTDiffusion/python/testidefix.py b/test/MHD/clessTDiffusion/python/testidefix.py index f37cf503..b68a52ba 100644 --- a/test/MHD/clessTDiffusion/python/testidefix.py +++ b/test/MHD/clessTDiffusion/python/testidefix.py @@ -11,26 +11,28 @@ gamma = conf["Hydro"]["gamma"] kappa = conf["Hydro"]["bragTDiffusion"][-1] -current_VTK = readVTK('../data.0003.vtk', geometry='spherical') +current_VTK = readVTK("../data.0003.vtk", geometry="spherical") -rho = current_VTK.data['RHO'].squeeze() -prs = current_VTK.data['PRS'].squeeze() -r = current_VTK.r +rho = current_VTK.data["RHO"].squeeze() +prs = current_VTK.data["PRS"].squeeze() +r = current_VTK.r + +idx1 = np.where(r > 20)[0][0] +idx2 = np.where(r > 30)[0][0] -idx1=np.where(r > 20)[0][0] -idx2=np.where(r > 30)[0][0] def gamma_prime(gamma, beta): - return (gamma+beta*(gamma-1))/(1+beta*(gamma-1)) + return (gamma + beta * (gamma - 1)) / (1 + beta * (gamma - 1)) + success = True eps = 3e-4 -gamma_eff = np.gradient(prs, r)/np.gradient(rho, r)*rho/prs +gamma_eff = np.gradient(prs, r) / np.gradient(rho, r) * rho / prs -Error=np.abs(gamma_eff[idx1:idx2]-gamma_prime(gamma, 1.5)).mean() +Error = np.abs(gamma_eff[idx1:idx2] - gamma_prime(gamma, 1.5)).mean() if Error > eps: - success=False + success = False if success: print("SUCCESS") diff --git a/test/MHD/clessTDiffusion/testme.py b/test/MHD/clessTDiffusion/testme.py index 98e3def9..0d1ca73c 100755 --- a/test/MHD/clessTDiffusion/testme.py +++ b/test/MHD/clessTDiffusion/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -13,31 +14,31 @@ def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - name = "dump.0001.dmp" - for ini in inifiles: - test.run(inputFile=ini) - test.standardTest() - if test.init: - test.makeReference(filename=name) - test.nonRegressionTest(filename=name, tolerance=2e-15) + # loop on all the ini files for this test + name = "dump.0001.dmp" + for ini in inifiles: + test.run(inputFile=ini) + test.standardTest() + if test.init: + test.makeReference(filename=name) + test.nonRegressionTest(filename=name, tolerance=2e-15) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp") - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp") + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) diff --git a/test/MHD/sod-iso/testme.py b/test/MHD/sod-iso/testme.py index 70c7daaa..f559eeb1 100755 --- a/test/MHD/sod-iso/testme.py +++ b/test/MHD/sod-iso/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,40 +12,41 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-hll.ini","idefix-hlld.ini","idefix-tvdlf.ini"] - if test.reconstruction==4: - inifiles=["idefix-rk3.ini","idefix-hlld-rk3.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-hll.ini", "idefix-hlld.ini", "idefix-tvdlf.ini"] + if test.reconstruction == 4: + inifiles = ["idefix-rk3.ini", "idefix-hlld-rk3.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init: - test.makeReference(filename=name) - test.nonRegressionTest(filename=name) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init: + test.makeReference(filename=name) + test.nonRegressionTest(filename=name) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename=name) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name) + else: + testMe(test) else: - test.noplot = True - for rec in range(2,5): - test.vectPot=False - test.single=False - test.reconstruction=rec - test.mpi=False + test.noplot = True + for rec in range(2, 5): + test.vectPot = False + test.single = False + test.reconstruction = rec + test.mpi = False + testMe(test) + + # test in single precision + test.reconstruction = 2 + test.single = True testMe(test) - - # test in single precision - test.reconstruction=2 - test.single=True - testMe(test) diff --git a/test/MHD/sod/testme.py b/test/MHD/sod/testme.py index 70c7daaa..f559eeb1 100755 --- a/test/MHD/sod/testme.py +++ b/test/MHD/sod/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,40 +12,41 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini","idefix-hll.ini","idefix-hlld.ini","idefix-tvdlf.ini"] - if test.reconstruction==4: - inifiles=["idefix-rk3.ini","idefix-hlld-rk3.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini", "idefix-hll.ini", "idefix-hlld.ini", "idefix-tvdlf.ini"] + if test.reconstruction == 4: + inifiles = ["idefix-rk3.ini", "idefix-hlld-rk3.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init: - test.makeReference(filename=name) - test.nonRegressionTest(filename=name) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init: + test.makeReference(filename=name) + test.nonRegressionTest(filename=name) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename=name) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name) + else: + testMe(test) else: - test.noplot = True - for rec in range(2,5): - test.vectPot=False - test.single=False - test.reconstruction=rec - test.mpi=False + test.noplot = True + for rec in range(2, 5): + test.vectPot = False + test.single = False + test.reconstruction = rec + test.mpi = False + testMe(test) + + # test in single precision + test.reconstruction = 2 + test.single = True testMe(test) - - # test in single precision - test.reconstruction=2 - test.single=True - testMe(test) diff --git a/test/MHD/sphBragTDiffusion/python/testidefix.py b/test/MHD/sphBragTDiffusion/python/testidefix.py index d42c8bc9..f96cbe6a 100644 --- a/test/MHD/sphBragTDiffusion/python/testidefix.py +++ b/test/MHD/sphBragTDiffusion/python/testidefix.py @@ -14,14 +14,14 @@ time_step = conf["Output"]["vtk"] nfile = 6 -list_time = [round(i*time_step, 3) for i in range(nfile)] +list_time = [round(i * time_step, 3) for i in range(nfile)] list_VTK = list() for k in range(nfile): if k >= 10: current_number = str(k) else: - current_number = '0' + str(k) - current_VTK = readVTK('../data.00' + current_number + '.vtk', geometry='spherical') + current_number = "0" + str(k) + current_VTK = readVTK("../data.00" + current_number + ".vtk", geometry="spherical") list_VTK.append(current_VTK) list_VX = list() @@ -30,55 +30,65 @@ list_K = list() list_T = list() for VTK in list_VTK: - list_VX.append(VTK.data['VX1']) - list_RHO.append(VTK.data['RHO']) - list_PRS.append(VTK.data['PRS']) - list_K.append(VTK.data['PRS']/(VTK.data['RHO']**gamma)) - list_T.append(VTK.data['PRS']/VTK.data['RHO']) + list_VX.append(VTK.data["VX1"]) + list_RHO.append(VTK.data["RHO"]) + list_PRS.append(VTK.data["PRS"]) + list_K.append(VTK.data["PRS"] / (VTK.data["RHO"] ** gamma)) + list_T.append(VTK.data["PRS"] / VTK.data["RHO"]) VTK = list_VTK[-1] R = VTK.r TH = VTK.theta PHI = VTK.phi -ir = len(R)//5 -ith = len(TH)//14 -iphi = len(PHI)*5//7 +ir = len(R) // 5 +ith = len(TH) // 14 +iphi = len(PHI) * 5 // 7 + +rho0 = 1.0 +c = 1.0 / (gamma - 1.0) +D = kappa / (rho0 * c) -rho0 = 1. -c = 1./(gamma - 1.) -D = kappa/(rho0*c) def analytic_sol(r, th, phi, t): - lambd2 = 1./r**2/np.sin(th)**2 - return amplitude*((3/r**2 - 1)*np.sin(r)/r - 3*np.cos(r)/r**2)*(-3*np.cos(th)*np.sin(th))*(np.sin(phi))*np.exp(-t*D*lambd2) + 1 + lambd2 = 1.0 / r**2 / np.sin(th) ** 2 + return ( + amplitude + * ((3 / r**2 - 1) * np.sin(r) / r - 3 * np.cos(r) / r**2) + * (-3 * np.cos(th) * np.sin(th)) + * (np.sin(phi)) + * np.exp(-t * D * lambd2) + + 1 + ) + def analytic_rad_sol(r, t): th = TH[ith] phi = PHI[iphi] - return analytic_sol(r,th,phi,t) + return analytic_sol(r, th, phi, t) + def analytic_th_sol(th, t): r = R[ir] phi = PHI[iphi] - return analytic_sol(r,th,phi,t) + return analytic_sol(r, th, phi, t) def analytic_phi_sol(phi, t): r = R[ir] th = TH[ith] - return analytic_sol(r,th,phi,t) + return analytic_sol(r, th, phi, t) success = True eps = 6.4e-8 -for it,t in enumerate(list_time): +for it, t in enumerate(list_time): analytic_T = np.zeros((R.shape[0], TH.shape[0], PHI.shape[0])) - for ir,r in enumerate(R): - for ith,th in enumerate(TH): - for iphi,phi in enumerate(PHI): - analytic_T[ir,ith,iphi] = analytic_sol(r,th,phi,t) + for ir, r in enumerate(R): + for ith, th in enumerate(TH): + for iphi, phi in enumerate(PHI): + analytic_T[ir, ith, iphi] = analytic_sol(r, th, phi, t) TEMP = list_T[it] if np.mean(np.fabs(TEMP - analytic_T)) > eps: - success = False + success = False if success: print("SUCCESS") diff --git a/test/MHD/sphBragTDiffusion/testme.py b/test/MHD/sphBragTDiffusion/testme.py index 98e3def9..0d1ca73c 100755 --- a/test/MHD/sphBragTDiffusion/testme.py +++ b/test/MHD/sphBragTDiffusion/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -13,31 +14,31 @@ def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - name = "dump.0001.dmp" - for ini in inifiles: - test.run(inputFile=ini) - test.standardTest() - if test.init: - test.makeReference(filename=name) - test.nonRegressionTest(filename=name, tolerance=2e-15) + # loop on all the ini files for this test + name = "dump.0001.dmp" + for ini in inifiles: + test.run(inputFile=ini) + test.standardTest() + if test.init: + test.makeReference(filename=name) + test.nonRegressionTest(filename=name, tolerance=2e-15) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp") - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp") + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) diff --git a/test/MHD/sphBragViscosity/python/testidefix.py b/test/MHD/sphBragViscosity/python/testidefix.py index 7369f636..84c60dfc 100644 --- a/test/MHD/sphBragViscosity/python/testidefix.py +++ b/test/MHD/sphBragViscosity/python/testidefix.py @@ -13,10 +13,10 @@ time_step = conf["Output"]["vtk"] nfile = 10 -list_time = [round(i*time_step, 7) for i in range(nfile)] +list_time = [round(i * time_step, 7) for i in range(nfile)] list_VTK = list() for k in range(nfile): - current_VTK = readVTK('../data.00{0:02d}.vtk'.format(k), 'spherical') + current_VTK = readVTK("../data.00{0:02d}.vtk".format(k), "spherical") list_VTK.append(current_VTK) list_VX = list() @@ -27,32 +27,33 @@ list_BY = list() list_BZ = list() for VTK in list_VTK: - list_VX.append(VTK.data['VX1']) - list_VY.append(VTK.data['VX2']) - list_VZ.append(VTK.data['VX3']) - list_BX.append(VTK.data['BX1']) - list_BY.append(VTK.data['BX2']) - list_BZ.append(VTK.data['BX3']) - list_RHO.append(VTK.data['RHO']) + list_VX.append(VTK.data["VX1"]) + list_VY.append(VTK.data["VX2"]) + list_VZ.append(VTK.data["VX3"]) + list_BX.append(VTK.data["BX1"]) + list_BY.append(VTK.data["BX2"]) + list_BZ.append(VTK.data["BX3"]) + list_RHO.append(VTK.data["RHO"]) X = VTK.r Y = VTK.theta Z = VTK.phi mu = conf["Hydro"]["bragViscosity"][-1] + def analytic_sol(r, th, phi, t): - lambd2 = 2. - return amplitude*jn(1,r)/r*np.exp(-t*lambd2*mu) + lambd2 = 2.0 + return amplitude * jn(1, r) / r * np.exp(-t * lambd2 * mu) eps = 7e-7 for k, Vx in enumerate(list_VX): t = list_time[k] - current_sol = np.zeros((len(X),len(Y),len(Z))) - for ix,x in enumerate(X): - for iy,y in enumerate(Y): - for iz,z in enumerate(Z): - current_sol[ix,iy,iz] = analytic_sol(x,y,z,t) + current_sol = np.zeros((len(X), len(Y), len(Z))) + for ix, x in enumerate(X): + for iy, y in enumerate(Y): + for iz, z in enumerate(Z): + current_sol[ix, iy, iz] = analytic_sol(x, y, z, t) if np.mean(abs(current_sol - Vx)) > eps: print("Failed") sys.exit(1) diff --git a/test/MHD/sphBragViscosity/testme.py b/test/MHD/sphBragViscosity/testme.py index 0ad8ef25..47c90c0e 100755 --- a/test/MHD/sphBragViscosity/testme.py +++ b/test/MHD/sphBragViscosity/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -13,34 +14,34 @@ def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - name = "dump.0001.dmp" - for ini in inifiles: - test.run(inputFile=ini) -## The current test is not physics, but merely a non-regression check -## The physical test needs to have the divergence of the velocity field artifiacially set to zero -# test.standardTest() - if test.init: - test.makeReference(filename=name) + # loop on all the ini files for this test + name = "dump.0001.dmp" + for ini in inifiles: + test.run(inputFile=ini) + ## The current test is not physics, but merely a non-regression check + ## The physical test needs to have the divergence of the velocity field artifiacially set to zero + # test.standardTest() + if test.init: + test.makeReference(filename=name) - test.nonRegressionTest(filename=name, tolerance=1e-15) + test.nonRegressionTest(filename=name, tolerance=1e-15) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename="dump.0001.dmp") - else: - testMe(test) + if test.check: + test.checkOnly(filename="dump.0001.dmp") + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) diff --git a/test/Planet/Planet3Body/python/testidefix.py b/test/Planet/Planet3Body/python/testidefix.py index 91a680d6..c66176f4 100755 --- a/test/Planet/Planet3Body/python/testidefix.py +++ b/test/Planet/Planet3Body/python/testidefix.py @@ -16,69 +16,91 @@ def separation(resonance): a = int(resonance[-1]) b = int(resonance[0]) - return pow(a/b,2/3) + return pow(a / b, 2 / 3) + def datafile(filename): return np.loadtxt(filename, dtype="float64").T + planet0 = datafile("../planet0.dat") planet1 = datafile("../planet1.dat") resonance = "3:2" plot_orbit = False -theta1 = np.arctan2(planet1[2],planet1[1]) +theta1 = np.arctan2(planet1[2], planet1[1]) apl1 = 1.0 ecc1 = 0.0 -rpl1 = (apl1*(1-ecc1**2)/(1-ecc1*np.cos(theta1))) -xpl1 = rpl1*np.cos(theta1) -ypl1 = rpl1*np.sin(theta1) +rpl1 = apl1 * (1 - ecc1**2) / (1 - ecc1 * np.cos(theta1)) +xpl1 = rpl1 * np.cos(theta1) +ypl1 = rpl1 * np.sin(theta1) -theta0 = np.arctan2(planet0[2],planet0[1]) +theta0 = np.arctan2(planet0[2], planet0[1]) apl0 = separation(resonance) ecc0 = 0.2 -rpl0 = (apl0*(1-ecc0**2)/(1+ecc0*np.cos(theta0))) -xpl0 = rpl0*np.cos(theta0) -ypl0 = rpl0*np.sin(theta0) +rpl0 = apl0 * (1 - ecc0**2) / (1 + ecc0 * np.cos(theta0)) +xpl0 = rpl0 * np.cos(theta0) +ypl0 = rpl0 * np.sin(theta0) if plot_orbit: import matplotlib.pyplot as plt + plt.close("all") - fig, ax = plt.subplots(figsize=(6,6)) + fig, ax = plt.subplots(figsize=(6, 6)) ax.scatter(planet0[1], planet0[2], ls="-", c="k") ax.scatter(planet1[1], planet1[2], ls="--", c="r") - ax.scatter(0,0,marker="+") + ax.scatter(0, 0, marker="+") ax.plot(xpl0, ypl0, c="b") ax.plot(xpl1, ypl1, c="g") - ax.set(xlabel="x", ylabel="y", aspect="equal", xlim=(-1.2,1.2), ylim=(-1.2,1.2)) + ax.set(xlabel="x", ylabel="y", aspect="equal", xlim=(-1.2, 1.2), ylim=(-1.2, 1.2)) # ax.set(xlabel="x", ylabel="y", aspect="equal", xlim=(-2.5,2.5), ylim=(-2.5,2.5)) - ax.set_title(r"m$_{\rm sat}=0$, e=%.1f, n/n$_s$=%s (a0=%.3f)"%(ecc0,resonance,separation(resonance=resonance)), pad=20) + ax.set_title( + r"m$_{\rm sat}=0$, e=%.1f, n/n$_s$=%s (a0=%.3f)" + % (ecc0, resonance, separation(resonance=resonance)), + pad=20, + ) # plt.savefig(f"resonance{resonance}_e{ecc}.png", dpi=200) fig.tight_layout() - fig2, ax2 = plt.subplots(figsize=(6,6)) - ax2.scatter(theta0, (rpl0-np.sqrt(planet0[1]**2 + planet0[2]**2))/(rpl0), c="b", s=200) - ax2.scatter(theta1, (rpl1-np.sqrt(planet1[1]**2 + planet1[2]**2))/(rpl1), c="g", s=200) + fig2, ax2 = plt.subplots(figsize=(6, 6)) + ax2.scatter( + theta0, + (rpl0 - np.sqrt(planet0[1] ** 2 + planet0[2] ** 2)) / (rpl0), + c="b", + s=200, + ) + ax2.scatter( + theta1, + (rpl1 - np.sqrt(planet1[1] ** 2 + planet1[2] ** 2)) / (rpl1), + c="g", + s=200, + ) plt.show() -rpl0_sim = np.sqrt(planet0[1]**2 + planet0[2]**2) -rpl1_sim = np.sqrt(planet1[1]**2 + planet1[2]**2) +rpl0_sim = np.sqrt(planet0[1] ** 2 + planet0[2] ** 2) +rpl1_sim = np.sqrt(planet1[1] ** 2 + planet1[2] ** 2) -mean0 = np.mean(100*(rpl0-rpl0_sim)/rpl0) -mean1 = np.mean(100*(rpl1-rpl1_sim)/rpl1) -std0 = np.std(100*(rpl0-rpl0_sim)/rpl0) -std1 = np.std(100*(rpl1-rpl1_sim)/rpl1) +mean0 = np.mean(100 * (rpl0 - rpl0_sim) / rpl0) +mean1 = np.mean(100 * (rpl1 - rpl1_sim) / rpl1) +std0 = np.std(100 * (rpl0 - rpl0_sim) / rpl0) +std1 = np.std(100 * (rpl1 - rpl1_sim) / rpl1) error_mean0 = mean0 error_mean1 = mean1 -error_dispersion0 = np.max([abs(mean0+std0),abs(mean0-std0)]) -error_dispersion1 = np.max([abs(mean1+std1),abs(mean1-std1)]) - -print("Mean Error (planet0)=%s"%error_mean0) -print("Mean Error (planet1)=%s"%error_mean1) -print("max(mean +/- std) (planet0)=%s"%error_dispersion0) -print("max(mean +/- std) (planet1)=%s"%error_dispersion1) -if error_mean0<1.0e-10 and error_dispersion0<1.0e-9 and error_mean1<1.0e-10 and error_dispersion1<1.0e-9: +error_dispersion0 = np.max([abs(mean0 + std0), abs(mean0 - std0)]) +error_dispersion1 = np.max([abs(mean1 + std1), abs(mean1 - std1)]) + +print("Mean Error (planet0)=%s" % error_mean0) +print("Mean Error (planet1)=%s" % error_mean1) +print("max(mean +/- std) (planet0)=%s" % error_dispersion0) +print("max(mean +/- std) (planet1)=%s" % error_dispersion1) +if ( + error_mean0 < 1.0e-10 + and error_dispersion0 < 1.0e-9 + and error_mean1 < 1.0e-10 + and error_dispersion1 < 1.0e-9 +): print("SUCCESS!") sys.exit(0) else: diff --git a/test/Planet/Planet3Body/testme.py b/test/Planet/Planet3Body/testme.py index bdb0e193..4a54d197 100755 --- a/test/Planet/Planet3Body/testme.py +++ b/test/Planet/Planet3Body/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,38 +12,39 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + +tolerance = 2e-11 -tolerance=2e-11 def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename=name) - test.standardTest() - test.nonRegressionTest(filename=name,tolerance=tolerance) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename=name) + test.standardTest() + test.nonRegressionTest(filename=name, tolerance=tolerance) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.checkOnly(filename=name,tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name, tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - test.mpi=True - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + test.mpi = True + testMe(test) diff --git a/test/Planet/PlanetMigration2D/python/testidefix.py b/test/Planet/PlanetMigration2D/python/testidefix.py index 0a47d24f..3414df98 100755 --- a/test/Planet/PlanetMigration2D/python/testidefix.py +++ b/test/Planet/PlanetMigration2D/python/testidefix.py @@ -16,36 +16,38 @@ def datafile(filename): return np.loadtxt(filename, dtype="float64").T + plot = False planet_ref = datafile("planet0.ref.dat") planet = datafile("../planet0.dat") -migration_ref = np.sqrt(planet_ref[1]**2+planet_ref[2]**2+planet_ref[3]**2) -migration = np.sqrt(planet[1]**2+planet[2]**2+planet[3]**2) +migration_ref = np.sqrt(planet_ref[1] ** 2 + planet_ref[2] ** 2 + planet_ref[3] ** 2) +migration = np.sqrt(planet[1] ** 2 + planet[2] ** 2 + planet[3] ** 2) # Compute the error on the planet distance -error_d=np.max(np.abs((migration-migration_ref)/migration_ref)) +error_d = np.max(np.abs((migration - migration_ref) / migration_ref)) torque_ref = datafile("tqwk0.ref.dat") torque = datafile("../tqwk0.dat") -tq_tot_ref = torque_ref[1]+torque_ref[2] -tq_tot = torque[1]+torque[2] +tq_tot_ref = torque_ref[1] + torque_ref[2] +tq_tot = torque[1] + torque[2] # Compute the error on the planet torque -error_t=np.max(np.abs((tq_tot-tq_tot_ref)/tq_tot_ref)) +error_t = np.max(np.abs((tq_tot - tq_tot_ref) / tq_tot_ref)) if plot: import matplotlib.pyplot as plt + plt.close("all") - fig, ax = plt.subplots(figsize=(6,6)) + fig, ax = plt.subplots(figsize=(6, 6)) ax.plot(planet_ref[-1], migration_ref, ls="-", c="k", label="dist ref") ax.plot(planet[-1], migration, ls="--", c="r", label="dist new") ax.set(xlabel="time") fig.tight_layout() ax.legend(frameon=False) - fig2, ax2 = plt.subplots(figsize=(6,6)) + fig2, ax2 = plt.subplots(figsize=(6, 6)) ax2.plot(torque_ref[-1], tq_tot_ref, ls="-", c="k", label="torque ref") ax2.plot(torque[-1], tq_tot, ls="--", c="r", label="torque new") ax2.set(xlabel="time") @@ -53,12 +55,12 @@ def datafile(filename): ax2.legend(frameon=False) plt.show() -diff_torque = np.max(np.abs(tq_tot-tq_tot_ref)) -print("diff_torque=%e"%diff_torque) -print("Error_dist=%e"%error_d) +diff_torque = np.max(np.abs(tq_tot - tq_tot_ref)) +print("diff_torque=%e" % diff_torque) +print("Error_dist=%e" % error_d) # print("Error_torque=%e"%error_t) # if error_d<5.0e-2 and error_t<5.0e-2 and error_rho<5.0e-2: -if error_d<5.0e-2 and diff_torque<1.0e-15: +if error_d < 5.0e-2 and diff_torque < 1.0e-15: print("SUCCESS!") sys.exit(0) else: diff --git a/test/Planet/PlanetMigration2D/testme.py b/test/Planet/PlanetMigration2D/testme.py index 7aff376d..94c27734 100755 --- a/test/Planet/PlanetMigration2D/testme.py +++ b/test/Planet/PlanetMigration2D/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,40 +12,41 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + +tolerance = 1e-13 -tolerance=1e-13 def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename=name) - test.standardTest() - test.nonRegressionTest(filename=name,tolerance=tolerance) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename=name) + test.standardTest() + test.nonRegressionTest(filename=name, tolerance=tolerance) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','2'] + test.dec = ["2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename=name,tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name, tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - test.mpi=True - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + test.mpi = True + testMe(test) diff --git a/test/Planet/PlanetPlanetRK42D/python/testidefix.py b/test/Planet/PlanetPlanetRK42D/python/testidefix.py index 656e55c2..fdb4aa2a 100755 --- a/test/Planet/PlanetPlanetRK42D/python/testidefix.py +++ b/test/Planet/PlanetPlanetRK42D/python/testidefix.py @@ -16,15 +16,16 @@ def filter_index(t): # get read of data previously generated (when the code was restarted) - tref=t[-1] - idx=np.array([2]) - idx[0]=t.size-1 - for i in range(t.size-2,0,-1): - if(t[i] 0: - integ = ( - np.nansum( - (data * dR)[k : kk + 1], - dtype="float64", - ) + integ = np.nansum( + (data * dR)[k : kk + 1], + dtype="float64", ) else: integ = -( @@ -97,8 +101,9 @@ def _integral( ) return integ + # Phase equation for m,n -#(phi-wavenumber,order of the spiral) +# (phi-wavenumber,order of the spiral) def phaseEquation( RR, *, @@ -108,17 +113,22 @@ def phaseEquation( Rp=Rp, info=False, ): - phieq = -np.sign(RR-Rp)*np.pi/4/m + 2*np.pi*n/m -_integral( - RR, - resonance=resonance, - m=mdom, - Rp=Rp, - h0=h0, - flaring=flaring, - info=info, - ) + phieq = ( + -np.sign(RR - Rp) * np.pi / 4 / m + + 2 * np.pi * n / m + - _integral( + RR, + resonance=resonance, + m=mdom, + Rp=Rp, + h0=h0, + flaring=flaring, + info=info, + ) + ) return phieq + plot = False RwkzMin = find_nearest(Rmed, 0.7) @@ -127,45 +137,45 @@ def phaseEquation( spiralR = [] spiralP = [] spiralP_theo = [] -for ir in range(RwkzMin, RwkzMax+1): - if (Rmed[ir] < Rmm): - rho = ds.data["RHO"][ir,:,0] +for ir in range(RwkzMin, RwkzMax + 1): + if Rmed[ir] < Rmm: + rho = ds.data["RHO"][ir, :, 0] spiralR.append(Rmed[ir]) - spiralP.append(phimed[find_nearest(rho,rho.max())]) - spiralP_theo.append(phaseEquation(Rmed[ir],resonance=-1)) - if (Rmed[ir] > Rmp): - rho = ds.data["RHO"][ir,:,0] + spiralP.append(phimed[find_nearest(rho, rho.max())]) + spiralP_theo.append(phaseEquation(Rmed[ir], resonance=-1)) + if Rmed[ir] > Rmp: + rho = ds.data["RHO"][ir, :, 0] spiralR.append(Rmed[ir]) - spiralP.append(phimed[find_nearest(rho,rho.max())]) - spiralP_theo.append(phaseEquation(Rmed[ir],resonance=+1)) + spiralP.append(phimed[find_nearest(rho, rho.max())]) + spiralP_theo.append(phaseEquation(Rmed[ir], resonance=+1)) spiralR_theo = spiralR spiralP = np.array(spiralP) -spiralP_theo = (np.array(spiralP_theo)/(np.pi))%2*np.pi-np.pi +spiralP_theo = (np.array(spiralP_theo) / (np.pi)) % 2 * np.pi - np.pi if plot: - fig, ax = plt.subplots(figsize=(5,4)) + fig, ax = plt.subplots(figsize=(5, 4)) ax.axvline(x=1.0, ls="--", c="k") - ax.scatter(spiralR_theo,spiralP_theo, c="r",marker="+", label="theoretical") - ax.scatter(spiralR,spiralP,c="k",marker="x", label=r"max $\rho$ simulation") + ax.scatter(spiralR_theo, spiralP_theo, c="r", marker="+", label="theoretical") + ax.scatter(spiralR, spiralP, c="k", marker="x", label=r"max $\rho$ simulation") fig.tight_layout() ax.legend(frameon=False) fig2, ax2 = plt.subplots() - ax2.scatter(spiralR, 100*(spiralP-spiralP_theo)/spiralP_theo) + ax2.scatter(spiralR, 100 * (spiralP - spiralP_theo) / spiralP_theo) fig2.tight_layout() plt.show() -mean = np.mean(100*(spiralP-spiralP_theo)/spiralP_theo) -std = np.std(100*(spiralP-spiralP_theo)/spiralP_theo) +mean = np.mean(100 * (spiralP - spiralP_theo) / spiralP_theo) +std = np.std(100 * (spiralP - spiralP_theo) / spiralP_theo) error_mean = mean -error_dispersion = np.max([abs(mean+std),abs(mean-std)]) +error_dispersion = np.max([abs(mean + std), abs(mean - std)]) -print("Mean Error=%.2f"%error_mean) -print("max(mean +/- std)=%.2f"%error_dispersion) -if error_mean<1.0 and error_dispersion<5.0: +print("Mean Error=%.2f" % error_mean) +print("max(mean +/- std)=%.2f" % error_dispersion) +if error_mean < 1.0 and error_dispersion < 5.0: print("SUCCESS!") sys.exit(0) else: diff --git a/test/Planet/PlanetSpiral2D/testme.py b/test/Planet/PlanetSpiral2D/testme.py index 7aff376d..94c27734 100755 --- a/test/Planet/PlanetSpiral2D/testme.py +++ b/test/Planet/PlanetSpiral2D/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,40 +12,41 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + +tolerance = 1e-13 -tolerance=1e-13 def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename=name) - test.standardTest() - test.nonRegressionTest(filename=name,tolerance=tolerance) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename=name) + test.standardTest() + test.nonRegressionTest(filename=name, tolerance=tolerance) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','2'] + test.dec = ["2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename=name,tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name, tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - test.mpi=True - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + test.mpi = True + testMe(test) diff --git a/test/Planet/PlanetTorque3D/python/testidefix.py b/test/Planet/PlanetTorque3D/python/testidefix.py index f51b4974..ce8d911d 100755 --- a/test/Planet/PlanetTorque3D/python/testidefix.py +++ b/test/Planet/PlanetTorque3D/python/testidefix.py @@ -16,6 +16,7 @@ def datafile(filename): return np.loadtxt(filename, dtype="float64").T + plot = False torque_ref = datafile("tqwk0.ref.dat") @@ -23,18 +24,19 @@ def datafile(filename): # Factor 2 here because the reference has been computed for half the disk # (hence the torque was half) -tq_tot_ref = 2*(torque_ref[1]+torque_ref[2]) -tq_tot = torque[1]+torque[2] -tqh_tot_ref = 2*(torque_ref[3]+torque_ref[4]) -tqh_tot = torque[3]+torque[4] +tq_tot_ref = 2 * (torque_ref[1] + torque_ref[2]) +tq_tot = torque[1] + torque[2] +tqh_tot_ref = 2 * (torque_ref[3] + torque_ref[4]) +tqh_tot = torque[3] + torque[4] # Compute the error on the planet torque -error_t=np.max(np.abs((tq_tot-tq_tot_ref)/tq_tot_ref)) -error_th=np.max(np.abs((tqh_tot-tqh_tot_ref)/tqh_tot_ref)) +error_t = np.max(np.abs((tq_tot - tq_tot_ref) / tq_tot_ref)) +error_th = np.max(np.abs((tqh_tot - tqh_tot_ref) / tqh_tot_ref)) if plot: import matplotlib.pyplot as plt - fig, ax = plt.subplots(figsize=(6,6)) + + fig, ax = plt.subplots(figsize=(6, 6)) ax.plot(torque_ref[-1], tq_tot_ref, ls="-", c="k", label="torque ref") ax.plot(torque[-1], tq_tot, ls="--", c="r", label="torque new") ax.plot(torque_ref[-1], tqh_tot_ref, ls="-", c="b", label="torque nohill ref") @@ -45,14 +47,14 @@ def datafile(filename): plt.show() -diff_torque = np.max(np.abs(tq_tot-tq_tot_ref)) -diff_torqueh = np.max(np.abs(tqh_tot-tqh_tot_ref)) -print("diff_torque=%e"%diff_torque) -print("diff_torque_nohill=%e"%diff_torqueh) +diff_torque = np.max(np.abs(tq_tot - tq_tot_ref)) +diff_torqueh = np.max(np.abs(tqh_tot - tqh_tot_ref)) +print("diff_torque=%e" % diff_torque) +print("diff_torque_nohill=%e" % diff_torqueh) # print("Error_torque=%e"%error_t) # print("Error_torque_nohill=%e"%error_th) # if error_t<5.0e-2 and error_th<5.0e-2 and error_rho<5.0e-2: -if diff_torque<1.0e-12 and diff_torqueh<1.0e-13: +if diff_torque < 1.0e-12 and diff_torqueh < 1.0e-13: print("SUCCESS!") sys.exit(0) else: diff --git a/test/Planet/PlanetTorque3D/testme.py b/test/Planet/PlanetTorque3D/testme.py index d68987df..ab5e8860 100755 --- a/test/Planet/PlanetTorque3D/testme.py +++ b/test/Planet/PlanetTorque3D/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,40 +12,41 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + +tolerance = 1e-13 -tolerance=1e-13 def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - if test.init and not test.mpi: - test.makeReference(filename=name) - test.standardTest() - test.nonRegressionTest(filename=name,tolerance=tolerance) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + if test.init and not test.mpi: + test.makeReference(filename=name) + test.standardTest() + test.nonRegressionTest(filename=name, tolerance=tolerance) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','2','2'] + test.dec = ["2", "2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename=name,tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name, tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - test.mpi=True - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + test.mpi = True + testMe(test) diff --git a/test/Planet/PlanetsIsActiveRK52D/python/testidefix.py b/test/Planet/PlanetsIsActiveRK52D/python/testidefix.py index c39f5178..0baa46a2 100755 --- a/test/Planet/PlanetsIsActiveRK52D/python/testidefix.py +++ b/test/Planet/PlanetsIsActiveRK52D/python/testidefix.py @@ -18,12 +18,15 @@ def find_nearest(array, value): idx = (np.abs(array - value)).argmin() return idx + def datafile(filename): return np.loadtxt(filename, dtype="float64").T + def time_dist(filename): columns = datafile(filename) - return(columns[-1], np.sqrt(columns[1]**2+columns[2]**2+columns[3]**2)) + return (columns[-1], np.sqrt(columns[1] ** 2 + columns[2] ** 2 + columns[3] ** 2)) + plot = False @@ -41,14 +44,15 @@ def time_dist(filename): time1, dist1 = time_dist("../planet1.dat") # Compute the error on the planet distance -error_d0=np.max(np.abs((dist0-dist0_rk5)/dist0_rk5)) -error_d1=np.max(np.abs((dist1-dist1_rk5)/dist1_rk5)) +error_d0 = np.max(np.abs((dist0 - dist0_rk5) / dist0_rk5)) +error_d1 = np.max(np.abs((dist1 - dist1_rk5) / dist1_rk5)) if plot: import matplotlib.pyplot as plt + plt.close("all") - fig, ax = plt.subplots(figsize=(6,6)) + fig, ax = plt.subplots(figsize=(6, 6)) ax.plot(time0_rk5, dist0_rk5, ls="-", c="k", label="dist 0 rk5") ax.plot(time0, dist0, ls="--", c="r", label="dist 0 new") ax.plot(time1_rk5, dist1_rk5, ls="-", c="b", label="dist 1 rk5") @@ -56,7 +60,7 @@ def time_dist(filename): ax.set(xlabel="time") fig.tight_layout() ax.legend(frameon=False) - fig1, ax1 = plt.subplots(figsize=(6,6)) + fig1, ax1 = plt.subplots(figsize=(6, 6)) ax1.plot(time0_active, active0, ls="-", c="k", label="isActive 0") ax1.plot(time1_active, active1, ls="--", c="r", label="isActive 1") ax1.axvline(x=time0_active[i0_active], c="k", ls=":") @@ -66,11 +70,19 @@ def time_dist(filename): ax1.legend(frameon=False) plt.show() -print("Error_dist0=%e"%error_d0) -print("Error_dist1=%e"%error_d1) -print("active0[ip0]=%e"%active0[i0_active]) -print("active1[ip1]=%e, active1[ip1-1]=%e"%(active1[i1_active],active1[i1_active-1])) -if error_d0<5.0e-10 and error_d1<5.0e-10 and active0[i0_active]==1.0 and active1[i1_active]==1.0 and active1[i1_active-1]==0.0: +print("Error_dist0=%e" % error_d0) +print("Error_dist1=%e" % error_d1) +print("active0[ip0]=%e" % active0[i0_active]) +print( + "active1[ip1]=%e, active1[ip1-1]=%e" % (active1[i1_active], active1[i1_active - 1]) +) +if ( + error_d0 < 5.0e-10 + and error_d1 < 5.0e-10 + and active0[i0_active] == 1.0 + and active1[i1_active] == 1.0 + and active1[i1_active - 1] == 0.0 +): print("SUCCESS!") sys.exit(0) else: diff --git a/test/Planet/PlanetsIsActiveRK52D/testme.py b/test/Planet/PlanetsIsActiveRK52D/testme.py index 562334d5..769be8fb 100755 --- a/test/Planet/PlanetsIsActiveRK52D/testme.py +++ b/test/Planet/PlanetsIsActiveRK52D/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,38 +12,39 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + +tolerance = 1e-13 -tolerance=1e-13 def testMe(test): - test.configure() - test.compile() - inifiles=["idefix-rk4.ini", "idefix-rk5.ini"] + test.configure() + test.compile() + inifiles = ["idefix-rk4.ini", "idefix-rk5.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - test.standardTest() - # test.nonRegressionTest(filename=name,tolerance=tolerance) + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + test.standardTest() + # test.nonRegressionTest(filename=name,tolerance=tolerance) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.dec: - test.dec=['2','2'] + test.dec = ["2", "2"] if not test.all: - if(test.check): - test.checkOnly(filename=name,tolerance=tolerance) - else: - testMe(test) + if test.check: + test.checkOnly(filename=name, tolerance=tolerance) + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) - - test.mpi=True - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) + + test.mpi = True + testMe(test) diff --git a/test/Pluto/HD/sod/python/idefixTools.py b/test/Pluto/HD/sod/python/idefixTools.py index 06e63afa..5dd73f7b 100644 --- a/test/Pluto/HD/sod/python/idefixTools.py +++ b/test/Pluto/HD/sod/python/idefixTools.py @@ -5,373 +5,385 @@ @author: glesur """ - import numpy as np class DataStructure: pass + # Read a vtk file def readVTKCart(filename): try: - fid=open(filename,"rb") + fid = open(filename, "rb") except: print("Can't open file") return 0 # define our datastructure - V=DataStructure() + V = DataStructure() # raw data which will be read from the file - V.data={} + V.data = {} - #print("Hello") + # print("Hello") # datatype we read - dt=np.dtype(">f") # Big endian single precision floats + dt = np.dtype(">f") # Big endian single precision floats - s=fid.readline() # VTK DataFile Version x.x - s=fid.readline() # Comments + s = fid.readline() # VTK DataFile Version x.x + s = fid.readline() # Comments - s=fid.readline() # BINARY - s=fid.readline() # DATASET RECTILINEAR_GRID - slist=s.split() - grid_type=str(slist[1],'utf-8') - if(grid_type != "RECTILINEAR_GRID"): + s = fid.readline() # BINARY + s = fid.readline() # DATASET RECTILINEAR_GRID + slist = s.split() + grid_type = str(slist[1], "utf-8") + if grid_type != "RECTILINEAR_GRID": print("ERROR: Wrong VTK file type.") print("This routine can only open Cartesian or Cylindrical VTK files.") fid.close() return 0 - s=fid.readline() # DIMENSIONS NX NY NZ - slist=s.split() - #s=fid.readline() # Extre line feed - V.nx=int(slist[1]) - V.ny=int(slist[2]) - V.nz=int(slist[3]) - - s=fid.readline() # X_COORDINATES NX float - - x=np.fromfile(fid,dt,V.nx) - s=fid.readline() # Extra line feed added by pluto + s = fid.readline() # DIMENSIONS NX NY NZ + slist = s.split() + # s=fid.readline() # Extre line feed + V.nx = int(slist[1]) + V.ny = int(slist[2]) + V.nz = int(slist[3]) + s = fid.readline() # X_COORDINATES NX float - s=fid.readline() # X_COORDINATES NX float + x = np.fromfile(fid, dt, V.nx) + s = fid.readline() # Extra line feed added by pluto - y=np.fromfile(fid,dt,V.ny) - s=fid.readline() # Extra line feed added by pluto + s = fid.readline() # X_COORDINATES NX float - s=fid.readline() # X_COORDINATES NX float + y = np.fromfile(fid, dt, V.ny) + s = fid.readline() # Extra line feed added by pluto - z=np.fromfile(fid,dt,V.nz) - s=fid.readline() # Extra line feed added by pluto + s = fid.readline() # X_COORDINATES NX float + z = np.fromfile(fid, dt, V.nz) + s = fid.readline() # Extra line feed added by pluto + s = fid.readline() # POINT_DATA NXNYNZ - s=fid.readline() # POINT_DATA NXNYNZ + slist = s.split() + point_type = str(slist[0], "utf-8") + npoints = int(slist[1]) + s = fid.readline() # EXTRA LINE FEED - slist=s.split() - point_type=str(slist[0],'utf-8') - npoints=int(slist[1]) - s=fid.readline() # EXTRA LINE FEED - - if(point_type == "CELL_DATA"): + if point_type == "CELL_DATA": # The file contains face coordinates, so we extrapolate to get the cell center coordinates. - if V.nx>1: - V.nx=V.nx-1 - V.x=0.5*(x[1:]+x[:-1]) + if V.nx > 1: + V.nx = V.nx - 1 + V.x = 0.5 * (x[1:] + x[:-1]) else: - V.x=x - if V.ny>1: - V.ny=V.ny-1 - V.y=0.5*(y[1:]+y[:-1]) + V.x = x + if V.ny > 1: + V.ny = V.ny - 1 + V.y = 0.5 * (y[1:] + y[:-1]) else: - V.y=y - if V.nz>1: - V.nz=V.nz-1 - V.z=0.5*(z[1:]+z[:-1]) + V.y = y + if V.nz > 1: + V.nz = V.nz - 1 + V.z = 0.5 * (z[1:] + z[:-1]) else: - V.z=z - elif(point_type == "POINT_DATA"): - V.x=x - V.y=y - V.z=z + V.z = z + elif point_type == "POINT_DATA": + V.x = x + V.y = y + V.z = z - if V.nx*V.ny*V.nz != npoints: + if V.nx * V.ny * V.nz != npoints: print("ERROR: Grid size incompatible with number of points in the data set") while 1: - s=fid.readline() # SCALARS/VECTORS name data_type (ex: SCALARS imagedata unsigned_char) - #print repr(s) - if len(s)<2: # leave if end of file + s = ( + fid.readline() + ) # SCALARS/VECTORS name data_type (ex: SCALARS imagedata unsigned_char) + # print repr(s) + if len(s) < 2: # leave if end of file break - slist=s.split() - datatype=str(slist[0],'utf-8') - varname=str(slist[1],'utf-8') + slist = s.split() + datatype = str(slist[0], "utf-8") + varname = str(slist[1], "utf-8") if datatype == "SCALARS": fid.readline() # LOOKUP TABLE - V.data[varname] = np.transpose(np.fromfile(fid,dt,V.nx*V.ny*V.nz).reshape(V.nz,V.ny,V.nx)) + V.data[varname] = np.transpose( + np.fromfile(fid, dt, V.nx * V.ny * V.nz).reshape(V.nz, V.ny, V.nx) + ) elif datatype == "VECTORS": - Q=np.fromfile(fid,dt,3*V.nx*V.ny*V.nz) + Q = np.fromfile(fid, dt, 3 * V.nx * V.ny * V.nz) - V.data[varname+'_X']=np.transpose(Q[::3].reshape(V.nz,V.ny,V.nx)) - V.data[varname+'_Y']=np.transpose(Q[1::3].reshape(V.nz,V.ny,V.nx)) - V.data[varname+'_Z']=np.transpose(Q[2::3].reshape(V.nz,V.ny,V.nx)) + V.data[varname + "_X"] = np.transpose(Q[::3].reshape(V.nz, V.ny, V.nx)) + V.data[varname + "_Y"] = np.transpose(Q[1::3].reshape(V.nz, V.ny, V.nx)) + V.data[varname + "_Z"] = np.transpose(Q[2::3].reshape(V.nz, V.ny, V.nx)) else: print("ERROR: Unknown datatype %s" % datatype) - break; + break - fid.readline() #extra line feed + fid.readline() # extra line feed fid.close() return V + # Read a vtk file def readVTKPolar(filename): try: - fid=open(filename,"rb") + fid = open(filename, "rb") except: print("Can't open file") return 0 # define our datastructure - V=DataStructure() + V = DataStructure() # raw data which will be read from the file - V.data={} + V.data = {} - #print("Hello") + # print("Hello") # datatype we read - dt=np.dtype(">f") # Big endian single precision floats + dt = np.dtype(">f") # Big endian single precision floats - s=fid.readline() # VTK DataFile Version x.x - s=fid.readline() # Comments + s = fid.readline() # VTK DataFile Version x.x + s = fid.readline() # Comments - s=fid.readline() # BINARY - s=fid.readline() # DATASET RECTILINEAR_GRID + s = fid.readline() # BINARY + s = fid.readline() # DATASET RECTILINEAR_GRID print(s) - slist=s.split() - grid_type=str(slist[1],'utf-8') - if(grid_type != "STRUCTURED_GRID"): + slist = s.split() + grid_type = str(slist[1], "utf-8") + if grid_type != "STRUCTURED_GRID": print("ERROR: Wrong VTK file type.") - print("Current type is: %s"%(grid_type)) + print("Current type is: %s" % (grid_type)) print("This routine can only open Polar VTK files.") fid.close() return 0 - s=fid.readline() # DIMENSIONS NX NY NZ - slist=s.split() - V.nx=int(slist[1]) - V.ny=int(slist[2]) - V.nz=int(slist[3]) + s = fid.readline() # DIMENSIONS NX NY NZ + slist = s.split() + V.nx = int(slist[1]) + V.ny = int(slist[2]) + V.nz = int(slist[3]) - print("nx=%d, ny=%d, nz=%d"%(V.nx,V.ny,V.nz)) + print("nx=%d, ny=%d, nz=%d" % (V.nx, V.ny, V.nz)) - s=fid.readline() # POINTS NXNYNZ float - slist=s.split() - npoints=int(slist[1]) - points=np.fromfile(fid,dt,3*npoints) - s=fid.readline() # EXTRA LINE FEED + s = fid.readline() # POINTS NXNYNZ float + slist = s.split() + npoints = int(slist[1]) + points = np.fromfile(fid, dt, 3 * npoints) + s = fid.readline() # EXTRA LINE FEED - V.points=points + V.points = points - if V.nx*V.ny*V.nz != npoints: + if V.nx * V.ny * V.nz != npoints: print("ERROR: Grid size incompatible with number of points in the data set") return 0 # Reconstruct the polar coordinate system - x1d=points[::3] - y1d=points[1::3] - z1d=points[2::3] - - xcart=np.transpose(x1d.reshape(V.nz,V.ny,V.nx)) - ycart=np.transpose(y1d.reshape(V.nz,V.ny,V.nx)) - zcart=np.transpose(z1d.reshape(V.nz,V.ny,V.nx)) - - r=np.sqrt(xcart[:,0,0]**2+ycart[:,0,0]**2) - theta=np.unwrap(np.arctan2(ycart[0,:,0],xcart[0,:,0])) - z=zcart[0,0,:] - - s=fid.readline() # CELL_DATA (NX-1)(NY-1)(NZ-1) - slist=s.split() - data_type=str(slist[0],'utf-8') - if(data_type != "CELL_DATA"): + x1d = points[::3] + y1d = points[1::3] + z1d = points[2::3] + + xcart = np.transpose(x1d.reshape(V.nz, V.ny, V.nx)) + ycart = np.transpose(y1d.reshape(V.nz, V.ny, V.nx)) + zcart = np.transpose(z1d.reshape(V.nz, V.ny, V.nx)) + + r = np.sqrt(xcart[:, 0, 0] ** 2 + ycart[:, 0, 0] ** 2) + theta = np.unwrap(np.arctan2(ycart[0, :, 0], xcart[0, :, 0])) + z = zcart[0, 0, :] + + s = fid.readline() # CELL_DATA (NX-1)(NY-1)(NZ-1) + slist = s.split() + data_type = str(slist[0], "utf-8") + if data_type != "CELL_DATA": print("ERROR: this routine expect CELL DATA as produced by PLUTO.") fid.close() return 0 - s=fid.readline() # Line feed + s = fid.readline() # Line feed # Perform averaging on coordinate system to get cell centers # The file contains face coordinates, so we extrapolate to get the cell center coordinates. - if V.nx>1: - V.nx=V.nx-1 - V.x=0.5*(r[1:]+r[:-1]) + if V.nx > 1: + V.nx = V.nx - 1 + V.x = 0.5 * (r[1:] + r[:-1]) else: - V.x=r - if V.ny>1: - V.ny=V.ny-1 - V.y=(0.5*(theta[1:]+theta[:-1])+np.pi)%(2.0*np.pi)-np.pi + V.x = r + if V.ny > 1: + V.ny = V.ny - 1 + V.y = (0.5 * (theta[1:] + theta[:-1]) + np.pi) % (2.0 * np.pi) - np.pi else: - V.y=theta - if V.nz>1: - V.nz=V.nz-1 - V.z=0.5*(z[1:]+z[:-1]) + V.y = theta + if V.nz > 1: + V.nz = V.nz - 1 + V.z = 0.5 * (z[1:] + z[:-1]) else: - V.z=z - + V.z = z while 1: - s=fid.readline() # SCALARS/VECTORS name data_type (ex: SCALARS imagedata unsigned_char) - #print repr(s) - if len(s)<2: # leave if end of file + s = ( + fid.readline() + ) # SCALARS/VECTORS name data_type (ex: SCALARS imagedata unsigned_char) + # print repr(s) + if len(s) < 2: # leave if end of file break - slist=s.split() - datatype=str(slist[0],'utf-8') - varname=str(slist[1],'utf-8') + slist = s.split() + datatype = str(slist[0], "utf-8") + varname = str(slist[1], "utf-8") if datatype == "SCALARS": fid.readline() # LOOKUP TABLE - V.data[varname] = np.transpose(np.fromfile(fid,dt,V.nx*V.ny*V.nz).reshape(V.nz,V.ny,V.nx)) + V.data[varname] = np.transpose( + np.fromfile(fid, dt, V.nx * V.ny * V.nz).reshape(V.nz, V.ny, V.nx) + ) elif datatype == "VECTORS": - Q=np.fromfile(fid,dt,3*V.nx*V.ny*V.nz) + Q = np.fromfile(fid, dt, 3 * V.nx * V.ny * V.nz) - V.data[varname+'_X']=np.transpose(Q[::3].reshape(V.nz,V.ny,V.nx)) - V.data[varname+'_Y']=np.transpose(Q[1::3].reshape(V.nz,V.ny,V.nx)) - V.data[varname+'_Z']=np.transpose(Q[2::3].reshape(V.nz,V.ny,V.nx)) + V.data[varname + "_X"] = np.transpose(Q[::3].reshape(V.nz, V.ny, V.nx)) + V.data[varname + "_Y"] = np.transpose(Q[1::3].reshape(V.nz, V.ny, V.nx)) + V.data[varname + "_Z"] = np.transpose(Q[2::3].reshape(V.nz, V.ny, V.nx)) else: print("ERROR: Unknown datatype %s" % datatype) - break; + break - fid.readline() #extra line feed + fid.readline() # extra line feed fid.close() return V + # Read a vtk file def readVTKSpherical(filename): try: - fid=open(filename,"rb") + fid = open(filename, "rb") except: print("Can't open file") return 0 # define our datastructure - V=DataStructure() + V = DataStructure() # raw data which will be read from the file - V.data={} + V.data = {} - #print("Hello") + # print("Hello") # datatype we read - dt=np.dtype(">f") # Big endian single precision floats + dt = np.dtype(">f") # Big endian single precision floats - s=fid.readline() # VTK DataFile Version x.x - s=fid.readline() # Comments + s = fid.readline() # VTK DataFile Version x.x + s = fid.readline() # Comments - s=fid.readline() # BINARY - s=fid.readline() # DATASET RECTILINEAR_GRID - slist=s.split() - grid_type=str(slist[1],'utf-8') - if(grid_type != "STRUCTURED_GRID"): + s = fid.readline() # BINARY + s = fid.readline() # DATASET RECTILINEAR_GRID + slist = s.split() + grid_type = str(slist[1], "utf-8") + if grid_type != "STRUCTURED_GRID": print("ERROR: Wrong VTK file type.") print("This routine can only open Spherical VTK files.") fid.close() return 0 - s=fid.readline() # DIMENSIONS NX NY NZ - slist=s.split() - V.nx=int(slist[1]) - V.ny=int(slist[2]) - V.nz=int(slist[3]) + s = fid.readline() # DIMENSIONS NX NY NZ + slist = s.split() + V.nx = int(slist[1]) + V.ny = int(slist[2]) + V.nz = int(slist[3]) - if(V.nz==1): - is2d=1 + if V.nz == 1: + is2d = 1 else: - is2d=0 + is2d = 0 - s=fid.readline() # POINTS NXNYNZ float - slist=s.split() - npoints=int(slist[1]) - points=np.fromfile(fid,dt,3*npoints) - s=fid.readline() # EXTRA LINE FEED + s = fid.readline() # POINTS NXNYNZ float + slist = s.split() + npoints = int(slist[1]) + points = np.fromfile(fid, dt, 3 * npoints) + s = fid.readline() # EXTRA LINE FEED - V.points=points + V.points = points - if V.nx*V.ny*V.nz != npoints: + if V.nx * V.ny * V.nz != npoints: print("ERROR: Grid size incompatible with number of points in the data set") return 0 # Reconstruct the spherical coordinate system + x1d = points[::3] + y1d = points[1::3] + z1d = points[2::3] - x1d=points[::3] - y1d=points[1::3] - z1d=points[2::3] + xcart = np.transpose(x1d.reshape(V.nz, V.ny, V.nx)) + ycart = np.transpose(y1d.reshape(V.nz, V.ny, V.nx)) + zcart = np.transpose(z1d.reshape(V.nz, V.ny, V.nx)) - xcart=np.transpose(x1d.reshape(V.nz,V.ny,V.nx)) - ycart=np.transpose(y1d.reshape(V.nz,V.ny,V.nx)) - zcart=np.transpose(z1d.reshape(V.nz,V.ny,V.nx)) - - if(is2d): - r=np.sqrt(xcart[:,0,0]**2+ycart[:,0,0]**2) - phi=np.unwrap(np.arctan2(zcart[0,0,:],xcart[0,0,:])) - theta=np.arccos(ycart[0,:,0]/np.sqrt(xcart[0,:,0]**2+ycart[0,:,0]**2)) + if is2d: + r = np.sqrt(xcart[:, 0, 0] ** 2 + ycart[:, 0, 0] ** 2) + phi = np.unwrap(np.arctan2(zcart[0, 0, :], xcart[0, 0, :])) + theta = np.arccos( + ycart[0, :, 0] / np.sqrt(xcart[0, :, 0] ** 2 + ycart[0, :, 0] ** 2) + ) else: - r=np.sqrt(xcart[:,0,0]**2+ycart[:,0,0]**2+zcart[:,0,0]**2) - phi=np.unwrap(np.arctan2(ycart[0,0,:],xcart[0,0,:])) - theta=np.arccos(zcart[0,:,0]/np.sqrt(xcart[0,:,0]**2+ycart[0,:,0]**2+zcart[0,:,0]**2)) - - - s=fid.readline() # CELL_DATA (NX-1)(NY-1)(NZ-1) - slist=s.split() - data_type=str(slist[0],'utf-8') - if(data_type != "CELL_DATA"): + r = np.sqrt(xcart[:, 0, 0] ** 2 + ycart[:, 0, 0] ** 2 + zcart[:, 0, 0] ** 2) + phi = np.unwrap(np.arctan2(ycart[0, 0, :], xcart[0, 0, :])) + theta = np.arccos( + zcart[0, :, 0] + / np.sqrt(xcart[0, :, 0] ** 2 + ycart[0, :, 0] ** 2 + zcart[0, :, 0] ** 2) + ) + + s = fid.readline() # CELL_DATA (NX-1)(NY-1)(NZ-1) + slist = s.split() + data_type = str(slist[0], "utf-8") + if data_type != "CELL_DATA": print("ERROR: this routine expect CELL DATA as produced by PLUTO.") fid.close() return 0 - s=fid.readline() # Line feed + s = fid.readline() # Line feed # Perform averaging on coordinate system to get cell centers # The file contains face coordinates, so we extrapolate to get the cell center coordinates. - if V.nx>1: - V.nx=V.nx-1 - V.r=0.5*(r[1:]+r[:-1]) + if V.nx > 1: + V.nx = V.nx - 1 + V.r = 0.5 * (r[1:] + r[:-1]) else: - V.x=r - if V.ny>1: - V.ny=V.ny-1 - V.theta=0.5*(theta[1:]+theta[:-1]) + V.x = r + if V.ny > 1: + V.ny = V.ny - 1 + V.theta = 0.5 * (theta[1:] + theta[:-1]) else: - V.y=theta - if V.nz>1: - V.nz=V.nz-1 - V.phi=(0.5*(phi[1:]+phi[:-1])+np.pi)%(2.0*np.pi)-np.pi + V.y = theta + if V.nz > 1: + V.nz = V.nz - 1 + V.phi = (0.5 * (phi[1:] + phi[:-1]) + np.pi) % (2.0 * np.pi) - np.pi else: - V.phi=phi - + V.phi = phi while 1: - s=fid.readline() # SCALARS/VECTORS name data_type (ex: SCALARS imagedata unsigned_char) - #print repr(s) - if len(s)<2: # leave if end of file + s = ( + fid.readline() + ) # SCALARS/VECTORS name data_type (ex: SCALARS imagedata unsigned_char) + # print repr(s) + if len(s) < 2: # leave if end of file break - slist=s.split() - datatype=str(slist[0],'utf-8') - varname=str(slist[1],'utf-8') + slist = s.split() + datatype = str(slist[0], "utf-8") + varname = str(slist[1], "utf-8") if datatype == "SCALARS": fid.readline() # LOOKUP TABLE - V.data[varname] = np.transpose(np.fromfile(fid,dt,V.nx*V.ny*V.nz).reshape(V.nz,V.ny,V.nx)) + V.data[varname] = np.transpose( + np.fromfile(fid, dt, V.nx * V.ny * V.nz).reshape(V.nz, V.ny, V.nx) + ) elif datatype == "VECTORS": - Q=np.fromfile(fid,dt,3*V.nx*V.ny*V.nz) + Q = np.fromfile(fid, dt, 3 * V.nx * V.ny * V.nz) - V.data[varname+'_X']=np.transpose(Q[::3].reshape(V.nz,V.ny,V.nx)) - V.data[varname+'_Y']=np.transpose(Q[1::3].reshape(V.nz,V.ny,V.nx)) - V.data[varname+'_Z']=np.transpose(Q[2::3].reshape(V.nz,V.ny,V.nx)) + V.data[varname + "_X"] = np.transpose(Q[::3].reshape(V.nz, V.ny, V.nx)) + V.data[varname + "_Y"] = np.transpose(Q[1::3].reshape(V.nz, V.ny, V.nx)) + V.data[varname + "_Z"] = np.transpose(Q[2::3].reshape(V.nz, V.ny, V.nx)) else: print("ERROR: Unknown datatype %s" % datatype) - break; + break - fid.readline() #extra line feed + fid.readline() # extra line feed fid.close() return V diff --git a/test/Pluto/HD/sod/python/sod.py b/test/Pluto/HD/sod/python/sod.py index 90ad62fe..99e55b7c 100755 --- a/test/Pluto/HD/sod/python/sod.py +++ b/test/Pluto/HD/sod/python/sod.py @@ -13,7 +13,7 @@ import scipy.optimize -def sound_speed(gamma, pressure, density, dustFrac=0.): +def sound_speed(gamma, pressure, density, dustFrac=0.0): """ Calculate sound speed, scaled by the dust fraction according to: @@ -23,27 +23,28 @@ def sound_speed(gamma, pressure, density, dustFrac=0.): Where :math:`\epsilon` is the dustFrac """ scale = np.sqrt(1 - dustFrac) - return np.sqrt(gamma * pressure/ density) * scale + return np.sqrt(gamma * pressure / density) * scale -def shock_tube_function(p4, p1, p5, rho1, rho5, gamma, dustFrac=0.): + +def shock_tube_function(p4, p1, p5, rho1, rho5, gamma, dustFrac=0.0): """ Shock tube equation """ - z = (p4 / p5 - 1.) + z = p4 / p5 - 1.0 c1 = sound_speed(gamma, p1, rho1, dustFrac) c5 = sound_speed(gamma, p5, rho5, dustFrac) - gm1 = gamma - 1. - gp1 = gamma + 1. - g2 = 2. * gamma + gm1 = gamma - 1.0 + gp1 = gamma + 1.0 + g2 = 2.0 * gamma - fact = gm1 / g2 * (c5 / c1) * z / np.sqrt(1. + gp1 / g2 * z) - fact = (1. - fact) ** (g2 / gm1) + fact = gm1 / g2 * (c5 / c1) * z / np.sqrt(1.0 + gp1 / g2 * z) + fact = (1.0 - fact) ** (g2 / gm1) return p1 * fact - p4 -def calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma=1.4, dustFrac=0.): +def calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma=1.4, dustFrac=0.0): """ Compute regions :rtype : tuple @@ -70,18 +71,18 @@ def calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma=1.4, dustFrac=0.): p4 = scipy.optimize.fsolve(shock_tube_function, p1, (p1, p5, rho1, rho5, gamma))[0] # compute post-shock density and velocity - z = (p4 / p5 - 1.) + z = p4 / p5 - 1.0 c5 = sound_speed(gamma, p5, rho5, dustFrac) - gm1 = gamma - 1. - gp1 = gamma + 1. + gm1 = gamma - 1.0 + gp1 = gamma + 1.0 gmfac1 = 0.5 * gm1 / gamma gmfac2 = 0.5 * gp1 / gamma - fact = np.sqrt(1. + gmfac2 * z) + fact = np.sqrt(1.0 + gmfac2 * z) u4 = c5 * z / (gamma * fact) - rho4 = rho5 * (1. + gmfac2 * z) / (1. + gmfac1 * z) + rho4 = rho5 * (1.0 + gmfac2 * z) / (1.0 + gmfac1 * z) # shock speed w = c5 * fact @@ -89,11 +90,11 @@ def calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma=1.4, dustFrac=0.): # compute values at foot of rarefaction p3 = p4 u3 = u4 - rho3 = rho1 * (p3 / p1)**(1. / gamma) + rho3 = rho1 * (p3 / p1) ** (1.0 / gamma) return (p1, rho1, u1), (p3, rho3, u3), (p4, rho4, u4), (p5, rho5, u5), w -def calc_positions(pl, pr, region1, region3, w, xi, t, gamma, dustFrac=0.): +def calc_positions(pl, pr, region1, region3, w, xi, t, gamma, dustFrac=0.0): """ :return: tuple of positions in the following order -> Head of Rarefaction: xhd, Foot of Rarefaction: xft, @@ -125,21 +126,39 @@ def region_states(pl, pr, region1, region3, region4, region5): where the value is a string, obviously """ if pl > pr: - return {'Region 1': region1, - 'Region 2': 'RAREFACTION', - 'Region 3': region3, - 'Region 4': region4, - 'Region 5': region5} + return { + "Region 1": region1, + "Region 2": "RAREFACTION", + "Region 3": region3, + "Region 4": region4, + "Region 5": region5, + } else: - return {'Region 1': region5, - 'Region 2': region4, - 'Region 3': region3, - 'Region 4': 'RAREFACTION', - 'Region 5': region1} - - -def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, - npts, gamma, t, xi, dustFrac=0.): + return { + "Region 1": region5, + "Region 2": region4, + "Region 3": region3, + "Region 4": "RAREFACTION", + "Region 5": region1, + } + + +def create_arrays( + pl, + pr, + xl, + xr, + positions, + state1, + state3, + state4, + state5, + npts, + gamma, + t, + xi, + dustFrac=0.0, +): """ :return: tuple of x, p, rho and u values across the domain of interest """ @@ -148,8 +167,8 @@ def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, p3, rho3, u3 = state3 p4, rho4, u4 = state4 p5, rho5, u5 = state5 - gm1 = gamma - 1. - gp1 = gamma + 1. + gm1 = gamma - 1.0 + gp1 = gamma + 1.0 x_arr = np.linspace(xl, xr, npts) rho = np.zeros(npts, dtype=float) @@ -163,10 +182,10 @@ def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, p[i] = p1 u[i] = u1 elif x < xft: - u[i] = 2. / gp1 * (c1 + (x - xi) / t) - fact = 1. - 0.5 * gm1 * u[i] / c1 - rho[i] = rho1 * fact ** (2. / gm1) - p[i] = p1 * fact ** (2. * gamma / gm1) + u[i] = 2.0 / gp1 * (c1 + (x - xi) / t) + fact = 1.0 - 0.5 * gm1 * u[i] / c1 + rho[i] = rho1 * fact ** (2.0 / gm1) + p[i] = p1 * fact ** (2.0 * gamma / gm1) elif x < xcd: rho[i] = rho3 p[i] = p3 @@ -194,10 +213,10 @@ def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, p[i] = p3 u[i] = -u3 elif x < xhd: - u[i] = -2. / gp1 * (c1 + (xi - x) / t) - fact = 1. + 0.5 * gm1 * u[i] / c1 - rho[i] = rho1 * fact ** (2. / gm1) - p[i] = p1 * fact ** (2. * gamma / gm1) + u[i] = -2.0 / gp1 * (c1 + (xi - x) / t) + fact = 1.0 + 0.5 * gm1 * u[i] / c1 + rho[i] = rho1 * fact ** (2.0 / gm1) + p[i] = p1 * fact ** (2.0 * gamma / gm1) else: rho[i] = rho1 p[i] = p1 @@ -206,8 +225,7 @@ def create_arrays(pl, pr, xl, xr, positions, state1, state3, state4, state5, return x_arr, p, rho, u -def solve(left_state, right_state, geometry, t, gamma=1.4, npts=500, - dustFrac=0.): +def solve(left_state, right_state, geometry, t, gamma=1.4, npts=500, dustFrac=0.0): """ Solves the Sod shock tube problem (i.e. riemann problem) of discontinuity across an interface. @@ -248,34 +266,57 @@ def solve(left_state, right_state, geometry, t, gamma=1.4, npts=500, # basic checking if xl >= xr: - print('xl has to be less than xr!') + print("xl has to be less than xr!") exit() if xi >= xr or xi <= xl: - print('xi has in between xl and xr!') + print("xi has in between xl and xr!") exit() # calculate regions - region1, region3, region4, region5, w = \ - calculate_regions(pl, ul, rhol, pr, ur, rhor, gamma, dustFrac) + region1, region3, region4, region5, w = calculate_regions( + pl, ul, rhol, pr, ur, rhor, gamma, dustFrac + ) regions = region_states(pl, pr, region1, region3, region4, region5) # calculate positions - x_positions = calc_positions(pl, pr, region1, region3, w, xi, t, gamma, - dustFrac) - - pos_description = ('Head of Rarefaction', 'Foot of Rarefaction', - 'Contact Discontinuity', 'Shock') + x_positions = calc_positions(pl, pr, region1, region3, w, xi, t, gamma, dustFrac) + + pos_description = ( + "Head of Rarefaction", + "Foot of Rarefaction", + "Contact Discontinuity", + "Shock", + ) positions = dict(zip(pos_description, x_positions, strict=True)) # create arrays - x, p, rho, u = create_arrays(pl, pr, xl, xr, x_positions, - region1, region3, region4, region5, - npts, gamma, t, xi, dustFrac) - - energy = p/(rho * (gamma - 1.0)) - rho_total = rho/(1.0 - dustFrac) - val_dict = {'x':x, 'p':p, 'rho':rho, 'u':u, 'energy':energy, - 'rho_total':rho_total} + x, p, rho, u = create_arrays( + pl, + pr, + xl, + xr, + x_positions, + region1, + region3, + region4, + region5, + npts, + gamma, + t, + xi, + dustFrac, + ) + + energy = p / (rho * (gamma - 1.0)) + rho_total = rho / (1.0 - dustFrac) + val_dict = { + "x": x, + "p": p, + "rho": rho, + "u": u, + "energy": energy, + "rho_total": rho_total, + } return positions, regions, val_dict diff --git a/test/Pluto/HD/sod/python/testidefix.py b/test/Pluto/HD/sod/python/testidefix.py index 239bb1c8..bdd6daf1 100755 --- a/test/Pluto/HD/sod/python/testidefix.py +++ b/test/Pluto/HD/sod/python/testidefix.py @@ -16,15 +16,14 @@ from scipy.interpolate import interp1d parser = argparse.ArgumentParser() -parser.add_argument("-noplot", - default=False, - help="disable plotting", - action="store_true") +parser.add_argument( + "-noplot", default=False, help="disable plotting", action="store_true" +) -args, unknown=parser.parse_known_args() +args, unknown = parser.parse_known_args() -V=idfx.readVTKCart('../data.0002.vtk') +V = idfx.readVTKCart("../data.0002.vtk") gamma = 1.4 npts = 5000 @@ -34,42 +33,48 @@ # t is the time evolution for which positions and states in tube should be calculated # gamma denotes specific heat # note that gamma and npts are default parameters (1.4 and 500) in solve function -positions, regions, values = sod.solve(left_state=(1, 1, 0), right_state=(0.1, 0.125, 0.), - geometry=(0., 1., 0.5), t=0.2, gamma=gamma, npts=npts) +positions, regions, values = sod.solve( + left_state=(1, 1, 0), + right_state=(0.1, 0.125, 0.0), + geometry=(0.0, 1.0, 0.5), + t=0.2, + gamma=gamma, + npts=npts, +) # Finally, let's plot solutions -p = values['p'] -rho = values['rho'] -u = values['u'] -x= values['x'] +p = values["p"] +rho = values["rho"] +u = values["u"] +x = values["x"] -solinterp=interp1d(x,p) +solinterp = interp1d(x, p) -if(not args.noplot): +if not args.noplot: plt.figure(1) - plt.plot(x,rho) - plt.plot(V.x,V.data['rho'][:,0,0],'+',markersize=2) - plt.title('Density') + plt.plot(x, rho) + plt.plot(V.x, V.data["rho"][:, 0, 0], "+", markersize=2) + plt.title("Density") plt.figure(2) - plt.plot(x,u) - plt.plot(V.x,V.data['vx1'][:,0,0],'+',markersize=2) - plt.title('Velocity') + plt.plot(x, u) + plt.plot(V.x, V.data["vx1"][:, 0, 0], "+", markersize=2) + plt.title("Velocity") plt.figure(3) - plt.plot(x,p) - plt.plot(V.x,V.data['prs'][:,0,0],'+',markersize=2) - plt.title('Pressure') + plt.plot(x, p) + plt.plot(V.x, V.data["prs"][:, 0, 0], "+", markersize=2) + plt.title("Pressure") plt.ioff() plt.show() -error=np.mean(np.fabs(V.data['prs'][:,0,0]-solinterp(V.x))) -print("Error=%e"%error) -if error<1e-2: +error = np.mean(np.fabs(V.data["prs"][:, 0, 0] - solinterp(V.x))) +print("Error=%e" % error) +if error < 1e-2: print("SUCCESS!") sys.exit(0) else: diff --git a/test/SelfGravity/DustyCollapse/python/testidefix.py b/test/SelfGravity/DustyCollapse/python/testidefix.py index 83ec380a..1ace55c5 100644 --- a/test/SelfGravity/DustyCollapse/python/testidefix.py +++ b/test/SelfGravity/DustyCollapse/python/testidefix.py @@ -6,18 +6,19 @@ sys.path.append(os.getenv("IDEFIX_DIR")) from pytools.vtk_io import readVTK -G=6.6743e-11 -M=1e16 #mass of spheric clump -r0=1e4 #radius of spheric clump -rmax = 3e4 #right edge of simulation box -d=M/(4*np.pi*r0**3/3) #density of spheric clump -phi0 = -22.247666667 #gravitationnal potential at rmax +G = 6.6743e-11 +M = 1e16 # mass of spheric clump +r0 = 1e4 # radius of spheric clump +rmax = 3e4 # right edge of simulation box +d = M / (4 * np.pi * r0**3 / 3) # density of spheric clump +phi0 = -22.247666667 # gravitationnal potential at rmax + def potential(r): retv = np.empty_like(r) - interior = (r0.01*plateau[-1]: + if np.abs(dsty - plateau[-1]) > 0.01 * plateau[-1]: return np.mean(plateau) else: plateau.append(dsty) return "Error : the whole density distribution is approximately constant !" + # Isolating density plateau -plateau=get_plateau(rho) +plateau = get_plateau(rho) # Calculating absolute error -r_th=time/tff # Theoretical ratio -eta = np.arccos((plateau/rho0)**(-1./6.)) -r_num = 2./np.pi*(eta + 1./2.*np.sin(2.*eta)) # Numerical ratio -error = np.abs(r_num-r_th) # Absolute error +r_th = time / tff # Theoretical ratio +eta = np.arccos((plateau / rho0) ** (-1.0 / 6.0)) +r_num = 2.0 / np.pi * (eta + 1.0 / 2.0 * np.sin(2.0 * eta)) # Numerical ratio +error = np.abs(r_num - r_th) # Absolute error # Print the result of the test -print("Error=%e"%error) -if error<2.0e-3: +print("Error=%e" % error) +if error < 2.0e-3: print("SUCCESS!") sys.exit(0) else: diff --git a/test/SelfGravity/UniformCollapse/testme.py b/test/SelfGravity/UniformCollapse/testme.py index f011cd21..123f417d 100755 --- a/test/SelfGravity/UniformCollapse/testme.py +++ b/test/SelfGravity/UniformCollapse/testme.py @@ -4,6 +4,7 @@ @author: glesur """ + import os import sys @@ -11,32 +12,33 @@ import pytools.idfx_test as tst -name="dump.0001.dmp" +name = "dump.0001.dmp" + +tolerance = 1e-13 -tolerance=1e-13 def testMe(test): - test.configure() - test.compile() - inifiles=["idefix.ini"] + test.configure() + test.compile() + inifiles = ["idefix.ini"] - # loop on all the ini files for this test - for ini in inifiles: - test.run(inputFile=ini) - test.standardTest() + # loop on all the ini files for this test + for ini in inifiles: + test.run(inputFile=ini) + test.standardTest() -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) if not test.all: - if(test.check): - test.standardTest() - else: - testMe(test) + if test.check: + test.standardTest() + else: + testMe(test) else: - test.noplot = True - test.vectPot=False - test.single=False - test.reconstruction=2 - test.mpi=False - testMe(test) + test.noplot = True + test.vectPot = False + test.single = False + test.reconstruction = 2 + test.mpi = False + testMe(test) diff --git a/test/utils/columnDensity/testme.py b/test/utils/columnDensity/testme.py index 4197a936..b0a6b573 100755 --- a/test/utils/columnDensity/testme.py +++ b/test/utils/columnDensity/testme.py @@ -4,13 +4,14 @@ @author: glesur """ + import os import sys sys.path.append(os.getenv("IDEFIX_DIR")) import pytools.idfx_test as tst -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) test.configure() test.compile() diff --git a/test/utils/dumpImage/testme.py b/test/utils/dumpImage/testme.py index 4197a936..b0a6b573 100755 --- a/test/utils/dumpImage/testme.py +++ b/test/utils/dumpImage/testme.py @@ -4,13 +4,14 @@ @author: glesur """ + import os import sys sys.path.append(os.getenv("IDEFIX_DIR")) import pytools.idfx_test as tst -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) test.configure() test.compile() diff --git a/test/utils/lookupTable/testme.py b/test/utils/lookupTable/testme.py index 1bc69a8f..c02f2b1c 100755 --- a/test/utils/lookupTable/testme.py +++ b/test/utils/lookupTable/testme.py @@ -4,36 +4,36 @@ @author: glesur """ + import os import sys sys.path.append(os.getenv("IDEFIX_DIR")) import numpy as np -#from scipy.interpolate import RegularGridInterpolator +# from scipy.interpolate import RegularGridInterpolator import pytools.idfx_test as tst def MakeNumpyFile(): - x=np.arange(1,10,1.0) - y=np.arange(5,10,1.0) - z=np.arange(2,5,1.0) - - xp, yp, zp = np.meshgrid(x,y,z,indexing='ij') + x = np.arange(1, 10, 1.0) + y = np.arange(5, 10, 1.0) + z = np.arange(2, 5, 1.0) - data=xp+2*yp-zp + xp, yp, zp = np.meshgrid(x, y, z, indexing="ij") - np.save("x.npy",x) - np.save("y.npy",y) - np.save("z.npy",z) - np.save("data.npy",data) - # show the expected result - #f=RegularGridInterpolator((x, y, z), data) - #print(f([2.7,7.4,3.9])) + data = xp + 2 * yp - zp + np.save("x.npy", x) + np.save("y.npy", y) + np.save("z.npy", z) + np.save("data.npy", data) + # show the expected result + # f=RegularGridInterpolator((x, y, z), data) + # print(f([2.7,7.4,3.9])) -test=tst.idfxTest(__file__) +test = tst.idfxTest(__file__) MakeNumpyFile() test.configure()