Test for https://github.com/Reference-LAPACK/lapack/pull/570

dlange-test-exact-float.cpp

// Copyright 2021 Weslley S Pereira

#include <cmath>
#include <limits>
#include <iostream>
#include <bitset>

using Integer = std::int32_t;

extern "C" {
    double dlange_(
        char* norm, Integer* m, Integer* n,
        const double* A, Integer* lda,
        double* work, std::size_t norm_len);
}

template<typename T>
void print_binary( const T& a ){
    const char* a_char = reinterpret_cast<const char*>(&a);
    for (size_t i = 0; i < sizeof(T); ++i)
    {
        std::bitset<8> x(a_char[sizeof(T)-i-1]);
        std::cout << x;
    }
    std::cout << std::endl;
}

int main() {
    using Real = double;

    auto normid = 'F';
    const Integer N = 3793;
    Integer m, n;
    Integer *lda = &m;
    auto work = std::numeric_limits<Real>::quiet_NaN();
    Real norm, expectedNorm, a, *A;

    a = pow(2, 1-1023);
    std::cout << "a            = " << a << " = (memory) ";
    print_binary( a );
    
    m = N;
    n = N;
    A = new Real[m*n];
    for (Integer i = 0; i < m*n; ++i) A[i] = a;

    expectedNorm = N*a;
    std::cout << "expectedNorm = " << expectedNorm << " = (memory) ";
    print_binary( expectedNorm );

    norm = dlange_(&normid, &m, &n, A, lda, &work, 1);
    std::cout << "dlange       = " << norm << " = (memory) ";
    print_binary( norm );

    if(!std::isfinite(norm) || norm != expectedNorm) {
        std::fprintf(
            stderr, "\nexpected norm %8.2e, got %e\n", expectedNorm, norm);
    }

    delete[] A;

    return 0;
}