slatec 0.1.0

Safe Rust interface to selected SLATEC numerical routines
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93

//! Manufactured constant solution for FISHPACK's unit-sphere surface driver.
// slatec-safe-example

use slatec::differential_equations::pde::{
    AxisymmetricSphericalHelmholtz2d, ColatitudeAxis, ColatitudeBoundary, CoordinateBoundary,
    FishpackGrid2, LongitudeAxis, RadialAxis, RadialBoundary, SphereSurfaceHelmholtz2d,
    StaggeredAxisymmetricSphericalHelmholtz2d, StaggeredColatitudeAxis, StaggeredLongitudeAxis,
    StaggeredRadialAxis, StaggeredSphereSurfaceHelmholtz2d,
};

fn main() -> Result<(), Box<dyn std::error::Error>> {
    let solution = SphereSurfaceHelmholtz2d::new(
        ColatitudeAxis::new(0.5, 2.5, 8)?,
        LongitudeAxis::full_circle(8)?,
        0.0,
        FishpackGrid2::zeros(9, 9)?,
        ColatitudeBoundary::Dirichlet {
            lower: vec![1.0; 9],
            upper: vec![1.0; 9],
        },
        CoordinateBoundary::Periodic,
    )?
    .solve()?;

    assert!(
        solution
            .values()
            .values()
            .iter()
            .all(|value| (*value - 1.0).abs() < 1.0e-4)
    );

    let staggered = StaggeredSphereSurfaceHelmholtz2d::new(
        StaggeredColatitudeAxis::new(0.5, 2.5, 8)?,
        StaggeredLongitudeAxis::full_circle(8)?,
        0.0,
        FishpackGrid2::zeros(8, 8)?,
        ColatitudeBoundary::Dirichlet {
            lower: vec![1.0; 8],
            upper: vec![1.0; 8],
        },
        CoordinateBoundary::Periodic,
    )?
    .solve()?;
    assert!(
        staggered
            .values()
            .values()
            .iter()
            .all(|value| (*value - 1.0).abs() < 1.0e-4)
    );

    let axisymmetric = AxisymmetricSphericalHelmholtz2d::new(
        ColatitudeAxis::full_sphere(8)?,
        RadialAxis::new(1.0, 2.0, 8)?,
        0.0,
        FishpackGrid2::zeros(9, 9)?,
        ColatitudeBoundary::BothPoles,
        RadialBoundary::Dirichlet {
            lower: vec![1.0; 9],
            upper: vec![1.0; 9],
        },
    )?
    .solve()?;
    assert!(
        axisymmetric
            .values()
            .values()
            .iter()
            .all(|value| (*value - 1.0).abs() < 5.0e-4)
    );

    let staggered_axisymmetric = StaggeredAxisymmetricSphericalHelmholtz2d::new(
        StaggeredColatitudeAxis::full_sphere(8)?,
        StaggeredRadialAxis::new(1.0, 2.0, 8)?,
        0.0,
        FishpackGrid2::zeros(8, 8)?,
        ColatitudeBoundary::BothPoles,
        RadialBoundary::Dirichlet {
            lower: vec![1.0; 8],
            upper: vec![1.0; 8],
        },
    )?
    .solve()?;
    assert!(
        staggered_axisymmetric
            .values()
            .values()
            .iter()
            .all(|value| (*value - 1.0).abs() < 5.0e-4)
    );
    Ok(())
}