Struct StationarySolver 

pub struct StationarySolver {
    pub sor_omega: f64,
    pub tol: f64,
    pub max_iter: usize,
    pub max_policy_iter: usize,
    pub policy_tol: f64,
}

Configuration for the stationary solver.

Fields

sor_omega: f64

Successive over-relaxation parameter.

tol: f64

SOR convergence tolerance.

max_iter: usize

Maximum SOR iterations.

max_policy_iter: usize

Maximum policy-iteration sweeps.

policy_tol: f64

Policy-iteration convergence tolerance on the value-function change.

Implementations

impl StationarySolver

pub fn new() -> Self

Creates a new stationary solver with default settings.

Examples

use solver::numeric::finite_difference::elliptic::StationarySolver;
let solver = StationarySolver::new();
assert_eq!(solver.sor_omega, 1.2);

pub fn with_sor_omega(self, omega: f64) -> Self

Selects the SOR relaxation parameter.

Examples

use solver::numeric::finite_difference::elliptic::StationarySolver;
let solver = StationarySolver::new().with_sor_omega(1.5);
assert_eq!(solver.sor_omega, 1.5);

pub fn with_tol(self, tol: f64) -> Self

Selects the SOR tolerance.

Examples

use solver::numeric::finite_difference::elliptic::StationarySolver;
let solver = StationarySolver::new().with_tol(1e-6);
assert_eq!(solver.tol, 1e-6);

pub fn with_max_policy_iter(self, max_policy_iter: usize) -> Self

Selects the maximum number of policy-iteration sweeps.

Examples

use solver::numeric::finite_difference::elliptic::StationarySolver;
let solver = StationarySolver::new().with_max_policy_iter(50);
assert_eq!(solver.max_policy_iter, 50);

pub fn with_policy_tol(self, policy_tol: f64) -> Self

Selects the policy-iteration convergence tolerance.

Examples

use solver::numeric::finite_difference::elliptic::StationarySolver;
let solver = StationarySolver::new().with_policy_tol(1e-6);
assert_eq!(solver.policy_tol, 1e-6);

pub fn solve<const N: usize, P: EllipticProblem<N>>( &self, grid: &Grid<N>, problem: &P, ) -> Vec<f64>

Solves -T u + reaction * u = f on grid.

Examples

use solver::core::grid::Grid;
use solver::models::control::StateDerivatives;
use solver::numeric::finite_difference::discretization::{
    BoundaryCondition, BoundaryConditions, DimensionKind, Transport,
};
use solver::numeric::finite_difference::elliptic::{
    EllipticProblem, StationarySolver,
};

struct Poisson { dx: f64 }
impl EllipticProblem<1> for Poisson {
    fn dimension_kind(&self, _dim: usize) -> DimensionKind {
        DimensionKind::Diffusion
    }
    fn transport(&self, _s: &[f64; 1], _d: &StateDerivatives<1>) -> Transport<1> {
        let c = 1.0 / (self.dx * self.dx);
        Transport::new([c], [c], 0.0)
    }
    fn rhs(&self, _s: &[f64; 1]) -> f64 { 0.0 }
    fn boundary_conditions(&self) -> BoundaryConditions<1> {
        BoundaryConditions::new(
            [BoundaryCondition::Dirichlet(0.0)],
            [BoundaryCondition::Dirichlet(1.0)],
        )
    }
}

let grid = Grid::<1>::new([11], [0.0], [1.0]);
let u = StationarySolver::new().solve(&grid, &Poisson { dx: grid.dx[0] });
assert!((u[0] - 0.0).abs() < 1e-12);
assert!((u[10] - 1.0).abs() < 1e-12);

pub fn solve_control<const N: usize, P: EllipticControlProblem<N> + Sync>( &self, grid: &Grid<N>, problem: &P, ) -> Vec<f64>

where P::Control: Send + Sync,

Solves a stationary stochastic optimal control problem by policy iteration.

The problem is 0 = sup_u { f + L^u V - r V }. At each sweep the solver optimizes the control at the current derivative bundle, assembles -T V + r V = source, solves the linear system, and repeats until the value function stops changing.

Examples

use solver::core::grid::Grid;
use solver::models::control::{ControlProblem, StateDerivatives};
use solver::numeric::finite_difference::discretization::{
    BoundaryConditions, DimensionKind, Transport,
};
use solver::numeric::finite_difference::elliptic::{
    EllipticControlProblem, StationarySolver,
};

struct Constant;
impl ControlProblem<1> for Constant {
    type Control = f64;
    fn optimize(&self, _t: f64, _s: &[f64; 1], _d: &StateDerivatives<1>) -> f64 { 0.0 }
    fn running_reward(&self, _t: f64, _s: &[f64; 1], _c: &f64) -> f64 { 0.0 }
    fn generator(&self, _t: f64, _s: &[f64; 1], _c: &f64, _d: &StateDerivatives<1>) -> f64 { 0.0 }
    fn terminal(&self, _s: &[f64; 1]) -> f64 { 0.0 }
    fn discount_rate(&self, _s: &[f64; 1]) -> f64 { 1.0 }
    fn next_step(&self, _t: f64, s: &[f64; 1], _dt: f64, _n: &[f64; 1]) -> [f64; 1] { *s }
}
impl EllipticControlProblem<1> for Constant {
    fn dimension_kind(&self, _dim: usize) -> DimensionKind { DimensionKind::Diffusion }
    fn transport(&self, _s: &[f64; 1], _c: &f64, _d: &StateDerivatives<1>) -> Transport<1> {
        Transport::new([0.0], [0.0], 0.0)
    }
}

let grid = Grid::<1>::new([11], [0.0], [1.0]);
let v = StationarySolver::new().solve_control(&grid, &Constant);
assert!(v.iter().all(|x| x.abs() < 1e-8));

Trait Implementations

impl Clone for StationarySolver

fn clone(&self) -> StationarySolver

Returns a duplicate of the value.

fn clone_from(&mut self, source: &Self)

Performs copy-assignment from source.

impl Copy for StationarySolver

impl Debug for StationarySolver

fn fmt(&self, f: &mut Formatter<'_>) -> Result

Formats the value using the given formatter.

impl Default for StationarySolver

fn default() -> Self

Returns the “default value” for a type.