Introduction to Computational Physics with CoCalc

Bridging Theory and Computation in Physics

Transforming Physics Understanding Through Computation

Physics—the study of matter, energy, and their interactions—has always been deeply mathematical. Today, computational methods have become essential tools for understanding complex physical phenomena that resist analytical solutions. CoCalc provides an ideal environment for learning physics through the powerful combination of theoretical understanding and computational exploration.

For information about available scientific computing tools and environments, see the CoCalc documentation.

Whether you're modeling planetary motion, analyzing quantum systems, or exploring electromagnetic fields, CoCalc's integrated tools help you visualize, simulate, and understand the physical world in ways that traditional methods alone cannot achieve.

Your Computational Physics Toolkit

Python for Physics: The Foundation

Python has become the lingua franca of computational physics, offering powerful libraries and intuitive syntax:

# Essential physics imports
import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate, optimize
import sympy as sp

# Welcome to computational physics!
print("Welcome to Computational Physics with CoCalc!")

# Physical constants (in SI units)
c = 299792458 # Speed of light (m/s)
h = 6.62607015e-34 # Planck constant (J⋅s)
hbar = h / (2 * np.pi) # Reduced Planck constant
k_B = 1.380649e-23 # Boltzmann constant (J/K)
e = 1.602176634e-19 # Elementary charge (C)
m_e = 9.1093837015e-31 # Electron mass (kg)
m_p = 1.67262192369e-27 # Proton mass (kg)

print(f"Speed of light: {c:.0e} m/s")
print(f"Planck constant: {h:.3e} J⋅s")
print(f"Electron mass: {m_e:.3e} kg")

Mechanics: Motion and Forces

Start your physics journey with classical mechanics:

# Classical mechanics: projectile motion
def projectile_motion():
"""
Simulate projectile motion with air resistance
"""

# Parameters
g = 9.81 # Gravitational acceleration (m/s²)
v0 = 50 # Initial velocity (m/s)
angle = 45 # Launch angle (degrees)
m = 1.0 # Mass (kg)
b = 0.1 # Air resistance coefficient

# Convert angle to radians
theta = np.radians(angle)

# Initial conditions
vx0 = v0 * np.cos(theta)
vy0 = v0 * np.sin(theta)

# Differential equation: F = ma = -mg - bv
def equations_of_motion(t, state):
x, y, vx, vy = state

# Air resistance force
v_magnitude = np.sqrt(vx**2 + vy**2)
if v_magnitude > 0:
ax = -b * vx * v_magnitude / m
ay = -g - b * vy * v_magnitude / m
else:
ax = 0
ay = -g

return [vx, vy, ax, ay]

# Solve the differential equation
t_span = (0, 10)
initial_state = [0, 0, vx0, vy0]

# Event function to stop when projectile hits ground
def hit_ground(t, state):
return state[1] # y-coordinate
hit_ground.terminal = True
hit_ground.direction = -1

solution = integrate.solve_ivp(
equations_of_motion, t_span, initial_state,
events=hit_ground, dense_output=True, rtol=1e-8
)

# Extract trajectory
t_flight = solution.t
trajectory = solution.y

# Plot trajectory
plt.figure(figsize=(10, 6))
plt.plot(trajectory[0], trajectory[1], 'b-', linewidth=2, label='With air resistance')

# Compare with no air resistance (analytical solution)
t_no_air = np.linspace(0, t_flight[-1], 100)
x_no_air = vx0 * t_no_air
y_no_air = vy0 * t_no_air - 0.5 * g * t_no_air**2

plt.plot(x_no_air, y_no_air, 'r--', linewidth=2, label='No air resistance')

plt.xlabel('Horizontal Distance (m)')
plt.ylabel('Height (m)')
plt.title(f'Projectile Motion (v₀={v0} m/s, θ={angle}°)')
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

# Calculate range
range_with_air = trajectory[0][-1]
range_no_air = v0**2 * np.sin(2*theta) / g

print(f"Range with air resistance: {range_with_air:.1f} m")
print(f"Range without air resistance: {range_no_air:.1f} m")
print(f"Air resistance reduces range by: {(1 - range_with_air/range_no_air)*100:.1f}%")

return t_flight, trajectory

# Execute projectile motion simulation
time, trajectory = projectile_motion()

Oscillations and Waves

Explore periodic motion and wave phenomena:

# Simple harmonic motion and damped oscillations
def harmonic_oscillator():
"""
Study simple harmonic motion and damping effects
"""

