animate.py

@@ -50,45 +50,19 @@ def read_output_file(filename):
                 
             name = data['name']
             pos = data['position']
-            
-            # Store position
             positions[name].append((pos[0], pos[1]))
             frame += 1
             times.append(frame)
     
     return positions, times
 
-def center_positions(positions, center_body):
+def create_animation(positions, times, output_file=None):
-    """Center all positions relative to the specified body."""
-    if center_body not in positions:
-        print(f"Warning: Center body '{center_body}' not found in simulation")
-        return positions
-        
-    centered_positions = defaultdict(list)
-    for frame in range(len(next(iter(positions.values())))):
-        # Get center body position for this frame
-        center_pos = positions[center_body][frame]
-        
-        # Shift all positions relative to center body
-        for name, pos_list in positions.items():
-            if frame < len(pos_list):
-                pos = pos_list[frame]
-                centered_pos = (pos[0] - center_pos[0], pos[1] - center_pos[1])
-                centered_positions[name].append(centered_pos)
-    
-    return centered_positions
-
-def create_animation(positions, times, output_file=None, center_body=None):
     """Create an animation of the bodies' orbits."""
     # Check if we have any data
     if not positions or not times:
         print("Error: No valid data found in the input file")
         return
         
-    # Center positions if requested
-    if center_body:
-        positions = center_positions(positions, center_body)
-        
     # Set up the figure and axis
     fig, ax = plt.subplots(figsize=(10, 10))
     
@@ -103,10 +77,7 @@ def create_animation(positions, times, output_file=None, center_body=None):
     # Set up the plot
     ax.set_xlabel('X (m)')
     ax.set_ylabel('Y (m)')
-    title = 'Orbital Simulation'
+    ax.set_title('Orbital Simulation')
-    if center_body:
-        title += f' (centered on {center_body})'
-    ax.set_title(title)
     ax.legend()
     
     # Find the bounds of the plot
@@ -155,11 +126,10 @@ def main():
     parser = argparse.ArgumentParser(description='Animate orbital simulation output')
     parser.add_argument('input_file', help='Input file from simulation')
     parser.add_argument('--output', '-o', help='Output video file (optional)')
-    parser.add_argument('--center', '-c', help='Center the animation on this body')
     args = parser.parse_args()
     
     positions, times = read_output_file(args.input_file)
-    create_animation(positions, times, args.output, args.center)
+    create_animation(positions, times, args.output)
 
 if __name__ == '__main__':
     main() 

config/planets.json

