Radical! Nearest neighbour tracing in nuclear automata. Below the program is one that incorporates general relativity and now gravity time-dialation ratio/scalar curvature as a k-body factor due to certain chemical stoichastic reporting systems you may choose to use either, this new update also explores a hypothetical relationship between the gravitational curvature parameter and a quantum mechanical property, this might elucidate how such a "bridge" might be conceptually built, even if not physically validated :-)
For the second program please search, I postulate
import numpy as np
import matplotlib.pyplot as plt
import pubchempy as pcp
import requests
import json
import time
import re
# --- Scientific Constants ---
# Nuclear Physics Constants
CURIE_TO_BQ = 3.7e10 # Conversion factor: 1 Curie = 3.7 x 10^10 Becquerels (disintegrations per second)
LN2 = np.log(2) # Natural logarithm of 2, used in half-life calculations
# Blackbody Radiation Constants
STEFAN_BOLTZMANN_CONSTANT = 5.670374419e-8 # Stefan-Boltzmann constant in W/(m^2·K^4)
PLANCK_CONSTANT = 6.62607015e-34 # Planck's constant in J·s
BOLTZMANN_CONSTANT = 1.380649e-23 # Boltzmann constant in J/K
SPEED_OF_LIGHT = 2.99792458e8 # Speed of light in a vacuum in m/s
WIEN_DISPLACEMENT_CONSTANT = 2.897771955e-3 # Wien's displacement law constant in m·K
# --- Manual Nuclear Data for Gamma Emission ---
# This dictionary stores known half-lives and gamma emission properties for isotopes.
# IMPORTANT: These values should be verified from reliable nuclear data sources (e.g., NNDC).
# PubChem provides general information, but specific gamma yields are best from dedicated databases.
MANUAL_GAMMA_DATA = {
"Cesium-137": {
"half_life_days": 30.08 * 365.25, # Half-life in days (approx 30.08 years)
"gamma_energies_MeV": [0.662],
"gamma_yields_per_decay": [0.851] # Gamma yield per decay (from Ba-137m isomer)
},
"Cobalt-60": {
"half_life_days": 1925.21, # Half-life in days
"gamma_energies_MeV": [1.173, 1.332],
"gamma_yields_per_decay": [0.999, 0.999] # Approximately 1 gamma per decay at each energy
},
"Iodine-131": {
"half_life_days": 8.02, # Half-life in days
"gamma_energies_MeV": [0.364],
"gamma_yields_per_decay": [0.817] # Gamma yield per decay
},
"Potassium-40": {
"half_life_days": 1.251e9 * 365.25, # Half-life in days (~1.251 billion years)
"gamma_energies_MeV": [1.4608],
"gamma_yields_per_decay": [0.1067] # ~10.67% of K-40 decays via electron capture
},
# Feel free to add more isotopes here with their respective data
}
# --- PubChem Helper Functions (for half-life lookup if not in manual data) ---
# These functions attempt to fetch half-life data from PubChem,
# but their parsing can be fragile due to variations in PubChem's data structure.
def _get_compound_cid_pubchem(compound_name):
"""
Searches PubChem for a compound's CID (Compound ID) by name.
Returns the CID if found, otherwise None.
"""
try:
compounds = pcp.get_compounds(compound_name, 'name')
return compounds[0].cid if compounds else None
except Exception:
return None
def _fetch_half_life_from_pubchem(cid):
"""
Fetches half-life data from PubChem using a compound's CID.
Attempts to parse numerical value and unit from the returned text.
Returns (value, unit) tuple if successful, otherwise (None, None).
"""
base_url = "https://pubchem.ncbi.nlm.nih.gov/rest/pug_view/data/compound/"
url = f"{base_url}cid/{cid}/JSON/"
try:
response = requests.get(url, timeout=5)
response.raise_for_status() # Raise an exception for HTTP errors
data = response.json()
search_sections = ['Atomic Mass, Half Life, and Decay', 'Chemical and Physical Properties', 'Depositor-Supplied']
for record in data.get('Record', {}).get('Section', []):
if record.get('TOCHeading') in search_sections:
for prop_section in record.get('Section', []):
if 'half-life' in prop_section.get('TOCHeading', '').lower():
for item in prop_section.get('Information', []):
for value_dict in item.get('Value', {}).get('StringWithMarkup', []):
text = value_dict.get('String')
if text:
# Use regex to find numerical value and unit
match = re.match(r"(\d+\.?\d*(?:e[+-]?\d+)?)\s*(seconds|minutes|hours|days|years|s|min|hr|d|y)", text.lower())
if match:
return float(match.group(1)), match.group(2)
return None, None
except Exception:
return None, None
def _convert_half_life_to_days(value, unit):
"""
Converts a half-life value from various units to days.
Returns the half-life in days, or None if the unit is unrecognized.
"""
if unit in ['s', 'seconds']: return value / (24 * 3600)
if unit in ['min', 'minutes']: return value / (24 * 60)
if unit in ['hr', 'hours']: return value / 24
if unit in ['d', 'days']: return value
if unit in ['y', 'years']: return value * 365.25
return None
# --- Core Radiation Calculation Functions ---
def calculate_decay_constant(half_life_days):
"""
Calculates the radioactive decay constant (lambda) from the half-life.
Returns the decay constant in inverse seconds (s^-1).
If half_life_days is non-positive, returns 0 (implying a stable isotope).
"""
if half_life_days is None or half_life_days <= 0:
return 0
half_life_seconds = half_life_days * 24 * 3600 # Convert days to seconds
return LN2 / half_life_seconds
def calculate_activity_at_time(initial_activity_Bq, decay_constant, time_elapsed_seconds):
"""
Calculates the activity (number of decays per second) of a radioactive source
at a given time, using the decay law (N = N0 * e^(-lambda*t)).
Returns activity in Becquerels (Bq).
"""
if decay_constant == 0:
return initial_activity_Bq # For stable isotopes or infinite half-life
return initial_activity_Bq * np.exp(-decay_constant * time_elapsed_seconds)
def get_gamma_curies_single_point(isotope_name, initial_activity_Bq, time_elapsed_days):
"""
Calculates the total gamma curies emitted by an isotope at a specific time.
It first tries to get half-life from the `MANUAL_GAMMA_DATA`, then falls back to PubChem.
Returns the gamma activity in Curies (Ci), or None if half-life data isn't found.
"""
half_life_days = None
if isotope_name in MANUAL_GAMMA_DATA:
half_life_days = MANUAL_GAMMA_DATA[isotope_name]["half_life_days"]
else:
# If not in manual data, try PubChem
cid = _get_compound_cid_pubchem(isotope_name)
if cid:
hl_val, hl_unit = _fetch_half_life_from_pubchem(cid)
if hl_val:
half_life_days = _convert_half_life_to_days(hl_val, hl_unit)
if half_life_days is None:
return None # Could not determine half-life
gamma_yields = MANUAL_GAMMA_DATA.get(isotope_name, {}).get("gamma_yields_per_decay", [0])
total_gamma_yield_per_decay = sum(gamma_yields)
decay_constant = calculate_decay_constant(half_life_days)
time_elapsed_seconds = time_elapsed_days * 24 * 3600
current_activity_Bq = calculate_activity_at_time(initial_activity_Bq, decay_constant, time_elapsed_seconds)
total_gamma_disintegrations_per_second = current_activity_Bq * total_gamma_yield_per_decay
gamma_curies = total_gamma_disintegrations_per_second / CURIE_TO_BQ
return gamma_curies
def calculate_blackbody_spectral_radiance(wavelength_m, temperature_K):
"""
Calculates the spectral radiance of a blackbody at a given wavelength and temperature,
using Planck's Law.
Returns spectral radiance in Watts per steradian per square meter per meter of wavelength (W/(m³·sr)).
Returns 0.0 if temperature or wavelength are non-positive or if the exponent is too large.
"""
if temperature_K <= 0 or wavelength_m <= 0:
return 0.0
term1 = (2 * PLANCK_CONSTANT * SPEED_OF_LIGHT**2) / (wavelength_m**5)
term2_exp_arg = (PLANCK_CONSTANT * SPEED_OF_LIGHT) / (wavelength_m * BOLTZMANN_CONSTANT * temperature_K)
# Check to prevent numerical overflow for very large exponent arguments
if term2_exp_arg > 700:
return 0.0 # Effectively zero emission at this wavelength/temperature
spectral_radiance = term1 / (np.exp(term2_exp_arg) - 1)
return spectral_radiance
# --- Utility Function for Generating Bell Curve Data ---
def generate_normal_distribution_series(mean, std_dev, num_points=200, spread_factor=3):
"""
Generates a series of values that follow a normal (bell curve) distribution.
This is used to create a distribution for a parameter, like measurement time.
Args:
mean (float): The mean (average) of the distribution.
std_dev (float): The standard deviation of the distribution.
num_points (int): The number of data points to generate.
spread_factor (int): How many standard deviations away from the mean
the series should span (e.g., 3 for 3-sigma range).
Returns:
tuple: A tuple containing:
- series_values (numpy.array): The generated data points.
- pdf_values (numpy.array): The probability density function (PDF) values
corresponding to `series_values`.
"""
start = mean - spread_factor * std_dev
end = mean + spread_factor * std_dev
if start < 0:
start = 0 # Ensure time or other physical quantities are non-negative
series_values = np.linspace(start, end, num_points)
pdf_values = (1 / (std_dev * np.sqrt(2 * np.pi))) * \
np.exp(-0.5 * ((series_values - mean) / std_dev)**2)
return series_values, pdf_values
# --- Analysis Functions for "Nearest Neighbor Boundaries" and "Displacement" ---
def analyze_gamma_curies_nearest_neighbor_boundaries(isotope_name, initial_activity_Bq,
mean_time_days, std_dev_time_days,
num_points=200, derivative_threshold_abs=None):
"""
Calculates gamma curies over a series of measurement times sampled from a bell curve.
It identifies "nearest neighbor boundaries" by highlighting regions where the
absolute rate of change (derivative) of gamma curies is above a given threshold.
This represents where the decay process is changing most rapidly.
Args:
isotope_name (str): The name of the radioactive isotope (e.g., "Cesium-137").
initial_activity_Bq (float): The initial activity of the isotope in Becquerels (Bq).
mean_time_days (float): The mean (average) of the bell curve distribution for measurement times in days.
std_dev_time_days (float): The standard deviation of the bell curve for measurement times in days.
num_points (int): The number of time points to simulate.
derivative_threshold_abs (float, optional): The absolute threshold for the derivative (Ci/day).
If None, the program will auto-calculate a threshold
(e.g., 10% of the maximum absolute derivative).
"""
print(f"\n--- Identifying 'Nearest Neighbor Boundaries' in Gamma Curies ---")
print(f" Isotope: {isotope_name}, Initial Activity: {initial_activity_Bq:.2e} Bq")
print(f" Time Distribution: Mean={mean_time_days:.2f} days, Std Dev={std_dev_time_days:.2f} days")
print(f" 'Nearest Neighbor Boundaries' are regions where the decay rate changes significantly.")
# 1. Generate time series based on a normal distribution
time_series_days, pdf_values = generate_normal_distribution_series(
mean_time_days, std_dev_time_days, num_points=num_points)
# 2. Calculate gamma curies for each time point
curies_series = []
for t_days in time_series_days:
curies = get_gamma_curies_single_point(isotope_name, initial_activity_Bq, t_days)
curies_series.append(curies if curies is not None else np.nan) # Store None as NaN
curies_series = np.array(curies_series) # Convert to NumPy array for numerical operations
# Handle cases where half-life was not found (curies_series will be all NaN)
if np.all(np.isnan(curies_series)):
print(f" Calculation failed. Please ensure '{isotope_name}' has half-life data.")
return
# 3. Calculate the numerical derivative (rate of change) of gamma curies with respect to time
# np.gradient provides a robust way to estimate the derivative for discrete data.
derivative_curies = np.gradient(curies_series, time_series_days)
# 4. Determine or use a threshold for identifying "boundaries"
if derivative_threshold_abs is None:
# If no threshold is provided, auto-calculate it as 10% of the maximum absolute derivative
# This highlights the top 10% steepest changes.
max_abs_deriv = np.max(np.abs(derivative_curies[np.isfinite(derivative_curies)]))
derivative_threshold_abs = max_abs_deriv * 0.1
print(f" Auto-calculated derivative threshold: {derivative_threshold_abs:.2e} Ci/day")
else:
print(f" Using provided derivative threshold: {derivative_threshold_abs:.2e} Ci/day")
# 5. Identify points where the absolute derivative exceeds the threshold
boundary_indices = np.where(np.abs(derivative_curies) > derivative_threshold_abs)
boundary_times = time_series_days[boundary_indices]
boundary_curies = curies_series[boundary_indices]
# --- Plotting Results ---
plt.figure(figsize=(12, 7))
plt.plot(time_series_days, curies_series, label='Gamma Curies Decay')
plt.scatter(boundary_times, boundary_curies, color='red', s=50,
label='Nearest Neighbor Boundaries (High Rate of Change)', zorder=5)
plt.title(f'Gamma Curies with "Nearest Neighbor Boundaries" for {isotope_name}')
plt.xlabel('Time (days)')
plt.ylabel('Gamma Curies (Ci)')
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()
# Optional: Plot the derivative itself for more insight
plt.figure(figsize=(10, 5))
plt.plot(time_series_days, derivative_curies, label='Rate of Change (dCuries/dt)')
plt.axhline(y=derivative_threshold_abs, color='grey', linestyle='--', label='Positive Threshold')
plt.axhline(y=-derivative_threshold_abs, color='grey', linestyle='--', label='Negative Threshold')
plt.title('Rate of Change of Gamma Curies over Time')
plt.xlabel('Time (days)')
plt.ylabel('dCuries/dt (Ci/day)')
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()
print(f"\n Identified {len(boundary_times)} data points in 'Nearest Neighbor Boundary' regions.")
if len(boundary_times) > 0:
print(" Example Boundary Points (Time, Curies, Derivative):")
# Print a few examples from the identified boundary points
for i in range(min(5, len(boundary_times))):
idx = boundary_indices[0][i] # Get the original index from the series
print(f" Time: {time_series_days[idx]:.2f} days, Curies: {curies_series[idx]:.2e} Ci, Derivative: {derivative_curies[idx]:.2e} Ci/day")
def simulate_blackbody_wien_displacement(temperatures_K, wavelength_range_nm=(1, 10000)):
"""
Simulates blackbody spectral radiance across a range of wavelengths for given temperatures,
and demonstrates Wien's Displacement Law by showing how the peak emission wavelength
("displaces" or shifts) with temperature.
Args:
temperatures_K (list): A list of temperatures in Kelvin (K) to simulate.
wavelength_range_nm (tuple): A tuple (min_wavelength, max_wavelength) in nanometers.
"""
print(f"\n--- Simulating Blackbody Spectral Radiance & Wien's Displacement Law ---")
print(f" Analyzing temperatures: {temperatures_K} K")
wavelengths_nm = np.linspace(wavelength_range_nm[0], wavelength_range_nm[1], 500)
wavelengths_m = wavelengths_nm * 1e-9 # Convert nanometers to meters
peak_wavelengths_data = [] # To store (Temperature, Peak Wavelength) for Wien's Law plot
plt.figure(figsize=(12, 7))
for T in temperatures_K:
if T <= 0:
print(f" Skipping temperature {T}K as it is non-positive.")
continue
radiance_series = [calculate_blackbody_spectral_radiance(wl_m, T) for wl_m in wavelengths_m]
# Find the peak wavelength for the current temperature
max_radiance_idx = np.argmax(radiance_series)
peak_wl_m = wavelengths_m[max_radiance_idx]
peak_wl_nm = wavelengths_nm[max_radiance_idx]
peak_wavelengths_data.append((T, peak_wl_nm))
plt.plot(wavelengths_nm, radiance_series, label=f'T = {T} K (Peak: {peak_wl_nm:.2f} nm)')
plt.title('Blackbody Spectral Radiance (Planck\'s Law)')
plt.xlabel('Wavelength (nm)')
plt.ylabel('Spectral Radiance (W/(m³·sr))')
plt.yscale('log') # Log scale helps visualize broad range of radiance values
plt.xscale('log') # Log scale for wavelength can also be useful for broad ranges
plt.grid(True, which="both", ls="--", alpha=0.7)
plt.legend()
plt.tight_layout()
plt.show()
print("\n--- Wien's Displacement Law Data Series ---")
print(" This shows how the peak emission wavelength 'displaces' (shifts) to shorter wavelengths as temperature increases.")
temp_list_sorted = [data[0] for data in sorted(peak_wavelengths_data)]
peak_wl_list = [data[1] for data in sorted(peak_wavelengths_data)]
for T, peak_wl_nm in sorted(peak_wavelengths_data):
print(f" Temperature: {T} K, Calculated Peak Wavelength: {peak_wl_nm:.2f} nm")
plt.figure(figsize=(8, 6))
plt.plot(temp_list_sorted, peak_wl_list, 'o-', label='Calculated Peak Wavelength')
# Calculate and plot the theoretical curve from Wien's Displacement Law
theoretical_peak_nm = (WIEN_DISPLACEMENT_CONSTANT / np.array(temp_list_sorted)) * 1e9
plt.plot(temp_list_sorted, theoretical_peak_nm, 'r--', label='Theoretical (Wien\'s Law)')
plt.title('Wien\'s Displacement Law: Peak Wavelength vs. Temperature')
plt.xlabel('Temperature (K)')
plt.ylabel('Peak Wavelength (nm)')
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()
# --- Main Program Execution ---
if __name__ == "__main__":
print("--- Radiation Analysis Program: Interpreting Data Behavior ---")
print("\nThis program helps you understand how different scientific concepts can be analyzed using statistical ideas.")
print("We'll explore two scenarios:")
print(" 1. **'Nearest Neighbor Boundaries' in Gamma Curies:** This identifies points where the rate of radioactive decay is changing most rapidly, relative to measurement times that follow a bell curve.")
print(" 2. **Blackbody Wavelength 'Displacement':** This demonstrates how the peak wavelength of light emitted by a hot object shifts as its temperature changes (Wien's Law).\n")
print("Please note: The terms 'Nearest Neighbor Boundaries' and 'Edge States' are statistical interpretations here, not direct physical phenomena like in quantum mechanics.\n")
while True:
main_choice = input("Select: (N)earest Neighbor Boundaries (Gamma Curies), (B)lackbody Wavelength Displacement, or (Q)uit? ").strip().lower()
if main_choice == 'q':
break
elif main_choice == 'n':
print("\n--- 'Nearest Neighbor Boundaries' for Gamma Curies ---")
isotope_name = input(" Enter radioactive isotope name (e.g., Cesium-137, Cobalt-60): ").strip()
# Check if isotope data is available either manually or via PubChem
if isotope_name not in MANUAL_GAMMA_DATA and _get_compound_cid_pubchem(isotope_name) is None:
print(f" '{isotope_name}' data not found. Please check spelling or add it to MANUAL_GAMMA_DATA.")
continue
try:
initial_activity = float(input(f" Enter initial activity of {isotope_name} in Becquerels (Bq): "))
mean_time = float(input(f" Enter mean measurement time for the time distribution (days): "))
std_dev_time = float(input(f" Enter standard deviation for the measurement time (days): "))
num_points = int(input(f" Enter number of data points for the simulation (e.g., 200): "))
deriv_thresh_input = input(f" Enter absolute derivative threshold for 'boundaries' (Ci/day, leave blank for auto-calc): ")
derivative_threshold_abs = float(deriv_thresh_input) if deriv_thresh_input else None
except ValueError:
print(" Invalid input. Please enter numerical values.")
continue
analyze_gamma_curies_nearest_neighbor_boundaries(
isotope_name, initial_activity, mean_time, std_dev_time, num_points, derivative_threshold_abs)
time.sleep(1) # Small pause for readability
elif main_choice == 'b':
print("\n--- Blackbody Wavelength 'Displacement' (Wien's Law) ---")
temp_input_str = input(" Enter temperatures (K) separated by commas (e.g., 300, 1000, 5000): ")
try:
temperatures_list = [float(t.strip()) for t in temp_input_str.split(',')]
except ValueError:
print(" Invalid temperature input. Please enter numbers separated by commas.")
continue
simulate_blackbody_wien_displacement(temperatures_list)
time.sleep(1) # Small pause for readability
else:
print("Invalid choice. Please enter 'N', 'B', or 'Q'.")
print("\nExiting program. Thank you for exploring these concepts!")
General Relativism as a K-Body factor for chemical study.
I postulate there is no Universe in a black hole :-o
import numpy as np
import matplotlib.pyplot as plt
import pubchempy as pcp
import requests
import json
import time
import re
# --- Scientific Constants ---
# Nuclear Physics Constants
CURIE_TO_BQ = 3.7e10 # Conversion factor: 1 Curie = 3.7 x 10^10 Becquerels (disintegrations per second)
LN2 = np.log(2) # Natural logarithm of 2, used in half-life calculations
# Blackbody Radiation Constants
STEFAN_BOLTZMANN_CONSTANT = 5.670374419e-8 # Stefan-Boltzmann constant in W/(m^2·K^4)
PLANCK_CONSTANT = 6.62607015e-34 # Planck's constant in J·s
BOLTZMANN_CONSTANT = 1.380649e-23 # Boltzmann constant in J/K
SPEED_OF_LIGHT = 2.99792458e8 # Speed of light in a vacuum in m/s
WIEN_DISPLACEMENT_CONSTANT = 2.897771955e-3 # Wien's displacement law constant in m·K
# General Relativity Constants (for the new functions)
GRAVITATIONAL_CONSTANT = 6.67430e-11 # G in m³·kg⁻¹·s⁻²
SPEED_OF_LIGHT_SQ = SPEED_OF_LIGHT**2
# --- Manual Nuclear Data for Gamma Emission ---
MANUAL_GAMMA_DATA = {
"Cesium-137": {
"half_life_days": 30.08 * 365.25,
"gamma_energies_MeV": [0.662],
"gamma_yields_per_decay": [0.851]
},
"Cobalt-60": {
"half_life_days": 1925.21,
"gamma_energies_MeV": [1.173, 1.332],
"gamma_yields_per_decay": [0.999, 0.999]
},
"Iodine-131": {
"half_life_days": 8.02,
"gamma_energies_MeV": [0.364],
"gamma_yields_per_decay": [0.817]
},
"Potassium-40": {
"half_life_days": 1.251e9 * 365.25,
"gamma_energies_MeV": [1.4608],
"gamma_yields_per_decay": [0.1067]
},
}
# --- PubChem Helper Functions (Unchanged) ---
def _get_compound_cid_pubchem(compound_name):
try:
compounds = pcp.get_compounds(compound_name, 'name')
return compounds[0].cid if compounds else None
except Exception:
return None
def _fetch_half_life_from_pubchem(cid):
base_url = "https://pubchem.ncbi.nlm.nih.gov/rest/pug_view/data/compound/"
url = f"{base_url}cid/{cid}/JSON/"
try:
response = requests.get(url, timeout=5)
response.raise_for_status()
data = response.json()
search_sections = ['Atomic Mass, Half Life, and Decay', 'Chemical and Physical Properties', 'Depositor-Supplied']
for record in data.get('Record', {}).get('Section', []):
if record.get('TOCHeading') in search_sections:
for prop_section in record.get('Section', []):
if 'half-life' in prop_section.get('TOCHeading', '').lower():
for item in prop_section.get('Information', []):
for value_dict in item.get('Value', {}).get('StringWithMarkup', []):
text = value_dict.get('String')
if text:
match = re.match(r"(\d+\.?\d*(?:e[+-]?\d+)?)\s*(seconds|minutes|hours|days|years|s|min|hr|d|y)", text.lower())
if match:
return float(match.group(1)), match.group(2)
return None, None
except Exception:
return None, None
def _convert_half_life_to_days(value, unit):
if unit in ['s', 'seconds']: return value / (24 * 3600)
if unit in ['min', 'minutes']: return value / (24 * 60)
if unit in ['hr', 'hours']: return value / 24
if unit in ['d', 'days']: return value
if unit in ['y', 'years']: return value * 365.25
return None
# --- Core Radiation Calculation Functions (Modified to accept dilation factor) ---
def calculate_decay_constant(half_life_days):
if half_life_days is None or half_life_days <= 0:
return 0
half_life_seconds = half_life_days * 24 * 3600
return LN2 / half_life_seconds
def calculate_activity_at_time(initial_activity_Bq, decay_constant, time_elapsed_seconds):
if decay_constant == 0:
return initial_activity_Bq
return initial_activity_Bq * np.exp(-decay_constant * time_elapsed_seconds)
def calculate_blackbody_spectral_radiance(wavelength_m, temperature_K):
if temperature_K <= 0 or wavelength_m <= 0:
return 0.0
term1 = (2 * PLANCK_CONSTANT * SPEED_OF_LIGHT_SQ) / (wavelength_m**5)
term2_exp_arg = (PLANCK_CONSTANT * SPEED_OF_LIGHT) / (wavelength_m * BOLTZMANN_CONSTANT * temperature_K)
if term2_exp_arg > 700:
return 0.0
spectral_radiance = term1 / (np.exp(term2_exp_arg) - 1)
return spectral_radiance
# --- UPDATED: General Relativity Functions ---
def calculate_schwarzschild_radius(mass_kg):
"""Calculates the Schwarzschild radius for a given mass."""
if mass_kg <= 0:
return 0.0
return (2 * GRAVITATIONAL_CONSTANT * mass_kg) / SPEED_OF_LIGHT_SQ
def calculate_dilation_factor(mass_kg, distance_m):
"""
Calculates the time dilation factor for an observer at infinity relative to an object in a gravitational field.
A factor of 1 means no dilation. Values < 1 mean time is dilated (slowed down).
"""
if mass_kg <= 0 or distance_m <= 0:
return 1.0 # No dilation for zero mass or distance
schwarzschild_radius = calculate_schwarzschild_radius(mass_kg)
if distance_m <= schwarzschild_radius:
return 0.0 # At or inside event horizon, time stops for an outside observer
return np.sqrt(1 - (schwarzschild_radius / distance_m))
def calculate_relativistic_activity(isotope_name, initial_activity_Bq, time_elapsed_days, mass_kg, distance_m):
"""
Calculates activity as observed by a distant observer, accounting for gravitational time dilation.
"""
half_life_days = MANUAL_GAMMA_DATA.get(isotope_name, {}).get("half_life_days")
if half_life_days is None:
cid = _get_compound_cid_pubchem(isotope_name)
if cid:
hl_val, hl_unit = _fetch_half_life_from_pubchem(cid)
if hl_val:
half_life_days = _convert_half_life_to_days(hl_val, hl_unit)
if half_life_days is None:
return None
decay_constant = calculate_decay_constant(half_life_days)
time_elapsed_seconds = time_elapsed_days * 24 * 3600
# Apply time dilation
dilation_factor = calculate_dilation_factor(mass_kg, distance_m)
# An observer far away sees time running slower for the object, so the decay process appears slowed.
effective_time_local_seconds = time_elapsed_seconds * dilation_factor
current_activity_Bq = calculate_activity_at_time(initial_activity_Bq, decay_constant, effective_time_local_seconds)
gamma_yields = MANUAL_GAMMA_DATA.get(isotope_name, {}).get("gamma_yields_per_decay", [0])
total_gamma_yield_per_decay = sum(gamma_yields)
total_gamma_disintegrations_per_second = current_activity_Bq * total_gamma_yield_per_decay
gamma_curies = total_gamma_disintegrations_per_second / CURIE_TO_BQ
return gamma_curies
def calculate_relativistic_peak_wavelength(local_temperature_K, mass_kg, distance_m):
"""
Calculates the observed peak wavelength of a blackbody, accounting for gravitational redshift.
"""
if local_temperature_K <= 0:
return None
# Gravitational redshift means a distant observer sees the blackbody as cooler.
dilation_factor = calculate_dilation_factor(mass_kg, distance_m)
observed_temperature_K = local_temperature_K * dilation_factor
# Use Wien's law with the observed temperature
if observed_temperature_K <= 0:
return np.inf # Effectively no peak wavelength observed
peak_wavelength_m = WIEN_DISPLACEMENT_CONSTANT / observed_temperature_K
return peak_wavelength_m * 1e9 # Convert to nm
# --- NEW: Quantum and Gravitational Bridge Functions ---
def calculate_gravitational_curvature_parameter(mass_kg, distance_m):
"""
A simplified, unitless parameter representing local gravitational curvature.
Higher values imply greater curvature.
"""
schwarzschild_radius = calculate_schwarzschild_radius(mass_kg)
if distance_m <= 0 or distance_m <= schwarzschild_radius:
return np.inf
return schwarzschild_radius / distance_m
def quantum_differentiation_hypothetical(curvature_parameter, base_quantum_value):
"""
Hypothetical function to illustrate 'quantum mechanical differentiation'.
This is NOT based on established physics but serves as an illustrative model.
It proposes that a quantum value (e.g., a spin state) changes with gravitational curvature.
For this example, we'll model it as a simple linear relationship for visualization.
"""
if curvature_parameter >= 1:
return 0 # At or inside the event horizon, the effect is undefined or at its maximum
# A hypothetical effect: the base value is 'differentiated' by the curvature.
# We use (1 - curvature_parameter) to get a value that goes from 1 to 0 as curvature increases.
differentiated_value = base_quantum_value * (1 - curvature_parameter)
return differentiated_value
# --- Utility Functions (Unchanged) ---
def generate_normal_distribution_series(mean, std_dev, num_points=200, spread_factor=3):
start = mean - spread_factor * std_dev
end = mean + spread_factor * std_dev
if start < 0:
start = 0
series_values = np.linspace(start, end, num_points)
pdf_values = (1 / (std_dev * np.sqrt(2 * np.pi))) * np.exp(-0.5 * ((series_values - mean) / std_dev)**2)
return series_values, pdf_values
def get_gamma_curies_single_point(isotope_name, initial_activity_Bq, time_elapsed_days):
half_life_days = None
if isotope_name in MANUAL_GAMMA_DATA:
half_life_days = MANUAL_GAMMA_DATA[isotope_name]["half_life_days"]
else:
cid = _get_compound_cid_pubchem(isotope_name)
if cid:
hl_val, hl_unit = _fetch_half_life_from_pubchem(cid)
if hl_val:
half_life_days = _convert_half_life_to_days(hl_val, hl_unit)
if half_life_days is None:
return None
gamma_yields = MANUAL_GAMMA_DATA.get(isotope_name, {}).get("gamma_yields_per_decay", [0])
total_gamma_yield_per_decay = sum(gamma_yields)
decay_constant = calculate_decay_constant(half_life_days)
time_elapsed_seconds = time_elapsed_days * 24 * 3600
current_activity_Bq = calculate_activity_at_time(initial_activity_Bq, decay_constant, time_elapsed_seconds)
total_gamma_disintegrations_per_second = current_activity_Bq * total_gamma_yield_per_decay
gamma_curies = total_gamma_disintegrations_per_second / CURIE_TO_BQ
return gamma_curies
# --- Analysis Functions (Modified to include new relativistic options) ---
def analyze_gamma_curies_comparison(isotope_name, initial_activity_Bq, mean_time_days, std_dev_time_days, mass_kg, distance_m, num_points=200):
"""
Compares gamma curies decay in flat spacetime vs. near a massive k-body.
"""
print(f"\n--- Comparing Gamma Curies: Flat Spacetime vs. Relativistic Environment ---")
print(f" Isotope: {isotope_name}, Initial Activity: {initial_activity_Bq:.2e} Bq")
print(f" Relativistic Environment: k-body Mass={mass_kg:.2e} kg, Distance={distance_m:.2e} m")
schwarzschild_radius = calculate_schwarzschild_radius(mass_kg)
dilation_factor = calculate_dilation_factor(mass_kg, distance_m)
print(f" Schwarzschild Radius: {schwarzschild_radius:.2e} m")
print(f" Gravitational Time Dilation Factor: {dilation_factor:.4f}")
if dilation_factor < 0.99:
print(f" (A dilation factor less than 1 means time runs slower in the gravitational field.)")
else:
print(f" (Dilation factor is near 1, indicating a weak gravitational effect.)")
time_series_days, _ = generate_normal_distribution_series(mean_time_days, std_dev_time_days, num_points=num_points)
# Flat Spacetime Calculation
flat_curies_series = [get_gamma_curies_single_point(isotope_name, initial_activity_Bq, t_days) for t_days in time_series_days]
# Relativistic Calculation
relativistic_curies_series = [calculate_relativistic_activity(isotope_name, initial_activity_Bq, t_days, mass_kg, distance_m) for t_days in time_series_days]
if all(c is None for c in flat_curies_series):
print(f" Calculation failed. Please ensure '{isotope_name}' has half-life data.")
return
# Plotting
plt.figure(figsize=(12, 7))
plt.plot(time_series_days, flat_curies_series, label='Flat Spacetime (Local Observer)')
plt.plot(time_series_days, relativistic_curies_series, label=f'Relativistic (Distant Observer)', linestyle='--')
plt.title(f'Gamma Curies Decay: {isotope_name} - Flat vs. Relativistic')
plt.xlabel('Time (days)')
plt.ylabel('Gamma Curies (Ci)')
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()
def simulate_blackbody_wien_comparison(temperatures_K, mass_kg, distance_m, wavelength_range_nm=(1, 10000)):
"""
Compares blackbody peak wavelength in flat spacetime vs. near a massive k-body.
"""
print(f"\n--- Comparing Blackbody Peak Wavelength: Flat vs. Relativistic ---")
print(f" Analyzing temperatures: {temperatures_K} K")
print(f" Relativistic Environment: k-body Mass={mass_kg:.2e} kg, Distance={distance_m:.2e} m")
schwarzschild_radius = calculate_schwarzschild_radius(mass_kg)
dilation_factor = calculate_dilation_factor(mass_kg, distance_m)
print(f" Schwarzschild Radius: {schwarzschild_radius:.2e} m")
print(f" Gravitational Time Dilation Factor: {dilation_factor:.4f}")
temp_list_sorted = sorted(temperatures_K)
flat_peak_wl_list_nm = [(WIEN_DISPLACEMENT_CONSTANT / T) * 1e9 for T in temp_list_sorted if T > 0]
relativistic_peak_wl_list_nm = [calculate_relativistic_peak_wavelength(T, mass_kg, distance_m) for T in temp_list_sorted if T > 0]
# Plotting
plt.figure(figsize=(8, 6))
plt.plot(temp_list_sorted, flat_peak_wl_list_nm, 'o-', label='Flat Spacetime (Local Observer)')
plt.plot(temp_list_sorted, relativistic_peak_wl_list_nm, 'r--', label='Relativistic (Distant Observer)')
plt.title('Wien\'s Law: Peak Wavelength vs. Temperature (Flat vs. Relativistic)')
plt.xlabel('Temperature (K)')
plt.ylabel('Peak Wavelength (nm)')
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()
def explore_quantum_gravitational_differentiation(mass_kg, distance_m_series):
"""
New function to explore the hypothetical relationship between gravity and a quantum state.
"""
print("\n--- Exploring a Hypothetical Gravitational-Quantum Bridge ---")
print(f" Analyzing Quantum State near a k-body with Mass={mass_kg:.2e} kg")
curvature_series = [calculate_gravitational_curvature_parameter(mass_kg, d) for d in distance_m_series]
# Let's use a base quantum value, for example, a spin state (up/down).
# We'll use 1.0 as a base value to represent an "undifferentiated" quantum state.
base_quantum_value = 1.0
differentiated_value_series = [quantum_differentiation_hypothetical(c, base_quantum_value) for c in curvature_series]
plt.figure(figsize=(10, 6))
plt.plot(distance_m_series, differentiated_value_series, 'b-')
plt.axvline(x=calculate_schwarzschild_radius(mass_kg), color='r', linestyle='--', label='Schwarzschild Radius')
plt.title('Hypothetical "Quantum Differentiation" vs. Distance from a k-body')
plt.xlabel('Distance from k-body (m)')
plt.ylabel('Hypothetical Differentiated Quantum Value (Unitless)')
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()
# --- Main Program Execution (Modified) ---
if __name__ == "__main__":
print("--- Radiation Analysis Program: Interpreting Data Behavior ---")
print("\nThis program explores two physics concepts in both a standard, flat-spacetime environment and a relativistic one.")
print("1. **Radioactive Decay:** Compares the decay rate of an isotope as seen by a local observer vs. a distant one in a gravitational field.")
print("2. **Blackbody Radiation:** Compares the peak emission wavelength of a blackbody as seen by a local observer vs. a distant one (gravitational redshift).")
print("3. **Hypothetical Quantum-Gravity Bridge:** A new, conceptual section to explore how a quantum state might be 'differentiated' by a gravitational field.\n")
while True:
main_choice = input("Select: (R)elativistic Comparison, (Q)uantum-Gravity, (E)xit? ").strip().lower()
if main_choice == 'e':
break
elif main_choice == 'r':
sub_choice = input("Select Relativistic Scenario: (N)uclear Decay or (B)lackbody Radiation? ").strip().lower()
try:
# Get common relativistic parameters
k_body_mass = float(input("\n Enter mass of the k-body (e.g., a star) in kg: "))
k_body_distance = float(input(" Enter distance from the center of the k-body in meters: "))
except ValueError:
print(" Invalid input. Please enter numerical values for mass and distance.")
continue
if sub_choice == 'n':
print("\n--- Relativistic Nuclear Decay Comparison ---")
isotope_name = input(" Enter radioactive isotope name (e.g., Cesium-137): ").strip()
if isotope_name not in MANUAL_GAMMA_DATA and _get_compound_cid_pubchem(isotope_name) is None:
print(f" '{isotope_name}' data not found.")
continue
try:
initial_activity = float(input(" Enter initial activity in Bq: "))
mean_time = float(input(" Enter mean measurement time (days): "))
std_dev_time = float(input(" Enter standard deviation for time (days): "))
except ValueError:
print(" Invalid input.")
continue
analyze_gamma_curies_comparison(isotope_name, initial_activity, mean_time, std_dev_time, k_body_mass, k_body_distance)
elif sub_choice == 'b':
print("\n--- Relativistic Blackbody Radiation Comparison ---")
temp_input_str = input(" Enter temperatures (K) separated by commas (e.g., 300, 1000, 5000): ")
try:
temperatures_list = [float(t.strip()) for t in temp_input_str.split(',')]
except ValueError:
print(" Invalid temperature input.")
continue
simulate_blackbody_wien_comparison(temperatures_list, k_body_mass, k_body_distance)
else:
print(" Invalid choice. Please enter 'N' or 'B'.")
elif main_choice == 'q':
print("\n--- Quantum-Gravitational Bridge Exploration ---")
try:
k_body_mass = float(input("\n Enter mass of the k-body (e.g., a black hole) in kg: "))
# We'll generate a range of distances to show the effect
max_distance_factor = float(input(" Enter max distance as a multiple of Schwarzschild Radius (e.g., 5): "))
schwarzschild_radius_val = calculate_schwarzschild_radius(k_body_mass)
if schwarzschild_radius_val == 0:
print(" Invalid mass. Cannot calculate Schwarzschild Radius.")
continue
distance_series = np.linspace(schwarzschild_radius_val + 1e-9, max_distance_factor * schwarzschild_radius_val, 200)
explore_quantum_gravitational_differentiation(k_body_mass, distance_series)
except ValueError:
print(" Invalid input. Please enter numerical values.")
continue
else:
print("Invalid choice. Please enter 'R', 'Q', or 'E'.")
print("\nExiting program. Thank you for exploring these concepts!")
