diff --git a/.gitignore b/.gitignore index 5056a889a..d9bee1a0a 100644 --- a/.gitignore +++ b/.gitignore @@ -4,7 +4,7 @@ !*/ # whitelist only the files that ever need to be tracked !*.f90 -!*.sh +swiftest_driver.sh !CMakeLists.txt !*.cmake !COPYING.txt @@ -23,6 +23,8 @@ dump* !examples/** *ipynb_checkpoints **/.DS_Store +!python/swiftest/tests/test_suite.py +!version.txt #Documentation !docs/* diff --git a/CMakeLists.txt b/CMakeLists.txt index 4d36d9612..558c93fab 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -13,11 +13,12 @@ # Define the project and the depencies that it has ################################################## -CMAKE_MINIMUM_REQUIRED(VERSION 2.8.5...3.20.1) -PROJECT(Swiftest Fortran) +CMAKE_MINIMUM_REQUIRED(VERSION 3.20.1) +# Get version stored in text file +FILE(READ "version.txt" VERSION) +PROJECT(Swiftest VERSION ${VERSION} LANGUAGES Fortran) -# Set the Swiftest version -SET(VERSION 1.0.0) +INCLUDE(CTest) # Add our local modlues to the module path SET(CMAKE_MODULE_PATH "${CMAKE_SOURCE_DIR}/cmake/Modules/") @@ -27,29 +28,32 @@ IF(NOT CMAKE_Fortran_COMPILER_SUPPORTS_F90) MESSAGE(FATAL_ERROR "Fortran compiler does not support F90") ENDIF(NOT CMAKE_Fortran_COMPILER_SUPPORTS_F90) + +IF (CMAKE_Fortran_COMPILER_ID MATCHES "^Intel") + SET(COMPILER_OPTIONS "Intel" CACHE STRING "Compiler identified as Intel") +ELSEIF (CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") + SET(COMPILER_OPTIONS "GNU" CACHE STRING "Compiler identified as gfortran") +ELSE () + MESSAGE(FATAL_ERROR "Compiler not recognized!") +ENDIF () + # Set some options the user may choose OPTION(USE_COARRAY "Use Coarray Fortran for parallelization of test particles" OFF) OPTION(USE_OPENMP "Use OpenMP for parallelization" ON) OPTION(USE_SIMD "Use SIMD vectorization" ON) -OPTION(CONTAINERIZE "Compiling for use in a Docker/Singularity container" OFF) OPTION(BUILD_SHARED_LIBS "Build using shared libraries" ON) INCLUDE(${CMAKE_MODULE_PATH}/SetParallelizationLibrary.cmake) INCLUDE(${CMAKE_MODULE_PATH}/SetUpNetCDF.cmake) -INCLUDE(${CMAKE_MODULE_PATH}/SetMKL.cmake) +IF (COMPILER_OPTIONS STREQUAL "Intel") + INCLUDE(${CMAKE_MODULE_PATH}/SetMKL.cmake) +ENDIF () # This INCLUDE statement executes code that sets the compile flags for DEBUG, # RELEASE, PROFILING, and TESTING. INCLUDE(${CMAKE_MODULE_PATH}/SetFortranFlags.cmake) INCLUDE_DIRECTORIES($ENV{NETCDF_FORTRAN_HOME}/include;$ENV{NETCDF_HOME}/include) - -# There is an error in CMAKE with this flag for pgf90. Unset it -GET_FILENAME_COMPONENT(FCNAME ${CMAKE_Fortran_COMPILER} NAME) -IF(FCNAME STREQUAL "pgf90") - UNSET(CMAKE_SHARED_LIBRARY_LINK_Fortran_FLAGS) -ENDIF(FCNAME STREQUAL "pgf90") - ############################################################ # Define the actual files and folders that make up the build ############################################################ @@ -62,6 +66,39 @@ SET(SRC ${CMAKE_SOURCE_DIR}/src) SET(LIB ${CMAKE_SOURCE_DIR}/lib) SET(BIN ${CMAKE_SOURCE_DIR}/bin) SET(MOD ${CMAKE_SOURCE_DIR}/include) +SET(TEST ${CMAKE_SOURCE_DIR}/test) +SET(PY ${CMAKE_SOURCE_DIR}/python/swiftest) + +FUNCTION(REPLACE_VERSION IN_FILE LANGUAGE) + # Make list of strings from file + FILE(STRINGS ${IN_FILE} LINES) + + # Replace file with new one + FILE(WRITE ${IN_FILE} "") + + # Look for the word "VERSION" and replace the line with the updated one + FOREACH(LINE IN LISTS LINES) + IF (LANGUAGE STREQUAL "Fortran") # This is the version found in the swiftest driver program + STRING(FIND "${LINE}" " VERSION =" LINE_HAS_VER) + IF (LINE_HAS_VER GREATER_EQUAL 0) # This is the version line + FILE(APPEND ${IN_FILE} " character(*), parameter :: VERSION = \"${CMAKE_PROJECT_VERSION}\" !! Swiftest version\n") + ELSE () + FILE(APPEND ${IN_FILE} "${LINE}\n") # No match. Put this line back like we found it + ENDIF () + ELSEIF (LANGUAGE STREQUAL "Python") # This is the version found in the Python package + STRING(FIND "${LINE}" "version=" LINE_HAS_VER) + IF (LINE_HAS_VER GREATER_EQUAL 0) # This is the version line + FILE(APPEND ${IN_FILE} " version='${CMAKE_PROJECT_VERSION}',\n") + ELSE () + FILE(APPEND ${IN_FILE} "${LINE}\n") # No match. Put this line back like we found it + ENDIF () + ENDIF () + ENDFOREACH () +ENDFUNCTION () + +REPLACE_VERSION(${SRC}/globals/globals_module.f90 "Fortran" ) +REPLACE_VERSION(${PY}/setup.py "Python") + # Have the .mod files placed in the lib folder SET(CMAKE_Fortran_MODULE_DIRECTORY ${MOD}) @@ -69,21 +106,11 @@ SET(CMAKE_Fortran_MODULE_DIRECTORY ${MOD}) # The source for the SWIFTEST binary and have it placed in the bin folder ADD_SUBDIRECTORY(${SRC} ${BIN}) +# Set up test directory +ENABLE_TESTING() +ADD_SUBDIRECTORY(${TEST}) + # Add a distclean target to the Makefile ADD_CUSTOM_TARGET(distclean COMMAND ${CMAKE_COMMAND} -P ${CMAKE_SOURCE_DIR}/distclean.cmake ) - -# Add an automatic test -ADD_CUSTOM_TARGET(Run_Test ALL - COMMAND echo "***Running Python script examples/Basic_Simulation/initial_conditions.py***" - COMMAND python ${CMAKE_SOURCE_DIR}/examples/Basic_Simulation/initial_conditions.py - COMMAND echo "***Saving data to directory examples/Basic_Simulation/simdata***" - COMMAND echo "***Running Python script examples/Basic_Simulation/output_reader.py***" - COMMAND python ${CMAKE_SOURCE_DIR}/examples/Basic_Simulation/output_reader.py - COMMAND echo "***Plotting results as examples/Basic_Simulation/output.eps***" - COMMAND echo "***Calculating errors with Python script examples/Basic_Simulation/errors.py***" - COMMAND python ${CMAKE_SOURCE_DIR}/examples/Basic_Simulation/errors.py - COMMAND echo "***Test Complete***" -) - diff --git a/Dockerfile b/Dockerfile index ce3e615c3..2c2b280ad 100644 --- a/Dockerfile +++ b/Dockerfile @@ -91,9 +91,11 @@ ENV F77="${FC}" # https://www.intel.com/content/www/us/en/docs/fortran-compiler/developer-guide-reference/2023-1/x-qx.html ARG MACHINE_CODE_VALUE="sse2" -# Build type options are DEBUG, RELEASE, PROFILE, or TESTING. ARG BUILD_TYPE="RELEASE" +# Additional CMAKE options: +ARG EXTRA_CMAKE_OPTIONS="" + # Swiftest ENV NETCDF_HOME=${INSTALL_DIR} ENV NETCDF_FORTRAN_HOME=${NETCDF_HOME} @@ -109,10 +111,11 @@ COPY ./src/ /swiftest/src/ COPY ./CMakeLists.txt /swiftest/ RUN cd swiftest && \ cmake -S . -B build -DCMAKE_INSTALL_PREFIX="${INSTALL_DIR}" \ - -DMACHINE_CODE_VALUE=${MACHINE_CODE} \ + -DMACHINE_CODE_VALUE=${MACHINE_CODE_VALUE} \ -DCMAKE_BUILD_TYPE=${BUILD_TYPE} \ -DUSE_COARRAY=OFF \ - -DBUILD_SHARED_LIBS=OFF && \ + -DBUILD_SHARED_LIBS=OFF \ + ${EXTRA_CMAKE_OPTIONS} && \ cmake --build build && \ cmake --install build @@ -133,20 +136,18 @@ ENV SHELL="/bin/bash" ENV PATH="/opt/conda/bin:${PATH}" ENV LD_LIBRARY_PATH="/usr/local/lib" +COPY --from=build_driver /usr/local/bin/swiftest_driver /opt/conda/bin/swiftest_driver +COPY ./python/. /opt/conda/pkgs/ +COPY environment.yml . + RUN conda update --all -y && \ conda install conda-libmamba-solver -y && \ - conda config --set solver libmamba - -COPY environment.yml . -RUN conda env create -f environment.yml && \ + conda config --set solver libmamba && \ + conda install conda-build -y && \ + conda env create -f environment.yml && \ + cd /opt/conda/pkgs/swiftest && conda develop . --name swiftest-env && \ conda init bash && \ - echo "conda activate swiftest-env" >> ~/.bashrc - -COPY ./python/. /opt/conda/pkgs/ -COPY --from=build_driver /usr/local/bin/swiftest_driver /opt/conda/bin/swiftest_driver - -# Start new shell to activate the environment and install Swiftest -RUN cd /opt/conda/pkgs/swiftest && conda develop . && \ + echo "conda activate swiftest-env" >> ~/.bashrc && \ conda clean --all -y && \ mkdir -p /.astropy && \ chmod -R 777 /.astropy && \ diff --git a/README.md b/README.md index 2de51325a..59a33b76b 100644 --- a/README.md +++ b/README.md @@ -63,55 +63,74 @@ $ git pull You now have a Swiftest repository on your personal machine that you may compile, edit, and run as you see fit. -**Compiling Swiftest** +**Compiling the Swiftest driver program** -Swiftest is written in modern Fortran and must be compiled before it can be run. After compilation, an executable, called the Swiftest driver, will have been created in the ```/swiftest/bin/``` directory. - -Swiftest is compiled through [CMake](https://cmake.org/). Compiling with CMake has a number of benefits that provide a streamlined experience for the Swiftest user and developer. At compilation, CMake will automatically select the set of flags that are compatible with the local compiler. CMake also allows a Swiftest developer to re-compile only the files that have been edited, instead of requiring the developer to re-compile the entire Swiftest program. Please visit the CMake website for more information on how to install CMake. - -Once CMake is installed, navigate to the topmost directory in your Swiftest repository. It is best practice to create a ```build``` directory in your topmost directory from which you will compile Swiftest. This way, temporary CMake files will not clutter up the ```swiftest/src/``` sub-directories. To create a new directory and then navigate into that directory, type the following: +***Compiling `swiftest_driver` using Docker*** +By far the simplest, most reliable way of compiling the driver program is via a Docker container. The Swiftest project contains a Dockerfile that may be used to generate an executable without needing to provide any external dependencies, other than the Docker engine itself (see [here](https://docs.docker.com/get-docker/) for instructions on obtaining Docker). Once Docker is installed and the Docker engine is running, execute: ``` -$ mkdir build -$ cd build +$ docker build --target=export_driver \ + --output=bin \ + --build-arg MACHINE_CODE_VALUE="Host" \ + [ --build-arg BUILD_TYPE="*RELEASE*|DEBUG|TESTING|PROFILE" ] \ + [ --build-arg EXTRA_CMAKE_OPTIONS="-D" ] ``` +The Docker build will download and compile all of the library dependencies (HDF5, NetCDF-C, and NetCDF-Fortran) as static libraries and the Swiftest driver using Intel compilers. Once completed, the Swiftest executable, called ```swiftest_driver```, should now be created in the ```bin/``` directory. + +> Note: The Dockerfile is designed to build an executable that is compatible with a broad range of CPU architectures by specifying the SSE2 instruction as a target for SIMD instructions using the `-x` compiler option. When compiling on the same CPU archictecture you plan to execute the driver program, for the highest possible SIMD performance, use `--build-arg MACHINE_CODE_VALUE="Host" to override the default `MACHINE_CODE_VALUE="SSE2"`. For additional options see [here](https://www.intel.com/content/www/us/en/docs/fortran-compiler/developer-guide-reference/2023-1/x-qx.html). + +The optional Docker argument `EXTRA_CMAKE_OPTIONS` is provided to pass any additional CMake arguments (see below). + + +***Compiling `swiftest_driver` using CMake*** + +The Swiftest driver program is written in modern Fortran and must be compiled before it can be run. After compilation, an executable, called the `swiftest_driver``, will have been created in the ```bin/``` directory. + +Swiftest is compiled through [CMake](https://cmake.org/). Compiling with CMake has a number of benefits that provide a streamlined experience for the Swiftest user and developer. At compilation, CMake will automatically select the set of flags that are compatible with the local compiler. CMake also allows a Swiftest developer to re-compile only the files that have been edited, instead of requiring the developer to re-compile the entire Swiftest program. Please visit the CMake website for more information on how to install CMake. + As mentioned in the **System Requirements** section, Swiftest requires the NetCDF and NetCDF Fortran libraries to be installed prior to compilation. If the libraries are installed in the standard library location on your machine, CMake should be able to find the libraries without specifying the path. However, if CMake struggles to find the NetCDF libraries, there are two ways to set the path to these libraries. 1. Create an environment variable called ```NETCDF_FORTRAN_HOME``` that contains the path to the location where the libraries are installed -2. Set the path at the time of compilation using ```-DCMAKE_PREFIX_PATH=/path/to/netcdf/``` +2. Set the path at the build step using ```-CMAKE_PREFIX_PATH=/path/to/netcdf/``` CMake allows the user to specify a set of compiler flags to use during compilation. We define five sets of compiler flags: release, testing, profile, math, and debug. To view and/or edit the flags included in each set, see ```swiftest/cmake/Modules/SetFortranFlags.cmake```. As a general rule, the release flags are fully optimized and best used when running Swiftest with the goal of generating results. This is the default set of flags. When making changes to the Swiftest source code, it is best to compile Swiftest using the debug set of flags. You may also define your own set of compiler flags. -To build Swiftest with the release flags (default) using the Intel fortran compiler (ifort), type the following: +Navigate to the topmost directory in your Swiftest repository. It is best practice to create a ```build``` directory in your topmost directory from which you will compile Swiftest. This way, temporary CMake files will not clutter up the ```swiftest/src/``` sub-directories. The commands to build the source code into a ```build``` directory and compile Swiftest are: + ``` -$ cmake .. +$ cmake -B build -S . +$ cmake --build build ``` -To build with the debug flags, type: +The [CMake Fortran template](https://github.com/SethMMorton/cmake_fortran_template) comes with a script that can be used to clean out any build artifacts and start from scratch: + ``` -$ cmake .. -DCMAKE_BUILD_TYPE=DEBUG +$ cmake -P distclean.cmake ``` -To build with another set of flags, simply replace ```DEBUG``` in the above line with the name of the flags you wish to use. -Add ```-CMAKE_PREFIX_PATH=/path/to/netcdf/``` to these commands as needed. +The Swiftest CMake configuration comes with several customization options: -If using the GCC fortran compiler (gfortran), add the following flags: -``` --DCMAKE_Fortran_FLAGS="-I/usr/lib64/gfortran/modules/ -ffree-line-length-512" -``` -You can manually specify the compiler you wish to use with the following flag: -``` -c-DCMAKE_Fortran_COMPILER=$(which ifort) -``` +| Option | CMake command | +| --------------------------------|------------------------------------------------------------| +| Build type | \-DCMAKE_BUILD_TYPE=[**RELEASE**\|DEBUG\|TESTING\|PROFILE] | +| Enable/Disable OpenMP support | \-DUSE_OPENMP=[**ON**\|OFF] | +| Enable/Disable SIMD directives | \-DUSE_SIMD=[**ON**\|OFF] | +| Enable/Disable Coarray support (experimental) | \-DUSE_COARRAY=[ON\|**OFF**] | +| Set Fortran compiler path | \-DCMAKE_Fortran_COMPILER=/path/to/fortran/compiler | +| Set path to make program | \-DCMAKE_MAKE_PROGRAM=/path/to/make | +| Enable/Disable shared libraries (Intel only) | \-DBUILD_SHARED_LIBS=[**ON\|OFF] | +| Add additional include path | \-DCMAKE_Fortran_FLAGS="-I/path/to/libries | -After building Swiftest, make the executable using: + +To see a list of all possible options available to CMake: ``` -$ make +$ cmake -B build -S . -LA ``` -The Swiftest executable, called ```swiftest_driver```, should now be created in the ```/swiftest/bin/``` directory. +The Swiftest executable, called ```swiftest_driver```, should now be created in the ```bin/``` directory. + **Download the `swiftest_driver` as a Docker or Singularity container.** @@ -493,7 +512,7 @@ This example walks through how to set up a standard solar system simulation. It - Semi-Interacting Massive Bodies - Gravitationally affect and are affected by fully-interacting massive bodies, do not gravitationally affect and are not affected by other semi-interacting massive bodies. - Test Particles - Gravitationally affected by fully-interacting massive bodies only. -To generate the initial conditions, run the Python script titled **initial_conditions.py**. This script also runs Swiftest SyMBA, generating output. To process the output file, run the script titled **output_reader.py**. +To generate the initial conditions, run the Python script titled **basic_simulation.py**. This script also runs Swiftest SyMBA, generating output. To process the output file, run the script titled **output_reader.py**. **Chambers2013** diff --git a/cmake/Modules/FindCoarray_Fortran.cmake b/cmake/Modules/FindCoarray_Fortran.cmake index 314371326..93f73f466 100644 --- a/cmake/Modules/FindCoarray_Fortran.cmake +++ b/cmake/Modules/FindCoarray_Fortran.cmake @@ -27,27 +27,40 @@ INCLUDE (${CMAKE_ROOT}/Modules/FindPackageHandleStandardArgs.cmake) STRING(TOUPPER "${CMAKE_BUILD_TYPE}" BT) IF(BT STREQUAL "DEBUG") - SET (Coarray_Fortran_FLAG_CANDIDATES - #Intel - "-coarray=single" - #Intel windows - "/Qcoarray:single" - #Gnu - "-fcoarray=single" - #Empty, if compiler automatically accepts coarray - " " + IF (COMPILER_OPTIONS STREQUAL "Intel") + SET (Coarray_Fortran_FLAG_CANDIDATES + #Intel + "-coarray=single" + #Intel windows + "/Qcoarray:single" + #Empty, if compiler automatically accepts coarray + " " + ) + ELSEIF (COMPILER_OPTIONS STREQUAL "GNU") + SET (Coarray_Fortran_FLAG_CANDIDATES + #Gnu + "-fcoarray=single" + #Empty, if compiler automatically accepts coarray + " " ) + ENDIF() ELSE() - SET (Coarray_Fortran_FLAG_CANDIDATES - #Intel - "-coarray=distributed" - #Intel windows - "/Qcoarray:distributed" - #Gnu - "-fcoarray=lib -lcaf_mpi" - #Empty, if compiler automatically accepts coarray - " " - ) + IF (COMPILER_OPTIONS STREQUAL "Intel") + SET (Coarray_Fortran_FLAG_CANDIDATES + #Intel + "-coarray=distributed" + #Intel windows + "/Qcoarray:distributed" + #Empty, if compiler automatically accepts coarray + " " + ) + ELSEIF (COMPILER_OPTIONS STREQUAL "GNU") + SET (Coarray_Fortran_FLAG_CANDIDATES + #Gnu + "-fcoarray=lib -lcaf_mpi" + #Empty, if compiler automatically accepts coarray + " " + ) ENDIF() diff --git a/cmake/Modules/FindOpenMP_Fortran.cmake b/cmake/Modules/FindOpenMP_Fortran.cmake index 06d679e7c..d3a0bc29b 100644 --- a/cmake/Modules/FindOpenMP_Fortran.cmake +++ b/cmake/Modules/FindOpenMP_Fortran.cmake @@ -25,19 +25,30 @@ INCLUDE (${CMAKE_ROOT}/Modules/FindPackageHandleStandardArgs.cmake) -IF (USE_SIMD) - SET (OpenMP_Fortran_FLAG_CANDIDATES - "-qopenmp" # Intel - "/Qopenmp" # Intel Windows - "-fopenmp" # GNU - ) -ELSE () - SET (OpenMP_Fortran_FLAG_CANDIDATES - "-qopenmp -qno-openmp-simd" # Intel - "/Qopenmp-simd-" # Intel Windows - "-fopenmp" # GNU - ) -ENDIF (USE_SIMD) +IF (COMPILER_OPTIONS STREQUAL "Intel") + IF (USE_SIMD) + SET (OpenMP_Fortran_FLAG_CANDIDATES + "-qopenmp" # Intel + "/Qopenmp" # Intel Windows + ) + ELSE () + SET (OpenMP_Fortran_FLAG_CANDIDATES + "-qopenmp -qno-openmp-simd" # Intel + "/Qopenmp-simd-" # Intel Windows + ) + ENDIF (USE_SIMD) +ELSEIF (COMPILER_OPTIONS STREQUAL "GNU") + IF (USE_SIMD) + SET (OpenMP_Fortran_FLAG_CANDIDATES + "-fopenmp" # GNU + ) + ELSE () + SET (OpenMP_Fortran_FLAG_CANDIDATES + "-fopenmp -fno-openmp-simd" # GNU + ) + ENDIF (USE_SIMD) + +ENDIF () IF (DEFINED OpenMP_Fortran_FLAGS) SET (OpenMP_Fortran_FLAG_CANDIDATES) diff --git a/cmake/Modules/SetFortranFlags.cmake b/cmake/Modules/SetFortranFlags.cmake index 4c8cc9b85..a7c46e10d 100644 --- a/cmake/Modules/SetFortranFlags.cmake +++ b/cmake/Modules/SetFortranFlags.cmake @@ -46,15 +46,6 @@ ELSE() MESSAGE(FATAL_ERROR "CMAKE_BUILD_TYPE not valid! ${BUILD_TYPE_MSG}") ENDIF(BT STREQUAL "RELEASE") -IF (CMAKE_Fortran_COMPILER_ID STREQUAL "GNU") - IF (APPLE) - SET(MACHINE_CODE_VALUE "tune=native" CACHE STRING "Tells the compiler which processor features it may target, including which instruction sets and optimizations it may generate.") - ELSE () - SET(MACHINE_CODE_VALUE "arch=native" CACHE STRING "Tells the compiler which processor features it may target, including which instruction sets and optimizations it may generate.") - ENDIF () -ELSE () - SET(MACHINE_CODE_VALUE "host" CACHE STRING "Tells the compiler which processor features it may target, including which instruction sets and optimizations it may generate.") -ENDIF () ######################################################### # If the compiler flags have already been set, return now @@ -77,338 +68,456 @@ ENDIF(CMAKE_Fortran_FLAGS_RELEASE AND CMAKE_Fortran_FLAGS_TESTING AND CMAKE_Fort ##################### # Free form -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-ffree-form" # GNU - ) - -# Don't add underscores in symbols for C-compatability -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-fno-underscoring" # GNU - ) -# Compile code assuming that IEEE signaling NaNs may generate user-visible traps during floating-point operations. -# Setting this option disables optimizations that may change the number of exceptions visible with signaling NaNs. -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-fsignaling-nans " # GNU - ) - -# Allows for lines longer than 80 characters without truncation -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-ffree-line-length-none" # GNU (gfortran) - ) +IF (COMPILER_OPTIONS STREQUAL "GNU") + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-ffree-form" # GNU + ) -# Disables right margin wrapping in list-directed output -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-no-wrap-margin" # Intel - "/wrap-margin-" # Intel Windows - ) + # Don't add underscores in symbols for C-compatability + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-fno-underscoring" # GNU + ) + # Compile code assuming that IEEE signaling NaNs may generate user-visible traps during floating-point operations. + # Setting this option disables optimizations that may change the number of exceptions visible with signaling NaNs. + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-fsignaling-nans " # GNU + ) + + # Allows for lines longer than 80 characters without truncation + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-ffree-line-length-512" # GNU (gfortran) + ) +ELSEIF (COMPILER_OPTIONS STREQUAL "Intel") + # Disables right margin wrapping in list-directed output + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-no-wrap-margin" # Intel + "/wrap-margin-" # Intel Windows + ) -# Aligns a variable to a specified boundary and offset -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-align all -align array64byte" # Intel - "/align:all /align:array64byte" # Intel Windows - ) + # Aligns a variable to a specified boundary and offset + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-align all -align array64byte" # Intel + "/align:all /align:array64byte" # Intel Windows + ) -# Enables changing the variable and array memory layout -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-pad" # Intel - "/Qpad" # Intel Windows + # Enables changing the variable and array memory layout + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-pad" # Intel + "/Qpad" # Intel Windows - ) + ) +ENDIF () IF (NOT BUILD_SHARED_LIBS) - # Use static Intel libraries - SET_COMPILE_FLAG(CMAKE_Fortran_LINK_FLAGS "${CMAKE_Fortran_LINK_FLAGS}" - Fortran "-static-intel" # Intel - ) - # Use static Intel MPI libraries - SET_COMPILE_FLAG(CMAKE_Fortran_LINK_FLAGS "${CMAKE_Fortran_LINK_FLAGS}" - Fortran "-static_mpi" # Intel - ) - - IF (USE_OPENMP) + + IF (COMPILER_OPTIONS STREQUAL "Intel") + # Use static Intel libraries SET_COMPILE_FLAG(CMAKE_Fortran_LINK_FLAGS "${CMAKE_Fortran_LINK_FLAGS}" - Fortran "-qopenmp-link=static" # Intel + Fortran "-static-intel" # Intel + ) + # Use static Intel MPI libraries + SET_COMPILE_FLAG(CMAKE_Fortran_LINK_FLAGS "${CMAKE_Fortran_LINK_FLAGS}" + Fortran "-static_mpi" # Intel ) - ENDIF (USE_OPENMP) -ENDIF () + + IF (USE_OPENMP) + SET_COMPILE_FLAG(CMAKE_Fortran_LINK_FLAGS "${CMAKE_Fortran_LINK_FLAGS}" + Fortran "-qopenmp-link=static" # Intel + ) + ENDIF (USE_OPENMP) + ELSEIF (COMPIELR_OPTIONS STREQUAL "GNU") + # Set GNU static libraries + SET_COMPILE_FLAG(CMAKE_Fortran_LINK_FLAGS "${CMAKE_Fortran_LINK_FLAGS}" + Fortran "-static-libgfortran" + ) + SET_COMPILE_FLAG(CMAKE_Fortran_LINK_FLAGS "${CMAKE_Fortran_LINK_FLAGS}" + Fortran "-static-libgcc" + ) + SET_COMPILE_FLAG(CMAKE_Fortran_LINK_FLAGS "${CMAKE_Fortran_LINK_FLAGS}" + Fortran "-static-libstdc++" + ) + SET_COMPILE_FLAG(CMAKE_Fortran_LINK_FLAGS "${CMAKE_Fortran_LINK_FLAGS}" + Fortran "-static-libquadmath" + ) + ENDIF () + +ENDIF (NOT BUILD_SHARED_LIBS) IF (USE_SIMD) - # Enables OpenMP SIMD compilation when OpenMP parallelization is disabled. - IF (NOT USE_OPENMP) + SET(MACHINE_CODE_VALUE "host" CACHE STRING "Tells the compiler which processor features it may target, including which instruction sets and optimizations it may generate.") + + IF (COMPILER_OPTIONS STREQUAL "Intel") + + # Enables OpenMP SIMD compilation when OpenMP parallelization is disabled. + IF (NOT USE_OPENMP) + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-qno-openmp -qopenmp-simd" # Intel + Fortran "/Qopenmp- /Qopenmp-simd" # Intel Windows + ) + ENDIF (NOT USE_OPENMP) + + # Optimize for an old enough processor that it should run on most computers + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-x${MACHINE_CODE_VALUE}" # Intel + "/Qx${MACHINE_CODE_VALUE}" # Intel Windows + ) + + # Generate an extended set of vector functions SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-qno-openmp -qopenmp-simd" # Intel - Fortran "/Qopenmp- /Qopenmp-simd" # Intel Windows + Fortran "-vecabi=cmdtarget" # Intel + "/Qvecabi:cmdtarget" # Intel Windows + ) + ELSEIF (COMPILER_OPTIONS STREQUAL "GNU") + + # Enables OpenMP SIMD compilation when OpenMP parallelization is disabled. + IF (NOT USE_OPENMP) + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-fno-openmp -fopenmp-simd" # GNU ) - ENDIF (NOT USE_OPENMP) + ENDIF (NOT USE_OPENMP) - # Optimize for an old enough processor that it should run on most computers - SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-x${MACHINE_CODE_VALUE}" # Intel - "/Qx${MACHINE_CODE_VALUE}" # Intel Windows - "-m${MACHINE_CODE_VALUE}" # GNU - ) + IF (MACHINE_CODE_VALUE STREQUAL "host") + SET(MACHINE_CODE_VALUE "native" CACHE STRING "Tells the compiler which processor features it may target, including which instruction sets and optimizations it may generate.") + ENDIF () - # Generate an extended set of vector functions - SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" - Fortran "-vecabi=cmdtarget" # Intel - Fortran "/Qvecabi:cmdtarget" # Intel Windows + IF (APPLE) + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-mtune=${MACHINE_CODE_VALUE}" ) -ENDIF (USE_SIMD) - + ELSE () + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS "${CMAKE_Fortran_FLAGS}" + Fortran "-march=${MACHINE_CODE_VALUE}" + ) + ENDIF () + ENDIF () +ENDIF (USE_SIMD) ################### ### DEBUG FLAGS ### ################### - -# Disable optimizations -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran REQUIRED "-O0" # All compilers not on Windows - "/Od" # Intel Windows - "-Og" # GNU (gfortran) - ) - -# Turn on all warnings -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-warn all" # Intel - "/warn:all" # Intel Windows - "-Wall" # GNU - ) -# This enables some extra warning flags that are not enabled by -Wall -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-Wextra" # GNU - ) - -# Disable the warning that arrays may be uninitialized, which comes up due to a known bug in gfortran -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-Wno-maybe-uninitialized" # GNU - ) -# Disable the warning about unused dummy arguments. These primarily occur due to interface rules for type-bound procedures used in extendable types. -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-Wno-unused-dummy-argument" # GNU - ) - -# Tells the compiler to issue compile-time messages for nonstandard language elements (Fortran 2018). -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-stand f18" # Intel - "/stand:f18" # Intel Windows - "-fstd=f2018" # GNU +IF (CMAKE_BUILD_TYPE STREQUAL "DEBUG" OR CMAKE_BUILD_TYPE STREQUAL "TESTING" ) + # Disable optimizations + IF (COMPILER_OPTIONS STREQUAL "Intel") + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran REQUIRED "-O0" # All compilers not on Windows + "/Od" # Intel Windows + ) + + # Turn on all warnings + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-warn all" # Intel + "/warn:all" # Intel Windows + ) + + # Tells the compiler to issue compile-time messages for nonstandard language elements (Fortran 2018). + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-stand f18" # Intel + "/stand:f18" # Intel Windows + ) + + # Traceback + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-traceback" # Intel Group + "/traceback" # Intel Windows + ) + + # Check everything + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-check all" # Intel + "/check:all" # Intel Windows + ) + + # Initializes matrices/arrays with NaN values + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-init=snan,arrays" # Intel + "/Qinit:snan,arrays" # Intel Windows + ) + + # Does not generate an interface block for each routine in a source file + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-nogen-interfaces" # Intel + "/nogen-interfaces" # Intel Windows + ) + + # Does not generate aposition independent executable + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-no-pie" # Intel + ) + + # Does not set denormal results from floating-point calculations to zero + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-no-ftz" # Intel + "/Qftz-" # Intel Windows + ) + + # Enables floating-point invalid, divide-by-zero, and overflow exceptions + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-fpe-all=0" # Intel + "/fpe-all:0" # Intel Windows + ) + + # Enables floating-point invalid, divide-by-zero, and overflow exceptions + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-fpe0" # Intel + "/fpe:0" # Intel Windows + ) + + # Enables debug info + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-debug all" # Intel + "/debug:all" # Intel Windows + ) + + # Disables additional interprocedural optimizations for a single file compilation + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-no-ip" # Intel + "/Qip-" # Intel Windows + ) + + # Disables prefetch insertion optimization + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-qno-opt-prefetch" # Intel + "/Qopt-prefetch-" # Intel Windows + ) + + ELSEIF (COMPILER_OPTIONS STREQUAL "GNU") + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran REQUIRED "-Og" # GNU (gfortran) + ) + + # Turn on all warnings + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-Wall" # GNU + ) + # This enables some extra warning flags that are not enabled by -Wall + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-Wextra" # GNU + ) + + # Disable the warning that arrays may be uninitialized, which comes up due to a known bug in gfortran + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-Wno-maybe-uninitialized" # GNU + ) + # Disable the warning about unused dummy arguments. These primarily occur due to interface rules for type-bound procedures used in extendable types. + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-Wno-unused-dummy-argument" # GNU + ) + + # Tells the compiler to issue compile-time messages for nonstandard language elements (Fortran 2018). + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-fstd=f2018" # GNU ) -# Traceback -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-traceback" # Intel Group - "/traceback" # Intel Windows - "-fbacktrace" # GNU (gfortran) + # Traceback + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-fbacktrace" # GNU (gfortran) ) -# Sanitize -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-fsanitize=address, undefined" # Gnu - ) - -# Check everything -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-check all" # Intel - "/check:all" # Intel Windows - "-fcheck=all" # GNU - ) -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-fstack-check" # GNU - ) - - -# Initializes matrices/arrays with NaN values -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-init=snan,arrays" # Intel - "/Qinit:snan,arrays" # Intel Windows - "-finit-real=snan" # GNU + # Sanitize + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-fsanitize=address, undefined" # Gnu ) -# Does not generate an interface block for each routine in a source file -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-nogen-interfaces" # Intel - "/nogen-interfaces" # Intel Windows + # Check everything + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-fcheck=all" # GNU ) - -# Does not generate aposition independent executable -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-no-pie" # Intel - ) - -# Does not set denormal results from floating-point calculations to zero -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-no-ftz" # Intel - "/Qftz-" # Intel Windows + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-fstack-check" # GNU ) -# Enables floating-point invalid, divide-by-zero, and overflow exceptions -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-fpe-all=0" # Intel - "/fpe-all:0" # Intel Windows - "-ffpe-trap=zero,overflow,underflow" # GNU + # Initializes matrices/arrays with NaN values + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-finit-real=snan" # GNU ) -# List of floating-point exceptions, whose flag status is printed to ERROR_UNIT when invoking STOP and ERROR STOP -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-ffpe-summary=all" # GNU + # Generates non position-independent code + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-fno-PIE" # GNU ) -# Enables floating-point invalid, divide-by-zero, and overflow exceptions -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-fpe0" # Intel - "/fpe:0" # Intel Windows + # Enables floating-point invalid, divide-by-zero, and overflow exceptions + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-ffpe-trap=zero,overflow,underflow" # GNU ) -# Enables debug info -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-debug all" # Intel - "/debug:all" # Intel Windows + # List of floating-point exceptions, whose flag status is printed to ERROR_UNIT when invoking STOP and ERROR STOP + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-ffpe-summary=all" # GNU ) -# Disables additional interprocedural optimizations for a single file compilation -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-no-ip" # Intel - "/Qip-" # Intel Windows + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran "-fstack-check" # GNU ) + ENDIF () -# Disables prefetch insertion optimization -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-qno-opt-prefetch" # Intel - "/Qopt-prefetch-" # Intel Windows - ) - -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran "-fstack-check" # GNU - ) +ENDIF () ##################### ### TESTING FLAGS ### ##################### -# Optimizations -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_TESTING "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran REQUIRED "-O3" # All compilers not on Windows - "/O3" # Intel Windows +IF (CMAKE_BUILD_TYPE STREQUAL "TESTING" ) + + # Optimizations + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_TESTING "${CMAKE_Fortran_FLAGS_DEBUG}" + Fortran REQUIRED "-O3" # All compilers not on Windows + "/O3" # Intel Windows ) +ENDIF () ##################### ### RELEASE FLAGS ### ##################### # NOTE: agressive optimizations (-O3) are already turned on by default -# Unroll loops -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-unroll" # Intel - "/unroll" # Intel Windows - "-funroll-loops" # GNU - ) -# Inline functions -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-inline" # Intel - "/Qinline" # Intel Windows - "-finline-functions" # GNU - ) - -# Calls the Matrix Multiply library -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-qopt-matmul" # Intel - "/Qopt-matmul" # Intel Windows - ) - -# Saves the compiler options and version number to the executable -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-sox" # Intel - ) - -# Aligns a variable to a specified boundary and offset -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-align all" # Intel - "/align:all" # Intel Windows - ) - -# No floating-point exceptions -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-fp-model no-except" # Intel - "/fp:no-except" # Intel Windows - ) - -# Generate fused multiply-add instructions - SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-fma" # Intel - "/Qfma" # Intel Windows - ) - -# Tells the compiler to link to certain libraries in the Intel oneAPI Math Kernel Library (oneMKL). - SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-qmkl=cluster" # Intel - "-qmkl" # Intel - "/Qmkl:cluster" # Intel Windows - "/Qmkl" # Intel Windows - ) - -# Enables additional interprocedural optimizations for a single file compilation -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-ip" # Intel - "/Qip" # Intel Windows - ) +IF (CMAKE_BUILD_TYPE STREQUAL "RELEASE" OR CMAKE_BUILD_TYPE STREQUAL "PROFILE") + + IF (COMPILER_OPTIONS STREQUAL "Intel") + # Unroll loops + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-unroll" # Intel + "/unroll" # Intel Windows + + ) + + # Inline functions + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-inline" # Intel + "/Qinline" # Intel Windows + ) + + # Calls the Matrix Multiply library + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-qopt-matmul" # Intel + "/Qopt-matmul" # Intel Windows + ) + + # Aligns a variable to a specified boundary and offset + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-align all" # Intel + "/align:all" # Intel Windows + ) + + # No floating-point exceptions + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-fp-model no-except" # Intel + "/fp:no-except" # Intel Windows + ) + + # Generate fused multiply-add instructions + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-fma" # Intel + "/Qfma" # Intel Windows + ) + + # Tells the compiler to link to certain libraries in the Intel oneAPI Math Kernel Library (oneMKL). + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-qmkl=cluster" # Intel + "-qmkl" # Intel + "/Qmkl:cluster" # Intel Windows + "/Qmkl" # Intel Windows + ) + + # Enables additional interprocedural optimizations for a single file compilation + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-ip" # Intel + "/Qip" # Intel Windows + ) + + ELSEIF(COMPILER_OPTIONS STREQUAL "GNU") + # Unroll loops + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-funroll-loops" # GNU + ) + + # Inline functions + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_RELEASE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-finline-functions" # GNU + ) + ENDIF () +ENDIF () ##################### ### MATH FLAGS ### ##################### -# Some subroutines require more strict floating point operation optimizations for repeatability -SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" - Fortran "-fp-model=precise" # Intel - "/fp:precise" # Intel Windows - "-fno-unsafe-math-optimizations" # GNU - ) -# Disable transformations and optimizations that assume default floating-point rounding behavior. -SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" - Fortran "-frounding-math" - ) +IF (COMPILER_OPTIONS STREQUAL "Intel") -SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" - Fortran "-prec-div" # Intel - "/Qprec-div" # Intel Windows - ) + # Some subroutines require more strict floating point operation optimizations for repeatability + SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" + Fortran "-fp-model=precise" # Intel + "/fp:precise" # Intel Windows + ) -SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" - Fortran "-prec-sqrt" # Intel - "/Qprec-sqrt" # Intel Windows - ) + SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" + Fortran "-prec-div" # Intel + "/Qprec-div" # Intel Windows + ) -SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" - Fortran "-assume protect-parens" # Intel - "/assume:protect-parens" # Intel Windows - ) - -# Improves floating-point precision and consistency -SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" - Fortran "-mp1" # Intel - "/Qprec" # Intel Windows - ) - -# Most subroutines can use aggressive optimization of floating point operations without problems. -SET_COMPILE_FLAG(FASTMATH_FLAGS "${FASTMATH_FLAGS}" - Fortran "-fp-model=fast" # Intel - "/fp:fast" # Intel Windows - "-ffast-math" # GNU - ) + SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" + Fortran "-prec-sqrt" # Intel + "/Qprec-sqrt" # Intel Windows + ) + SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" + Fortran "-assume protect-parens" # Intel + "/assume:protect-parens" # Intel Windows + ) + + # Improves floating-point precision and consistency + SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" + Fortran "-mp1" # Intel + "/Qprec" # Intel Windows + ) + + # Most subroutines can use aggressive optimization of floating point operations without problems. + SET_COMPILE_FLAG(FASTMATH_FLAGS "${FASTMATH_FLAGS}" + Fortran "-fp-model=fast" # Intel + "/fp:fast" # Intel Windows + ) + + +ELSEIF (COMPILER_OPTIONS STREQUAL "GNU") + + # Some subroutines require more strict floating point operation optimizations for repeatability + SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" + Fortran "-fno-unsafe-math-optimizations" # GNU + ) + # Disable transformations and optimizations that assume default floating-point rounding behavior. + SET_COMPILE_FLAG(STRICTMATH_FLAGS "${STRICTMATH_FLAGS}" + Fortran "-frounding-math" + ) + # Most subroutines can use aggressive optimization of floating point operations without problems. + SET_COMPILE_FLAG(FASTMATH_FLAGS "${FASTMATH_FLAGS}" + Fortran "-ffast-math" # GNU + ) +ENDIF () + +# Debug mode always uses strict math SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_DEBUG "${CMAKE_Fortran_FLAGS_DEBUG}" - Fortran ${STRICTMATH_FLAGS} - ) + Fortran ${STRICTMATH_FLAGS} +) ##################### ### PROFILE FLAGS ### ##################### -# Enables the optimization reports to be generated -SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_PROFILE "${CMAKE_Fortran_FLAGS_RELEASE}" - Fortran "-O2 -pg -qopt-report=5 -traceback -p -g3" # Intel - "/O2 /Qopt-report:5 /traceback /Z7" # Intel Windows - "-O2 -pg -fbacktrace" # GNU - ) +IF (CMAKE_BUILD_TYPE STREQUAL "PROFILE") + + IF (COMPILER_OPTIONS STREQUAL "Intel") + # Enables the optimization reports to be generated + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_PROFILE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-O2 -pg -qopt-report=5 -traceback -p -g3" # Intel + "/O2 /Qopt-report:5 /traceback /Z7" # Intel Windows + ) + ELSEIF (COMPILER_OPTIONS STREQUAL "GNU") + # Enables the optimization reports to be generated + SET_COMPILE_FLAG(CMAKE_Fortran_FLAGS_PROFILE "${CMAKE_Fortran_FLAGS_RELEASE}" + Fortran "-O2 -pg -fbacktrace" # GNU + ) + ENDIF () +ENDIF () diff --git a/distclean.cmake b/distclean.cmake index 00d2e454c..a06162074 100644 --- a/distclean.cmake +++ b/distclean.cmake @@ -44,7 +44,6 @@ FILE(GLOB_RECURSE CMAKEINSTALL "${TOPDIR}/*cmake_install.cmake" FILE(GLOB_RECURSE MAKEFILE "${TOPDIR}/*Makefile") FILE(GLOB_RECURSE CMAKETESTFILES "${TOPDIR}/*CTestTestfile.cmake") SET(TOPDIRECTORIES "${TOPDIR}/lib" - "${TOPDIR}/test" "${TOPDIR}/bin" ) diff --git a/environment.yml b/environment.yml index c2bca5877..cb49c2304 100644 --- a/environment.yml +++ b/environment.yml @@ -19,5 +19,4 @@ dependencies: - astroquery>=0.4.6 - tqdm>=4.65.0 - x264>=1!157.20191217 - - ffmpeg>=4.3.2 - - conda-build \ No newline at end of file + - ffmpeg>=4.3.2 \ No newline at end of file diff --git a/examples/Basic_Simulation/.gitignore b/examples/Basic_Simulation/.gitignore index 9e4d02aad..eb544ec84 100644 --- a/examples/Basic_Simulation/.gitignore +++ b/examples/Basic_Simulation/.gitignore @@ -1,6 +1,6 @@ * !.gitignore -!initial_conditions.py +!basic_simulation.py !output_reader.py !errors.py !README.txt diff --git a/examples/Basic_Simulation/README.txt b/examples/Basic_Simulation/README.txt index ae626a9c5..968415176 100644 --- a/examples/Basic_Simulation/README.txt +++ b/examples/Basic_Simulation/README.txt @@ -16,7 +16,7 @@ Date : June 27, 2023 Included in the Basic_Simulation example directory are the following files: - README.txt : This file - - initial_conditions.py : A Python Script that generates and runs a set of initial conditions. + - basic_simulation.py : A Python Script that generates and runs a set of initial conditions. - output_reader.py : A Python Script that processes out.nc and generates output.eps - errors.py : A Python Script that processes out.nc and reports the simulation errors to the terminal diff --git a/examples/Basic_Simulation/initial_conditions.py b/examples/Basic_Simulation/basic_simulation.py similarity index 92% rename from examples/Basic_Simulation/initial_conditions.py rename to examples/Basic_Simulation/basic_simulation.py index 83614fb07..bee49ac8b 100755 --- a/examples/Basic_Simulation/initial_conditions.py +++ b/examples/Basic_Simulation/basic_simulation.py @@ -42,7 +42,13 @@ rng = default_rng(seed=123) # Add the modern planets and the Sun using the JPL Horizons Database. -sim.add_solar_system_body(["Sun","Mercury","Venus","Earth","Mars","Jupiter","Saturn","Uranus","Neptune","Pluto"]) +sim.add_solar_system_body(["Sun","Mercury","Venus","Earth","Mars","Jupiter","Saturn","Uranus","Neptune"]) + +# Add in some main belt asteroids +sim.add_solar_system_body(name=["Ceres","Vesta","Pallas","Hygiea"],id_type="smallbody") + +# Add in some big KBOs +sim.add_solar_system_body(name=["Pluto","Eris","Sedna","Haumea","Makemake","Quaoar","Orcus","Gonggong","Salacia"]) # Add 5 user-defined massive bodies. npl = 5 diff --git a/examples/Basic_Simulation/errors.py b/examples/Basic_Simulation/errors.py deleted file mode 100644 index 5cbdfaa55..000000000 --- a/examples/Basic_Simulation/errors.py +++ /dev/null @@ -1,69 +0,0 @@ -#!/usr/bin/env python3 - -""" - Copyright 2023 - David Minton, Carlisle Wishard, Jennifer Pouplin, Jake Elliott, & Dana Singh - This file is part of Swiftest. - Swiftest is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License - as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. - Swiftest is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty - of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. - You should have received a copy of the GNU General Public License along with Swiftest. - If not, see: https://www.gnu.org/licenses. -""" - -""" -Reads and processes a Swiftest output file. - -Input ------- -data.nc : A NetCDF file containing the simulation output. - -Output ------- -None -""" - -import swiftest -import xarray as xr -import matplotlib.pyplot as plt - -# Read in the simulation output and store it as an Xarray dataset. -ds = swiftest.Simulation(read_data=True) - -# Calculate the angular momentum error -ds.data['L_tot'] = ds.data['L_orbit'] + ds.data['L_spin'] + ds.data['L_escape'] -ds.data['DL'] = ds.data['L_tot'] - ds.data['L_tot'].isel(time=0) -ds.data['L_error'] = swiftest.tool.magnitude(ds.data,'DL') / swiftest.tool.magnitude(ds.data.isel(time=0), 'L_tot') -L_final = ds.data['L_error'][-1].values - -# Calculate the energy error -ds.data['E_error'] = (ds.data['TE'] - ds.data['TE'].isel(time=0)) / ds.data['TE'].isel(time=0) -E_final = ds.data['E_error'][-1].values - -# Calculate the mass error -ds.data['GMtot'] = ds.data['Gmass'].sum(dim='name',skipna=True) + ds.data['GMescape'] -ds.data['GM_error'] = (ds.data['GMtot'] - ds.data['GMtot'].isel(time=0)) / ds.data['GMtot'].isel(time=0) -GM_final = ds.data['GM_error'][-1].values - -# Print the final errors -print("Final Angular Momentum Error: ", L_final) -print("Final Energy Error: ", E_final) -print("Final Mass Error: ", GM_final) - -# Determine if the errors are within bounds -L_limit = 1e-10 -E_limit = 1e-5 -GM_limit = 1e-14 - -lerror = 0 -if (L_final > L_limit): - lerror = 1 - print("Angular Momentum Error of", L_final, " higher than threshold value of", L_limit, ". Test failed.") -if (E_final > E_limit): - lerror = 1 - print("Energy Error of", E_final, " higher than threshold value of", E_limit, ". Test failed.") -if (GM_final > GM_limit): - lerror = 1 - print("Mass Error of", GM_final, " higher than threshold value of", GM_limit, ". Test failed.") -if (lerror == 0): - print("Test successful.") diff --git a/python/swiftest/setup.py b/python/swiftest/setup.py index 917ad64e4..3a3c2a674 100644 --- a/python/swiftest/setup.py +++ b/python/swiftest/setup.py @@ -12,6 +12,6 @@ from setuptools import setup, find_packages setup(name='swiftest', - version='0.2', + version='2023.08.00', author='David A. Minton', - packages=find_packages()) \ No newline at end of file + packages=find_packages()) diff --git a/python/swiftest/swiftest/constants.py b/python/swiftest/swiftest/constants.py index 2d3f89f7c..6fa87c888 100644 --- a/python/swiftest/swiftest/constants.py +++ b/python/swiftest/swiftest/constants.py @@ -11,6 +11,8 @@ import numpy as np import astropy.constants as const +import astropy.units as u +from astropy.coordinates import SkyCoord # Constants in SI units GC = const.G.value[()] @@ -27,4 +29,6 @@ # Solar oblatenes values: From Mecheri et al. (2004), using Corbard (b) 2002 values (Table II) J2Sun = 2.198e-7 J4Sun = -4.805e-9 +rotpoleSun = SkyCoord(ra=286.13 * u.degree, dec=63.87 * u.degree).cartesian +rotSun = (360.0 / 25.05) / JD2S * rotpoleSun diff --git a/python/swiftest/swiftest/init_cond.py b/python/swiftest/swiftest/init_cond.py index 1b140f119..978fe9808 100644 --- a/python/swiftest/swiftest/init_cond.py +++ b/python/swiftest/swiftest/init_cond.py @@ -24,90 +24,270 @@ List, Any ) + +def horizons_get_physical_properties(altid,**kwargs): + """ + Parses the raw output from JPL Horizons in order to extract physical properties of a body if they exist + + + Parameters + ---------- + altid : list of str + List of ids to use for Horizons query + **kwargs: Any + Additional keyword arguments to pass to the query method (see https://astroquery.readthedocs.io/en/latest/jplhorizons/jplhorizons.html) + + + Returns + ------- + MSun_over_Mpl : float + The ratio of MSun/M of the body + radius : float + The radius of the body in m + rot: (3) float vector + The rotation rate vector oriented toward the north pole + """ + + def get_Gmass(raw_response): + GM = [s for s in raw_response.split('\n') if 'GM' in s + and 'Rel' not in s + and 'GMT' not in s + and 'ANGMOM' not in s] + if len(GM) == 0: + return None + GM = GM[0] + if len(GM) > 1: + GM = GM.split('GM')[1] + if len(GM) > 1: + GM = GM.split('=') + if len(GM) > 1: + GM = GM[1].strip().split(' ')[0].split('+')[0] + try: + GM = float(GM) + except: + GM = None + return GM + + def get_radius(raw_response): + + radius = [s for s in raw_response.split('\n') if 'mean radius' in s.lower() or 'RAD' in s or 'Radius (km)' in s] + if len(radius) == 0: + return None + radius = radius[0] + if "Radius (km)" in radius: + radius = radius.split("Radius (km)")[1].strip(' =').split('+')[0].split()[0].strip() + elif "R_eff" in radius: # Some small bodies list an effective radius + radius = radius.split('R_eff')[1].strip(' =').split('+')[0].split()[0].strip() + elif "RAD" in radius: # Some small bodies list the radius like this + radius = radius.split('RAD')[1].strip(' =').split('+')[0].split()[0].strip() + if 'x' in radius: # Triaxial ellipsoid bodies like Haumea may have multiple dimensions which need to be averaged + radius = radius.split('x') + try: + for i,v in enumerate(radius): + radius[i] = float(v.split()[0].strip()) + radius = np.average(radius) + except: + radius = None + else: # Handles most major bodies + radius = radius.split('=')[1].strip().split('+')[0].split()[0].strip() + try: + radius = float(radius) + except: + radius = None + return radius + + def get_rotrate(raw_response): + raw_response=jpl.raw_response + rotrate = [s for s in raw_response.split('\n') if 'rot. rat' in s.lower()] + if len(rotrate) == 0: + rotrate = [s for s in raw_response.split('\n') if 'ROTPER' in s.upper()] # Try the small body version + if len(rotrate) > 0: + rotrate = rotrate[0].split('ROTPER=')[1].strip() + try: + rotrate = 2*np.pi / (float(rotrate) * 3600) + except: + rotrate = None + else: + if "Synchronous" in raw_response: # Satellites have this: + rotrate = [s for s in raw_response.split('\n') if 'Orbital period' in s][0] + rotrate = rotrate.split('Orbital period')[1].replace('~',' ').replace('d',' ').replace('=',' ').strip() + rotrate = 2*np.pi / (float(rotrate) * swiftest.JD2S) + else: + rotrate = None + else: + rotrate = rotrate[0].lower().split('rot. rat')[1].split('=')[1].strip().split(' ')[0].strip() + try: + rotrate = float(rotrate) + except: + rotrate = None + return rotrate + + def get_rotpole(jpl): + RA = jpl.ephemerides()['NPole_RA'][0] + DEC = jpl.ephemerides()['NPole_DEC'][0] + + if np.ma.is_masked(RA) or np.ma.is_masked(DEC): + return np.array([0.0,0.0,1.0]) + + rotpole = SkyCoord(ra=RA * u.degree, dec=DEC * u.degree).cartesian + return np.array([rotpole.x.value, rotpole.y.value, rotpole.z.value]) + + if type(altid) != list: + altid = [altid] + + for id in altid: + jpl,idlist,namelist = horizons_query(id=id,ephemerides_start_date='2023-07-26',verbose=False,**kwargs) + Rpl = get_radius(jpl.raw_response) + if Rpl is not None: + Rpl *= 1e3 + break + + Gmass = get_Gmass(jpl.raw_response) + if Rpl is None or Gmass is None: + rot = np.full(3,np.nan) + else: + print(f"Physical properties found for {namelist[0]}") + Gmass *= 1e-9 # JPL returns GM in units of km**3 / s**2 + + rotrate = get_rotrate(jpl.raw_response) + if rotrate is None: + rotrate = 0.0 + else: + rotrate = np.rad2deg(rotrate) + + rotpole = get_rotpole(jpl) + rot = rotpole*rotrate + + return Gmass,Rpl,rot + + +def horizons_query(id, ephemerides_start_date, exclude_spacecraft=True, verbose=False,**kwargs): + """ + Queries JPL/Horizons for a body matching the id. If one is found, a HorizonsClass object is returned for the first object that + matches the passed id string. If more than one match is found, a list of alternate ids is also returned. If no object is found + then None is returned. + + Parameters + ---------- + id : string + A string identifying which body is requested from JPL/Horizons + ephemerides_start_date : string + Date to use when obtaining the ephemerides in the format YYYY-MM-DD. + exclude_spacecraft: bool (optional) - Default True + Indicate whether spacecraft ids should be exluded from the alternate id list + verbose: bool (optional) - Default True + Indicate whether to print messages about the query or not + + Returns + ------- + jpl: HorizonsClass | None + An astroquery.jplhorizons HorizonsClass object. Or None if no match was found. + altid: string list | None + A list of alternate ids if more than one object matches the list + + """ + + def get_altid(errstr,exclude_spacecraft=True): + """ + Parses the error message returned from a failed Horizons query. If the query failed because of an ambiguous id, then it will + return a list of id values that could possibly match the query. + not found + + Parameters + ---------- + raw_response : string + Raw response from the JPL Horizons query + + Returns + ------- + MSun_over_Mpl : float + """ + if "ID" in errstr: + altid = errstr.split('ID')[1] + altid = altid.split('\n')[2:-1] + altname = [] + for n,v in enumerate(altid): + altid[n] = v.strip().split(' ')[0] + altname.append(' '.join(v.strip().split(' ')[1:]).strip().split(' ')[0]) + if exclude_spacecraft: + altname = [altname[i] for i,n in enumerate(altid) if int(n) > 0] + altid = [n for n in altid if int(n) > 0] + + return altid,altname + else: + return None,None + + + # Horizons date time internal variables + tstart = datetime.date.fromisoformat(ephemerides_start_date) + tstep = datetime.timedelta(days=1) + tend = tstart + tstep + ephemerides_end_date = tend.isoformat() + ephemerides_step = '1d' + + try: + jpl = Horizons(id=id, location='@sun', + epochs={'start': ephemerides_start_date, 'stop': ephemerides_end_date, + 'step': ephemerides_step},**kwargs) + eph=jpl.ephemerides() + altid = [id] + altname =[jpl.table['targetname'][0]] + except Exception as e: + altid,altname = get_altid(str(e)) + if altid is not None and len(altid) >0: # Return the first matching id + id = altid[0] + jpl = Horizons(id=id, location='@sun', + epochs={'start': ephemerides_start_date, 'stop': ephemerides_end_date, + 'step': ephemerides_step}) + eph=jpl.ephemerides() + else: + print(f"Could not find {id} in the JPL/Horizons system") + return None,None,None + if verbose: + print(f"Found matching body: {altname[0]} ({altid[0]})") + if len(altid) > 1: + print(" Alternate matching bodies:") + for i,v in enumerate(altid): + if i > 0: + print(f" {altname[i]} ({v})") + + return jpl,altid,altname + + def solar_system_horizons(name: str, param: Dict, ephemerides_start_date: str, - id: int | None = None): + ephemeris_id: str | None = None, + **kwargs: Any): """ Initializes a Swiftest dataset containing the major planets of the Solar System at a particular data from JPL/Horizons Parameters ---------- + name : string + Name of body to add to Dataset. If `id` is not supplied, this is also what will be searche for in the JPL/Horizon's database. + The first matching body is found (for major planets, this is the barycenter of a planet-satellite system) param : dict Swiftest paramuration parameters. This method uses the unit conversion factors to convert from JPL's AU-day system into the system specified in the param file ephemerides_start_date : string Date to use when obtaining the ephemerides in the format YYYY-MM-DD. + ephemeris_id : string (optional) + If passed, this is passed to Horizons instead of `name`. This can be used to find a more precise body than given by `name`. + **kwargs: Any + Additional keyword arguments to pass to the query method (see https://astroquery.readthedocs.io/en/latest/jplhorizons/jplhorizons.html) Returns ------- ds : xarray dataset - """ - # Planet ids - planetid = { - 'Sun': '0', - 'Mercury': '1', - 'Venus': '2', - 'Earth': '3', - 'Mars': '4', - 'Jupiter': '5', - 'Saturn': '6', - 'Uranus': '7', - 'Neptune': '8', - 'Pluto': '9' - } + Initial conditions of body formatted for Swiftest - if name in planetid: - ispl = True - id = planetid[name] - else: - ispl = False - print(f"\nMassive body {name} not found or not yet supported") - print("This will be created as a massless test particle") - if id is None: - print("ID value required for this input type") - return - - # Planet MSun/M ratio - MSun_over_Mpl = { - 'Sun' : np.longdouble(1.0), - 'Mercury': np.longdouble(6023600.0), - 'Venus': np.longdouble(408523.71), - 'Earth': np.longdouble(328900.56), - 'Mars': np.longdouble(3098708.), - 'Jupiter': np.longdouble(1047.3486), - 'Saturn': np.longdouble(3497.898), - 'Uranus': np.longdouble(22902.98), - 'Neptune': np.longdouble(19412.24), - 'Pluto': np.longdouble(1.35e8) - } - - # Planet radii in AU - planetradius = { - 'Sun' : swiftest.RSun / swiftest.AU2M, - 'Mercury': np.longdouble(2439.4e3) / swiftest.AU2M, - 'Venus': np.longdouble(6051.8e3) / swiftest.AU2M, - 'Earth': np.longdouble(6371.0084e3) / swiftest.AU2M, # Earth only for radius - 'Mars': np.longdouble(3389.50e3) / swiftest.AU2M, - 'Jupiter': np.longdouble(69911e3) / swiftest.AU2M, - 'Saturn': np.longdouble(58232.0e3) / swiftest.AU2M, - 'Uranus': np.longdouble(25362.e3) / swiftest.AU2M, - 'Neptune': np.longdouble(24622.e3) / swiftest.AU2M, - 'Pluto': np.longdouble(1188.3e3 / swiftest.AU2M) - } + Notes + -------- + When passing `name` == "Earth" or `name` == "Pluto", it a body is generated that has initial conditions matching the system + barycenter and mass equal to the sum of Earth+Moon or Pluto+Charon. To obtain initial conditions for either Earth or Pluto alone, + pass `ephemeris_id` == "399" for Earth or `id` == "999" for Pluto. + """ - planetrot = { - 'Sun' : np.longdouble(360.0 / 25.05) / swiftest.JD2S, # Approximate - 'Mercury': np.longdouble(360.0 / 58.646) / swiftest.JD2S, - 'Venus': np.longdouble(360.0 / 243.0226 ) / swiftest.JD2S, - 'Earth': np.longdouble(360.0 / 0.99726968) / swiftest.JD2S, - 'Mars': np.longdouble(360.0 / 1.025957) / swiftest.JD2S, - 'Jupiter': np.longdouble(360.0 / (9.9250 / 24.0) ) / swiftest.JD2S, - 'Saturn': np.longdouble(360.0 / (10.656 / 24.0) ) / swiftest.JD2S, - 'Uranus': np.longdouble(360.0 / 0.71833) / swiftest.JD2S, - 'Neptune': np.longdouble(360.0 / 0.6713) / swiftest.JD2S, - 'Pluto': np.longdouble(360.0 / 6.387230) / swiftest.JD2S - } - planetIpz = { # Only the polar moments of inertia are used currently. Where the quantity is unkown, we just use the value of a sphere = 0.4 'Sun' : np.longdouble(0.070), 'Mercury' : np.longdouble(0.346), @@ -132,9 +312,7 @@ def solar_system_horizons(name: str, J2RP2 = swiftest.J2Sun * (swiftest.RSun / param['DU2M']) ** 2 J4RP4 = swiftest.J4Sun * (swiftest.RSun / param['DU2M']) ** 4 - solarpole = SkyCoord(ra=286.13 * u.degree, dec=63.87 * u.degree) - solarrot = planetrot['Sun'] * param['TU2S'] - rotcb = solarpole.cartesian * solarrot + rotcb = swiftest.rotSun * param['TU2S'] rotcb = np.array([rotcb.x.value, rotcb.y.value, rotcb.z.value]) Ipsun = np.array([0.0, 0.0, planetIpz['Sun']]) @@ -157,7 +335,7 @@ def solar_system_horizons(name: str, J2 = None J4 = None - if name == "Sun" : # Create central body + if name == "Sun" or ephemeris_id == "0": # Create central body print("Creating the Sun as a central body") Gmass = GMcb Rpl = Rcb @@ -167,62 +345,71 @@ def solar_system_horizons(name: str, Ip = Ipsun rot = rotcb else: # Fetch solar system ephemerides from Horizons - print(f"Fetching ephemerides data for {name} from JPL/Horizons") - - # Horizons date time internal variables - tstart = datetime.date.fromisoformat(ephemerides_start_date) - tstep = datetime.timedelta(days=1) - tend = tstart + tstep - ephemerides_end_date = tend.isoformat() - ephemerides_step = '1d' - - pldata = {} - pldata[name] = Horizons(id=id, location='@sun', - epochs={'start': ephemerides_start_date, 'stop': ephemerides_end_date, - 'step': ephemerides_step}) + if ephemeris_id is None: + ephemeris_id = name + + print(f"Fetching ephemerides data for {ephemeris_id} from JPL/Horizons") + + jpl,altid,altname = horizons_query(ephemeris_id,ephemerides_start_date,**kwargs) + if jpl is not None: + print(f"Found ephemerides data for {altname[0]} ({altid[0]}) from JPL/Horizons") + if name == None: + name = altname[0] + else: + return None if param['IN_FORM'] == 'XV': - rx = pldata[name].vectors()['x'][0] * DCONV - ry = pldata[name].vectors()['y'][0] * DCONV - rz = pldata[name].vectors()['z'][0] * DCONV - vx = pldata[name].vectors()['vx'][0] * VCONV - vy = pldata[name].vectors()['vy'][0] * VCONV - vz = pldata[name].vectors()['vz'][0] * VCONV + rx = jpl.vectors()['x'][0] * DCONV + ry = jpl.vectors()['y'][0] * DCONV + rz = jpl.vectors()['z'][0] * DCONV + vx = jpl.vectors()['vx'][0] * VCONV + vy = jpl.vectors()['vy'][0] * VCONV + vz = jpl.vectors()['vz'][0] * VCONV rh = np.array([rx,ry,rz]) vh = np.array([vx,vy,vz]) elif param['IN_FORM'] == 'EL': - a = pldata[name].elements()['a'][0] * DCONV - e = pldata[name].elements()['e'][0] - inc = pldata[name].elements()['incl'][0] - capom = pldata[name].elements()['Omega'][0] - omega = pldata[name].elements()['w'][0] - capm = pldata[name].elements()['M'][0] - - if ispl: - Gmass = GMcb / MSun_over_Mpl[name] + a = jpl.elements()['a'][0] * DCONV + e = jpl.elements()['e'][0] + inc = jpl.elements()['incl'][0] + capom = jpl.elements()['Omega'][0] + omega = jpl.elements()['w'][0] + capm = jpl.elements()['M'][0] + + Gmass,Rpl,rot = horizons_get_physical_properties(altid,**kwargs) + # If the user inputs "Earth" or Pluto, then the Earth-Moon or Pluto-Charon barycenter and combined mass is used. + # To use the Earth or Pluto alone, simply pass "399" or "999", respectively to name + if ephemeris_id == "Earth": + print("Combining mass of Earth and the Moon") + Gmass_moon,tmp,tmp = horizons_get_physical_properties(["301"],**kwargs) + Gmass += Gmass_moon + elif ephemeris_id == "Pluto": + print("Combining mass of Pluto and Charon") + Gmass_charon,tmp,tmp = horizons_get_physical_properties(["901"],**kwargs) + Gmass += Gmass_charon + + if Gmass is not None: if param['CHK_CLOSE']: - Rpl = planetradius[name] * DCONV + Rpl /= param['DU2M'] # Generate planet value vectors if (param['RHILL_PRESENT']): - rhill = pldata[name].elements()['a'][0] * DCONV * (3 * MSun_over_Mpl[name]) ** (-THIRDLONG) + rhill = jpl.elements()['a'][0] * DCONV * (3 * Gmass / GMcb) ** (-THIRDLONG) if (param['ROTATION']): - RA = pldata[name].ephemerides()['NPole_RA'][0] - DEC = pldata[name].ephemerides()['NPole_DEC'][0] - - rotpole = SkyCoord(ra=RA * u.degree, dec=DEC * u.degree) - rotrate = planetrot[name] * param['TU2S'] - rot = rotpole.cartesian * rotrate - rot = np.array([rot.x.value, rot.y.value, rot.z.value]) - Ip = np.array([0.0, 0.0, planetIpz[name]]) - + rot *= param['TU2S'] + if name in planetIpz: + Ip = np.array([0.0, 0.0, planetIpz[name]]) + else: + Ip = np.array([0.4, 0.4, 0.4]) else: Gmass = None - if id is None: - id = planetid[name] + # Only the Sun gets assigned its own special id for now. All other ids will be sorted later + if name == "Sun": + id = 0 + else: + id = -1 return id,name,a,e,inc,capom,omega,capm,rh,vh,Gmass,Rpl,rhill,Ip,rot,J2,J4 diff --git a/python/swiftest/swiftest/simulation_class.py b/python/swiftest/swiftest/simulation_class.py index 399ab0d38..c9f125173 100644 --- a/python/swiftest/swiftest/simulation_class.py +++ b/python/swiftest/swiftest/simulation_class.py @@ -113,7 +113,7 @@ def __init__(self,read_param: bool = False, Parameter input file equivalent: `ISTEP_OUT` tstep_out : float, optional The approximate time between when outputs are written to file. Passing this computes - `istep_out = floor(tstep_out/dt)`. *Note*: only `istep_out` or `tstep_out` can be set. + `istep_out = floor(tstep_out/dt)`. *Note*: only `istep_out` or `tstep_out` can be set. `tstep_out` must be less than `tstop` Parameter input file equivalent: None nstep_out : int, optional The total number of times that outputs are written to file. Passing this allows for a geometric progression of output steps: @@ -405,7 +405,7 @@ def _run_swiftest_driver(self): # Get current environment variables env = os.environ.copy() - cmd = f"{env['SHELL']} -l {self.driver_script}" + cmd = f"{env['SHELL']} {self.driver_script}" def _type_scrub(output_data): @@ -676,6 +676,10 @@ def set_simulation_time(self, return {} else: update_list.append("istep_out") + + if tstep_out is not None and tstep_out > tstop: + warnings.warn("tstep_out must be less than tstop. Setting tstep_out=tstop",stacklevel=2) + tstep_out = tstop if tstep_out is not None and dt is not None: istep_out = int(tstep_out / dt) @@ -683,6 +687,7 @@ def set_simulation_time(self, if istep_out is not None: self.param['ISTEP_OUT'] = int(istep_out) + if nstep_out is not None: if istep_out is None: warnings.warn("nstep_out requires either istep_out or tstep_out to also be set", stacklevel=2) @@ -2157,57 +2162,61 @@ def get_distance_range(self, arg_list: str | List[str] | None = None, verbose: b return range_dict def add_solar_system_body(self, - name: str | List[str], + name: str | List[str] | None = None, ephemeris_id: int | List[int] | None = None, date: str | None = None, - source: str = "HORIZONS"): + source: str = "HORIZONS", + **kwargs: Any): """ - Adds a solar system body to an existing simulation Dataset from the JPL Horizons ephemeris service. - - The following are name/ephemeris_id pairs that are currently known to Swiftest, and therefore have - physical properties that can be used to make massive bodies. + Adds a solar system body to an existing simulation Dataset from the JPL Horizons ephemeris service. The JPL Horizons service + will be searched for a body matching the string passed by `name`, or alternatively `ephemeris_id` if passed. Bodies will be + named in the Swiftest initial conditions Dataset using `name`. Use `ephemeris_id` to have finer control over which body is + searched in Horizons while using a custom name. + + If `name` is not passed, then the target name property is used as the name. You must pass either `name` and/or `ephemeris_id` + + When passing `name` == "Earth" or `name` == "Pluto", it a body is generated that has initial conditions matching the system + barycenter and mass equal to the sum of Earth+Moon or Pluto+Charon. + + To obtain initial conditions for either Earth or Pluto alone, pass `ephemeris_id` == "399" for Earth or + `ephemeris_id` == "999" for Pluto. - Sun : 0 - Mercury : 1 - Venus : 2 - Earth : 3 - Mars : 4 - Jupiter : 5 - Saturn : 6 - Uranus : 7 - Neptune : 8 - Pluto : 9 Parameters ---------- - name : str | List[str] - Add solar system body by name. - Bodies not on this list will be added as test particles, but additional properties can be added later if - desired. + name : str, optional | List[str] + Add solar system body by name. This will be the name used in the Swiftest initial conditions Dataset unless not supplied ephemeris_id : int | List[int], optional but must be the same length as `name` if passed. - Use id if the body you wish to add is recognized by Swiftest. In that case, the id is passed to the - ephemeris service and the name is used. The body specified by `id` supercedes that given by `name`. + In that case, the id is passed to the ephemeris service and the name is used. date : str, optional ISO-formatted date sto use when obtaining the ephemerides in the format YYYY-MM-DD. Defaults to value set by `set_ephemeris_date`. source : str, default "Horizons" The source of the ephemerides. >*Note.* Currently only the JPL Horizons ephemeris is implemented, so this is ignored. + **kwargs: Any + Additional keyword arguments to pass to the query method (i.e. astroquery.Horizons) Returns ------- None initial conditions data stored as an Xarray Dataset in the init_cond instance variable """ + + if name == None and ephemeris_id == None: + warnings.warn("Either `name` and/or `ephemeris_id` must be supplied to add_solar_system_body") + return None - if type(name) is str: - name = [name] if ephemeris_id is not None: - if type(ephemeris_id) is int: + if type(ephemeris_id) is int or type(ephemeris_id) is str: ephemeris_id = [ephemeris_id] - if len(ephemeris_id) != len(name): + if name is None: + name = [None] * len(ephemeris_id) + elif len(ephemeris_id) != len(name): warnings.warn(f"The length of ephemeris_id ({len(ephemeris_id)}) does not match the length of name ({len(name)})",stacklevel=2) return None else: + if type(name) is str: + name = [name] ephemeris_id = [None] * len(name) if self.ephemeris_date is None: @@ -2226,10 +2235,15 @@ def add_solar_system_body(self, body_list = [] for i,n in enumerate(name): - body_list.append(init_cond.solar_system_horizons(n, self.param, date, id=ephemeris_id[i])) + body = init_cond.solar_system_horizons(n, self.param, date, ephemeris_id=ephemeris_id[i],**kwargs) + if body is not None: + body_list.append(body) #Convert the list receieved from the solar_system_horizons output and turn it into arguments to vec2xr - if len(body_list) == 1: + if len(body_list) == 0: + print("No valid bodies found") + return + elif len(body_list) == 1: values = list(np.hsplit(np.array(body_list[0],dtype=np.dtype(object)),17)) else: values = list(np.squeeze(np.hsplit(np.array(body_list,np.dtype(object)),17))) @@ -2241,6 +2255,7 @@ def add_solar_system_body(self, for k,v in kwargs.items(): if k in scalar_ints: + v[v == None] = -1 kwargs[k] = v.astype(int) elif k in scalar_floats: kwargs[k] = v.astype(np.float64) @@ -2254,11 +2269,21 @@ def add_solar_system_body(self, kwargs['time'] = np.array([self.param['TSTART']]) + if len(self.data) == 0: + maxid = -1 + else: + maxid = self.data.id.max().values[()] + + nbodies = kwargs["name"].size + kwargs['id'] = np.where(kwargs['id'] < 0,np.arange(start=maxid+1,stop=maxid+1+nbodies,dtype=int),kwargs['id']) + dsnew = init_cond.vec2xr(self.param,**kwargs) dsnew = self._combine_and_fix_dsnew(dsnew) if dsnew['npl'] > 0 or dsnew['ntp'] > 0: self.save(verbose=False) + + self.init_cond = self.data.copy(deep=True) diff --git a/python/swiftest/swiftest/tool.py b/python/swiftest/swiftest/tool.py index fdcd5585f..e81194d37 100644 --- a/python/swiftest/swiftest/tool.py +++ b/python/swiftest/swiftest/tool.py @@ -359,7 +359,7 @@ def el2xv_vec(mu, a, ecc, inc, Omega, omega, M): longitude of ascending node (degrees) omega : (n) float array argument of periapsis (degrees) - M : (n) flaot array + M : (n) float array Mean anomaly (degrees) Returns diff --git a/python/swiftest/tests/test_suite.py b/python/swiftest/tests/test_suite.py new file mode 100755 index 000000000..d4ddc7c50 --- /dev/null +++ b/python/swiftest/tests/test_suite.py @@ -0,0 +1,166 @@ +#!/usr/bin/env python3 + +""" + Copyright 2023 - David Minton, Carlisle Wishard, Jennifer Pouplin, Jake Elliott, & Dana Singh + This file is part of Swiftest. + Swiftest is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License + as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. + Swiftest is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty + of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. + You should have received a copy of the GNU General Public License along with Swiftest. + If not, see: https://www.gnu.org/licenses. +""" + +""" +Tests that energy and momentum errors are within tolerances in a Swiftest simulation + +Input +------ + +Output +------ +None +""" + +import swiftest +import unittest +import os +import shutil +import numpy as np +from numpy.random import default_rng + +rng = default_rng(seed=123) + +major_bodies = ["Sun","Mercury","Venus","Earth","Mars","Jupiter","Saturn","Uranus","Neptune"] +param = {} + + +class TestSwiftest(unittest.TestCase): + + def test01_gen_ic(self): + """ + Tests that Swiftest is able to successfully generate a set of initial conditions in a file without any exceptions being raised + """ + print("\ntest_gen_ic: Test whether we can generate simulation initial conditions test") + # Files that are expected to be generated: + simdir = "simdata" + file_list = [simdir, os.path.join(simdir,"param.in"), os.path.join(simdir,"init_cond.nc")] + + sim = swiftest.Simulation() + sim.clean() + + # Add the modern planets and the Sun using the JPL Horizons Database. + # Add the modern planets and the Sun using the JPL Horizons Database. + sim.add_solar_system_body(major_bodies) + + # Display the run configuration parameters. + param = sim.get_parameter(verbose=False) + sim.save() + + for f in file_list: + self.assertTrue(os.path.exists(f)) + + def test02_read_ic(self): + """ + Tests that Swiftest is able to read a set of pre-existing initial conditions files and that they contain the correct data + """ + print("\ntest_read_ic: Test whether we can read back initial conditions files created by test_gen_ic") + sim = swiftest.Simulation(read_param=True) + # Check if all names in Dataset read in from file match the expected list of names + self.assertTrue((major_bodies == sim.init_cond['name']).all(), msg="Name mismatch in Dataset") + + # Check to see if all parameter values read in from file match the expected parameters saved when generating the file + self.assertTrue(all([v == param[k] for k,v in sim.param.items() if k in param])) + + def test03_integrators(self): + """ + Tests that Swiftest is able to integrate a collection of massive bodies and test particles with all available integrators + """ + print("\ntest_integrators: Tests that Swiftest is able to integrate a collection of massive bodies and test particles with all available integrators") + sim = swiftest.Simulation(read_param=True) + + # Add 10 user-defined test particles. + ntp = 10 + + name_tp = [f"TestParticle_{i:02}" for i in range(1,ntp+1)] + a_tp = rng.uniform(0.3, 1.5, ntp) + e_tp = rng.uniform(0.0, 0.2, ntp) + inc_tp = rng.uniform(0.0, 10, ntp) + capom_tp = rng.uniform(0.0, 360.0, ntp) + omega_tp = rng.uniform(0.0, 360.0, ntp) + capm_tp = rng.uniform(0.0, 360.0, ntp) + + integrators= ["whm","helio","rmvs","symba"] + sim.add_body(name=name_tp, a=a_tp, e=e_tp, inc=inc_tp, capom=capom_tp, omega=omega_tp, capm=capm_tp) + sim.set_parameter(tstart=0.0, tstop=0.02, dt=0.01, istep_out=1, dump_cadence=0) + integrators= ["whm","helio","rmvs","symba"] + for i in integrators: + try: + sim.run(integrator=i) + except: + self.fail(f"Failed with integrator {i}") + + + def test04_conservation(self): + """ + Tests that Swiftest conserves mass, energy, and momentum to within acceptable tolerances. + """ + print("\ntest_conservation: Tests that Swiftest conserves mass, energy, and momentum to within acceptable tolerances.") + + # Error limits + L_slope_limit = 1e-10 + E_slope_limit = 1e-8 + GM_limit = 1e-14 + + sim = swiftest.Simulation() + sim.clean() + + sim.add_solar_system_body(major_bodies) + + dt = 0.01 + nout = 1000 + tstop = 1e4 + tstep_out = tstop / nout + + sim.run(tstart=0.0, tstop=tstop, dt=dt, tstep_out=tstep_out, dump_cadence=0, compute_conservation_values=True, integrator="symba") + + def fit_func(x,slope,b): + """ + Linear function for curve fitting + """ + return slope * x + b + + # Calculate the angular momentum error + sim.data['L_tot'] = sim.data['L_orbit'] + sim.data['L_spin'] + sim.data['L_escape'] + sim.data['DL'] = sim.data['L_tot'] - sim.data['L_tot'].isel(time=0) + L_error = swiftest.tool.magnitude(sim.data,'DL') / swiftest.tool.magnitude(sim.data.isel(time=0), 'L_tot') + + # Calculate the energy error + E_error = (sim.data['TE'] - sim.data['TE'].isel(time=0)) / sim.data['TE'].isel(time=0) + + # Calculate the mass error + sim.data['GMtot'] = sim.data['Gmass'].sum(dim='name',skipna=True) + sim.data['GMescape'] + GM_error = (sim.data['GMtot'] - sim.data['GMtot'].isel(time=0)) / sim.data['GMtot'].isel(time=0) + GM_final = GM_error.isel(time=-1).values + + # Compute the slope of the error vs time fit + E_fit_result = E_error.curvefit("time",fit_func) + L_fit_result = L_error.curvefit("time",fit_func) + + E_slope = E_fit_result['curvefit_coefficients'].sel(param='slope').values + L_slope = L_fit_result['curvefit_coefficients'].sel(param='slope').values + + # Print the final errors + print("\n") + print(f"Angular momentum error slope: {L_slope:.2e}/{sim.TU_name}") + print(f"Energy error slope: {E_slope:.2e}/{sim.TU_name}") + print(f"Final Mass Error: {GM_final:.2e}") + + self.assertLess(np.abs(L_slope),L_slope_limit, msg=f"Angular Momentum Error of {L_slope:.2e}/{sim.TU_name} higher than threshold value of {L_slope_limit:.2e}/{sim.TU_name}") + self.assertLess(np.abs(E_slope),E_slope_limit, msg=f"Energy Error of {E_slope:.2e}/{sim.TU_name} higher than threshold value of {E_slope_limit:.2e}/{sim.TU_name}") + self.assertLess(np.abs(GM_final),GM_limit, msg=f"Mass Error of {GM_final:.2e} higher than threshold value of {GM_limit:.2e}") + +if __name__ == '__main__': + unittest.main() + if os.path.exists("simdir"): + shutil.rmtree("simdir") diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index f360bcdbc..2899ae49f 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -167,3 +167,5 @@ end program TestDoConcurrentLoc TARGET_COMPILE_DEFINITIONS(${SWIFTEST_DRIVER} PRIVATE -DDOCONLOC) ENDIF (DOCONLOC_WORKS) ENDIF (USE_OPENMP) + + diff --git a/src/base/base_module.f90 b/src/base/base_module.f90 index 51266238f..f47bb9f7a 100644 --- a/src/base/base_module.f90 +++ b/src/base/base_module.f90 @@ -586,7 +586,7 @@ subroutine base_util_dealloc_storage(self) end subroutine base_util_dealloc_storage - subroutine base_util_exit(code) + subroutine base_util_exit(code,unit) !! author: David A. Minton !! !! Print termination message and exit program @@ -596,25 +596,33 @@ subroutine base_util_exit(code) implicit none ! Arguments integer(I4B), intent(in) :: code + integer(I4B), intent(in), optional :: unit ! Internals - character(*), parameter :: BAR = '("------------------------------------------------")' - character(*), parameter :: SUCCESS_MSG = '(/, "Normal termination of Swiftest (version ", f3.1, ")")' - character(*), parameter :: FAIL_MSG = '(/, "Terminating Swiftest (version ", f3.1, ") due to error!!")' + character(*), parameter :: BAR = '("---------------------------------------------------")' + character(*), parameter :: SUCCESS_MSG = '(/, "Normal termination of Swiftest (version ", A, ")")' + character(*), parameter :: FAIL_MSG = '(/, "Terminating Swiftest (version ", A, ") due to error!!")' character(*), parameter :: USAGE_MSG = '("Usage: swiftest ' // & '[{standard}|compact|progress]")' character(*), parameter :: HELP_MSG = USAGE_MSG + integer(I4B) :: iu + + if (present(unit)) then + iu = unit + else + iu = OUTPUT_UNIT + end if select case(code) case(SUCCESS) - write(*, SUCCESS_MSG) VERSION_NUMBER - write(*, BAR) + write(iu, SUCCESS_MSG) VERSION + write(iu, BAR) case(USAGE) - write(*, USAGE_MSG) + write(iu, USAGE_MSG) case(HELP) - write(*, HELP_MSG) + write(iu, HELP_MSG) case default - write(*, FAIL_MSG) VERSION_NUMBER - write(*, BAR) + write(iu, FAIL_MSG) VERSION + write(iu, BAR) error stop end select diff --git a/src/collision/collision_generate.f90 b/src/collision/collision_generate.f90 index 7f529ba02..fc3b7ed99 100644 --- a/src/collision/collision_generate.f90 +++ b/src/collision/collision_generate.f90 @@ -88,7 +88,7 @@ module subroutine collision_generate_bounce(self, nbody_system, param, t) call self%merge(nbody_system, param, t) ! Use the default collision model, which is merge case default call swiftest_io_log_one_message(COLLISION_LOG_OUT,"Error in swiftest_collision, unrecognized collision regime") - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,unit=param%display_unit) end select end associate end select diff --git a/src/collision/collision_io.f90 b/src/collision/collision_io.f90 index cb4a544a0..c1b5015de 100644 --- a/src/collision/collision_io.f90 +++ b/src/collision/collision_io.f90 @@ -267,7 +267,7 @@ module subroutine collision_io_netcdf_initialize_output(self, param) 667 continue write(*,*) "Error creating fragmentation output file. " // trim(adjustl(errmsg)) - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,unit=param%display_unit) end subroutine collision_io_netcdf_initialize_output diff --git a/src/collision/collision_resolve.f90 b/src/collision/collision_resolve.f90 index 0166711ab..72e6ad445 100644 --- a/src/collision/collision_resolve.f90 +++ b/src/collision/collision_resolve.f90 @@ -648,7 +648,7 @@ module subroutine collision_resolve_plpl(self, nbody_system, param, t, dt, irec) if (loop == MAXCASCADE) then call swiftest_io_log_one_message(COLLISION_LOG_OUT,"A runaway collisional cascade has been detected in collision_resolve_plpl.") call swiftest_io_log_one_message(COLLISION_LOG_OUT,"Consider reducing the step size or changing the parameters in the collisional model to reduce the number of fragments.") - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,unit=param%display_unit) end if end associate end do diff --git a/src/encounter/encounter_io.f90 b/src/encounter/encounter_io.f90 index fb5de048f..171191b39 100644 --- a/src/encounter/encounter_io.f90 +++ b/src/encounter/encounter_io.f90 @@ -145,7 +145,7 @@ module subroutine encounter_io_netcdf_initialize_output(self, param) 667 continue write(*,*) "Error creating encounter output file. " // trim(adjustl(errmsg)) - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) end subroutine encounter_io_netcdf_initialize_output diff --git a/src/fraggle/fraggle_generate.f90 b/src/fraggle/fraggle_generate.f90 index cd1df7033..aeeae27d9 100644 --- a/src/fraggle/fraggle_generate.f90 +++ b/src/fraggle/fraggle_generate.f90 @@ -52,7 +52,7 @@ module subroutine fraggle_generate(self, nbody_system, param, t) message = "Supercatastrophic disruption between" case default write(*,*) "Error in swiftest_collision, unrecognized collision regime" - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) end select call collision_io_collider_message(pl, impactors%id, message) call swiftest_io_log_one_message(COLLISION_LOG_OUT, trim(adjustl(message))) diff --git a/src/globals/globals_module.f90 b/src/globals/globals_module.f90 index dd58d6dae..a708b080d 100644 --- a/src/globals/globals_module.f90 +++ b/src/globals/globals_module.f90 @@ -44,7 +44,7 @@ module globals integer(I4B), parameter :: UPPERCASE_OFFSET = iachar('A') - iachar('a') !! ASCII character set parameter for lower to upper !! conversion - offset between upper and lower - real(SP), parameter :: VERSION_NUMBER = 1.0_SP !! Swiftest version + character(*), parameter :: VERSION = "2023.08.00" !! Swiftest version !> Symbolic name for integrator types character(*), parameter :: UNKNOWN_INTEGRATOR = "UKNOWN INTEGRATOR" diff --git a/src/rmvs/rmvs_step.f90 b/src/rmvs/rmvs_step.f90 index 4a5dcf3af..098622037 100644 --- a/src/rmvs/rmvs_step.f90 +++ b/src/rmvs/rmvs_step.f90 @@ -293,7 +293,7 @@ subroutine rmvs_interp_in(cb, pl, nbody_system, param, dt, outer_index) vtmp(:,i), new_line('a'), & " STOPPING " call swiftest_io_log_one_message(COLLISION_LOG_OUT, message) - call base_util_exit(failure) + call base_util_exit(FAILURE,param%display_unit) end if end do end if @@ -317,7 +317,7 @@ subroutine rmvs_interp_in(cb, pl, nbody_system, param, dt, outer_index) write(*, *) xtmp(:,i) write(*, *) vtmp(:,i) write(*, *) " STOPPING " - call base_util_exit(failure) + call base_util_exit(FAILURE,param%display_unit) end if end do end if diff --git a/src/swiftest/swiftest_driver.f90 b/src/swiftest/swiftest_driver.f90 index 6cec24feb..3aea4aab7 100644 --- a/src/swiftest/swiftest_driver.f90 +++ b/src/swiftest/swiftest_driver.f90 @@ -184,7 +184,7 @@ program swiftest_driver #ifdef COARRAY if (this_image() == 1) then #endif - call base_util_exit(SUCCESS) + call base_util_exit(SUCCESS,unit=param%display_unit) #ifdef COARRAY end if ! (this_image() == 1) #endif diff --git a/src/swiftest/swiftest_io.f90 b/src/swiftest/swiftest_io.f90 index 70ddac5f7..64731a7eb 100644 --- a/src/swiftest/swiftest_io.f90 +++ b/src/swiftest/swiftest_io.f90 @@ -196,7 +196,7 @@ module subroutine swiftest_io_conservation_report(self, param, lterminal) if (abs(nbody_system%Mtot_error) > 100 * epsilon(nbody_system%Mtot_error)) then write(*,*) "Severe error! Mass not conserved! Halting!" ! Save the frame of data to the bin file in the slot just after the present one for diagnostics - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) end if end if end associate @@ -344,7 +344,7 @@ module subroutine swiftest_io_dump_param(self, param_file_name) 667 continue write(*,*) "Error opening parameter dump file " // trim(adjustl(errmsg)) - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,self%display_unit) end subroutine swiftest_io_dump_param @@ -896,7 +896,7 @@ module subroutine swiftest_io_netcdf_initialize_output(self, param) 667 continue write(*,*) "Error creating NetCDF output file. " // trim(adjustl(errmsg)) - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) end subroutine swiftest_io_netcdf_initialize_output @@ -1161,7 +1161,7 @@ module function swiftest_io_netcdf_read_frame_system(self, nc, param) result(ier write(*,*) "Error reading in NetCDF file: The recorded value of ntp does not match the number of active test particles" write(*,*) "Recorded: ",ntp write(*,*) "Active : ",ntp_check - call base_util_exit(failure) + call base_util_exit(FAILURE,param%display_unit) end if ! Now read in each variable and split the outputs by body type @@ -2879,7 +2879,7 @@ module subroutine swiftest_io_read_in_cb(self, param) 667 continue write(*,*) "Error reading central body file: " // trim(adjustl(errmsg)) - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) end subroutine swiftest_io_read_in_cb @@ -2922,7 +2922,7 @@ module subroutine swiftest_io_read_in_system(self, nc, param) end if ierr = self%read_frame(nc, tmp_param) deallocate(tmp_param) - if (ierr /=0) call base_util_exit(FAILURE) + if (ierr /=0) call base_util_exit(FAILURE,param%display_unit) end if param%loblatecb = ((abs(self%cb%j2rp2) > 0.0_DP) .or. (abs(self%cb%j4rp4) > 0.0_DP)) @@ -3044,7 +3044,7 @@ module function swiftest_io_read_frame_body(self, iu, param) result(ierr) class default write(*,*) "Error reading body file: " // trim(adjustl(errmsg)) end select - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) end function swiftest_io_read_frame_body @@ -3122,7 +3122,7 @@ module subroutine swiftest_io_set_display_param(self, display_style) 667 continue write(*,*) "Error opening swiftest log file: " // trim(adjustl(errmsg)) - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,self%display_unit) end subroutine swiftest_io_set_display_param @@ -3203,7 +3203,7 @@ module subroutine swiftest_io_initialize_output_file_system(self, nc, param) 667 continue write(*,*) "Error writing nbody_system frame: " // trim(adjustl(errmsg)) - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) end subroutine swiftest_io_initialize_output_file_system end submodule s_swiftest_io diff --git a/src/swiftest/swiftest_util.f90 b/src/swiftest/swiftest_util.f90 index f3ce3082f..c390be826 100644 --- a/src/swiftest/swiftest_util.f90 +++ b/src/swiftest/swiftest_util.f90 @@ -2313,6 +2313,7 @@ module subroutine swiftest_util_setup_construct_system(nbody_system, param) case (INT_BS) write(*,*) 'Bulirsch-Stoer integrator not yet enabled' case (INT_HELIO) + allocate(helio_nbody_system :: nbody_system) select type(nbody_system) class is (helio_nbody_system) allocate(helio_cb :: nbody_system%cb) @@ -2366,7 +2367,7 @@ module subroutine swiftest_util_setup_construct_system(nbody_system, param) write(*,*) 'RINGMOONS-SyMBA integrator not yet enabled' case default write(*,*) 'Unkown integrator',param%integrator - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) end select allocate(swiftest_particle_info :: nbody_system%cb%info) @@ -3328,7 +3329,7 @@ module subroutine swiftest_util_version() !! !! Adapted from David E. Kaufmann's Swifter routine: util_version.f90 implicit none - write(*, 200) VERSION_NUMBER + write(*, 200) VERSION 200 format(/, "************* Swiftest: Version ", f3.1, " *************", //, & "Based off of Swifter:", //, & "Authors:", //, & diff --git a/src/symba/symba_step.f90 b/src/symba/symba_step.f90 index 22cdddeb6..9da52440a 100644 --- a/src/symba/symba_step.f90 +++ b/src/symba/symba_step.f90 @@ -197,7 +197,7 @@ recursive module subroutine symba_step_recur_system(self, param, t, ireci) write(*, *) "SWIFTEST Warning:" write(*, *) " In symba_step_recur_system, local time step is too small" write(*, *) " Roundoff error will be important!" - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) END IF irecp = ireci + 1 if (ireci == 0) then diff --git a/src/whm/whm_drift.f90 b/src/whm/whm_drift.f90 index 1161919f2..fa0dcdd8f 100644 --- a/src/whm/whm_drift.f90 +++ b/src/whm/whm_drift.f90 @@ -50,7 +50,7 @@ module subroutine whm_drift_pl(self, nbody_system, param, dt) call swiftest_io_log_one_message(COLLISION_LOG_OUT, message) end if end do - call base_util_exit(FAILURE) + call base_util_exit(FAILURE,param%display_unit) end if end associate diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt new file mode 100644 index 000000000..17e6a6b62 --- /dev/null +++ b/test/CMakeLists.txt @@ -0,0 +1,17 @@ +ENABLE_TESTING() +FIND_PACKAGE(Python COMPONENTS Interpreter REQUIRED) + +SET(UNIT_TEST swiftest_test_suite) +ADD_CUSTOM_TARGET(${UNIT_TEST} ALL) +ADD_DEPENDENCIES(${UNIT_TEST} ${SWIFTEST_DRIVER}) + +ADD_TEST(NAME ${UNIT_TEST} COMMAND python ${CMAKE_SOURCE_DIR}/python/swiftest/tests/test_suite.py) + +ADD_CUSTOM_COMMAND( + TARGET ${UNIT_TEST} + COMMENT "Run tests" + POST_BUILD + WORKING_DIRECTORY ${TEST} + COMMAND ${CMAKE_CTEST_COMMAND} "${UNIT_TEST}" --output-on-failures +) + diff --git a/version.txt b/version.txt new file mode 100644 index 000000000..f2bfdb11a --- /dev/null +++ b/version.txt @@ -0,0 +1 @@ +2023.08.00 \ No newline at end of file