# Parameters
omega_0 = 2.0 # Natural frequency (rad/s)
gamma = 0.1 # Damping coefficient
F0 = 1.0 # Driving force amplitude
omega_d = 1.8 # Driving frequency

# Equation of motion: m*x'' + gamma*x' + k*x = F*cos(omega_d*t)
# Using omega_0^2 = k/m, we get: x'' + gamma*x' + omega_0^2*x = (F0/m)*cos(omega_d*t)

def oscillator_equation(t, state):
x, v = state

# Driving force
driving_force = F0 * np.cos(omega_d * t)

# Acceleration
a = -2*gamma*v - omega_0**2*x + driving_force

return [v, a]

# Time array
t = np.linspace(0, 20, 1000)

# Solve for different initial conditions
solutions = {}
initial_conditions = [
(1.0, 0.0, "x₀=1, v₀=0"),
(0.0, 2.0, "x₀=0, v₀=2"),
(0.5, 1.0, "x₀=0.5, v₀=1")
]

plt.figure(figsize=(12, 8))

for i, (x0, v0, label) in enumerate(initial_conditions):
solution = integrate.solve_ivp(
oscillator_equation, (0, 20), [x0, v0],
t_eval=t, rtol=1e-8
)

solutions[label] = solution

plt.subplot(2, 2, i+1)
plt.plot(solution.t, solution.y[0], 'b-', linewidth=1.5)
plt.xlabel('Time (s)')
plt.ylabel('Position (m)')
plt.title(f'Damped Oscillator: {label}')
plt.grid(True, alpha=0.3)

# Phase space plot
plt.subplot(2, 2, 4)
for label, solution in solutions.items():
plt.plot(solution.y[0], solution.y[1], linewidth=1.5, label=label)

plt.xlabel('Position (m)')
plt.ylabel('Velocity (m/s)')
plt.title('Phase Space')
plt.legend()
plt.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

# Analyze frequency response
frequencies = np.linspace(0.1, 4.0, 100)
amplitude_response = []

for omega in frequencies:
# Steady-state amplitude for driven oscillator
amplitude = F0 / np.sqrt((omega_0**2 - omega**2)**2 + (2*gamma*omega)**2)
amplitude_response.append(amplitude)

plt.figure(figsize=(10, 6))
plt.plot(frequencies, amplitude_response, 'r-', linewidth=2)
plt.axvline(omega_0, color='b', linestyle='--', alpha=0.7, label=f'ω₀ = {omega_0}')
plt.xlabel('Driving Frequency (rad/s)')
plt.ylabel('Steady-State Amplitude')
plt.title('Frequency Response of Damped Harmonic Oscillator')
plt.legend()
plt.grid(True, alpha=0.3)
plt.show()

return solutions, amplitude_response

# Execute oscillator analysis
oscillator_solutions, frequency_response = harmonic_oscillator()

Electromagnetism: Fields and Forces

Explore electric and magnetic phenomena:

# Electromagnetic field visualization
def electromagnetic_fields():
"""
Visualize electric and magnetic fields
"""

# Electric field of point charges
def electric_field_point_charges():
"""Visualize electric field of multiple point charges"""

# Define point charges: (x, y, charge)
charges = [
(2, 2, 1.0), # Positive charge
(-2, -2, -1.0), # Negative charge
(0, 3, 0.5), # Smaller positive charge
]

# Create grid
x = np.linspace(-4, 4, 20)
y = np.linspace(-4, 4, 20)
X, Y = np.meshgrid(x, y)

# Calculate electric field
Ex = np.zeros_like(X)
Ey = np.zeros_like(Y)
V = np.zeros_like(X) # Electric potential

k = 8.99e9 # Coulomb constant (N⋅m²/C²)

for charge_x, charge_y, q in charges:
# Distance from each grid point to the charge
dx = X - charge_x
dy = Y - charge_y
r = np.sqrt(dx**2 + dy**2)

# Avoid division by zero
r = np.where(r < 0.1, 0.1, r)

# Electric field components
Ex += k * q * dx / r**3
Ey += k * q * dy / r**3

# Electric potential
V += k * q / r

# Plot electric field and equipotential lines
plt.figure(figsize=(12, 5))

# Electric field vectors
plt.subplot(1, 2, 1)
plt.quiver(X, Y, Ex, Ey, np.sqrt(Ex**2 + Ey**2), cmap='viridis', alpha=0.7)