@@ -1,28 +1,16 @@
 {
-    "bodies": [
+    "planets": [
-        {
-            "name": "Mercury",
-            "mass": 3.30104e23,
-            "position": [4.6000e10, 0, 0],
-            "velocity": [0, 58970, 0]
-        },
-        {
-            "name": "Venus",
-            "mass": 4.867e24,
-            "position": [1.08941e11, 0, 0],
-            "velocity": [0, 34780, 0]
-        },
         {
             "name": "Earth",
             "mass": 5.972e24,
-            "position": [1.47095e11, 0, 0],
+            "position": [149597870700, 0, 0],
-            "velocity": [0, 29290, 0]
+            "velocity": [0, 29780, 0]
         },
         {
             "name": "Moon",
             "mass": 7.34767309e22,
             "position": [149982270700, 0, 0],
-            "velocity": [0, 30822, 0]
+            "velocity": [0, 30802, 0]
         },
         {
             "name": "Sun",

src/body.cpp

@@ -3,6 +3,7 @@
 #include <cmath>
 #include <ranges>
 
+
 Body::Body(const Position& X, const Velocity& V, const Mass& m, const std::string& name)
     : X(X), V(V), m(m), name(name) {
     A = Acceleration{Decimal(0), Decimal(0), Decimal(0)};
@@ -14,12 +15,6 @@ void Body::addAcceleration(const Acceleration& new_A) {
     }
 }
 
-void Body::subAcceleration(const Acceleration& new_A) {
-    for (const auto& [new_a, cur_a]: std::views::zip(new_A, this->A)){
-        cur_a -= new_a;
-    }
-}
-
 std::tuple<Position, Velocity, Mass> Body::save() const {
     return std::make_tuple(X, V, m);
 }
@@ -71,8 +66,8 @@ std::string Body::speed() const {
 std::string Body::toString() const {
     std::string pos_str, vel_str;
     for (int i = 0; i < 3; ++i) {
-        pos_str += format_sig_figs(real_pos(X[i]), 10) + "m ";
+        pos_str += format_sig_figs(real_pos(X[i]), 3) + "m ";
-        vel_str += format_sig_figs(real_vel(V[i]), 10) + "m/s ";
+        vel_str += format_sig_figs(real_vel(V[i]), 3) + "m/s ";
     }
     return name + ": X = " + pos_str + ", V = " + vel_str;
 } 

src/body.hpp

@@ -3,12 +3,10 @@
 #include "units.hpp"
 #include <string>
 #include <tuple>
-#include <memory>
 
 class Body {
 public:
-    Body(const Position& X, const Velocity& V, const Mass& m, 
+    Body(const Position& X, const Velocity& V, const Mass& m, const std::string& name = "");
-         const std::string& name = "");
     
     // Save and load state
     std::tuple<Position, Velocity, Mass> save() const;
@@ -33,7 +31,7 @@ public:
     void setAcceleration(const Acceleration& new_A) { A = new_A; }
     
     void addAcceleration(const Acceleration& new_A);
-    void subAcceleration(const Acceleration& new_A);
+    
     // String representation
     std::string toString() const;
 

src/main.cpp

@@ -52,7 +52,7 @@ SimulationConfig parse_command_line(int argc, char* argv[]) {
     desc.add_options()
         ("help,h", "Show help message")
         ("config,c", po::value<std::string>(&config.config_file)->required(), 
-         "Path to body configuration file (JSON)")
+         "Path to planet configuration file (JSON)")
         ("output,o", po::value<std::string>(&config.output_file)->required(), 
          "Path to output file")
         ("time,t", po::value<std::string>()->required(),
@@ -97,7 +97,7 @@ SimulationConfig parse_command_line(int argc, char* argv[]) {
     return config;
 }
 
-std::vector<Body> load_bodies(const std::string& config_file) {
+std::vector<Body> load_planets(const std::string& config_file) {
     std::ifstream f(config_file);
     if (!f.is_open()) {
         throw std::runtime_error("Could not open config file: " + config_file);
@@ -107,12 +107,12 @@ std::vector<Body> load_bodies(const std::string& config_file) {
     f >> j;
     
     std::vector<Body> bodies;
-    for (const auto& body : j["bodies"]) {
+    for (const auto& planet : j["planets"]) {
-        std::string name = body["name"];
+        std::string name = planet["name"];
-        double mass = body["mass"];
+        double mass = planet["mass"];
         
-        std::vector<double> pos = body["position"];
+        std::vector<double> pos = planet["position"];
-        std::vector<double> vel = body["velocity"];
+        std::vector<double> vel = planet["velocity"];
         
         if (pos.size() != 3 || vel.size() != 3) {
             throw std::runtime_error("Position and velocity must be 3D vectors");
@@ -124,7 +124,6 @@ std::vector<Body> load_bodies(const std::string& config_file) {
         Mass normalized_mass = norm_mass(mass);
         
         bodies.emplace_back(position, velocity, normalized_mass, name);
-        
         std::cout << "Loaded " << name << " with mass " << mass << " kg\n";
     }
     
@@ -136,8 +135,8 @@ int main(int argc, char* argv[]) {
         // Parse command line options
         auto config = parse_command_line(argc, argv);
         
-        // Load bodies from config file
+        // Load planets from config file
-        auto bodies = load_bodies(config.config_file);
+        auto bodies = load_planets(config.config_file);
         
         // Create and run simulator
         Simulator simulator(

src/simulator.cpp

@@ -11,11 +11,11 @@ Simulator::Simulator(const std::vector<Body>& bodies,
                      const std::filesystem::path& output_file,
                      int current_step,
                      bool overwrite_output)
-    : bodies(bodies),
+    : bodies(bodies)
-      step_size(step_size),
+    , step_size(step_size)
-      steps_per_save(steps_per_save),
+    , steps_per_save(steps_per_save)
-      output_file(output_file),
+    , output_file(output_file)
-      current_step(current_step) {
+    , current_step(current_step) {
     
     if (std::filesystem::exists(output_file) && !overwrite_output) {
         throw std::runtime_error("File " + output_file.string() + " exists and overwrite flag not given.");
@@ -28,14 +28,14 @@ Simulator::Simulator(const std::vector<Body>& bodies,
         clear.close();
     }
     
+    // Write initial header with masses
     out.open(output_file, std::ios::app);
     if (!out) {
         throw std::runtime_error("Failed to open output file: " + output_file.string());
     }
     
     for (const auto& body : bodies) {
-        out << body.getName() << ": mass=" << body.getMass();
+        out << body.getName() << ": " << body.getMass() << "\n";
-        out << "\n";
     }
     out << "\n";
     
@@ -80,9 +80,10 @@ void Simulator::calculate_forces() {
         for (size_t j = i + 1; j < bodies.size(); j++) {
             Position vec;
             for (int k = 0; k < 3; ++k) {
-                vec[k] = positions[i][k] - positions[j][k];
+                vec[k] = bodies[i].getPosition()[k] - bodies[j].getPosition()[k];
             }
             
+            // Calculate distance
             Decimal dist = std::sqrt(vec[0] * vec[0] + vec[1] * vec[1] + vec[2] * vec[2]);
             
             // Calculate force magnitude using Newton's law of gravitation
@@ -90,6 +91,7 @@ void Simulator::calculate_forces() {
             // for the acceleration. Use r^3 to avoid normalizing vec when multiplying later
             Decimal force_magnitude = 1 / (dist * dist * dist); 
             
+            // Calculate acceleration for both bodies
             Decimal acc_magnitude_i = force_magnitude * bodies[j].getMass();
             Decimal acc_magnitude_j = force_magnitude * bodies[i].getMass();
 
@@ -100,11 +102,11 @@ void Simulator::calculate_forces() {
                 acc_j[k] = vec[k] * acc_magnitude_j;   
             }
             
+            // Add this component to the two bodies
             bodies[i].addAcceleration(acc_i);
             bodies[j].addAcceleration(acc_j);
         }
     }
-
 }
 
 void Simulator::move_bodies() {
@@ -117,5 +119,6 @@ void Simulator::checkpoint() {
     for (const auto& body : bodies) {
         out << body.toString() << "\n";
     }
+    out << "\n";
     out.flush();
 } 

tests/simulator_test.cpp

@@ -1,54 +1,54 @@
-// #include <gtest/gtest.h>
+#include <gtest/gtest.h>
-// #include "simulator.hpp"
+#include "simulator.hpp"
-// #include <vector>
+#include <vector>
-// #include <filesystem>
+#include <filesystem>
-// #include <memory>
+#include <memory>
 
-// class SimulatorTest : public testing::Test {
+class SimulatorTest : public testing::Test {
-// protected:
+protected:
-//     SimulatorTest() : 
+    SimulatorTest() : 
-//         bodies({
+        bodies({
-//             Body(
+            Body(
-//                 Position{Decimal(0), Decimal(0), Decimal(0)},
+                Position{Decimal(0), Decimal(0), Decimal(0)},
-//                 Velocity{Decimal(0), Decimal(0), Decimal(0)},
+                Velocity{Decimal(0), Decimal(0), Decimal(0)},
-//                 Mass(1.0),
+                Mass(1.0),
-//                 "body1"
+                "body1"
-//             ),
+            ),
-//             Body(
+            Body(
-//                 Position{Decimal(1), Decimal(0), Decimal(0)},
+                Position{Decimal(1), Decimal(0), Decimal(0)},
-//                 Velocity{Decimal(0), Decimal(1), Decimal(0)},
+                Velocity{Decimal(0), Decimal(1), Decimal(0)},
-//                 Mass(1.0),
+                Mass(1.0),
-//                 "body2"
+                "body2"
-//             )
+            )
-//         })
+        })
-//     {
+    {
-//         simulator = std::make_unique<Simulator>(
+        simulator = std::make_unique<Simulator>(
-//             bodies,
+            bodies,
-//             Decimal(0.1),  // step_size
+            Decimal(0.1),  // step_size
-//             1,            // steps_per_save
+            1,            // steps_per_save
-//             std::filesystem::path("test_output.json"),  // output_file
+            std::filesystem::path("test_output.json"),  // output_file
-//             0,            // current_step
+            0,            // current_step
-//             true          // overwrite_output
+            true          // overwrite_output
-//         );
+        );
-//     }
+    }
 
-//     std::vector<Body> bodies;
+    std::vector<Body> bodies;
-//     std::unique_ptr<Simulator> simulator;
+    std::unique_ptr<Simulator> simulator;
-// };
+};
 
-// TEST_F(SimulatorTest, Step) {
+TEST_F(SimulatorTest, Step) {
-//     // Take a step
+    // Take a step
-//     //TODO
+    //TODO
-// }
+}
 
-// TEST_F(SimulatorTest, TotalEnergy) {
+TEST_F(SimulatorTest, TotalEnergy) {
-//     //TODO
+    //TODO
-//     // Decimal initial_energy = simulator->total_energy();
+    // Decimal initial_energy = simulator->total_energy();
 
-//     // simulator->step(Decimal(0.1));
+    // simulator->step(Decimal(0.1));
 
-//     // Decimal new_energy = simulator->total_energy();
+    // Decimal new_energy = simulator->total_energy();
 
-//     // EXPECT_NEAR(initial_energy, new_energy, 1e-10);
+    // EXPECT_NEAR(initial_energy, new_energy, 1e-10);
-// }
+}