Implemented Op::{Ln, Exp} added precsion arg to raphson

This commit is contained in:
Shimun 2019-11-22 15:46:01 +01:00
parent 8d19b238f1
commit da47b6baa2
Signed by: shimun
GPG Key ID: E81D8382DC2F971B

View File

@ -4,16 +4,29 @@ use std::collections::HashMap;
use std::fmt; use std::fmt;
use std::iter::successors; use std::iter::successors;
#[derive(Clone)] #[derive(Clone)]
struct Context { struct MultiContext {
vars: HashMap<String, f64>, vars: HashMap<String, f64>,
} }
impl Context { trait Context {
fn get(&self, name: &str) -> Option<f64>;
}
/// Simple Context will return it's value for any variable name
impl Context for f64 {
fn get(&self, name: &str) -> Option<f64> {
Some(*self)
}
}
impl Context for MultiContext {
fn get(&self, name: &str) -> Option<f64> { fn get(&self, name: &str) -> Option<f64> {
return self.vars.get(name).cloned(); return self.vars.get(name).cloned();
} }
}
fn with(&self, name: String, value: f64) -> Context { impl MultiContext {
fn with(&self, name: String, value: f64) -> Self {
let mut n = self.clone(); let mut n = self.clone();
n.vars.insert(name, value); n.vars.insert(name, value);
n n
@ -21,7 +34,7 @@ impl Context {
} }
trait Eval: Derive + fmt::Display { trait Eval: Derive + fmt::Display {
fn eval(&self, ctx: &Context) -> f64; fn eval(&self, ctx: &dyn Context) -> f64;
fn boxed_clone(&self) -> Box<dyn Eval>; fn boxed_clone(&self) -> Box<dyn Eval>;
} }
@ -56,10 +69,10 @@ impl Derive for Value {
} }
impl Eval for Value { impl Eval for Value {
fn eval(&self, ctx: &Context) -> f64 { fn eval(&self, ctx: &dyn Context) -> f64 {
match self { match self {
Self::Const(c) => *c, Self::Const(c) => *c,
Self::Eval(e) => e.eval(&ctx), Self::Eval(e) => e.eval(ctx),
} }
} }
fn boxed_clone(&self) -> Box<dyn Eval> { fn boxed_clone(&self) -> Box<dyn Eval> {
@ -72,7 +85,7 @@ impl Eval for Value {
struct Var(String); struct Var(String);
impl Eval for Var { impl Eval for Var {
fn eval(&self, ctx: &Context) -> f64 { fn eval(&self, ctx: &dyn Context) -> f64 {
ctx.get(&self.0).unwrap() ctx.get(&self.0).unwrap()
} }
fn boxed_clone(&self) -> Box<dyn Eval> { fn boxed_clone(&self) -> Box<dyn Eval> {
@ -112,16 +125,22 @@ enum Op {
Sub(Value, Value), Sub(Value, Value),
#[display(fmt = "{}**({})", _0, _1)] #[display(fmt = "{}**({})", _0, _1)]
Pow(Value, Value), Pow(Value, Value),
#[display(fmt = "ln({})", _0)]
Ln(Value),
#[display(fmt = "exp({})", _0)]
Exp(Value),
} }
impl Eval for Op { impl Eval for Op {
fn eval(&self, ctx: &Context) -> f64 { fn eval(&self, ctx: &dyn Context) -> f64 {
match self { match self {
Self::Mul(a, b) => a.eval(ctx) * b.eval(ctx), Self::Mul(a, b) => a.eval(ctx) * b.eval(ctx),
Self::Div(a, b) => a.eval(ctx) / b.eval(ctx), Self::Div(a, b) => a.eval(ctx) / b.eval(ctx),
Self::Add(a, b) => a.eval(ctx) + b.eval(ctx), Self::Add(a, b) => a.eval(ctx) + b.eval(ctx),
Self::Sub(a, b) => a.eval(ctx) - b.eval(ctx), Self::Sub(a, b) => a.eval(ctx) - b.eval(ctx),
Self::Pow(a, b) => a.eval(ctx).powf(b.eval(ctx)), Self::Pow(a, b) => a.eval(ctx).powf(b.eval(ctx)),
Self::Ln(x) => x.eval(ctx).ln(),
Self::Exp(x) => x.eval(ctx).exp(),
} }
} }
fn boxed_clone(&self) -> Box<dyn Eval> { fn boxed_clone(&self) -> Box<dyn Eval> {
@ -148,7 +167,7 @@ impl Derive for Op {
Self::Sub(a, b) => Op::Sub(a.derive(), b.derive()), Self::Sub(a, b) => Op::Sub(a.derive(), b.derive()),
Self::Pow(a, b) => match b { Self::Pow(a, b) => match b {
Value::Eval(e) => Op::Mul( Value::Eval(e) => Op::Mul(
unimplemented!("Ln(x)"), Op::Ln(a.clone()).into(),
Op::Pow(a.clone(), b.clone()).into(), Op::Pow(a.clone(), b.clone()).into(),
), ),
b => Op::Pow( b => Op::Pow(
@ -156,6 +175,8 @@ impl Derive for Op {
Op::Sub(b.clone(), Value::Const(1.0)).into(), Op::Sub(b.clone(), Value::Const(1.0)).into(),
), ),
}, },
Self::Ln(x) => Op::Div(Value::Const(1.0), x.clone()),
Self::Exp(x) => Op::Mul(x.derive(), Op::Exp(x.clone()).into()),
} }
.into() .into()
} }
@ -167,22 +188,30 @@ impl Into<Value> for Op {
} }
} }
fn newton_raphson<'a>(eq: &'a dyn Eval, x0: f64) -> impl Iterator<Item = f64> + 'a { fn newton_raphson<'a>(eq: &'a dyn Eval, x0: f64, precision: i32) -> impl Iterator<Item = f64> + 'a {
//xp - f(x)/f'(x) //xp - f(x)/f'(x)
let d = eq.derive(); let d = eq.derive();
let eq = eq.clone(); let eq = eq.clone();
let mut i = 0; let mut i = 0;
let mut same_counter = 0usize;
successors(Some(x0), move |xp: &f64| -> Option<f64> { successors(Some(x0), move |xp: &f64| -> Option<f64> {
let ctx = Context { let ctx = xp;
vars: vec![("x".to_string(), *xp)].into_iter().collect(),
};
let calc = Op::Sub( let calc = Op::Sub(
Var::val("x"), Var::val("x"),
Op::Div(Value::Eval(eq.boxed_clone()), d.clone()).into(), Op::Div(Value::Eval(eq.boxed_clone()), d.clone()).into(),
); );
println!("( {}.) {} -> {}\n", &i, &calc, &xp); println!("( {}.) {} -> {}\n", &i, &calc, &xp);
i += 1; i += 1;
Some(calc.eval(&ctx)).filter(|x| (x - xp).abs() > (10.0f64).powf(-8.0)) let x = calc.eval(ctx);
let equal = (x - xp).abs() < (10.0f64).powi(-1 * precision as i32);
if equal {
dbg!(("Reached precision after ", i));
}
if same_counter > 3 {
dbg!(("Aborted after", i));
}
same_counter = (same_counter + equal as usize) * (equal as usize);
Some(x).filter(|_| !equal && same_counter < 3)
}) })
} }
@ -192,11 +221,11 @@ fn n_root(number: f64, n: f64) -> Option<f64> {
Value::Const(number), Value::Const(number),
) )
.into(); .into();
newton_raphson(&eq, 1.0).last() newton_raphson(&eq, 1.0, 10).last()
} }
fn main() { fn main() {
let ctx = Context { let ctx = MultiContext {
vars: vec![("x".to_string(), 5.0)].into_iter().collect(), vars: vec![("x".to_string(), 5.0)].into_iter().collect(),
}; };
let eq = Op::Sub( let eq = Op::Sub(
@ -216,7 +245,7 @@ fn main() {
Op::Pow(Value::Const(b), Var::val("x")).into() Op::Pow(Value::Const(b), Var::val("x")).into()
} }
let eq: Value = Op::Sub(Op::Add(pow_to_x(4.0), pow_to_x(6.0)).into(), pow_to_x(9.0)).into(); let eq: Value = Op::Sub(Op::Add(pow_to_x(4.0), pow_to_x(6.0)).into(), pow_to_x(9.0)).into();
//let res = n_root(9.0, 2.0); let res = n_root(9.0, 2.0);
let res = newton_raphson(&eq, 1.0).last(); let res = newton_raphson(&eq, 1.0, 10).last();
println!("1.0 - {}/{} ==> {}", &eq, eq.derive(), res.unwrap()); println!("1.0 - {}/{} ==> {}", &eq, "", res.unwrap());
} }