# Mark charge locations
for charge_x, charge_y, q in charges:
color = 'red' if q > 0 else 'blue'
size = abs(q) * 100
plt.scatter(charge_x, charge_y, c=color, s=size, edgecolors='black')

plt.xlabel('x (m)')
plt.ylabel('y (m)')
plt.title('Electric Field')
plt.axis('equal')
plt.colorbar(label='Field Magnitude')

# Equipotential lines
plt.subplot(1, 2, 2)
contour = plt.contour(X, Y, V, levels=20, colors='blue', alpha=0.6)
plt.clabel(contour, inline=True, fontsize=8)

# Mark charge locations
for charge_x, charge_y, q in charges:
color = 'red' if q > 0 else 'blue'
size = abs(q) * 100
plt.scatter(charge_x, charge_y, c=color, s=size, edgecolors='black')

plt.xlabel('x (m)')
plt.ylabel('y (m)')
plt.title('Electric Potential')
plt.axis('equal')

plt.tight_layout()
plt.show()

return X, Y, Ex, Ey, V

# Magnetic field of current loop
def magnetic_field_current_loop():
"""Magnetic field of a circular current loop"""

# Parameters
I = 1.0 # Current (A)
R = 1.0 # Loop radius (m)
mu_0 = 4*np.pi*1e-7 # Permeability of free space

# Calculate field along axis
z = np.linspace(-3, 3, 100)
Bz = mu_0 * I * R**2 / (2 * (R**2 + z**2)**(3/2))

plt.figure(figsize=(10, 6))

plt.subplot(1, 2, 1)
plt.plot(z, Bz*1e6, 'b-', linewidth=2)
plt.xlabel('Distance along axis (m)')
plt.ylabel('Magnetic Field (μT)')
plt.title(f'Magnetic Field of Current Loop (I={I}A, R={R}m)')
plt.grid(True, alpha=0.3)

# Field lines visualization (simplified)
plt.subplot(1, 2, 2)
theta = np.linspace(0, 2*np.pi, 100)
loop_x = R * np.cos(theta)
loop_y = R * np.sin(theta)

# Draw current loop
plt.plot(loop_x, loop_y, 'r-', linewidth=3, label='Current Loop')

# Simplified field line representation
r_field = np.linspace(0.2, 3, 10)
for r in r_field:
# Approximate field lines
phi = np.linspace(0, 2*np.pi, 50)
field_x = r * np.cos(phi)
field_y = 0.5 * r * np.sin(phi)
plt.plot(field_x, field_y, 'b-', alpha=0.5, linewidth=1)
plt.plot(field_x, -field_y, 'b-', alpha=0.5, linewidth=1)

plt.xlabel('x (m)')
plt.ylabel('y (m)')
plt.title('Magnetic Field Lines')
plt.axis('equal')
plt.legend()

plt.tight_layout()
plt.show()

return z, Bz

# Execute electromagnetic field calculations
electric_results = electric_field_point_charges()
magnetic_results = magnetic_field_current_loop()

return electric_results, magnetic_results

# Execute electromagnetic field analysis
em_results = electromagnetic_fields()

Thermodynamics and Statistical Mechanics

Explore thermal phenomena and statistical behavior:

# Thermodynamics and kinetic theory
def thermal_physics():
"""
Explore thermodynamics and statistical mechanics
"""

# Ideal gas law and kinetic theory
def kinetic_theory_simulation():
"""Simulate kinetic theory of gases"""

# Parameters
N = 1000 # Number of particles
T = 300 # Temperature (K)
m = 4.65e-26 # Mass of N2 molecule (kg)

# Maxwell-Boltzmann distribution
def maxwell_boltzmann(v, T, m):
"""Maxwell-Boltzmann speed distribution"""
k_B = 1.38e-23
factor = 4 * np.pi * (m / (2 * np.pi * k_B * T))**(3/2)
return factor * v**2 * np.exp(-m * v**2 / (2 * k_B * T))

# Generate random velocities according to Maxwell-Boltzmann
# (Simplified: using normal distribution for each component)
sigma = np.sqrt(k_B * T / m) # Standard deviation for each component

vx = np.random.normal(0, sigma, N)
vy = np.random.normal(0, sigma, N)
vz = np.random.normal(0, sigma, N)

# Calculate speeds
speeds = np.sqrt(vx**2 + vy**2 + vz**2)

