#!/usr/bin/perl
package Math::Interval;

use warnings;
use strict;
use List::Util qw( max min reduce );

sub _eps_for {
    my ($num, $epsilon) = (shift) x 2;   # copy arg to both vars
    $epsilon /= 2.0 while $num + $epsilon / 2.0 != $num;
    return $epsilon;
}

sub _interval {
    my ($min, $max) = ( min(@_), max(@_) );
    return bless [$min - _eps_for($min), $max + _eps_for($max)], __PACKAGE__;
}

sub lb {
    my ($x) = @_;
    return $x->[0];
}

sub ub {
    my ($x) = @_;
    return $x->[1];
}

sub abs {
    my ($x) = @_;
    return max( abs($x->lb()), abs($x->ub()) );
}


sub matrix_norm {
    my ($mat_ref) = @_;
    my @row_sums;
    foreach my $row_ref (@$mat_ref) {
        my $sum = reduce { $a->abs() + $b->abs() } @$row_ref;
        push @row_sums, $sum;
    }
    return max @row_sums;
}


package main;

# getting the absolute value of a single interval works
my $a = Math::Interval::_interval(0.333333, 0.333334);
print 'abs([0.333333; 0.333334] is ', $a->abs(), "\n";

# getting the matrix norm for a matrix of intervals doesn't
my $x = [
    [
        Math::Interval::_interval(3),
        Math::Interval::_interval(2),
        Math::Interval::_interval(3)
    ],
    [
        Math::Interval::_interval(5),
        Math::Interval::_interval(9),
        Math::Interval::_interval(8)
    ],
];

print 'mat norm ', Math::Interval::matrix_norm($x), "\n";