diff --git a/.gitignore b/.gitignore index 1121265..6727398 100644 --- a/.gitignore +++ b/.gitignore @@ -17,3 +17,5 @@ coverage/ .claude/settings.local.json +ffi/medicallib.h +todo.md diff --git a/CLAUDE.md b/CLAUDE.md new file mode 100644 index 0000000..6ab4cfa --- /dev/null +++ b/CLAUDE.md @@ -0,0 +1,99 @@ +# CLAUDE.md + +This file provides guidance to Claude Code (claude.ai/code) when working with code in this repository. + +## Build, Test, and Development Commands + +### Core workflows +- **Build library**: `cargo build` (debug) or `cargo build --release` (optimized) +- **Build FFI shared library**: `cargo build --release --features ffi` (generates C header via cbindgen) +- **Run tests**: `cargo test` (all tests) or `cargo test --test ` (single integration test file) +- **Run examples**: `cargo run --example usage`, `cargo run --example patient`, `cargo run --example tracing_demo` +- **Run demo monitor**: `cargo run --example demo_app --features demo-monitor` +- **Benchmarks**: `cargo bench` (runs Criterion benchmarks; results in `target/criterion/`) +- **Format code**: `cargo fmt --all` (run before committing) +- **Lint**: `cargo clippy --all-targets -- -D warnings` (enforce before committing) + +### Packaging binary distributions +- **Linux/macOS**: Run `../scripts/package.sh [target-triple]` from `medicallib_rust/` +- **Windows**: Run `..\scripts\package.ps1 [target-triple]` from `medicallib_rust\` +- Artifacts placed in `medicallib_rust/dist/` as both `.tar.gz` and `.zip` + +## Architecture Overview + +### Core structure +- **`src/lib.rs`**: Crate facade re-exporting public API (`Patient`, `Organ`, types, errors) +- **`src/types.rs`**: Domain types (`Blood`, `BloodPressure`, `OrganType`, `VitalSign`) +- **`src/error.rs`**: `MedicalError` enum using `thiserror` +- **`src/patient.rs`**: `Patient` struct orchestrating `Vec>` with inter-organ signal routing +- **`src/organs/mod.rs`**: `Organ` trait and `OrganInfo`; per-organ modules (heart, lungs, brain, etc.) +- **`src/ekg/mod.rs`**: EKG monitoring system with configurable leads +- **`src/ffi.rs`**: C-compatible FFI layer (feature `ffi`); header auto-generated to `ffi/medicallib.h` by `build.rs` using cbindgen + +### Organ system +- Each organ module implements the `Organ` trait (`id()`, `organ_type()`, `update(dt)`, `summary()`, `as_any()`, `as_any_mut()`) +- `Patient` owns organs as trait objects, calls `update()` each tick, and routes signals between organs +- Inter-organ communication examples: + - Lungs SpO2 → Heart rate adjustment + - Kidneys GFR → Bladder volume accumulation + - Brain autonomic signals → Heart rate variability and respiratory drive + - Bloodstream perfusion → multiple organ metabolic states + +### FFI boundary +- **Opaque handle**: `MLPatient` pointer to boxed `Patient` +- **Functions**: `ml_patient_new/free/update`, `ml_patient_summary`, `ml_patient_organ_summary`, `medicallib_bmi` +- **Memory**: All returned strings allocated by Rust must be freed via `ml_string_free` +- **Build script**: `build.rs` runs cbindgen when `--features ffi` is enabled, updates `ffi/medicallib.h` if changed +- **Validation**: Any change to `src/ffi.rs` must be tested against `examples/c/ffi_example.c` + +### Testing and benchmarks +- **Integration tests**: `tests/` directory (e.g., `tests/patient.rs`, `tests/ffi.rs`) +- **Benchmarks**: `benches/heart.rs` uses Criterion; compare results across runs for performance regressions +- **Unit tests**: In-module `#[cfg(test)]` blocks for focused invariants + +## Important Development Notes + +### When modifying FFI +1. Edit `src/ffi.rs` +2. Run `cargo build --features ffi` to regenerate `ffi/medicallib.h` via cbindgen +3. Verify C example still compiles: `gcc -o ffi_example examples/c/ffi_example.c -L target/release -lmedicallib_rust` +4. Update `cbindgen.toml` only if header generation settings need adjustment + +### Adding new organs +1. Create module in `src/organs/.rs` +2. Define struct implementing `Organ` trait +3. Add module declaration and public re-export in `src/organs/mod.rs` +4. Update `Patient` signal structs and `couple_organs()` method if inter-organ communication needed +5. Add corresponding `OrganType` variant in `src/types.rs` +6. Update FFI constants in `src/ffi.rs` if organ needs FFI exposure + +### Cargo features +- **`serde`**: Enables serialization on public types +- **`ffi`**: Exposes C ABI; triggers cbindgen in `build.rs` +- **`demo-monitor`**: Required for `demo_app` example with colorized dashboard + +### CI workflows +- GitHub Actions: `.github/workflows/ci.yml` +- Gitea Actions: `.gitea/workflows/ci.yml` +- Both run tests, clippy, formatting checks, and cross-platform builds + +## Key Patterns + +### Patient simulation loop +```rust +let mut patient = Patient::new("case-01")?.initialize_default(); +loop { + patient.update(dt_seconds); + let summary = patient.patient_summary(); + // Process vitals... +} +``` + +### Accessing specific organ +```rust +let heart = patient.find_organ_typed::()?; +println!("HR: {}", heart.heart_rate_bpm()); +``` + +### Inter-organ signals +Patient's `couple_organs()` method reads state from one organ (e.g., `Lungs::spo2_pct()`) and injects it into another (e.g., `Heart::receive_spo2_signal()`). Add new signal fields to the internal `*Signals` structs in `patient.rs` when coupling new organ interactions. \ No newline at end of file diff --git a/Cargo.toml b/Cargo.toml index b135259..b303014 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -1,6 +1,6 @@ [package] name = "medicallib_rust" -version = "0.2.0" +version = "0.3.0" edition = "2021" description = "MedicalSim core library rewrite in Rust: basic clinical calculations and types." authors = ["MedicalSim Team"] diff --git a/docs/research_summary.md b/docs/research_summary.md new file mode 100644 index 0000000..ce71699 --- /dev/null +++ b/docs/research_summary.md @@ -0,0 +1,424 @@ +# Research Summary + +## Heart + +**Current simulation coverage** +- Lumped ventricular pump tracks heart rate, arterial pressures, stroke volume, end-diastolic/systolic volumes, ejection fraction, contractility, preload/afterload and systemic vascular resistance, all coupled to cardiac output control (`src/organs/heart.rs:20-218`). +- Baroreflex-like autonomic loop adjusts sinoatrial pacing, AV delay, vascular resistance and contractility while classifying rhythm states and arrhythmia burden (`src/organs/heart.rs:125-260`). +- Coronary perfusion, myocardial oxygen supply/demand balance and stroke work are approximated through simple pressure-driven formulas (`src/organs/heart.rs:199-247`). + +**Physiology findings and observed gaps** +- Model lacks explicit atrial chambers and valve mechanics, yet real hearts rely on four coordinated valves guiding flow through right/left atria and ventricles for one-way circulation and chamber filling.[1][2] +- The cardiac cycle here omits staged systole/diastole events (atrial kick, isovolumic contraction/relaxation, rapid and reduced filling) that shape physiologic pressure-volume loops and heart sounds.[6] +- Electrical conduction is reduced to SA node rate and an AV delay; physiological activation propagates through the bundle of His, bundle branches and Purkinje fibers to synchronize ventricular contraction and prevent dys-synchrony.[3] +- Coronary supply is approximated by linear clamps, whereas left ventricular perfusion primarily occurs during diastole and depends on the aortic diastolic pressure minus LV end-diastolic pressure (coronary perfusion pressure).[4] +- No endocrine or renal modulation is represented even though the renin-angiotensin-aldosterone system governs long-term blood pressure, volume status and sympathetic tone that in turn alter preload/afterload.[5] + +**Opportunities for improvement** +- Introduce discrete atrial compartments and valve state logic (open/closed/regurgitant) to capture atrial kick, regurgitation, shunts and valve pathology scenarios.[1][2][6] +- Extend conduction modeling with explicit His-Purkinje pathways, refractory periods and conduction delays to support bundle branch blocks, paced rhythms and ventricular dyssynchrony cases.[3] +- Replace fixed coronary supply proxies with a diastolic perfusion model that references coronary perfusion pressure, heart rate-dependent diastolic time fraction and metabolic autoregulation.[4] +- Layer hormonal regulation (e.g., RAAS-driven blood volume and vascular tone adjustments) atop existing autonomic controls to reflect chronic adaptations and pharmacologic interventions.[5] +- Incorporate phase-based cardiac cycle state machine (atrial systole, isovolumic contraction, ejection, isovolumic relaxation, filling) to align simulated pressures/volumes with physiologic waveforms and auscultation cues.[6] + +**Sources** +1. Cleveland Clinic - Heart Valves: What They Are and How They Work.[1] +2. Cleveland Clinic - Chambers of the Heart.[2] +3. Cleveland Clinic - Heart Conduction System.[3] +4. StatPearls/NCBI - Coronary Perfusion Pressure.[4] +5. Cleveland Clinic - Renin-Angiotensin-Aldosterone System (RAAS).[5] +6. Merck Manual - Diagram of the Cardiac Cycle.[6] + +## Brain + +**Current simulation coverage** +- Sleep-wake regulation couples circadian phase, homeostatic sleep pressure, and a stage state machine to drive cortical arousal, sleep stage transitions, and EEG proxies (`src/organs/brain.rs:54-219`). +- Brainstem autonomic drive aggregates respiratory, pain, thirst, hunger, and thermoregulatory inputs to modulate sympathetic tone and variability (`src/organs/brain.rs:220-298`). +- Cerebral perfusion, intracranial pressure, oxygenation, and cerebral blood flow are updated each tick via simplified clamp models tied to metabolic demand and autonomic output (`src/organs/brain.rs:299-373`). +- Neurotransmitter proxies for glutamate, GABA, and dopamine modulate metabolic demand, arousal, seizure risk, and cognitive load (`src/organs/brain.rs:374-445`). +- Consciousness index, seizure risk, and syncope propensity summarize global cortical state for downstream consumers (`src/organs/brain.rs:446-480`). + +**Physiology findings and observed gaps** +- The sleep model advances through staged NREM and REM sequences but lacks representation of REM atonia, stage-specific EEG waveforms, and variable cycle length that characterize physiologic sleep architecture.[7] +- Brainstem autonomic control is condensed into a single drive, omitting the reflex circuitry across nucleus tractus solitarii, ventrolateral medulla, and dorsal vagal outputs that regulate cardiovascular and respiratory coupling.[8] +- Cerebral perfusion is linearized around a fixed CPP target, whereas real brains maintain ~50 mL/100 g/min flow using autoregulation across 60-160 mmHg CPP and react strongly to CO2 shifts, ischemic thresholds, and gray/white matter differences.[9] +- Neurotransmitter variables float within heuristically clamped ranges, yet physiological excitatory-inhibitory balance depends on compartmentalized glutamatergic and GABAergic signaling, receptor kinetics, and astrocytic clearance that drive excitotoxic risk.[10][11] +- Hunger, thirst, and thermoregulatory drives evolve independently of endocrine and hypothalamic feedback, diverging from integrated osmo- and volumetric sensing networks that coordinate angiotensin, vasopressin, and circadian cues.[12] + +**Opportunities for improvement** +- Extend the sleep state machine with polysomnographic markers (e.g., REM atonia, spindle counts) and adaptive cycle timing informed by age or prior sleep debt to better match clinical sleep staging.[7] +- Model key brainstem nuclei, allowing baroreflex, chemoreflex, and vagal pathways to feed back into cardiovascular and respiratory organs with latency and gain parameters derived from physiology texts.[8] +- Replace static CPP/CBF clamps with an autoregulatory module that tracks vessel resistance, CO2 reactivity, and ischemic thresholds to capture plateau and failure zones of cerebral flow.[9] +- Introduce neurotransmitter pools and receptor-specific dynamics (AMPA/NMDA, GABA_A/GABA_B) with astrocytic buffering to simulate excitotoxic cascades and pharmacologic interventions.[10][11] +- Tie hunger, thirst, and thermoregulation drives to hypothalamic and endocrine mediators (e.g., ghrelin, vasopressin, angiotensin II) so fluid and energy balance respond to hormonal and circadian signals.[12] + +**Sources** +7. StatPearls - Physiology of Sleep.[7] +8. Frontiers in Physiology - Synaptic Mechanisms Underlying Elevated Sympathetic Outflow.[8] +9. Stroke Manual - Regulation of Cerebral Blood Flow.[9] +10. NCBI Bookshelf - Glutamate and Aspartate Are the Major Excitatory Transmitters in the Brain.[10] +11. StatPearls - GABA Receptor.[11] +12. Springer Review - Thirst: Neuroendocrine Regulation in Mammals.[12] + +## Bladder + +**Current simulation coverage** +- Three-phase state machine (`Filling`, `Voiding`, `PostVoidRefractory`) governs storage dynamics, reflex triggers, and refractory timing to avoid immediate reactivation (`src/organs/bladder.rs:5-133`). +- Afferent stretch and urgency perception derive from volume thresholds and compliance normalization so filling maps to sensation (`src/organs/bladder.rs:58-96`). +- Autonomic (parasympathetic/sympathetic) and somatic drives converge toward phase-specific targets to coordinate detrusor activation with internal and external sphincter tone (`src/organs/bladder.rs:34-88`). +- Pressure dynamics blend passive compliance, abdominal baseline, and detrusor contraction to clamp intravesical pressure during filling and voiding (`src/organs/bladder.rs:97-147`). +- Cortical inhibition and pontine guarding/void loops drive hypogastric, pudendal, and detrusor outputs with metrics exposed through `Bladder::metrics` and FFI for voluntary continence modeling (`src/organs/bladder.rs:170-452`; `src/ffi.rs`).[14][15] + +**Physiology findings and observed gaps** +- Compliance and capacity stay fixed, yet healthy bladders accommodate volume with minimal pressure rise and elevate detrusor pressure only near voiding; sustained storage pressures above safety limits threaten upper tract health.[13][17] +- Guarding-loop magnitudes still use heuristic gains; calibrating cortical-pontine gating and pudendal discharge against human EMG and urodynamic datasets is needed to capture pathology-specific continence changes.[14][15] +- Urgency is volume-only, whereas mature continence uses cortical oversight of the pontine micturition center to suppress reflex voiding until socially appropriate.[14][15] +- Internal and external sphincters are merged, despite smooth-muscle alpha-adrenergic tone and striated, pudendal-innervated control failing independently in disease.[16] +- Urge and micturition thresholds remain static, even though typical reflex activation spans roughly 250-400 mL and shifts with age, hydration, and neurologic status.[14][18] + +**Opportunities for improvement** +- Replace static compliance with a pressure-volume curve that adapts to bladder history, hydration, or pathology while tracking detrusor pressure against safety limits.[13][17] +- Model explicit guarding and voiding reflex pathways (pontine storage/micturition centers, hypogastric, pelvic, pudendal nerves) so autonomic and somatic loops respond to systemic inputs.[14][15] +- Separate internal and external sphincter models with receptor-specific pharmacology to simulate outlet obstruction, pelvic floor dysfunction, or targeted therapies.[16] +- Couple urgency and continence to cortical and behavioral inputs so developmental milestones, stress, or voluntary suppression can modulate voiding thresholds.[14] +- Parameterize urge and micturition thresholds by age, renal output, or neurologic status to span pediatric, neurogenic, and overactive bladder scenarios.[18] + +**Sources** +13. StatPearls - Urodynamic Testing and Interpretation.[13] +14. StatPearls - Physiology, Urination.[14] +15. Nature Reviews Neuroscience - The Neural Control of Micturition.[15] +16. StatPearls - Anatomy of the bladder wall and sphincter receptor distributions.[16] +17. Physiological Reviews - Detrusor mechanics and compliance regulation in health and disease.[17] +18. Indiana University Pressbooks - Typical micturition reflex thresholds and nerve pathways.[18] + +## Bloodstream + +**Current simulation coverage** +- Maintains plasma volume, red cell volume, total circulating volume, hematocrit, and hemoglobin while syncing cardiac output, mean arterial pressure, and oxygen saturation targets (`src/organs/bloodstream.rs:18-170`). +- Computes arterial/venous oxygen content, delivery, consumption, supply-demand ratio, and circulation time, feeding perfusion and metabolic state classifiers (`src/organs/bloodstream.rs:70-210`). +- Aggregates metabolic waste load, lactate, pH proxy, temperature, glucose, and clearance indices for renal and hepatic pathways (`src/organs/bloodstream.rs:120-270`). +- Maintains albumin and globulin pools with hepatic synthesis targets, Starling oncotic pressure terms, lymphatic return modulation, and edema risk scoring in the plasma volume controller (`src/organs/bloodstream.rs:170-230`).[19][20] +- Tracks erythrocyte age cohorts with spleen-mediated clearance, platelet tagging, and reticuloendothelial iron recycling feeding hepatic stores and transferrin saturation metrics (`src/organs/bloodstream.rs:240-320`).[21][22] +- Couples hypoxia-inducible EPO drive with micronutrient sufficiency (iron saturation, folate, cobalamin) to gate erythropoiesis and reticulocyte surges (`src/organs/bloodstream.rs:320-360`).[23][24] +- Integrates platelet mass, coagulation factor activity, fibrinogen dynamics, fibrinolysis feedback, and thrombosis/bleeding risk indices modulated by splenic platelet reservoirs (`src/organs/bloodstream.rs:360-820`).[25][26] +- Drives leukocyte, differential counts, complement activation, and inflammation indices from spleen immune activity for downstream organs and telemetry surfaces (`src/organs/bloodstream.rs:720-860`).[27][28] +- Replaces the static pH clamp with bicarbonate, base excess, anion gap, arterial PCO2, and lactate guided respiratory/renal compensation to update systemic pH targets (`src/organs/bloodstream.rs:851-910`).[29][30] +- Sets ventilation-perfusion, pulmonary gas exchange targets, and renal/hepatic clearance goals to coordinate with lung and kidney controllers (`src/organs/bloodstream.rs:150-310`). + +**Physiology findings and observed gaps** +- Acute-phase shifts, capillary leak syndromes, and protein-losing pathologies are not yet parameterized, limiting how the new oncotic module responds to inflammation or liver dysfunction.[19][20] +- Cohort-based erythrocyte turnover lacks explicit macrophage phenotypes, hemolysis triggers, or disease-specific remodeling compared with observed splenic clearance pathways.[21][22] +- Hemostasis dynamics still rely on generic activation/fibrinolysis curves and omit factor-specific deficiencies, platelet granule secretion, and transfusion or antithrombotic therapy responses.[25][26] +- Leukocyte modeling excludes adaptive lymphocyte subsets, cytokine networks, and pathogen-specific complement cascades necessary for sepsis or immunosuppression case studies.[27][28] +- Acid-base controller does not yet simulate strong ion difference, renal ammoniagenesis, or mixed metabolic-respiratory disorders beyond linear compensation curves.[29][30] + +**Opportunities for improvement** +- Calibrate acute-phase protein responses and capillary permeability effects within the new oncotic and lymphatic model to capture edema and hypoalbuminemia scenarios.[19][20] +- Add macrophage/monocyte phenotypes, hemolysis triggers, and disease-specific erythrocyte remodeling pathways to the cohort turnover logic.[21][22] +- Expand coagulation modeling with von Willebrand interactions, factor-specific deficits, platelet granule release, and transfusion protocols to reflect trauma and anticoagulation management.[25][26] +- Introduce cytokine-mediated leukocyte recruitment, adaptive immune compartments, and pathogen load feedback to enrich inflammation signaling.[27][28] +- Extend acid-base buffering with strong ion difference, renal ammoniagenesis, and ventilatory control loops tied to organ dysfunction scenarios.[29][30] + +**Sources** +19. Hahn RG. Plasma Volume Oscillations Induced by Hyperoncotic Albumin Infusion. Life (Basel). 2025;15(1):111.[19] +20. Wu JW, Mack GW. Effect of lymphatic outflow on albumin flux from exercising skeletal muscle. J Appl Physiol. 2001;90(5):1912-1918.[20] +21. Vautrinot J, Poole AW. Platelets mediate the clearance of senescent red blood cells. Blood. 2024;143(7):800-812.[21] +22. Mohandas N, Gallagher PG. Accelerated aging of red blood cells in pathologic states. Blood. 2021;137(18):2429-2437.[22] +23. Peng W, Zhan Y, Yu T, et al. Regulation of erythropoiesis by hypoxia-inducible factors and nutrient availability. BMC Med. 2024;22(1):194.[23] +24. Coneyworth LJ, Ford D, Mathers JC. Vitamin B12 and folate interactions in erythropoiesis and neurological function. Nutrients. 2023;15(5):1120.[24] +25. Real-time imaging of platelet-initiated plasma clot formation and lysis unveils distinct impacts of anticoagulants. Res Pract Thromb Haemost. 2024.[25] +26. Thrombopoiesis regulation by hepatic thrombopoietin and splenic clearance ensures platelet homeostasis. J Thromb Thrombolysis. 2025.[26] +27. The role of complement in thromboinflammation. J Trauma Acute Care Surg. 2024.[27] +28. Complement orchestrates innate immune cell differentiation in sepsis. iScience. 2025.[28] +29. Limitations of serum bicarbonate in the ED for diagnosing acid-base disorders. J Emerg Med. 2024.[29] +30. A physiology-based approach to acid-base disorders. BJA Educ. 2024.[30] + +## Esophagus + +**Current simulation coverage** +- State machine cycles swallow initiation, primary/secondary peristalsis, clearing, and reflux exposure while tracking bolus volume and peristaltic progress (`src/organs/esophagus.rs:105-218`). +- Swallow drive integrates oral dryness and mucosal irritation to retune swallow intervals, vagal tone, and peristaltic strength (`src/organs/esophagus.rs:113-131`). +- Lower and upper esophageal sphincter tones adapt via stage-specific modifiers and approach dynamics (`src/organs/esophagus.rs:133-150`). +- Hiatal pressure gradient and reflux propensity are blended with sphincter tone to govern reflux transitions and event rates (`src/organs/esophagus.rs:153-244`). +- Acid balance updates saliva buffering, mucosal integrity, luminal pH, and estimated reflux frequency each tick (`src/organs/esophagus.rs:221-249`). + +**Physiology findings and observed gaps** +- Physiologic swallowing relies on proximal striated and distal smooth muscle with deglutitive inhibition orchestrated by nucleus ambiguus and enteric circuits, but the model uses a single peristaltic strength scalar without segmental timing.[26][27] +- Primary peristaltic waves in healthy adults traverse ~3-6 cm/s with durations modulated by bolus consistency and posture, whereas wave speeds and bolus emptying here stay fixed regardless of load or position.[27][28] +- Esophageal acid clearance depends on saliva-stimulated secondary swallows, gravity assistance, and esophageal shortening; current logic omits clearance latency, posture effects, and bicarbonate secretion variability.[29][30] +- The anti-reflux barrier combines intrinsic LES tone with diaphragmatic crural pinch and transient LES relaxations triggered by gastric distention, yet the simulation only scales a hiatal gradient without diaphragmatic coupling or TLESR triggers.[31] + +**Opportunities for improvement** +- Split the tube into proximal striated and distal smooth segments with deglutitive inhibition timing and enteric reflex loops to cover neurogenic dysphagia and achalasia variants.[26][27] +- Parameterize peristaltic velocity and bolus transport against meal consistency, volume, and posture, and allow failed primary waves to spawn variable secondary peristalsis.[27][28] +- Extend acid clearance to track saliva flow, sequential swallows, gravitational drainage, and bicarbonate buffering so supine, xerostomia, and nocturnal reflux scenarios emerge.[29][30] +- Model diaphragmatic contributions and transient LES relaxation triggers tied to gastric load, belching, and vagal reflexes to enable GERD, hiatal hernia, and fundoplication training cases.[31] + +**Sources** +26. StatPearls - Physiology, Esophagus.[26] +27. GI Motility Online - Physiology of esophageal motility.[27] +28. TSRA Primer - Esophageal Motility & Function Testing.[28] +29. PubMed - Esophageal acid clearance testing and clinical significance.[29] +30. PubMed - Salivary bicarbonate secretion in gastroesophageal reflux disease.[30] +31. Gastroenterology & Hepatology - Gastroesophageal Reflux Disease: Pathophysiology.[31] +## Gallbladder + +**Current simulation coverage** +- Tracks bile reservoir dynamics (volume, acid concentration, hepatic inflow, bile acid pool, recycling efficiency, mucosal absorption, and gallstone index) within the gallbladder struct (`src/organs/gallbladder.rs:24-102`). +- Meal-drive controller sequences fasting clock, meal signal decay, CCK level targeting, and vagal tone adjustments to gate activation cues (`src/organs/gallbladder.rs:103-145`). +- Phase state machine (Filling -> Primed -> Contraction -> Expulsion -> Recovery) tunes sphincter of Oddi tone and bile outflow during each stage (`src/organs/gallbladder.rs:146-208`). +- Bile pool updater concentrates or dilutes bile, clamps cholesterol saturation, and updates the gallstone nucleation index for downstream risk reporting (`src/organs/gallbladder.rs:209-233`). + +**Physiology findings and observed gaps** +- Real gallbladders absorb about 90% of bile water during fasting and eject 50-75% of their contents when CCK triggers contraction alongside sphincter of Oddi relaxation; the model uses fixed constants rather than hormone- and pressure-driven coupling.[32] +- Interdigestive motility features motilin-driven emptying pulses and sphincter phasic contractions preceding migrating motor complex phase III, but the simulation lacks motilin signaling and oscillatory sphincter behavior.[33][37] +- Enterohepatic circulation recycles a roughly 3 g bile acid pool 4-12 times per day with 95% ileal reabsorption; the current model does not exchange bile acids with intestinal or hepatic compartments, obscuring pool depletion or malabsorption states.[34] +- Gallstone formation requires cholesterol supersaturation, mucin-mediated nucleation, and gallbladder stasis, whereas the simulation reduces risk to a single scalar that omits mucin dynamics, bile composition shifts, and sludge progression.[35][36] + +**Opportunities for improvement** +- Replace fixed absorption and outflow clamps with transport models that concentrate bile via electrolyte exchange, allow fractional emptying, and coordinate sphincter relaxation with CCK levels to match post-prandial kinetics.[32] +- Introduce motilin and vagovagal reflex inputs that trigger intermittent fasting-phase contractions and modulate sphincter of Oddi tone, enabling biliary dyskinesia and post-cholecystectomy motility scenarios.[33][37] +- Link the gallbladder bile acid pool to a shared enterohepatic circuit with hepatic synthesis, intestinal reuptake, and fecal losses so bile acid sequestrants or ileal disease deplete the pool and alter digestion.[34] +- Expand the gallstone framework to track bile lipid composition, mucin secretion, sludge accumulation, and stasis duration to differentiate cholesterol versus pigment stone risks and evaluate preventive therapies.[35][36] + +**Sources** +32. Merck Manual Professional Edition - Overview of Biliary Function.[32] +33. PubMed - Cyclic motility of the sphincter of Oddi.[33] +34. PMC - Nuclear receptor control of enterohepatic circulation.[34] +35. StatPearls - Gallstones (Cholelithiasis).[35] +36. PubMed - Role of gallbladder mucin in pathophysiology of gallstones.[36] +37. PubMed - Differential effects of motilin on interdigestive motility.[37] + + + +## Intestines + +**Current simulation coverage** +- Maintains macronutrient absorption, electrolyte reclamation, water reuptake, bile acid recycling, microbiome balance, and inflammatory tone fields within the intestinal struct (`src/organs/intestines.rs:15-113`). +- Internal feeding clock injects nutrient loads, updates motilin and GLP-1 proxies, and drives phase transitions among interdigestive, fed, MMC, ileal brake, and dysmotility states (`src/organs/intestines.rs:111-173`). +- Motility updater blends peristaltic and segmentation indices, lumen volume, and enteric tone to set motility_index, segmentation_index, and mmc_activity scaling (`src/organs/intestines.rs:174-197`). +- Absorption, microbiome, and mucosal routines convert nutrient energy into carbohydrate/fat/protein uptake, adjust electrolyte-water handling, SCFA production, pH, mucosal integrity, inflammation, and bile acid recirculation (`src/organs/intestines.rs:198-287`). + +**Physiology findings and observed gaps** +- Physiologically the small intestine absorbs ~95% of carbohydrates and proteins, 90% of water, and recovers bile salts in the ileum while the colon reclaims the last 1-2 L with electrolyte-coupled transport; the model uses fixed rates without segment-specific transporters or fluid budgets.[38][43] +- MMC cycles repeat every ~90-120 minutes with distinct Phase I-IV patterns governed by motilin pulses that prevent bacterial overgrowth; simulation only tracks a scalar mmc_activity and phase enum without propagating cyclical motor waves or bacterial clearing.[39] +- Ileal brake hormones GLP-1 and PYY respond to distal nutrient exposure to slow gastric emptying and upper gut motility, yet the model lacks nutrient-specific triggers or feedback to upstream organs despite tracking hormone_glp1.[40][41] +- Terminal ileum bile-salt reabsorption and vitamin B12 uptake depend on mucosal transporters and flow, whereas the simulation clamps bile acid recirculation to a motility-dependent target without hepatic pool linkage or malabsorption states.[38] +- Microbiota ferment 30 g/day of carbohydrates to produce ~300 mmol SCFAs that drive sodium/water absorption and fuel colonocytes; current logic approximates SCFA generation linearly from fiber load and does not expose ion transport coupling or microbial community shifts.[42][43] + +**Opportunities for improvement** +- Partition the intestine into duodenal, jejunal, ileal, and colonic segments with transporter-limited nutrient and water absorption, dynamic bile acid pools, and luminal fluid accounting to capture malabsorption and diarrhea phenotypes.[38][43] +- Implement MMC phase cycling driven by motilin bursts with spatial propagation, allowing fasting length, vagal tone, or opioids to disrupt waves and precipitate SIBO scenarios.[39] +- Couple nutrient sensing to GLP-1/PYY secretion, feeding-clock intervals, and feedback onto stomach, pancreas, and gallbladder controllers so fat vs carbohydrate loads elicit distinct ileal brake responses.[40][41] +- Expand microbiome modeling to track fiber species, SCFA spectra, lumen pH, and epithelial fuel utilization, enabling dysbiosis, antibiotic, or prebiotic interventions to shift absorption and mucosal integrity.[42] + +**Sources** +38. StatPearls - Physiology, Small Bowel.[38] +39. Gastroenterology & Hepatology Board Review - Small Intestinal Motility Disorders.[39] +40. PMC - Effects of GLP-1 and incretin-based therapies on gastrointestinal motor function.[40] +41. PubMed - PYY and GLP-1 contribute to feedback inhibition from the canine ileum and colon.[41] +42. PubMed - Colonic health: fermentation and short chain fatty acids.[42] +43. PubMed - Colonic absorption: the importance of short chain fatty acids in man.[43] +## Kidneys + +**Current simulation coverage** +- Tracks glomerular filtration rate, renal plasma flow, filtration fraction, osmolality metrics, and endocrine proxies within the kidney struct (`src/organs/kidneys.rs:14-94`). +- Autoregulation routine classifies perfusion into Autoregulated, Hypoperfused, Hyperperfused, or Obstructed states based on renal plasma flow and obstruction heuristics (`src/organs/kidneys.rs:100-151`). +- Perfusion and hormonal controllers adjust sympathetic tone, renin release, aldosterone drive, and ADH sensitivity in response to plasma volume and osmolality (`src/organs/kidneys.rs:115-151`). +- Tubular handling, acid-base, and erythropoietin updates set sodium reabsorption, potassium and urea excretion, urine flow/osmolality, serum osmolality, plasma volume, and EPO secretion each tick (`src/organs/kidneys.rs:152-231`). + +**Physiology findings and observed gaps** +- Healthy kidneys filter ~120 mL/min from ~600-720 mL/min renal plasma flow and hold GFR constant across 80-180 mmHg via coupled myogenic and macula densa feedback; the model uses fixed thresholds without afferent/efferent resistance dynamics or nephron flow sensing.[44][45] +- Segmental transport normally reclaims ~65-70% of sodium and water in the proximal tubule, ~25% in Henle segments, and fine-tunes electrolytes distally; the simulation collapses everything into a single tubular reabsorption fraction, limiting malabsorption or diuretic scenarios.[46] +- Urine concentration depends on the countercurrent multiplier, medullary gradient maintenance, and ADH-gated aquaporins to span ~50-1200 mOsm; current clamps ignore gradient washout and aquaporin trafficking.[47] +- Renal acid-base balance requires bicarbonate reclamation plus ammonium and titratable acid excretion across segments, whereas the model condenses compensation into a single acid_base_balance scalar.[48] +- Juxtaglomerular renin release integrates baroreceptor, macula densa, and sympathetic inputs, yet renin and aldosterone here follow simplified algebra that cannot capture RAAS pharmacology or dysregulation.[49] +- Erythropoietin secretion arises from hypoxia-responsive peritubular interstitial cells and intercalated cells, but the simulation ties EPO to a coarse renal oxygenation clamp only.[50] + +**Opportunities for improvement** +- Add afferent/efferent arteriole resistance modeling with myogenic and tubuloglomerular feedback loops plus RAAS modulation to reproduce autoregulatory plateaus and pressure-natriuresis shifts.[44][45][49] +- Break the nephron into proximal, loop, distal, and collecting segments with transporter-limited sodium/water reabsorption and endocrine regulation, enabling segment-specific injuries and diuretic effects.[46][47] +- Implement countercurrent multiplier/ exchanger dynamics with medullary gradient washout, aquaporin trafficking, and osmolality feedback to simulate diabetes insipidus, SIADH, or osmotic diuresis.[47] +- Expand acid-base handling to compute bicarbonate reclamation, titratable acid, and ammonium excretion, linking to respiratory compensation and chronic kidney disease buffering limits.[48] +- Drive erythropoietin output from local oxygen tension, fibrosis, and inflammatory cues to support anemia-of-CKD progression and ESA therapy responses.[50] + +**Sources** +44. StatPearls - Physiology, Glomerular Filtration Rate.[44] +45. PubMed - Renal autoregulation in health and disease.[45] +46. PMC - Mechanistic insights into renal ion and water transport in the distal nephron.[46] +47. Kidney: Physiology of the Tubular Reabsorption.[47] +48. PubMed - Acid-Base Homeostasis.[48] +49. StatPearls - Physiology, Renin Angiotensin System.[49] +50. PubMed - Renal epithelium regulates erythropoiesis via HIF-dependent suppression of erythropoietin.[50] + +## Liver + +**Current simulation coverage** +- Multi-state hepatic controller transitions between postabsorptive, fed, fasting, acute phase, and regenerating modes while updating glycogen, lipids, ammonia clearance, and hormone signals each tick (`src/organs/liver.rs:1`). +- Meal-driven hormone routine modulates insulin, glucagon, and cortisol proxies alongside glycogenolysis, gluconeogenesis, and lipogenesis rates to shape fuel handling (`src/organs/liver.rs:73`). +- Bile synthesis, secretion, detox capacity, Kupffer activation, and portal hemodynamics are adjusted through bile/enzymatic update loops (`src/organs/liver.rs:176`). +- Protein synthesis block tracks albumin, clotting factor output, and hepatic fat fraction to summarize synthetic function (`src/organs/liver.rs:214`). + +**Physiology findings and observed gaps** +- The simulation uses static portal and arterial flow clamps, but healthy livers receive ~80% portal venous inflow with a hepatic arterial buffer response that compensates dynamically for portal changes.[51][52] +- Hepatic lobules exhibit periportal-pericentral zonation governing carbohydrate, lipid, and ammonia metabolism, whereas the model treats hepatocytes as a single compartment.[55] +- Enterohepatic bile acid cycling involves secretion, ileal reabsorption, and hepatic reconjugation; current logic generates bile locally without coupling to intestinal pools or transporter limits.[54] +- Fasting gluconeogenesis in humans scales from ~1 to 2 mg·kg⁻¹·min⁻¹ with prolonged fasts driving Cori-cycle recycling; fixed-rate clamps in the model underrepresent these range shifts.[56] +- Kupffer cell activation drives cytokine and acute phase cascades based on pattern-recognition signaling; present implementation raises a scalar without linking to inflammatory inputs or downstream APR protein switches.[53] + +**Opportunities for improvement** +- Implement dual-inflow hemodynamics with arterial buffer feedback and sinusoidal resistance modulation to recreate portal hypertension and ischemia scenarios.[51][52] +- Add zonated hepatocyte segments with zone-specific enzyme sets (urea cycle periportal, glycolysis/pericentral lipogenesis) to capture differential injury and drug metabolism.[55] +- Couple bile acid production to an enterohepatic pool shared with the intestines, including transporter saturation and fecal loss terms for cholestasis modeling.[54] +- Replace fixed gluconeogenesis/glycogenolysis clamps with hormone- and substrate-responsive pathways calibrated to fasting studies, enabling hypoglycemia and stress testing.[56] +- Expand Kupffer signaling to accept pathogen/toxin inputs and drive acute phase protein synthesis, oxidative stress, and stellate activation responses.[53] + +**Sources** +51. World Journal of Gastroenterology - Liver hemodynamics reference values.[51] +52. American Society of Anesthesiologists - Hepatic arterial buffer response overview.[52] +53. StatPearls - Physiology, Liver.[53] +54. StatPearls - Physiology, Bile Secretion.[54] +55. Elsevier - Zonation of hepatic fatty acid metabolism review.[55] +56. American Journal of Physiology - Gluconeogenesis and the Cori cycle in prolonged fasting.[56] + +## Lungs + +**Current simulation coverage** +- Breathing phase loop advances inhalation, exhalation, and pause states while updating diaphragm kinematics, tidal volume, and respiratory rate targets (`src/organs/lungs.rs:37`). +- Ventilatory state machine shifts between resting, hypercapnic, hypoxic, exercise, and distress modes with chemoreceptor and muscle drive scalars (`src/organs/lungs.rs:108`). +- Gas exchange routine adjusts alveolar PO₂/PCO₂, shunt fraction, SpO₂, and CO₂ elimination based on ventilation and metabolic demand (`src/organs/lungs.rs:246`). +- Pulmonary vascular block tunes pulmonary artery and wedge pressures along with V/Q ratio and dead space fractions (`src/organs/lungs.rs:292`). + +**Physiology findings and observed gaps** +- The model collapses regional ventilation/perfusion into single ratios even though human lungs show gravity-dependent gradients (low V/Q at bases, higher at apices) crucial for hypoxemia phenotyping.[57] +- Lung compliance is treated as a scalar, yet in vivo compliance reflects surfactant dynamics, chest wall interaction, and disease-specific hysteresis that shift the pressure–volume curve.[58] +- Diffusing capacity is estimated by linear clamps, whereas normal DL(O₂) ≈ 20–25 ml·min⁻¹·mmHg⁻¹ and increases threefold with exercise due to capillary recruitment—effects absent in the current abstraction.[59][60] +- Chemoreceptor control is condensed into a single drive, but central and peripheral chemoreceptors differ in latency and CO₂/O₂ sensitivity; blended logic masks disorders like carotid body failure.[61] +- Distress flag drives shunt fraction heuristics without modeling airway resistance, surfactant loss, or V/Q scatter that define ARDS phenotypes.[57] + +**Opportunities for improvement** +- Introduce multi-compartment V/Q modeling (apex/mid/base) with gravity and posture modifiers to support shunt-vs-dead-space diagnostics.[57] +- Track static and dynamic compliance separately with surfactant depletion, fibrosis, and chest wall modifiers to capture recruitment and hysteresis behavior.[58] +- Add diffusing-capacity calculations tied to capillary blood volume and exercise state so DL and end-tidal gradients respond to flow changes.[59][60] +- Split chemoreceptor control into central CO₂/pH and peripheral O₂/CO₂ pathways with time constants and hypoxic potentiation to emulate ventilatory drive disorders.[61] +- Expand distress modeling to include airway resistance, alveolar flooding, and recruitable units rather than a single shunt scalar.[57] + +**Sources** +57. StatPearls - Physiology, Pulmonary Ventilation and Perfusion.[57] +58. StatPearls - Physiology, Pulmonary Compliance.[58] +59. MedMuv - Diffusion capacity of the lungs for oxygen.[59] +60. European Respiratory Journal - Reference values for alveolar membrane diffusion capacity.[60] +61. NCBI Bookshelf - Chemical Regulation of Respiration.[61] + +## Pancreas + +**Current simulation coverage** +- State machine toggles basal, postprandial anabolic, hypoglycemic counterregulation, and beta-cell exhaustion modes while updating endocrine outputs (`src/organs/pancreas.rs:15`). +- Meal simulator modulates glucose, incretin, and autonomic tone inputs to drive insulin, glucagon, somatostatin, and pancreatic polypeptide responses (`src/organs/pancreas.rs:55`). +- Exocrine routines adjust enzyme secretion, bicarbonate output, acinar flow, and ductal pressure against incretin and autonomic cues (`src/organs/pancreas.rs:146`). +- Chronic stress logic tracks beta-cell mass fraction and islet stress index to approximate long-term endocrine reserve (`src/organs/pancreas.rs:119`). + +**Physiology findings and observed gaps** +- Ductal secretion depends on CFTR-mediated chloride/bicarbonate exchange achieving ~140 mM bicarbonate, yet the model lacks CFTR or flow-dependent coupling to maintain alkaline secretion.[62][63] +- Incretin physiology (GLP-1/GIP) enhances glucose-stimulated insulin secretion and suppresses glucagon in a glucose-dependent fashion; current logic uses fixed incretin boosts without receptor kinetics.[64] +- Chronic ER stress triggers reversible beta-cell de-differentiation before apoptosis, implying nonlinear mass dynamics beyond the linear decay implemented here.[65] +- Enzyme composition adapts to macronutrient content (e.g., fat increases lipase output), but the simulation scales enzymes uniformly with incretin tone.[66] +- Autonomic tone integrates vagal and sympathetic signaling that co-modulate endocrine and exocrine release; the single autonomic scalar omits frequency-specific vagal bursts and adrenergic suppression.[62] + +**Opportunities for improvement** +- Add CFTR and SLC26 exchanger models with flow-dependent secretion and sensitivity to ductal pressure to emulate cystic fibrosis and pancreatitis.[62][63] +- Implement incretin receptor kinetics with glucose thresholds and pharmacologic agonist profiles to study GLP-1 therapies.[64] +- Model beta-cell stress with reversible identity states and thresholds for apoptosis, capturing adaptation vs. failure under chronic load.[65] +- Differentiate enzyme synthesis pathways for amylase, proteases, and lipase responding to nutrient sensing and CCK feedback.[66] +- Split autonomic inputs into vagal (bursting) and sympathetic (tonic) components to simulate stress-induced endocrine shifts.[62] + +**Sources** +62. Cells - Bicarbonate Transport in the Exocrine Pancreas.[62] +63. NCBI Bookshelf - Water and Ion Secretion from the Pancreatic Ductal System.[63] +64. Diabetes - Multiple actions of GLP-1 on glucose-stimulated insulin secretion.[64] +65. Cell Reports - Adaptation to chronic ER stress enforces beta-cell plasticity.[65] +66. Pancreapedia - Secretion of the human exocrine pancreas in health and disease.[66] + +## Spleen + +**Current simulation coverage** +- Splenic controller tracks immune activity, red pulp volume, platelet reservoir, sympathetic tone, cytokine output, and contraction fraction states (`src/organs/spleen.rs:13`). +- State machine transitions among homeostatic, sympathetic contraction, hyperimmune activation, sequestration, and hypofunction modes driving pulp volumes and cytokines (`src/organs/spleen.rs:45`). +- Contraction and immune routines update platelet release, erythrocyte culling, IgM production, and cytokine signals each tick (`src/organs/spleen.rs:72`). + +**Physiology findings and observed gaps** +- Real spleens hold ~1/3 of circulating platelets and mobilize contracted red pulp during sympathetic surges; model contraction lacks integration with circulating platelet counts or venous return.[69][71] +- White pulp architecture (periarteriolar lymphoid sheaths, marginal zone B cells) drives rapid responses to blood-borne antigens, whereas simulation aggregates immune activity into a single scalar.[67][70] +- Splenic contraction alters hemoglobin and hematocrit during apnea/diving; current logic adjusts reservoir without systemic hematologic feedbacks.[69] +- Marginal zone B cells and macrophages maintain humoral defense; absence of compartment-specific responses limits modeling of asplenia risks.[70] +- Chronic splenomegaly and hypersplenism alter sequestration thresholds; static volume clamps miss disease-dependent compliance changes.[71] + +**Opportunities for improvement** +- Couple splenic reservoir to systemic platelet/erythrocyte pools and venous return to reflect contraction-induced hematologic shifts.[69][71] +- Represent white pulp, marginal zone, and red pulp compartments with discrete immune cell populations and activation kinetics.[67][70] +- Add sympathetic burst inputs tied to cardiovascular simulations to trigger dynamic spleen contraction during exercise or hemorrhage.[69] +- Model splenic compliance changes in infiltrative disease to capture hypersplenism and cytopenia risks.[71] +- Track antigen capture and antibody production timelines to assess vaccine efficacy in asplenic states.[67][70] + +**Sources** +67. StatPearls - Physiology, Spleen.[67] +68. Kenhub - Microscopic anatomy of the spleen.[68] +69. Journal of Applied Physiology - Spleen as an erythrocyte reservoir during diving responses.[69] +70. Journal of Immunology - B cells are crucial for splenic marginal zone development.[70] +71. StatPearls - Splenomegaly.[71] + +## Spinal Cord + +**Current simulation coverage** +- Spinal cord organ tracks signal integrity, ascending/descending conduction velocities, reflex gain, autonomic outflows, locomotor CPG tone, nociceptive facilitation, and perfusion metrics (`src/organs/spinal_cord.rs:14`). +- State machine differentiates intact, concussed, inflammatory, ischemic, and neurogenic shock presentations with corresponding autonomic targets (`src/organs/spinal_cord.rs:48`). +- Integrity and perfusion routines update glial scar index, inflammation, sympathetic/parasympathetic outputs, and locomotor CPG tone over time (`src/organs/spinal_cord.rs:71`). + +**Physiology findings and observed gaps** +- Acute SCI guidelines recommend maintaining MAP 75–95 mmHg for 3–7 days and targeting spinal cord perfusion pressure ≥60–65 mmHg; model perfusion responds to heuristics without explicit SCPP control.[72][73] +- Locomotor central pattern generators reside in segmental networks coordinating limb pairs; current CPG tone scalar lacks limb-specific oscillators and interlimb coordination.[74] +- Human corticospinal conduction velocities average ~67 m/s with disease-dependent slowing; model uses unvalidated values without temperature or demyelination effects.[75] +- Glial scar formation involves astrocyte proliferation over 1–2 weeks with STAT3 signaling, producing barriers and cytokine gradients; simulation increments a scar index without temporal staging or astrocyte roles.[76] +- Neurogenic shock features sympathetic failure, hypotension, and bradycardia; present autonomic outputs shift but are not coupled to cardiovascular modules for systemic responses.[72][73] + +**Opportunities for improvement** +- Add SCPP calculations (MAP − intrathecal pressure) with targets from acute SCI guidelines and allow vasopressor/CSF drainage interventions.[72][73] +- Build bilateral limb CPG models with flexor/extensor half-centers and commissural coupling to study gait adaptations.[74] +- Parameterize conduction velocities with temperature, demyelination, and injury length to match clinical evoked potential data.[75] +- Stage glial scar development (hours–weeks) with astrocyte, fibroblast, and inflammatory cell modules influencing regeneration and cytokines.[76] +- Link autonomic outputs to cardiovascular simulations to reproduce neurogenic shock hemodynamics and vasopressor responses.[72][73] + +**Sources** +72. Neurology Practice Guideline - Hemodynamic management of acute spinal cord injury.[72] +73. Journal of Anesthesia, Analgesia and Critical Care - Hemodynamic management narrative review.[73] +74. Wikipedia - Central pattern generator physiology summary with vertebrate locomotion focus.[74] +75. Journal of Neurology, Neurosurgery & Psychiatry - Motor conduction velocity in the human spinal cord.[75] +76. Cells - Current advancements in spinal cord injury glial scar research.[76] + +## Stomach + +**Current simulation coverage** +- Gastric organ tracks phase (fasting, cephalic, gastric, intestinal, delayed emptying) with vagal tone, hormonal outputs, acid level, motility indices, and gastric volume (`src/organs/stomach.rs:7`). +- Meal routine adjusts ghrelin, gastrin, vagal drive, and nutrient load timing to shift phases and target meal intervals (`src/organs/stomach.rs:63`). +- Secretory update sets acid output, histamine, somatostatin, mucus, and intrinsic factor proxies while motility block governs antral pump strength and emptying rate (`src/organs/stomach.rs:118`). + +**Physiology findings and observed gaps** +- Gastrin stimulates ECL-cell histamine release while somatostatin provides paracrine inhibition; current model applies direct scalar adjustments without receptor-mediated feedback loops.[77][78] +- Ghrelin rises pre-meal and is suppressed by nutrient load, integrating with hypothalamic circuits; simulation lowers ghrelin via simple volume clamps lacking macronutrient and circadian modulation.[79] +- MMC phases originate in stomach/duodenum via motilin and 5-HT feedback; model phases do not enforce interdigestive MMC cycling or motilin triggers.[80] +- Human gastric emptying delivers ~2–3 kcal·min⁻¹ to duodenum with slowing as energy density rises; the current emptying rate is set by motility index without energy-density feedback.[81] +- Gastric emptying kinetics differ for liquids vs solids and respond to osmolarity; single clamp cannot represent nutrient-specific delays.[81] + +**Opportunities for improvement** +- Implement receptor-level gastrin→ECL histamine→parietal pathways with somatostatin inhibition and cholinergic disinhibition loops.[77][78] +- Add ghrelin dynamics tied to macronutrient sensing, sleep state, and leptin/IL-1 signaling to integrate appetite control.[79] +- Introduce interdigestive MMC oscillator driven by motilin and 5-HT with suppression during fed state to coordinate gastric clearing.[80] +- Couple gastric emptying to meal energy density, volume, and macronutrient type to reproduce caloric delivery constraints.[81] +- Differentiate liquid versus solid emptying curves with sieving function and osmolar feedback for hypertonic loads.[81] + +**Sources** +77. Comprehensive Physiology - Gastric peptides gastrin and somatostatin.[77] +78. American Journal of Physiology - Role of histamine in control of gastric acid secretion.[78] +79. Physiological Reviews - Regulation of ghrelin secretion.[79] +80. Neurogastroenterology & Motility - Migrating motor complex control mechanisms.[80] +81. Gastroenterology - Effect of meal volume and energy density on gastric emptying of carbohydrates.[81] diff --git a/examples/c/ffi_example.c b/examples/c/ffi_example.c index 4ab1a7e..39d33d4 100644 --- a/examples/c/ffi_example.c +++ b/examples/c/ffi_example.c @@ -17,6 +17,29 @@ int main(void) { char* sum = ml_patient_summary(p); if (sum) { printf("%s\n", sum); ml_string_free(sum); } + MLBloodstreamMetrics blood = {0}; + if (ml_patient_bloodstream_metrics(p, &blood) == ML_OK) { + printf("Bloodstream: perfusion=%u metabolic=%u oncotic=%.1f mmHg lymph=%.2f mL/min RBC young/mature/senescent=%.0f/%.0f/%.0f%% hif=%.2f iron=%.0f mg\n", + (unsigned)blood.perfusion_state, + (unsigned)blood.metabolic_state, + blood.oncotic_pressure_mm_hg, + blood.lymphatic_return_ml_min, + blood.rbc_young_fraction * 100.0f, + blood.rbc_mature_fraction * 100.0f, + blood.rbc_senescent_fraction * 100.0f, + blood.hif_activation, + blood.iron_store_mg); + } + + MLBladderMetrics bladder = {0}; + if (ml_patient_bladder_metrics(p, &bladder) == ML_OK) { + printf("Bladder: phase=%u urgency=%.0f%% guard=%.0f%% hold=%.0f%%\n", + (unsigned)bladder.phase, + bladder.urgency * 100.0f, + bladder.guarding_reflex_gain * 100.0f, + bladder.voluntary_hold_fraction * 100.0f); + } + char* h = ml_patient_organ_summary(p, ML_ORGAN_HEART); if (h) { printf("%s\n", h); ml_string_free(h); } @@ -24,3 +47,4 @@ int main(void) { return 0; } + diff --git a/ffi/medicallib.h b/ffi/medicallib.h index 6199501..b2fbc77 100644 --- a/ffi/medicallib.h +++ b/ffi/medicallib.h @@ -96,6 +96,46 @@ */ #define ML_ORGAN_BLOODSTREAM 13 +enum MLBladderPhase +#ifdef __cplusplus + : uint32_t +#endif // __cplusplus + { + Filling = 0, + Voiding = 1, + PostVoidRefractory = 2, +}; +#ifndef __cplusplus +typedef uint32_t MLBladderPhase; +#endif // __cplusplus + +enum MLMetabolicState +#ifdef __cplusplus + : uint32_t +#endif // __cplusplus + { + Aerobic = 0, + CompensatedAnaerobic = 1, + AnaerobicCrisis = 2, +}; +#ifndef __cplusplus +typedef uint32_t MLMetabolicState; +#endif // __cplusplus + +enum MLPerfusionState +#ifdef __cplusplus + : uint32_t +#endif // __cplusplus + { + Balanced = 0, + Compensated = 1, + Hypovolemic = 2, + Shock = 3, +}; +#ifndef __cplusplus +typedef uint32_t MLPerfusionState; +#endif // __cplusplus + /** * Opaque patient handle type for C consumers. Wraps a heap-allocated `Patient`. */ @@ -103,6 +143,54 @@ typedef struct MLPatient { void *inner; } MLPatient; +typedef struct MLBloodstreamMetrics { + MLPerfusionState perfusion_state; + MLMetabolicState metabolic_state; + float total_volume_l; + float plasma_volume_l; + float red_cell_volume_l; + float plasma_albumin_g_dl; + float plasma_globulin_g_dl; + float plasma_fibrinogen_g_dl; + float oncotic_pressure_mm_hg; + float lymphatic_return_ml_min; + float rbc_young_fraction; + float rbc_mature_fraction; + float rbc_senescent_fraction; + float erythropoiesis_rate_ml_per_day; + float erythrocyte_clearance_ml_per_day; + float iron_store_mg; + float folate_store_mg; + float b12_store_mcg; + float hif_activation; + float oxygen_supply_demand_ratio; +} MLBloodstreamMetrics; + +typedef struct MLBladderMetrics { + MLBladderPhase phase; + float volume_ml; + float capacity_ml; + float pressure_cm_h2o; + float detrusor_pressure_cm_h2o; + float pressure_safety_limit_cm_h2o; + float afferent_signal; + float urgency; + float cortical_inhibition; + float voluntary_hold_command; + float cortical_gate_fraction; + float voluntary_hold_fraction; + float guarding_reflex_gain; + float pontine_storage_signal; + float pontine_micturition_signal; + float hypogastric_efferent; + float pudendal_efferent; + float parasympathetic_drive; + float sympathetic_drive; + float somatic_drive; + float pressure_stress_index; + float time_above_safety_limit_s; +} MLBladderMetrics; + #ifdef __cplusplus extern "C" { #endif // __cplusplus @@ -139,6 +227,17 @@ void ml_string_free(char *s); */ int32_t ml_patient_update(struct MLPatient *p, float dt_seconds); +/** + * Populate bloodstream metrics for the patient. Returns ML_OK on success. + */ +int32_t ml_patient_bloodstream_metrics(const struct MLPatient *p, + struct MLBloodstreamMetrics *out_metrics); + +/** + * Populate bladder metrics for the patient. Returns ML_OK on success. + */ +int32_t ml_patient_bladder_metrics(const struct MLPatient *p, struct MLBladderMetrics *out_metrics); + /** * Return organ summary by type code. See header for codes. Caller frees string. */ @@ -149,3 +248,4 @@ char *ml_patient_organ_summary(const struct MLPatient *p, uint32_t organ_code); #endif // __cplusplus #endif /* MEDICALLIB_RUST_MEDICALLIB_H */ + diff --git a/src/ffi.rs b/src/ffi.rs index 7526a62..7f8fabc 100644 --- a/src/ffi.rs +++ b/src/ffi.rs @@ -11,6 +11,9 @@ use core::ffi::{c_char, c_void, CStr}; use std::ffi::CString; use std::ptr; +use crate::organs::{ + BladderMetrics, BladderPhase, BloodstreamMetrics, MetabolicState, PerfusionState, +}; use crate::patient::Patient; use crate::types::OrganType; @@ -51,6 +54,202 @@ pub const ML_ORGAN_SPLEEN: u32 = 12; /// Organ code for `OrganType::Bloodstream`. pub const ML_ORGAN_BLOODSTREAM: u32 = 13; +#[repr(u32)] +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum MLBladderPhase { + Filling = 0, + Voiding = 1, + PostVoidRefractory = 2, +} + +impl From for MLBladderPhase { + fn from(phase: BladderPhase) -> Self { + match phase { + BladderPhase::Filling => MLBladderPhase::Filling, + BladderPhase::Voiding => MLBladderPhase::Voiding, + BladderPhase::PostVoidRefractory => MLBladderPhase::PostVoidRefractory, + } + } +} + +#[repr(u32)] +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum MLPerfusionState { + Balanced = 0, + Compensated = 1, + Hypovolemic = 2, + Shock = 3, +} + +impl From for MLPerfusionState { + fn from(state: PerfusionState) -> Self { + match state { + PerfusionState::Balanced => MLPerfusionState::Balanced, + PerfusionState::Compensated => MLPerfusionState::Compensated, + PerfusionState::Hypovolemic => MLPerfusionState::Hypovolemic, + PerfusionState::Shock => MLPerfusionState::Shock, + } + } +} + +#[repr(u32)] +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum MLMetabolicState { + Aerobic = 0, + CompensatedAnaerobic = 1, + AnaerobicCrisis = 2, +} + +impl From for MLMetabolicState { + fn from(state: MetabolicState) -> Self { + match state { + MetabolicState::Aerobic => MLMetabolicState::Aerobic, + MetabolicState::CompensatedAnaerobic => MLMetabolicState::CompensatedAnaerobic, + MetabolicState::AnaerobicCrisis => MLMetabolicState::AnaerobicCrisis, + } + } +} + +#[repr(C)] +#[derive(Debug, Clone, Copy)] +pub struct MLBladderMetrics { + pub phase: MLBladderPhase, + pub volume_ml: f32, + pub capacity_ml: f32, + pub pressure_cm_h2o: f32, + pub detrusor_pressure_cm_h2o: f32, + pub pressure_safety_limit_cm_h2o: f32, + pub afferent_signal: f32, + pub urgency: f32, + pub cortical_inhibition: f32, + pub voluntary_hold_command: f32, + pub cortical_gate_fraction: f32, + pub voluntary_hold_fraction: f32, + pub guarding_reflex_gain: f32, + pub pontine_storage_signal: f32, + pub pontine_micturition_signal: f32, + pub hypogastric_efferent: f32, + pub pudendal_efferent: f32, + pub parasympathetic_drive: f32, + pub sympathetic_drive: f32, + pub somatic_drive: f32, + pub pressure_stress_index: f32, + pub time_above_safety_limit_s: f32, +} + +impl From for MLBladderMetrics { + fn from(metrics: BladderMetrics) -> Self { + Self { + phase: MLBladderPhase::from(metrics.phase), + volume_ml: metrics.volume_ml, + capacity_ml: metrics.capacity_ml, + pressure_cm_h2o: metrics.pressure_cm_h2o, + detrusor_pressure_cm_h2o: metrics.detrusor_pressure_cm_h2o, + pressure_safety_limit_cm_h2o: metrics.pressure_safety_limit_cm_h2o, + afferent_signal: metrics.afferent_signal, + urgency: metrics.urgency, + cortical_inhibition: metrics.cortical_inhibition, + voluntary_hold_command: metrics.voluntary_hold_command, + cortical_gate_fraction: metrics.cortical_gate_fraction, + voluntary_hold_fraction: metrics.voluntary_hold_fraction, + guarding_reflex_gain: metrics.guarding_reflex_gain, + pontine_storage_signal: metrics.pontine_storage_signal, + pontine_micturition_signal: metrics.pontine_micturition_signal, + hypogastric_efferent: metrics.hypogastric_efferent, + pudendal_efferent: metrics.pudendal_efferent, + parasympathetic_drive: metrics.parasympathetic_drive, + sympathetic_drive: metrics.sympathetic_drive, + somatic_drive: metrics.somatic_drive, + pressure_stress_index: metrics.pressure_stress_index, + time_above_safety_limit_s: metrics.time_above_safety_limit_s, + } + } +} + +#[repr(C)] +#[derive(Debug, Clone, Copy)] +pub struct MLBloodstreamMetrics { + pub perfusion_state: MLPerfusionState, + pub metabolic_state: MLMetabolicState, + pub total_volume_l: f32, + pub plasma_volume_l: f32, + pub red_cell_volume_l: f32, + pub plasma_albumin_g_dl: f32, + pub plasma_globulin_g_dl: f32, + pub plasma_fibrinogen_g_dl: f32, + pub oncotic_pressure_mm_hg: f32, + pub lymphatic_return_ml_min: f32, + pub rbc_young_fraction: f32, + pub rbc_mature_fraction: f32, + pub rbc_senescent_fraction: f32, + pub platelet_count_1e3_per_ul: f32, + pub platelet_activation: f32, + pub coagulation_factor_activity: f32, + pub fibrinolysis_activity: f32, + pub thrombosis_risk_index: f32, + pub bleeding_risk_index: f32, + pub wbc_count_giga_per_l: f32, + pub neutrophil_fraction: f32, + pub lymphocyte_fraction: f32, + pub monocyte_fraction: f32, + pub complement_activity: f32, + pub inflammation_index: f32, + pub bicarbonate_mmol_l: f32, + pub base_excess_mmol_l: f32, + pub anion_gap_mmol_l: f32, + pub arterial_pco2_mm_hg: f32, + pub erythropoiesis_rate_ml_per_day: f32, + pub erythrocyte_clearance_ml_per_day: f32, + pub iron_store_mg: f32, + pub folate_store_mg: f32, + pub b12_store_mcg: f32, + pub hif_activation: f32, + pub oxygen_supply_demand_ratio: f32, +} + +impl From for MLBloodstreamMetrics { + fn from(metrics: BloodstreamMetrics) -> Self { + Self { + perfusion_state: MLPerfusionState::from(metrics.perfusion_state), + metabolic_state: MLMetabolicState::from(metrics.metabolic_state), + total_volume_l: metrics.total_volume_l, + plasma_volume_l: metrics.plasma_volume_l, + red_cell_volume_l: metrics.red_cell_volume_l, + plasma_albumin_g_dl: metrics.plasma_albumin_g_dl, + plasma_globulin_g_dl: metrics.plasma_globulin_g_dl, + plasma_fibrinogen_g_dl: metrics.plasma_fibrinogen_g_dl, + oncotic_pressure_mm_hg: metrics.oncotic_pressure_mm_hg, + lymphatic_return_ml_min: metrics.lymphatic_return_ml_min, + rbc_young_fraction: metrics.rbc_young_fraction, + rbc_mature_fraction: metrics.rbc_mature_fraction, + rbc_senescent_fraction: metrics.rbc_senescent_fraction, + platelet_count_1e3_per_ul: metrics.platelet_count_1e3_per_ul, + platelet_activation: metrics.platelet_activation, + coagulation_factor_activity: metrics.coagulation_factor_activity, + fibrinolysis_activity: metrics.fibrinolysis_activity, + thrombosis_risk_index: metrics.thrombosis_risk_index, + bleeding_risk_index: metrics.bleeding_risk_index, + wbc_count_giga_per_l: metrics.wbc_count_giga_per_l, + neutrophil_fraction: metrics.neutrophil_fraction, + lymphocyte_fraction: metrics.lymphocyte_fraction, + monocyte_fraction: metrics.monocyte_fraction, + complement_activity: metrics.complement_activity, + inflammation_index: metrics.inflammation_index, + bicarbonate_mmol_l: metrics.bicarbonate_mmol_l, + base_excess_mmol_l: metrics.base_excess_mmol_l, + anion_gap_mmol_l: metrics.anion_gap_mmol_l, + arterial_pco2_mm_hg: metrics.arterial_pco2_mm_hg, + erythropoiesis_rate_ml_per_day: metrics.erythropoiesis_rate_ml_per_day, + erythrocyte_clearance_ml_per_day: metrics.erythrocyte_clearance_ml_per_day, + iron_store_mg: metrics.iron_store_mg, + folate_store_mg: metrics.folate_store_mg, + b12_store_mcg: metrics.b12_store_mcg, + hif_activation: metrics.hif_activation, + oxygen_supply_demand_ratio: metrics.oxygen_supply_demand_ratio, + } + } +} + /// Opaque patient handle type for C consumers. Wraps a heap-allocated `Patient`. #[repr(C)] #[derive(Debug)] @@ -155,6 +354,56 @@ pub extern "C" fn ml_patient_update(p: *mut MLPatient, dt_seconds: f32) -> i32 { ML_OK } +/// Populate bloodstream metrics for the patient. Returns ML_OK on success. +#[no_mangle] +pub extern "C" fn ml_patient_bloodstream_metrics( + p: *const MLPatient, + out_metrics: *mut MLBloodstreamMetrics, +) -> i32 { + if p.is_null() || out_metrics.is_null() { + return ML_EINVAL; + } + let patient_ptr = unsafe { (*p).inner as *const Patient }; + if patient_ptr.is_null() { + return ML_EINVAL; + } + let patient = unsafe { &*patient_ptr }; + match patient.bloodstream_metrics() { + Some(metrics) => { + unsafe { + *out_metrics = MLBloodstreamMetrics::from(metrics); + } + ML_OK + } + None => ML_ERR, + } +} + +/// Populate bladder metrics for the patient. Returns ML_OK on success. +#[no_mangle] +pub extern "C" fn ml_patient_bladder_metrics( + p: *const MLPatient, + out_metrics: *mut MLBladderMetrics, +) -> i32 { + if p.is_null() || out_metrics.is_null() { + return ML_EINVAL; + } + let patient_ptr = unsafe { (*p).inner as *const Patient }; + if patient_ptr.is_null() { + return ML_EINVAL; + } + let patient = unsafe { &*patient_ptr }; + match patient.bladder_metrics() { + Some(metrics) => { + unsafe { + *out_metrics = MLBladderMetrics::from(metrics); + } + ML_OK + } + None => ML_ERR, + } +} + /// Return organ summary by type code. See header for codes. Caller frees string. #[no_mangle] pub extern "C" fn ml_patient_organ_summary(p: *const MLPatient, organ_code: u32) -> *mut c_char { @@ -191,3 +440,5 @@ pub extern "C" fn ml_patient_organ_summary(p: *const MLPatient, organ_code: u32) Err(_) => ptr::null_mut(), } } + + diff --git a/src/lib.rs b/src/lib.rs index e0f1bef..ba53c72 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -38,7 +38,8 @@ pub mod ffi; pub use crate::ekg::{EkgLead, EkgMonitor, EkgSnapshot, HeartElectricalState}; pub use crate::error::MedicalError; pub use crate::organs::{ - Bloodstream, BreathingPhase, Heart, Lungs, MetabolicState, Organ, PerfusionState, + Bladder, BladderMetrics, BladderPhase, Bloodstream, BreathingPhase, Heart, Lungs, + MetabolicState, Organ, PerfusionState, }; pub use crate::patient::Patient; pub use crate::types::{Blood, BloodPressure, OrganType}; diff --git a/src/organs/bladder.rs b/src/organs/bladder.rs index 79470e2..5242463 100644 --- a/src/organs/bladder.rs +++ b/src/organs/bladder.rs @@ -8,6 +8,124 @@ pub enum BladderPhase { PostVoidRefractory, } +#[derive(Debug, Clone)] +struct ComplianceModel { + low_volume_compliance_ml_per_cm_h2o: f32, + mid_volume_compliance_ml_per_cm_h2o: f32, + high_volume_compliance_ml_per_cm_h2o: f32, + detrusor_safety_limit_cm_h2o: f32, + viscoelastic_relaxation: f32, + stiffness_factor: f32, + stress_adaptation_rate: f32, + recovery_rate: f32, +} + +impl ComplianceModel { + fn new() -> Self { + Self { + low_volume_compliance_ml_per_cm_h2o: 32.0, + mid_volume_compliance_ml_per_cm_h2o: 22.0, + high_volume_compliance_ml_per_cm_h2o: 8.5, + detrusor_safety_limit_cm_h2o: 40.0, + viscoelastic_relaxation: 1.0, + stiffness_factor: 1.0, + stress_adaptation_rate: 0.08, + recovery_rate: 0.015, + } + } + + fn safety_limit(&self) -> f32 { + self.detrusor_safety_limit_cm_h2o + } + + fn effective_compliance(&self, normalized_volume: f32, filling_rate_ml_per_min: f32) -> f32 { + let nv = normalized_volume.clamp(0.0, 2.0); + let logistic_mid = 1.0 / (1.0 + (-10.0 * (nv - 0.55)).exp()); + let logistic_high = 1.0 / (1.0 + (-16.0 * (nv - 0.9)).exp()); + + let base_mid = self.low_volume_compliance_ml_per_cm_h2o + - (self.low_volume_compliance_ml_per_cm_h2o - self.mid_volume_compliance_ml_per_cm_h2o) + * logistic_mid; + let base = + base_mid - (base_mid - self.high_volume_compliance_ml_per_cm_h2o) * logistic_high; + + let rate_penalty = (1.0 - filling_rate_ml_per_min / 240.0).clamp(0.65, 1.0); + let viscoelastic = base * self.viscoelastic_relaxation * rate_penalty; + + (viscoelastic / self.stiffness_factor).clamp(3.0, 80.0) + } + + fn update_state(&mut self, detrusor_pressure: f32, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + let limit = self.detrusor_safety_limit_cm_h2o; + if detrusor_pressure > limit { + let overload = detrusor_pressure - limit; + self.viscoelastic_relaxation = (self.viscoelastic_relaxation + + self.stress_adaptation_rate * dt_seconds * (1.0 + 0.02 * overload)) + .clamp(0.85, 1.35); + self.stiffness_factor = (self.stiffness_factor + + 0.15 * self.stress_adaptation_rate * dt_seconds * (1.0 + overload / limit)) + .clamp(1.0, 1.8); + } else { + self.viscoelastic_relaxation = + (self.viscoelastic_relaxation - self.recovery_rate * dt_seconds).clamp(0.85, 1.35); + self.stiffness_factor = + (self.stiffness_factor - 0.25 * self.recovery_rate * dt_seconds).clamp(1.0, 1.8); + } + } +} + +#[derive(Debug, Clone)] +struct ReflexState { + pelvic_afferent_rate: f32, + guarding_reflex_gain: f32, + pontine_storage_signal: f32, + pontine_micturition_signal: f32, + hypogastric_efferent: f32, + pudendal_efferent: f32, +} + +impl ReflexState { + fn new() -> Self { + Self { + pelvic_afferent_rate: 0.2, + guarding_reflex_gain: 0.4, + pontine_storage_signal: 0.5, + pontine_micturition_signal: 0.1, + hypogastric_efferent: 0.6, + pudendal_efferent: 0.6, + } + } +} + +#[derive(Debug, Clone, Copy)] +pub struct BladderMetrics { + pub phase: BladderPhase, + pub volume_ml: f32, + pub capacity_ml: f32, + pub pressure_cm_h2o: f32, + pub detrusor_pressure_cm_h2o: f32, + pub pressure_safety_limit_cm_h2o: f32, + pub afferent_signal: f32, + pub urgency: f32, + pub cortical_inhibition: f32, + pub voluntary_hold_command: f32, + pub cortical_gate_fraction: f32, + pub voluntary_hold_fraction: f32, + pub guarding_reflex_gain: f32, + pub pontine_storage_signal: f32, + pub pontine_micturition_signal: f32, + pub hypogastric_efferent: f32, + pub pudendal_efferent: f32, + pub parasympathetic_drive: f32, + pub sympathetic_drive: f32, + pub somatic_drive: f32, + pub pressure_stress_index: f32, + pub time_above_safety_limit_s: f32, +} + #[derive(Debug, Clone)] pub struct Bladder { info: OrganInfo, @@ -15,6 +133,8 @@ pub struct Bladder { pub volume_ml: f32, /// Intraluminal/detrusor pressure (cmH2O). pub pressure: f32, + /// Detected detrusor component of intravesical pressure (cmH2O). + pub detrusor_pressure_cm_h2o: f32, /// Normalized afferent firing (0..=1) representing stretch receptor activity. pub afferent_signal: f32, /// Normalized urgency perception (0..=1). @@ -25,8 +145,26 @@ pub struct Bladder { pub capacity_ml: f32, /// Residual volume expected after a complete void (ml). pub residual_volume_ml: f32, - /// Compliance relating volume change to pressure (ml per cmH2O). - pub compliance_ml_per_cm_h2o: f32, + /// Baseline cystometric capacity for demographic scaling (ml). + pub baseline_capacity_ml: f32, + /// Chronological age in years influencing thresholds. + pub age_years: f32, + /// Integrity of supraspinal and spinal pathways (0..=1). + pub neurologic_integrity: f32, + /// Habitual continence training or pelvic floor conditioning (0..=1). + pub continence_training: f32, + /// Descending cortical inhibition of the micturition reflex (0..=1). + pub cortical_inhibition: f32, + /// Voluntary hold motor command (0..=1) driving pelvic floor recruitment. + pub voluntary_hold_command: f32, + /// Filtered cortical gate strength (0..=1) applied to reflex loops. + pub cortical_gate_fraction: f32, + /// Filtered voluntary hold strength (0..=1) applied to sphincter control. + pub voluntary_hold_fraction: f32, + /// Adaptive compliance model governing pressure-volume behavior. + compliance: ComplianceModel, + /// Integrated reflex pathways coordinating storage and voiding. + reflex: ReflexState, /// Baseline pressure generated by abdominal cavity (cmH2O). pub baseline_pressure_cm_h2o: f32, /// Volume threshold at which a full micturition reflex is triggered (ml). @@ -47,6 +185,10 @@ pub struct Bladder { pub sympathetic_drive: f32, /// Somatic drive through the pudendal nerve to the external sphincter (0..=1). pub somatic_drive: f32, + /// Exponentially-weighted overload index tracking detrusor pressure burden (cmH2O·s). + pub pressure_stress_index: f32, + /// Cumulative time the detrusor pressure exceeded the safety limit (seconds). + pub time_above_safety_limit_s: f32, time_since_last_void_s: f32, refractory_seconds: f32, } @@ -56,13 +198,23 @@ impl Bladder { Self { info: OrganInfo::new(id, OrganType::Bladder), volume_ml: 120.0, - pressure: 5.0, + pressure: 8.0, + detrusor_pressure_cm_h2o: 5.0, afferent_signal: 0.1, urgency: 0.1, phase: BladderPhase::Filling, capacity_ml: 500.0, residual_volume_ml: 30.0, - compliance_ml_per_cm_h2o: 18.0, + baseline_capacity_ml: 500.0, + age_years: 38.0, + neurologic_integrity: 1.0, + continence_training: 0.65, + cortical_inhibition: 0.6, + voluntary_hold_command: 0.45, + cortical_gate_fraction: 0.0, + voluntary_hold_fraction: 0.0, + compliance: ComplianceModel::new(), + reflex: ReflexState::new(), baseline_pressure_cm_h2o: 5.0, micturition_threshold_ml: 350.0, urge_threshold_ml: 200.0, @@ -73,7 +225,9 @@ impl Bladder { parasympathetic_drive: 0.05, sympathetic_drive: 0.8, somatic_drive: 0.8, + pressure_stress_index: 0.0, time_since_last_void_s: 0.0, + time_above_safety_limit_s: 0.0, refractory_seconds: 15.0, } } @@ -94,35 +248,237 @@ impl Bladder { } } - fn update_drives(&mut self, dt_seconds: f32) { - let urgency = self.urgency; - let (parasym_target, sym_target, somatic_target) = match self.phase { - BladderPhase::Filling => { - let parasym = urgency.powf(1.2).clamp(0.0, 0.85); - let sym = (1.0 - 0.6 * urgency).clamp(0.2, 1.0); - let somatic = (0.95 - 0.5 * urgency).clamp(0.3, 0.95); - (parasym, sym, somatic) - } - BladderPhase::Voiding => (1.0, 0.1, 0.2), - BladderPhase::PostVoidRefractory => (0.2, 0.6, 0.7), + fn cortical_gate_strength(&self) -> f32 { + let training = 0.4 + 0.6 * self.continence_training; + (self.cortical_inhibition * self.neurologic_integrity * training).clamp(0.0, 1.0) + } + + fn voluntary_hold_strength(&self) -> f32 { + let training = 0.35 + 0.65 * self.continence_training; + (self.voluntary_hold_command * training * self.neurologic_integrity).clamp(0.0, 1.0) + } + + fn update_voluntary_modulators(&mut self, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + let cortical_target = self.cortical_gate_strength(); + self.cortical_gate_fraction = Self::approach( + self.cortical_gate_fraction, + cortical_target, + 1.2, + dt_seconds, + ) + .clamp(0.0, 1.0); + let hold_target = self.voluntary_hold_strength(); + self.voluntary_hold_fraction = + Self::approach(self.voluntary_hold_fraction, hold_target, 1.6, dt_seconds) + .clamp(0.0, 1.0); + } + + fn update_threshold_targets(&mut self, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + let baseline_capacity = self.baseline_capacity_ml.max(320.0); + let age_factor = if self.age_years <= 40.0 { + 1.0 + (40.0 - self.age_years).max(0.0) * 0.0035 + } else { + (1.0 - (self.age_years - 40.0) * 0.003).clamp(0.74, 1.08) + }; + let neurologic_factor = (0.55 + 0.45 * self.neurologic_integrity).clamp(0.45, 1.1); + let capacity_target = + (baseline_capacity * age_factor * neurologic_factor).clamp(240.0, 720.0); + self.capacity_ml = + Self::approach(self.capacity_ml, capacity_target, 0.5, dt_seconds).clamp(160.0, 900.0); + + let diuresis_ratio = (self.filling_rate_ml_per_min / 60.0).clamp(0.25, 2.5); + let diuresis_factor = diuresis_ratio.powf(-0.18).clamp(0.82, 1.18); + let voluntary_buffer = (1.0 + + 0.22 * self.cortical_gate_fraction + + 0.25 * self.voluntary_hold_fraction + + 0.12 * self.continence_training) + .clamp(0.8, 1.5); + let neurologic_drop = (1.0 - 0.3 * (1.0 - self.neurologic_integrity)).clamp(0.55, 1.05); + + let urge_fraction = (0.38 + 0.08 * self.voluntary_hold_fraction + - 0.1 * (1.0 - self.neurologic_integrity)) + .clamp(0.25, 0.6); + let mut urge_target = + capacity_target * urge_fraction * diuresis_factor * voluntary_buffer * neurologic_drop; + urge_target = urge_target.clamp(capacity_target * 0.22, capacity_target * 0.75); + + let mict_fraction = (0.68 + 0.1 * self.voluntary_hold_fraction + - 0.15 * (1.0 - self.neurologic_integrity)) + .clamp(0.5, 0.92); + let mut mict_target = capacity_target + * mict_fraction + * diuresis_factor.powf(0.75) + * voluntary_buffer + * neurologic_drop; + mict_target = mict_target.clamp(urge_target + 25.0, capacity_target * 0.95); + + self.urge_threshold_ml = + Self::approach(self.urge_threshold_ml, urge_target, 0.8, dt_seconds).clamp(50.0, 800.0); + self.micturition_threshold_ml = + Self::approach(self.micturition_threshold_ml, mict_target, 0.7, dt_seconds) + .clamp(60.0, 850.0); + } + + fn update_reflex_gate(&mut self, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + let normalized_volume = (self.volume_ml / self.capacity_ml).clamp(0.0, 1.7); + let detrusor_ratio = + (self.detrusor_pressure_cm_h2o / self.compliance.safety_limit()).clamp(0.0, 2.0); + let cortical_damping = (1.0 - 0.35 * self.cortical_gate_fraction).clamp(0.55, 1.0); + let neuro_sensitivity = (1.0 + 0.25 * (1.0 - self.neurologic_integrity)).clamp(0.8, 1.3); + let pelvic_drive = + (0.32 * normalized_volume + 0.45 * self.afferent_signal + 0.23 * detrusor_ratio) + * cortical_damping + * neuro_sensitivity; + let pelvic_target = pelvic_drive.clamp(0.0, 1.2); + self.reflex.pelvic_afferent_rate = Self::approach( + self.reflex.pelvic_afferent_rate, + pelvic_target.min(1.0), + 1.6, + dt_seconds, + ); + + let voluntary_buffer = (self.voluntary_hold_fraction * 0.45).clamp(0.0, 0.45); + let storage_modulator = + (1.0 - self.urgency * (1.0 - 0.35 * self.cortical_gate_fraction)).clamp(0.0, 1.0); + let guarding_target = if matches!(self.phase, BladderPhase::Voiding) { + 0.05 + } else { + ((self.reflex.pelvic_afferent_rate + 0.12) * storage_modulator).powf(0.7) + + voluntary_buffer + }; + self.reflex.guarding_reflex_gain = Self::approach( + self.reflex.guarding_reflex_gain, + guarding_target.clamp(0.0, 1.0), + 1.4, + dt_seconds, + ); + + let mut supraspinal_trigger = ((self.reflex.pelvic_afferent_rate - 0.55).max(0.0) + * self.urgency.powf(1.2) + + (self.detrusor_pressure_cm_h2o - self.compliance.safety_limit()).max(0.0) / 40.0) + .clamp(0.0, 1.2); + supraspinal_trigger = + (supraspinal_trigger * (1.0 - 0.5 * self.cortical_gate_fraction)).clamp(0.0, 1.2); + + let pontine_void_base = match self.phase { + BladderPhase::Voiding => (0.5 + 0.6 * supraspinal_trigger).clamp(0.0, 1.0), + BladderPhase::PostVoidRefractory => (0.2 * supraspinal_trigger).clamp(0.0, 0.4), + BladderPhase::Filling => (0.35 * supraspinal_trigger).clamp(0.0, 0.85), + }; + let pontine_void_target = + (pontine_void_base * (1.0 - 0.55 * self.cortical_gate_fraction)).clamp(0.0, 1.0); + let pontine_storage_target = if matches!(self.phase, BladderPhase::Voiding) { + (0.1 + 0.2 * (1.0 - supraspinal_trigger) + 0.15 * self.voluntary_hold_fraction) + .clamp(0.0, 0.45) + } else { + (0.45 + 0.5 * self.reflex.guarding_reflex_gain - 0.4 * supraspinal_trigger + + 0.35 * self.cortical_gate_fraction + + 0.2 * self.voluntary_hold_fraction) + .clamp(0.0, 1.0) }; - self.parasympathetic_drive = - Self::approach(self.parasympathetic_drive, parasym_target, 0.8, dt_seconds); - self.sympathetic_drive = - Self::approach(self.sympathetic_drive, sym_target, 0.6, dt_seconds); - self.somatic_drive = Self::approach(self.somatic_drive, somatic_target, 0.9, dt_seconds); + self.reflex.pontine_micturition_signal = Self::approach( + self.reflex.pontine_micturition_signal, + pontine_void_target, + 1.0, + dt_seconds, + ) + .clamp(0.0, 1.0); + self.reflex.pontine_storage_signal = Self::approach( + self.reflex.pontine_storage_signal, + pontine_storage_target, + 1.0, + dt_seconds, + ) + .clamp(0.0, 1.0); + } - self.parasympathetic_drive = self.parasympathetic_drive.clamp(0.0, 1.0); - self.sympathetic_drive = self.sympathetic_drive.clamp(0.0, 1.0); - self.somatic_drive = self.somatic_drive.clamp(0.0, 1.0); + fn update_reflex_efferents(&mut self, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + let sympathetic_bias = (0.6 + 0.2 * self.neurologic_integrity).clamp(0.4, 0.9); + let hypogastric_target = (sympathetic_bias * self.sympathetic_drive + + 0.35 * self.reflex.pontine_storage_signal + - 0.45 * self.reflex.pontine_micturition_signal + + 0.15 * self.cortical_gate_fraction) + .clamp(0.0, 1.0); + let pudendal_target = ((0.55 + 0.15 * self.continence_training) * self.somatic_drive + + 0.25 * self.reflex.guarding_reflex_gain + + 0.35 * self.voluntary_hold_fraction + - 0.55 * self.reflex.pontine_micturition_signal) + * (0.7 + 0.3 * self.neurologic_integrity); + let pudendal_target = pudendal_target.clamp(0.0, 1.0); + self.reflex.hypogastric_efferent = Self::approach( + self.reflex.hypogastric_efferent, + hypogastric_target, + 1.8, + dt_seconds, + ); + self.reflex.pudendal_efferent = Self::approach( + self.reflex.pudendal_efferent, + pudendal_target, + 1.8, + dt_seconds, + ); + } + + fn update_drives(&mut self, dt_seconds: f32) { + let urgency = self.urgency; + let safety_ratio = + (self.detrusor_pressure_cm_h2o / self.compliance.safety_limit()).clamp(0.0, 2.0); + let safety_bias = (safety_ratio.powf(1.2) * 0.35).clamp(0.0, 0.4); + let storage_bias = (self.reflex.pontine_storage_signal * 0.6).clamp(0.0, 0.6); + let void_bias = self.reflex.pontine_micturition_signal; + let cortical_filter = (1.0 - 0.4 * self.cortical_gate_fraction).clamp(0.35, 1.0); + let parasym_target = + (urgency.powf(1.2) * cortical_filter + safety_bias + 0.7 * void_bias).clamp(0.0, 1.0); + let sym_target = (1.0 + - 0.6 * urgency * (1.0 - 0.3 * self.cortical_gate_fraction) + - 0.3 * safety_bias + - 0.6 * void_bias + + storage_bias + + 0.15 * self.cortical_gate_fraction) + .clamp(0.15, 1.0); + let voluntary_boost = + (0.35 * self.voluntary_hold_fraction + 0.2 * self.continence_training).clamp(0.0, 0.45); + let somatic_target = ((0.55 + 0.25 * self.continence_training) * (1.0 - urgency * 0.55) + + storage_bias * 0.45 + + voluntary_boost + - void_bias * 0.5) + * (0.75 + 0.25 * self.neurologic_integrity); + + self.parasympathetic_drive = + Self::approach(self.parasympathetic_drive, parasym_target, 0.8, dt_seconds) + .clamp(0.0, 1.0); + self.sympathetic_drive = + Self::approach(self.sympathetic_drive, sym_target, 0.6, dt_seconds).clamp(0.0, 1.0); + self.somatic_drive = Self::approach( + self.somatic_drive, + somatic_target.clamp(0.2, 1.05), + 0.9, + dt_seconds, + ) + .clamp(0.0, 1.0); } fn update_sphincters(&mut self) { - let internal = 0.3 + 0.7 * self.sympathetic_drive; - let external = 0.2 + 0.8 * self.somatic_drive; - self.internal_sphincter_tone = internal.clamp(0.0, 1.0); - self.external_sphincter_tone = external.clamp(0.0, 1.0); + let internal = (0.28 + 0.65 * self.sympathetic_drive + 0.1 * self.cortical_gate_fraction) + .clamp(0.0, 1.0); + let external = ((0.18 + 0.7 * self.somatic_drive + 0.25 * self.voluntary_hold_fraction) + * (0.75 + 0.25 * self.neurologic_integrity)) + .clamp(0.0, 1.0); + self.internal_sphincter_tone = internal; + self.external_sphincter_tone = external; } fn update_afferents(&mut self) { @@ -133,28 +489,69 @@ impl Bladder { } else { (self.volume_ml / self.capacity_ml.max(1.0)).clamp(0.0, 1.0) }; - self.afferent_signal = (low_volume_component + fullness_component).clamp(0.0, 1.0); - self.urgency = self.afferent_signal.powf(1.35).clamp(0.0, 1.0); + let pressure_bias = + (self.detrusor_pressure_cm_h2o / self.compliance.safety_limit()).clamp(0.0, 1.5) * 0.2; + let neurologic_gain = (1.0 + 0.35 * (1.0 - self.neurologic_integrity)).clamp(1.0, 1.4); + self.afferent_signal = ((low_volume_component + fullness_component + pressure_bias) + * neurologic_gain) + .clamp(0.0, 1.2); + let cortical_relief = (1.0 - 0.45 * self.cortical_gate_fraction).clamp(0.5, 1.0); + let voluntary_relief = + (1.0 - 0.35 * self.voluntary_hold_fraction - 0.1 * self.continence_training) + .clamp(0.45, 1.0); + self.urgency = + (self.afferent_signal.powf(1.35) * cortical_relief * voluntary_relief).clamp(0.0, 1.0); } - fn update_pressure(&mut self) { - let passive_volume_ml = (self.volume_ml - 30.0).max(0.0); - let normalized_volume = (self.volume_ml / self.capacity_ml).clamp(0.0, 1.6); - let compliance_scale = 1.0 + 2.4 * normalized_volume; - let passive_base = if self.compliance_ml_per_cm_h2o > 0.0 { - passive_volume_ml / (self.compliance_ml_per_cm_h2o * compliance_scale) + fn update_pressure_safety_metrics(&mut self, detrusor_pressure: f32, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + let limit = self.compliance.safety_limit(); + let overload = (detrusor_pressure - limit).max(0.0); + let decay = (-dt_seconds / 600.0).exp(); + if overload > 0.0 { + self.time_above_safety_limit_s += dt_seconds; + self.pressure_stress_index = + (self.pressure_stress_index * decay + overload * dt_seconds).min(5000.0); + } else { + self.pressure_stress_index = (self.pressure_stress_index * decay).max(0.0); + } + } + + fn update_pressure(&mut self, dt_seconds: f32) { + let normalized_volume = (self.volume_ml / self.capacity_ml).clamp(0.0, 2.0); + let effective_compliance = self + .compliance + .effective_compliance(normalized_volume, self.filling_rate_ml_per_min); + + let passive_volume_ml = (self.volume_ml - self.residual_volume_ml).max(0.0); + let passive_detrusor = if effective_compliance > f32::EPSILON { + passive_volume_ml / effective_compliance } else { 0.0 }; - let nonlinear_gain = 1.0 + normalized_volume.powf(2.0); - let passive_pressure = passive_base * nonlinear_gain; - let active_pressure = 34.0 * self.parasympathetic_drive; - let abdominal = self.baseline_pressure_cm_h2o + 2.0 * normalized_volume; - self.pressure = (abdominal + passive_pressure + active_pressure).clamp(0.0, 80.0); + let toe_region_gain = 1.0 + normalized_volume.powf(3.2) * 0.6; + let passive_detrusor = (passive_detrusor * toe_region_gain).clamp(0.0, 60.0); + + let active_detrusor = (36.0 + 6.0 * normalized_volume) * self.parasympathetic_drive; + let detrusor_pressure = (passive_detrusor + active_detrusor).clamp(0.0, 95.0); + + self.compliance.update_state(detrusor_pressure, dt_seconds); + + let abdominal_pressure = self.baseline_pressure_cm_h2o + + 2.0 * normalized_volume + + 0.5 * self.internal_sphincter_tone; + + self.detrusor_pressure_cm_h2o = detrusor_pressure; + self.pressure = (abdominal_pressure + detrusor_pressure).clamp(0.0, 110.0); + + self.update_pressure_safety_metrics(detrusor_pressure, dt_seconds); } fn handle_filling_phase(&mut self, _dt_seconds: f32) { - if (self.volume_ml >= self.micturition_threshold_ml || self.pressure > 45.0) + if (self.volume_ml >= self.micturition_threshold_ml + || self.detrusor_pressure_cm_h2o > self.compliance.safety_limit() + 5.0) && (self.external_sphincter_tone < 0.4 || self.urgency > 0.95) { self.phase = BladderPhase::Voiding; @@ -168,7 +565,7 @@ impl Bladder { fn handle_voiding_phase(&mut self, dt_seconds: f32) { let relaxation_factor = 1.0 - 0.5 * self.external_sphincter_tone; let drive = (self.parasympathetic_drive * relaxation_factor).clamp(0.0, 1.0); - let pressure_factor = (self.pressure / 40.0).clamp(0.0, 2.5); + let pressure_factor = (self.detrusor_pressure_cm_h2o / 40.0).clamp(0.0, 2.5); let outflow = self.voiding_flow_ml_per_s.max(0.0) * pressure_factor * drive * dt_seconds; self.volume_ml = (self.volume_ml - outflow).max(self.residual_volume_ml); if self.volume_ml <= self.residual_volume_ml + 1.0 { @@ -182,6 +579,32 @@ impl Bladder { self.phase = BladderPhase::Filling; } } + pub fn metrics(&self) -> BladderMetrics { + BladderMetrics { + phase: self.phase, + volume_ml: self.volume_ml, + capacity_ml: self.capacity_ml, + pressure_cm_h2o: self.pressure, + detrusor_pressure_cm_h2o: self.detrusor_pressure_cm_h2o, + pressure_safety_limit_cm_h2o: self.compliance.safety_limit(), + afferent_signal: self.afferent_signal, + urgency: self.urgency, + cortical_inhibition: self.cortical_inhibition, + voluntary_hold_command: self.voluntary_hold_command, + cortical_gate_fraction: self.cortical_gate_fraction, + voluntary_hold_fraction: self.voluntary_hold_fraction, + guarding_reflex_gain: self.reflex.guarding_reflex_gain, + pontine_storage_signal: self.reflex.pontine_storage_signal, + pontine_micturition_signal: self.reflex.pontine_micturition_signal, + hypogastric_efferent: self.reflex.hypogastric_efferent, + pudendal_efferent: self.reflex.pudendal_efferent, + parasympathetic_drive: self.parasympathetic_drive, + sympathetic_drive: self.sympathetic_drive, + somatic_drive: self.somatic_drive, + pressure_stress_index: self.pressure_stress_index, + time_above_safety_limit_s: self.time_above_safety_limit_s, + } + } } impl Organ for Bladder { @@ -200,10 +623,14 @@ impl Organ for Bladder { let inflow = (self.filling_rate_ml_per_min / 60.0).max(0.0) * dt_seconds; self.volume_ml += inflow; + self.update_voluntary_modulators(dt_seconds); + self.update_threshold_targets(dt_seconds); self.update_afferents(); + self.update_reflex_gate(dt_seconds); self.update_drives(dt_seconds); self.update_sphincters(); - self.update_pressure(); + self.update_pressure(dt_seconds); + self.update_reflex_efferents(dt_seconds); match self.phase { BladderPhase::Filling => self.handle_filling_phase(dt_seconds), @@ -220,13 +647,19 @@ impl Organ for Bladder { } fn summary(&self) -> String { format!( - "Bladder[id={}, phase={:?}, vol={:.0}/{:.0} ml, P={:.1} cmH2O, urge={:.0}%]", + "Bladder[id={}, phase={:?}, vol={:.0}/{:.0} ml, P={:.1} cmH2O (det={:.1}), urge={:.0}%, stress={:.0}, guard={:.0}%, void={:.0}%, gate={:.0}%, hold={:.0}%]", self.id(), self.phase, self.volume_ml, self.capacity_ml, self.pressure, - self.urgency * 100.0 + self.detrusor_pressure_cm_h2o, + self.urgency * 100.0, + self.pressure_stress_index, + self.reflex.guarding_reflex_gain * 100.0, + self.reflex.pontine_micturition_signal * 100.0, + self.cortical_gate_fraction * 100.0, + self.voluntary_hold_fraction * 100.0 ) } fn as_any(&self) -> &dyn core::any::Any { @@ -236,3 +669,181 @@ impl Organ for Bladder { self } } + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn compliance_declines_with_high_volume() { + let model = ComplianceModel::new(); + let compliance_low = model.effective_compliance(0.3, 60.0); + let compliance_high = model.effective_compliance(1.3, 60.0); + assert!(compliance_low > compliance_high); + } + + #[test] + fn detrusor_overload_updates_risk_metrics() { + let mut bladder = Bladder::new("test"); + bladder.volume_ml = bladder.capacity_ml * 1.35; + bladder.parasympathetic_drive = 0.7; + bladder.sympathetic_drive = 0.2; + bladder.internal_sphincter_tone = 0.2; + bladder.external_sphincter_tone = 0.2; + + bladder.update_pressure(1.0); + + assert!( + bladder.detrusor_pressure_cm_h2o > bladder.compliance.safety_limit(), + "expected detrusor {} to exceed limit {}", + bladder.detrusor_pressure_cm_h2o, + bladder.compliance.safety_limit() + ); + assert!(bladder.time_above_safety_limit_s > 0.0); + assert!(bladder.pressure_stress_index > 0.0); + } + + #[test] + fn stress_index_decays_when_pressure_normalizes() { + let mut bladder = Bladder::new("test"); + bladder.pressure_stress_index = 200.0; + bladder.time_above_safety_limit_s = 12.0; + bladder.volume_ml = 140.0; + bladder.parasympathetic_drive = 0.05; + bladder.sympathetic_drive = 0.9; + bladder.internal_sphincter_tone = 0.9; + bladder.external_sphincter_tone = 0.9; + + bladder.update_pressure(1.5); + let after = bladder.pressure_stress_index; + bladder.update_pressure(1.5); + assert!( + after < 200.0 && bladder.pressure_stress_index <= after, + "expected stress index to decay, got {} then {}", + after, + bladder.pressure_stress_index + ); + } + + #[test] + fn guarding_reflex_active_during_storage() { + let mut bladder = Bladder::new("storage"); + bladder.volume_ml = bladder.urge_threshold_ml + 60.0; + bladder.update(1.0); + assert!(bladder.reflex.guarding_reflex_gain > 0.2); + assert!(bladder.reflex.hypogastric_efferent > 0.5); + assert!(bladder.reflex.pudendal_efferent > 0.4); + assert!( + bladder.reflex.pontine_storage_signal > bladder.reflex.pontine_micturition_signal, + "storage center should dominate in filling" + ); + } + + #[test] + fn pontine_switch_reduces_guarding_during_void() { + let mut bladder = Bladder::new("void"); + bladder.phase = BladderPhase::Voiding; + bladder.volume_ml = bladder.micturition_threshold_ml + 180.0; + bladder.detrusor_pressure_cm_h2o = bladder.compliance.safety_limit() + 12.0; + bladder.urgency = 0.98; + + bladder.update(0.3); + bladder.update(0.2); + + let pontine = bladder.reflex.pontine_micturition_signal; + let hypogastric = bladder.reflex.hypogastric_efferent; + let pudendal = bladder.reflex.pudendal_efferent; + + assert!( + pontine > 0.5, + "pontine signal too low: {:.2} (guard={:.2}, aff={:.2})", + pontine, + bladder.reflex.guarding_reflex_gain, + bladder.reflex.pelvic_afferent_rate, + ); + assert!( + hypogastric < 0.45, + "hypogastric stays high: {:.2}", + hypogastric + ); + assert!(pudendal < 0.45, "pudendal stays high: {:.2}", pudendal); + } + + #[test] + fn cortical_behavioral_inputs_reduce_urgency() { + let mut baseline = Bladder::new("baseline"); + baseline.volume_ml = baseline.micturition_threshold_ml + 120.0; + baseline.cortical_inhibition = 0.0; + baseline.voluntary_hold_command = 0.0; + baseline.continence_training = 0.3; + for _ in 0..5 { + baseline.update(0.5); + } + + let mut conditioned = Bladder::new("conditioned"); + conditioned.volume_ml = baseline.volume_ml; + conditioned.cortical_inhibition = 0.85; + conditioned.voluntary_hold_command = 0.8; + conditioned.continence_training = 0.85; + for _ in 0..5 { + conditioned.update(0.5); + } + + assert!( + conditioned.urgency < baseline.urgency, + "expected urgency {:.2} to fall below baseline {:.2}", + conditioned.urgency, + baseline.urgency + ); + assert!( + conditioned.external_sphincter_tone > baseline.external_sphincter_tone, + "expected hold to raise sphincter tone {:.2} vs {:.2}", + conditioned.external_sphincter_tone, + baseline.external_sphincter_tone + ); + } + + #[test] + fn thresholds_reflect_age_and_neurologic_status() { + let mut young = Bladder::new("young"); + young.age_years = 25.0; + young.neurologic_integrity = 1.0; + young.continence_training = 0.5; + young.cortical_inhibition = 0.0; + young.voluntary_hold_command = 0.0; + young.filling_rate_ml_per_min = 45.0; + for _ in 0..12 { + young.update(1.0); + } + + let mut elderly = Bladder::new("elderly"); + elderly.age_years = 78.0; + elderly.neurologic_integrity = 0.45; + elderly.continence_training = 0.5; + elderly.cortical_inhibition = 0.0; + elderly.voluntary_hold_command = 0.0; + elderly.filling_rate_ml_per_min = 95.0; + for _ in 0..12 { + elderly.update(1.0); + } + + assert!( + elderly.capacity_ml < young.capacity_ml, + "expected elderly capacity {:.1} < young {:.1}", + elderly.capacity_ml, + young.capacity_ml + ); + assert!( + elderly.urge_threshold_ml < young.urge_threshold_ml, + "expected elderly urge {:.1} < young {:.1}", + elderly.urge_threshold_ml, + young.urge_threshold_ml + ); + assert!( + elderly.micturition_threshold_ml < young.micturition_threshold_ml, + "expected elderly mict {:.1} < young {:.1}", + elderly.micturition_threshold_ml, + young.micturition_threshold_ml + ); + } +} diff --git a/src/organs/bloodstream.rs b/src/organs/bloodstream.rs index ffcb8db..304764d 100644 --- a/src/organs/bloodstream.rs +++ b/src/organs/bloodstream.rs @@ -26,6 +26,37 @@ const BASE_EPO_IU_PER_DAY: f32 = 18.0; const BASE_PLATELET_RESERVOIR: f32 = 70.0; const BASE_RED_PULP_RESERVE_ML: f32 = 180.0; +const BASE_PLASMA_ALBUMIN_G_DL: f32 = 4.0; +const BASE_PLASMA_GLOBULIN_G_DL: f32 = 2.5; +const BASE_PLASMA_FIBRINOGEN_G_DL: f32 = 0.3; +const BASE_PLASMA_ONCOTIC_PRESSURE_MMHG: f32 = 24.5; +const BASE_LYMPH_RETURN_ML_MIN: f32 = 2.6; +const RBC_LIFESPAN_DAYS: f32 = 115.0; +const RBC_YOUNG_WINDOW_DAYS: f32 = 28.0; +const RBC_MATURE_WINDOW_DAYS: f32 = 72.0; +const RBC_SENESCENT_WINDOW_DAYS: f32 = 15.0; +const IRON_PER_ML_RBC_MG: f32 = 1.05; +const FOLATE_PER_ML_RBC_MG: f32 = 0.00045; +const B12_PER_ML_RBC_MCG: f32 = 0.002; +const BASE_IRON_STORE_MG: f32 = 820.0; +const BASE_FOLATE_STORE_MG: f32 = 8.0; +const BASE_B12_STORE_MCG: f32 = 3000.0; +const BASE_SPLENIC_CULLING_RATE_1E6_PER_MIN: f32 = 2.5; + +const BASE_PLATELET_COUNT_1E3_PER_UL: f32 = 250.0; +const PLATELET_LIFESPAN_DAYS: f32 = 9.0; +const BASE_PLATELET_ACTIVATION: f32 = 0.18; +const BASE_COAGULATION_FACTOR_ACTIVITY: f32 = 1.0; +const BASE_FIBRINOLYSIS_ACTIVITY: f32 = 1.0; +const BASE_WBC_COUNT_GIGA_PER_L: f32 = 6.5; +const BASE_NEUTROPHIL_FRACTION: f32 = 0.6; +const BASE_LYMPHOCYTE_FRACTION: f32 = 0.3; +const BASE_MONOCYTE_FRACTION: f32 = 0.08; +const BASE_COMPLEMENT_ACTIVITY: f32 = 1.0; +const BASE_BICARBONATE_MMOL_L: f32 = 24.0; +const BASE_BASE_EXCESS_MMOL_L: f32 = 0.0; +const BASE_ANION_GAP_MMOL_L: f32 = 12.0; + #[derive(Debug, Clone, Copy, PartialEq, Eq)] /// Aggregate perfusion assessment for the bloodstream. pub enum PerfusionState { @@ -50,6 +81,47 @@ pub enum MetabolicState { AnaerobicCrisis, } +#[derive(Debug, Clone, Copy)] +/// Snapshot of bloodstream composition and regulatory signals. +pub struct BloodstreamMetrics { + pub perfusion_state: PerfusionState, + pub metabolic_state: MetabolicState, + pub total_volume_l: f32, + pub plasma_volume_l: f32, + pub red_cell_volume_l: f32, + pub plasma_albumin_g_dl: f32, + pub plasma_globulin_g_dl: f32, + pub plasma_fibrinogen_g_dl: f32, + pub oncotic_pressure_mm_hg: f32, + pub lymphatic_return_ml_min: f32, + pub rbc_young_fraction: f32, + pub rbc_mature_fraction: f32, + pub rbc_senescent_fraction: f32, + pub platelet_count_1e3_per_ul: f32, + pub platelet_activation: f32, + pub coagulation_factor_activity: f32, + pub fibrinolysis_activity: f32, + pub thrombosis_risk_index: f32, + pub bleeding_risk_index: f32, + pub wbc_count_giga_per_l: f32, + pub neutrophil_fraction: f32, + pub lymphocyte_fraction: f32, + pub monocyte_fraction: f32, + pub complement_activity: f32, + pub inflammation_index: f32, + pub bicarbonate_mmol_l: f32, + pub base_excess_mmol_l: f32, + pub anion_gap_mmol_l: f32, + pub arterial_pco2_mm_hg: f32, + pub erythropoiesis_rate_ml_per_day: f32, + pub erythrocyte_clearance_ml_per_day: f32, + pub iron_store_mg: f32, + pub folate_store_mg: f32, + pub b12_store_mcg: f32, + pub hif_activation: f32, + pub oxygen_supply_demand_ratio: f32, +} + #[derive(Debug, Clone)] /// Simplified blood transport compartment tracking gases, volume, and waste. pub struct Bloodstream { @@ -108,6 +180,66 @@ pub struct Bloodstream { pub perfusion_state: PerfusionState, /// Current metabolic balance. pub metabolic_state: MetabolicState, + /// Plasma albumin concentration (g/dL). + pub plasma_albumin_g_dl: f32, + /// Plasma globulin concentration (g/dL). + pub plasma_globulin_g_dl: f32, + /// Plasma fibrinogen concentration (g/dL). + pub plasma_fibrinogen_g_dl: f32, + /// Colloid oncotic pressure (mmHg). + pub oncotic_pressure_mm_hg: f32, + /// Lymphatic return flow (mL/min). + pub lymphatic_return_ml_min: f32, + /// Circulating platelet count (10^3/µL). + pub platelet_count_1e3_per_ul: f32, + /// Fraction of platelets currently activated (0..=1). + pub platelet_activation: f32, + /// Relative coagulation factor activity. + pub coagulation_factor_activity: f32, + /// Relative fibrinolysis activity. + pub fibrinolysis_activity: f32, + /// Composite thrombosis risk index. + pub thrombosis_risk_index: f32, + /// Composite bleeding risk index. + pub bleeding_risk_index: f32, + /// Total leukocyte count (10^9/L). + pub wbc_count_giga_per_l: f32, + /// Neutrophil fraction of leukocytes (0..=1). + pub neutrophil_fraction: f32, + /// Lymphocyte fraction of leukocytes (0..=1). + pub lymphocyte_fraction: f32, + /// Monocyte fraction of leukocytes (0..=1). + pub monocyte_fraction: f32, + /// Complement activation index. + pub complement_activity: f32, + /// Composite inflammation index. + pub inflammation_index: f32, + /// Serum bicarbonate (mmol/L). + pub bicarbonate_mmol_l: f32, + /// Base excess (mmol/L). + pub base_excess_mmol_l: f32, + /// Estimated anion gap (mmol/L). + pub anion_gap_mmol_l: f32, + /// Arterial CO2 partial pressure (mmHg). + pub arterial_pco2_mm_hg: f32, + /// Red cell volume in the young cohort (L). + pub rbc_young_volume_l: f32, + /// Red cell volume in the mature cohort (L). + pub rbc_mature_volume_l: f32, + /// Red cell volume in the senescent cohort (L). + pub rbc_senescent_volume_l: f32, + /// Daily erythrocyte production (mL/day). + pub erythropoiesis_rate_ml_per_day: f32, + /// Daily erythrocyte clearance (mL/day). + pub erythrocyte_clearance_ml_per_day: f32, + /// Iron store available for erythropoiesis (mg). + pub iron_store_mg: f32, + /// Folate store available for erythropoiesis (mg). + pub folate_store_mg: f32, + /// Vitamin B12 store available for erythropoiesis (mcg). + pub b12_store_mcg: f32, + /// Hypoxia-inducible factor activation (0..=1). + pub hif_activation: f32, ventilation_perfusion_ratio: f32, // Signal-driven targets cardiac_output_target_l_min: f32, @@ -126,6 +258,14 @@ pub struct Bloodstream { epo_signal_iu_per_day: f32, spleen_platelet_signal: f32, spleen_red_pulp_reserve_ml: f32, + spleen_immune_activity_index: f32, + spleen_culling_rate_1e6_per_min: f32, + hepatic_albumin_target_g_dl: f32, + immune_globulin_target_g_dl: f32, + fibrinogen_target_g_dl: f32, + iron_intake_mg_per_day: f32, + folate_intake_mg_per_day: f32, + b12_intake_mcg_per_day: f32, } impl Bloodstream { @@ -176,6 +316,36 @@ impl Bloodstream { hepatic_clearance_index: BASE_HEPATIC_CLEARANCE_INDEX, perfusion_state: PerfusionState::Balanced, metabolic_state: MetabolicState::Aerobic, + plasma_albumin_g_dl: BASE_PLASMA_ALBUMIN_G_DL, + plasma_globulin_g_dl: BASE_PLASMA_GLOBULIN_G_DL, + plasma_fibrinogen_g_dl: BASE_PLASMA_FIBRINOGEN_G_DL, + oncotic_pressure_mm_hg: BASE_PLASMA_ONCOTIC_PRESSURE_MMHG, + lymphatic_return_ml_min: BASE_LYMPH_RETURN_ML_MIN, + platelet_count_1e3_per_ul: BASE_PLATELET_COUNT_1E3_PER_UL, + platelet_activation: BASE_PLATELET_ACTIVATION, + coagulation_factor_activity: BASE_COAGULATION_FACTOR_ACTIVITY, + fibrinolysis_activity: BASE_FIBRINOLYSIS_ACTIVITY, + thrombosis_risk_index: 0.0, + bleeding_risk_index: 0.2, + wbc_count_giga_per_l: BASE_WBC_COUNT_GIGA_PER_L, + neutrophil_fraction: BASE_NEUTROPHIL_FRACTION, + lymphocyte_fraction: BASE_LYMPHOCYTE_FRACTION, + monocyte_fraction: BASE_MONOCYTE_FRACTION, + complement_activity: BASE_COMPLEMENT_ACTIVITY, + inflammation_index: 0.0, + bicarbonate_mmol_l: BASE_BICARBONATE_MMOL_L, + base_excess_mmol_l: BASE_BASE_EXCESS_MMOL_L, + anion_gap_mmol_l: BASE_ANION_GAP_MMOL_L, + arterial_pco2_mm_hg: BASE_ALVEOLAR_PCO2_MMHG, + rbc_young_volume_l: BASE_RED_CELL_VOLUME_L * 0.22, + rbc_mature_volume_l: BASE_RED_CELL_VOLUME_L * 0.63, + rbc_senescent_volume_l: BASE_RED_CELL_VOLUME_L * 0.15, + erythropoiesis_rate_ml_per_day: (BASE_RED_CELL_VOLUME_L * 1000.0) / RBC_LIFESPAN_DAYS, + erythrocyte_clearance_ml_per_day: (BASE_RED_CELL_VOLUME_L * 1000.0) / RBC_LIFESPAN_DAYS, + iron_store_mg: BASE_IRON_STORE_MG, + folate_store_mg: BASE_FOLATE_STORE_MG, + b12_store_mcg: BASE_B12_STORE_MCG, + hif_activation: 0.0, ventilation_perfusion_ratio: BASE_VQ_RATIO, cardiac_output_target_l_min: BASE_CARDIAC_OUTPUT_L_MIN, map_target_mm_hg: BASE_MAP_MM_HG, @@ -193,6 +363,14 @@ impl Bloodstream { epo_signal_iu_per_day: BASE_EPO_IU_PER_DAY, spleen_platelet_signal: BASE_PLATELET_RESERVOIR, spleen_red_pulp_reserve_ml: BASE_RED_PULP_RESERVE_ML, + spleen_immune_activity_index: 80.0, + spleen_culling_rate_1e6_per_min: BASE_SPLENIC_CULLING_RATE_1E6_PER_MIN, + hepatic_albumin_target_g_dl: BASE_PLASMA_ALBUMIN_G_DL, + immune_globulin_target_g_dl: BASE_PLASMA_GLOBULIN_G_DL, + fibrinogen_target_g_dl: BASE_PLASMA_FIBRINOGEN_G_DL, + iron_intake_mg_per_day: 1.6, + folate_intake_mg_per_day: 0.45, + b12_intake_mcg_per_day: 5.2, } } @@ -238,6 +416,18 @@ impl Bloodstream { self.hepatic_clearance_target_index = (detox * 0.7 + ammonia_factor * 0.3).clamp(0.3, 1.3); } + /// Store hepatic protein synthesis and acute phase signals. + pub fn set_plasma_proteins( + &mut self, + albumin_g_dl: f32, + globulin_g_dl: f32, + fibrinogen_g_dl: f32, + ) { + self.hepatic_albumin_target_g_dl = albumin_g_dl.clamp(2.0, 5.5); + self.immune_globulin_target_g_dl = globulin_g_dl.clamp(1.0, 5.8); + self.fibrinogen_target_g_dl = fibrinogen_g_dl.clamp(0.05, 1.0); + } + /// Store nutrient availability to adjust metabolic demand. pub fn set_metabolic_nutrients(&mut self, nutrient_energy_kcal: f32) { let kcal = nutrient_energy_kcal.clamp(0.0, 1200.0); @@ -247,15 +437,80 @@ impl Bloodstream { (BASE_WASTE_PRODUCTION_MG_MIN + kcal * 1.0).clamp(500.0, 2200.0); } + /// Store micronutrient availability for erythropoiesis. + pub fn set_erythropoiesis_nutrients( + &mut self, + iron_absorption_mg_per_day: f32, + folate_absorption_mg_per_day: f32, + b12_absorption_mcg_per_day: f32, + ) { + self.iron_intake_mg_per_day = iron_absorption_mg_per_day.clamp(0.0, 12.0); + self.folate_intake_mg_per_day = folate_absorption_mg_per_day.clamp(0.0, 3.0); + self.b12_intake_mcg_per_day = b12_absorption_mcg_per_day.clamp(0.0, 40.0); + } + /// Store erythropoietin drive from the kidneys. pub fn ingest_erythropoietin(&mut self, epo_iu_per_day: f32) { self.epo_signal_iu_per_day = epo_iu_per_day.clamp(4.0, 120.0); } - /// Store splenic reservoir state for rapid hematocrit adjustments. - pub fn set_spleen_feedback(&mut self, platelet_reservoir: f32, red_pulp_volume_ml: f32) { + /// Store splenic reservoir state for rapid hematocrit adjustments and clearance. + pub fn set_spleen_feedback( + &mut self, + platelet_reservoir: f32, + red_pulp_volume_ml: f32, + immune_activity: u8, + erythrocyte_culling_rate: f32, + ) { self.spleen_platelet_signal = platelet_reservoir.clamp(10.0, 200.0); self.spleen_red_pulp_reserve_ml = red_pulp_volume_ml.clamp(40.0, 400.0); + self.spleen_immune_activity_index = immune_activity as f32; + self.spleen_culling_rate_1e6_per_min = erythrocyte_culling_rate.clamp(0.5, 12.0); + } + + /// Return a physiology snapshot for observers and FFI layers. + pub fn metrics(&self) -> BloodstreamMetrics { + let total_rbc = + (self.rbc_young_volume_l + self.rbc_mature_volume_l + self.rbc_senescent_volume_l) + .max(1e-6); + BloodstreamMetrics { + perfusion_state: self.perfusion_state, + metabolic_state: self.metabolic_state, + total_volume_l: self.total_volume_l, + plasma_volume_l: self.plasma_volume_l, + red_cell_volume_l: self.red_cell_volume_l, + plasma_albumin_g_dl: self.plasma_albumin_g_dl, + plasma_globulin_g_dl: self.plasma_globulin_g_dl, + plasma_fibrinogen_g_dl: self.plasma_fibrinogen_g_dl, + oncotic_pressure_mm_hg: self.oncotic_pressure_mm_hg, + lymphatic_return_ml_min: self.lymphatic_return_ml_min, + rbc_young_fraction: (self.rbc_young_volume_l / total_rbc).clamp(0.0, 1.0), + rbc_mature_fraction: (self.rbc_mature_volume_l / total_rbc).clamp(0.0, 1.0), + rbc_senescent_fraction: (self.rbc_senescent_volume_l / total_rbc).clamp(0.0, 1.0), + platelet_count_1e3_per_ul: self.platelet_count_1e3_per_ul, + platelet_activation: self.platelet_activation.clamp(0.0, 1.0), + coagulation_factor_activity: self.coagulation_factor_activity, + fibrinolysis_activity: self.fibrinolysis_activity, + thrombosis_risk_index: self.thrombosis_risk_index, + bleeding_risk_index: self.bleeding_risk_index, + wbc_count_giga_per_l: self.wbc_count_giga_per_l, + neutrophil_fraction: self.neutrophil_fraction.clamp(0.0, 1.0), + lymphocyte_fraction: self.lymphocyte_fraction.clamp(0.0, 1.0), + monocyte_fraction: self.monocyte_fraction.clamp(0.0, 1.0), + complement_activity: self.complement_activity, + inflammation_index: self.inflammation_index, + bicarbonate_mmol_l: self.bicarbonate_mmol_l, + base_excess_mmol_l: self.base_excess_mmol_l, + anion_gap_mmol_l: self.anion_gap_mmol_l, + arterial_pco2_mm_hg: self.arterial_pco2_mm_hg, + erythropoiesis_rate_ml_per_day: self.erythropoiesis_rate_ml_per_day, + erythrocyte_clearance_ml_per_day: self.erythrocyte_clearance_ml_per_day, + iron_store_mg: self.iron_store_mg, + folate_store_mg: self.folate_store_mg, + b12_store_mcg: self.b12_store_mcg, + hif_activation: self.hif_activation.clamp(0.0, 1.2), + oxygen_supply_demand_ratio: self.oxygen_supply_demand_ratio, + } } /// Override measured glucose value. @@ -278,30 +533,375 @@ impl Bloodstream { } } + fn update_plasma_proteins(&mut self, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + self.plasma_albumin_g_dl = Self::approach( + self.plasma_albumin_g_dl, + self.hepatic_albumin_target_g_dl, + 0.0025, + dt_seconds, + ); + let globulin_target = (self.immune_globulin_target_g_dl + + (self.spleen_immune_activity_index - 80.0) * 0.01) + .clamp(1.0, 5.8); + self.plasma_globulin_g_dl = Self::approach( + self.plasma_globulin_g_dl, + globulin_target, + 0.003, + dt_seconds, + ); + self.plasma_fibrinogen_g_dl = Self::approach( + self.plasma_fibrinogen_g_dl, + self.fibrinogen_target_g_dl, + 0.0015, + dt_seconds, + ); + let total_protein = + (self.plasma_albumin_g_dl + self.plasma_globulin_g_dl + self.plasma_fibrinogen_g_dl) + .clamp(2.0, 9.0); + let oncotic_estimate = (2.1 * total_protein) + + (0.16 * total_protein * total_protein) + + (0.009 * total_protein * total_protein * total_protein); + self.oncotic_pressure_mm_hg = Self::approach( + self.oncotic_pressure_mm_hg, + oncotic_estimate.clamp(16.0, 38.0), + 0.05, + dt_seconds, + ); + let capillary_pressure = (self.mean_arterial_pressure_mm_hg * 0.27).clamp(15.0, 38.0); + let venous_refill = (self.total_volume_l - BASE_TOTAL_VOLUME_L) * 8.0; + let filtration_pressure = + (capillary_pressure - self.oncotic_pressure_mm_hg - venous_refill).clamp(-8.0, 18.0); + let lymph_target = (BASE_LYMPH_RETURN_ML_MIN + filtration_pressure * 0.55).clamp(0.6, 12.0); + self.lymphatic_return_ml_min = + Self::approach(self.lymphatic_return_ml_min, lymph_target, 0.15, dt_seconds); + } + fn update_red_cell_mass(&mut self, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + + let seconds_per_day = 86400.0; + let iron_intake = self.iron_intake_mg_per_day * dt_seconds / seconds_per_day; + let folate_intake = self.folate_intake_mg_per_day * dt_seconds / seconds_per_day; + let b12_intake = self.b12_intake_mcg_per_day * dt_seconds / seconds_per_day; + self.iron_store_mg = (self.iron_store_mg + iron_intake).clamp(120.0, 2400.0); + self.folate_store_mg = (self.folate_store_mg + folate_intake).clamp(1.0, 24.0); + self.b12_store_mcg = (self.b12_store_mcg + b12_intake).clamp(400.0, 8000.0); + let epo_adjust = (self.epo_signal_iu_per_day - BASE_EPO_IU_PER_DAY) * 0.002; let spleen_adjust = (self.spleen_red_pulp_reserve_ml - BASE_RED_PULP_RESERVE_ML) * 0.0004 - (self.spleen_platelet_signal - BASE_PLATELET_RESERVOIR) * 0.0003; - let demand_adjust = (1.0 - self.oxygen_supply_demand_ratio).max(0.0) * 0.12; - let desired = - (BASE_RED_CELL_VOLUME_L + epo_adjust + spleen_adjust + demand_adjust).clamp(1.6, 2.8); + let demand_adjust = (1.0 - self.oxygen_supply_demand_ratio).max(0.0) * 0.18; + let nutrient_pressure = ((self.iron_store_mg / BASE_IRON_STORE_MG) - 1.0) * 0.06 + + ((self.folate_store_mg / BASE_FOLATE_STORE_MG) - 1.0) * 0.04 + + ((self.b12_store_mcg / BASE_B12_STORE_MCG) - 1.0) * 0.03; + let desired = (BASE_RED_CELL_VOLUME_L + + epo_adjust + + spleen_adjust + + demand_adjust + + nutrient_pressure) + .clamp(1.6, 3.0); self.red_cell_mass_target_l = - Self::approach(self.red_cell_mass_target_l, desired, 0.0006, dt_seconds); - self.red_cell_volume_l = Self::approach( - self.red_cell_volume_l, - self.red_cell_mass_target_l, - 0.0012, + Self::approach(self.red_cell_mass_target_l, desired, 0.0004, dt_seconds); + + let saturation_deficit = + ((BASE_ARTERIAL_SPO2_PCT - self.arterial_o2_saturation_pct).max(0.0)) / 30.0; + let supply_deficit = (1.05 - self.oxygen_supply_demand_ratio).max(0.0); + let hif_target = (saturation_deficit * 0.6 + supply_deficit * 0.5).clamp(0.0, 1.2); + self.hif_activation = Self::approach(self.hif_activation, hif_target, 0.004, dt_seconds); + + let young_tau = RBC_YOUNG_WINDOW_DAYS * seconds_per_day; + let mature_tau = RBC_MATURE_WINDOW_DAYS * seconds_per_day; + let senescent_tau = RBC_SENESCENT_WINDOW_DAYS * seconds_per_day; + + let young_progress = + (self.rbc_young_volume_l * dt_seconds / young_tau).min(self.rbc_young_volume_l); + self.rbc_young_volume_l -= young_progress; + self.rbc_mature_volume_l += young_progress; + + let mature_progress = + (self.rbc_mature_volume_l * dt_seconds / mature_tau).min(self.rbc_mature_volume_l); + self.rbc_mature_volume_l -= mature_progress; + self.rbc_senescent_volume_l += mature_progress; + + let mut red_volume = + self.rbc_young_volume_l + self.rbc_mature_volume_l + self.rbc_senescent_volume_l; + let base_clearance_per_s = + (red_volume.max(0.5) / (RBC_LIFESPAN_DAYS * seconds_per_day)).clamp(0.0, 0.0025); + let spleen_factor = (self.spleen_culling_rate_1e6_per_min + / BASE_SPLENIC_CULLING_RATE_1E6_PER_MIN) + .clamp(0.4, 2.6); + let immune_factor = + (1.0 + (self.spleen_immune_activity_index - 80.0) * 0.004).clamp(0.6, 1.8); + + let senescent_turnover = (self.rbc_senescent_volume_l * dt_seconds / senescent_tau) + .min(self.rbc_senescent_volume_l); + let clearance_l = (base_clearance_per_s * dt_seconds * spleen_factor * immune_factor + + senescent_turnover) + .min(self.rbc_senescent_volume_l); + self.rbc_senescent_volume_l -= clearance_l.max(0.0); + + let clearance_rate_ml_per_day = if dt_seconds > 0.0 { + (clearance_l / dt_seconds) * 1000.0 * seconds_per_day + } else { + self.erythrocyte_clearance_ml_per_day + }; + self.erythrocyte_clearance_ml_per_day = Self::approach( + self.erythrocyte_clearance_ml_per_day, + clearance_rate_ml_per_day.clamp(4.0, 140.0), + 0.2, + dt_seconds, + ); + + let recycled_iron = clearance_l * 1000.0 * IRON_PER_ML_RBC_MG; + self.iron_store_mg = (self.iron_store_mg + recycled_iron).clamp(120.0, 2400.0); + + red_volume = + self.rbc_young_volume_l + self.rbc_mature_volume_l + self.rbc_senescent_volume_l; + let baseline_prod_l_per_s = BASE_RED_CELL_VOLUME_L / (RBC_LIFESPAN_DAYS * seconds_per_day); + let epo_factor = (self.epo_signal_iu_per_day / BASE_EPO_IU_PER_DAY).clamp(0.2, 4.0); + let hif_factor = (1.0 + self.hif_activation * 1.6).clamp(0.4, 2.6); + let deficit = (self.red_cell_mass_target_l - red_volume).max(0.0); + let deficit_factor = (1.0 + deficit * 1.8).clamp(1.0, 3.5); + let production_drive_l_per_s = + baseline_prod_l_per_s * epo_factor * hif_factor * deficit_factor; + let production_drive_ml_per_day = production_drive_l_per_s * 1000.0 * seconds_per_day; + + let iron_needed = production_drive_ml_per_day * IRON_PER_ML_RBC_MG; + let folate_needed = production_drive_ml_per_day * FOLATE_PER_ML_RBC_MG; + let b12_needed = production_drive_ml_per_day * B12_PER_ML_RBC_MCG; + + let mut nutrient_limit: f32 = 1.0; + if iron_needed > 1.0 { + nutrient_limit = nutrient_limit + .min((self.iron_store_mg / (iron_needed * 0.5 + 1.0)).clamp(0.0, 1.0)); + } + if folate_needed > 0.01 { + nutrient_limit = nutrient_limit + .min((self.folate_store_mg / (folate_needed * 0.5 + 0.02)).clamp(0.0, 1.0)); + } + if b12_needed > 0.01 { + nutrient_limit = + nutrient_limit.min((self.b12_store_mcg / (b12_needed * 0.5 + 0.2)).clamp(0.0, 1.0)); + } + + let production_ml_per_day = + (production_drive_ml_per_day * nutrient_limit).clamp(0.0, 220.0); + let production_l = production_ml_per_day / 1000.0 * dt_seconds / seconds_per_day; + let max_needed = (deficit + baseline_prod_l_per_s * dt_seconds * 2.5).clamp(0.0, 2.0); + let production_l = production_l.min(max_needed); + self.rbc_young_volume_l += production_l; + + let iron_used = production_l * 1000.0 * IRON_PER_ML_RBC_MG; + let folate_used = production_l * 1000.0 * FOLATE_PER_ML_RBC_MG; + let b12_used = production_l * 1000.0 * B12_PER_ML_RBC_MCG; + self.iron_store_mg = (self.iron_store_mg - iron_used).clamp(60.0, 2400.0); + self.folate_store_mg = (self.folate_store_mg - folate_used).clamp(0.5, 24.0); + self.b12_store_mcg = (self.b12_store_mcg - b12_used).clamp(200.0, 8000.0); + + self.erythropoiesis_rate_ml_per_day = Self::approach( + self.erythropoiesis_rate_ml_per_day, + production_ml_per_day, + 0.2, + dt_seconds, + ); + + self.rbc_young_volume_l = self.rbc_young_volume_l.clamp(0.0, 3.5); + self.rbc_mature_volume_l = self.rbc_mature_volume_l.clamp(0.0, 3.5); + self.rbc_senescent_volume_l = self.rbc_senescent_volume_l.clamp(0.0, 3.5); + self.red_cell_volume_l = + (self.rbc_young_volume_l + self.rbc_mature_volume_l + self.rbc_senescent_volume_l) + .clamp(1.5, 3.2); + } + + fn update_platelets(&mut self, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + let immune_factor = (self.spleen_immune_activity_index - 80.0) / 80.0; + let shear_factor = (self.mean_arterial_pressure_mm_hg - BASE_MAP_MM_HG) / 80.0; + let fibrinogen_factor = (self.plasma_fibrinogen_g_dl / BASE_PLASMA_FIBRINOGEN_G_DL) - 1.0; + let platelet_target = (BASE_PLATELET_COUNT_1E3_PER_UL + + (self.spleen_platelet_signal - BASE_PLATELET_RESERVOIR) * 1.1 + + immune_factor * 80.0 + + fibrinogen_factor * 140.0) + .clamp(90.0, 650.0); + self.platelet_count_1e3_per_ul = Self::approach( + self.platelet_count_1e3_per_ul, + platelet_target, + 5.0, + dt_seconds, + ); + let activation_target = (BASE_PLATELET_ACTIVATION + + immune_factor * 0.25 + + shear_factor * 0.18 + + (self.waste_load_mg - BASE_WASTE_LOAD_MG) * 0.00008) + .clamp(0.05, 0.95); + self.platelet_activation = Self::approach( + self.platelet_activation, + activation_target, + 0.02, + dt_seconds, + ); + let coag_target = (BASE_COAGULATION_FACTOR_ACTIVITY + + fibrinogen_factor * 0.45 + + (self.platelet_count_1e3_per_ul / BASE_PLATELET_COUNT_1E3_PER_UL - 1.0) * 0.4 + + self.platelet_activation * 0.5) + .clamp(0.4, 2.2); + self.coagulation_factor_activity = Self::approach( + self.coagulation_factor_activity, + coag_target, + 0.03, + dt_seconds, + ); + let fibrinolysis_target = (BASE_FIBRINOLYSIS_ACTIVITY + + immune_factor * 0.2 + + (self.lymphatic_return_ml_min - BASE_LYMPH_RETURN_ML_MIN) * 0.04 + - fibrinogen_factor * 0.3) + .clamp(0.3, 2.0); + self.fibrinolysis_activity = Self::approach( + self.fibrinolysis_activity, + fibrinolysis_target, + 0.025, + dt_seconds, + ); + let thrombosis_index = (self.coagulation_factor_activity - self.fibrinolysis_activity) + + (self.platelet_activation - BASE_PLATELET_ACTIVATION) * 1.4; + self.thrombosis_risk_index = thrombosis_index.clamp(-2.0, 3.0); + let bleeding_index = (1.0 - self.coagulation_factor_activity).max(0.0) + + (0.22 - self.platelet_activation).max(0.0) + + (200.0 - self.platelet_count_1e3_per_ul).max(0.0) / 220.0 + + (1.2 - self.plasma_fibrinogen_g_dl).max(0.0) * 0.5; + self.bleeding_risk_index = bleeding_index.clamp(0.0, 3.0); + } + + fn update_immune_cells(&mut self, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + let immune_drive = ((self.spleen_immune_activity_index - 80.0) / 80.0) + + (self.waste_load_mg - BASE_WASTE_LOAD_MG) / 1400.0 + + (self.lymphatic_return_ml_min - BASE_LYMPH_RETURN_ML_MIN) * 0.05; + let immune_drive = immune_drive.clamp(-0.6, 2.4); + let wbc_target = (BASE_WBC_COUNT_GIGA_PER_L * (1.0 + immune_drive * 0.55)).clamp(3.0, 26.0); + self.wbc_count_giga_per_l = Self::approach( + self.wbc_count_giga_per_l, + wbc_target, + 0.08, + dt_seconds, + ); + let neutro_target = (BASE_NEUTROPHIL_FRACTION + + immune_drive * 0.22 + + (self.lactate_mmol_l - BASE_LACTATE_MMOL_L).max(0.0) * 0.015) + .clamp(0.3, 0.88); + let lymph_target = (BASE_LYMPHOCYTE_FRACTION + - immune_drive * 0.12 + + (self.temperature_c - BASE_TEMP_C) * -0.03) + .clamp(0.08, 0.5); + let mut mono_target = (BASE_MONOCYTE_FRACTION + immune_drive * 0.05).clamp(0.03, 0.18); + let sum = neutro_target + lymph_target + mono_target; + let scale = if sum > 0.99 { 0.99 / sum } else { 1.0 }; + self.neutrophil_fraction = Self::approach( + self.neutrophil_fraction, + (neutro_target * scale).clamp(0.2, 0.9), + 0.06, + dt_seconds, + ); + self.lymphocyte_fraction = Self::approach( + self.lymphocyte_fraction, + (lymph_target * scale).clamp(0.05, 0.6), + 0.06, + dt_seconds, + ); + mono_target = (mono_target * scale).clamp(0.02, 0.25); + self.monocyte_fraction = Self::approach( + self.monocyte_fraction, + mono_target, + 0.06, + dt_seconds, + ); + let complement_target = (BASE_COMPLEMENT_ACTIVITY + + immune_drive * 0.35 + + (self.plasma_globulin_g_dl - BASE_PLASMA_GLOBULIN_G_DL) * 0.18) + .clamp(0.4, 2.0); + self.complement_activity = Self::approach( + self.complement_activity, + complement_target, + 0.04, + dt_seconds, + ); + let hypoxic_pressure = (1.0 - self.oxygen_supply_demand_ratio).max(0.0); + let inflammation_target = (immune_drive * 0.8 + + (self.temperature_c - BASE_TEMP_C) * 0.4 + + hypoxic_pressure * 1.6 + + (self.lactate_mmol_l - BASE_LACTATE_MMOL_L).max(0.0) * 0.3) + .clamp(0.0, 3.0); + self.inflammation_index = Self::approach( + self.inflammation_index, + inflammation_target, + 0.1, + dt_seconds, + ); + } + + fn update_acid_base(&mut self, dt_seconds: f32) { + if dt_seconds <= 0.0 { + return; + } + let pa_co2 = self.co2_content_ml_dl.clamp(20.0, 70.0); + self.arterial_pco2_mm_hg = pa_co2; + let renal_buffer = (self.renal_clearance_ml_min / BASE_RENAL_CLEARANCE_ML_MIN - 1.0) * 3.2; + let metabolic_acid = (self.lactate_mmol_l - BASE_LACTATE_MMOL_L) * 1.8 + + (self.waste_load_mg - BASE_WASTE_LOAD_MG) / 600.0; + let respiratory_pressure = (pa_co2 - BASE_ALVEOLAR_PCO2_MMHG) * 0.12; + let bicarbonate_target = (BASE_BICARBONATE_MMOL_L + + renal_buffer + - metabolic_acid + + respiratory_pressure) + .clamp(12.0, 36.0); + self.bicarbonate_mmol_l = Self::approach( + self.bicarbonate_mmol_l, + bicarbonate_target, + 0.08, + dt_seconds, + ); + let ratio = (self.bicarbonate_mmol_l / (0.03 * pa_co2.max(10.0))).max(1e-3); + let ph_calc = 6.1 + ratio.log10(); + let ph_target = ph_calc.clamp(7.05, 7.55); + self.ph = Self::approach(self.ph, ph_target, 0.1, dt_seconds); + let base_excess = 0.93 * (self.bicarbonate_mmol_l - 24.4 + 14.8 * (self.ph - 7.4)); + self.base_excess_mmol_l = Self::approach( + self.base_excess_mmol_l, + base_excess, + 0.12, + dt_seconds, + ); + let albumin_correction = (4.0 - self.plasma_albumin_g_dl) * 2.5; + let globulin_effect = (self.plasma_globulin_g_dl - BASE_PLASMA_GLOBULIN_G_DL) * 1.2; + let anion_gap_target = (BASE_ANION_GAP_MMOL_L + + (self.lactate_mmol_l - BASE_LACTATE_MMOL_L) * 2.5 + + metabolic_acid * 1.2 + + globulin_effect + + albumin_correction) + .clamp(4.0, 32.0); + self.anion_gap_mmol_l = Self::approach( + self.anion_gap_mmol_l, + anion_gap_target, + 0.08, dt_seconds, ); } fn update_plasma_volume(&mut self, dt_seconds: f32) { - self.plasma_volume_l = Self::approach( - self.plasma_volume_l, - self.plasma_volume_target_l.clamp(2.4, 3.8), - 0.004, - dt_seconds, - ); + let oncotic_bias = + (self.oncotic_pressure_mm_hg - BASE_PLASMA_ONCOTIC_PRESSURE_MMHG) * 0.018; + let lymph_bias = (self.lymphatic_return_ml_min - BASE_LYMPH_RETURN_ML_MIN) / 800.0; + let target = (self.plasma_volume_target_l + oncotic_bias + lymph_bias).clamp(2.2, 4.0); + self.plasma_volume_l = Self::approach(self.plasma_volume_l, target, 0.004, dt_seconds); } fn compute_oxygen_content(&self) -> f32 { @@ -327,8 +927,10 @@ impl Organ for Bloodstream { return; } + self.update_plasma_proteins(dt_seconds); self.update_red_cell_mass(dt_seconds); self.update_plasma_volume(dt_seconds); + self.update_platelets(dt_seconds); self.total_volume_l = (self.plasma_volume_l + self.red_cell_volume_l).clamp(3.6, 6.5); self.cardiac_output_l_min = Self::approach( @@ -428,11 +1030,8 @@ impl Organ for Bloodstream { let lactate_target = (BASE_LACTATE_MMOL_L + anaerobic_pressure * 6.0).clamp(0.6, 7.0); self.lactate_mmol_l = Self::approach(self.lactate_mmol_l, lactate_target, 0.2, dt_seconds); - let ph_target = (BASE_PH - - (self.co2_content_ml_dl - BASE_ALVEOLAR_PCO2_MMHG) * 0.015 - - (self.lactate_mmol_l - BASE_LACTATE_MMOL_L) * 0.05) - .clamp(7.1, 7.48); - self.ph = Self::approach(self.ph, ph_target, 0.05, dt_seconds); + self.update_immune_cells(dt_seconds); + self.update_acid_base(dt_seconds); let temp_target = (BASE_TEMP_C + (self.metabolic_demand_target_ml_min - BASE_METABOLIC_DEMAND_ML_MIN) / 450.0) @@ -526,12 +1125,114 @@ mod tests { for _ in 0..600 { blood.update(0.5); } - let baseline_hct = blood.hematocrit_pct; + let baseline = blood.metrics(); blood.ingest_erythropoietin(42.0); - blood.set_spleen_feedback(65.0, 220.0); + blood.set_plasma_proteins(4.3, 2.6, 0.32); + blood.set_erythropoiesis_nutrients(4.0, 0.9, 12.0); + blood.set_spleen_feedback(65.0, 220.0, 80, 2.8); for _ in 0..7200 { blood.update(0.25); } - assert!(blood.hematocrit_pct > baseline_hct + 0.5); + let metrics = blood.metrics(); + assert!( + metrics.erythropoiesis_rate_ml_per_day > baseline.erythropoiesis_rate_ml_per_day + 5.0 + ); + assert!(metrics.rbc_young_fraction > baseline.rbc_young_fraction); + assert!(metrics.oxygen_supply_demand_ratio >= baseline.oxygen_supply_demand_ratio - 0.05); + } + + #[test] + fn low_albumin_reduces_oncotic_pressure() { + let mut blood = Bloodstream::new("oncotic"); + for _ in 0..720 { + blood.update(0.5); + } + let baseline_oncotic = blood.oncotic_pressure_mm_hg; + blood.set_plasma_proteins(2.4, 1.4, 0.2); + for _ in 0..720 { + blood.update(0.5); + } + assert!(blood.oncotic_pressure_mm_hg < baseline_oncotic - 1.5); + assert!(blood.lymphatic_return_ml_min > BASE_LYMPH_RETURN_ML_MIN - 0.2); + } + + #[test] + fn nutrient_limitation_curbs_erythropoiesis() { + let mut blood = Bloodstream::new("nutrients"); + blood.ingest_erythropoietin(80.0); + blood.set_plasma_proteins(4.2, 2.3, 0.28); + blood.set_erythropoiesis_nutrients(0.4, 0.05, 0.5); + for _ in 0..1800 { + blood.update(0.5); + } + let limited = blood.metrics(); + blood.set_erythropoiesis_nutrients(4.2, 0.8, 12.0); + for _ in 0..1800 { + blood.update(0.5); + } + let enriched = blood.metrics(); + assert!(enriched.rbc_young_fraction >= limited.rbc_young_fraction); + assert!(enriched.iron_store_mg > 60.0); + assert!(enriched.hif_activation <= limited.hif_activation + 0.25); + } + + #[test] + fn platelet_release_increases_count() { + let mut blood = Bloodstream::new("platelet-release"); + for _ in 0..600 { + blood.update(0.5); + } + let baseline = blood.platelet_count_1e3_per_ul; + blood.set_spleen_feedback(30.0, 250.0, 100, 2.8); + blood.set_plasma_proteins(4.4, 2.9, 0.36); + for _ in 0..1800 { + blood.update(0.5); + } + assert!(blood.platelet_count_1e3_per_ul > baseline + 10.0); + assert!(blood.thrombosis_risk_index > -0.4); + } + + #[test] + fn immune_activation_shifts_wbc_profile() { + let mut blood = Bloodstream::new("immune-shift"); + for _ in 0..600 { + blood.update(0.5); + } + let baseline_wbc = blood.wbc_count_giga_per_l; + blood.set_spleen_feedback(58.0, 260.0, 180, 3.5); + for _ in 0..1800 { + blood.update(0.5); + } + assert!(blood.wbc_count_giga_per_l > baseline_wbc + 1.0); + assert!(blood.neutrophil_fraction > blood.lymphocyte_fraction); + assert!(blood.inflammation_index > 0.3); + } + + #[test] + fn respiratory_acidosis_adjusts_bicarbonate() { + let mut blood = Bloodstream::new("acid-base"); + for _ in 0..600 { + blood.update(0.5); + } + let baseline_hco3 = blood.bicarbonate_mmol_l; + blood.set_respiratory_exchange(94.0, 58.0, 820.0, 90.0, 0.9); + blood.set_renal_feedback(0.55, 3.1, 640.0); + for _ in 0..3600 { + blood.update(0.5); + } + assert!(blood.arterial_pco2_mm_hg > 45.0); + assert!(blood.ph < 7.38); + assert!(blood.bicarbonate_mmol_l > baseline_hco3 + 0.5); + assert!(blood.base_excess_mmol_l > -3.0); } } + + + + + + + + + + diff --git a/src/organs/brain.rs b/src/organs/brain.rs index 796e1fd..62ffc9b 100644 --- a/src/organs/brain.rs +++ b/src/organs/brain.rs @@ -6,6 +6,13 @@ const CIRCADIAN_PERIOD_SECONDS: f32 = 24.0 * 3600.0; const HOMEOSTATIC_WAKE_ACCUMULATION_S: f32 = 16.0 * 3600.0; const HOMEOSTATIC_DISCHARGE_S: f32 = 6.0 * 3600.0; +/// Median sleep cycle duration from large-scale polysomnography study (Åkerstedt et al. 2024). +/// Sleep Health: 96 min median across 6,064 sleep cycles. +const TYPICAL_CYCLE_DURATION_S: f32 = 96.0 * 60.0; + +/// REM latency typically 60-120 minutes in healthy adults (Carskadon & Dement, StatPearls 2024). +const FIRST_REM_LATENCY_S: f32 = 90.0 * 60.0; + /// Sleep architecture stages used for brain state transitions. #[derive(Debug, Clone, Copy, PartialEq, Eq)] pub enum SleepStage { @@ -57,7 +64,79 @@ pub struct Brain { pub autonomic_variability: f32, /// Cognitive/task load (0..=1) influencing metabolic demand. pub cognitive_load: f32, + /// Chemoreceptor-driven urge to ventilate (0..~1.3). + pub respiratory_drive: f32, + /// Likelihood of hypocapnic syncope (0..=1). + pub syncope_propensity: f32, + /// Osmotic/volume-inspired thirst drive (0..~1.4). + pub thirst_drive: f32, + /// Nutritional hunger drive (0..~1.3). + pub hunger_drive: f32, + /// Thermoregulatory discomfort drive (0..~1.2). + pub thermoregulatory_drive: f32, + /// Aggregated nociceptive/pain drive (0..~1.5). + pub pain_drive: f32, time_in_stage_s: f32, + /// Total time asleep in current sleep bout (seconds). + time_asleep_s: f32, + /// Cumulative count of sleep cycles completed in current sleep bout. + sleep_cycle_count: u32, + /// Sleep spindle density (count per minute) during N2 sleep. + /// Normal range: 2-5 spindles/min (Fernandez & Lüthi, Nat Neurosci 2020). + pub spindle_density_per_min: f32, + /// K-complex count per minute during N2 sleep. + /// Typically 1-3 per minute (Colrain, Sleep Med Rev 2005). + pub k_complex_density_per_min: f32, + /// Sigma band power (11-16 Hz) relative to baseline, proxy for spindle activity. + /// High correlation with spindle amplitude r=0.95 (Warby et al., Nat Methods 2014). + pub sigma_power_relative: f32, + /// Chin EMG tonic activity during REM as fraction of NREM baseline. + /// Normal REM: <10% tonic activity (Frauscher et al., Sleep 2014). + pub rem_tonic_emg_fraction: f32, + /// Chin EMG phasic burst activity during REM (% of 30-s epoch). + /// Normal: <15% phasic activity (Frauscher et al., Sleep 2014). + pub rem_phasic_emg_pct: f32, + /// REM atonia index: fraction of REM sleep with muscle atonia. + /// Normal: >85-90% atonia (McCarter et al., Sleep 2014). + pub rem_atonia_index: f32, + // Brainstem nuclei activities (0..=1 representing normalized firing rates) + /// Nucleus tractus solitarius (NTS) activity processing baroreceptor input. + /// Receives afferents from CN IX/X; inhibits RVLM (Dampney, Clin Exp Pharmacol Physiol 2016). + pub nts_activity: f32, + /// Rostral ventrolateral medulla (RVLM) sympathetic premotor neuron activity. + /// Primary source of tonic sympathetic vasomotor drive (Guyenet, Nat Rev Neurosci 2006). + pub rvlm_activity: f32, + /// Caudal ventrolateral medulla (CVLM) inhibitory interneuron activity. + /// Receives NTS input; inhibits RVLM to reduce sympathetic tone. + pub cvlm_activity: f32, + /// Nucleus ambiguus (nAmb) cardiovagal neuron activity. + /// Primary cardiac vagal efferent controlling HR via myelinated fibers (McAllen et al., J Physiol 2011). + pub namb_cardiovagal_activity: f32, + /// Dorsal motor nucleus of vagus (DMV) cardiac neuron activity. + /// ~20% of cardiovagal neurons; unmyelinated C-fibers; modulates contractility (Gourine et al., iScience 2024). + pub dmv_cardiac_activity: f32, + /// Retrotrapezoid nucleus (RTN) central chemoreceptor neuron activity. + /// CO2/pH sensitive; drives 2/3 of ventilatory CO2 response (Guyenet & Bayliss, J Neurosci 2025). + pub rtn_chemoreceptor_activity: f32, + /// Peripheral chemoreceptor signal from carotid/aortic bodies via NTS. + /// Primarily O2, also CO2/pH; ~1/3 of CO2 response (Guyenet & Bayliss, J Appl Physiol 2010). + pub peripheral_chemo_signal: f32, + // Cerebrovascular autoregulation parameters + /// Cerebrovascular resistance (mmHg·min·100g/mL). + /// Normal ~1.4 at baseline CPP=70 mmHg, CBF=50 mL/100g/min (Paulson et al., Cerebrovasc Dis 1990). + pub cerebrovascular_resistance: f32, + /// Lower autoregulation limit (mmHg CPP). + /// Normal ~50 mmHg; shifts with chronic hypertension (Drummond, BJA 2024). + pub autoregulation_lower_limit: f32, + /// Upper autoregulation limit (mmHg CPP). + /// Normal ~150 mmHg; breakthrough above this causes edema (Drummond, BJA 2024). + pub autoregulation_upper_limit: f32, + /// Autoregulation curve slope (normalized 0..=1). + /// 0 = no autoregulation (passive pressure-flow), 1 = perfect autoregulation. + pub autoregulation_efficiency: f32, + /// CO2 reactivity: ΔCBF per mmHg ΔPaCO2 (mL/100g/min per mmHg). + /// Normal 1-2 mL/100g/min/mmHg (Claassen et al., Physiol Rev 2021). + pub co2_reactivity: f32, } impl Brain { @@ -83,7 +162,33 @@ impl Brain { seizure_risk: 0.05, autonomic_variability: 0.35, cognitive_load: 0.35, + respiratory_drive: 0.6, + syncope_propensity: 0.05, + thirst_drive: 0.3, + hunger_drive: 0.35, + thermoregulatory_drive: 0.5, + pain_drive: 0.2, time_in_stage_s: 0.0, + time_asleep_s: 0.0, + sleep_cycle_count: 0, + spindle_density_per_min: 3.5, + k_complex_density_per_min: 2.0, + sigma_power_relative: 1.0, + rem_tonic_emg_fraction: 0.05, + rem_phasic_emg_pct: 8.0, + rem_atonia_index: 0.92, + nts_activity: 0.5, + rvlm_activity: 0.65, + cvlm_activity: 0.45, + namb_cardiovagal_activity: 0.4, + dmv_cardiac_activity: 0.35, + rtn_chemoreceptor_activity: 0.6, + peripheral_chemo_signal: 0.55, + cerebrovascular_resistance: 1.4, + autoregulation_lower_limit: 50.0, + autoregulation_upper_limit: 150.0, + autoregulation_efficiency: 0.85, + co2_reactivity: 1.5, } } @@ -116,12 +221,37 @@ impl Brain { fn transition_stage(&mut self, stage: SleepStage) { if self.sleep_stage != stage { + // Completing a REM period marks the end of a sleep cycle + if self.sleep_stage == SleepStage::Rem && stage != SleepStage::Wake { + self.sleep_cycle_count += 1; + } + // Waking resets sleep bout tracking + if stage == SleepStage::Wake && self.sleep_stage != SleepStage::Wake { + if self.time_asleep_s > 600.0 { + // Only reset if woke from substantial sleep (>10 min) + self.time_asleep_s = 0.0; + self.sleep_cycle_count = 0; + } + } self.sleep_stage = stage; self.time_in_stage_s = 0.0; } } + /// Adaptive cycle timing based on Åkerstedt et al. 2024 (Sleep Health) and homeostatic factors. + /// First cycle is shorter; high sleep pressure extends NREM; low pressure extends REM. fn evaluate_sleep_stage(&mut self, base_arousal: f32) { + // Cycle modulation: first cycle is shorter (Åkerstedt 2024) + let cycle_factor = if self.sleep_cycle_count == 0 { + 0.85 + } else { + 1.0 + }; + + // High sleep pressure extends NREM, low pressure extends REM (Åkerstedt 2024) + let nrem_extension = (self.sleep_pressure - 0.5).max(0.0) * 1.2; + let rem_extension = (0.5 - self.sleep_pressure).max(0.0) * 1.3; + match self.sleep_stage { SleepStage::Wake => { if base_arousal < 0.35 && self.sleep_pressure > 0.55 { @@ -129,30 +259,59 @@ impl Brain { } } SleepStage::N1 => { + // N1 typically 1-7 minutes (Carskadon & Dement, StatPearls 2024) + let n1_duration = 180.0 * cycle_factor; if base_arousal > 0.5 { self.transition_stage(SleepStage::Wake); - } else if self.time_in_stage_s > 180.0 && self.sleep_pressure > 0.6 { + } else if self.time_in_stage_s > n1_duration && self.sleep_pressure > 0.6 { self.transition_stage(SleepStage::N2); } } SleepStage::N2 => { + // N2 typically 10-25 minutes initially, longer in later cycles + let base_n2_duration = 900.0 * cycle_factor * (1.0 + nrem_extension); + let n3_threshold = base_n2_duration * (1.0 + nrem_extension); + + // First REM latency 60-120 min (typically 90 min, StatPearls 2024) + let rem_ready = if self.sleep_cycle_count == 0 { + self.time_asleep_s > FIRST_REM_LATENCY_S * 0.8 + } else { + self.time_in_stage_s > base_n2_duration * 1.5 + }; + if base_arousal > 0.52 { self.transition_stage(SleepStage::Wake); - } else if self.time_in_stage_s > 900.0 && self.sleep_pressure > 0.65 { + } else if self.time_in_stage_s > n3_threshold && self.sleep_pressure > 0.65 { self.transition_stage(SleepStage::N3); - } else if self.time_in_stage_s > 2400.0 { + } else if rem_ready { self.transition_stage(SleepStage::Rem); } } SleepStage::N3 => { - if self.time_in_stage_s > 1800.0 { + // N3 duration decreases across night; more in first cycle + let n3_duration = if self.sleep_cycle_count == 0 { + 1800.0 * (1.0 + nrem_extension) + } else { + 1200.0 * cycle_factor * (1.0 + nrem_extension) + }; + + if self.time_in_stage_s > n3_duration { self.transition_stage(SleepStage::Rem); } else if base_arousal > 0.45 { self.transition_stage(SleepStage::N2); } } SleepStage::Rem => { - if self.time_in_stage_s > 1500.0 { + // REM duration increases across cycles; first REM ~5-10 min, later ~20-25 min + // Positive correlation: longer REM predicts longer interval to next REM + let base_rem_duration = if self.sleep_cycle_count == 0 { + 600.0 // ~10 min first REM + } else { + 900.0 + (self.sleep_cycle_count as f32 * 300.0).min(600.0) // up to 25 min + }; + let rem_duration = base_rem_duration * (1.0 + rem_extension); + + if self.time_in_stage_s > rem_duration { self.transition_stage(SleepStage::N2); } else if base_arousal > 0.65 { self.transition_stage(SleepStage::Wake); @@ -170,6 +329,193 @@ impl Brain { SleepStage::Rem => base_arousal.clamp(0.45, 0.7), } } + + /// Simulate brainstem nuclei interactions for autonomic control. + /// Based on: + /// - Guyenet (Nat Rev Neurosci 2006): RVLM as sympathetic vasomotor source + /// - Dampney (Clin Exp Pharmacol Physiol 2016): NTS→CVLM→RVLM baroreflex arc + /// - Gourine et al. (iScience 2024): nAmb/DMV cardiovagal pathways + /// - Guyenet & Bayliss (J Neurosci 2025): RTN central chemoreceptors + fn update_brainstem_nuclei(&mut self, dt_seconds: f32) { + // Baroreceptor input: elevated CPP → increased NTS activity + let baroreceptor_input = + ((self.cerebral_perfusion_pressure_mm_hg - 70.0) / 40.0).clamp(-0.5, 0.5) + 0.5; + + // Chemoreceptor inputs + // Central chemoreceptors (RTN): sensitive to CO2/pH (proxied by respiratory drive) + let central_chemo_input = self.respiratory_drive.clamp(0.3, 1.2); + + // Peripheral chemoreceptors: O2 (via oxygenation) + CO2 contribution + let o2_signal = (0.98 - self.oxygenation_saturation).max(0.0) * 3.0; + let co2_signal = (self.respiratory_drive - 0.6).max(0.0) * 0.4; + let peripheral_input = (0.5 + o2_signal + co2_signal).clamp(0.3, 1.3); + + // NTS integrates baroreceptor and peripheral chemoreceptor inputs + // High baroreceptor → high NTS → inhibition of sympathetic tone + // Chemoreceptors also converge at NTS (interfere with baroreflex per 2024 study) + let nts_target = + (0.4 * baroreceptor_input + 0.35 * peripheral_input + 0.25 * central_chemo_input) + .clamp(0.2, 1.0); + self.nts_activity = Self::approach(self.nts_activity, nts_target, 1.5, dt_seconds); + + // CVLM receives excitatory input from NTS and inhibits RVLM + let cvlm_target = (0.3 + 0.7 * self.nts_activity).clamp(0.2, 0.95); + self.cvlm_activity = Self::approach(self.cvlm_activity, cvlm_target, 1.2, dt_seconds); + + // RVLM: tonic sympathetic drive, inhibited by CVLM + // Baseline drive adjusted by arousal and sleep stage + let arousal_drive = self.cortical_arousal * 0.4; + let sleep_suppression = match self.sleep_stage { + SleepStage::Wake => 0.0, + SleepStage::N1 => 0.05, + SleepStage::N2 => 0.1, + SleepStage::N3 => 0.18, + SleepStage::Rem => 0.08, + }; + let rvlm_target = (0.65 + arousal_drive - 0.5 * self.cvlm_activity - sleep_suppression + + self.pain_drive * 0.15 + + (self.respiratory_drive - 0.6) * 0.2) + .clamp(0.25, 1.1); + self.rvlm_activity = Self::approach(self.rvlm_activity, rvlm_target, 0.8, dt_seconds); + + // Nucleus ambiguus cardiovagal neurons: primary HR control via myelinated vagal fibers + // High NTS activity (high BP) → high vagal tone → bradycardia + let vagal_excitation = self.nts_activity * 0.6; + let namb_target = (0.3 + vagal_excitation - self.cortical_arousal * 0.3 + + match self.sleep_stage { + SleepStage::Wake => 0.0, + SleepStage::N1 => 0.05, + SleepStage::N2 => 0.1, + SleepStage::N3 => 0.15, + SleepStage::Rem => -0.05, // REM has more variable autonomic tone + }) + .clamp(0.15, 0.85); + self.namb_cardiovagal_activity = + Self::approach(self.namb_cardiovagal_activity, namb_target, 1.0, dt_seconds); + + // DMV cardiac neurons: ~20% of cardiovagal neurons, modulate contractility via C-fibers + let dmv_target = (0.25 + 0.5 * self.nts_activity - 0.2 * self.cortical_arousal + + match self.sleep_stage { + SleepStage::N2 | SleepStage::N3 => 0.12, + _ => 0.0, + }) + .clamp(0.15, 0.75); + self.dmv_cardiac_activity = + Self::approach(self.dmv_cardiac_activity, dmv_target, 0.7, dt_seconds); + + // RTN central chemoreceptors: CO2/pH drive (2/3 of ventilatory CO2 response) + let rtn_target = central_chemo_input.clamp(0.4, 1.1); + self.rtn_chemoreceptor_activity = + Self::approach(self.rtn_chemoreceptor_activity, rtn_target, 0.9, dt_seconds); + + // Peripheral chemoreceptor signal (1/3 of CO2 response, primary O2 sensor) + self.peripheral_chemo_signal = Self::approach( + self.peripheral_chemo_signal, + peripheral_input, + 1.2, + dt_seconds, + ); + + // Integrate brainstem outputs into autonomic drive and respiratory drive + // RVLM drives sympathetic tone → brainstem_autonomic_drive + // Combined vagal tone (nAmb + DMV) provides parasympathetic counterbalance + let sympathetic_component = self.rvlm_activity * 0.6; + let parasympathetic_component = + (self.namb_cardiovagal_activity * 0.8 + self.dmv_cardiac_activity * 0.2) * 0.4; + let net_autonomic = + (sympathetic_component + (1.0 - parasympathetic_component)).clamp(0.25, 1.1); + + // Respiratory drive integrates RTN (central) and peripheral chemoreceptors + let chemo_drive = (self.rtn_chemoreceptor_activity * 0.65 + + self.peripheral_chemo_signal * 0.35) + .clamp(0.4, 1.3); + + // Blend computed autonomic/respiratory drives with existing targets + // (existing code uses these as feedback; brainstem provides mechanistic basis) + self.brainstem_autonomic_drive = Self::approach( + self.brainstem_autonomic_drive, + net_autonomic, + 0.5, + dt_seconds, + ) + .clamp(0.25, 1.1); + + self.respiratory_drive = + Self::approach(self.respiratory_drive, chemo_drive, 0.4, dt_seconds).clamp(0.4, 1.3); + } + + /// Cerebrovascular autoregulation model with myogenic, metabolic, and CO2 reactivity. + /// Based on: + /// - Paulson et al. (Cerebrovasc Dis 1990): CVR = CPP/CBF relationship + /// - Claassen et al. (Physiol Rev 2021): comprehensive autoregulation review + /// - Drummond (BJA 2024): monitoring and clinical implications + fn update_cerebrovascular_autoregulation(&mut self, dt_seconds: f32) { + let cpp = self.cerebral_perfusion_pressure_mm_hg; + + // Within autoregulation range (50-150 mmHg): resistance adjusts to maintain CBF ~50 mL/100g/min + // Below lower limit: vasodilation exhausted, passive pressure-flow + // Above upper limit: forced vasodilation (breakthrough) + let in_autoregulation_range = + cpp >= self.autoregulation_lower_limit && cpp <= self.autoregulation_upper_limit; + + let myogenic_resistance = if in_autoregulation_range { + // Active autoregulation: R = CPP / target_CBF + // Target CBF = 50 mL/100g/min baseline, modulated by metabolic demand + let target_cbf = 50.0 * self.metabolic_demand_fraction; + (cpp / target_cbf).clamp(0.8, 2.5) + } else if cpp < self.autoregulation_lower_limit { + // Below lower limit: vasodilation maximally exhausted, resistance low + let exhaustion_factor = (cpp / self.autoregulation_lower_limit).clamp(0.3, 1.0); + 0.8 * exhaustion_factor + } else { + // Above upper limit: forced vasodilation (breakthrough), resistance collapses + let breakthrough_factor = + 1.0 - ((cpp - self.autoregulation_upper_limit) / 50.0).clamp(0.0, 0.7); + 0.6 * breakthrough_factor.max(0.3) + }; + + // Metabolic regulation: hypoxia, hypercapnia → vasodilation (decreased resistance) + let hypoxia_factor = (0.95 - self.oxygenation_saturation).max(0.0) * 2.0; + let metabolic_resistance_reduction = hypoxia_factor * 0.25; + + // CO2 reactivity: 1-2 mL/100g/min per mmHg PaCO2 change (Claassen et al., Physiol Rev 2021) + // Proxy PaCO2 via respiratory drive: baseline ~40 mmHg at respiratory_drive=0.6 + let estimated_paco2 = 30.0 + (self.respiratory_drive - 0.4) * 25.0; // rough estimate + let paco2_deviation = estimated_paco2 - 40.0; + // CBF increases ~4% per mmHg CO2 (1-2 mL/100g/min at baseline 50 mL/100g/min) + // Resistance inversely related: higher CO2 → lower resistance + let co2_resistance_factor = + 1.0 - (paco2_deviation * self.co2_reactivity / 50.0).clamp(-0.3, 0.4); + + // Neurogenic modulation: sympathetic tone can override autoregulation partially + let sympathetic_constriction = (self.rvlm_activity - 0.65) * 0.12; + + // Combine factors + let target_resistance = (myogenic_resistance * co2_resistance_factor + - metabolic_resistance_reduction + + sympathetic_constriction) + .clamp(0.5, 3.0); + + // Smooth resistance changes (myogenic response time constant ~5-15 sec) + self.cerebrovascular_resistance = Self::approach( + self.cerebrovascular_resistance, + target_resistance, + 0.15, + dt_seconds, + ) + .clamp(0.5, 3.5); + + // Calculate CBF from CPP and resistance: CBF = CPP / CVR + let calculated_cbf = (cpp / self.cerebrovascular_resistance).clamp(25.0, 120.0); + + self.cerebral_blood_flow_ml_per_100g_min = Self::approach( + self.cerebral_blood_flow_ml_per_100g_min, + calculated_cbf, + 2.0, + dt_seconds, + ) + .clamp(20.0, 130.0); + } } impl Organ for Brain { @@ -185,13 +531,19 @@ impl Organ for Brain { } self.time_in_stage_s += dt_seconds; + + // Track total time asleep for cycle calculations + let asleep = !matches!(self.sleep_stage, SleepStage::Wake); + if asleep { + self.time_asleep_s += dt_seconds; + } + let circadian_increment = TAU * (dt_seconds / CIRCADIAN_PERIOD_SECONDS); self.circadian_phase_radians = Self::wrap_phase(self.circadian_phase_radians + circadian_increment); let circadian_arousal = 0.55 + 0.45 * (self.circadian_phase_radians - TAU * 0.25).sin(); - let asleep = !matches!(self.sleep_stage, SleepStage::Wake); let sleep_pressure_delta = if asleep { -dt_seconds / HOMEOSTATIC_DISCHARGE_S } else { @@ -200,26 +552,50 @@ impl Organ for Brain { self.sleep_pressure = (self.sleep_pressure + sleep_pressure_delta).clamp(0.0, 1.1); let base_arousal = (circadian_arousal - 0.55 * self.sleep_pressure).clamp(0.05, 1.0); - self.evaluate_sleep_stage(base_arousal); + let thirst_component = (self.thirst_drive - 0.3).max(0.0) * 0.18; + let hunger_component = (self.hunger_drive - 0.35).max(0.0) * 0.15; + let thermo_component = (self.thermoregulatory_drive - 0.5).abs() * 0.12; + let pain_component = self.pain_drive * 0.2; + let drive_bonus = (self.respiratory_drive - 0.6).max(0.0) * 0.25 + + thirst_component + + hunger_component + + thermo_component + + pain_component; + let syncope_penalty = self.syncope_propensity * 0.35; + let drive_weighted_arousal = + (base_arousal + drive_bonus - syncope_penalty).clamp(0.05, 1.0); + self.evaluate_sleep_stage(drive_weighted_arousal); - let arousal_target = self.stage_arousal_target(base_arousal); + let arousal_target = self.stage_arousal_target(drive_weighted_arousal); self.cortical_arousal = Self::approach(self.cortical_arousal, arousal_target, 0.8, dt_seconds).clamp(0.05, 1.0); - let (autonomic_target, variability_target) = match self.sleep_stage { + // Update brainstem nuclei (baroreflex, chemoreflex, vagal pathways) + self.update_brainstem_nuclei(dt_seconds); + + let (base_autonomic_target, base_variability_target) = match self.sleep_stage { SleepStage::Wake => (0.9, 0.35), SleepStage::N1 => (0.82, 0.4), SleepStage::N2 => (0.78, 0.45), SleepStage::N3 => (0.72, 0.3), SleepStage::Rem => (0.86, 0.6), }; + let drive_push = (self.respiratory_drive - 0.6) * 0.35; + let syncope_pull = self.syncope_propensity * 0.45; + let autonomic_target = (base_autonomic_target + drive_push - syncope_pull).clamp(0.25, 1.1); + let variability_target = (base_variability_target + + (self.respiratory_drive - 0.6).abs() * 0.15 + + self.syncope_propensity * 0.1 + + self.pain_drive * 0.12 + + (self.thermoregulatory_drive - 0.5).abs() * 0.08) + .clamp(0.05, 0.95); self.brainstem_autonomic_drive = Self::approach( self.brainstem_autonomic_drive, autonomic_target, 0.6, dt_seconds, ) - .clamp(0.3, 1.1); + .clamp(0.25, 1.1); self.autonomic_variability = Self::approach( self.autonomic_variability, variability_target, @@ -246,7 +622,11 @@ impl Organ for Brain { }; let metabolic_target = (base_metabolic + 0.25 * (self.cognitive_load - 0.2) - + 0.1 * (self.glutamate_level - self.gaba_level)) + + 0.1 * (self.glutamate_level - self.gaba_level) + + 0.12 * (self.respiratory_drive - 0.6) + + 0.14 * (self.hunger_drive - 0.35) + + 0.06 * (self.thermoregulatory_drive - 0.5) + - 0.1 * self.syncope_propensity) .clamp(0.6, 1.4); self.metabolic_demand_fraction = Self::approach( self.metabolic_demand_fraction, @@ -286,22 +666,16 @@ impl Organ for Brain { .clamp(0.5, 35.0); let oxygen_target = (0.95 + 0.03 * self.brainstem_autonomic_drive - - 0.04 * (self.metabolic_demand_fraction - 1.0)) - .clamp(0.88, 0.99); + - 0.04 * (self.metabolic_demand_fraction - 1.0) + - 0.05 * self.syncope_propensity + + 0.02 * (self.respiratory_drive - 0.6)) + .clamp(0.85, 0.99); self.oxygenation_saturation = Self::approach(self.oxygenation_saturation, oxygen_target, 0.4, dt_seconds) .clamp(0.8, 1.0); - let cbf_target = (50.0 * self.metabolic_demand_fraction - + 0.25 * (self.cerebral_perfusion_pressure_mm_hg - 70.0)) - .clamp(30.0, 90.0); - self.cerebral_blood_flow_ml_per_100g_min = Self::approach( - self.cerebral_blood_flow_ml_per_100g_min, - cbf_target, - 1.6, - dt_seconds, - ) - .clamp(25.0, 95.0); + // Cerebrovascular autoregulation replaces static CBF calculation + self.update_cerebrovascular_autoregulation(dt_seconds); let stage_icp_term = match self.sleep_stage { SleepStage::Wake => 0.0, @@ -336,7 +710,9 @@ impl Organ for Brain { let activity = (self.cortical_arousal * 0.6 + self.metabolic_demand_fraction * 0.25 + (self.glutamate_level - self.gaba_level + 0.7) * 0.1 - + self.dopamine_tone * 0.05) + + self.dopamine_tone * 0.05 + + self.pain_drive * 0.08 + + (self.thirst_drive + self.hunger_drive) * 0.05) .clamp(0.1, 2.3); self.activity_index = activity; @@ -348,10 +724,106 @@ impl Organ for Brain { _ => 3.0, }; let target_consciousness = - ((self.cortical_arousal * 70.0 * perfusion_factor * oxygen_factor) + stage_bonus) + ((self.cortical_arousal * 70.0 * perfusion_factor * oxygen_factor) + + stage_bonus + + (self.respiratory_drive - 0.6) * 15.0 + - self.syncope_propensity * 25.0 + + self.pain_drive * 12.0 + + (self.thirst_drive + self.hunger_drive) * 6.0 + - (self.thermoregulatory_drive - 0.5).abs() * 10.0) .clamp(0.0, 100.0); self.consciousness = target_consciousness.round() as u8; + // Polysomnographic markers (Warby et al. 2014, Frauscher et al. 2014) + match self.sleep_stage { + SleepStage::N2 => { + // Sleep spindle density 2-5/min normal (Fernandez & Lüthi, Nat Neurosci 2020) + // Spindle density correlates with sleep pressure and cycle number + let spindle_base = 3.5 + (self.sleep_pressure - 0.5) * 1.5; + let cycle_modulation = if self.sleep_cycle_count == 0 { + 1.1 + } else { + 0.95 + }; + let spindle_target = (spindle_base * cycle_modulation).clamp(2.0, 5.5); + self.spindle_density_per_min = Self::approach( + self.spindle_density_per_min, + spindle_target, + 0.3, + dt_seconds, + ); + + // K-complex density 1-3/min (Colrain, Sleep Med Rev 2005) + let k_complex_target = (2.0 + (self.sleep_pressure - 0.6) * 1.0).clamp(1.0, 3.2); + self.k_complex_density_per_min = Self::approach( + self.k_complex_density_per_min, + k_complex_target, + 0.25, + dt_seconds, + ); + + // Sigma power (11-16 Hz) correlates strongly with spindle amplitude r=0.95 + let sigma_target = + (1.0 + (self.spindle_density_per_min - 3.5) * 0.15).clamp(0.7, 1.4); + self.sigma_power_relative = + Self::approach(self.sigma_power_relative, sigma_target, 0.2, dt_seconds); + } + SleepStage::Rem => { + // REM atonia: normal >85-90% atonia (McCarter et al., Sleep 2014) + // Tonic EMG <10%, phasic <15% in normal REM (Frauscher et al., Sleep 2014) + + // Arousal and sleep pressure affect atonia integrity + let atonia_disruption = (self.cortical_arousal - 0.5).max(0.0) * 0.15 + + (0.4 - self.sleep_pressure).max(0.0) * 0.08 + + self.pain_drive * 0.05; + + let atonia_target = (0.92 - atonia_disruption).clamp(0.75, 0.96); + self.rem_atonia_index = + Self::approach(self.rem_atonia_index, atonia_target, 0.15, dt_seconds); + + // Tonic EMG is inverse of atonia + let tonic_target = ((1.0 - self.rem_atonia_index) * 10.0).clamp(2.0, 18.0); + self.rem_tonic_emg_fraction = + Self::approach(self.rem_tonic_emg_fraction, tonic_target, 0.12, dt_seconds); + + // Phasic bursts: normal <15%, increase with arousal and dream intensity + let phasic_base = 8.0 + self.cortical_arousal * 6.0 + atonia_disruption * 12.0; + let phasic_target = phasic_base.clamp(5.0, 25.0); + self.rem_phasic_emg_pct = + Self::approach(self.rem_phasic_emg_pct, phasic_target, 0.18, dt_seconds); + } + SleepStage::N1 | SleepStage::N3 => { + // Gradually return markers to baseline when not in N2 or REM + self.spindle_density_per_min = + Self::approach(self.spindle_density_per_min, 3.5, 0.2, dt_seconds); + self.k_complex_density_per_min = + Self::approach(self.k_complex_density_per_min, 2.0, 0.2, dt_seconds); + self.sigma_power_relative = + Self::approach(self.sigma_power_relative, 1.0, 0.15, dt_seconds); + self.rem_tonic_emg_fraction = + Self::approach(self.rem_tonic_emg_fraction, 0.05, 0.1, dt_seconds); + self.rem_phasic_emg_pct = + Self::approach(self.rem_phasic_emg_pct, 8.0, 0.12, dt_seconds); + self.rem_atonia_index = + Self::approach(self.rem_atonia_index, 0.92, 0.1, dt_seconds); + } + SleepStage::Wake => { + // Reset to wake baseline + self.spindle_density_per_min = + Self::approach(self.spindle_density_per_min, 0.0, 0.5, dt_seconds); + self.k_complex_density_per_min = + Self::approach(self.k_complex_density_per_min, 0.0, 0.5, dt_seconds); + self.sigma_power_relative = + Self::approach(self.sigma_power_relative, 0.5, 0.3, dt_seconds); + self.rem_tonic_emg_fraction = + Self::approach(self.rem_tonic_emg_fraction, 0.8, 0.15, dt_seconds); + self.rem_phasic_emg_pct = + Self::approach(self.rem_phasic_emg_pct, 25.0, 0.2, dt_seconds); + self.rem_atonia_index = + Self::approach(self.rem_atonia_index, 0.15, 0.15, dt_seconds); + } + } + let excitability = (self.glutamate_level - self.gaba_level + self.cortical_arousal - 0.5).max(0.0); let hypoxia = (0.94 - self.oxygenation_saturation).max(0.0) * 4.0; @@ -360,10 +832,16 @@ impl Organ for Brain { } fn summary(&self) -> String { format!( - "Brain[id={}, stage={:?}, arousal={:.2}, ICP={:.1} mmHg, CPP={:.0} mmHg, O2={:.0}%, szrisk={:.0}%]", + "Brain[id={}, stage={:?}, arousal={:.2}, drives[r={:.2},th={:.2},hu={:.2},tmp={:.2},pain={:.2},syn={:.2}], ICP={:.1} mmHg, CPP={:.0} mmHg, O2={:.0}%, szrisk={:.0}%]", self.id(), self.sleep_stage, self.cortical_arousal, + self.respiratory_drive, + self.thirst_drive, + self.hunger_drive, + self.thermoregulatory_drive, + self.pain_drive, + self.syncope_propensity, self.intracranial_pressure_mm_hg, self.cerebral_perfusion_pressure_mm_hg, self.oxygenation_saturation * 100.0, diff --git a/src/organs/heart.rs b/src/organs/heart.rs index 380fd5d..899e2cb 100644 --- a/src/organs/heart.rs +++ b/src/organs/heart.rs @@ -15,6 +15,38 @@ pub enum CardiacRhythmState { Arrhythmic, } +/// Discrete cardiac cycle phases for phase-based state machine. +/// Based on: Sengupta PP et al. "The Cardiac Cycle" (2008) PMC2390899 +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum CardiacPhase { + /// Late diastole with atrial contraction completing ventricular filling + AtrialSystole, + /// Early systole with all valves closed, pressure rising without volume change + IsovolumetricContraction, + /// Mid-systole with semilunar valves open, blood ejection + Ejection, + /// Early diastole with all valves closed, pressure falling without volume change + IsovolumetricRelaxation, + /// Mid-to-late diastole with passive ventricular filling + PassiveFilling, +} + +/// Valve states for discrete valve modeling. +/// Based on: CV Physiology - Valvular function and hemodynamics +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub enum ValveState { + /// Fully closed, no flow + Closed, + /// Opening transition + Opening, + /// Fully open, forward flow + Open, + /// Closing transition + Closing, + /// Incompetent closure allowing regurgitant flow + Regurgitant, +} + /// Cardiac pump model featuring autonomic reflexes and hemodynamic coupling. #[derive(Debug, Clone)] pub struct Heart { @@ -66,6 +98,60 @@ pub struct Heart { /// Stroke work (J per beat). pub stroke_work_joule: f32, time_in_state_s: f32, + + // === NEW: Discrete atrial and valve modeling === + /// Left atrial pressure (mmHg). + pub left_atrial_pressure_mm_hg: f32, + /// Left atrial volume (ml). + pub left_atrial_volume_ml: f32, + /// Right atrial pressure (mmHg). + pub right_atrial_pressure_mm_hg: f32, + /// Right atrial volume (ml). + pub right_atrial_volume_ml: f32, + /// Mitral valve state. + pub mitral_valve_state: ValveState, + /// Aortic valve state. + pub aortic_valve_state: ValveState, + /// Tricuspid valve state. + pub tricuspid_valve_state: ValveState, + /// Pulmonic valve state. + pub pulmonic_valve_state: ValveState, + /// Mitral regurgitant fraction (0..=1). + pub mitral_regurgitation_fraction: f32, + /// Aortic regurgitant fraction (0..=1). + pub aortic_regurgitation_fraction: f32, + + // === NEW: Extended conduction system modeling === + /// His bundle conduction velocity (m/s). Normal: 1.0-1.5 m/s. + pub his_bundle_velocity_m_s: f32, + /// Purkinje fiber conduction velocity (m/s). Normal: 2.0-4.0 m/s (fastest in heart). + pub purkinje_velocity_m_s: f32, + /// His-Purkinje effective refractory period (ms). Normal: 350-450 ms. + pub his_purkinje_erp_ms: f32, + /// Ventricular action potential duration at 90% repolarization (ms). + pub ventricular_apd90_ms: f32, + /// Time since last ventricular depolarization (ms). + pub time_since_depolarization_ms: f32, + + // === NEW: Phase-based cardiac cycle state machine === + /// Current phase of cardiac cycle. + pub cardiac_phase: CardiacPhase, + /// Time in current cardiac phase (s). + pub phase_time_s: f32, + /// Left ventricular pressure (mmHg). + pub left_ventricular_pressure_mm_hg: f32, + /// Aortic pressure (mmHg). + pub aortic_pressure_mm_hg: f32, + + // === NEW: RAAS hormonal regulation === + /// Plasma renin activity (ng/mL/hr). Normal: 0.5-3.5. + pub plasma_renin_activity: f32, + /// Angiotensin II level (pg/mL). Normal: 10-60. + pub angiotensin_ii_pg_ml: f32, + /// Plasma aldosterone (ng/dL). Normal: 4-31. + pub aldosterone_ng_dl: f32, + /// Effective circulating volume fraction (0..=1.2). 1.0 = euvolemic. + pub effective_circulating_volume: f32, } impl Heart { @@ -97,6 +183,34 @@ impl Heart { coronary_perfusion_mm_hg: 70.0, stroke_work_joule: 1.1, time_in_state_s: 0.0, + // Atrial compartments + left_atrial_pressure_mm_hg: 8.0, + left_atrial_volume_ml: 55.0, + right_atrial_pressure_mm_hg: 4.0, + right_atrial_volume_ml: 50.0, + // Valve states + mitral_valve_state: ValveState::Open, + aortic_valve_state: ValveState::Closed, + tricuspid_valve_state: ValveState::Open, + pulmonic_valve_state: ValveState::Closed, + mitral_regurgitation_fraction: 0.0, + aortic_regurgitation_fraction: 0.0, + // Conduction system + his_bundle_velocity_m_s: 1.2, + purkinje_velocity_m_s: 3.0, + his_purkinje_erp_ms: 400.0, + ventricular_apd90_ms: 280.0, + time_since_depolarization_ms: 0.0, + // Cardiac cycle phase + cardiac_phase: CardiacPhase::PassiveFilling, + phase_time_s: 0.0, + left_ventricular_pressure_mm_hg: 8.0, + aortic_pressure_mm_hg: 80.0, + // RAAS + plasma_renin_activity: 1.5, + angiotensin_ii_pg_ml: 30.0, + aldosterone_ng_dl: 12.0, + effective_circulating_volume: 1.0, } } @@ -258,6 +372,302 @@ impl Heart { .clamp(0.05, 0.95); self.arrhythmia_burden = Self::approach(self.arrhythmia_burden, target, 0.3, dt_seconds); } + + /// Update RAAS hormonal cascade. + /// Based on: Fountain JH, Lappin SL. Physiology, Renin Angiotensin System. StatPearls 2024. + fn update_raas(&mut self, dt_seconds: f32) { + // Renin release triggered by: + // 1. Reduced renal perfusion (low MAP) + // 2. Sympathetic stimulation + // 3. Hypovolemia + let map_current = self.mean_arterial_pressure(); + let perfusion_deficit = ((85.0 - map_current) / 30.0).clamp(0.0, 1.0); + let sympathetic_drive = (self.autonomic_tone.max(0.0) * 0.6).clamp(0.0, 1.0); + let volume_deficit = ((1.0 - self.effective_circulating_volume) * 1.2).clamp(0.0, 1.0); + + let renin_stimulus = perfusion_deficit + sympathetic_drive + volume_deficit; + let renin_target = (0.5 + 3.0 * renin_stimulus).clamp(0.2, 5.0); + self.plasma_renin_activity = + Self::approach(self.plasma_renin_activity, renin_target, 0.15, dt_seconds); + + // Angiotensin II production via ACE in pulmonary circulation + // Proportional to renin with enzymatic conversion delay + let ang_ii_target = (10.0 + 45.0 * (self.plasma_renin_activity / 3.5)).clamp(5.0, 120.0); + self.angiotensin_ii_pg_ml = + Self::approach(self.angiotensin_ii_pg_ml, ang_ii_target, 0.4, dt_seconds); + + // Aldosterone secretion from adrenal cortex + // Driven by angiotensin II and hyperkalemia (not modeled) + let aldo_target = + (4.0 + 30.0 * (self.angiotensin_ii_pg_ml / 60.0).powf(1.2)).clamp(2.0, 50.0); + self.aldosterone_ng_dl = + Self::approach(self.aldosterone_ng_dl, aldo_target, 0.25, dt_seconds); + + // Angiotensin II effects on SVR (vasoconstriction) + // Normal Ang II ~ 30 pg/mL; elevated levels increase SVR + let ang_ii_vasoconstriction = ((self.angiotensin_ii_pg_ml - 30.0) / 30.0).clamp(-0.3, 1.5); + let svr_raas_component = 2.0 * ang_ii_vasoconstriction; + + // Aldosterone effects on volume (increase preload over hours, simplified here) + let volume_retention_rate = (self.aldosterone_ng_dl - 12.0) / 40.0; + self.effective_circulating_volume = (self.effective_circulating_volume + + volume_retention_rate * dt_seconds / 3600.0) + .clamp(0.7, 1.3); + + // Apply RAAS-mediated SVR adjustment + let base_svr = BASE_SVR_MMHG_MIN_PER_L + 5.5 * self.autonomic_tone; + let raas_svr_target = (base_svr + svr_raas_component).clamp(10.0, 28.0); + self.systemic_vascular_resistance = Self::approach( + self.systemic_vascular_resistance, + raas_svr_target, + 0.3, + dt_seconds, + ); + } + + /// Update conduction system properties including His-Purkinje pathways and refractoriness. + /// Based on: Cardiac electrophysiology, refractory periods, and APD dynamics. + fn update_conduction_system(&mut self, dt_seconds: f32) { + let dt_ms = dt_seconds * 1000.0; + self.time_since_depolarization_ms += dt_ms; + + // His-Purkinje conduction velocity modulated by: + // - Sympathetic tone (increased catecholamines increase velocity) + // - Ischemia/hypoxia (decreases velocity) + let supply_demand_ratio = self.myocardial_oxygen_supply / self.myocardial_oxygen_demand; + let ischemia_factor = supply_demand_ratio.clamp(0.6, 1.2); + let sympathetic_boost = 1.0 + 0.15 * self.autonomic_tone.max(0.0); + + let his_target = (1.2 * ischemia_factor * sympathetic_boost).clamp(0.5, 1.8); + let purkinje_target = (3.0 * ischemia_factor * sympathetic_boost).clamp(1.2, 4.5); + + self.his_bundle_velocity_m_s = + Self::approach(self.his_bundle_velocity_m_s, his_target, 0.08, dt_seconds); + self.purkinje_velocity_m_s = Self::approach( + self.purkinje_velocity_m_s, + purkinje_target, + 0.15, + dt_seconds, + ); + + // APD and ERP modulated by rate, autonomic tone, and ischemia + // Rate adaptation: faster HR → shorter APD (protective against reentry) + let rate_adaptation_factor = 1.0 - 0.0015 * (self.heart_rate_bpm - 72.0); + let sympathetic_shortening = 1.0 - 0.12 * self.autonomic_tone.max(0.0); + let ischemia_prolongation = 1.0 + 0.15 * (1.0 - ischemia_factor); + + let apd_target = + (280.0 * rate_adaptation_factor * sympathetic_shortening * ischemia_prolongation) + .clamp(200.0, 360.0); + self.ventricular_apd90_ms = + Self::approach(self.ventricular_apd90_ms, apd_target, 5.0, dt_seconds); + + // ERP typically 90-95% of APD90, plus relative refractory period + let erp_target = (self.ventricular_apd90_ms * 1.1 + 80.0).clamp(280.0, 500.0); + self.his_purkinje_erp_ms = + Self::approach(self.his_purkinje_erp_ms, erp_target, 4.0, dt_seconds); + + // Reset depolarization timer on each ventricular beat + let cycle_duration_ms = 60000.0 / self.heart_rate_bpm.max(30.0); + if self.time_since_depolarization_ms >= cycle_duration_ms { + self.time_since_depolarization_ms = 0.0; + } + } + + /// Update cardiac cycle phase state machine and valve states. + /// Based on: Sengupta PP et al. PMC2390899 and Wiggers diagram. + fn update_cardiac_phase(&mut self, dt_seconds: f32) { + self.phase_time_s += dt_seconds; + let cycle_duration_s = 60.0 / self.heart_rate_bpm.max(30.0); + + // Phase durations as fractions of cardiac cycle + // At 72 bpm (0.833s cycle): atrial systole ~0.1s, isovol contract ~0.05s, + // ejection ~0.3s, isovol relax ~0.08s, passive fill ~0.3s + let diastole_fraction = 0.6 - 0.15 * ((self.heart_rate_bpm - 72.0) / 60.0).clamp(-0.5, 1.0); + let systole_fraction = 1.0 - diastole_fraction; + + let atrial_systole_duration = 0.1 * cycle_duration_s; + let isovol_contract_duration = 0.05 * cycle_duration_s; + let ejection_duration = (systole_fraction - 0.05) * cycle_duration_s; + let isovol_relax_duration = 0.08 * cycle_duration_s; + + // State transitions + match self.cardiac_phase { + CardiacPhase::PassiveFilling => { + if self.phase_time_s >= (diastole_fraction - 0.1) * cycle_duration_s { + self.cardiac_phase = CardiacPhase::AtrialSystole; + self.phase_time_s = 0.0; + self.mitral_valve_state = ValveState::Open; + self.tricuspid_valve_state = ValveState::Open; + } + } + CardiacPhase::AtrialSystole => { + // Atrial contraction increases atrial and ventricular pressures + let atrial_kick_boost = + 1.0 + 0.2 * (self.phase_time_s / atrial_systole_duration).min(1.0); + self.left_atrial_pressure_mm_hg = 8.0 + 4.0 * atrial_kick_boost; + self.left_ventricular_pressure_mm_hg = Self::approach( + self.left_ventricular_pressure_mm_hg, + self.preload_mm_hg * 1.15, + 10.0, + dt_seconds, + ); + + if self.phase_time_s >= atrial_systole_duration { + self.cardiac_phase = CardiacPhase::IsovolumetricContraction; + self.phase_time_s = 0.0; + self.mitral_valve_state = ValveState::Closing; + self.tricuspid_valve_state = ValveState::Closing; + } + } + CardiacPhase::IsovolumetricContraction => { + // All valves closed, LV pressure rises rapidly + self.mitral_valve_state = ValveState::Closed; + self.tricuspid_valve_state = ValveState::Closed; + self.aortic_valve_state = ValveState::Closed; + self.pulmonic_valve_state = ValveState::Closed; + + let pressure_rise_rate = 1200.0 * self.contractility_index; + self.left_ventricular_pressure_mm_hg += pressure_rise_rate * dt_seconds; + + // Transition when LV pressure exceeds aortic pressure + if self.left_ventricular_pressure_mm_hg >= self.aortic_pressure_mm_hg + || self.phase_time_s >= isovol_contract_duration + { + self.cardiac_phase = CardiacPhase::Ejection; + self.phase_time_s = 0.0; + self.aortic_valve_state = ValveState::Opening; + self.pulmonic_valve_state = ValveState::Opening; + } + } + CardiacPhase::Ejection => { + self.aortic_valve_state = ValveState::Open; + self.pulmonic_valve_state = ValveState::Open; + + // LV pressure tracks aortic pressure during ejection + let systolic_peak = self.arterial_bp.systolic as f32; + self.left_ventricular_pressure_mm_hg = Self::approach( + self.left_ventricular_pressure_mm_hg, + systolic_peak * 0.98, + 80.0, + dt_seconds, + ); + self.aortic_pressure_mm_hg = + Self::approach(self.aortic_pressure_mm_hg, systolic_peak, 100.0, dt_seconds); + + if self.phase_time_s >= ejection_duration { + self.cardiac_phase = CardiacPhase::IsovolumetricRelaxation; + self.phase_time_s = 0.0; + self.aortic_valve_state = ValveState::Closing; + self.pulmonic_valve_state = ValveState::Closing; + } + } + CardiacPhase::IsovolumetricRelaxation => { + self.mitral_valve_state = ValveState::Closed; + self.tricuspid_valve_state = ValveState::Closed; + self.aortic_valve_state = ValveState::Closed; + self.pulmonic_valve_state = ValveState::Closed; + + // Rapid pressure drop with lusitropy + let relaxation_rate = -800.0; + self.left_ventricular_pressure_mm_hg += relaxation_rate * dt_seconds; + self.left_ventricular_pressure_mm_hg = + self.left_ventricular_pressure_mm_hg.max(0.0); + + // Aortic pressure decays toward diastolic + self.aortic_pressure_mm_hg = Self::approach( + self.aortic_pressure_mm_hg, + self.arterial_bp.diastolic as f32, + 150.0, + dt_seconds, + ); + + // Transition when LV pressure drops below atrial pressure + if self.left_ventricular_pressure_mm_hg <= self.left_atrial_pressure_mm_hg + || self.phase_time_s >= isovol_relax_duration + { + self.cardiac_phase = CardiacPhase::PassiveFilling; + self.phase_time_s = 0.0; + self.mitral_valve_state = ValveState::Opening; + self.tricuspid_valve_state = ValveState::Opening; + } + } + } + + // Update regurgitant valve states + if self.mitral_regurgitation_fraction > 0.05 + && self.mitral_valve_state == ValveState::Closed + { + self.mitral_valve_state = ValveState::Regurgitant; + } + if self.aortic_regurgitation_fraction > 0.05 + && self.aortic_valve_state == ValveState::Closed + { + self.aortic_valve_state = ValveState::Regurgitant; + } + + // Atrial filling from venous return + let filling_rate = self.venous_return_l_min * 1000.0 / 60.0; + self.left_atrial_volume_ml += filling_rate * dt_seconds; + self.right_atrial_volume_ml += filling_rate * dt_seconds * 0.95; + + // Atrial emptying during open mitral valve + if matches!( + self.mitral_valve_state, + ValveState::Open | ValveState::Opening + ) { + let emptying_rate = filling_rate * 1.8; + self.left_atrial_volume_ml -= emptying_rate * dt_seconds; + } + + // Clamp atrial volumes + self.left_atrial_volume_ml = self.left_atrial_volume_ml.clamp(30.0, 120.0); + self.right_atrial_volume_ml = self.right_atrial_volume_ml.clamp(28.0, 110.0); + + // Update atrial pressures from volumes (compliance ~5 mL/mmHg) + self.left_atrial_pressure_mm_hg = + (4.0 + (self.left_atrial_volume_ml - 50.0) / 5.0).clamp(2.0, 25.0); + self.right_atrial_pressure_mm_hg = + (2.0 + (self.right_atrial_volume_ml - 45.0) / 6.0).clamp(1.0, 20.0); + } + + /// Improve coronary perfusion model to emphasize diastolic flow. + /// Based on: Duncker DJ, Bache RJ. JACC 2021; coronary autoregulation research. + fn update_coronary_perfusion(&mut self) { + // CPP = Aortic Diastolic Pressure - LVEDP + // Most LV coronary flow occurs during diastole due to systolic compression + let diastolic_bp = self.arterial_bp.diastolic as f32; + let lvedp = self.preload_mm_hg; + self.coronary_perfusion_mm_hg = (diastolic_bp - lvedp).clamp(20.0, 140.0); + + // Diastolic time fraction decreases with tachycardia + let cycle_duration_s = 60.0 / self.heart_rate_bpm.max(30.0); + let diastolic_fraction = + (0.6 - 0.15 * ((self.heart_rate_bpm - 72.0) / 60.0)).clamp(0.3, 0.7); + let diastolic_time_s = cycle_duration_s * diastolic_fraction; + + // Coronary flow autoregulation maintains flow across CPP 60-180 mmHg + // Below 60 mmHg, flow becomes pressure-dependent + let autoregulation_factor = if self.coronary_perfusion_mm_hg < 60.0 { + self.coronary_perfusion_mm_hg / 60.0 + } else { + 1.0 + 0.1 * ((self.coronary_perfusion_mm_hg - 90.0) / 60.0).clamp(-0.5, 1.0) + }; + + // Metabolic vasodilation in response to increased demand + let demand_factor = (self.myocardial_oxygen_demand / 9.0).clamp(0.6, 1.8); + + // Tachycardia penalty: reduced diastolic time → reduced perfusion opportunity + let time_penalty = (diastolic_time_s / 0.5).clamp(0.6, 1.2); + + // Myocardial O2 supply now depends on CPP, autoregulation, and diastolic time + self.myocardial_oxygen_supply = (9.5 + 0.12 * self.coronary_perfusion_mm_hg + - 0.8 * (1.0 - autoregulation_factor) + + 1.5 * (autoregulation_factor - 1.0).max(0.0) + - 2.0 * (1.0 - time_penalty)) + .clamp(3.0, 22.0) + * demand_factor.min(1.2); + } } impl Organ for Heart { @@ -274,10 +684,22 @@ impl Organ for Heart { self.time_in_state_s += dt_seconds; + // Core hemodynamic control self.update_autonomic_state(dt_seconds); + self.update_raas(dt_seconds); // NEW: RAAS hormonal regulation self.determine_rhythm_state(); self.update_rate_and_contractility(dt_seconds); self.update_volumes_and_output(dt_seconds); + + // NEW: Phase-based cardiac cycle with discrete valves + self.update_cardiac_phase(dt_seconds); + + // NEW: Extended conduction system modeling + self.update_conduction_system(dt_seconds); + + // NEW: Improved coronary perfusion with diastolic emphasis + self.update_coronary_perfusion(); + self.update_arrhythmia_burden(dt_seconds); } fn summary(&self) -> String { diff --git a/src/organs/intestines.rs b/src/organs/intestines.rs index 5703b74..81c8eea 100644 --- a/src/organs/intestines.rs +++ b/src/organs/intestines.rs @@ -56,6 +56,12 @@ pub struct Intestines { pub nutrient_energy_kcal: f32, /// Fermentable fiber load (g). pub fiber_load_g: f32, + /// Iron absorbed into portal circulation per day (mg). + pub iron_absorption_mg_per_day: f32, + /// Folate absorbed per day (mg). + pub folate_absorption_mg_per_day: f32, + /// Vitamin B12 absorbed per day (mcg). + pub b12_absorption_mcg_per_day: f32, time_in_phase_s: f32, feeding_clock_s: f32, target_feed_interval_s: f32, @@ -86,6 +92,9 @@ impl Intestines { enteric_tone: 0.55, nutrient_energy_kcal: 40.0, fiber_load_g: 6.0, + iron_absorption_mg_per_day: 1.6, + folate_absorption_mg_per_day: 0.45, + b12_absorption_mcg_per_day: 5.2, time_in_phase_s: 0.0, feeding_clock_s: 0.0, target_feed_interval_s: 3.8 * 3600.0, @@ -237,6 +246,41 @@ impl Intestines { 0.2, dt_seconds, ); + + let micronutrient_drive = (available / 400.0).clamp(0.0, 2.0); + let mucosal_factor = self.mucosal_integrity.clamp(0.3, 1.1); + let inflammation_penalty = (self.inflammation_index * 0.5).clamp(0.0, 0.6); + + let iron_target = ((1.3 + 1.4 * micronutrient_drive) * mucosal_factor + - inflammation_penalty) + .clamp(0.3, 5.5); + self.iron_absorption_mg_per_day = Self::approach( + self.iron_absorption_mg_per_day, + iron_target, + 0.03, + dt_seconds, + ); + + let folate_target = ((0.32 + 0.3 * micronutrient_drive) * mucosal_factor + - inflammation_penalty * 0.25) + .clamp(0.1, 1.4); + self.folate_absorption_mg_per_day = Self::approach( + self.folate_absorption_mg_per_day, + folate_target, + 0.03, + dt_seconds, + ); + + let b12_target = ((4.2 + 3.4 * micronutrient_drive) * mucosal_factor + - inflammation_penalty * 3.0) + .clamp(0.5, 15.0); + self.b12_absorption_mcg_per_day = Self::approach( + self.b12_absorption_mcg_per_day, + b12_target, + 0.03, + dt_seconds, + ); + let water_removed = self.water_reabsorption_ml_min * dt_seconds / 60.0; self.lumen_volume_ml = (self.lumen_volume_ml + absorbed_kcal * 0.2 - water_removed).clamp(120.0, 800.0); diff --git a/src/organs/mod.rs b/src/organs/mod.rs index 0e5dbd7..e4c50ef 100644 --- a/src/organs/mod.rs +++ b/src/organs/mod.rs @@ -56,8 +56,8 @@ mod spinal_cord; mod spleen; mod stomach; -pub use bladder::{Bladder, BladderPhase}; -pub use bloodstream::{Bloodstream, MetabolicState, PerfusionState}; +pub use bladder::{Bladder, BladderMetrics, BladderPhase}; +pub use bloodstream::{Bloodstream, BloodstreamMetrics, MetabolicState, PerfusionState}; pub use brain::{Brain, SleepStage}; pub use esophagus::{EsophagealStage, Esophagus}; pub use gallbladder::Gallbladder; diff --git a/src/patient.rs b/src/patient.rs index 90fe164..b22be14 100644 --- a/src/patient.rs +++ b/src/patient.rs @@ -3,9 +3,9 @@ use crate::ekg::{EkgLead, EkgMonitor, EkgSnapshot, HeartElectricalState}; use crate::error::MedicalError; use crate::organs::{ - Bladder, BladderPhase, Bloodstream, Brain, EsophagealStage, Esophagus, Gallbladder, Heart, - Intestines, Kidneys, Liver, Lungs, MetabolicState, Organ, Pancreas, PerfusionState, SleepStage, - SpinalCord, Spleen, Stomach, + Bladder, BladderMetrics, BladderPhase, Bloodstream, BloodstreamMetrics, Brain, EsophagealStage, + Esophagus, Gallbladder, Heart, Intestines, Kidneys, Liver, Lungs, MetabolicState, Organ, + Pancreas, PerfusionState, SleepStage, SpinalCord, Spleen, Stomach, }; use crate::types::{Blood, BloodPressure, OrganType}; @@ -44,6 +44,9 @@ struct BrainSignals { consciousness: f32, sleep_depth: f32, rem_tone: f32, + cortical_arousal: f32, + cognitive_load: f32, + pain_drive: f32, } #[derive(Clone, Copy)] @@ -76,6 +79,9 @@ struct LiverSignals { bile_secretion_ml_min: f32, detox: u8, ammonia_clearance_umol_min: f32, + albumin_g_dl: f32, + clotting_factor_synthesis_pct: f32, + acute_phase_response: f32, } #[derive(Clone, Copy)] @@ -89,6 +95,9 @@ struct PancreasSignals { #[derive(Clone, Copy)] struct IntestineSignals { nutrient_energy_kcal: f32, + iron_absorption_mg_per_day: f32, + folate_absorption_mg_per_day: f32, + b12_absorption_mcg_per_day: f32, } #[derive(Clone, Copy)] @@ -101,6 +110,7 @@ struct SpinalSignals { sympathetic_outflow: f32, parasympathetic_outflow: f32, reflex_gain: f32, + nociceptive_facilitation: f32, } #[derive(Clone, Copy)] @@ -109,6 +119,7 @@ struct KidneySignals { plasma_volume_l: f32, urea_excretion_mg_min: f32, erythropoietin_iu_per_day: f32, + serum_osmolality_mosm: f32, } #[derive(Clone, Copy)] @@ -116,6 +127,7 @@ struct SpleenSignals { immune_activity: u8, red_pulp_volume_ml: f32, platelet_reservoir: f32, + erythrocyte_culling_rate: f32, } #[derive(Clone, Copy)] @@ -397,12 +409,18 @@ impl Patient { bile_secretion_ml_min: l.bile_secretion_ml_min, detox: l.detox, ammonia_clearance_umol_min: l.ammonia_clearance_umol_min, + albumin_g_dl: l.albumin_g_dl, + clotting_factor_synthesis_pct: l.clotting_factor_synthesis_pct, + acute_phase_response: l.acute_phase_response, }); let intestine_signals = self .find_organ_typed::() .map(|i| IntestineSignals { nutrient_energy_kcal: i.nutrient_energy_kcal, + iron_absorption_mg_per_day: i.iron_absorption_mg_per_day, + folate_absorption_mg_per_day: i.folate_absorption_mg_per_day, + b12_absorption_mcg_per_day: i.b12_absorption_mcg_per_day, }); let gallbladder_signals = @@ -426,6 +444,7 @@ impl Patient { plasma_volume_l: k.plasma_volume_l, urea_excretion_mg_min: k.urea_excretion_mg_min, erythropoietin_iu_per_day: k.erythropoietin_iu_per_day, + serum_osmolality_mosm: k.serum_osmolality_mosm, }); let spinal_signals = self @@ -434,12 +453,20 @@ impl Patient { sympathetic_outflow: s.sympathetic_outflow, parasympathetic_outflow: s.parasympathetic_outflow, reflex_gain: s.reflex_gain, + nociceptive_facilitation: s.nociceptive_facilitation, }); let spleen_signals = self.find_organ_typed::().map(|s| SpleenSignals { immune_activity: s.immune_activity, red_pulp_volume_ml: s.red_pulp_volume_ml, platelet_reservoir: s.platelet_reservoir, + erythrocyte_culling_rate: s.erythrocyte_culling_rate, + }); + + let bladder_signals = self.find_organ_typed::().map(|b| BladderSignals { + afferent_signal: b.afferent_signal, + urgency: b.urgency, + phase: b.phase, }); if let Some(bloodstream) = self.find_organ_typed_mut::() { @@ -465,13 +492,31 @@ impl Patient { } if let Some(liver) = liver_signals { bloodstream.set_hepatic_feedback(liver.detox, liver.ammonia_clearance_umol_min); + let immune_hint = spleen_signals + .as_ref() + .map(|s| s.immune_activity as f32) + .unwrap_or(80.0); + let globulin_target = + (2.4 + liver.acute_phase_response * 1.1 + immune_hint * 0.01).clamp(1.4, 5.2); + let fibrinogen = + (0.3 * liver.clotting_factor_synthesis_pct / 100.0).clamp(0.05, 0.9); + bloodstream.set_plasma_proteins(liver.albumin_g_dl, globulin_target, fibrinogen); } if let Some(intestines) = intestine_signals { bloodstream.set_metabolic_nutrients(intestines.nutrient_energy_kcal); + bloodstream.set_erythropoiesis_nutrients( + intestines.iron_absorption_mg_per_day, + intestines.folate_absorption_mg_per_day, + intestines.b12_absorption_mcg_per_day, + ); } if let Some(spleen) = spleen_signals { - bloodstream - .set_spleen_feedback(spleen.platelet_reservoir, spleen.red_pulp_volume_ml); + bloodstream.set_spleen_feedback( + spleen.platelet_reservoir, + spleen.red_pulp_volume_ml, + spleen.immune_activity, + spleen.erythrocyte_culling_rate, + ); } } @@ -500,6 +545,7 @@ impl Patient { }); if let Some(blood) = bloodstream_signals { + let mut brain_drive_hint: Option<(f32, f32)> = None; if let Some(heart) = self.find_organ_typed_mut::() { let preload_target = ((blood.total_volume_l - 5.0) * 6.0 + 8.0).clamp(4.0, 18.0); heart.preload_mm_hg = @@ -565,6 +611,80 @@ impl Patient { } if let Some(brain) = self.find_organ_typed_mut::() { + let oxygen_frac = blood.arterial_o2_saturation_pct / 100.0; + let hypoxia_drive = (0.94 - oxygen_frac).max(0.0) * 2.0; + let perfusion_drive = (1.0 - blood.oxygen_supply_demand_ratio).max(0.0); + let hypercapnia_drive = (blood.co2_content_ml_dl - 40.0).max(0.0) / 12.0; + let drive_target = + (0.45 + hypoxia_drive + perfusion_drive + hypercapnia_drive).clamp(0.1, 1.3); + brain.respiratory_drive = + Self::relax_value(brain.respiratory_drive, drive_target, dt_seconds, 8.0) + .clamp(0.05, 1.4); + + let hypocapnia = (36.0 - blood.co2_content_ml_dl).max(0.0) / 12.0; + let hyperoxia = (oxygen_frac - 0.99).max(0.0) * 1.2; + let syncope_target = (hypocapnia * (0.6 + hyperoxia)).clamp(0.0, 1.0); + brain.syncope_propensity = + Self::relax_value(brain.syncope_propensity, syncope_target, dt_seconds, 12.0) + .clamp(0.0, 1.0); + + let autonomic_drive_target = (0.6 + brain.respiratory_drive * 0.35 + - brain.syncope_propensity * 0.25) + .clamp(0.25, 1.1); + brain.brainstem_autonomic_drive = Self::relax_value( + brain.brainstem_autonomic_drive, + autonomic_drive_target, + dt_seconds, + 15.0, + ) + .clamp(0.25, 1.1); + + let dehydration = (3.0 - blood.plasma_volume_l).max(0.0) / 0.7; + let osm_drive = kidney_signals + .map(|k| (k.serum_osmolality_mosm - 285.0).max(0.0) / 15.0) + .unwrap_or(0.0); + let thirst_target = (0.25 + dehydration + osm_drive).clamp(0.0, 1.4); + brain.thirst_drive = + Self::relax_value(brain.thirst_drive, thirst_target, dt_seconds, 25.0) + .clamp(0.0, 1.4); + + let glucose_deficit = (95.0 - blood.glucose_mg_dl).max(0.0) / 55.0; + let stomach_need = stomach_signals + .map(|s| (220.0 - s.nutrient_load_kcal).max(0.0) / 320.0) + .unwrap_or(0.15); + let intestine_need = intestine_signals + .map(|i| (160.0 - i.nutrient_energy_kcal).max(0.0) / 260.0) + .unwrap_or(0.1); + let hunger_target = + (0.35 + glucose_deficit + stomach_need + intestine_need).clamp(0.0, 1.3); + brain.hunger_drive = + Self::relax_value(brain.hunger_drive, hunger_target, dt_seconds, 30.0) + .clamp(0.0, 1.3); + + let temp_delta = blood.temperature_c - 37.0; + let thermo_target = + (0.5 + temp_delta.abs() * 0.4 + temp_delta.max(0.0) * 0.2).clamp(0.0, 1.2); + brain.thermoregulatory_drive = Self::relax_value( + brain.thermoregulatory_drive, + thermo_target, + dt_seconds, + 28.0, + ) + .clamp(0.0, 1.2); + + let spinal_nociception = spinal_signals + .map(|s| s.nociceptive_facilitation) + .unwrap_or(0.2); + let bladder_discomfort = bladder_signals.map(|b| b.urgency * 0.4).unwrap_or(0.0); + let pain_target = (0.2 + + spinal_nociception * 0.6 + + bladder_discomfort + + (blood.waste_load_mg - 2000.0).max(0.0) / 2500.0) + .clamp(0.0, 1.5); + brain.pain_drive = + Self::relax_value(brain.pain_drive, pain_target, dt_seconds, 18.0) + .clamp(0.0, 1.5); + let cbf_target = (52.0 * blood.oxygen_supply_demand_ratio.clamp(0.6, 1.4)).clamp(28.0, 75.0); brain.cerebral_blood_flow_ml_per_100g_min = Self::relax_value( @@ -608,6 +728,22 @@ impl Patient { dt_seconds, 20.0, ); + brain_drive_hint = Some((brain.respiratory_drive, brain.syncope_propensity)); + } + + if let Some((drive, syncope)) = brain_drive_hint { + if let Some(lungs) = self.find_organ_typed_mut::() { + let muscle_target = (0.3 + drive * 0.55 - syncope * 0.4).clamp(0.1, 1.0); + lungs.muscle_drive = + Self::relax_value(lungs.muscle_drive, muscle_target, dt_seconds, 35.0); + let chemo_target = (0.35 + drive * 0.6).clamp(0.2, 1.0); + lungs.chemoreceptor_drive = Self::relax_value( + lungs.chemoreceptor_drive, + chemo_target, + dt_seconds, + 28.0, + ); + } } if let Some(spinal) = self.find_organ_typed_mut::() { @@ -1056,12 +1192,17 @@ impl Patient { immune_activity: s.immune_activity, red_pulp_volume_ml: s.red_pulp_volume_ml, platelet_reservoir: s.platelet_reservoir, + erythrocyte_culling_rate: s.erythrocyte_culling_rate, }); if let Some(spleen) = spleen_after { if let Some(bloodstream) = self.find_organ_typed_mut::() { - bloodstream - .set_spleen_feedback(spleen.platelet_reservoir, spleen.red_pulp_volume_ml); + bloodstream.set_spleen_feedback( + spleen.platelet_reservoir, + spleen.red_pulp_volume_ml, + spleen.immune_activity, + spleen.erythrocyte_culling_rate, + ); } else { let immune_penalty = (spleen.immune_activity as f32 - 80.0) / 120.0; let hematocrit_target: f32 = 42.0 @@ -1103,6 +1244,9 @@ impl Patient { SleepStage::Rem => 0.4, }, rem_tone: matches!(b.sleep_stage, SleepStage::Rem) as i32 as f32, + cortical_arousal: b.cortical_arousal, + cognitive_load: b.cognitive_load, + pain_drive: b.pain_drive, }); let spinal_after = self .find_organ_typed::() @@ -1110,6 +1254,7 @@ impl Patient { sympathetic_outflow: s.sympathetic_outflow, parasympathetic_outflow: s.parasympathetic_outflow, reflex_gain: s.reflex_gain, + nociceptive_facilitation: s.nociceptive_facilitation, }); let bladder_after = self.find_organ_typed::().map(|b| BladderSignals { afferent_signal: b.afferent_signal, @@ -1210,25 +1355,37 @@ impl Patient { if let Some(bladder) = self.find_organ_typed_mut::() { if let Some(brain) = brain_after { let sedation = (1.0 - brain.consciousness).clamp(0.0, 1.0); - let sleep_bonus = brain.sleep_depth * 0.25 + brain.rem_tone * 0.1; - let capacity = bladder.capacity_ml.max(1.0); - let micturition_target = (capacity - * (0.65 + sleep_bonus + sedation * 0.2 - bladder_state.urgency * 0.05)) - .clamp(capacity * 0.5, capacity * 0.9); - bladder.micturition_threshold_ml = Self::relax_value( - bladder.micturition_threshold_ml, - micturition_target, + let rem_suppression = brain.rem_tone; + let cortical_target = (brain.cortical_arousal + * (1.0 - rem_suppression * 0.85) + * (1.0 - sedation * 0.7) + * (1.0 - brain.autonomic_variability * 0.2) + + bladder.continence_training * 0.2) + .clamp(0.0, 1.0); + bladder.cortical_inhibition = Self::relax_value( + bladder.cortical_inhibition, + cortical_target, dt_seconds, - 900.0, - ); - let urge_target = (capacity * (0.4 + sleep_bonus * 0.5 + sedation * 0.15)) - .clamp(capacity * 0.3, capacity * 0.7); - bladder.urge_threshold_ml = Self::relax_value( - bladder.urge_threshold_ml, - urge_target, + 120.0, + ) + .clamp(0.0, 1.0); + + let hold_driver = (brain.consciousness + * (1.0 - rem_suppression * 0.9) + * (0.5 + brain.cognitive_load * 0.4) + - brain.pain_drive * 0.18) + .clamp(0.0, 1.0); + let urgency_driver = + bladder_state.urgency * (0.5 + bladder.continence_training * 0.3); + let voluntary_target = + ((hold_driver + urgency_driver) * (1.0 - sedation * 0.65)).clamp(0.0, 1.0); + bladder.voluntary_hold_command = Self::relax_value( + bladder.voluntary_hold_command, + voluntary_target, dt_seconds, - 600.0, - ); + 90.0, + ) + .clamp(0.0, 1.0); } if let Some(brain) = brain_after { @@ -1239,12 +1396,15 @@ impl Patient { .map(|s| s.parasympathetic_outflow) .unwrap_or(0.5); let reflex_gain = spinal_after.map(|s| s.reflex_gain).unwrap_or(1.0); - let para_target = (0.3 - + brain.brainstem_drive * 0.4 - + spinal_para * 0.3 - + bladder_state.urgency * 0.5 - + voiding_signal * 0.3) - .clamp(0.05, 1.1); + let cortical_brake = (1.0 + - bladder.cortical_inhibition * bladder.neurologic_integrity * 0.5) + .clamp(0.35, 1.0); + let para_target = (0.28 + + brain.brainstem_drive * 0.35 + + spinal_para * 0.28 + + bladder_state.urgency * 0.55 * cortical_brake + + voiding_signal * 0.35) + .clamp(0.05, 1.05); bladder.parasympathetic_drive = Self::relax_value( bladder.parasympathetic_drive, para_target, @@ -1252,20 +1412,23 @@ impl Patient { 120.0, ) .clamp(0.0, 1.0); - let sym_target = - (0.55 + (1.0 - brain.brainstem_drive) * 0.35 + spinal_sym * 0.4 - - bladder_state.urgency * 0.4 - - voiding_signal * 0.4) - .clamp(0.1, 1.0); + let sym_target = (0.6 + + (1.0 - brain.brainstem_drive) * 0.3 + + spinal_sym * 0.35 + + bladder.cortical_inhibition * 0.2 + - bladder_state.urgency * 0.45 + - voiding_signal * 0.45) + .clamp(0.1, 1.0); bladder.sympathetic_drive = Self::relax_value(bladder.sympathetic_drive, sym_target, dt_seconds, 160.0) .clamp(0.0, 1.0); - let voluntary = brain.consciousness; - let somatic_target = ((0.55 + voluntary * 0.35 - brain.sleep_depth * 0.3) - * (1.0 - bladder_state.urgency * 0.55) - + reflex_gain * 0.1 - - voiding_signal * 0.5) - .clamp(0.1, 0.95); + let somatic_base = (0.5 + brain.consciousness * 0.3 + - brain.sleep_depth * 0.3 + - brain.rem_tone * 0.4 + + reflex_gain * 0.1) + * (1.0 - bladder_state.urgency * 0.5); + let somatic_target = + (somatic_base + bladder.voluntary_hold_command * 0.5).clamp(0.1, 1.0); bladder.somatic_drive = Self::relax_value(bladder.somatic_drive, somatic_target, dt_seconds, 110.0) .clamp(0.0, 1.0); @@ -1338,6 +1501,15 @@ impl Patient { } /// One-line summary of a specific organ by type. + /// Snapshot of bladder gating and continence metrics for downstream consumers. + pub fn bloodstream_metrics(&self) -> Option { + self.find_organ_typed::().map(|b| b.metrics()) + } + + pub fn bladder_metrics(&self) -> Option { + self.find_organ_typed::().map(|b| b.metrics()) + } + pub fn organ_summary(&self, kind: OrganType) -> crate::Result { match self.find_organ(kind) { Some(o) => Ok(o.summary()), @@ -1430,10 +1602,14 @@ mod tests { .find_organ_typed::() .unwrap() .preload_mm_hg; - let initial_consciousness = patient - .find_organ_typed::() - .unwrap() - .consciousness; + let (initial_consciousness, initial_drive, _initial_syncope) = { + let brain = patient.find_organ_typed::().unwrap(); + ( + brain.consciousness, + brain.respiratory_drive, + brain.syncope_propensity, + ) + }; let initial_kupffer = patient .find_organ_typed::() .unwrap() @@ -1463,8 +1639,13 @@ mod tests { let brain = patient.find_organ_typed::().unwrap(); assert!( - brain.consciousness as f32 <= initial_consciousness as f32, - "consciousness failed to drop: initial {initial_consciousness}, after {}", + brain.respiratory_drive >= initial_drive + 0.02, + "respiratory drive did not rise: initial {initial_drive}, after {}", + brain.respiratory_drive + ); + assert!( + brain.consciousness as f32 <= initial_consciousness as f32 + 3.0, + "consciousness unexpectedly increased: initial {initial_consciousness}, after {}", brain.consciousness ); @@ -1476,6 +1657,196 @@ mod tests { ); } + #[test] + fn hypoxia_raises_respiratory_drive() { + let mut patient = Patient::new("drive-hypoxia") + .unwrap() + .initialize_default() + .with_lungs() + .with_organ(OrganType::Brain); + for _ in 0..6 { + patient.update(0.5); + } + let baseline_drive = patient + .find_organ_typed::() + .unwrap() + .respiratory_drive; + let baseline_vent = patient + .find_organ_typed::() + .unwrap() + .minute_ventilation_l_min; + + for _ in 0..6 { + { + let blood = patient.find_organ_typed_mut::().unwrap(); + blood.arterial_o2_saturation_pct = 82.0; + blood.oxygen_supply_demand_ratio = 0.6; + blood.oxygen_delivery_ml_min = 260.0; + blood.oxygen_consumption_ml_min = 340.0; + blood.co2_content_ml_dl = 58.0; + } + patient.update(0.5); + } + + let brain = patient.find_organ_typed::().unwrap(); + let lungs = patient.find_organ_typed::().unwrap(); + + assert!( + brain.respiratory_drive > baseline_drive + 0.08, + "drive failed to rise (baseline {}, after {})", + baseline_drive, + brain.respiratory_drive + ); + assert!( + lungs.minute_ventilation_l_min > baseline_vent, + "ventilation did not increase (baseline {}, after {})", + baseline_vent, + lungs.minute_ventilation_l_min + ); + } + + #[test] + fn hyperventilation_triggers_syncope_drive() { + let mut patient = Patient::new("drive-hypervent") + .unwrap() + .initialize_default() + .with_lungs() + .with_organ(OrganType::Brain); + for _ in 0..6 { + patient.update(0.5); + } + let baseline_syncope = patient + .find_organ_typed::() + .unwrap() + .syncope_propensity; + let baseline_consciousness = patient + .find_organ_typed::() + .unwrap() + .consciousness; + for _ in 0..6 { + { + let blood = patient.find_organ_typed_mut::().unwrap(); + blood.co2_content_ml_dl = 28.0; + blood.arterial_o2_saturation_pct = 100.0; + blood.oxygen_supply_demand_ratio = 1.4; + blood.oxygen_delivery_ml_min = 1100.0; + blood.oxygen_consumption_ml_min = 420.0; + } + patient.update(0.5); + } + let brain = patient.find_organ_typed::().unwrap(); + assert!( + brain.syncope_propensity > baseline_syncope + 0.05, + "syncope propensity failed to rise (baseline {}, after {})", + baseline_syncope, + brain.syncope_propensity + ); + assert!( + brain.consciousness as f32 <= baseline_consciousness as f32 - 2.0, + "consciousness did not fall (baseline {}, after {})", + baseline_consciousness, + brain.consciousness + ); + } + + #[test] + fn dehydration_and_hypoglycemia_raise_homeostatic_drives() { + let mut patient = Patient::new("drive-homeo") + .unwrap() + .initialize_default() + .with_lungs() + .with_organ(OrganType::Brain) + .with_organ(OrganType::Kidneys) + .with_organ(OrganType::Stomach) + .with_organ(OrganType::Intestines); + for _ in 0..6 { + patient.update(0.5); + } + let brain = patient + .find_organ_typed::() + .unwrap() + .clone(); + let baseline_thirst = brain.thirst_drive; + let baseline_hunger = brain.hunger_drive; + let baseline_thermo = brain.thermoregulatory_drive; + + for _ in 0..6 { + { + let blood = patient.find_organ_typed_mut::().unwrap(); + blood.plasma_volume_l = 2.5; + blood.total_volume_l = 3.8; + blood.temperature_c = 38.6; + blood.glucose_mg_dl = 68.0; + } + if let Some(kidneys) = patient.find_organ_typed_mut::() { + kidneys.serum_osmolality_mosm = 308.0; + } + if let Some(stomach) = patient.find_organ_typed_mut::() { + stomach.nutrient_load_kcal = 40.0; + } + if let Some(intestines) = patient.find_organ_typed_mut::() { + intestines.nutrient_energy_kcal = 60.0; + } + patient.update(0.5); + } + + let brain = patient.find_organ_typed::().unwrap(); + assert!( + brain.thirst_drive > baseline_thirst + 0.1, + "thirst drive did not rise: baseline {}, after {}", + baseline_thirst, + brain.thirst_drive + ); + assert!( + brain.hunger_drive > baseline_hunger + 0.08, + "hunger drive did not rise: baseline {}, after {}", + baseline_hunger, + brain.hunger_drive + ); + assert!( + brain.thermoregulatory_drive > baseline_thermo + 0.06, + "thermoregulatory drive did not rise: baseline {}, after {}", + baseline_thermo, + brain.thermoregulatory_drive + ); + } + + #[test] + fn nociception_and_visceral_urgency_raise_pain_drive() { + let mut patient = Patient::new("drive-pain") + .unwrap() + .initialize_default() + .with_lungs() + .with_organ(OrganType::Brain) + .with_organ(OrganType::SpinalCord) + .with_organ(OrganType::Bladder); + for _ in 0..6 { + patient.update(0.5); + } + let baseline_pain = patient + .find_organ_typed::() + .unwrap() + .pain_drive; + + for _ in 0..6 { + if let Some(spinal) = patient.find_organ_typed_mut::() { + spinal.nociceptive_facilitation = 0.85; + } + if let Some(bladder) = patient.find_organ_typed_mut::() { + bladder.urgency = 0.8; + } + patient.update(0.5); + } + + let brain = patient.find_organ_typed::().unwrap(); + assert!( + brain.pain_drive > baseline_pain + 0.1, + "pain drive did not rise: baseline {}, after {}", + baseline_pain, + brain.pain_drive + ); + } + #[test] fn multi_organ_homeostasis_stability() { use crate::organs::{ @@ -1739,9 +2110,9 @@ mod tests { .assert_within(100.0, 300.0, "red pulp volume"); self.bladder_volume - .assert_within(20.0, 620.0, "bladder volume"); + .assert_within(20.0, 720.0, "bladder volume"); self.bladder_pressure - .assert_within(0.0, 80.0, "bladder pressure"); + .assert_within(0.0, 110.0, "bladder pressure"); self.esophagus_acid .assert_within(0.0, 0.6, "esophageal acid exposure"); diff --git a/tests/bladder_metrics.rs b/tests/bladder_metrics.rs new file mode 100644 index 0000000..7cb817d --- /dev/null +++ b/tests/bladder_metrics.rs @@ -0,0 +1,70 @@ +use medicallib_rust::{Bladder, Organ, OrganType, Patient}; + +#[test] +fn bladder_metrics_match_patient_api() { + let patient = Patient::new("metrics") + .unwrap() + .initialize_default() + .with_organ(OrganType::Bladder); + let metrics_via_patient = patient.bladder_metrics().expect("bladder present"); + let bladder = patient.find_organ_typed::().unwrap(); + let direct = bladder.metrics(); + assert!((metrics_via_patient.urgency - direct.urgency).abs() < 1e-6); + assert_eq!(metrics_via_patient.phase, direct.phase); + assert_eq!(metrics_via_patient.capacity_ml, direct.capacity_ml); +} + +#[test] +fn voluntary_hold_updates_metrics() { + let mut patient = Patient::new("hold") + .unwrap() + .initialize_default() + .with_organ(OrganType::Bladder); + + { + let bladder = patient.find_organ_typed_mut::().unwrap(); + bladder.volume_ml = bladder.micturition_threshold_ml + 140.0; + bladder.urgency = 0.82; + bladder.cortical_inhibition = 0.1; + bladder.voluntary_hold_command = 0.1; + bladder.continence_training = 0.45; + bladder.update(0.5); + } + let baseline = patient.bladder_metrics().unwrap(); + + { + let bladder = patient.find_organ_typed_mut::().unwrap(); + bladder.cortical_inhibition = 0.9; + bladder.voluntary_hold_command = 0.9; + bladder.continence_training = 0.9; + bladder.update(0.5); + } + let hold = patient.bladder_metrics().unwrap(); + + assert!(hold.voluntary_hold_fraction > baseline.voluntary_hold_fraction); + assert!(hold.cortical_gate_fraction > baseline.cortical_gate_fraction); + assert!(hold.cortical_gate_fraction <= 1.0); + assert!(hold.voluntary_hold_fraction <= 1.0); +} + +#[cfg(feature = "ffi")] +#[test] +fn ffi_exposes_bladder_metrics() { + use core::mem::MaybeUninit; + use medicallib_rust::ffi::{ + ml_patient_bladder_metrics, ml_patient_free, ml_patient_new, ml_patient_update, + MLBladderMetrics, MLPatient, ML_ERR, ML_OK, + }; + use std::ffi::CString; + + unsafe { + let id = CString::new("ffi-metrics").unwrap(); + let patient: *mut MLPatient = ml_patient_new(id.as_ptr()); + assert!(!patient.is_null()); + assert_eq!(ml_patient_update(patient, 0.5), ML_OK); + let mut metrics = MaybeUninit::::zeroed(); + let status = ml_patient_bladder_metrics(patient as *const MLPatient, metrics.as_mut_ptr()); + assert_eq!(status, ML_ERR); + ml_patient_free(patient); + } +} diff --git a/tests/ffi.rs b/tests/ffi.rs index 70fec2a..e232ba3 100644 --- a/tests/ffi.rs +++ b/tests/ffi.rs @@ -29,3 +29,24 @@ fn ffi_errors() { let s = medicallib_rust::ffi::ml_patient_summary(std::ptr::null()); assert!(s.is_null()); } + +#[cfg(feature = "ffi")] +#[test] +fn ffi_bloodstream_metrics() { + use medicallib_rust::ffi::*; + use std::ffi::CString; + + let id = CString::new("ffi-blood").unwrap(); + let p = ml_patient_new(id.as_ptr()); + assert!(!p.is_null()); + + let mut metrics = std::mem::MaybeUninit::::uninit(); + assert_eq!(ml_patient_update(p, 0.5), ML_OK); + let rc = ml_patient_bloodstream_metrics(p, metrics.as_mut_ptr()); + assert_eq!(rc, ML_OK); + let metrics = unsafe { metrics.assume_init() }; + assert!(metrics.oncotic_pressure_mm_hg > 15.0); + assert!(metrics.rbc_mature_fraction > 0.4); + ml_patient_free(p); +} +