# Plot speed distribution
plt.figure(figsize=(12, 8))

plt.subplot(2, 2, 1)
plt.hist(speeds, bins=50, density=True, alpha=0.7, label='Simulation')

# Theoretical Maxwell-Boltzmann distribution
v_theory = np.linspace(0, np.max(speeds), 200)
mb_theory = maxwell_boltzmann(v_theory, T, m)
plt.plot(v_theory, mb_theory, 'r-', linewidth=2, label='Theory')

plt.xlabel('Speed (m/s)')
plt.ylabel('Probability Density')
plt.title(f'Maxwell-Boltzmann Distribution (T={T}K)')
plt.legend()
plt.grid(True, alpha=0.3)

# Calculate characteristic speeds
v_avg = np.mean(speeds)
v_rms = np.sqrt(np.mean(speeds**2))
v_mp = np.sqrt(2 * k_B * T / m) # Most probable speed

plt.axvline(v_avg, color='blue', linestyle='--', label=f'Average: {v_avg:.0f} m/s')
plt.axvline(v_rms, color='green', linestyle='--', label=f'RMS: {v_rms:.0f} m/s')
plt.axvline(v_mp, color='red', linestyle='--', label=f'Most probable: {v_mp:.0f} m/s')
plt.legend()

# Energy distribution
plt.subplot(2, 2, 2)
kinetic_energies = 0.5 * m * speeds**2
plt.hist(kinetic_energies / k_B, bins=50, density=True, alpha=0.7)
plt.xlabel('Kinetic Energy / k_B (K)')
plt.ylabel('Probability Density')
plt.title('Kinetic Energy Distribution')
plt.grid(True, alpha=0.3)

# Temperature dependence
plt.subplot(2, 2, 3)
temperatures = [200, 300, 400, 500]
v_range = np.linspace(0, 1500, 200)

for T_temp in temperatures:
mb_dist = maxwell_boltzmann(v_range, T_temp, m)
plt.plot(v_range, mb_dist, linewidth=2, label=f'T = {T_temp}K')

plt.xlabel('Speed (m/s)')
plt.ylabel('Probability Density')
plt.title('Temperature Dependence')
plt.legend()
plt.grid(True, alpha=0.3)

# Pressure calculation from kinetic theory
plt.subplot(2, 2, 4)
volumes = np.linspace(0.01, 0.1, 100) # m³
n_moles = 1.0 # mol
R = 8.314 # J/(mol⋅K)

# Ideal gas law: PV = nRT
pressures_ideal = n_moles * R * T / volumes

plt.plot(volumes*1000, pressures_ideal/1000, 'b-', linewidth=2, label='Ideal Gas Law')
plt.xlabel('Volume (L)')
plt.ylabel('Pressure (kPa)')
plt.title(f'P-V Diagram (T={T}K)')
plt.grid(True, alpha=0.3)
plt.legend()

plt.tight_layout()
plt.show()

print(f"Theoretical values at T = {T}K:")
print(f"Average speed: {np.sqrt(8*k_B*T/(np.pi*m)):.0f} m/s")
print(f"RMS speed: {np.sqrt(3*k_B*T/m):.0f} m/s")
print(f"Most probable speed: {np.sqrt(2*k_B*T/m):.0f} m/s")

return speeds, kinetic_energies

# Heat transfer
def heat_transfer_simulation():
"""Simulate heat conduction"""

# One-dimensional heat equation: ∂T/∂t = α ∂²T/∂x²
# where α is thermal diffusivity

# Parameters
L = 1.0 # Length (m)
alpha = 1e-4 # Thermal diffusivity (m²/s)
dx = 0.01 # Spatial step
dt = 0.1 # Time step
t_max = 3600 # Maximum time (s)

# Grid
x = np.arange(0, L + dx, dx)
nx = len(x)

# Check stability condition
stability = alpha * dt / dx**2
print(f"Stability parameter: {stability:.3f} (should be ≤ 0.5)")

# Initial condition: hot at one end, cold at the other
T = np.zeros(nx)
T[0] = 100 # Hot end (°C)
T[-1] = 0 # Cold end (°C)

# Time stepping
times = np.arange(0, t_max + dt, dt)
T_history = []

for t in times[::int(len(times)/10)]: # Store every 10th time step
T_history.append(T.copy())

# Update temperature using finite difference
T_new = T.copy()
for i in range(1, nx-1):
T_new[i] = T[i] + alpha * dt / dx**2 * (T[i+1] - 2*T[i] + T[i-1])

