From aebd7702b5f34e3f9dccbbe687549ea6179c0557 Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Tue, 27 Jun 2023 10:31:35 -0400 Subject: [PATCH 01/11] Update README.md fixed link to clone repo, also added in cmake flags for gfortran --- README.md | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/README.md b/README.md index 3671832cb..e011b8430 100644 --- a/README.md +++ b/README.md @@ -46,7 +46,7 @@ Parallelization in Swiftest is done with OpenMP. Version 3.1.4 or higher is nece The easiest way to get Swiftest on your machine is to clone the GitHub repository. To do so, open a terminal window and type the following: ``` -$ git clone https://github.itap.purdue.edu/MintonGroup/swiftest.git +$ git clone https://github.com/carlislewishard/swiftest.git ``` If your cloned version is not already set to the master branch: @@ -85,7 +85,7 @@ CMake allows the user to specify a set of compiler flags to use during compilati 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), type the following: +To build Swiftest with the release flags (default) using the Intel fortran compiler (ifort), type the following: ``` $ cmake .. ``` @@ -97,6 +97,11 @@ To build with another set of flags, simply replace ```DEBUG``` in the above line Add ```-CMAKE_PREFIX_PATH=/path/to/netcdf/``` to these commands as needed. +If using the GCC fortran compiler (gfortran), add the following flags: +``` +cmake -DCMAKE_Fortran_FLAGS="-I/usr/lib64/gfortran/modules/ -ffree-line-length-512" ... +``` + After building Swiftest, make the executable using: ``` From f1bfc2fb48090f09b79a039db13f17ea624fdbfe Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Tue, 27 Jun 2023 10:37:31 -0400 Subject: [PATCH 02/11] Update README.md added flag to specify fortran compiler --- README.md | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/README.md b/README.md index e011b8430..f45a812c3 100644 --- a/README.md +++ b/README.md @@ -99,11 +99,14 @@ Add ```-CMAKE_PREFIX_PATH=/path/to/netcdf/``` to these commands as needed. If using the GCC fortran compiler (gfortran), add the following flags: ``` -cmake -DCMAKE_Fortran_FLAGS="-I/usr/lib64/gfortran/modules/ -ffree-line-length-512" ... +-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) ``` After building Swiftest, make the executable using: - ``` $ make ``` From 8feca172adcfe01395ebb49a4a993bc0d9249eb5 Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Tue, 27 Jun 2023 13:52:46 -0400 Subject: [PATCH 03/11] updated readme to add errors.py --- examples/Basic_Simulation/README.txt | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/examples/Basic_Simulation/README.txt b/examples/Basic_Simulation/README.txt index 9db2d5f54..ae626a9c5 100644 --- a/examples/Basic_Simulation/README.txt +++ b/examples/Basic_Simulation/README.txt @@ -11,13 +11,14 @@ README.txt Swiftest Example : Basic_Simulation Author : Carlisle Wishard and David Minton -Date : December 6, 2022 +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. - 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 This example is intended to be run with Swiftest SyMBA. For details on how to generate, run, and analyze this example, -see the Swiftest User Manual. \ No newline at end of file +see the Swiftest User Manual. From 47b0ad55bd3f7228dc58c3b1c15e6ee3411023bc Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Tue, 27 Jun 2023 13:53:39 -0400 Subject: [PATCH 04/11] added error analysis python script --- examples/Basic_Simulation/errors.py | 69 +++++++++++++++++++++++++++++ 1 file changed, 69 insertions(+) create mode 100644 examples/Basic_Simulation/errors.py diff --git a/examples/Basic_Simulation/errors.py b/examples/Basic_Simulation/errors.py new file mode 100644 index 000000000..581677949 --- /dev/null +++ b/examples/Basic_Simulation/errors.py @@ -0,0 +1,69 @@ +#!/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-4 +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.") \ No newline at end of file From 1fc4bef397867cc348b0fb832cc093c7629eb56b Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Tue, 27 Jun 2023 13:54:14 -0400 Subject: [PATCH 05/11] Update .gitignore to add errors.py --- examples/Basic_Simulation/.gitignore | 1 + 1 file changed, 1 insertion(+) diff --git a/examples/Basic_Simulation/.gitignore b/examples/Basic_Simulation/.gitignore index 728fe8873..9e4d02aad 100644 --- a/examples/Basic_Simulation/.gitignore +++ b/examples/Basic_Simulation/.gitignore @@ -2,4 +2,5 @@ !.gitignore !initial_conditions.py !output_reader.py +!errors.py !README.txt From e5181a6890a9e8f32e145d92a7df08a0e7d08e47 Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Tue, 27 Jun 2023 13:54:56 -0400 Subject: [PATCH 06/11] Updated output_reader python script --- examples/Basic_Simulation/output_reader.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/examples/Basic_Simulation/output_reader.py b/examples/Basic_Simulation/output_reader.py index 851384bf3..b8cf36a70 100644 --- a/examples/Basic_Simulation/output_reader.py +++ b/examples/Basic_Simulation/output_reader.py @@ -43,4 +43,4 @@ ax.text(1.5, 0.35, f"t = {sim.data['time'].isel(time=0).values} years", size=10, ha="left") plt.tight_layout() plt.show() -fig.savefig("output.eps", dpi=300, facecolor='white', transparent=False, bbox_inches="tight") +fig.savefig("../examples/Basic_Simulation/output.eps", dpi=300, facecolor='white', transparent=False, bbox_inches="tight") From 1551a512e86f119ef3795e145098548fbe5b225b Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Tue, 27 Jun 2023 13:55:49 -0400 Subject: [PATCH 07/11] added automatic test that runs the Basic_Simulation example and generates output --- CMakeLists.txt | 14 ++++++++++++++ 1 file changed, 14 insertions(+) diff --git a/CMakeLists.txt b/CMakeLists.txt index 961683f3e..4d36d9612 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -73,3 +73,17 @@ ADD_SUBDIRECTORY(${SRC} ${BIN}) 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***" +) + From 99aff8e3c98bdc3c5ea915162d987e3bb87a8e11 Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Tue, 27 Jun 2023 13:59:52 -0400 Subject: [PATCH 08/11] prevented output_reader from showing plot --- examples/Basic_Simulation/output_reader.py | 1 - 1 file changed, 1 deletion(-) diff --git a/examples/Basic_Simulation/output_reader.py b/examples/Basic_Simulation/output_reader.py index b8cf36a70..f0b08c98c 100644 --- a/examples/Basic_Simulation/output_reader.py +++ b/examples/Basic_Simulation/output_reader.py @@ -42,5 +42,4 @@ ax.set_ylim(0, 0.4) ax.text(1.5, 0.35, f"t = {sim.data['time'].isel(time=0).values} years", size=10, ha="left") plt.tight_layout() -plt.show() fig.savefig("../examples/Basic_Simulation/output.eps", dpi=300, facecolor='white', transparent=False, bbox_inches="tight") From 8814dffd7e4b076bebcbe679ac4be690619d687d Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Thu, 29 Jun 2023 11:07:42 -0400 Subject: [PATCH 09/11] Update initial_conditions.py set a seed for the random number generator, clarified names of bodies, and adjusted eccentricity and inclination of bodies to be less energetic --- .../Basic_Simulation/initial_conditions.py | 35 ++++++++++--------- 1 file changed, 18 insertions(+), 17 deletions(-) diff --git a/examples/Basic_Simulation/initial_conditions.py b/examples/Basic_Simulation/initial_conditions.py index 813bd8dfb..83614fb07 100755 --- a/examples/Basic_Simulation/initial_conditions.py +++ b/examples/Basic_Simulation/initial_conditions.py @@ -39,6 +39,7 @@ #sim = swiftest.Simulation(fragmentation=True, minimum_fragment_mass = 2.5e-11, mtiny=2.5e-8) sim = swiftest.Simulation() sim.clean() +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"]) @@ -47,18 +48,18 @@ npl = 5 density_pl = 3000.0 / (sim.param['MU2KG'] / sim.param['DU2M'] ** 3) -name_pl = ["MassiveBody_01", "MassiveBody_02", "MassiveBody_03", "MassiveBody_04", "MassiveBody_05"] -a_pl = default_rng().uniform(0.3, 1.5, npl) -e_pl = default_rng().uniform(0.0, 0.3, npl) -inc_pl = default_rng().uniform(0.0, 90, npl) -capom_pl = default_rng().uniform(0.0, 360.0, npl) -omega_pl = default_rng().uniform(0.0, 360.0, npl) -capm_pl = default_rng().uniform(0.0, 360.0, npl) -M_pl = np.array([6e23, 8e23, 1e24, 3e24, 5e24]) * sim.KG2MU +name_pl = ["SemiBody_01", "SemiBody_02", "SemiBody_03", "SemiBody_04", "SemiBody_05"] +a_pl = rng.uniform(0.3, 1.5, npl) +e_pl = rng.uniform(0.0, 0.2, npl) +inc_pl = rng.uniform(0.0, 10, npl) +capom_pl = rng.uniform(0.0, 360.0, npl) +omega_pl = rng.uniform(0.0, 360.0, npl) +capm_pl = rng.uniform(0.0, 360.0, npl) +M_pl = np.array([6e20, 8e20, 1e21, 3e21, 5e21]) * sim.KG2MU R_pl = np.full(npl, (3 * M_pl/ (4 * np.pi * density_pl)) ** (1.0 / 3.0)) -Ip_pl = np.full((npl,3),0.4,) -rot_pl = np.zeros((npl,3)) -mtiny = 1.01 * np.max(M_pl) +Ip_pl = np.full((npl,3),0.4,) +rot_pl = np.zeros((npl,3)) +mtiny = 1.1 * np.max(M_pl) sim.add_body(name=name_pl, a=a_pl, e=e_pl, inc=inc_pl, capom=capom_pl, omega=omega_pl, capm=capm_pl, mass=M_pl, radius=R_pl, Ip=Ip_pl, rot=rot_pl) @@ -66,12 +67,12 @@ ntp = 10 name_tp = ["TestParticle_01", "TestParticle_02", "TestParticle_03", "TestParticle_04", "TestParticle_05", "TestParticle_06", "TestParticle_07", "TestParticle_08", "TestParticle_09", "TestParticle_10"] -a_tp = default_rng().uniform(0.3, 1.5, ntp) -e_tp = default_rng().uniform(0.0, 0.3, ntp) -inc_tp = default_rng().uniform(0.0, 90, ntp) -capom_tp = default_rng().uniform(0.0, 360.0, ntp) -omega_tp = default_rng().uniform(0.0, 360.0, ntp) -capm_tp = default_rng().uniform(0.0, 360.0, ntp) +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) 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=1.0e3, dt=0.01, istep_out=100, dump_cadence=0, compute_conservation_values=True, mtiny=mtiny) From 6181c8867ab330fb5331cbaebea081f28940312c Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Thu, 29 Jun 2023 11:09:53 -0400 Subject: [PATCH 10/11] Update output_reader.py adjusted bounds of plot to account for less excited system --- examples/Basic_Simulation/output_reader.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/examples/Basic_Simulation/output_reader.py b/examples/Basic_Simulation/output_reader.py index f0b08c98c..a46b6cefc 100644 --- a/examples/Basic_Simulation/output_reader.py +++ b/examples/Basic_Simulation/output_reader.py @@ -39,7 +39,7 @@ ax.set(xlabel="Semimajor Axis (AU)", ylabel="Eccentricity", title="Simulation Initial Conditions (t=0)") ax.scatter(sim.data['a'].isel(time=0), sim.data['e'].isel(time=0), c=colors, s=sizes, edgecolor='black') ax.set_xlim(0, 2.0) -ax.set_ylim(0, 0.4) -ax.text(1.5, 0.35, f"t = {sim.data['time'].isel(time=0).values} years", size=10, ha="left") +ax.set_ylim(0, 0.25) +ax.text(1.5, 0.2, f"t = {sim.data['time'].isel(time=0).values} years", size=10, ha="left") plt.tight_layout() fig.savefig("../examples/Basic_Simulation/output.eps", dpi=300, facecolor='white', transparent=False, bbox_inches="tight") From 6a82e69cc2e50cc376bcf50897b060a6be58157e Mon Sep 17 00:00:00 2001 From: carlislewishard <70146819+carlislewishard@users.noreply.github.com> Date: Thu, 29 Jun 2023 11:11:10 -0400 Subject: [PATCH 11/11] Update errors.py adjusted energy error limit slightly to account for less energetic system --- examples/Basic_Simulation/errors.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/examples/Basic_Simulation/errors.py b/examples/Basic_Simulation/errors.py index 581677949..5cbdfaa55 100644 --- a/examples/Basic_Simulation/errors.py +++ b/examples/Basic_Simulation/errors.py @@ -52,7 +52,7 @@ # Determine if the errors are within bounds L_limit = 1e-10 -E_limit = 1e-4 +E_limit = 1e-5 GM_limit = 1e-14 lerror = 0 @@ -66,4 +66,4 @@ lerror = 1 print("Mass Error of", GM_final, " higher than threshold value of", GM_limit, ". Test failed.") if (lerror == 0): - print("Test successful.") \ No newline at end of file + print("Test successful.")