diff --git a/Binary/README.md b/Binary/README.md index d69a6ab..229d0d8 100644 --- a/Binary/README.md +++ b/Binary/README.md @@ -1,7 +1,25 @@ +# Binary star example + A simple simulation of a binary star system, to exercise the code and give some experience with the file formats. The example takes two stars and evolves them under direct interaction. The initial conditions for the stars are specified in `binary.bods`, which follows standard `EXP` input format: the first line is the number of bodies (followed by two integers that specify the number of extra integer and double fields -- here set to be zero). Lines after the first line specify, for each body, the mass, xyz position, and xyz velocity. +## How to run + Run as follows. The YAML configuration file is set for the Docker container. If you are using a native build and installation, change the ldlibdir parameter to point at your library install directory. Then, execute the command: `mpirun -np 1 exp binary.yml`. +## A fixed potential example + There is also an example where the stars no longer feel mutual self-gravity, but instead are evolved in an external logarithmic potential. This may be run using `mpirun -np 1 exp logpot.yml`. This example demonstrates `noforce`. +## Viewing results + +Each example creates a trajectory file named `ORBTRACE.runX` where `X` +is either `1` or `2`, for `binary.yml` or `logpot.yml` +respectively. These may be plotted using the command: + +``` bash +python3 plotOrbit.py --input=ORBTRACE.runX +``` + +## Additional exploration + A first exercise to extend these examples would be to instead use the `usermw` module to simulate the orbit of the sun in the MW. diff --git a/Binary/logpot.yml b/Binary/logpot.yml index 3af6f2e..80ac0f0 100755 --- a/Binary/logpot.yml +++ b/Binary/logpot.yml @@ -15,7 +15,7 @@ Global: multistep : 12 dynfracA : 0.01 dynfracV : 0.01 - VERBOSE : 4 + VERBOSE : 1 # ------------------------------------------------------------------------ # This is a sequence of components. The parameters for the force are # now included as a parameter map, rather than a separate file. diff --git a/Binary/plotOrbit.py b/Binary/plotOrbit.py new file mode 100644 index 0000000..3f2c7e9 --- /dev/null +++ b/Binary/plotOrbit.py @@ -0,0 +1,57 @@ +import matplotlib.pyplot as plt +import numpy as np +import argparse + +def plot_orbit(positions, title="Orbital Trajectory", save_path=None): + """ + Plots two orbital trajectories given a list of 3D positions. + + Parameters: + positions: The x, y, z positions for N times both trajectories as an ndarray of shape (N, 6). + title (str): The title of the plot. + save_path (str): If provided, the plot will be saved to this path. + """ + + # Extract x, y, z coordinates + x1 = positions[:, 0] + y1 = positions[:, 1] + z1 = positions[:, 2] + x2 = positions[:, 3] + y2 = positions[:, 4] + z2 = positions[:, 5] + + # Create a 3D plot + fig = plt.figure() + ax = fig.add_subplot(111, projection='3d') + + # Plot the trajectory + ax.plot(x1, y1, z1, '-') + ax.plot(x2, y2, z2, '-') + + # Set labels and title + ax.set_xlabel('X Position') + ax.set_ylabel('Y Position') + ax.set_zlabel('Z Position') + ax.set_title(title) + + # Show grid + ax.grid(True) + + # Save or show the plot + if save_path: + plt.savefig(save_path) + else: + plt.show() + + +if __name__ == "__main__": + parser = argparse.ArgumentParser(description="Plot orbital trajectory from 3D positions.") + parser.add_argument('--input', type=str, + default='ORBTRACE.run0', + help="Path to the input file containing 3D positions (one per line as x,y,z).") + + args = parser.parse_args() + + data = np.loadtxt(args.input) + + plot_orbit(data[:, [1,2,3,11,12,13]], title="Binary orbit") diff --git a/DiskHalo/config.yml b/DiskHalo/config.yml index 0c11780..598dafc 100644 --- a/DiskHalo/config.yml +++ b/DiskHalo/config.yml @@ -39,7 +39,7 @@ Components: bodyfile : disk.bod force : id : cylinder - parameters : {acyl: 1.0, nmaxfid: 32, mmax: 2, hcyl: 0.05, nmax: 10, ncylodd: 2, ncylnx: 128, ncylny: 64, cachename: disk.cache} + parameters : {acyl: 1.0, lmaxfid: 64, mmax: 2, hcyl: 0.05, sech2: false, nmax: 10, ncylodd: 2, ncylnx: 200, ncylny: 100, cachename: disk.cache} # # The parameters could be expressed like this, equivalently: # diff --git a/DiskHalo/disk.cache b/DiskHalo/disk.cache index 1c6f134..d681cf5 100644 Binary files a/DiskHalo/disk.cache and b/DiskHalo/disk.cache differ diff --git a/Nbody/README.md b/Nbody/README.md index 927beec..6537a90 100644 --- a/Nbody/README.md +++ b/Nbody/README.md @@ -48,3 +48,13 @@ testing. The directory `data` contains some snapshots and coefficients for use with the pyEXP tutorials. + +## Viewing results + +This example creates a trajectory file with selected disk orbits named +`ORBTRACE.run0`. These orbits may be plotted using the command: + +``` bash +python3 plotOrbit.py +``` + diff --git a/Nbody/plotOrbit.py b/Nbody/plotOrbit.py new file mode 100644 index 0000000..0647a97 --- /dev/null +++ b/Nbody/plotOrbit.py @@ -0,0 +1,56 @@ +import matplotlib.pyplot as plt +import numpy as np +import argparse + +def plot_orbit(positions, title="Orbital Trajectory", save_path=None): + """ + Plots two orbital trajectories given a list of 3D positions. + + Parameters: + positions: A list of 3d particle trajectories (x, y, z) positions for N times, as a list of M ndarrays of shape (N, 3). + title (str): The title of the plot. + save_path (str): If provided, the plot will be saved to this path. + """ + + # Create a 3D plot + fig = plt.figure() + ax = fig.add_subplot(111, projection='3d') + + # Plot the trajectories + for n in range(len(positions)): + ax.plot(positions[n][:, 0], positions[n][:, 1], positions[n][:, 2], '-') + + # Set labels and title + ax.set_xlabel('X Position') + ax.set_ylabel('Y Position') + ax.set_zlabel('Z Position') + ax.set_title(title) + + # Show grid + ax.grid(True) + + # Save or show the plot + if save_path: + plt.savefig(save_path) + else: + plt.show() + + +if __name__ == "__main__": + parser = argparse.ArgumentParser(description="Plot orbital trajectory from 3D positions.") + parser.add_argument('--input', type=str, + default='ORBTRACE.run0', + help="Path to the input file containing 3D positions (one per line as x,y,z).") + + args = parser.parse_args() + + data = np.loadtxt(args.input) + + pos = [] + pos.append(data[:, 1:4]) # First object positions + pos.append(data[:, 7:10]) # Second object positions + pos.append(data[:, 13:16]) # Third object positions + pos.append(data[:, 19:22]) # Fourth object positions + pos.append(data[:, 25:28]) # Fifth object positions + + plot_orbit(pos, title="Selected disk orbits")