T = T_new

# Plot temperature evolution
plt.figure(figsize=(10, 6))

plt.subplot(1, 2, 1)
for i, T_snapshot in enumerate(T_history):
time_label = f't = {i * t_max / 10 / 60:.0f} min'
plt.plot(x, T_snapshot, linewidth=2, label=time_label)

plt.xlabel('Position (m)')
plt.ylabel('Temperature (°C)')
plt.title('Heat Conduction Over Time')
plt.legend()
plt.grid(True, alpha=0.3)

# Steady-state solution (linear)
plt.subplot(1, 2, 2)
T_steady = 100 * (1 - x / L) # Linear temperature profile
plt.plot(x, T_history[-1], 'b-', linewidth=2, label='Numerical solution')
plt.plot(x, T_steady, 'r--', linewidth=2, label='Analytical steady-state')
plt.xlabel('Position (m)')
plt.ylabel('Temperature (°C)')
plt.title('Final Temperature Profile')
plt.legend()
plt.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

return x, T_history

# Execute thermal physics simulations
speeds, energies = kinetic_theory_simulation()
x_heat, T_heat = heat_transfer_simulation()

return speeds, energies, x_heat, T_heat

# Execute thermal physics analysis
thermal_results = thermal_physics()

Building Physics Intuition

Dimensional Analysis and Scaling

Physics understanding begins with dimensional analysis:

# Dimensional analysis tools
def dimensional_analysis():
"""
Explore dimensional analysis and scaling laws
"""

# Fundamental dimensions
dimensions = {
'length': 'L',
'mass': 'M',
'time': 'T',
'electric_current': 'I',
'temperature': 'Θ',
'amount': 'N',
'luminous_intensity': 'J'
}

# Common physical quantities and their dimensions
quantities = {
'velocity': 'L T⁻¹',
'acceleration': 'L T⁻²',
'force': 'M L T⁻²',
'energy': 'M L² T⁻²',
'power': 'M L² T⁻³',
'pressure': 'M L⁻¹ T⁻²',
'electric_field': 'M L T⁻³ I⁻¹',
'magnetic_field': 'M T⁻² I⁻¹'
}

print("DIMENSIONAL ANALYSIS EXAMPLES")
print("="*40)

# Example 1: Pendulum period
print("Example 1: Simple Pendulum Period")
print("Variables: length L, gravity g, mass m")
print("Dimensions: [L] = L, [g] = L T⁻², [m] = M")
print("Proposed form: T = L^a × g^b × m^c")
print("Dimension equation: T = L^a × (L T⁻²)^b × M^c")
print(" T = L^(a+b) × T^(-2b) × M^c")
print("Matching dimensions:")
print(" Time: 1 = -2b → b = -1/2")
print(" Length: 0 = a + b → a = 1/2")
print(" Mass: 0 = c → c = 0")
print("Result: T ∝ √(L/g)")

# Numerical verification
L_values = np.array([0.5, 1.0, 1.5, 2.0]) # lengths in meters
g = 9.81
T_predicted = 2 * np.pi * np.sqrt(L_values / g)

plt.figure(figsize=(10, 6))
plt.subplot(1, 2, 1)
plt.plot(L_values, T_predicted, 'bo-', linewidth=2, markersize=8)
plt.xlabel('Length (m)')
plt.ylabel('Period (s)')
plt.title('Pendulum Period vs Length')
plt.grid(True, alpha=0.3)

# Example 2: Drag force scaling
print("\nExample 2: Drag Force at High Speed")
print("Variables: velocity v, density ρ, area A, coefficient C_d")
print("Expected form: F_drag = ½ × C_d × ρ × A × v²")

# Show velocity scaling
velocities = np.linspace(10, 100, 20)
rho = 1.225 # Air density (kg/m³)
A = 2.0 # Cross-sectional area (m²)
C_d = 0.3 # Drag coefficient

F_drag = 0.5 * C_d * rho * A * velocities**2

plt.subplot(1, 2, 2)
plt.plot(velocities, F_drag, 'r-', linewidth=2)
plt.xlabel('Velocity (m/s)')
plt.ylabel('Drag Force (N)')
plt.title('Drag Force vs Velocity (∝ v²)')
plt.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

return quantities, T_predicted, F_drag

# Execute dimensional analysis
dim_analysis = dimensional_analysis()

Error Analysis and Uncertainty

