Entropyk/docs/tutorial/05-solver-configuration.md

Solver Configuration

Entropyk provides multiple solver strategies for thermodynamic system simulation.

Overview

Solver	Convergence	Use Case
Newton-Raphson	Quadratic	Well-conditioned systems with good initial guess
Picard	Linear	Strongly nonlinear systems, more robust
Fallback	Adaptive	Automatic Newton→Picard fallback on divergence

Newton-Raphson Solver

The Newton-Raphson method solves F(x) = 0 by iteratively computing:

x_{n+1} = x_n - J^{-1} * F(x_n)

Where J is the Jacobian matrix of partial derivatives.

Configuration

use entropyk_solver::{NewtonConfig, NewtonSolver, ConvergenceCriteria};

let config = NewtonConfig {
    max_iterations: 100,
    tolerance: 1e-6,
    relaxation: 1.0,  // Damping factor (0 < α ≤ 1)
    ..Default::default()
};

let solver = NewtonSolver::new(config);

With Convergence Criteria

let criteria = ConvergenceCriteria {
    pressure_tolerance: 100.0,     // Pa
    enthalpy_tolerance: 1000.0,    // J/kg
    residual_tolerance: 1e-6,
    ..Default::default()
};

let solver = NewtonSolver::new(NewtonConfig::default())
    .with_convergence_criteria(criteria);

With Initial State

// Provide initial guess for state vector
let initial_state = vec![
    1e6,    // P_edge0 (Pa)
    400e3,  // h_edge0 (J/kg)
    0.9e6,  // P_edge1 (Pa)
    250e3,  // h_edge1 (J/kg)
];

let solver = NewtonSolver::new(NewtonConfig::default())
    .with_initial_state(initial_state);

Jacobian Freezing

For performance, the Jacobian can be reused across iterations:

use entropyk_solver::JacobianFreezingConfig;

let freeze_config = JacobianFreezingConfig {
    max_frozen_iterations: 5,  // Reuse J for up to 5 iterations
    recompute_on_divergence: true,
};

let solver = NewtonSolver::new(NewtonConfig {
    jacobian_freezing: Some(freeze_config),
    ..Default::default()
});

Picard Solver

The Picard (fixed-point) iteration is more robust for strongly nonlinear systems:

x_{n+1} = G(x_n)

Configuration

use entropyk_solver::{PicardConfig, PicardSolver};

let config = PicardConfig {
    max_iterations: 200,
    tolerance: 1e-5,
    relaxation: 0.8,  // Under-relaxation for stability
    ..Default::default()
};

let solver = PicardSolver::new(config);

When to Use Picard

• Phase change problems: Evaporation/condensation with sharp property changes
• Poor initial guess: When Newton diverges
• Near singularities: Close to critical points

Fallback Solver

The Fallback solver automatically switches from Newton to Picard on divergence:

use entropyk_solver::{FallbackConfig, FallbackSolver};

let config = FallbackConfig {
    newton_max_iterations: 50,
    picard_max_iterations: 100,
    switch_on_newton_divergence: true,
    ..Default::default()
};

let solver = FallbackSolver::new(config);

Fallback Behavior

1. Start with Newton-Raphson
2. If Newton diverges (residual increases), switch to Picard
3. If Picard converges, optionally try Newton again for faster convergence

Convergence Criteria

Basic Tolerance

let criteria = ConvergenceCriteria {
    pressure_tolerance: 100.0,     // |ΔP| < 100 Pa
    enthalpy_tolerance: 1000.0,    // |Δh| < 1000 J/kg
    residual_tolerance: 1e-6,      // ||F(x)|| < 1e-6
    max_iterations: 100,
};

Per-Circuit Convergence

For multi-circuit systems, convergence is tracked per circuit:

// After solving, check the convergence report
if let Some(report) = converged_state.convergence_report() {
    for circuit in &report.circuits {
        println!(
            "Circuit {}: converged={}, max_ΔP={:.1} Pa, max_Δh={:.1} J/kg",
            circuit.circuit_id,
            circuit.converged,
            circuit.max_pressure_change,
            circuit.max_enthalpy_change
        );
    }
}

Smart Initialization

Good initial guess is critical for convergence. Use SmartInitializer:

use entropyk_solver::{SmartInitializer, InitializerConfig};

let config = InitializerConfig {
    fluid: "R134a".to_string(),
    source_temperature_celsius: 40.0,  // Condensing temperature
    sink_temperature_celsius: 5.0,     // Evaporating temperature
    superheat_kelvin: 5.0,
    subcool_kelvin: 3.0,
};

let initializer = SmartInitializer::new(config);

// Estimate saturation pressures
let (p_high, p_low) = initializer.estimate_pressures(
    40.0,  // Source temperature (°C)
    5.0,   // Sink temperature (°C)
);

// Populate state vector
let state = initializer.populate_state(&system, 40.0, 5.0)?;

Running the Solver

Basic Solve

use entropyk_solver::{System, NewtonSolver, NewtonConfig};

fn main() -> Result<(), Box<dyn std::error::Error>> {
    let mut system = build_system();
    system.finalize()?;
    
    let solver = NewtonSolver::new(NewtonConfig::default());
    let result = solver.solve(&system)?;
    
    println!("Converged in {} iterations", result.iterations);
    println!("Final state: {:?}", result.state);
    
    Ok(())
}

With Timeout

use entropyk_solver::TimeoutConfig;

let config = NewtonConfig {
    timeout: Some(TimeoutConfig {
        max_time_seconds: 10.0,
        return_best_on_timeout: true,
    }),
    ..Default::default()
};

Error Handling

use entropyk_solver::SolverError;

match solver.solve(&system) {
    Ok(result) => {
        println!("Converged: {:?}", result.status);
    }
    Err(SolverError::MaxIterationsExceeded { iterations, residual }) => {
        eprintln!("Failed to converge after {} iterations", iterations);
        eprintln!("Final residual: {}", residual);
    }
    Err(SolverError::JacobianSingular { condition_number }) => {
        eprintln!("Jacobian is singular (condition number: {})", condition_number);
        eprintln!("Check system topology and initial guess");
    }
    Err(e) => {
        eprintln!("Solver error: {:?}", e);
    }
}

Solver Strategy Selection

use entropyk_solver::SolverStrategy;

let strategy = SolverStrategy::Fallback;  // Recommended for production

match strategy {
    SolverStrategy::Newton => {
        let solver = NewtonSolver::new(NewtonConfig::default());
    }
    SolverStrategy::Picard => {
        let solver = PicardSolver::new(PicardConfig::default());
    }
    SolverStrategy::Fallback => {
        let solver = FallbackSolver::new(FallbackConfig::default());
    }
}

Best Practices

1. Always use SmartInitializer - Good initial guess is critical
2. Use Fallback solver - More robust for production use
3. Set appropriate tolerances - Balance accuracy vs. convergence time
4. Monitor convergence - Log iteration progress for debugging
5. Handle errors gracefully - Provide meaningful feedback on failure

Python Bindings

import entropyk

# Note: Python bindings use placeholder adapters
# The solver won't converge because there are no real physics equations

system = entropyk.System()
# ... build system ...

solver = entropyk.NewtonSolver(
    max_iterations=100,
    tolerance=1e-6
)

result = solver.solve(system)
print(f"Status: {result.status}")
print(f"Iterations: {result.iterations}")

Next Steps

• Components Reference - Detailed component documentation
• Building Systems - Create system topology