Project Goal

This project aims to model the Hydrogen atom using two primary physical mechanisms: the Coulomb force to represent the electrostatic potential between the proton and electron, and Larmor radiation to account for energy loss due to electromagnetic radiation. By constructing this simulation, I hope to compare the results to those of previous contributors who used different software platforms for similar analyses.

Getting Started: The Two-Body Problem

To build foundational understanding, I began by modeling the classic two-body problem. Rather than jumping directly into the hydrogen atom, I chose to first simulate a gravitational system involving two general masses. This served as a conceptual stepping stone and helped me get familiar with COMSOL’s capabilities.

While researching, I came across an online example simulating a satellite orbiting Earth—complete with an mph file. This example, created by Walter Frei (a well-known contributor to COMSOL’s learning resources), became a valuable reference for me. It was also my first experience working with the Global ODEs and DAEs interface in COMSOL.

Tutorial Resource:
Walter Frei’s COMSOL tutorial – a helpful guide for beginners and advanced users alike.
(You can download the model from the bottom of his webpage.)

As a test of my understanding, I created a simple orbital simulation: a randomly shaped object orbiting Earth at a radius roughly equivalent to the Moon’s distance. The force governing the motion was gravitational, and the model was built using with the Global ODEs and DAEs physics. This exercise helped me get familiar with setting up the mechanical equations, defining initial conditions, and choosing appropriate solver configurations. It served as an early confidence to building step before transitioning to the more complex task of modeling atomic-scale systems.

This write up was created to test my knowledge and confirm I knew the parts to the COMSOL

These were 2 orbits that I created

From Orbital Mechanics to Atomic Modeling

After refining the initial setup, I began modifying the conditions to more closely reflect a true two-body system. These adjustments included simplifying the model to remove unnecessary complexities, such as matrix definitions used in satellite-orbit problems. Unlike Walter Frei’s satellite model—which required transformation into the satellite’s reference frame—this version remains entirely in an inertial frame, which is sufficient for modeling central forces like gravity or Coulomb attraction.

The updated COMSOL model, available in the accompanying .mph file, reflects these changes and serves as a clean base from which to pivot toward atomic-scale physics.

Hydrogen Atom

The next step was straightforward: Adjust the masses to match those of the proton and electron. Replace the gravitational force with the Coulomb force (electrostatic attraction). Set an appropriate initial radius and velocity to correspond to a known energy level of the hydrogen atom.

With these changes, I created a new simulation that approximates the hydrogen atom classically. The model shown here corresponds to quantum level n = 3, and at this stage, radiation loss is not included.

The goal for this setup was to test energy conservation. Ideally, since radiation is disabled, the total energy should remain constant over time. Any small drift in energy is due to numerical integration error, which I minimized by using the Classical Runge-Kutta (RK4) Integrator with a time-step of 1e-3 and of order 4, as in this mph file. This configuration provides a stable baseline for testing the effects of adding radiation later.

Physical Potential Well Approach

Another way to visualize the hydrogen atom is to represent the Coulomb potential, V(r)= -1/r, as a physical surface on which a sphere can roll. With carefully chosen initial conditions, the motion of the sphere can resemble classical orbits of the electron. To reach this point, I first explored a simpler system to better understand how COMSOL handles surface-to-surface contact mechanics.

I started with a 2D parabolic well. This setup was simpler to construct and served as a great testing ground. To make the simulation work, I included elements like contact pairs, gravity, and material properties. This was my first time applying COMSOL’s multibody dynamics features to surface motion, and it helped me get familiar with the settings required for proper physical interaction.

A key test for this setup was checking how much energy the system lost when no energy-loss physics (like friction or damping) were present. In theory, the ball should return to the same height after each swing, since no forces were acting to remove energy from the system. As seen in the accompanying GIF, the ball doesn’t return to exactly the same height, but it comes close. The small energy drift reflects the accuracy limits of the integrator. After some experimentation, we found that using a manual BDF integrator with a time step of 1e-3 and order 1 gave us consistent and stable results.

After figuring this out, we quickly moved onto the -1/r potential well in 3D. This process took longer, because everything 3D takes longer, as well as the mesh needed to be fine in order to stop any bouncing of the ball, but that also increased computation time.

3d -1/r Well

 

Radiation

Now we wanted to figure out how to model the radiation that occurs when an electron moves between energy levels. There were many different ideas and methods, we will include write ups about them all.

 

A Force-Based Approach

To introduce energy loss into the hydrogen atom model, we turned to a paper that presents a classical framework for photon emission:

Marian Kowalski,
Classical Description of Photon Emission from Atomic Hydrogen,
Physics Essays, Volume 12, Number 2, 1999.

Kowalski’s approach expresses radiation loss as a force term of the form: Frad=μγv 

where μ is the reduced mass of the system, γ is an empirical parameter associated with radiation, and r is the position vector. This term acts as a kind of “damping” that pulls the electron inward, simulating the radiation-induced loss of energy over time.

ODE for the x position, with corresponding ones in the y and z position.

mu * rxtt + (rx)/ rmag^3 + M2 * G1* rxt = 0

I added this force term directly to the ODEs governing the Coulomb interaction in the hydrogen atom model. Once implemented, the simulation produced behavior that aligned fairly well with Kowalski’s results.

The model shows expected behavior: both energy and angular momentum decrease gradually over time. For example, when starting at n=2, the orbit transitions toward $1$, and the angular momentum (represented by ℓ) drops from 2 to 1, consistent with quantum predictions in a classical analog framework.

However, the decay isn’t a perfect match. The simulated trajectory tends to overshoot the final energy level slightly. We’re currently unsure whether this is due to differences in the radiation parameter γ, or if it stems from numerical inaccuracies in the integrator. To better control this, I implemented a stop condition in the solver that halts the simulation once a target energy is reached:

E_goal=−1/2  * 1/ ((n−1)^ 2), comp1.minop1(comp1.E) < E_Goal

This is what stops the computation in the video above, and as you see, the energy finally reaches the goal slightly past the (1,0) as excepted. That’s where the issue lies.