diff --git a/bachelorarbeit/src/calc.rs b/bachelorarbeit/src/calc.rs index 08d9e53..fdddbef 100644 --- a/bachelorarbeit/src/calc.rs +++ b/bachelorarbeit/src/calc.rs @@ -4,8 +4,6 @@ use std::collections::HashMap; use std::io; use std::rc::Rc; -use ahash::{AHasher, RandomState}; - use std::{fs::File, io::BufReader}; use crate::gq_storage; @@ -217,7 +215,6 @@ pub fn omnes(cache: &mut Cache, s: f64, n: u32, use_tan: bool) -> Complex { return *res; } } - if !cache.omnes_lut.contains_key(&n) { cache .omnes_lut @@ -275,6 +272,10 @@ pub fn phi0_integrand( use_tan_phi0: bool, use_xsub: bool, _use_reduced_integrand: bool, + use_cutoff: bool, + cutoff_power: f64, + cutoff_s: f64, + cutoff_factor: Complex, ) -> f64 { let sub = match (use_tan_phi0, use_xsub) { (true, true) => s_tick.tan().powi(2) + S0, @@ -282,8 +283,18 @@ pub fn phi0_integrand( (false, true) => s_tick.powi(2) + S0, (false, false) => s_tick, }; + // let omnes_s = if use_cutoff && s >= cutoff_s { + // cutoff_factor * s.powf(cutoff_power) + // } else { + // omnes(cache, s, n_omnes, use_tan_omnes) + // }; + let omnes_sub = if use_cutoff && sub >= cutoff_s { + cutoff_factor * sub.powf(cutoff_power) + // Complex::new((0.00002f64).powf(cutoff_power), 0.0) + } else { + omnes(cache, sub, n_omnes, use_tan_omnes) + }; let omnes_s = omnes(cache, s, n_omnes, use_tan_omnes); - let omnes_sub = omnes(cache, sub, n_omnes, use_tan_omnes); let mut res = (omnes_sub / omnes_s).norm().ln() / (sub * (sub - s)); res *= match (use_tan_phi0, use_xsub) { (true, true) => 2.0 / s_tick.cos().powi(2), @@ -304,6 +315,9 @@ pub fn phi0( use_xsub: bool, use_reduced_integrand: bool, print_prog: bool, + use_cutoff: bool, + cutoff_power: f64, + cutoff_s: f64 ) -> f64 { if !cache.gauss_quad_lut.contains_key(&n_phi0) { let gq_values = GaussLegendre::nodes_and_weights(n_phi0.try_into().unwrap()); @@ -323,6 +337,9 @@ pub fn phi0( .ln() / (s * (s - S0).powf(1.5)) }; + + let cutoff_factor = omnes(cache, cutoff_s, n_omnes, use_tan_omnes) * cutoff_s.powf(-cutoff_power); + let intgrl = integrate( &*roots, &*weights, @@ -344,6 +361,10 @@ pub fn phi0( use_tan_phi0, use_xsub, use_reduced_integrand, + use_cutoff, + cutoff_power, + cutoff_s, + cutoff_factor, ) * s / PI }, @@ -352,8 +373,6 @@ pub fn phi0( -(s - S0).sqrt() * (intgrl - analyt) } - - #[allow(unused)] pub fn bit_pattern_randomness_tester(a: f64, b: f64, n_points: u64) -> Vec { let mut modulos = vec![0.0; n_points as usize]; diff --git a/bachelorarbeit/src/gq_storage.rs b/bachelorarbeit/src/gq_storage.rs index 75c8182..40ebafa 100644 --- a/bachelorarbeit/src/gq_storage.rs +++ b/bachelorarbeit/src/gq_storage.rs @@ -7,7 +7,7 @@ use std::{ fs::File, - io::{self, BufReader, BufWriter, Read, Write, Seek, SeekFrom}, + io::{self, BufReader, BufWriter, Read, Seek, SeekFrom, Write}, }; use gauss_quad::GaussLegendre; @@ -28,7 +28,10 @@ pub fn serialize(filepath: &str, params: Vec<(Vec, Vec)>) -> io::Resul Ok(()) } -pub fn deserialize(reader: &mut R, filter: impl Fn(u64) -> bool) -> io::Result, Vec)>> { +pub fn deserialize( + reader: &mut R, + filter: impl Fn(u64) -> bool, +) -> io::Result, Vec)>> { let mut buf = [0u8; 8]; let mut values = Vec::new(); @@ -41,10 +44,12 @@ pub fn deserialize(reader: &mut R, filter: impl Fn(u64) -> bool) // skips stored floats for n that should not be stored if !filter(n) { - reader.seek(SeekFrom::Current(i64::try_from(2*n*8).unwrap())).unwrap(); + reader + .seek(SeekFrom::Current(i64::try_from(2 * n * 8).unwrap())) + .unwrap(); continue; } - + // Stores roots and weights into values Vec let mut roots = Vec::new(); @@ -81,4 +86,4 @@ pub fn fill_with_gauss_quad(n_values: Vec) -> std::io::Result<()> { let res = deserialize(&mut file, |_| true); println!("{:?}", res); Ok(()) -} \ No newline at end of file +} diff --git a/bachelorarbeit/src/main.rs b/bachelorarbeit/src/main.rs index 22c61e6..6612609 100644 --- a/bachelorarbeit/src/main.rs +++ b/bachelorarbeit/src/main.rs @@ -6,14 +6,14 @@ use indicatif; use std::{collections::HashMap, time::Instant}; use egui_plot::{log_grid_spacer, AxisHints, GridInput, GridMark, Legend, Line, Plot}; -use num::complex::ComplexFloat; +use num::{complex::ComplexFloat, Float}; mod calc; mod gq_storage; struct App { plot_data: Vec, - plots_available: [String; 4], + plots_available: [String; 7], calc_cache: calc::Cache, scaling_type: ScalingType, a: f64, @@ -26,8 +26,8 @@ struct App { use_xsub_phi0: bool, use_reduced_integrand: bool, use_cutoff: bool, + cutoff_power: f64, cutoff_s: f64, - cont_power: f64, } #[derive(Clone)] @@ -45,6 +45,9 @@ impl Default for App { "Omnes".to_string(), "Phi0".to_string(), "Phi0 / delta".to_string(), + "s^cutoff_power".to_string(), + "phi0 integrand".to_string(), + "Phi0_cut / Phi0".to_string(), ], calc_cache: calc::Cache::default(), scaling_type: ScalingType::default(), @@ -62,8 +65,8 @@ impl Default for App { // Cutoff params use_cutoff: false, + cutoff_power: -1.0, cutoff_s: 100.0, - cont_power: -1.0, } } } @@ -82,42 +85,47 @@ impl eframe::App for App { ui.vertical(|ui| { ui.horizontal(|ui| { - egui::Grid::new("grid1") - .min_col_width(125.0) - .show(ui, |ui| { - ui.label("calc interval start:").on_hover_text( - "the lower bound of the interval that functions are calculated in", - ); - float_text_edit_singleline(ui, &mut self.a, true); - ui.end_row(); + ui.vertical(|ui| { + ui.label(egui::RichText::new("Boundary values")); + ui.add_space(10.0); + egui::Grid::new("grid1") + .min_col_width(125.0) + .show(ui, |ui| { + ui.label("calc interval start:").on_hover_text( + "the lower bound of the interval that functions are calculated in", + ); + float_text_edit_singleline(ui, &mut self.a, true); + ui.end_row(); - ui.label("calc interval end:").on_hover_text( - "the upper bound of the interval that functions are calculated in", - ); - float_text_edit_singleline(ui, &mut self.b, true); - ui.end_row(); + ui.label("calc interval end:").on_hover_text( + "the upper bound of the interval that functions are calculated in", + ); + float_text_edit_singleline(ui, &mut self.b, true); + ui.end_row(); - ui.label("N for omnes:").on_hover_text( - "the number of points used in calculation of omnes function", - ); - int_text_edit_singleline(ui, &mut self.n_omnes, true); - ui.end_row(); + ui.label("N for omnes:").on_hover_text( + "the number of points used in calculation of omnes function", + ); + int_text_edit_singleline(ui, &mut self.n_omnes, true); + ui.end_row(); - ui.label("N for phi0:").on_hover_text( - "the number of points used in calculation of phi0 function", - ); - int_text_edit_singleline(ui, &mut self.n_phi0, true); - ui.end_row(); - - ui.label("number of plotpoints:") - .on_hover_text("the number of points shown in the plot"); - int_text_edit_singleline(ui, &mut self.n_points, true); - ui.end_row(); - }); + ui.label("N for phi0:").on_hover_text( + "the number of points used in calculation of phi0 function", + ); + int_text_edit_singleline(ui, &mut self.n_phi0, true); + ui.end_row(); + ui.label("number of plotpoints:") + .on_hover_text("the number of points shown in the plot"); + int_text_edit_singleline(ui, &mut self.n_points, true); + ui.end_row(); + }); + }); ui.separator(); ui.vertical(|ui| { + ui.label(egui::RichText::new("Use substitutions")); + ui.add_space(10.0); let omnes_tan_cb = ui.checkbox(&mut self.use_tan_omnes, "Use tan()-subst. for Omnes ('ot')"); if omnes_tan_cb.clicked() { @@ -155,18 +163,19 @@ impl eframe::App for App { ui.separator(); - ui.vertical(|ui| { - ui.label("Omnes Continuation with polynomial").on_hover_text( + ui.with_layout(egui::Layout::top_down(egui::Align::Min),|ui| { + ui.label(egui::RichText::new("Omnes Continuation with polynomial")).on_hover_text( "If activated, from a user chosen cutoff value, the Omnes function will be replaced with s to a user chosen power, starting at the last Omnes value of the numerically calculated values", ); - ui.checkbox(&mut self.use_cutoff, "Use cutoff"); + ui.checkbox(&mut self.use_cutoff, "Use cutoff for phi0"); + ui.add_space(10.0); egui::Grid::new("grid2").min_col_width(125.0).show(ui, |ui| { - ui.label("Power for continuation").on_hover_text("Omnes will be replaced with s^{this power}."); - float_text_edit_singleline(ui, &mut self.cont_power, self.use_cutoff); + ui.label("Power of cont. polynomial").on_hover_text("Omnes will be replaced with s^{this power}."); + float_text_edit_singleline(ui, &mut self.cutoff_power, true); ui.end_row(); ui.label("Cutoff point").on_hover_text("Omnes will be replaced from this s-value onward."); - float_text_edit_singleline(ui, &mut self.cutoff_s, self.use_cutoff); + float_text_edit_singleline(ui, &mut self.cutoff_s, true); ui.end_row(); }); }) @@ -268,49 +277,41 @@ impl eframe::App for App { }) } 2 => { - let y_values_phi0: Vec = - if self.scaling_type == ScalingType::Linear { - let bar = - indicatif::ProgressBar::new(u64::from(self.n_points)); - let t0 = Instant::now(); - let values = x_values - .iter() - .map(|&x| { - bar.inc(1); - calc::phi0( - &mut self.calc_cache, - x, - self.n_omnes, - self.n_phi0, - self.use_tan_omnes, - self.use_tan_phi0, - self.use_xsub_phi0, - self.use_reduced_integrand, - x == self.a, - ) - }) - .collect(); - println!("{:?}", t0.elapsed()); - values - } else { - x_values - .iter() - .map(|&x| { - calc::phi0( - &mut self.calc_cache, - x, - self.n_omnes, - self.n_phi0, - self.use_tan_omnes, - self.use_tan_phi0, - self.use_xsub_phi0, - self.use_reduced_integrand, - x == self.a, - ) - .log10() - }) - .collect() - }; + let y_values_phi0: Vec = { + let bar = + indicatif::ProgressBar::new(u64::from(self.n_points)); + let t0 = Instant::now(); + let values = x_values + .iter() + .map(|&x| { + bar.inc(1); + let val = calc::phi0( + &mut self.calc_cache, + x, + self.n_omnes, + self.n_phi0, + self.use_tan_omnes, + self.use_tan_phi0, + self.use_xsub_phi0, + self.use_reduced_integrand, + x == self.a, + self.use_cutoff, + self.cutoff_power, + self.cutoff_s, + ); + match self.scaling_type { + ScalingType::Linear => { + val + } + ScalingType::Logarithmic => { + val.log10() + } + } + } + ).collect(); + println!("{:?}", t0.elapsed()); + values + }; self.plot_data.push(PlotCurve { points: x_values .iter() @@ -326,52 +327,40 @@ impl eframe::App for App { self.n_omnes, self.n_phi0 ), - }) + }); } 3 => { - let y_values_phi0: Vec = - if self.scaling_type == ScalingType::Linear { - x_values - .iter() - .map(|&x| { - let val = calc::phi0( - &mut self.calc_cache, - x, - self.n_omnes, - self.n_phi0, - self.use_tan_omnes, - self.use_tan_phi0, - self.use_xsub_phi0, - self.use_reduced_integrand, - x == self.a, - ) / calc::delta_with_lut( - &mut self.calc_cache, - x, - ); - val - }) - .collect() - } else { - x_values - .iter() - .map(|&x| { - let val = calc::phi0( - &mut self.calc_cache, - x, - self.n_omnes, - self.n_phi0, - self.use_tan_omnes, - self.use_tan_phi0, - self.use_xsub_phi0, - self.use_reduced_integrand, - x == self.a, - ) / calc::delta_with_lut( - &mut self.calc_cache, - x, - ); - val.log10() - }) - .collect() + let y_values_phi0: Vec = { + x_values + .iter() + .map(|&x| { + let val = calc::phi0( + &mut self.calc_cache, + x, + self.n_omnes, + self.n_phi0, + self.use_tan_omnes, + self.use_tan_phi0, + self.use_xsub_phi0, + self.use_reduced_integrand, + x == self.a, + self.use_cutoff, + self.cutoff_power, + self.cutoff_s, + ) / calc::delta_with_lut( + &mut self.calc_cache, + x, + ); + match self.scaling_type { + ScalingType::Linear => { + val + } + ScalingType::Logarithmic => { + val.log10() + } + } + }) + .collect() }; self.plot_data.push(PlotCurve { points: x_values @@ -388,10 +377,130 @@ impl eframe::App for App { self.n_omnes, self.n_phi0 ), + + }); + } + 4 => { + let cutoff_factor = calc::omnes(&mut self.calc_cache, self.cutoff_s, self.n_omnes, self.use_tan_omnes).abs() * self.cutoff_s; + println!("{:?}", self.scaling_type); + let y_values_func: Vec = x_values.iter().map(|&x| { + let val = (cutoff_factor*x.powf(self.cutoff_power)).powi(2); + match self.scaling_type { + ScalingType::Linear => { + val + } + ScalingType::Logarithmic => { + val.log10() + } + } + }).collect(); + self.plot_data.push(PlotCurve { + points: x_values + .iter() + .zip(y_values_func.iter()) + .map(|(&x, &y)| [x.powf(0.5), y]) + .collect(), + name:"cont. polynomial".to_string() + }); + } + 5 => { + let cutoff_factor = calc::omnes(&mut self.calc_cache, self.cutoff_s, self.n_omnes, self.use_tan_omnes) * self.cutoff_s; + println!("{:?}", cutoff_factor); + let y_values_func: Vec = x_values + .iter() + .map(|&x| { + let val = calc::phi0_integrand( + &mut self.calc_cache, + x, + 2.0, + self.n_omnes, + self.use_tan_omnes, + self.use_tan_phi0, + self.use_xsub_phi0, + self.use_reduced_integrand, + self.use_cutoff, + self.cutoff_power, + self.cutoff_s, + cutoff_factor); + match self.scaling_type { + ScalingType::Linear => { + val + } + ScalingType::Logarithmic => { + val.log10() + } + } + } + ).collect(); + self.plot_data.push(PlotCurve { + points: x_values + .iter() + .zip(y_values_func.iter()) + .map(|(&x, &y)| [x.powf(0.5), y]) + .collect(), + name: format!( + "(#{}) phi0 integrand", + self.plot_data.len() + 1, + ), }) } - 4_usize.. => { - panic!("Not all buttons have been assigned data!"); + 6 => { + let bar = indicatif::ProgressBar::new(u64::from(self.n_points)); + let y_values: Vec = + x_values.iter().map(|&x| { + bar.inc(1); + let phi0_nocut = calc::phi0( + &mut self.calc_cache, + 3.0, + self.n_omnes, + self.n_phi0, + self.use_tan_omnes, + self.use_tan_phi0, + self.use_xsub_phi0, + self.use_reduced_integrand, + false, + false, + self.cutoff_power, + x + ); + let phi0_cut = calc::phi0( + &mut self.calc_cache, + 3.0, + self.n_omnes, + self.n_phi0, + self.use_tan_omnes, + self.use_tan_phi0, + self.use_xsub_phi0, + self.use_reduced_integrand, + false, + true, + self.cutoff_power, + x + ); + match self.scaling_type { + ScalingType::Linear => { + (phi0_cut / phi0_nocut - 1.0).abs() + } + ScalingType::Logarithmic => { + (phi0_cut / phi0_nocut - 1.0).abs().log10() + } + } + }).collect(); + self.plot_data.push(PlotCurve { + points: x_values + .iter() + .zip(y_values.iter()) + .map(|(&x, &y)| [x.powf(0.5), y]) + .collect(), + name: format!( + "(#{}) phi0_cut / phi0, cutoff power: {}", + self.plot_data.len() + 1, + self.cutoff_power, + ), + }); + } + 7_usize.. => { + panic!("This button has no assigned purpose!"); } } } @@ -464,30 +573,33 @@ fn main() { let mut app = App::default(); app.calc_cache = calc::Cache::from_file("./gauss_quad_lut.morello").unwrap(); - let x_values: Vec = (0..app.n_points) - .map(|x| -> f64 { f64::from(x) * ((app.b - app.a) / f64::from(app.n_points)) + app.a }) - .collect(); - let bar = indicatif::ProgressBar::new(u64::from(app.n_points)); - let values: Vec = x_values - .iter() - .map(|&x| { - bar.inc(1); - calc::phi0( - &mut app.calc_cache, - x, - 1000, - 1000, - app.use_tan_omnes, - app.use_tan_phi0, - app.use_xsub_phi0, - app.use_reduced_integrand, - x == app.a, - ) - }) - .collect(); + // let x_values: Vec = (0..app.n_points) + // .map(|x| -> f64 { f64::from(x) * ((app.b - app.a) / f64::from(app.n_points)) + app.a }) + // .collect(); + // let bar = indicatif::ProgressBar::new(u64::from(app.n_points)); + // let values: Vec = x_values + // .iter() + // .map(|&x| { + // bar.inc(1); + // calc::phi0( + // &mut app.calc_cache, + // x, + // 1000, + // 1000, + // app.use_tan_omnes, + // app.use_tan_phi0, + // app.use_xsub_phi0, + // app.use_reduced_integrand, + // x == app.a, + // ) + // }) + // .collect(); + + let phi0 = calc::phi0(&mut app.calc_cache, 10000.0, 10000, 10000, true, true, true, true, true, false, 0.0, 0.0); + println!("relative difference: {}", phi0 / calc::delta(10000.0) - 1.0); // eframe::run_native( - // "Titel", + // "Omnes Calculator", // eframe::NativeOptions::default(), // Box::new(|_cc| Box::new(app)), // )