Understanding measurement uncertainty is crucial in physics:

# Error analysis and uncertainty propagation
def error_analysis():
"""
Demonstrate error analysis and uncertainty propagation
"""

print("ERROR ANALYSIS IN PHYSICS")
print("="*30)

# Example: Measuring gravitational acceleration
# Using pendulum: g = 4π²L/T²

# Measurements with uncertainties
L = 1.000 # Length (m)
dL = 0.001 # Uncertainty in length (m)
T = 2.006 # Period (s)
dT = 0.010 # Uncertainty in period (s)

# Calculate g
g = 4 * np.pi**2 * L / T**2

# Uncertainty propagation: g = 4π²L/T²
# Relative uncertainty: δg/g = δL/L + 2δT/T
relative_uncertainty = dL/L + 2*dT/T
dg = g * relative_uncertainty

print(f"Measurement results:")
print(f"Length: L = {L:.3f} ± {dL:.3f} m")
print(f"Period: T = {T:.3f} ± {dT:.3f} s")
print(f"Calculated g = {g:.2f} ± {dg:.2f} m/s²")
print(f"Accepted value: 9.81 m/s²")
print(f"Percent error: {abs(g - 9.81)/9.81 * 100:.1f}%")

# Monte Carlo error propagation
def monte_carlo_errors(n_samples=10000):
"""Use Monte Carlo method for error propagation"""

# Generate random samples
L_samples = np.random.normal(L, dL, n_samples)
T_samples = np.random.normal(T, dT, n_samples)

# Calculate g for each sample
g_samples = 4 * np.pi**2 * L_samples / T_samples**2

# Statistics
g_mean = np.mean(g_samples)
g_std = np.std(g_samples)

print(f"\nMonte Carlo results ({n_samples} samples):")
print(f"g = {g_mean:.2f} ± {g_std:.2f} m/s²")

# Plot distribution
plt.figure(figsize=(10, 6))

plt.subplot(1, 2, 1)
plt.hist(g_samples, bins=50, density=True, alpha=0.7, edgecolor='black')
plt.axvline(g_mean, color='red', linestyle='--', linewidth=2, label=f'Mean: {g_mean:.2f}')
plt.axvline(9.81, color='green', linestyle='--', linewidth=2, label='True value: 9.81')
plt.xlabel('g (m/s²)')
plt.ylabel('Probability Density')
plt.title('Distribution of g Values')
plt.legend()
plt.grid(True, alpha=0.3)

# Error vs number of measurements
plt.subplot(1, 2, 2)
n_measurements = np.logspace(1, 4, 20).astype(int)
standard_errors = []

for n in n_measurements:
# Standard error decreases as 1/√n
standard_error = g_std / np.sqrt(n)
standard_errors.append(standard_error)

plt.loglog(n_measurements, standard_errors, 'bo-', linewidth=2)
plt.loglog(n_measurements, 1/np.sqrt(n_measurements), 'r--',
linewidth=2, label='∝ 1/√n')
plt.xlabel('Number of Measurements')
plt.ylabel('Standard Error (m/s²)')
plt.title('Error Reduction with More Measurements')
plt.legend()
plt.grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

return g_samples, g_mean, g_std

# Execute Monte Carlo analysis
g_samples, g_mean, g_std = monte_carlo_errors()

return g, dg, g_samples

# Execute error analysis
error_results = error_analysis()

Next Steps in Your Physics Journey

Immediate Explorations

Planetary Motion: Simulate orbital mechanics and Kepler's laws

Wave Interference: Explore wave superposition and interference patterns

Quantum Basics: Introduction to wave-particle duality

Thermal Equilibrium: Statistical mechanics fundamentals

Building Toward Advanced Physics

Classical Field Theory: Maxwell's equations and electromagnetic waves

Quantum Mechanics: Schrödinger equation and quantum systems

Statistical Physics: Phase transitions and critical phenomena

Relativity: Special and general relativistic effects

Research Skills Development

Experimental Design: Planning and analyzing physics experiments

Data Analysis: Statistical methods for physics data

Modeling: Creating and validating physical models

Communication: Presenting physics results effectively

Computational physics in CoCalc opens new ways to understand the physical world. Start with fundamental concepts, build your computational skills, and gradually explore more sophisticated phenomena. The combination of theory, computation, and visualization makes complex physics accessible and engaging.

Begin your computational physics journey. Access physics simulations and start exploring at cocalc.com