diff --git a/crates/atomic-math/src/lib.rs b/crates/atomic-math/src/lib.rs index 9b1b87f..ae92b9e 100644 --- a/crates/atomic-math/src/lib.rs +++ b/crates/atomic-math/src/lib.rs @@ -90,7 +90,7 @@ pub enum OrbitalMode { RealChemist(RealOrbitalKind), } -pub fn probability_density( +pub fn wavefunction_value( qn: &QuantumNumbers, mode: &OrbitalMode, z_eff: f64, @@ -101,15 +101,36 @@ pub fn probability_density( qn.validate()?; let r_part = wavefunctions::r_nl(qn.n, qn.l, z_eff, r)?; - let angular_part_sq = match mode { - OrbitalMode::PureEigenstate => spherical_harmonics::y_lm_real_squared(qn.l, qn.m, theta)?, + let angular_part = match mode { + OrbitalMode::PureEigenstate => { + let m_abs = qn.m.unsigned_abs(); + let x = theta.cos(); + let plm = math_utils::associated_legendre(qn.l, m_abs as i32, x)?; + let l_f = qn.l as f64; + let num_fact = math_utils::factorial(qn.l - m_abs)?; + let den_fact = math_utils::factorial(qn.l + m_abs)?; + let norm = (((2.0 * l_f + 1.0) / (4.0 * std::f64::consts::PI)) * (num_fact / den_fact)).sqrt(); + let trig = if qn.m >= 0 { (qn.m as f64 * phi).cos() } else { (qn.m.abs() as f64 * phi).sin() }; + norm * plm * trig + } OrbitalMode::RealChemist(kind) => { - let y = spherical_harmonics::real_orbital_angular(kind, theta, phi); - y * y + spherical_harmonics::real_orbital_angular(kind, theta, phi) } }; - Ok(r_part * r_part * angular_part_sq) + Ok(r_part * angular_part) +} + +pub fn probability_density( + qn: &QuantumNumbers, + mode: &OrbitalMode, + z_eff: f64, + r: f64, + theta: f64, + phi: f64, +) -> Result { + let psi = wavefunction_value(qn, mode, z_eff, r, theta, phi)?; + Ok(psi * psi) } pub fn sample_points( @@ -118,7 +139,7 @@ pub fn sample_points( z_eff: f64, n_points: usize, seed: u64, -) -> Result, String> { +) -> Result, String> { sampling::sample_points_internal(qn, mode, z_eff, n_points, seed) } @@ -145,11 +166,12 @@ pub fn sample_orbital_points( let points = sample_points(&qn, &mode, z_eff, n_points, seed)?; - let mut flat = Vec::with_capacity(points.len() * 3); - for p in points { + let mut flat = Vec::with_capacity(points.len() * 4); + for (p, sign) in points { flat.push(p[0]); flat.push(p[1]); flat.push(p[2]); + flat.push(sign); } Ok(js_sys::Float32Array::from(flat.as_slice())) diff --git a/crates/atomic-math/src/sampling.rs b/crates/atomic-math/src/sampling.rs index 626343b..a63e2dd 100644 --- a/crates/atomic-math/src/sampling.rs +++ b/crates/atomic-math/src/sampling.rs @@ -1,4 +1,4 @@ -use crate::{probability_density, OrbitalMode, QuantumNumbers}; +use crate::{wavefunction_value, OrbitalMode, QuantumNumbers}; struct Lcg { state: u64, @@ -24,7 +24,7 @@ pub fn sample_points_internal( z_eff: f64, n_points: usize, seed: u64, -) -> Result, String> { +) -> Result, String> { if n_points == 0 { return Ok(Vec::new()); } @@ -55,9 +55,9 @@ pub fn sample_points_internal( for &r in &grid_radii { for &theta in &grid_thetas { for &phi in &grid_phis { - let p = probability_density(qn, mode, z_eff, r, theta, phi)?; + let psi = wavefunction_value(qn, mode, z_eff, r, theta, phi)?; let weight = r * r * theta.sin(); - let density = p * weight; + let density = psi * psi * weight; if density > p_max { p_max = density; } @@ -71,9 +71,9 @@ pub fn sample_points_internal( let theta = rng.next_f64() * std::f64::consts::PI; let phi = rng.next_f64() * 2.0 * std::f64::consts::PI; - let p = probability_density(qn, mode, z_eff, r, theta, phi)?; + let psi = wavefunction_value(qn, mode, z_eff, r, theta, phi)?; let weight = r * r * theta.sin(); - let density = p * weight; + let density = psi * psi * weight; if density > p_max { p_max = density; } @@ -97,9 +97,9 @@ pub fn sample_points_internal( let theta = rng.next_f64() * std::f64::consts::PI; let phi = rng.next_f64() * 2.0 * std::f64::consts::PI; - let p = probability_density(qn, mode, z_eff, r, theta, phi)?; + let psi = wavefunction_value(qn, mode, z_eff, r, theta, phi)?; let weight = r * r * theta.sin(); - let density = p * weight; + let density = psi * psi * weight; if density > p_max { p_max = density * 1.2; } @@ -109,7 +109,8 @@ pub fn sample_points_internal( let x = r * theta.sin() * phi.cos(); let y = r * theta.sin() * phi.sin(); let z = r * theta.cos(); - points.push([x as f32, y as f32, z as f32]); + let sign = if psi >= 0.0 { 1.0f32 } else { -1.0f32 }; + points.push(([x as f32, y as f32, z as f32], sign)); } } diff --git a/crates/atomic-math/src/spherical_harmonics.rs b/crates/atomic-math/src/spherical_harmonics.rs index 76b3e5b..f3ef433 100644 --- a/crates/atomic-math/src/spherical_harmonics.rs +++ b/crates/atomic-math/src/spherical_harmonics.rs @@ -1,7 +1,7 @@ use crate::math_utils::{associated_legendre, factorial}; use crate::RealOrbitalKind; -pub fn y_lm_real_squared(l: u32, m: i32, theta: f64) -> Result { +pub fn y_lm_real(l: u32, m: i32, theta: f64, phi: f64) -> Result { let m_abs = m.unsigned_abs(); if m_abs > l { return Err(format!( @@ -15,9 +15,22 @@ pub fn y_lm_real_squared(l: u32, m: i32, theta: f64) -> Result { let l_f = l as f64; let num_fact = factorial(l - m_abs)?; let den_fact = factorial(l + m_abs)?; - let prefactor = ((2.0 * l_f + 1.0) / (4.0 * std::f64::consts::PI)) * (num_fact / den_fact); + let prefactor = (((2.0 * l_f + 1.0) / (4.0 * std::f64::consts::PI)) * (num_fact / den_fact)).sqrt(); + + let phi_part = if m == 0 { + 1.0 + } else if m > 0 { + std::f64::consts::SQRT_2 * (m as f64 * phi).cos() + } else { + std::f64::consts::SQRT_2 * (m.abs() as f64 * phi).sin() + }; - Ok(prefactor * plm * plm) + Ok(prefactor * plm * phi_part) +} + +pub fn y_lm_real_squared(l: u32, m: i32, theta: f64) -> Result { + let y = y_lm_real(l, m, theta, 0.0)?; + Ok(y * y) } pub fn real_orbital_angular(kind: &RealOrbitalKind, theta: f64, phi: f64) -> f64 { diff --git a/package.json b/package.json index 438fd45..e287801 100644 --- a/package.json +++ b/package.json @@ -1,7 +1,7 @@ { "name": "atomic-explorer", "private": true, - "version": "0.2.4", + "version": "0.2.5", "type": "module", "scripts": { "dev": "vite", diff --git a/src-tauri/tauri.conf.json b/src-tauri/tauri.conf.json index 95e0f56..a16faa3 100644 --- a/src-tauri/tauri.conf.json +++ b/src-tauri/tauri.conf.json @@ -1,7 +1,7 @@ { "$schema": "https://schema.tauri.app/config/2", "productName": "Atomic Explorer", - "version": "0.2.4", + "version": "0.2.5", "identifier": "com.arrase.atomic-explorer", "build": { "beforeDevCommand": "npm run dev", diff --git a/src/main.ts b/src/main.ts index 1a9a388..b02ad60 100644 --- a/src/main.ts +++ b/src/main.ts @@ -106,23 +106,34 @@ async function init() { return orbitalRenderer.captureSnapshot(options); }); + let currentLoadRequestId = 0; + const loadOrbital = async (params: ExtendedOrbitalParams) => { + const requestId = ++currentLoadRequestId; try { document.body.classList.add('loading'); orbitalRenderer.updateParams(params); if (params.mode === 'points') { const points = await sampleOrbitalPoints(params); - orbitalRenderer.setPointCloud(points); + if (requestId === currentLoadRequestId) { + orbitalRenderer.setPointCloud(points); + } } else if (params.mode === 'isosurface') { - orbitalRenderer.updateIsosurface(params); + if (requestId === currentLoadRequestId) { + orbitalRenderer.updateIsosurface(params); + } } else if (params.mode === 'raymarching') { - orbitalRenderer.updateRaymarching(params); + if (requestId === currentLoadRequestId) { + orbitalRenderer.updateRaymarching(params); + } } } catch (err) { console.error('Failed to load orbital:', err); } finally { - document.body.classList.remove('loading'); + if (requestId === currentLoadRequestId) { + document.body.classList.remove('loading'); + } } }; diff --git a/src/render/orbital-renderer.ts b/src/render/orbital-renderer.ts index 747c851..83ae469 100644 --- a/src/render/orbital-renderer.ts +++ b/src/render/orbital-renderer.ts @@ -357,29 +357,21 @@ export class OrbitalRenderer { this.scene.add(dirLight2); } - - public setPointCloud(positions: Float32Array): void { + public setPointCloud(buffer: Float32Array): void { this.clearCurrentMesh(); this.currentMode = 'points'; - const count = positions.length / 3; + const count = Math.floor(buffer.length / 4); + if (count === 0) return; + + const positions = new Float32Array(count * 3); const signs = new Float32Array(count); - const { n, l, m, useRealOrbital, zEff } = this.currentParams; for (let i = 0; i < count; i++) { - const px = positions[i * 3]; - const py = positions[i * 3 + 1]; - const pz = positions[i * 3 + 2]; - const r = Math.hypot(px, py, pz); - if (r < 1e-4) { - signs[i] = 1.0; - continue; - } - const theta = Math.acos(Math.max(-1, Math.min(1, pz / r))); - const phi = Math.atan2(py, px); - const R = this.evalRadial(n, l, zEff, r); - const Y = this.evalAngular(l, m, useRealOrbital, theta, phi); - signs[i] = (R * Y >= 0) ? 1.0 : -1.0; + positions[i * 3] = buffer[i * 4]; + positions[i * 3 + 1] = buffer[i * 4 + 1]; + positions[i * 3 + 2] = buffer[i * 4 + 2]; + signs[i] = buffer[i * 4 + 3]; } const geometry = new THREE.BufferGeometry(); @@ -692,6 +684,8 @@ export class OrbitalRenderer { this.currentMode = mode; this.currentParams = mergedParams; + this.clearCurrentMesh(); + if (mode === 'isosurface') { this.updateIsosurface(mergedParams); } else if (mode === 'raymarching') {